数值分析计算方法实验报告.docx
课程名称计算方法实验项目名 称线性方程组的数值解法实验项目类型验证演示综合设计其他指导教师胡小兵成绩V1. 实验目的：(1) 高斯列主元消去法求解线性方程组的过程(2) 熟悉用迭代法求解线性方程组的过程(3) 设计出相应的算法，编制相应的函数子程序2. 实验内容分别用高斯列主元消去法Jacobi迭代法,Gauss-Saidel迭代法，超松弛迭代法求解线性方程组_ 2100-3 -尤110 -3-4-1213=5123-4尤3-24149-13人1- 473. 实验环境及实验文件存档名写出实验环境及实验文件存档名4. 实验结果及分析输出计算结果，结果分析和小结等。解：1,高斯列主元消去法：%用高斯列主元消去法解实验1%高斯列主元消元法求解线性方程组Ax=b%A为输入矩阵系数,b为方程组右端系数%方程组的解保存在x变量中format long;A=2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13b=10,5,-2,7'm,n=size(A);if m=nerror('矩阵A的行数和列数必须相同')；return;endif m=size(b)error('b的大小必须和A的行数或A的列数相同')；return;endif rank(A)=rank(A,b)error('A矩阵的秩和增广矩阵的秩不相同，方程不存在唯一解');return;endc=n+1;A(:,c)=b;for k=1:n-1r,m=max(abs(A(k:n,k);m=m+k-1;if(A(m,k)=0)if(m=k)A(k m,:)=A(m k,:);endA(k+1:n, k:c)=A(k+1:n, k:c)-(A(k+1:n,k)/ A(k,k)*A(k, k:c); endendx=zeros(length(b),1);x(n)=A(n,c)/A(n,n);for k=n-1:-1:1x(k)=(A(k,c)-A(k,k+1:n)*x(k+1:n)/A(k,k);enddisp('X=');disp(x);format short;输出结果：2. Jacobi迭代法:%Jacobi迭代法求解实验1%A为方程组的增广矩阵clc;A=2 10 0 -3 10;-3 -4 -12 13 5;1 2 3 -4 -2;4 14 9 -13 7 MAXTIME=50;eps=1e-5;n,m=size(A);x=zeros(n,1);y=zeros(n,1);k=0;disp('迭代过程X的值情况如下:')disp('X=');while 1disp(x');for i=1:1:ns=0.0;for j=1:1:n if j=i s=s+A(i,j)*x(j);endy(i)=(A(i,n+1)-s)/A(i,i);endendfor i=1:1:nmaxeps=max(0,abs(x(i)-y(i);endif maxeps<=epsfor i=1:1:nx(i)=y(i);endreturn;endfor i=1:1:nx(i)=y(i);y(i)=0.0;endk=k+1;if k>MAXTIMEerror('超过最大迭代次数,退出')； return;endend由于不收敛，故出现上述情况。3. Gauss-Saidel 迭代法%Gauss_Seidel迭代法求解实验1%A为方程组的增广矩阵clc;format long;A=2 10 0 -3 10;-3 -4 -12 13 5;1 2 3 -4 -2;4 14 9 -13 7 n,m=size(A);Maxtime=50;Eps=10E-5;x=zeros(1,n);disp('x=');for k=1:Maxtimedisp(x);for i=1:ns=0.0;for j=1:n if i=j s=s+A(i,j)*x(j); end end x(i)=(A(i,n+1)-s)/A(i,i);endif sum(x-floor(x).A2)<Eps break;end;end;X=x;disp('迭代结果:');Xformat short;输出结果：因为不收敛，故出现上述情况。4. 超松弛迭代法：%SOR法求解实验1%w=1.45%方程组系数矩阵clc;A=2,10,0,-3;-3,-4,-12,13;1,2,3,-4;4,14,9,-13b=10,5,-2,7'w=1.45;Maxtime=100;Eps=1E-5;format long;n=length(A);k=0;x=ones(n,1);y=x;disp('迭代过程:')；disp('x=');while 1y=x;disp(x');for i=1:ns=b(i);for j=1:nifj=is=s-A(i,j)*x(j);endendif abs(A(i,i)<1E-10 I k>=Maxtimeerror(已达最大迭代次数或矩阵系数近似为0,无法进行迭代');return;ends=s/A(i,i);x(i)=(1-w)*x(i)+w*s;endif norm(y-x,inf)<Epsbreak;endk=k+1;enddisp('最后迭代结果:')；X=x'format short；由于不收敛，故出现上述情况。课程 名称计算方法实验项目名 称实验项目类型插值方法验证演示综合设计其他指导教师胡小兵成绩V1. 实验目的：(1) 学会拉格朗日插值、牛顿插值等基本方法(2) 设计出相应的算法，编制相应的函数子程序(3) 会用这些函数解决实际问题2. 实验内容(1) 设计拉格朗日插值算法，编制并调试相应的函数子程序(2) 设计牛顿插值算法，编制并调试相应的函数子程序(3) 给定函数四个点的数据如下：X1.12.33.95.1Y3.8874.2764.6512.117试用拉格朗日插值确定函数在x=2.101，4.234处的函数值。4)已知切=L，4 = 2,9 = 3,用牛顿插值公式求15的近似值。3. 实验原理写出本次实验所用算法的算法步骤叙述或画出算法程序框图4. 实验环境及实验文件存档名写出实验环境及实验文件存档名5. 实验结果及分析输出计算结果，结果分析和小结等。(1)拉格朗日法程序：function lagrange(x)format long;x0=1.1 2.3 3.9 5.1;y0=3.887 4.276 4.651 2.117;n=length(x0);s=0;for j=0:(n-1)t=1;for i=0:(n-1) if i=j t=t*(x-x0(i+1)/(x0(j+1)-x0(i+1);endends=s+t*y0(j+1);endsformat short;运行结果：CoBBajid TindowFile Edit H ebug Ilesktcip Wirnlow Help>> lagrange (2. 101)4. 14568709626432>> lagrange (4. 234)4.30074399459970>>1I(2)牛顿差值法：function Newton(x)format long;% 显示 15位x0=1 4 9;%x 的值y0=1 2 3;%y 的值x=5;%插值点n=max(size(x0);y=y0(1);%迭代初始值disp(y);s = 1;dx=y0;for i=1:n-1%构造差商表dx0=dx;for j=1:n-idx(j) = (dx0(j+1)-dx0(j)/(x0(i+j)-x0(j);enddf=dx(1);s=s*(x-x0(i);y=y+s*df;%算disp(y);end运行结果：|_lnl x|File Edi t Debug De zktop iVindow Help>> Nevrton(512.333333333333332.266666666666677、 1d实际结果是2.236066797749979,有一定的误差。总体来说还是很准确的2.实验内容(1)(2)(4)课程 名称计算方法实验项目名 称数值微积分实验项目类型验证演示综合设计其他指导教师胡小兵成绩V1.实验目的：(1)学会复化梯形、复化辛浦生求积公式的应用(3) 设计出相应的算法，编制相应的函数子程序(4) 会用这些函数解决实际问题设计复化梯形公式求积算法，编制并调试相应的函数子程序 设计复化辛浦生求积算法，编制并调试相应的函数子程序 分别用复化梯形公式和复化辛浦生公式计算定积分.sin x ,1 dx取 n=2, 4, 8, 16,精确解为 0.94608313、实验原理写出本次实验所用算法的算法步骤叙述或画出算法程序框图4.实验环境及实验文件存档名写出实验环境及实验文件存档名5. 实验结果及分析输出计算结果，结果分析和小结等。本实验是用MATLAB来完成。梯形法存档名为trap.m辛浦生法存档名为simpson.m 复化梯形公式求积算法：function T=trap(n)f=sym('sin(x)/x');a=0;b=1;h=(b-a)/n;T=0;for k=1:(n-1)x0=a+h*k;T=T+limit(f,x0);endT=h*(limit(f,a)+limit(f,b)/2+h*T;T=double(T);复化辛浦生法：function S=simpson(n) f=sym('sin(x)/x');a=0;b=1;h=(b-a)/(2火n);s1 = 0;s2 = 0;for k=1:nx0=a+h*(2*k-1);s1=s1+limit(f,x0);endfor k=1:(n-1)x0=a+h*2*k;s2=s2+limit(f,x0);endS=h*(limit(f,a)+limit(f,b)+4火s1+2火s2)/3;S=double(S);运行结果Jt： Cobb and Tindov_|nl x|1- CoBBand. Tindowa xiFile Edi + D ebug Desktop WiridowHelpFile Edi t Debug Desktop 里indo.Help>> trap (2)ans =0.93979328480618>> trap (4)ans =0.94451352166539>> trap (8)ans =0.94569086358270>> trap(16)ans =0.945985029934391>> simpson(2)ans =0.94608693395179>> simpson(4)ans =0.94608331088847>> simpson(8)ans 二0.94608308538495>> simpson(16)ans =0.946083071305561课程 名称计算方法实验项目名 称常微分方程的数值解法实验项目类型验证演示综合设计其他指导教师胡小兵成绩V1. 实验目的：(1) 学会四阶龙格-库塔方法的使用(2) 设计出相应的算法，编制相应的函数子程序(3) 会用这些函数解决实际问题2. 实验内容(1) 分别取h=0.05, N=10； h=0.025, N=20； h=0.01, N=50,用四阶龙格-库塔方法求解微分 方程初值 问题：y=-50y, y(0)=10(2) 某跳伞者在t=0时刻从飞机上跳出，假设初始时刻的垂直速度为0,且跳伞者垂直下落。 已知空气阻力为F=cv2，其中c为常数，v为垂直速度，向下方方向为正。写出此跳伞者的速度满足 的微分方程；若此跳伞者的质量为M=70kg,且已知c=0.27kg/m，利用四阶龙格-库塔公式计算t<=20s 的速度(取h=0.1s)3. 实验原理写出本次实验所用算法的算法步骤叙述或画出算法程序框图3、实验结果及分析输出计算结果，结果分析和小结等。(1)程序：function RK(h,n)F='-50火y'a=0;b=a+h*n;X1=a:h:b;Y1=zeros(1,n+1);Y1(1)=10;for i=1:nx=X1(i);y=Y1(i);K1=h*eval(F);x=x+h/2;y=y+K1/2;K2=h*eval(F);x=x;y=Y1(i)+K2/2;K3=h*eval(F);x=X1(i)+h;y=Y1(i)+K3;K4=h*eval(F);Y1(i+1)=Y1(i)+(K1+2火K2+2火K3+K4)/6;endX1disp(Y1);endh=0.05, N=10时，运行结果：h=0.025,n=20时，运行结果:(2)微分方程为：v=9.8-0.003857vA2,v(0)=0程序：function velosity(h,n)F='9.8-0.003857火y*y'a=0;b=a+h*n;t=a:h:b;v=zeros(1,n+1);v(1)=0;for i=1:nx=t(i);y=v(i);K1=h*eval(F);x=x+h/2;y=y+K1/2;K2=h*eval(F);x=x;y=v(i)+K2/2;K3=h*eval(F);x=t(i)+h;y=v(i)+K3;K4=h*eval(F);v(i+1)=v(i)+(K1+2火K2+2火K3+K4)/6;enddisp(v);end运行结果
