数值代数实验报告1

合集下载

数值计算方法实验报告(含所有)

数值计算方法实验报告(含所有)

本科实验报告课程名称:计算机数值方法实验项目:计算机数值方法实验实验地点:专业班级:学号:学生姓名:xxx指导教师:xxx太原理工大学学生实验报告学院名称软件学院专业班级1217班学号201200xxxx 学生姓名xx 实验日期2014.05.21 成绩课程名称数值计算方法实验题目实验一方程求解一、实验目的和要求熟悉使用、迭代法、牛顿法、割线法等方法对给定的方程进行根的求解。

选择上述方法中的两种方法求方程:二分法f(x)=x3+4x2-10=0在[1,2]内的一个实根,且要求满足精度|x*-x n|<0.5×10-5二、主要设备笔记本 HP ProBook 6470b 一台编译软件:VC++6.0三、实验内容和原理函数f(x)在区间(x,y)上连续,先在区间(x,y)确定a与b,若f(a),f(b)异号,说明在区间(a,b)内存在零点,然后求f[(a+b)/2]。

假设F(a)<0,F(b)>0,a<b,①如果f[(a+b)/2]=0,该点即为零点;②如果f[(a+b)/2]<0,则区间((a+b)/2,b)内存在零点,(a+b)/2≥a;③如果f[(a+b)/2]>0,则区间(a,(a+b)/2)内存在零点,(a+b)/2≤b;返回①重新循环,不断接近零点。

通过每次把f(x)的零点所在区间收缩一半的方法,使区间内的两个端点逐步逼近函数零点,最终求得零点近似值。

四、操作方法与实验步骤1. 二分法:#include<stdio.h>#include<stdlib.h>#include<math.h>int main(){double a=1.0, b=2.0;double x,s;printf(" An\t\tBn\t\tF(Xn)\n");while(1){x=(a+b)/2;s=pow(x,3)+4*x*x-10;if (-0.000005 < s && s < 0.000005){break;}else if(s < 0){a=x;}else if(s > 0){b=x;}printf("%f\t%f\t%f\n",a,b,s);}printf("X的值为:%f\n",x);printf("误差:\t%f\n",s);return 0;}2. 割线法:#include"stdio.h"#include"math.h"int main(){float c,a=1.0,b=2.0;printf("每次得到的X的近似值:\n");while(1){c=b-(b*b*b+4*b*b-10)*(b-a)/(b*b*b+4*b*b-(a*a*a+4*a*a));if(fabs(b-c)<0.5*0.00001)break;b=c;printf("%f\n",b);}printf("X的值为:%f\n",c);}五、实验结果与分析二分法割线法分析:由程序知,使用二分法和割线法均能计算出方程的根,但利用割线法要比二分法计算的次数少,并且能够较早的达到精度要求。

数值代数实验报告

数值代数实验报告

数值代数实验报告数值代数实验报告引言:数值代数是一门研究数值计算方法和算法的学科,它在科学计算和工程应用中起着重要的作用。

本实验报告旨在通过实际的数值计算问题,探讨数值代数的应用和效果。

实验一:线性方程组求解线性方程组求解是数值代数中的一个重要问题。

在实验中,我们使用了高斯消元法和LU分解法两种求解线性方程组的方法,并对比了它们的效果。

首先,我们考虑一个3×3的线性方程组:2x + 3y - z = 54x - 2y + 2z = 1x + y + z = 3通过高斯消元法,我们将该方程组转化为上三角形式,并得到解x=1, y=2, z=0。

而通过LU分解法,我们将该方程组分解为LU两个矩阵的乘积,并得到相同的解。

接下来,我们考虑一个更大的线性方程组,例如10×10的方程组。

通过比较高斯消元法和LU分解法的运行时间,我们可以发现LU分解法在处理大规模方程组时更加高效。

实验二:特征值与特征向量计算特征值与特征向量计算是数值代数中的另一个重要问题。

在实验中,我们使用了幂法和QR方法两种求解特征值与特征向量的方法,并对比了它们的效果。

首先,我们考虑一个3×3的矩阵:1 2 34 5 67 8 9通过幂法,我们可以得到该矩阵的最大特征值为15.372,对应的特征向量为[0.384, 0.707, 0.577]。

而通过QR方法,我们也可以得到相同的结果。

接下来,我们考虑一个更大的矩阵,例如10×10的矩阵。

通过比较幂法和QR 方法的运行时间,我们可以发现QR方法在处理大规模矩阵时更加高效。

实验三:奇异值分解奇异值分解是数值代数中的一种重要技术,它可以将一个矩阵分解为三个矩阵的乘积,从而实现数据降维和信息提取的目的。

在实验中,我们使用了奇异值分解方法,并通过实际的数据集进行了验证。

我们选取了一个包含1000个样本和20个特征的数据集,通过奇异值分解,我们将该数据集分解为三个矩阵U、S和V的乘积。

数值代数上机实验报告

数值代数上机实验报告

数值代数上机实验报告试验项目名称:平方根法与改进平方根法实验内容:先用你熟悉的计算机语言将平方根法和改进平方根法编写成通用的子程序,然后用你编写的程序求解对称正定方程组Ax=b,其中,A=[101 10 1…1 10 11 10]100*100b随机生成,比较计算结果,评论方法优劣。

实验要求:平方根法与改进的平方根的解法步骤;存储单元,变量名称说明;系数矩阵与右端项的生成;结果分析。

