《数值分析》实验报告书.

1、数值分析实验报告实验一、误差分析误差问题是数值分析的基础，又是数值分析中一个困难的课题。在实际计算中，如 果选用了不同的算法，由于舍入误差的影响，将会得到截然不同的结果。因此，选取算 法时注重分析舍入误差的影响，在实际计算中是十分重要的。同时，由于在数值求解过 程中用有限的过程代替无限的过程会产生截断误差，因此算法的好坏会影响到数值结果的精度。一、实验目的1、通过上机编程，复习巩固以前所学程序设计语言及上机操作指令；2、通过上机计算，了解误差、绝对误差、误差界、相对误差界的有关概念;3、通过上机计算，了解舍入误差所引起的数值不稳定性。、实验任务对n =0,1,2, ,20,计算定积分yn1 n

2、0dx 算法1：利用递推公式1取 y0 =0n =1,2,20,=1 n6-l n 5 0.182322. 0 x 5算法2:利用递推公式1 1y5T5yn= 20,19,1.注意到1 1x20X0x 5x20dx11 1取y20(八 0.008730.20 105126思考：从计算结果看，哪个算法是不稳定的，哪个算法是稳定的算法1：t=log(6.0)-log(5.0); n=0;y=zeros(1,21); y(1)=t;for k=2:21y(k)=1/k-5*y(k-1); n=n+1;endy(1:6) y(7:21)运行结果：ans = 0.1823 -0.4116 2.3914

3、-11.7069 58.7343 -293.5049算法2:y=zeros(21,1);n=1;y仁(1/105+1/126)/20;for k=21:-1:2y(k-1)=1/(5*k)-y(k)/5;n=n+1;end运行结果：y =0.0884 0.0580 0.0431 0.0343 0.0285 0.0243 0.0212 0.0188 0.01690.00890.01540.01410.01300.01200.01120.01050.00990.00930.0081 0.0095 0由数据对比可知，算法2较为稳定。实验二、插值法插值法是函数逼近的一种重要方法，它是数值积分、微分方程

4、数值解等数值计算的 基础与工具，其中多项式插值是最常用和最基本的方法。拉格朗日插值多项式的优点是 表达式简单明确，形式对称，便于记忆，它的缺点是如果想要增加插值节点，公式必须 整个改变，这就增加了计算工作量。而牛顿插值多项式对此做了改进，当增加一个节点 时只需在原牛顿插值多项式基础上增加一项，此时原有的项无需改变，从而达到节省计 算次数、节约存储单元、应用较少节点达到应有精度的目的。一、实验目的1、理解插值的基本概念，掌握各种插值方法，包括拉格朗日插值和牛顿插值等，注意 其不同特点；2、通过实验进一步理解并掌握各种插值的基本算法。二、实验任务1、已知函数表xi 0.56160 0.56280

5、0.56401 0.56521y 0.82741 0.82659 0.82577 0.82495用二次拉格朗日插值多项式求 x二0.5635时的函数近似值。2、已知函数表xi 0.40.550.650.80.9y 0.41075 0.578150.69675 0.888111.02652用牛顿插值多项式求 2(0.596)和N4(0.895)。1.fun ctio n y,R=lagra nzi(X,Y,x,M) x=0.5635;M=2;X=0.56160,0.56280,0.56401,0.56521;Y=0.82741,0.82659,0.82577,0.82495; n=length(

6、X); m=length(x);for i=1:mz=x(i);s=0.0;for k=1:n p=1.0;q1=1.0; c1=1.0;for j=1:n if j=k p=p*(z-X(j)/(X(k)-X(j);end q1=abs(q1*(z-X(j); c1=c1*j; end s=p*Y(k)+s;end y(i)=s;end R=M.*q1./c1; 运行结果： ans = 0.82612. N3(0.596) function y,R= newcz(X,Y,x,M) x=0.596;M=3;X=0.4,0.65,0.9; Y=0.41075,0.69675,1.02652; n

7、=length(X); m=length(x);for t=1:mz=x(t); A=zeros(n,n);A(:,1)=Y;s=0.0; p=1.0; q1=1.0; c1=1.0;for j=2:nfor i=j:n A(i,j)=(A(i,j-1)-A(i-1,j-1)/(X(i)-X(i-j+1);end q1=abs(q1*(z-X(j-1);c1=c1*j;endC=A(n,n);q1=abs(q1*(z-X(n);for k=(n-1):-1:1C=conv(C,poly(X(k);d=length(C); C(d)=C(d)+A(k,k); endy(k)= polyval(C

8、, z);endR=M*q1/c1;运行结果：ans = 0.6313N4(0.895 ) function y,R= newcz(X,Y,x,M)x=0.895;M=4;X=0.4,0.55,0.65,0.8,0.9;Y=0.41075,0.57815,0.69675,0.88811,1.02652; n=length(X); m=length(x);for t=1:mz=x(t); A=zeros(n,n);A(:,1)=Y;s=0.0; p=1.0; q1=1.0; c1=1.0;for j=2:nfor i=j:nA(i,j)=(A(i,j-1)-A(i-1,j-1)/(X(i)-X(

9、i-j+1);end q1=abs(q1*(z-X(j-1);c1=c1*j; endC=A(n,n);q1=abs(q1*(z-X(n);for k=(n-1):-1:1C=conv(C,poly(X(k);d=length(C); C(d)=C(d)+A(k,k); endy(k)= polyval(C, z);endR=M*q1/c1; 运行结果： ans = 1.0194实验三、解线性方程组的直接法解线性方程组的直接法是指经过有限步运算后能求得方程组精确解的方法。但由于实际计算中舍入误差是客观存在的，因而使用这类方法也只能得到近似解。目前较实用 的直接法是古老的高斯消去法的变形，即主元

10、素消去法及矩阵的三角分解法。引进选主 元的技巧是为了控制计算过程中舍入误差的增长，减少舍入误差的影响。一般说来，列 主元消去法及列主元三角分解法是数值稳定的算法，它具有精确度较高、计算量不大和 算法组织容易等优点，是目前计算机上解中、小型稠密矩阵方程组可靠而有效的常用方 法。、实验目的1、了解求线性方程组的直接法的有关理论和方法；2、会编制列主元消去法、LU分解法的程序；3、通过实际计算，进一步了解各种方法的优缺点，选择合适的数值方法。、实验任务1、 用列主元高斯消去法求解方程组0.101X1 +2.304X2 +3.555X3 =1.183 1.347X1 +3.712X2 +4.623X3

11、 = 2.137. 2.835xj +1.072X2 +5.643X3 =3.0352、用矩阵直接三角分解法求解方程组 Ax二b，其中一12-128 1547-2A =-37956-1283 一2714b =111.主程序:fun ctio n RA,RB, n,X=liezhu(A,b)B=A b; n=length(b); RA=rank(A);(求矩阵的秩)RB=ra nk(B);zhica=RB-RA; if zhica0, disp（请注意：因为RA=RB所以此方程组无解.）return endif RA=RBif RA=ndisp( 请注意：因为 RA=RB=，n 所以此方程组有唯

12、一解 .) X=zeros(n,1); C=zeros(1,n+1);for p= 1:n-1 Y,j=max(abs(B(p:n,p); C=B(p,:); B(p,:)= B(j+p-1,:); B(j+p-1,:)=C;for k=p+1:nm= B(k,p)/ B(p,p);B(k,p:n+1)= B(k,p:n+1)-m* B(p,p:n+1); end endb=B(1:n,n+1);A=B(1:n,1:n); X(n)=b(n)/A(n,n);for q=n-1:-1:1X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)/A(q,q);endelsedisp(

13、请注意：因为 RA=RB，n 所以此方程组有无穷多解 .) end end计算程序：A=0.101 2.304 3.555;-1.347 3.712 4.623;-2.835 1.072 5.643; b=1.183;2.137;3.035; RA,RB,n,X=liezhu(A,b) 运行结果： ans = 求矩阵的秩请注意：因为 RA=RB=，n 所以此方程组有唯一解 .RA = 3RB = 3X = -0.39820.01380.33512.程序：function X=LUjfcz(A,b) n,n =size(A);X=zeros(n,1);Y=zeros(n,1);C=zeros(1

14、,n);r=1:n;for p=1:n-1max1,j=max(abs(A(p:n,p);C=A(p,:);A(p,:)= A(j+p-1,:);A(j+p-1,:)=C;g=r(p);r(p)= r(j+p-1);r(j+p-1)=g;if A(p,p)=0disp(A 是奇异阵 , 方程组无唯一解 ); break;endfor k=p+1:nH= A(k,p)/A(p,p);A(k,p) = H;A(k,p+1:n)=A(k,p+1:n)- H* A(p,p+1:n);end endY(1)=b(r(1);for k=2:nY(k)= b(r(k)- A(k,1:k-1)* Y(1:k-

15、1);endX(n)= Y(n)/ A(n,n);for i=n-1:-1:1X(i)= (Y(i)- A(i, i+1:n) * X (i+1:n)/ A(i,i); endEnd 计算程序：A=1,2,-12,8;5,4,7,-2;-3,7,9,5;6,-12,-8,3; b=27;4;11;49;X=LUjfcz(A,b) 运行结果：X = 3.0000-2.00001.00005.0000实验四、解线性方程组的迭代法解线性方程组的迭代法是用某种极限过程去逐步逼近线性方程组精确解的方法，即是从一个初始向量x(0)出发，按照一定的迭代格式产生一个向量序列x(k)，使其收敛到 方程组Ax=b

16、的解。迭代法的优点是所需计算机存储单元少，程序设计简单，原始系数 矩阵在计算过程中始终不变等。但迭代法存在收敛性及收敛速度问题。迭代法是解大型 稀疏矩阵方程组的重要方法。一、实验目的1、熟悉迭代法的有关理论和方法；2、会编制雅可比迭代法、高斯-塞德尔迭代法的程序；3、注意所用方法的收敛性及其收敛速度问题。二、实验任务1、用雅可比迭代法解方程组X 十 2x2 - 2x3 = 7* Xi * X? + X3 = 2.2x 2x2 +x3 = 5注意：若用高斯-塞德尔迭代法则发散。2、用高斯-塞德尔迭代法解方程组x-i 0.9x20.9x3 二 1.9* 0.9x x2 +0.9x3 = 2.0.0

17、9捲 + 0.9x2 + x3 =1.7注意：若用雅可比迭代法则发散。1.主程序：fun ction X=jacdd(A,b,X0,P,wucha,max1)n m=size(A);for k=1:max1for j=1:mX(j)=(b(j)-A(j,1:j-1,j+1:m)*X0(1:j-1,j+1:m)/A(j,j); endXdjwcX=norm(X-X0,P); xdwcX=djwcX/(norm(X,P)+eps); X0=X;if (djwcXwucha)&(xdwcXwucha)&(xdwcXwucha)disp( 请注意：雅可比迭代次数已经超过最大迭代次数 max1 )End

18、计算程序：A=1 2 -2;1 1 1;2 2 1;b=7;2;5; X0=0 0 0;X=jacdd(A,b,X0,inf,0.01,100) 运行结果：k = 1X = 7 2 5 k = 2X = 13 -10 -13 k = 3X = 12 -1k = 4X = 12 -1主程序：function X=gsdddy(A,b,X0,P,wucha,max1)D=diag(diag(A);U=-triu(A,1);L=-tril(A,-1); dD=det(D);if dD=0disp( 请注意：因为对角矩阵 D 奇异，所以此方程组无解 .)elsedisp( 请注意：因为对角矩阵 D 非

19、奇异，所以此方程组有解 .) iD=inv(D-L); B2=iD*U;f2=iD*b;jX=Ab; X=X0;n m=size(A);for k=1:max1X1= B2*X+f2; djwcX=norm(X1-X,P); xdwcX=djwcX/(norm(X,P)+eps);if (djwcXwucha)|(xdwcXwucha)returnelsek,X1,k=k+1;X=X1;endendif (djwcXwucha)|(xdwcXwucha)disp( 请注意：高斯 - 塞德尔迭代收敛 , 此 A 的分解矩阵 D,U,L 和方 程组的精确解 jX 和近似解 X 如下： )elsed

20、isp( 请注意：高斯 - 塞德尔迭代的结果没有达到给定的精度，并且迭代次数已经超过 最大迭代次数maxi,方程组的精确解jX和迭代向量X如下：)X=X;jX=jXendendX=X;D,U,L,jX=jX计算程序：A=1 0.9 0.9;0.9 1 0.9;0.9 0.9 1;b=1.9;2.0;1.7;X0=0 0 0;X=gsdddy(A,b,X0,inf, 0.001,100)运行结果：k =1 ans =1.9000k =2ans =1.8829k = 3ans =1.8457k =4ans =1.7949k =5 ans =1.7356k =6ans =1.6718k = 7ans =1.6065k = 8ans =1.5419k =9ans =1.4796k =10ans =1.4208k =11ans =1.3661k =12ans =1.3160k =13ans =1.2706k =14ans =1.2300k =15ans =1.19390.2900 -0.27100.5493 -0.48900.7789 -0.66220.9805