东南大学数值分析上机

精选优质文档-倾情为你奉上第一章一、题目设，其精确值为。（1）编制按从大到小的顺序，计算SN的通用程序。（2）编制按从小到大的顺序，计算SN的通用程序。（3）按两种顺序分别计算,并指出有效位数。（编制程序时用单精度）（4）通过本次上机题，你明白了什么？二、MATLAB程序N=input(请输入N(N1)：);AccurateValue=single(0-1/(N+1)-1/N+3/2)/2); %single使其为单精度Sn1=single(0); %从小到大的顺序 for a=2:N; Sn1=Sn1+1/(a2-1);endSn2=single(0); %从大到小的顺序 for a=2:N; Sn2=Sn2+1/(N-a+2)2-1);endfprintf(Sn的值 (N=%d)n,N);disp(_)fprintf(精确值 %fn,AccurateValue);fprintf(从大到小计算的结果 %fn,Sn1);fprintf(从小到大计算的结果 %fn,Sn2);disp(_)三、结果请输入N(N1)：100Sn的值 (N=100)_精确值 0.从大到小计算的结果 0.从小到大计算的结果 0._请输入N(N1)：10000Sn的值 (N=10000)_精确值 0.从大到小计算的结果 0.从小到大计算的结果 0._请输入N(N1)：Sn的值 (N=)_精确值 0.从大到小计算的结果 0.从小到大计算的结果 0._四、结果分析可以得出，算法对误差的传播又一定的影响，在计算时选一种好的算法可以使结果更为精确。从以上的结果可以看到从大到小的顺序导致大数吃小数的现象，容易产生较大的误差，求和运算从小数到大数算所得到的结果才比较准确。第二章一、题目（1）给定初值及容许误差，编制牛顿法解方程f(x)=0的通用程序。（2）给定方程,易知其有三个根由牛顿方法的局部收敛性可知存在当时，Newton迭代序列收敛于根x2*。试确定尽可能大的。试取若干初始值，观察当时Newton序列的收敛性以及收敛于哪一个根。（3）通过本上机题，你明白了什么？二、MATLAB程序 文件fx.mfunction Fx=fx(x) % 定义函数f(x)Fx=x3/3-x;文件dfx.mfunction Fx=dfx(x) % 定义导函数df(x)Fx=x2-1;% Newton法求方程的根%clearef=10-6; %这里取容许误差10-6k=0;x0=input(请输入 Xo的值:);disp(k Xk); %使用空格将其分隔开fprintf(0 %fn,x0); flag=1;while flag=1 & k=103 x1=x0-fx(x0)/dfx(x0); if abs(x1-x0)ef flag=0; end k=k+1; x0=x1;fprintf(%d %fn,k,x0); end %寻找最大的delta值%clear%flag=1;k=1;x0=0;while flag=1 delta=k*10-6; %delta与k有关 x0=delta; k=k+1; m=0; flag1=1; while flag1=1 & m=103 x1=x0-fx(x0)/dfx(x0); if abs(x1-x0)=10-6 %未小于给定误差时停止循环 flag=0; endendfprintf( delta 的最大值是 %fn,delta);三、结果1.运行search.m文件结果为： delta 的最大值为 0.即得最大的为0.，Newton迭代序列收敛于根=0的最大区间为（-0.，0.）。2. （1）区间上取-1000，-100，-50，-30，-10，-8，-7，-5，-3k Xk0 -10000.1 -6666.2 -4444.3 -2962.4 -1975.5 -1316.6 -877.7 -585.8 -390.9 -260.10 -173.11 -115.12 -77.13 -51.14 -34.15 -22.16 -15.17 -10.18 -6.19 -4.20 -3.21 -2.22 -1.23 -1.24 -1.25 -1.26 -1.k Xk0 -100.1 -66.2 -44.3 -29.4 -19.5 -13.6 -8.7 -5.8 -4.9 -2.10 -2.11 -1.12 -1.13 -1.14 -1.15 -1.k Xk0 -50.1 -33.2 -22.3 -14.4 -9.5 -6.6 -4.7 -3.8 -2.9 -1.10 -1.11 -1.12 -1.k Xk0 -30.1 -20.2 -13.3 -8.4 -6.5 -4.6 -2.7 -2.8 -1.9 -1.10 -1.11 -1.12 -1.k Xk0 -10.1 -6.2 -4.3 -3.4 -2.5 -1.6 -1.7 -1.8 -1.9 -1.k Xk0 -8.1 -5.2 -3.3 -2.4 -2.5 -1.6 -1.7 -1.8 -1.9 -1.k Xk0 -7.1 -4.2 -3.3 -2.4 -1.5 -1.6 -1.7 -1.8 -1.k Xk0 -5.1 -3.2 -2.3 -1.4 -1.5 -1.6 -1.7 -1.8 -1.k Xk0 -3.1 -2.2 -1.3 -1.4 -1.5 -1.6 -1.结果显示，以上初值迭代序列均收敛于-1.，即根。(2) 在区间即区间（-1，-0.）上取-0.，-0.8，-0.85，-0.9，-0.99，计算结果如下：k Xk0 -0.1 0.2 -0.3 0.4 -0.5 0.6 -0.7 1.8 1.9 1.10 1.11 1.12 1.13 1.k Xk0 -0.1 0.2 -5.3 -3.4 -2.5 -2.6 -1.7 -1.8 -1.9 -1.10 -1.k Xk0 0.1 -1.2 -1.3 -1.4 -1.5 -1.6 -1.k Xk0 -0.1 2.2 2.3 1.4 1.5 1.6 1.7 1.k Xk0 -0.1 32.2 21.3 14.4 9.5 6.6 4.7 3.8 2.9 1.10 1.11 1.12 1.13 1.计算结果显示，迭代序列局部收敛于-1.，即根，局部收敛于1.，即根。(3) 有上题可知，在区间（-0.，0.）上，在整个区间上均收敛于0，即根。(4) 在区间即区间（0.，1）上取0.，0.8，0.85，0.9，0.99，计算结果如下：k Xk0 0.1 -0.2 0.3 -0.4 0.5 -0.6 0.7 -1.8 -1.9 -1.10 -1.11 -1.12 -1.13 -1.k Xk0 0.1 -0.2 5.3 3.4 2.5 2.6 1.7 1.8 1.9 1.10 1.k Xk0 0.1 -1.2 -1.3 -1.4 -1.5 -1.6 -1.k Xk0 0.1 -2.2 -2.3 -1.4 -1.5 -1.6 -1.7 -1.k Xk0 0.1 -32.2 -21.3 -14.4 -9.5 -6.6 -4.7 -3.8 -2.9 -1.10 -1.11 -1.12 -1.13 -1.计算结果显示，迭代序列局部收敛于-1.，即根，局部收敛于1.，即根。(5) 区间上取100，60，20，10，7，6，4，3，1.5，计算结果如下：k Xk0 100.1 66.2 44.3 29.4 19.5 13.6 8.7 5.8 4.9 2.10 2.11 1.12 1.13 1.14 1.15 1.k Xk0 60.1 40.2 26.3 17.4 11.5 8.6 5.7 3.8 2.9 2.10 1.11 1.12 1.13 1.14 1.k Xk0 20.1 13.2 8.3 6.4 4.6 2.7 1.8 1.9 1.10 1.11 1.k Xk0 10.1 6.2 4.3 3.4 2.5 1.6 1.7 1.8 1.9 1.k Xk0 7.1 4.2 3.3 2.4 1.5 1.6 1.7 1.8 1.k Xk0 6.1 4.2 2.3 2.4 1.5 1.6 1.7 1.8 1.k Xk0 4.1 2.2 2.3 1.4 1.5 1.6 1.7 1.k Xk0 3.1 2.2 1.3 1.4 1.5 1.6 1.k Xk0 1.1 1.2 1.3 1.4 1.5 1.结果显示，以上初值迭代序列均收敛于1.，即根。四、结果分析综上所述：区间收敛于-1.73205， 区间局部收敛于1.73205,局部收敛于-1.73205，区间收敛于0，区间类似于区间，收敛于1.73205。通过本上机题，明白了对于多根方程，Newton法求方程根时，迭代序列收敛于某一个根有一定的区间限制，在一个区间上，可能会局部收敛于不同的根。第三章 一、题目列主元Gauss消去法对于某电路的分析，归结为求解线性方程组。其中 （1）编制解n阶线性方程组的列主元高斯消去法的通用程序；（2）用所编程序线性方程组，并打印出解向量，保留5位有效数；二、MATLAB程序% 列主元Gauss消去法求解线性方程组%参数输入n=input(请输入矩阵 A的阶数: n=); %输入线性方程组阶数nb=zeros(1,n);A=input(请输入矩阵 A:);b(1,:)=input(请输入行向量 b:); %输入行向量bb=b; %得到列向量b C=A,b; %得到增广矩阵%列主元消去得上三角矩阵for i=1:n-1 maximum,index=max(abs(C(i:n,i); %将最大元素位置放在index行中 index=index+i-1; T=C(index,:); %T作为一个中转站，交换两行 C(index,:)=C(i,:); C(i,:)=T; for k=i+1:n %列主元消去 if C(k,i)=0 C(k,:)=C(k,:)-C(k,i)/C(i,i)*C(i,:); end endend% 回代求解 %x=zeros(n,1);x(n)=C(n,n+1)/C(n,n);for i=n-1:-1:1 x(i)=(C(i,n+1)-C(i,i+1:n)*x(i+1:n,1)/C(i,i);endA=C(1:n,1:n); disp(上三角矩阵为：)for k=1:n fprintf(%f ,A(k,:); fprintf(n);enddisp(方程的解为：);fprintf(%.5gn,x); %以5位有效数字输出结果三、结果请输入矩阵 A的阶数: n=9请输入矩阵 A:31 -13 0 0 0 -10 0 0 0;-13 35 -9 0 -11 0 0 0 0;0 -9 31 -10 0 0 0 0 0;0 0 -10 79 -30 0 0 0 -9;0 0 0 -30 57 -7 0 -5 0;0 0 0 0 -7 47 -30 0 0;0 0 0 0 0 -30 41 0 0;0 0 0 0 -5 0 0 27 -2;0 0 0 -9 0 0 0 -2 9请输入行向量 b: -15 27 -23 0 -20 12 -7 7 10上三角矩阵为：31. -13. 0. 0. 0. -10. 0. 0. 0. 0. 29. -9. 0. -11. -4. 0. 0. 0. 0. 0. 28. -10. -3. -1. 0. 0. 0. 0. 0. 0. 75. -31. -0. 0. 0. -9. 0. 0. 0. 0. 44. -7. 0. -5. -3. 0. 0. 0. 0. -0. 45. -30. -0. -0. 0. 0. 0. 0. -0. -0. 21. -0. -0. 0. 0. 0. 0. -0. -0. 0. 26. -2. 0. 0. 0. 0. -0. -0. 0. 0. 27.方程的解为：-0.289230.34544-0.71281-0.22061-0.43040.15431-0.0.20105第四章：多项式插值与函数最佳逼近一、题目：（1）编制求第一型3次样条插值函数的通用程序; （2） 已知汽车曲线型值点的数据如下：0123456789102.513.304.044.705.225.545.785.405.575.705.80端点条件为，。用所编制程序求车门的3次样条插值函数,并打印出。二、MATLAB程序n=input(Input n: n=); %n=10n=n+1;x=zeros(1,n); %x用来存储xy=zeros(1,n); %y用来存储yx(1,:)=input(输入x:); %分别输入x,yy(1,:)=input(输入y:);dy0=input(输入y(0)的一阶导数:); %输入边界条件dyn=input(输入y(n)的一阶导数:);d=zeros(n,1);h=zeros(1,n-1); %h为两个点之间的距离f1=zeros(1,n-1);f2=zeros(1,n-2);for i=1:n-1 h(i)=x(i+1)-x(i); f1(i)=(y(i+1)-y(i)/h(i); %一阶差商endfor i=2:n-1 f2(i)=(f1(i)-f1(i-1)/(x(i+1)-x(i-1); %求二阶差商 d(i)=6*f2(i) ;endd(1)=6*(f1(1)-dy0)/h(1); %补充d1和dnd(n)=6*(dy-f1(n-1)/h(n-1);A=zeros(n);miu=zeros(1,n-2);lamda=zeros(1,n-2);for i=1:n-2 miu(i)=h(i)/(h(i)+h(i+1); lamda(i)=1-miu(i);endA(1,2)=1;A(n,n-1)=1;for i=1:n A(i,i)=2;endfor i=2:n-1 A(i,i-1)=miu(i-1); A(i,i+1)=lamda(i-1);endM=Ad;syms x;for i=1:n-1 Sx(i)=collect(y(i)+(f1(i)-(M(i)/3+M(i+1)/6)*h(i)*(x-x(i)+M(i)/2*(x-x(i)2+(M(i+1)-M(i)/(6*h(i)*(x-x(i)3); Sx(i)=vpa(Sx(i),4);endS=zeros(1,n-1);for i=1:n-1 x=x(i)+0.5;S(i)=y(i)+(f1(i)-(M(i)/3+M(i+1)/6)*h(i)*(x-x(i)+M(i)/2*(x-x(i)2+、(M(i+1)-M(i)/(6*h(i)*(x-x(i)3;end%输出结果disp(S(x)=);for i=1:n-1 fprintf( %s (%d,%d)n,char(Sx(i),x(i),x(i+1) disp( )enddisp(S(i+0.5)disp( i x(i+0.5) S(i+0.5)for i=1:n-1 fprintf( %d %0.4f %0.4fn,x(i)+0.5,S(i)end三、结果1.数据输入输入n:n=10输入x:0 1 2 3 4 5 6 7 8 9 10输入y:2.51 3.30 4.04 4.70 5.22 5.54 5.78 5.40 5.57 5.70 5.80输入y(0)的一阶导数:0.8输入y(n)的一阶导数:0.22.计算结果x属于区间0,1时；S(x)=2.51+0.8(x)-0.(x)(x)-0.(x)(x)(x)x属于区间1,2时；S(x)=3.3+0.(x-1)-0.(x-1)(x-1)-0.(x-1)(x-1)(x-1)x属于区间2,3时；S(x)=4.04+0.(x-2)-0.(x-2)(x-2)-0.(x-2)(x-2)(x-2)x属于区间3,4时；S(x)=4.7+0.(x-3)-0.(x-3)(x-3)-0.(x-3)(x-3)(x-3)x属于区间4,5时；S(x)=5.22+0.(x-4)-0.(x-4)(x-4)+0.(x-4)(x-4)(x-4)x属于区间5,6时；S(x)=5.54+0.(x-5)+0.(x-5)(x-5)-0.(x-5)(x-5)(x-5)x属于区间6,7时；S(x)=5.78-0.(x-6)-0.(x-6)(x-6)+0.(x-6)(x-6)(x-6)x属于区间7,8时；S(x)=5.4-0.(x-7)+0.(x-7)(x-7)-0.(x-7)(x-7)(x-7)x属于区间8,9时；S(x)=5.57+0.(x-8)-0.(x-8)(x-8)+0.(x-8)(x-8)(x-8)x属于区间9,10时；S(x)=5.7+0.(x-9)-0.(x-9)(x-9)+0.(x-9)(x-9)(x-9)ix(i+0.5)S(i+0.5)10.52.21.53.32.54.43.54.54.55.65.55.76.55.87.55.98.55.109.55.四、体会通过这次对三次样条插值函数的编程，对三次样条插值的计算有了更加深刻的理解。对复杂公式编程实现能力有了一定的提高。专心-专注-专业第五章 数值积分与数值微分一、题目：（1）给定积分。取初始步长和，及精度。应用复化梯形公式，采用逐次二分步长的方法，并应用外推思想编制计算通用程序。计算至相邻两次近似值之差的绝对值不超过为止。（2）用所编程序计算积分，取。二、MATLAB程序functionBEST,R=GET(f,a,b,c,d,m,n,epsilon)R=;COUNT=1;R(COUNT,1)=getT(f,a,b,c,d,2COUNT*m,2COUNT*n);COUNT=2;R(COUNT,1)=getT(f,a,b,c,d,2COUNT*m,2COUNT*n);R(COUNT-1,2)=4/3*R(COUNT,1)-1/3*R(COUNT-1,1);R1=R(1,1);R2=R(2,1);co=4/3 -1/3;16/15 -1/15;64/63 -1/63;while abs(R1-R2)epsilon COUNT=COUNT+1; R(COUNT,1)=getT(f,a,b,c,d,2COUNT*m,2COUNT*n); for i=1:size(co,1) if COUNTi R(COUNT-i,i+1)=co(i,1)*R(COUNT-i+1,i)+co(i,2)*R(COUNT-i,i); end end i=size(co,1); while size(R,1)-i-10 i=i-1; end R1=R(COUNT-i,i); R2=R(COUNT-i+1,i); BEST=R2;endendfunction R=getT(f,a,b,c,d,m,n)h=(b-a)/m; %h为步长k=(d-c)/n;x=a+(0:m)*h;y=c+(0:n)*k;R=0;for i=1:m for j=1:n R=R+h*k/4*(f(x(i),y(j)+f(x(i),y(j+1)+f(x(i+1),y(j)+f(x(i+1),y(j+1); endendendfunction exp5f=myfun; %调用被积函数a=0; %abcd分别为积分上下限b=pi/6;c=0;d=pi/3; m=1; n=1;epsilon=0.5e-5; %精度BEST,R=GET(f,a,b,c,d,m,n,epsilon);fprintf(Rult is %0.7f,detail:n,BEST)fprintf(kttt2ktttT(2k)ttts(2k)tttC(2k)tttR(2k)n)for k=1:size(R,1) fprintf(it%10.0f,k,2k); str=; for i=1:size(R,1)-k+1 if i=size(R,2) str=sprintf(%stt %0.7f,str,R(k,i); end end fprintf(%sn,str);endendfunction R=myfun(x,y) %被积函数R=tan(x2+y2);end三、结果Rult is 0.,detail:k2kT(2k)S(2k)C(2k)R(2k)120.0.0.0.240.0.0.0.380.0.0.0.4160.0.0.5320.0.6640.四、体会通过求解梯形公式、Simpson公式、科特斯公式、Romberg公式，可以在达到同样精度的前提下大大节省计算量。不过此题的编程还是相对复杂的。