实验报告姓名:罗胜利班级:信息与计算科学0802 学号:u200810087 实验一、平方根法与改进平方根法先用你所熟悉的计算机语言将平方根法和改进的平方根法编成通用的子程序,然后用你编写的程序求解对称正定方程组AX=b,其中系数矩阵为40阶Hilbert矩阵,即系数矩阵A的第i行第j列元素为=,向量b的第i个分量为=.平方根法函数程序如下:function [x,b]=pingfanggenfa(A,b)n=size(A);n=n(1);x=A^-1*b; %矩阵求解disp('Matlab自带解即为x');for k=1:nA(k,k)=sqrt(A(k,k));A(k+1:n,k)=A(k+1:n,k)/A(k,k);for j=k+1:n;A(j:n,j)=A(j:n,j)-A(j:n,k)*A(j,k);endend %Cholesky分解for j=1:n-1b(j)=b(j)/A(j,j);b(j+1:n)=b(j+1:n)-b(j)*A(j+1:n,j);endb(n)=b(n)/A(n,n); %前代法A=A';for j=n:-1:2b(j)=b(j)/A(j,j);b(1:j-1)=b(1:j-1)-b(j)*A(1:j-1,j);endb(1)=b(1)/A(1,1); %回代法disp('平方根法的解即为b');endfunction [x]=ave(A,b,n) %用改进平方根法求解Ax=b L=zeros(n,n); %L为n*n矩阵D=diag(n,0); %D为n*n的主对角矩阵S=L*D;for i=1:n %L的主对角元素均为1L(i,i)=1;for i=1:nfor j=1:n %验证A是否为对称正定矩阵if (eig(A)<=0)|(A(i,j)~=A(j,i)) %A的特征值小于0或A非对称时,输出wrong disp('wrong');break;endendendD(1,1)=A(1,1); %将A分解使得A=LDL Tfor i=2:nfor j=1:i-1S(i,j)=A(i,j)-sum(S(i,1:j-1)*L(j,1:j-1)');L(i,1:i-1)=S(i,1:i-1)/D(1:i-1,1:i-1);endD(i,i)=A(i,i)-sum(S(i,1:i-1)*L(i,1:i-1)');endy=zeros(n,1); % x,y为n*1阶矩阵x=zeros(n,1);for i=1:ny(i)=(b(i)-sum(L(i,1:i-1)*D(1:i-1,1:i-1)*y(1:i-1)))/D(i,i); %通过LDy=b解得y的值endfor i=n:-1:1x(i)=y(i)-sum(L(i+1:n,i)'*x(i+1:n)); %通过L T x=y解得x的值改进平方根法函数程序如下:function b=gaijinpinfanggenfa(A,b)n=size(A);n=n(1);v=zeros(n,1);for j=1:nfor i=1:j-1v(i)=A(j,i)*A(i,i);endA(j,j)=A(j,j)-A(j,1:j-1)*v(1:j-1);A(j+1:n,j)=(A(j+1:n,j)-A(j+1:n,1:j-1)*v(1:j-1))/A(j,j);end %LDL'分解B=diag(A);D=zeros(n);for i=1:nD(i,i)=B(i);A(i,i)=1;EndA=tril(A); %得到L和Dfor j=1:n-1b(j)=b(j)/A(j,j);b(j+1:n)=b(j+1:n)-b(j)*A(j+1:n,j);endb(n)=b(n)/A(n,n); %前代法A=D*(A');for j=n:-1:2b(j)=b(j)/A(j,j);b(1:j-1)=b(1:j-1)-b(j)*A(1:j-1,j);endb(1)=b(1)/A(1,1); %回代法disp('改进平方根法解得的解即为b');end调用函数解题:clear;clc;n=input('请输入矩阵维数:');b=zeros(n,1);A=zeros(n);for i=1:nfor j=1:nA(i,j)=1/(i+j-1);b(i)=b(i)+1/(i+j-1);endend %生成hilbert矩阵[x,b]=pingfanggenfa(A,b) b=gaijinpinfanggenfa(A,b)运行结果:请输入矩阵维数:40Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND = 6.570692e-020. > In pingfanggenfa at 4In qiujie at 10Matlab自带解即为x平方根法的解即为bx =1.60358.96850.85621.01950.9375-50.2500-3.0000-16.000024.0000-49.5000-30.000039.000022.0000-64.0000-12.00002.000010.2500-10.5000-1.0000-10.875083.000046.0000-98.0000-69.000068.000021.0000-50.7188-8.7500-8.0000 112.0000 6.0000 -68.7500 22.000044.0000 -28.0000 8.0000 -44.000012.0000b =1.0e+007 *0.0000-0.00000.0001-0.0004-0.00140.0424-0.29801.1419-2.73354.2539-4.30182.7733-1.19890.5406-0.36880.32850.4621-0.25130.05650.0000-0.00510.0071-0.0027-0.0031-0.00190.00090.0002-0.0002-0.00060.00040.0001-0.00020.00010.0000-0.00000.0000-0.0000-0.0000改进平方根法解得的解即为bb =1.0e+024 *0.0000-0.00000.0001-0.0012-0.0954 0.4208 -1.2101 2.0624 -1.0394 -3.3343 6.2567 -0.2463 -7.45942.80303.6990 0.7277 -1.7484 -0.4854 -3.6010 0.2532 5.1862 1.4410 0.8738 -4.5654 1.0422 4.0920 -2.7764 -2.2148 -0.8953 0.3665 4.8967 1.0416 0.1281-1.1902-2.83348.4610-3.6008实验二、利用QR分解解线性方程组:利用QR分解解线性方程组Ax=b,其中A=[16 4 8 4;4 10 8 4;8 8 12 10;4 4 10 12];b=[32 26 38 30];求解程序如下:定义house函数:function [v,B]=house(x)n=length(x);y=norm(x,inf);x=x/y;Q=x(2:n)'*x(2:n);v(1)=1;v(2:n)=x(2:n);if n==1B=0;elsea=sqrt(x(1)^2+Q);if x(1)<=0v(1)=x(1)-a;elsev(1)=-Q/(x(1)+a);endB=2*v(1)^2/(Q+v(1)^2);endend进行QR分解:clear;clc;A=[16 4 8 4;4 10 8 4;8 8 12 10;4 4 10 12]; b=[32 26 38 30];b=b';x=size(A);m=x(1);n=x(2);d=zeros(n,1);for j=1:n[v,B]=house(A(j:m,j));A(j:m,j:n)=(eye(m-j+1)-B*(v')*v)*A(j:m,j:n); d(j)=B;if j<m< p="">A(j+1:m,j)=v(2:m-j+1);endend %QR分解R=triu(A); %得到R D=A;I=eye(m,n);Q=I;for i=1:nD(i,i)=1;endH=tril(D);M=H';for i=1:nN=I-d(i)*H(1:m,i)*M(i,1:m);Q=Q*N;end %得到Qb=(Q')*b; %Q是正交阵for j=n:-1:2b(j)=b(j)/R(j,j);b(1:j-1)=b(1:j-1)-b(j)*R(1:j-1,j);endb(1)=b(1)/R(1,1); %回带法运行结果如下:R =18.7617 9.8072 15.7769 11.08640 9.9909 9.3358 7.53410 0 5.9945 9.80130 0 0 -0.5126Q =0.8528 -0.4368 -0.2297 -0.17090.2132 0.7916 -0.4594 -0.34170.4264 0.3822 0.2844 0.76890.2132 0.1911 0.8095 -0.5126b=1.000000000000001.000000000000010.9999999999999881.00000000000001实验三、Newton下山法解非线性方程组:3x-cos(yz)-=0,-81+sinz+1.06=0,exp(-xy)+20z+=0;要求满足数值解=满足或.定义所求方程组的函数:Newtonfun.mfunction F = Newtonfun(X)F(1,1)=3*X(1)-cos(X(2)*X(3))-1/2;F(2,1)=X(1)^2-81*(X(2)+0.1)^2+sin(X(3))+1.06;F(3,1)=exp(-X(1)*X(2))+20*X(3)+(10*pi-3)/3;End向量求导:Xiangliangqiudao.mfunction J=xiangliangqiudao()syms x y zX=[x,y,z];F=[3*X(1)-cos(X(2)*X(3))-1/2;X(1)^2-81*(X(2)+0.1)^2+sin(X(3))+1.06;exp(-X(1)*X(2))+20*X(3)+(10*pi-3)/3];J=jacobian(F,[x y z]);End代值函数:Jacobi.mfunction F=Jacobi(x)F=[ 3,x(3)*sin(x(2)*x(3)), x(2)*sin(x(2)*x(3));2*x(1), -162*x(2)-81/5,cos(x(3));-x(2)/exp(x(1)*x(2)),-x(1)/exp(x(1)*x(2)),20];End方程组求解:format long; %数据表示为双精度型X1=[0,0,0]';eps=10^(-8);k=1;i=1;X2=X1-Jacobi(X1)^(-1)*Newtonfun(X1);while (norm(subs(X2-X1,pi,3.1415926),2)>=eps)&&(norm(Newtonfun(X1),2)>=eps) if norm(Newtonfun(X2),2)<="" p="">X1=X2;B=inv(Jacobi(X2));C=Newtonfun(X2);X2=X2-B*C;i=i+1;elsev=1/(2^k); %引入下山因子X1=X2;B=inv(Jacobi(X2));C=Newtonfun(X2);X2=X2-v*B*C;k=k+1;endendj=i+k-1 %迭代次数X=X2 %输出结果运行结果如下:j =5X =0.500000000000000 -0.000000000000000 -0.523598775598299</m<>。

数值计算基础实验报告(3篇)

数值计算基础实验报告(3篇)

第1篇一、实验目的1. 理解数值计算的基本概念和常用算法;2. 掌握Python编程语言进行数值计算的基本操作;3. 熟悉科学计算库NumPy和SciPy的使用;4. 分析算法的数值稳定性和误差分析。

二、实验内容1. 实验环境操作系统:Windows 10编程语言:Python 3.8科学计算库:NumPy 1.19.2,SciPy 1.5.02. 实验步骤(1)Python编程基础1)变量与数据类型2)运算符与表达式3)控制流4)函数与模块(2)NumPy库1)数组的创建与操作2)数组运算3)矩阵运算(3)SciPy库1)求解线性方程组2)插值与拟合3)数值积分(4)误差分析1)舍入误差2)截断误差3)数值稳定性三、实验结果与分析1. 实验一:Python编程基础(1)变量与数据类型通过实验,掌握了Python中变量与数据类型的定义方法,包括整数、浮点数、字符串、列表、元组、字典和集合等。

