西南交大MATLAB上机题解Word文件下载.docx
《西南交大MATLAB上机题解Word文件下载.docx》由会员分享，可在线阅读，更多相关《西南交大MATLAB上机题解Word文件下载.docx（19页珍藏版）》请在冰豆网上搜索。
1,4,-2;
-3,1,4];
b=[-3;
2;
4];
x=Gauss（a,b）
运行结果：
x=[;
;
]
题目2在CommandWindow中的程序：
a=[1,,;
1,;
,1];
b=[3;
1];
x=Gauss（a,b）
题目3在CommandWindow中的程序：
a=[1,3;
-7,1];
b=[4;
6];
2.Gauss-Seidel迭代法
Gauss-Seidel迭代算法类似Jacobi迭代算法，Jacobi迭代算法的一般迭代格式为：
式
（1）
在Jacobi迭代算法收敛时，式
（1）中第一个方程计算得到的
马上投入到第二方程中使用，把前两个方程已经得到的
，
马上投入到第三个方程中使用，如此循环，通过这种计算线性方程组解的迭代格式称为Gauss-Seidel迭代，其迭代格式为：
式
（2）
矩阵形式如式（3）所示：
式（3）
解得
，得到的结果如式（4）所示：
式（4）
一般线性方程组通用的Gauss-Seidel迭代法的迭代程序见附录二。
由于Gauss-Seidel迭代法计算时要首先考虑迭代收敛与否的问题，所以在M文件中，通过语言控制，首先对求解的线性方程组的矩阵进行相关计算，得到迭代阵B的谱半径值，如果谱半径值R<
1，则进行Gauss-Seidel迭代，否则不进行迭代。
A=[6,2,-1;
b=[-3;
x0=[0;
0;
0];
x=GaussSeidel（A,b,x0）
经过计算，得到题目1的迭代阵B的谱半径R=<
1，所以能够得到收敛解。
经过41次迭代后，得到该线性方程组的解为：
A=[1,,;
经过计算，得到题目2的迭代阵B的谱半径R=<
经过12次迭代后，得到该线性方程组的解为：
A=[1,3;
经过计算，得到题目3的迭代阵B的谱半径R=21>
1，所以该线性方程组不能得到用Gauss-Seidel迭代方法求得收敛解。
通过使用上面三个例题求解的过程发现，在MATLAB中编写的M文件能够很好地求解线性方程组，Gauss消元法能够得到上面三个例题的全部收敛解，而Gauss-Seidel虽然无法求解到题目3的收敛解，但是整个M文件程序的使用性能还是非常的不错，能够真实反映出线性方程组系数矩阵内在的一些性质。
第二题
问题：
给定初值问题
（精确解为
），用Runge-Kutta4阶算法按步长
求解，分析其中遇到的现象及问题。
龙格-库塔法（Runge-Kutta）经常用来解决常微分方程数值问题，与欧拉（Euler）法、线性多步法、Adams内插公式等其他几种常微分方程数值求解方法在应用上各有千秋。
龙格-库塔法的推导是基于泰勒（Taylor）展开方法，故而此法要求求解的微分方程解具有较好的光滑性质。
根据推导得到四阶龙格-库塔法的具体格式如式（5）所示：
其中：
式（5）
求解本题微分方程的程序保存在M文件中，具体的内容见附录三。
当步长h=时，运行程序能够得到相应的结果，将该结果绘制成图形，如图
（1）所示：
图
（1）
当h=时，结果对应的图形如图
（2）所示：
图
（2）
当h=时，结果对应的图形如图（3）所示：
图（3）
第三题
给定函数
，及节点
，利用高次多项式、分段线性多项式进行插值，画图比较几种插值多项式与
的图形，并说明其中的数学现象。
利用Newton插值法对该函数进行插值，通用的Newton插值法程序保存在M文件中，具体程序清单见附录四。
实际插值时，先将对应的M文件打开，在MATLAB的CommandWindow中输入下面的程序：
x=[-5:
1:
5];
y=1./（1+x.^2）;
x0=[-5:
:
y0=Lagrange（x,y,x0）;
y1=1./（1+x0.^2）;
gridon;
plot（x0,y0,'
--r'
）;
holdon;
plot（x0,y1,'
--b'
xlabel（'
自变量--x'
ylabel（'
因变量--y'
legend（'
Lagrange插值图形'
'
函数准确图形'
title（'
Lagrange插值图形和函数准确图形'
gridon
最终得到如图（4）的图形：
图（4）
利用分段线性插值多项式进行插值计算，程序保存在M文件中，具体程序清单见附录四。
u=[-5:
v=pline（x,y,u）;
plot（x,y,'
plot（u,v,'
pline插值图形'
pline插值图形和函数准确图形'
最终得到如图（5）的图形：
图（5）
理论上而言，对于一个给定区间[a,b]，根据节点分布采用插值多项式得到函数的近似值时，插值多项式的次数越高，逼近
值的精度就越高。
Newton法插值得到的多项式
次数明显高于分段线性插值，所以Newton插值法逼近
值的精度要高于分段线性插值法，但是通过上面两个实验获得的图形可以发现，事实并非如此，图（4）和图（5）展示出了Newton插值法在区间的首位处逼近
值的精度明显低于分段线性插值，从而得出高次多项式插值的精度并不是在区间上任何地方都高于分段线性的结论。
造成本题出现该问题产生的原因是：
函数
在区间[-5,5]上的各阶导数存在，但是由n个节点所构成的Lagrange插值多项式在全区间内并非都收敛。
第四题
利用牛顿法求方程
于区间
的根，考虑不同初值下牛顿法的收敛情况。
迭代法是求解非线性方程（组）的主要方法，而Newton法是最有效的迭代方法之一，它在单根附近具有较高阶的收敛速度，统统是还能够求解到方程的复根（此时初值需要时复数）。
通常是在方程
的根
附近，用
的某个近似值去代替
，再令这个近似式为零，解出的根就作为根
的新的近似值。
基于这样的一种思想，可以推导得到Newton法的计算公式，首先设
的近似根为
，函数
可以在
附近作Taylor展开，得到如式（6）：
式（6）
利用式（6）的线性部分作为函数
的近似，得到如式（7）的公式：
式（7）
求解式（7），得到根值为
，即（8）：
式（8）
由此推断得到Newton迭代公式为式（9）：
式（9）
针对此题的MATLAB程序见附录六，由于题中给定根值区间为
，故首先取初值
在MATLAB中CommandWindow窗口中输入以下应用程序：
f=inline（'
x-log（x）-2'
df=inline（'
1-1./x'
Newton（f,df,,
运行得到结果为：
第0次迭代（初值）x[0]=
第1次迭代x[1]=
第2次迭代x[2]=
第3次迭代x[3]=
第4次迭代x[4]=
ans=
改变初值
，运行得到结果为：
第0次迭代（初值）x[0]=
第5次迭代x[5]=
根据不同初值的迭代过程可以发现，不同的初值对应的迭代次数并不完全相同，距离最终根的距离越近，迭代的次数相对越少。
经过几次取值，计算结果都没有出现根不收敛的情况，所以无法分析初值对根最终收敛的影响。
附录
第一题程序清单
Gauss迭代的M文件源程序
function[x]=Gauss（a,b）
n=length（a）;
x=zeros（n,1）;
a=[ab];
fork=1:
n-1
max=k;
fori=k+1:
n
ifa（i,k）>
a（max,k）
max=i;
end
temp=a（k,k:
n+1）;
a（k,k:
n+1）=a（max,k:
a（max,k:
n+1）=temp;
n
a（i,k）=-a（i,k）/a（k,k）;
a（i,k+1:
n+1）=a（i,k+1:
n+1）+a（i,k）*a（k,k+1:
x（n,1）=a（n,n+1）/a（n,n）;
fori=n-1:
-1:
1
sum=0;
forj=i+1:
sum=sum+x（j,1）*a（i,j）;
x（i,1）=（a（i,n+1）-sum）/a（i,i）;
Gauss-Seidel迭代的M文件源程序：
functionx=GaussSeidel（A,b,x0）
D=diag（diag（A））;
L=-tril（A,-1）;
U=-triu（A,1）;
B=（D-L）\U;
f=（D-L）\b;
x=zeros（length（x0）,1）;
i=0
x0=x;
x=B*x+f;
R=max（abs（eig（B）））
ifR<
fori=1:
100
ifx0==x;
i=i+1
break;
else
i=i+1;
x0=x;
x=B*x0+f;
end
fprintf（'
不能够用Gauss-Seidel迭代法计算’）
end
第二题程序清单
Runge-Kutta4阶算法程序清单
clear;
F='
-20*x'
a=0;
b=2;
h=;
n=（b-a）/;
X=a:
b;
Y=zeros（1,n+1）;
Y
（1）=1;
fori=1:
x=X（i）;
y=Y（i）;
K1=eval（F）;
x=x+h/2;
y=y+h*K1/2;
K2=eval（F）;
y=Y（i）+h*K2/2;
K3=eval（F）;
x=X（i）+h;
y=Y（i）+h*K3;
K4=eval（F）;
Y（i+1）=Y（i）+h*（K1+2*K2+2*K3+K4）/6;
%准确值
temp=[];
f=dsolve（'
Dy=-20*x'
y（0）=1'
x'
df=zeros（1,n+1）;
n+1
temp=subs（f,'
X（i））;
df（i）=double（vpa（temp））;
disp（'
步长四阶经典R-K法准确值'
disp（[X'
Y'
df'
]）;
figure;
plot（X,df,'
k*'
X,Y,'
四阶经典R-K法解常微分方程'
准确值'
四阶经典R-K法'
第三题程序清单
高次多项式插值法（Newton法）程序清单
functionyy=Lagrange（x,y,xi）
m=length（x）;
n=length（y）;
ifm~=n,error（'
向量x和y的长度必须一致'
s=0;
z=ones（1,length（xi））;
forj=1:
ifj~=i
z=z.*（xi-x（j））/（x（i）-x（j））;
s=s+z*y（i）;
yy=s;
MATLAB中CommandWindow输入的应用程序清单
分段线性多项式插值法程序清单
functionv=pline（x,y,u）
delta=diff（y）./diff（x）;
n=length（x）;
k=ones（size（u））;
forj=2:
k（x（j）<
=u）=j;
end
s=u-x（k）;
v=y（k）+s.*delta（k）;
第四题程序清单
Newton法求解非线性方程组程序清单
functionx=Newton（f,df,x0,e,N）
ifnargin<
5
N=500;
4
e=1e-4;
x=x0;
x0=x+2*e;
k=0;
fprintf（'
第%2d次迭代（初值）x[%2d]=%\n'
k,k,x）
whileabs（x0-x）>
e&
k<
N,
k=k+1;
x=x0-feval（f,x0）/feval（df,x0）;
第%2d次迭代x[%2d]=%\n'
ifk==N,fprintf（'
迭代次数已经达到上限值'
end;
升级会员 