1、实验项目类型课程名称计算方法 实验项目名 称 线性方程组的数值解法 验证 演示 综合 设计 其他指导教师胡小兵 成 绩 1实验目的：（1）高斯列主元消去法求解线性方程组的过程（2）熟悉用迭代法求解线性方程组的过程（3）设计出相应的算法，编制相应的函数子程序2实验内容分别用高斯列主元消去法 ,Jacobi 迭代法,Gauss-Saidel 迭代法,超松弛迭代法求解线性方程组7251039142304x3实验环境及实验文件存档名写出实验环境及实验文件存档名4.实验结果及分析输出计算结果，结果分析和小结等。解：1. 高斯列主元消去法：%用高斯列主元消去法解实验 1%高斯列主元消元法求解线性方程组 A
2、x=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,7m,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; 第 1 页，共 14
3、 页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;输出结果：2Jacobi 迭代法：第 2 页
4、，共 14 页%Jacobi 迭代法求解实验 1% A 为方程组的增广矩阵clc;A=2 10 0 -3 10;-3 -4 -12 13 5;1 2 3 -4 -2;4 14 9 -13 7MAXTIME=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:nif j=is=s+A(i,j)*x(j); endy(i)=(A(i,n+1)-s)/A(i,i);endendfor i=1:1:nmax
5、eps=max(0,abs(x(i)-y(i); endif maxepsMAXTIMEerror(超过最大迭代次数,退出);return;endend第 3 页，共 14 页输出结果：由于不收敛，故出现上述情况。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 7n,m=size(A);Maxtime=50;Eps=10E-5;x=zeros(1,n);disp(x=);for k=1:Maxtim
6、e disp(x); for i=1:ns=0.0;第 4 页，共 14 页for j=1:nif i=js=s+A(i,j)*x(j);endendx(i)=(A(i,n+1)-s)/A(i,i);endif sum(x-floor(x).2)=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=xformat short；第 6 页，共 14 页输出结果：由于不收敛，故出现
7、上述情况。第 7 页，共 14 页实验项目类型课程名称计算方法 实验项目名 称 插值方法 验证 演示 综合 设计 其他指导教师胡小兵 成 绩 1.实验目的：（1） 学会拉格朗日插值、牛顿插值等基本方法（2） 设计出相应的算法，编制相应的函数子程序（3） 会用这些函数解决实际问题2实验内容（1）设计拉格朗日插值算法，编制并调试相应的函数子程序（2）设计牛顿插值算法，编制并调试相应的函数子程序（3）给定函数四个点的数据如下：X 1.1 2.3 3.9 5.1Y 3.887 4.276 4.651 2.117试用拉格朗日插值确定函数在 x=2.101，4.234 处的函数值。4）已知 用牛顿插值公式
8、求 的近似值。， 3924153实验原理写出本次实验所用算法的算法步骤叙述或画出算法程序框图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=jt=t*(x-x0(i+1)/(x0(j+1)-x0(i+1);endends=s+t*y0(j+1);end第 1 页
9、，共 14 页s format short;运行结果：（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运行结果：第 2 页，共 14 页实际结果是 2.2
10、36066797749979.有一定的误差。总体来说还是很准确的第 3 页，共 14 页实验项目类型课程名称计算方法 实验项目名 称 数值微积分 验证 演示 综合 设计 其他指导教师胡小兵 成 绩 1实验目的：（1）学会复化梯形、复化辛浦生求积公式的应用（3）设计出相应的算法，编制相应的函数子程序（4）会用这些函数解决实际问题2实验内容（1）设计复化梯形公式求积算法，编制并调试相应的函数子程序（2）设计复化辛浦生求积算法，编制并调试相应的函数子程序（4）分别用复化梯形公式和复化辛浦生公式计算定积分 10sindx取 n=2，4，8，16，精确解为 0.94608313、 实验原理写出本次实验所
11、用算法的算法步骤叙述或画出算法程序框图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(
12、n)f=sym(sin(x)/x);第 1 页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);运行结果第 2 页实验项目类型课程名称计算方法 实验项目名 称 常微分方程的数值解法 验证 演示 综合 设计 其他指导教师胡小兵 成 绩 1实验目的：（1）学会四阶龙格-库塔方法的使用（2）设计出相应的算法
13、，编制相应的函数子程序（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.实验原理写出本次实验所用算法的算法步
14、骤叙述或画出算法程序框图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;第 2 页K4=h*eval(F);Y1(i+1)=Y1(i)+(K1+2*K2+2*K3+K4)/6;endX1disp(Y1);e
15、ndh=0.05，N=10时，运行结果：h=0.025,n=20 时，运行结果：第 3 页h=0.01，N=50 时，输出结果：第 4 页第 5 页第 6 页（2）微分方程为：v=9.8-0.003857v2,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;endtdisp(v);end运行结果：第 7 页第 8 页第 9 页第 10 页第 11 页第 12 页