(2)运算符与表达式实验验证了Python中的算术运算、关系运算、逻辑运算等运算符,并学习了如何使用表达式进行计算。

(3)控制流实验学习了if-else、for、while等控制流语句,掌握了条件判断、循环控制等编程技巧。

(4)函数与模块实验介绍了Python中函数的定义、调用、参数传递和返回值,并学习了如何使用模块进行代码复用。

2. 实验二:NumPy库(1)数组的创建与操作通过实验,掌握了NumPy数组的基本操作,包括创建数组、索引、切片、排序等。

(2)数组运算实验验证了NumPy数组在数学运算方面的优势,包括加、减、乘、除、幂运算等。

(3)矩阵运算实验学习了NumPy中矩阵的创建、操作和运算,包括矩阵乘法、求逆、行列式等。

3. 实验三:SciPy库(1)求解线性方程组实验使用了SciPy库中的线性代数模块,通过高斯消元法、LU分解等方法求解线性方程组。

(2)插值与拟合实验使用了SciPy库中的插值和拟合模块,实现了对数据的插值和拟合,并分析了拟合效果。

数值计算方法实验报告

数值计算方法实验报告

数值计算方法实验报告实验目的:通过实验验证不同数值计算方法在求解数学问题时的精度和效率,并分析其优缺点。

实验原理:实验内容:本实验选取了三个典型的数值计算问题,并分别采用了二分法、牛顿迭代法和梯度下降法进行求解。

具体问题和求解方法如下:1. 问题一:求解方程sin(x)=0的解。

-二分法:利用函数值的符号变化将解空间不断缩小,直到找到满足精度要求的解。

-牛顿迭代法:通过使用函数的斜率来逼近方程的解,并不断逼近真实解。

-梯度下降法:将方程转化为一个极小化问题,并利用梯度下降的方式逼近极小值点,进而找到方程的解。

2.问题二:求解函数f(x)=x^2-3x+2的极小值点。

-二分法:通过确定函数在一个区间内的变化趋势,将极小值所在的区间不断缩小,从而找到极小值点。

-牛顿迭代法:通过使用函数的导数和二阶导数来逼近极小值点,并不断逼近真实解。

-梯度下降法:将函数转化为一个极小化问题,并利用梯度下降的方式逼近极小值点,进而找到函数的极小值点。

3. 问题三:求解微分方程dy/dx = -0.1*y的解。

-二分法:通过离散化微分方程,将微分方程转化为一个差分方程,然后通过迭代计算不同点的函数值,从而得到函数的近似解。

-牛顿迭代法:将微分方程转化为一个积分方程,并通过迭代计算得到不同点的函数值,从而得到函数的近似解。

-梯度下降法:将微分方程转化为一个极小化问题,并利用梯度下降的方式逼近极小值点,从而得到函数的近似解。

实验步骤:1.编写代码实现各个数值计算方法的求解过程。

2.对每个数值计算问题,设置合适的初始值和终止条件。

3.运行程序,记录求解过程中的迭代次数和每次迭代的结果。

4.比较不同数值计算方法的精度和效率,并分析其优缺点。

实验结果:经过实验测试,得到了如下结果:-问题一的二分法迭代次数为10次,求解结果为x=0;牛顿迭代法迭代次数为4次,求解结果为x=0;梯度下降法迭代次数为6次,求解结果为x=0。

-问题二的二分法迭代次数为10次,求解结果为x=1;牛顿迭代法迭代次数为3次,求解结果为x=1;梯度下降法迭代次数为4次,求解结果为x=1-问题三的二分法迭代次数为100次,求解结果为y=e^(-0.1x);牛顿迭代法迭代次数为5次,求解结果为y=e^(-0.1x);梯度下降法迭代次数为10次,求解结果为y=e^(-0.1x)。

数值分析第一次实验报告

数值分析第一次实验报告

数值分析实验报告(一)2016级数学基地班尹烁翔320160928411一、问题重述:hamming级数求和二、问题分析级数为∑1k(k+x)∞k=1易知当X=1时,φ(1)=1我们可以考虑这个新级数:φ(x)−φ(1)用这个级数可以使精度更高,误差更小且迭代次数变少。

通分易得:φ(x)−φ(1)=1k(k+x)−1k(k+1)=1−xk(k+x)(k+1)我们还可以继续算得φ(2)及φ(x)−φ(2)这样精度会继续提高,且迭代次数也会减少。

下面考虑误差:由公式可得∑1−xk(k+x)(k+1)∞k=1<1k3<∫1k3∞n−1<10−10要把误差控制在范围内,需要k即迭代次数至少70001次。

