数学实验报告 2.docx
《数学实验报告 2.docx》由会员分享，可在线阅读，更多相关《数学实验报告 2.docx（54页珍藏版）》请在冰豆网上搜索。
数学实验报告2
数学实验报告
院系计算机学院
班级2004级二学位
指导教师张兴永
姓名田琳
实验目录
1．计算
程序内容：
源程序：
\s.m；
程序清单：
functions（x）
s=1;
for（i=1:
8）;
s=s*cos（2^i*x）;
end
s
实验结果：
S
（1）=3.4160e-004;s
（2）=8.1827e-004;s（pi）=1;s（pi/2）=-1;s（pi/3）=0.0039
2.求
的fourier、Laplace、ztrans变换
程序内容：
fourier、Laplace、ztrans变换
源程序：
\bianhuan.m；
程序清单：
symst;
f=5*sin（2*t）-3*cos（2*t）
fourier_f=fourier（f）
laplace_f=laplace（f）
ztrans_f=ztrans（f）
实验结果：
f=5*sin（2*t）-3*cos（2*t）
fourier_f=pi*（5*i*Dirac（w+2）-3*Dirac（w+2）-5*i*Dirac（w-2）-3*Dirac（w-2））
laplace_f=10/（s^2+4）-3*s/（s^2+4）
ztrans_f=
10*z*cos
（1）*sin
（1）/（-4*z*cos
（1）^2+z^2+2*z+1）-3*（z+1-2*cos
（1）^2）*z/（-4*z*cos
（1）^2+z^2+2*z+1）
3．矩阵高斯消去法
对矩阵
进行高斯消去法变换
求三次初等变换的总的等价乘子,并用MATLAB求原来A的行列式、秩和迹
源程序：
\gaosi.m
程序清单：
A=[107;415;2-19];A0=A;%输入A，并保留一个备份
A（2,:
）=-4*A（1,:
）+A（2,:
）
A1=A,
B1=A1/A0%消去A（2,1）,求B1
A（3,:
）=-2*A（1,:
）+A（3,:
）
A2=A,
B2=A2/A1%消去A（3,1）
A（3,:
）=-A（3,2）/A（2,2）*A（2,:
）+A（3,:
）
A3=A,
B3=A3/A2%消去A（3,2）
B0=A3/A0%求三次初等变换的总的等价乘子
det_A=det（A0）,
rank_a=rank（A0）,
tr_A=trace（A0）,%求原来A的行列式、秩和迹
实验结果：
，
det_A=-28rank_a=3tr_A=11
4．画出心型线、星型线、双纽线和四叶玫瑰线的图形
心形线绘制：
Clear
t=0:
0.001:
2*pi;
subplot（2,2,1）;
polar（a,1+cos（t））
subplot（2,2,2）;
plot（cos（t）.^3,sin（t）.^3）
subplot（2,2,3）;
polar（t,abs（sin（t）.*cos（t）））
subplot（2,2,4）;
polar（t,（cos（2*t））.^0.5）
图像：
图表1四曲线
5．求fibonacci数
图表2fib曲线
实验结果：
100以内的fibonacci数：
f=
1123581321345589
实验结果：
1000以内的fibonacci数：
f=
1123581321345589
144233377610987
6．数据三次拟合曲线
数据：
1234567899123165198243277353345303288275的三次拟合
程序：
\nihe3.m
图表3三次拟合曲线
7．递推公式的稳定性
实验内容：
[1]P11页试验课题1
实验程序:
clear;
clc;
symsx;
resualt=zeros（4）;
N_R=1;
a=input（'输入a的值:
'）;
%------------------------------------方案1
I1=log（a+1）-log（a）;
forn=1:
10
I1=-1*a*I1+1/n;
f=x^n/（a+x）;
I0=int（f,'x',0,1）;
I0=vpa（I0,500）;
I0=vpa（I0,8）;
I=vpa（I1,6）;
wucha1=abs（（I0-I1）/I0）;
wucha1=vpa（wucha1,8）;
db=[nI0Iwucha1];
ifn==1
resualt=db;
else
resualt=[resualt;db];
end
end
resualt
%------------------------------------方案2
resualt=zeros（4）;
N=13;
ifa>=N/（N+1）
I2=（2*a+1）/（2*a*（a+1）*（N+1））;
else
I2=0.5*（1/（（a+1）*（N+1））+1/N）;
end
forn=N:
-1:
1
f=x^n/（a+x）;
I0=int（f,'x',0,1）;
I0=vpa（I0,500）;
I0=vpa（I0,8）;
I2=（-1*I2+1/n）/a;
I=vpa（I2,6）;
wucha2=abs（（I0-I2）/I0）;
wucha2=vpa（wucha2,8）;
db=[nI0Iwucha2];
ifn==N
resualt=db;
else
resualt=[resualt;db];
end
end
resualt
实验结果:
输入a的值:
0.05
次数精确值迭代值相对误差
[1,.84777388,.847774,.22248517e-8]
[2,.45761131,.457611,.85349539e-8]
[3,.31045277,.310453,.63500227e-8]
[4,.23447736,.234477,.68175840e-8]
[5,.18827613,.188276,.10198168e-7]
[6,.15725286,.157253,.44935959e-9]
[7,.13499450,.134994,.10844168e-8]
[8,.11825028,.118250,.42221299e-7]
[9,.10519860,.105199,.25088308e-7]
[10,.94740070e-1,.947401e-1,.13928927e-8]
次数精确值迭代值相对误差
[13,.72971840e-1,.889587e-1,.21908205]
[12,.79024729e-1,-.112507,2.4236878]
[11,.86172087e-1,4.06831,46.211490]
[10,.94740070e-1,-79.3663,838.72635]
[9,.10519860,1589.55,15108.966]
[8,.11825028,-31788.4,268824.43]
[7,.13499450,635772.,4709611.4]
[6,.15725286,-.127154e8,80859783.]
[5,.18827613,.254309e9,.13507216e10]
[4,.23447736,-.508617e10,.21691531e11]
[3,.31045277,.101723e12,.32766162e12]
[2,.45761131,-.203447e13,.44458454e13]
[1,.84777388,.406894e14,.47995561e14]
输入a的值:
15
次数精确值迭代值相对误差
[1,.31922183e-1,.319222e-1,.19912866e-8]
[2,.21167256e-1,.211673e-1,.21971093e-8]
[3,.15824494e-1,.158245e-1,.19548773e-8]
[4,.12632590e-1,.126326e-1,.36732305e-7]
[5,.10511157e-1,.105112e-1,.37710110e-8]
[6,.89993133e-2,.899931e-2,.11542747e-6]
[7,.78674434e-2,.786746e-2,.19750588e-5]
[8,.69883483e-2,.698812e-2,.33252485e-4]
[9,.62858867e-2,.628937e-2,.55451370e-3]
[10,.57116994e-2,.565942e-2,.91538526e-2]
次数精确值迭代值相对误差
[13,.44830339e-2,.482067e-2,.75313181e-1]
[12,.48293362e-2,.523418e-2,.83829670e-1]
[11,.52335998e-2,.571166e-2,.91344598e-1]
[10,.57116994e-2,.628589e-2,.10052873]
[9,.62858867e-2,.698835e-2,.11175216]
[8,.69883483e-2,.786744e-2,.12579441]
[7,.78674434e-2,.899931e-2,.14386756]
[6,.89993133e-2,.105112e-1,.16799544]
[5,.10511157e-1,.126326e-1,.20182674]
[4,.12632590e-1,.158245e-1,.25267218]
[3,.15824494e-1,.211673e-1,.33762608]
[2,.21167256e-1,.319222e-1,.50809264]
[1,.31922183e-1,.645385e-1,1.0217452]
明显,对方案2计算结果都不可靠.
8．迭代法的收敛性与收敛速度的比较
实验内容：
[1]P37页试验课题二
实验程序
clear;
clc;
symsx;
f=x^3-sin（x）-12*x+1;
resualt=solve（'x^3-sin（x）-12*x+1','x'）;
disp（'Matlab求根结果:
'）
r=vpa（resualt,8）
pause;
df=diff（f,'x'）;
x0=input（'输入迭代初值:
'）;
N=input（'输入最多迭代次数:
'）;
e=input（'输入迭代精度:
'）;
fork=1:
N
f0=subs（f,x0）;
df0=subs（df,x0）;
ifdf0==0
disp（'导数为0,停止计算'）
break;
else
xx=x0-f0/df0;
ifabs（xx）<=1
E=abs（xx-x0）;
else
E=abs（（x0-xx）/xx）;
end
ifEdisp（'牛顿法__计算结果:
'）
x=vpa（xx,8）
disp（'计算误差为:
'）
error=vpa（abs（（xx-r）/r）,5）
f_x=vpa（subs（f,xx）,8）
disp（'最终迭代次数:
'）
k
break;
else
x0=xx;
end
end
end
ifk==N
disp（'经过设置的迭代次数没有收敛,计算失败'）
end
pause;
disp（'普通迭代法:
'）
disp（'迭代式1:
'）
x0=rand（-4,-3）;
f1=（12*x+sin（x）-1）^（1/3）;
fork=1:
N
xx=subs（f1,x0）;
ifabs（xx）<=1
E=abs（xx-x0）;
else
E=abs（（x0-xx）/xx）;
end
ifEdisp（'[-4,-3]计算结果:
'）
x=vpa（xx,8）
disp（'计算误差为:
'）
error=vpa（abs（（xx-r）/r）,5）
f_x=vpa（subs（f,xx）,8）
disp（'最终迭代次数:
'）
k
break;
else
x0=xx;
end
end
ifk==N
disp（'经过设置的迭代次数没有收敛,计算失败'）
end
pause;
x0=rand（3,4）;
f1=（12*x+sin（x）-1）^（1/3）;
fork=1:
N
xx=subs（f1,x0）;
ifabs（xx）<=1
E=abs（xx-x0）;
else
E=abs（（x0-xx）/xx）;
end
ifEdisp（'[3,4]牛顿法__计算结果:
'）
x=vpa（xx,8）
disp（'计算误差为:
'）
error=vpa（abs（（xx-r）/r）,5）
f_x=vpa（subs（f,xx）,8）
disp（'最终迭代次数:
'）
k
break;
else
x0=xx;
end
end
ifk==N
disp（'经过设置的迭代次数没有收敛,计算失败'）
end
pause;
disp（'普通迭代法:
'）
disp（'迭代式2:
'）
x0=rand（0,0.2）;
f1=（12*x+sin（x）-1）^（1/3）;
fork=1:
N
xx=subs（f1,x0）;
ifabs（xx）<=1
E=abs（xx-x0）;
else
E=abs（（x0-xx）/xx）;
end
ifEdisp（'[0,0.2]计算结果:
'）
x=vpa（xx,8）
disp（'计算误差为:
'）
error=vpa（abs（（xx-r）/r）,5）
f_x=vpa（subs（f,xx）,8）
disp（'最终迭代次数:
'）
k
break;
else
x0=xx;
end
end
ifk==N
disp（'经过设置的迭代次数没有收敛,计算失败'）
end
实验结果
函数
的图像为：
图表4
方案
初值
迭代次数
迭代结果
误差
牛顿法
-3
5
-3.4911788
0
普通迭代
（1）
-3
发散
普通迭代
（2）
-3
发散
牛顿法
3
5
3.4101251
0
普通迭代
（1）
3
发散
普通迭代
（2）
3
发散
牛顿法
-1
4
.76963989e-1
.74718010e-13
牛顿法
0
4
.76963989e-1
0
普通迭代
（1）
0
2
.76964248e-1
-.33644973e-5
[34]
0
2
.77153212e-1
-.24559693e-2
普通迭代
（2）
0
2
.14285360
-.85369616
图表5
9．雅可比迭代法与高斯塞德尔法的收敛性与收敛速度
实验内容：
[1]P71页试验课题（四），雅可比迭代法与高斯——塞德尔迭代法的收敛性与收敛速度。
实验程序：
clear;
clc;
A=input（'输入系数矩阵A:
'）;
b1=input（'输入矩阵b1:
'）;
b2=input（'输入矩阵b2:
'）;
disp（'Matlab计算结果:
'）
x1=inv（A）*b1;
X1=vpa（x1,6）
x2=inv（A）*b2;
X2=vpa（x2,6）
L=tril（A）;
U0=triu（A）;
U=L-A;
L=U0-A;
D=A+L+U;
pause;
N=input（'输入最多迭代次数:
'）;
disp（'雅可比迭代计算结果:
'）
MJ=inv（D）*（L+U）;
disp（'计算b1情况'）;
X0=zeros（sqrt（numel（A））,1）;
%与输入矩阵配置相同的初始值
fork=1:
N
X=MJ*X0+D^-1*b1;
if（norm（X-X0）<1e-6）
b1
disp（'计算结果:
'）
X1=vpa（X,6）
disp（'迭代次数'）
k
AX=vpa（A*X,6）
break;
else
X0=X;
end
end
ifk==N
disp（'迭代没有收敛'）
end
pause;
disp（'计算b2情况'）;
X0=zeros（sqrt（numel（A））,1）;
%与输入矩阵配置相同的初始值
fork=1:
N
X=MJ*X0+D^-1*b2;
if（norm（X-X0）<1e-6）
b2
disp（'计算结果:
'）
X2=vpa（X,6）
disp（'迭代次数'）
k
AX=vpa（A*X,6）
break;
else
X0=X;
end
end
ifk==N
disp（'迭代没有收敛'）
end
pause
disp（'高斯--塞德尔迭代计算结果:
'）
MG=（D-L）^-1*U;
disp（'计算b1情况'）;
X0=zeros（sqrt（numel（A））,1）;
%与输入矩阵配置相同的初始值
fork=1:
N
X=MG*X0+（D-L）^-1*b1;
if（norm（X-X0）<1e-6）
b1
disp（'计算结果:
'）
X1=vpa（X,6）
disp（'迭代次数'）
k
AX=vpa（A*X,6）
break;
else
X0=X;
end
end
ifk==N
disp（'迭代没有收敛'）
end
pause;
disp（'计算b2情况'）;
X0=zeros（sqrt（numel（A））,1）;
%与输入矩阵配置相同的初始值
fork=1:
N
X=MG*X0+（D-L）^-1*b2;
if（norm（X-X0）<1e-6）
b2
disp（'计算结果:
'）
X2=vpa（X,6）
disp（'迭代次数'）
k
AX=vpa（A*X,6）
break;
else
X0=X;
end
end
ifk==N
disp（'迭代没有收敛'）
end
实验结果
1．
，
，
结果：
，
雅可比迭代法计算b1迭代24次，计算b2迭代30次。
高斯塞德尔法计算b1迭代15次，计算b2迭代20次。
2．
，
，
结果：
雅可比迭代法计算计算不收敛。
高斯塞德尔法计算b1迭代45次，计算b2迭代52次。
3．
，
结果为
但使用两种方法都不收敛。
10．龙格现象的发生、防止和插值效果的比较
实验内容：
[1]P102页实验课题
（一）龙格现象的发生、防止和插值效果的比较
实验程序：
clear;
clc;
f='x/（1+x^4）';
N=input（'输入插值次数'）;
xi=zeros（1,N+1）;
h=10/N;
fork=1:
N+1
xi（k）=-5+（k-1）*h;
end
yi=subs（f,xi）;
xk=zeros（1,41）;
fork=1:
41
xk（k）=-5+0.25*（k-1）;
end;
%拉格朗日插值算法
ykl=zeros（1,41）;
fori=1:
41
ykl（i）=0;
fork=1:
N+1
t=1;
forj=1:
N+1
ifj~=k
t=（xk（i）-xi（j））/（xi（k）-xi（j））*t;
ykl（i）=ykl（i）+t*yi（k）;
end
end
end
end
subplot（2,2,1）;
plot（xk,ykl,'b'）
title（'拉格朗日插值曲线'）
%分段线性插值
ykx=zeros（1,41）;
fori=1:
41
k=1;
k=k+1;
whilek=xi（k）
k=k+1;
end
ykx（i）=yi（k-1）*（xk（i）-xi（k））/（xi（k-1）-xi（k））+yi（k）*（xk（i）-xi（k-1））/（xi（k）-xi（k-1））;
end
subplot（2,2,2）;
plot（xk,ykx,'g'）
title（'分段线性插值曲线'）
%
yky=zeros（1,41）;
yky=interp1（xi,yi,xk,'spline'）;
subplot（2,2,3）;
plot（xk,yky,'r'）
title（'三次样条插值曲线'）
disp（'计算结果输出:
'）
z=[xk'ykl'ykx'yky']
结算结果
1．
的十次插值图像为：
二十次插值图像：
插值顺序：
，拉格朗日插值、分段线形插值、I型三次样条插值：
十次插值：
-5.0000-0.0799-0.0080-0.0080
-4.75002.7400-0.0099-0.0007
-4.50002.3631-0.0118-0.0011
-4.25001.0023-0.0137-0.0068
-4.0000-0.2309-0.0156-0.0156
-3.7500-0.9262-0.0208-0.0252
-3.5000-1.0791-0.0261-0.0334
-3.2500-0.8821-0.0313-0.0380
-3.0000-0.5846-0.0366-0.0366
-2.7500-0.4067-0.0569-0.0300
-2.5000-0.4943-0.0771-0.0309
-2.2500-0.9042-0.0974-0.0549
-2.0000-1.6096-0.1176-0.1176
-1.7500-2.5185-0.2132-0.2259
-1.5000-3.4987-
升级会员 