中南大学数值分析计算方法实验报告完整版Word格式.docx
《中南大学数值分析计算方法实验报告完整版Word格式.docx》由会员分享，可在线阅读，更多相关《中南大学数值分析计算方法实验报告完整版Word格式.docx（30页珍藏版）》请在冰点文库上搜索。
fori=1:
k
x1（i+1）=1+1/（x1（i））^2
x2（i+1）=（1+x2（i）^2）^（1/3）
x3（i+1）=1/sqrt（x3（i）-1）
4、实验结果：
>
>
diedai1
y1是收敛的
y2是收敛的
y3是发散的
输入迭代次数1
k=
1
x1=
1.50001.4444
x2=
1.50001.4812
x3=
1.50001.4142
>
5、实验总结：
实验二
非线性方程与方程组的数值解法
（2）
1、设计目的：
（1）symsx
y=exp（x）+10*x-2;
f=solve（'
exp（x）+10*x-2=0'
a=0;
b=1;
设置检验次数'
l（i）=subs（y,0）%l
（1）为负值
k（i）=subs（y,1）%k（i）为正值
m=subs（y,（a+b）/2）;
ifm<
a=（a+b）/2
else
b=（a+b）/2
end
（2）clearall
formatlong
请输入迭代次数'
x
（1）=0
x（i+1）=（2-exp（x（i）））/10
4、实验结果：
（1）>
erfenfa
设置检验次数3
k=
3
l=
-1
10.718281828459045
b=
0.500000000000000
-1-1
10.71828182845904510.718281828459045
0.250000000000000
-1-1-1
10.71828182845904510.71828182845904510.718281828459045
0.125000000000000
（2）>
diedai2
请输入迭代次数2
x=
0
00.100000000000000
00.1000000000000000.089482908192435
实验三
牛顿下山法
用牛顿下山法
0.6附近的根。
此程序有一个主函数，储存在一个M文件中，实验时调用函数完成。
主函数一（niudunxiashan1）
clearall
x
（1）=0.6;
inf
x（i+1）=x（i）-（x（i）^3-x（i）^2-1）/（3*x（i）^2-2*x（i））;
y（i）=x（i）^3-x（i）^2-1;
y（i+1）=x（i+1）^3-x（i+1）^2-1;
ifabs（y（i+1））>
abs（y（i））
fork=1:
50
x（i+1）=x（i）-（1/2^（k-1））*（x（i）^3-x（i）^2-1）/（3*x（i）^2-2*x（i））;
ifabs（y（i+1））<
break
iffloor（y（i+1）*10000）==0
j=i+1;
disp（['
迭代次数为:
'
]）
disp（j）
x
输入
niudunxiashan1
输出
8
0.60000.0042-0.00321.21031.56091.47341.46561.4656
实验四
史蒂芬孙迭代法
求
在0.6附近的根，迭代公式为
（1）
（2）
（3）
此程序有两个主函数，储存在两个M文件中，实验时分别调用函数完成。
主函数一（steffensen1）
y（i）=1+1/x（i）^2;
z（i）=1+1/y（i）^2;
x（i+1）=x（i）-（y（i）-x（i））^2/（z（i）-2*y（i）+x（i））;
iffloor（x（i+1）*1000）==floor（x（i）*1000）
y（j）=0;
z（j）=0;
A=['
k'
'
x'
y'
z'
];
M=zeros（j,4）;
j
M（i,1）=i-1;
M（i,2）=x（i）;
M（i,3）=y（i）;
M（i,4）=z（i）;
xlswrite（'
E:
\ProgramFiles\MATLAB\bin\迭代加速\steffensen.xls'
A,'
sheet1'
A1:
D1'
M,'
A2:
D7'
主函数二（diedai）
y（i）=（1+x（i）^2）^（1/3）;
z（i）=（1+y（i）^2）^（1/3）;
sheet2'
D5'
输入
steffensen1
steffensen2
实验五
迭代法的计算
用列主元高斯法解下列线性方程组：
3.01006.03001.9900x11
1.27004.1600-1.2300x2=1
0.9870-4.81009.3400x31
3.00006.03001.9900x11
0.9900-4.81009.3400x31
用雅可比迭代法和高斯-塞得尔迭代法解下列线性方程组并检查收敛性：
5x1+2x2+x3=-12
-x1+4x2+2x3=20
+2x1-3x2+10x3=3
3、详细设计：
（1）functionx=gauss_lie（A,b）
n=length（b）;
p=1:
n;
lu=A;
y=[];
fork=1:
n
[c,i]=max（abs（lu（k:
n,k）））;
ik=i+k-1;
ifik~=k
m=p（k）;
p（k）=p（ik）;
p（ik）=m;
ck=lu（k,:
lu（k,:
）=lu（ik,:
lu（ik,:
）=ck;
ifk==n
break;
lu（k+1:
n,k）=lu（k+1:
n,k）/lu（k,k）;
n,k+1:
n）=lu（k+1:
n）-lu（k+1:
n,k）*lu（k,k+1:
n）;
l=diag（ones（n,1））+tril（lu,-1）;
u=triu（lu）;
y
（1）=b（p
（1））;
fori=2:
y（i）=b（p（i））-l（i,1:
i-1）*y（1:
i-1）'
;
x（n）=y（n）/u（n,n）;
fori=n-1:
-1:
x（i）=（y（i）-u（i,i+1:
n）*x（i+1:
n）'
）/u（i,i）;
A
b
detA=det（A）
x=x'
（2）
functionx=gauss_seidel（A,b,x0,tol）
if（nargin==2）
x0=ones（size（b））;
tol=1e-6;
elseif（nargin==3）
sprintf（'
USAGE:
Gauss_Seidel（A,b,x0,tol）'
D=diag（diag（A））;
U=triu（A,1）;
L=tril（A,-1）;
G=-（D+L）\U;
d1=（D+L）\b;
x=G*x0+d1;
n=1;
whilenorm（x-x0）>
=tol
x0=x;
x=G*x0+d1;
n=n+1;
functionx=Jacobi（A,b,x0）
n=10000;
%最大迭代次数
ee=1e-4;
%精度
n1=length（b）;
x1=zeros（n1,1）;
p=0;
p=p+1;
forj=1:
n1
s=0;
ifk~=j
s=s+A（j,k）*x0（k）;
x1（j）=（b（j）-s）/A（j,j）;
end
ifnorm（x1-x0）<
ee
x0=x1;
n=p
x=x1;
（1）将代码保存为gauss_lie.m，在Matlab主界面输入下述代码，得到结果：
A=[3.01,6.03,1.99;
1.27,4.16,-1.23;
0.987,-4.81,9.34];
b=[1,1,1]'
gauss_lie（A,b）
A=[3.00,6.03,1.99;
0.990,-4.81,9.34];
（2）将代码保存为Jacobi.m和Gauss_Seidel.m，在Matlab主界面输入下述代码，得到结果：
A=[5,-1,2;
2,4,-3;
1,2,10]'
b=[-12,20,3]'
x0=[0,0,0]'
Jacobi（A,b,x0）
gauss_seidel（A,b）
实验六
列主元三角分解
编写列主元矩阵分解法的程序，并以实例验证。
A=
.
3.详细设计：
主函数一（liezhuyuanxiaoqu（A）
functionoutput=liezhuyuanxiaoqu（A）
n=size（A）;
forr=2:
fori=r:
n-1
b=0;
t=0;
r-1
b=A（r,k）*A（k,i）+b;
A（r,i）=A（r,i）-b;
s（i）=A（r,i）;
[a,m]=max（abs（s））;
B=A（r,:
A（r,:
）=A（r+m-1,:
A（r+m-1,:
）=B;
t=A（i+1,k）*A（k,r）+t;
A（i+1,r）=（A（i+1,r）-t）/A（r,r）;
b=A（r,k）*A（k,n）+b;
A（r,n）=A（r,n）-b;
L=eye（n）;
L（i+1:
n,i）=A（i+1:
n,i）;
U=zeros（n）;
U（1:
i,i）=A（1:
i,i）;
L
U
A=[1/41/51/6
1/31/41/5
1/212]
liezhuyuanxiaoqu（A）
A=
0.25000.20000.1667
0.50001.00001.9167
0.33330.1167-0.0792
L=
1.000000
0.50001.00000
0.33330.11671.0000
U=
01.00001.9167
00-0.0792
5、实验总结
实验七
牛顿插值
编写牛顿插值的计算程序，并以实例验证。
t=1.7
主函数一（niudunchazhi（x,y,t））
functions=niudunchazhi（x,y,t）
symsp;
n=length（x）;
s=y
（1）;
dd=0;
dxs=1;
forj=i+1:
dd（j）=（y（j）-y（i））/（x（j）-x（i））;
temp1（i）=dd（i+1）;
dxs=dxs*（p-x（i））;
s=s+temp1（i）*dxs;
y=dd;
simplify（s）;
ifnargin==2
s=subs（s,'
p'
s=collect（s）;
s=vpa（s,4）;
m=length（t）;
fori=1:
m
temp（i）=subs（s,'
t（i））;
s=temp;
x=1:
7；
y=[00.67890.79631.03751.29361.37451.7954]；
niudunchazhi（x,y）
niudunchazhi（x,y,1.7）
ans=
0.0001021*x^6+0.003798*x^5-0.1044*x^4+0.8751*x^3-3.295*x^2+5.88*x-3.359
0.5975
实验八
书上习题
分别用梯形公式和辛普森公式计算下列积分
（1）
，n=8；
（2）
，n=6.
此程序有两个主函数，储存在两个M文件中，实验时调分别用函数完成。
主函数一（twoone）
n=8;
symsx
y=x/（4+x^2）;
fa=subs（y,'
a）;
fb=subs（y,'
b）;
T=fa+fb;
S=fa+fb;
h=（b-a）/n;
c=a+h/2;
fc=subs（y,'
c）;
S=S+4*fc;
（n-1）
xk=a+k*h;
fk=subs（y,'
xk）;
T=T+2*fk;
T=T*h/2;
xkb=xk+h/2;
fkb=subs（y,'
xkb）;
S=S+4*fkb+2*fk;
S=S*h/6;
disp（'
用梯形公式结果为:
disp（T）;
用辛普森公式结果为:
disp（S）;
主函数一（twothree）
clc
n=6;
y=sqrt（4-（sin（x））^2）;
b=pi/6;
twoone
0.111402354529548
0.111571813252631
twothree
.0356********
.0357********
实验九
梯形公式辛普森公式
用龙贝格求积方法计算下列积分，使误差不超过
；
（2）
主函数一（eightone）
clc
T=zeros（10,10）;
y=2*exp（-x）/sqrt（pi）;
T（1,1）=（fa+fb）/2;
9
g=0;
forn=1:
（2^（k-1））
xk=a+（2*n-1）*（b-a）/（2^k）;
g=g+fk;
T（k+1,1）=T（k,1）/2+（b-a）*g/（2^k）;
10
forj=2:
i
m=j-1;
l=4^m;
T（i,j）=（l*T（i,j-1）-T（i-1,j-1））/（l-1）;
ifceil（T（i,j）*100000）==ceil（T（i,j-1）*100000）
c=i;
d=j;
T
T（i,j）
eightone
0.713271669814180
eightthree
10.162443613495572
实验十
龙格贝塔迭代原理
用辛普森公式的自适应积分计算
，允许误差
主函数一（nine）
y=（x^2）*log（x）;
a=1;
b=1.5;
c=（a+b）/2;
S=（b-a）*（fa+4*fc+fb）/6;
d=（a+c）/2;
e=（c+b）/2;
fd=subs（y,'
d）;
fe=subs（y,'
e）;
s2=（c-a）*（fa+4*fd+fc）/6;
s3=（b-c）*（fc+4*fe+fb）/6
升级会员 