三、算法实现:#include<iostream>#include<iomanip>>using namespace std;int main(){double sum;//sum为级数和double x;//x为代入的自变量int k=1;//k为迭代次数for (x=0; x<=10; x=x+0.1)//对0到10以内进行迭代运算,每次加0.1{sum=0;//每迭代完一个x,级数归零for (k=1; k<=70001; k++)//固定x并对k进行运算{sum=sum+1/(k*(k+x)*(k+1));}sum=(1-x)*sum+1.0;cout<<setiosflags(ios::fixed)<<" "<<setprecision(1)<<x;cout<<setiosflags(ios::fixed)<<" "<<setprecision(10)<<sum<<endl;}for (x=11; x<=290; x++)//对11到290以内进行迭代运算,每次加1{sum=0;for (k=1; k<=70001; k++)//固定x{sum=sum+1/(k*(k+x)*(k+1));}sum=(1-x)*sum+1.0;cout<<setiosflags(ios::fixed)<<" "<<setprecision(1)<<x;cout<<setiosflags(ios::fixed)<<" "<<setprecision(10)<<sum<<endl;}for (x=290; x<=300; x=x+0.1)//对290.1到300以内进行迭代运算,每次加0.1 {sum=0;for (k=1; k<=70001; k++)//固定x{sum=sum+1/(k*(k+x)*(k+1));}sum=(1-x)*sum+1.0;cout<<setiosflags(ios::fixed)<<" "<<setprecision(1)<<x;cout<<setiosflags(ios::fixed)<<" "<<setprecision(10)<<sum<<endl;}return 0;}四、数据结果:0.0 1.6449340667 0.1 1.5346072448 0.2 1.4408788415 0.3 1.3600825867 0.4 1.2895778007 0.5 1.2274112777 0.6 1.1721051961 0.7 1.1225193425 0.8 1.07775887270.9 1.03711091781.0 1.0000000000 1.1 0.9659560305 1.2 0.9345909181 1.3 0.9055811887 1.4 0.8786548819 1.5 0.853******* 1.6 0.8301644486 1.7 0.8082346082 1.8 0.78764591881.9 0.76827137672.0 0.7500000000 2.1 0.7327343381 2.2 0.7163884348 2.3 0.7008861540 2.4 0.6861597923 2.5 0.6721489224 2.6 0.6587994241 2.7 0.6460626684 2.8 0.63389482552.9 0.62225627673.0 0.6111111113 3.1 0.6004266954 3.2 0.5901732990 3.3 0.5803237751 3.4 0.5708532792 3.5 0.5617390263 3.6 0.5529600781 3.7 0.5444971556 3.8 0.53633247553.9 0.52844960504.0 0.5208333336 4.1 0.5134695598 4.2 0.5063451894 4.3 0.49944804604.4 0.49276679034.5 0.48629084784.6 0.48001034484.7 0.47391604974.8 0.46799932104.9 0.46225205975.0 0.45666666715.1 0.45123600545.2 0.44595336325.3 0.44081242345.4 0.43580723395.5 0.43093218145.6 0.42618196715.7 0.42155158445.8 0.41703629915.9 0.41263163046.0 0.40833333386.1 0.40413738606.2 0.40003996986.3 0.39603746096.4 0.39212641636.5 0.38830356206.6 0.38456578316.7 0.38091011406.8 0.37733372946.9 0.37383393577.0 0.37040816397.1 0.36705396157.2 0.36376898657.3 0.36055100097.4 0.35739786507.5 0.35430753177.6 0.35127804177.7 0.34830751887.8 0.34539416537.9 0.34253625788.0 0.33973214368.1 0.33698023688.2 0.33427901518.3 0.33162701648.4 0.32902283598.5 0.32646512338.6 0.32395258008.7 0.32148395698.8 0.31905805168.9 0.31667370669.0 0.31432980689.1 0.31202527809.2 0.30975908459.3 0.30753022799.4 0.30533774499.5 0.30318070609.6 0.30105821429.7 0.29896940319.8 0.29691343609.9 0.294889504210.0 0.292896826311.0 0.274534305112.0 0.258600891013.0 0.244625674714.0 0.232254453215.0 0.221215267616.0 0.211295563617.0 0.202326620618.0 0.194172672719.0 0.186723141720.0 0.179886984821.0 0.173588511822.0 0.167764240823.0 0.162360502724.0 0.157331593125.0 0.152638329626.0 0.148246914727.0 0.144128030628.0 0.140256111329.0 0.136608754530.0 0.133166240731.0 0.129911138432.0 0.126827978033.0 0.123902979834.0 0.121123826635.0 0.118479472636.0 0.115959981337.0 0.113556388138.0 0.111260583139.0 0.109065210040.0 0.106963580041.0 0.104949596342.0 0.103017690143.0 0.101162762944.0 0.099380138345.0 0.097665518246.0 0.096014944747.0 0.094424767348.0 0.092891612649.0 0.091412358750.0 0.089984111851.0 0.088604185152.0 0.087270081253.0 0.085979474654.0 0.084730197955.0 0.083520227556.0 0.082347672757.0 0.081210763958.0 0.080107843659.0 0.079037357560.0 0.077997846261.0 0.076987938262.0 0.076006343163.0 0.075051846164.0 0.074123301865.0 0.073219629966.0 0.072339810267.0 0.071482878568.0 0.070647922969.0 0.069834080070.0 0.069040532171.0 0.068266503872.0 0.067511259473.0 0.066774100374.0 0.066054362875.0 0.065351416076.0 0.064664659377.0 0.063993521278.0 0.063337457279.0 0.062695948280.0 0.062068499081.0 0.061454637382.0 0.0608539117 83.0 0.060265891284.0 0.059690163685.0 0.059126334986.0 0.058574027887.0 0.058032881288.0 0.057502549189.0 0.056982699990.0 0.056473015891.0 0.055973191792.0 0.055482935193.0 0.055001964994.0 0.054530011295.0 0.054066814696.0 0.053612125897.0 0.053165704998.0 0.052727321299.0 0.0522967526100.0 0.0518737853101.0 0.0514582132102.0 0.0510498380103.0 0.0506484683104.0 0.0502539197105.0 0.0498660140106.0 0.0494845798107.0 0.0491094512108.0 0.0487404681109.0 0.0483774760110.0 0.0480203256111.0 0.0476688725112.0 0.0473229772113.0 0.0469825047114.0 0.0466473244115.0 0.0463173100116.0 0.0459923394117.0 0.0456722940118.0 0.0453570593119.0 0.0450465242120.0 0.0447405812121.0 0.0444391259122.0 0.0441420572123.0 0.0438492771124.0 0.0435606905125.0 0.0432762052126.0 0.0429957316127.0 0.0427191829128.0 0.0424464746129.0 0.0421775249130.0 0.0419122542131.0 0.0416505852132.0 0.0413924428133.0 0.0411377539134.0 0.0408864476135.0 0.0406384549136.0 0.0403937087137.0 0.0401521437138.0 0.0399136963139.0 0.0396783048140.0 0.0394459089141.0 0.0392164502142.0 0.0389898715143.0 0.0387661174144.0 0.0385451338145.0 0.0383268679146.0 0.0381112684147.0 0.0378982853148.0 0.0376878698149.0 0.0374799743150.0 0.0372745524151.0 0.0370715590152.0 0.0368709499153.0 0.0366726822154.0 0.0364767137155.0 0.0362830036156.0 0.0360915118157.0 0.0359021994158.0 0.0357150281159.0 0.0355299609160.0 0.0353469614161.0 0.0351659940162.0 0.0349870241163.0 0.0348100178164.0 0.0346349421165.0 0.0344617645166.0 0.0342904534167.0 0.0341209780168.0 0.0339533080169.0 0.0337874138170.0 0.0336232666171.0 0.0334608381 172.0 0.0333001006 173.0 0.0331410270 174.0 0.0329835910 175.0 0.0328277666 176.0 0.0326735285 177.0 0.0325208518 178.0 0.0323697123 179.0 0.0322200861 180.0 0.0320719500 181.0 0.0319252812 182.0 0.0317800574 183.0 0.0316362566 184.0 0.0314938575 185.0 0.0313528391 186.0 0.0312131807 187.0 0.0310748622 188.0 0.0309378640 189.0 0.0308021665 190.0 0.0306677509 191.0 0.0305345985 192.0 0.0304026910 193.0 0.0302720107 194.0 0.0301425399 195.0 0.0300142615 196.0 0.029******* 197.0 0.029******* 198.0 0.029******* 199.0 0.029******* 200.0 0.029******* 201.0 0.029******* 202.0 0.029******* 203.0 0.029******* 204.0 0.028******* 205.0 0.028******* 206.0 0.028******* 207.0 0.028******* 208.0 0.028******* 209.0 0.028******* 210.0 0.028******* 211.0 0.028******* 212.0 0.028******* 213.0 0.027******* 214.0 0.027******* 215.0 0.027*******216.0 0.027*******217.0 0.027*******218.0 0.027*******219.0 0.027*******220.0 0.027*******221.0 0.027*******222.0 0.0269466153223.0 0.0268458877224.0 0.0267459700225.0 0.0266468523226.0 0.0265485248227.0 0.0264509777228.0 0.0263542015229.0 0.0262581869230.0 0.0261629247231.0 0.0260684057232.0 0.025*******233.0 0.025*******234.0 0.025*******235.0 0.025*******236.0 0.025*******237.0 0.025*******238.0 0.025*******239.0 0.025*******240.0 0.025*******241.0 0.025*******242.0 0.025*******243.0 0.024*******244.0 0.024*******245.0 0.024*******246.0 0.024*******247.0 0.024*******248.0 0.024*******249.0 0.024*******250.0 0.024*******251.0 0.024*******252.0 0.024*******253.0 0.024*******254.0 0.024*******255.0 0.024*******256.0 0.023*******257.0 0.023*******258.0 0.023*******259.0 0.023*******260.0 0.023*******261.0 0.023*******262.0 0.023*******263.0 0.023*******264.0 0.023*******265.0 0.023*******266.0 0.023*******267.0 0.023*******268.0 0.023*******269.0 0.022*******270.0 0.022*******271.0 0.022*******272.0 0.022*******273.0 0.022*******274.0 0.022*******275.0 0.022*******276.0 0.022*******277.0 0.022*******278.0 0.022*******279.0 0.022*******280.0 0.022*******281.0 0.022*******282.0 0.022*******283.0 0.021*******284.0 0.021*******285.0 0.021*******286.0 0.021*******287.0 0.021*******288.0 0.021*******289.0 0.021*******290.0 0.021*******290.1 0.021*******290.2 0.021*******290.3 0.021*******290.4 0.021*******290.5 0.021*******290.6 0.021*******290.7 0.021*******290.8 0.021*******290.9 0.021*******291.0 0.021*******291.1 0.021*******291.2 0.021*******291.3 0.021******* 291.4 0.021******* 291.5 0.021******* 291.6 0.021******* 291.7 0.021******* 291.8 0.021******* 291.9 0.021******* 292.0 0.021******* 292.1 0.021******* 292.2 0.021******* 292.3 0.021******* 292.4 0.021******* 292.5 0.021******* 292.6 0.021******* 292.7 0.021******* 292.8 0.021******* 292.9 0.021******* 293.0 0.021******* 293.1 0.021******* 293.2 0.021******* 293.3 0.021******* 293.4 0.021******* 293.5 0.021******* 293.6 0.021******* 293.7 0.021******* 293.8 0.021******* 293.9 0.021******* 294.0 0.021******* 294.1 0.021******* 294.2 0.021******* 294.3 0.021******* 294.4 0.021******* 294.5 0.021******* 294.6 0.021******* 294.7 0.021******* 294.8 0.021******* 294.9 0.021******* 295.0 0.021******* 295.1 0.021******* 295.2 0.021******* 295.3 0.021******* 295.4 0.021******* 295.5 0.021******* 295.6 0.021******* 295.7 0.021******* 295.8 0.021******* 295.9 0.021******* 296.0 0.021******* 296.1 0.021******* 296.2 0.021******* 296.3 0.021******* 296.4 0.021******* 296.5 0.021******* 296.6 0.021******* 296.7 0.021******* 296.8 0.021******* 296.9 0.021******* 297.0 0.021******* 297.1 0.021******* 297.2 0.021******* 297.3 0.021******* 297.4 0.021******* 297.5 0.021******* 297.6 0.021******* 297.7 0.021******* 297.8 0.021******* 297.9 0.021******* 298.0 0.021******* 298.1 0.021******* 298.2 0.021******* 298.3 0.021******* 298.4 0.021******* 298.5 0.021******* 298.6 0.021******* 298.7 0.021******* 298.8 0.021******* 298.9 0.021******* 299.0 0.021******* 299.1 0.020******* 299.2 0.020******* 299.3 0.020******* 299.4 0.020******* 299.5 0.020******* 299.6 0.020******* 299.7 0.020******* 299.8 0.020******* 299.9 0.020******* 300.0 0.020*******。

