数值分析计算方法实验报告Word下载.docx
《数值分析计算方法实验报告Word下载.docx》由会员分享，可在线阅读，更多相关《数值分析计算方法实验报告Word下载.docx（29页珍藏版）》请在冰点文库上搜索。
%高斯列主元消元法求解线性程组Ax=b
%A为输入矩阵系数,b为程组右端系数
%程组的解保存在x变量中
formatlong;
A=[2,10,0,-3;
-3,-4,-12,13;
1,2,3,-4;
4,14,9,-13]
b=[10,5,-2,7]'
[m,n]=size（A）;
ifm~=n
error（'
矩阵A的行数和列数必须相同'
）;
return;
end
ifm~=size（b）
b的大小必须和A的行数或A的列数相同'
ifrank（A）~=rank（[A,b]）
A矩阵的秩和增广矩阵的秩不相同,程不存在唯一解'
c=n+1;
A（:
c）=b;
fork=1:
n-1
[r,m]=max（abs（A（k:
n,k）））;
m=m+k-1;
if（A（m,k）~=0）
if（m~=k）
A（[km],:
）=A（[mk],:
end
A（k+1:
n,k:
c）=A（k+1:
c）-（A（k+1:
n,k）/A（k,k））*A（k,k:
c）;
x=zeros（length（b）,1）;
x（n）=A（n,c）/A（n,n）;
fork=n-1:
-1:
1
x（k）=（A（k,c）-A（k,k+1:
n）*x（k+1:
n））/A（k,k）;
disp（'
X='
disp（x）;
formatshort;
输出结果：
2．Jacobi迭代法：
%Jacobi迭代法求解实验1
%A为程组的增广矩阵
clc;
A=[2100-310;
-3-4-12135;
123-4-2;
4149-137]
MAXTIME=50;
eps=1e-5;
[n,m]=size（A）;
x=zeros（n,1）;
y=zeros（n,1）;
k=0;
disp（'
迭代过程X的值情况如下:
'
）
while1
disp（x'
fori=1:
1:
n
s=0.0;
forj=1:
ifj~=i
s=s+A（i,j）*x（j）;
y（i）=（A（i,n+1）-s）/A（i,i）;
fori=1:
maxeps=max（0,abs（x（i）-y（i）））;
ifmaxeps<
=eps
x（i）=y（i）;
y（i）=0.0;
k=k+1;
ifk>
MAXTIME
超过最大迭代次数,退出'
由于不收敛，故出现上述情况。
3.Gauss--Saidel迭代法
%Gauss_Seidel迭代法求解实验1
Maxtime=50;
Eps=10E-5;
x=zeros（1,n）;
x='
Maxtime
disp（x）;
ifi~=j
x（i）=（A（i,n+1）-s）/A（i,i）;
ifsum（（x-floor（x））.^2）<
Eps
break;
end;
end;
X=x;
迭代结果:
X
因为不收敛，故出现上述情况。
4.超松弛迭代法：
%SOR法求解实验1
%w=1.45
%程组系数矩阵
w=1.45;
Maxtime=100;
Eps=1E-5;
n=length（A）;
k=0;
x=ones（n,1）;
y=x;
迭代过程:
while1
y=x;
s=b（i）;
s=s-A（i,j）*x（j）;
ifabs（A（i,i））<
1E-10|k>
=Maxtime
已达最大迭代次数或矩阵系数近似为0,无法进行迭代'
s=s/A（i,i）;
x（i）=（1-w）*x（i）+w*s;
ifnorm（y-x,inf）<
最后迭代结果:
X=x'
formatshort；
插值法
1.实验目的：
（1）学会拉格朗日插值、牛顿插值等基本法
（2）设计出相应的算法，编制相应的函数子程序
（3）会用这些函数解决实际问题
（1）设计拉格朗日插值算法，编制并调试相应的函数子程序
（2）设计牛顿插值算法，编制并调试相应的函数子程序
（3）给定函数四个点的数据如下：
1.1
2.3
3.9
5.1
Y
3.887
4.276
4.651
2.117
试用拉格朗日插值确定函数在x=2.101，4.234处的函数值。
4）已知
用牛顿插值公式求
的近似值。
3．实验原理
写出本次实验所用算法的算法步骤叙述或画出算法程序框图
4．实验环境及实验文件存档名
5.实验结果及分析
（1）拉格朗日法
程序：
functionlagrange（x）
x0=[1.12.33.95.1];
y0=[3.8874.2764.6512.117];
n=length（x0）;
s=0;
forj=0:
（n-1）
t=1;
fori=0:
t=t*（x-x0（i+1））/（x0（j+1）-x0（i+1））;
s=s+t*y0（j+1）;
s
运行结果：
（2）牛顿差值法：
functionNewton（x）
%显示15位
x0=[149];
%x的值
y0=[123];
%y的值
x=5;
%插值点
n=max（size（x0））;
y=y0
（1）;
%迭代初始值
disp（y）;
s=1;
dx=y0;
fori=1:
n-1%构造差商表
dx0=dx;
n-i
dx（j）=（dx0（j+1）-dx0（j））/（x0（i+j）-x0（j））;
df=dx
（1）;
s=s*（x-x0（i））;
y=y+s*df;
%计算
disp（y）;
实际结果是2.9979.有一定的误差。
总体来说还是很准确的
数值微积分
（1）学会复化梯形、复化辛浦生求积公式的应用
（4）会用这些函数解决实际问题
（1）设计复化梯形公式求积算法，编制并调试相应的函数子程序
（2）设计复化辛浦生求积算法，编制并调试相应的函数子程序
（4）分别用复化梯形公式和复化辛浦生公式计算定积分
取n=2，4，8，16，精确解为0.9460831
3、实验原理
本实验是用MATLAB来完成。
梯形法存档名为trap.m辛浦生法存档名为simpson.m
复化梯形公式求积算法：
functionT=trap（n）
f=sym（'
sin（x）/x'
a=0;
b=1;
h=（b-a）/n;
T=0;
x0=a+h*k;
T=T+limit（f,x0）;
T=h*（limit（f,a）+limit（f,b））/2+h*T;
T=double（T）;
复化辛浦生法：
functionS=simpson（n）
h=（b-a）/（2*n）;
s1=0;
s2=0;
x0=a+h*（2*k-1）;
s1=s1+limit（f,x0）;
x0=a+h*2*k;
s2=s2+limit（f,x0）;
S=h*（limit（f,a）+limit（f,b）+4*s1+2*s2）/3;
S=double（S）;
运行结果
常微分程的数值解法
（1）学会四阶龙格-库塔法的使用
（2）设计出相应的算法，编制相应的函数子程序
（3）会用这些函数解决实际问题
（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）
functionRK（h,n）
F='
-50*y'
;
b=a+h*n;
X1=a:
h:
b;
Y1=zeros（1,n+1）;
Y1
（1）=10;
x=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;
X1
disp（Y1）;
h=0.05，N=10时，运行结果：
h=0.025,n=20时，运行结果：
h=0.01，N=50时，输出结果：
（2）微分程为：
v’=9.8-0.003857v^2,v（0）=0
functionvelosity（h,n）
9.8-0.003857*y*y'
t=a:
v=zeros（1,n+1）;
v
（1）=0;
x=t（i）;
y=v（i）;
y=v（i）+K2/2;
x=t（i）+h;
y=v（i）+K3;
v（i+1）=v（i）+（K1+2*K2+2*K3+K4）/6;
t
disp（v）;
升级会员 