数值代数qr方法实验报告

数值代数qr方法实验报告

一、引言很多情况下在工程中抽象出来的数学方程组是超定的,没有精确解,这样就需要找一个在某种意义下最接近精确解的解。

设A 是m*n 的实矩阵,22A (A )T x bQ x b -=-,这样2min A y b -就等价与2min (A )T Q x b -,为方便求解,需要A T Q 是上三角矩阵,这样引入QR 分解就比较求解方便。

二、豪斯霍尔德变换豪斯霍尔德变换(Householder transformation )又称初等反射(Elementary reflection ),最初由A.C Aitken 在1932年提出,Alston Scott Householder 在1958年指出了这一变换在数值线性代数上的意义。

这一变换将一个向量变换为由一个超平面反射的镜像,是一种线性变换。

其变换矩阵被称作豪斯霍尔德矩阵,在一般内积空间中的类比被称作豪斯霍尔德算子,超平面的法向量被称作豪斯霍尔德向量。

三、理论依据householder 变换:任意的一个向量x ,经过一个正交变换2**T H I W W =- 后总可以变为一个与之范数相等的另一个向量Hx 。

如上图中所示,记v Hx x =-,2/w v v =,2**TH I W W =-上述H 既为所要求的householder 变换。

具体操作时将需要变化的矩阵的每一列当做一个向量,第一列变为除了第一个元素都为0的向量…第i 列变为除了前i 个元素都为0的向量...为了保证在每次变化时不改变已经变好的0元素,第i 次只变化每列第i 个元素到第n 个元素,每个变化矩阵的形式是[I 0;0 H]。

算法如下:qr_Houserhoilder.m function [Q,R]=qr_Houserhoilder(A) [m,n]=size(A); Q=eye(n); for i=1:nu=house(A(i:m,i));P=eye(m-i+1)-2/(norm(u)^2)*u*u'; A(i:m,i:n)=P*A(i:m,i:n);PP=blkdiag(eye(i-1),P);Q=Q*PP;endR=triu(A);程序house.mfunction u=house(Ai)Ai(1)=Ai(1)+sign(Ai(1))*norm(Ai);u=Ai/norm(Ai);四、数值实验结果:为了保证在一定阶数下,所做的数值试验具有代表性,我们随机产生3个n阶矩阵,随着n的不断变化,得到如下一系列图:(1)A=rand(4,4)A =0.3027 0.7894 0.1793 0.39220.2364 0.7134 0.4350 0.39360.4690 0.0742 0.9747 0.11150.5381 0.6768 0.7245 0.8650[Q,R]=qr_Houserhoilder(A)Q =-0.3735 -0.5369 0.4602 0.6003-0.2916 -0.5439 -0.7840 -0.0668-0.5786 0.6445 -0.2679 0.4219-0.6639 -0.0208 0.3190 -0.6761R =-0.8106 -0.9951 -1.2387 -0.90010 -0.7781 0.2803 -0.37070 0 -0.2886 0.11800 0 0 -0.3286test_qr3err =1.0e-10 *0.1950 0.0120 0.0119 0.3029-111 1.52 2.53 3.54(2)A=rand(5,5)A =0.7575 0.7005 0.6420 0.8233 0.88770.6373 0.4698 0.5815 0.5290 0.37700.6643 0.3036 0.9458 0.7210 0.12950.1056 0.0254 0.3903 0.1048 0.30270.2813 0.3787 0.5302 0.0627 0.1156[Q,R]=qr_Houserhoilder(A)Q =-0.6161 0.4378 0.2968 -0.5805 -0.0598-0.5183 -0.0381 0.2932 0.7098 -0.3743-0.5403 -0.6729 -0.2440 -0.1030 0.4303-0.0859 -0.1847 -0.5513 -0.2507 -0.7693-0.2288 0.5657 -0.6801 0.2928 0.2817R =-1.2295 -0.9280 -1.3629 -1.1944 -0.86480 0.2940 -0.1496 -0.1288 0.29660 0 -0.4455 0.1231 0.09690 0 0 -0.1847 -0.30310 0 0 0 -0.3388test_qr3err =1.0e-10 *0.1722 0.0119 0.0120 0.4239-111 1.52 2.53 3.54(3)A=rand(6,6)A =0.5961 0.8422 0.4279 0.3473 0.2731 0.34860.1335 0.2860 0.2157 0.6181 0.6306 0.19300.4479 0.0467 0.3284 0.4588 0.4228 0.95480.9430 0.7296 0.7951 0.0920 0.0322 0.31140.7581 0.7811 0.4338 0.1512 0.6185 0.12710.8647 0.1359 0.1735 0.9036 0.4393 0.8964 [Q,R]=qr_Houserhoilder(A)Q =-0.3571 -0.5571 0.2402 0.2038 0.4576 0.5035-0.0800 -0.2530 -0.1880 0.8253 -0.1645 -0.4314-0.2684 0.3532 -0.4622 0.2419 -0.3487 0.6399-0.5650 -0.0874 -0.6012 -0.3372 0.2898 -0.3376-0.4542 -0.3266 0.3151 -0.2560 -0.7185 -0.0773-0.5181 0.6217 0.4823 0.1991 0.2045 -0.1760R =-1.6690 -1.1737 -0.9944 -0.8854 -0.7882 -1.09430 -0.7594 -0.2803 0.3166 -0.0940 0.58290 0 -0.3472 0.1832 0.1390 -0.10880 0 0 0.8020 0.5966 0.50220 0 0 0 -0.4714 -0.02290 0 0 0 0 0.4305test_qr3err =1.0e-10 *0.1518 0.0119 0.0120 0.2845-111 1.52 2.53 3.54五、总结与分析由上面的数值试验结果和绘图可知:比较之下,我自己写的houserhoilder程序误差是最大的,CGS与MGS的误差相差不大,系统自带的matlab软件速度是最快的,我自己写的houserhoilder程序还是最慢的。

数值代数实验报告(1)绝对经典

数值代数实验报告(1)绝对经典

数值代数实验报告Numerical Linear Algebra And ItsApplications学生所在学院:理学院学生所在班级:计算数学10-1学生姓名:戈东潮指导教师:于春肖教务处2012年12月实验一实验名称: Poisson 方程边值问题的五点差分格式 实验时间: 2012年12月13日 星期四 实验成绩: 一、实验目的:通过上机利用Matlab 数学软件实现Poisson 方程的边值问题的五点差分格式的线性代数方程组来认真解读Poisson 方程边值问题的具体思想与方法,使我们掌握得更加深刻,并且做到学习与使用并存,增加学习的实际动手性,不再让学习局限于书本和纸上,而是利用计算机学习来增加我们的学习兴趣。

二、实验内容:利用Poisson 方程来解下列问题:⎪⎩⎪⎨⎧∂∈=∈=∂∂+∂∂-ΩΩ),(,),(),(,sin sin )(y x y x u y x y x y u x u 0222222πππ 其中}{的边界是ΩΩ∂<<∈Ω,1,0),(y x y x 。

边值问题的解释是y x y x u ππsin sin ),(=(1)用例题中模型问题的方法取N=5,10(也可以取N=20)列出五点差分格式的线性代数方程组三、实验过程:系数矩阵A 的Matlab 实现函数如下:%文件名:poisson.m function T=poisson(N) for i=1:N a(i)=4; end b=diag(a); for j=1:N-1b(j,j+1)=-1;endfor j=2:Nb(j,j-1)=-1;endfor i=1:NI(i)=-1;endI=diag(I);for i=1:N:N*Nfor j=1:NT(i+j-1,i:i+N-1)=b(j,:);endendfor i=1:N:N*N-Nfor j=1:NT(i+j-1,i+N:i+N-1+N)=I(j,:);endendfor i=1+N:N:N*Nfor j=1:NT(i+j-1,i-N:i-1)=I(j,:);endend求解常数项b的Matlab实现函数如下:%函数名:constant.mfunction b=constant(N)n=N+1;h=1/n;i=1;for y=1/n:1/n:N/nfor x=1/n:1/n:N/nf(i)=2*pi*pi*sin(pi*x)*sin(pi*y);i=i+1;endendb=h*h*f;b=b';四、实验结果(总结/方案)在Matlab运行窗口输入>>A=poisson(5) Enter输出结果A =Columns 1 through 114 -1 0 0 0 -1 0 0 0 0 0-1 4 -1 0 0 0 -1 0 0 0 00 -1 4 -1 0 0 0 -1 0 0 00 0 -1 4 -1 0 0 0 -1 0 00 0 0 -1 4 0 0 0 0 -1 0-1 0 0 0 0 4 -1 0 0 0 -10 -1 0 0 0 -1 4 -1 0 0 00 0 -1 0 0 0 -1 4 -1 0 00 0 0 -1 0 0 0 -1 4 -1 00 0 0 0 -1 0 0 0 -1 4 00 0 0 0 0 -1 0 0 0 0 40 0 0 0 0 0 -1 0 0 0 -10 0 0 0 0 0 0 -1 0 0 00 0 0 0 0 0 0 0 -1 0 00 0 0 0 0 0 0 0 0 -1 00 0 0 0 0 0 0 0 0 0 -10 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 0 Columns 12 through 220 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 00 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 0 0 0 0 0 0 00 -1 0 0 0 0 0 0 0 0 00 0 -1 0 0 0 0 0 0 0 00 0 0 -1 0 0 0 0 0 0 0 -1 0 0 0 -1 0 0 0 0 0 04 -1 0 0 0 -1 0 0 0 0 0 -1 4 -1 0 0 0 -1 0 0 0 00 -1 4 -1 0 0 0 -1 0 0 00 0 -1 4 0 0 0 0 -1 0 00 0 0 0 4 -1 0 0 0 -1 0 -1 0 0 0 -1 4 -1 0 0 0 -10 -1 0 0 0 -1 4 -1 0 0 00 0 -1 0 0 0 -1 4 -1 0 00 0 0 -1 0 0 0 -1 4 0 00 0 0 0 -1 0 0 0 0 4 -10 0 0 0 0 -1 0 0 0 -1 40 0 0 0 0 0 -1 0 0 0 -10 0 0 0 0 0 0 -1 0 0 00 0 0 0 0 0 0 0 -1 0 0 Columns 23 through 250 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 00 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -1 0 0 0 -1 0 0 0 -1 0 0 0 -1 0 0 4 -1 0 -1 4 -10 -1 4 在运行窗口输入 >>b=constant(5) Enter 常数项的输出结果b =[ 0.1371 0.2374 0.2742 0.2374 0.1371 0.2374 0.4112 0.4749 0.4112 0.2374 0.2742 0.4749 0.5483 0.4749 0.2742 0.2374 0.4112 0.4749 0.4112 0.2374 0.1371 0.2374 0.2742 0.2374 0.1371]T所以A*u=b 的五点差分格式为j i j i j i j i j i j i f h u u u u u ,,,,,,211114=-----+-+ (i,j=1,2……)实验二实验名称: 用Jacobi 迭代法和SOR 迭代法求解方程组 实验时间: 2012年12月13日 星期四 实验成绩: 一、实验目的:通过上机利用Matlab 数学软件采用Jacobi 迭代法和SOR 迭代法来求解方程组,求解过程中了解两种迭代方法的基本思想与迭代过程,并分析两种迭代方法的收敛性和实用用以及他们的异同点。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

暨南大学本科实验报告专用纸课程名称数值代数成绩评定实验项目名称(填写实验所属章节名称) 指导教师刘娟实验项目编号0701******* 实验项目类型验证实验地点机房学生姓名王伟文学号2010051727学院信息科学技术学院数学系信息与计算科学专业2010级实验时间2011年9月 1 日~12月30日温度24℃第一章1.实验选题:(写出你在该项目所选的实验题目)第一章上机练习12.谈谈你对该算法的理解:(简单谈一下你是如何理解该算法的?)对算法的理解:先将84阶的矩阵A分解为一个下三角矩阵L和上三角矩阵U,先考虑下三角形方程组Ly=b,利用前代法解出y,再考虑上三角形方程组Ux=y;利用回代法解出x。

列主元gauss消去法Ax b PA LU Ly Pb Ux y=⇔===,,;3.实验内容(将实验程序及其实验结果粘贴,最好对程序各部分注释清楚,比如设置了哪些函数,这些函数的输入输出是什么,具有什么功能?)实验程序对矩阵A进行LU分解的程序function [ L,U ] = LUfac( A )for k=1:n-1A(k+1:n,k)=A(k+1:n,k)/A(k,k);A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-A(k+1:n,k)*A(k,k+1:n);endL=tril(A,0);for i=1:nL(i,i)=1;endU=triu(A,0);end利用前代法解出y值的程序function [ b ] = TSL( L,b )n=size(L,1);for j=1:n-1b(j)=b(j)/L(j,j);b(j+1:n)=b(j+1:n)-b(j)*L(j+1:n,j);endb(n)=b(n)/L(n,n);end利用回代法解出x值的程序function [ b ] = TSU( U,b )n=size(U,1);for j=n:-1:2b(j)=b(j)/U(j,j);b(1:j-1)=b(1:j-1)-b(j)*U(1:j-1,j);endb(1)=b(1)/U(1,1);end主函数程序(‘生成84阶的矩阵A’)A=eye(84);A=6*A;for i=2:84A(i,i-1)=8;A(i-1,i)=1;end(‘生成84乘1的矩阵b’)b=ones(84,1);b=b*15;b(1)=7;b(84)=14;[L,U]=LUfac(A);(‘调用函数LUfac对矩阵A进行分解’)y=TSL(L,b);(‘调用函数TSL求解Ly=b方程’)x=TSU(U,y);(‘调用函数TSU求解UX=b方程’)用MATLAB求解得到的结果x’ans =1.0e+008 *Columns 1 through 70.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 8 through 140.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 15 through 210.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 22 through 280.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 29 through 350.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 36 through 420.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 43 through 490.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 50 through 560.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 0.0000 Columns 57 through 63-0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 Columns 64 through 700.0000 -0.0000 0.0000 -0.0001 0.0002 -0.0003 0.0007 Columns 71 through 77-0.0013 0.0026 -0.0052 0.0105 -0.0209 0.0419 -0.0836 Columns 78 through 840.1665 -0.3303 0.6501 -1.2582 2.3487 -4.0263 5.3684列主元gauss消去法function [L,U,P]=Lufac(A)n=size(A,1);P=eye(n,n);for i=1:n-1[r,m]=max(abs(A(i:n,i)));m=m+i-1;A([i,m],:)=A([m,i],:);P([i,m])=P([m,i]);if A(i,i)~=0A(i+1:n,i)=A(i+1:n,i)/ A(i,i);A(i+1:n,i+1:n)= A(i+1:n,i+1:n)-A(i+1:n,i)*A(i,i+1:n); endendU=triu(A);L=tril(A,-1)+eye(n,n);endA=eye(84);A=6*A;for i=2:84A(i,i-1)=8; A(i-1,i)=1; endb=ones(84,1);b=b*15;b(1)=7;b(84)=14;[L,U,P]=Lufac(A); b=P*b;y=TSL(L,b);x=TSU(U,y);求解结果为x =1.0e+025 *0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000 -0.0000 0.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00000.0000-0.00010.0002-0.00040.0007-0.00150.0030-0.00600.0119-0.02380.0475-0.09460.1878-0.36960.7154-1.33552.2894-3.05251.教师评语、评分:(请认真对待每次数值代数实验,期末交实验报告,将计算实验成绩)2.实验选题:第三章上机习题1(3)。

3.谈谈你对该算法的理解:该问题可转化为解方程组Ax=b的最小二乘解可运用正则化方法求解(1)计算T=;C A Ad A b=,T(2)用平方根方法计算C的Cholesky分解:T=;C L L(3)求解三角方程组Ly d==和T L x y4.实验内容A=[1 4.9176 1.0 3.4720 0.9980 1.0 7 4 42 3 1 0;1 5.0208 1.0 3.5310 1.5000 2.0 7 4 62 1 1 0;1 4.5429 1.0 2.2750 1.1750 1.0 6 3 402 1 0;1 4.5573 1.0 4.0500 1.2320 1.0 6 3 54 4 1 0;1 5.0597 1.0 4.4550 1.1210 1.0 6 3 423 1 0;1 3.8910 1.0 4.4550 0.9880 1.0 6 3 562 1 0;1 5.8980 1.0 5.8500 1.2400 1.0 7 3 512 1 1;1 5.6039 1.0 9.5200 1.5000 0.0 6 3 32 1 1 0;1 15.4202 2.5 9.800 3.4200 2.0 10 5 42 2 1 1;1 14.4598 2.5 12.80 3.000 2.0 8 5 14 4 1 1;1 5.8282 1.0 6.4350 1.2250 2.0 63 32 1 1 0;1 5.3003 1.0 4.9883 1.5520 1.0 6 3 30 12 0;1 6.2712 1.0 5.5200 0.9750 1.0 5 2 30 1 2 0;1 5.9592 1.0 6.6660 1.1210 2.0 63 32 2 1 0;1 5.0500 1.0 5.0000 1.0200 0.0 52 46 4 1 1;1 5.6039 1.0 9.5200 1.5010 0.0 6 3 32 1 1 0;1 8.2464 1.5 5.1500 1.6640 2.0 8 4 50 4 1 0;1 6.6969 1.5 6.0920 1.4880 1.5 7 3 22 1 1 1;1 7.7841 1.5 7.1020 1.3760 1.0 6 3 172 1 0;1 9.0384 1.0 7.8000 1.5000 1.5 7 3 23 3 3 0;1 5.9894 1.0 5.5200 1.2560 2.0 6 3 40 4 1 1;1 7.5422 1.5 4.0000 1.6900 1.0 63 22 1 1 0;1 8.7951 1.5 9.8900 1.8200 2.0 8 4 50 1 1 1;1 6.0931 1.5 6.7265 1.6520 1.0 6 3 44 4 1 0;1 8.3607 1.5 9.1500 1.7770 2.0 8 4 48 1 1 1;1 8.1400 1.0 8.0000 1.5040 2.0 7 3 3 1 3 0;1 9.1416 1.5 7.3262 1.8310 1.5 8 4 31 4 1 01 12.000 1.5 5.0000 1.2000 2.0 6 3 30 3 1 1;];b=[25.9 29.5 27.9 25.9 29.9 29.9 30.9 28.9 84.9 82.9 35.9 31.5 31.0 30.9 30.0 28.9 36.9 41.9 40.5 43.9 37.5 37.9 44.5 37.9 38.9 36.9 45.8 41.0]'; C=A'*A;C=A'*A;b=[25.9 29.5 27.9 25.9 29.9 29.9 30.9 28.9 84.9 82.9 35.9 31.5 31.0 30.9 30.0 28.9 36.9 41.9 40.5 43.9 37.5 37.9 44.5 37.9 38.9 36.9 45.8 41.0]'; d=A'*b;[L]=Cholesky(C);//用平方根方法求Cholesky分解U=L';[y]=TSL(L,d);//用前代法求解Ly d=[x]=TSU(U,y);//用回代法求解T L x y=x =1.0e+002 *0.3307 + 0.0000i0.0082 - 0.0000i-0.0001 - 0.0000i-0.0031 + 0.0000i0.0199 + 0.0000i-0.0002 - 0.0000i-0.0036 + 0.0000i0.0163 - 0.0000i-0.0024 + 0.0000i0.0274 + 0.0000i-0.3980 + 0.0000i1.4399用此可知模型中参数的最小二乘结果为x =1.0e+002 *0.3307 + 0.0000i0.0082 - 0.0000i-0.0001 - 0.0000i-0.0031 + 0.0000i0.0199 + 0.0000i-0.0002 - 0.0000i-0.0036 + 0.0000i0.0163 - 0.0000i-0.0024 + 0.0000i0.0274 + 0.0000i-0.3980 + 0.0000i1.4399Cholesky分解函数function [ L ] = Cholesky( A )n=size(A,1);for k=1:n;A(k,k)=sqrt(A(k,k));A(k+1:n,k)=A(k+1:n,k)/A(k,k);for j=k+1:n;A(j:n,j)=A(j:n,k)-A(j:n,k)*A(j,k);endendL=tril(A);end前代法求解函数function [ b ] = TSL( L,b )n=size(L,1);for j=1:n-1b(j)=b(j)/L(j,j);b(j+1:n)=b(j+1:n)-b(j)*L(j+1:n,j);endb(n)=b(n)/L(n,n);end回代法求解函数function [ b ] = TSU( U,b )n=size(U,1);for j=n:-1:2b(j)=b(j)/U(j,j);b(1:j-1)=b(1:j-1)-b(j)*U(1:j-1,j);endb(1)=b(1)/U(1,1);end(将实验程序及其实验结果粘贴,最好对程序各部分注释清楚,比如设置了哪些函数,这些函数的输入输出是什么,具有什么功能?)暨南大学本科实验报告专用纸(附页)5.教师评语、评分:(请认真对待每次数值代数实验,期末交实验报告,将计算实验成绩)6. 实验选题:(写出你在该项目所选的实验题目) 第四章上机习题7. 谈谈你对该算法的理解:(简单谈一下你是如何理解该算法的?) (1)Jacobi 迭代11()Ax bA D L Ux D L U x D b --==--=++111111()k k k k x D L U x D b D Lx D Ux D b -----+=++=++ 1111()/i nk k k ii ij jijj ii j j i xb a x a xa -+==+=--∑∑(i=1:n )(2)Guass_Seidel 迭代11()()Ax bA D L Ux D L Ux D L b --==--=-+-1111111()()k k k k x D L Ux D L b D Lx D Ux D b -----++=-+-=++ 11111()/i nk k k ii ij jijj ii j j i xb a xa xa -++==+=--∑∑(i=1:n)(3)SOR 迭代111111()(1)k k k k k x L x D L b x D Lx D Ux D b ωωωω----++=+-=-+++ 11111(1)()/i nk k k k iii ij jijjii j j i xx b a xa xa ωω-++==+=-+--∑∑(i=1:n )8. 实验内容(将实验程序及其实验结果粘贴,最好对程序各部分注释清楚,比如设置了哪些函数,这些函数的输入输出是什么,具有什么功能?) (1)矩阵生成命令 A=eye(100); b=ones(100,1); b=b*(1/20000); A=(-2.01)*A; for i=2:100 A(i,i-1)=1; A(i-1,i)=1.01; endD=eye(100); D=-2.01*D; L=zeros(100); U=zeros(100); for i=2:100L(i,i-1)=-1;L(i-1,i)=0;endfor i=2:100U(i,i-1)=0;U(i-1,i)=-1.01;EndB=inv(D)*(L+U);[h]=eig(B);% 求矩阵B的特征值矩阵B的谱半径为0.9995(2)Jacobi迭代函数function [x]=Jacobi(A,b)n=size(A,1);x=zeros(n,1);y=zeros(n,1);for k=1:100 %迭代次数为100for i=1:ny(i)=(b(i)-A(1:i-1)*x(1:i-1)-A(i+1:n)*x(i+1:n))/A(i,i);endx=y;endend求解结果为x'ans =1.0e-004 *Columns 1 through 7-0.2488 0 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 8 through 140.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 15 through 210.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 22 through 280.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 29 through 350.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 36 through 420.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 43 through 490.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 50 through 560.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 57 through 630.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 64 through 700.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 71 through 770.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 78 through 840.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 85 through 910.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 92 through 980.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 99 through 1000.0000 0.0000(2)Guass-Seidel迭代function [x]=GuassSeidel(A,b)n=size(A,1);x=zeros(n,1);for k=1:100for i=1:nx(i)=(b(i)-A(1:i-1)*x(1:i-1)-A(i+1:n)*x(i+1:n))/A(i,i);endendend求解结果为x'ans =1.0e-004 *Columns 1 through 7-0.2488 0 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 8 through 140.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 15 through 210.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 22 through 280.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 29 through 350.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 36 through 420.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 43 through 490.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 50 through 560.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.00000.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 64 through 700.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 71 through 770.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 78 through 840.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 85 through 910.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 92 through 980.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 Columns 99 through 100(3)SOR 迭代function [x]=SOR(A,b)n=size(A,1);x=zeros(n,1);for k=1:100%迭代次数为100for i=1:nx(i)=(1-ω)*x(i)+ *(b(i)-A(1:i-1)*x(1:i-1)-A(i+1:n)*x(i+1:n))/A(i,i);1())B D L U ω-==+endendend求解结果为x'ans =1.0e+049 *Columns 1 through 71.1723 -3.1072 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 8 through 14-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 15 through 21-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047Columns 22 through 28-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 29 through 35-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 36 through 42-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 50 through 56-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 57 through 63-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 64 through 70-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 71 through 77-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 78 through 84-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 85 through 91-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 92 through 98-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 99 through 100-7.2047 -7.2047考虑ε=0.1的情况Jacobi和GuassSeidel迭代求结果为x'ans =1.0e-003 *Columns 1 through 7-0.2381 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 8 through 14-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 15 through 21-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 22 through 28-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 29 through 35-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 36 through 42-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 43 through 49-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000Columns 50 through 56-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 57 through 63-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 64 through 70-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 71 through 77-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 78 through 84-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 85 through 91-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 92 through 98-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 99 through 100-0.0000 -0.0000SOR迭代求解结果为x'ans =1.0e+049 *Columns 1 through 71.1723 -3.1072 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 8 through 14-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 15 through 21-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 22 through 28-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 29 through 35-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 36 through 42-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 43 through 49-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 50 through 56-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 57 through 63-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 64 through 70-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 71 through 77-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 78 through 84-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 85 through 91-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 92 through 98-7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 -7.2047 Columns 99 through 100-7.2047 -7.2047考虑ε=0.01的情况Jacobi和GuassSeidel迭代求结果为x'ans =Columns 1 through 7-0.0017 0 0 0 0 0 0 Columns 8 through 140 0 0 0 0 0 0 Columns 15 through 210 0 0 0 0 0 0 Columns 22 through 280 0 0 0 0 0 0 Columns 29 through 350 0 0 0 0 0 0 Columns 36 through 420 0 0 0 0 0 0 Columns 43 through 490 0 0 0 0 0 0 Columns 50 through 560 0 0 0 0 0 0 Columns 57 through 630 0 0 0 0 0 0 Columns 64 through 700 0 0 0 0 0 0 Columns 71 through 770 0 0 0 0 0 0 Columns 78 through 840 0 0 0 0 0 0 Columns 85 through 910 0 0 0 0 0 0 Columns 92 through 980 0 0 0 0 0 0 Columns 99 through 1000 0SOR迭代求解结果为x'ans =1.0e+017 *Columns 1 through 70.5371 -1.1729 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 15 through 21-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 22 through 28-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 29 through 35-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 36 through 42-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 43 through 49-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 50 through 56-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 57 through 63-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 64 through 70-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 71 through 77-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 78 through 84-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 85 through 91-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 92 through 98-2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 -2.0266 Columns 99 through 100-2.0266 -2.0266考虑ε=0.01的情况Jacobi和GuassSeidel迭代求结果为x'ans =Columns 1 through 7-0.0049 0 0 0 0 0 0 Columns 8 through 140 0 0 0 0 0 0 Columns 15 through 210 0 0 0 0 0 0 Columns 22 through 280 0 0 0 0 0 0 Columns 29 through 350 0 0 0 0 0 00 0 0 0 0 0 0 Columns 43 through 490 0 0 0 0 0 0 Columns 50 through 560 0 0 0 0 0 0 Columns 57 through 630 0 0 0 0 0 0 Columns 64 through 700 0 0 0 0 0 0 Columns 71 through 770 0 0 0 0 0 0 Columns 78 through 840 0 0 0 0 0 0 Columns 85 through 910 0 0 0 0 0 0 Columns 92 through 980 0 0 0 0 0 0 Columns 99 through 1000 0SOR迭代求解结果为x'ans =Columns 1 through 7-0.0049 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 8 through 14-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 15 through 21-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 22 through 28-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 29 through 35-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 36 through 42-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 43 through 49-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 50 through 56-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 57 through 63-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 64 through 70-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 71 through 77-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 78 through 84-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 85 through 91-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 92 through 98-0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 -0.0000 Columns 99 through 100-0.0000 -0.00009.教师评语、评分:(请认真对待每次数值代数实验,期末交实验报告,将计算实验成绩)。

相关文档
最新文档