矩阵实习作业Word格式文档下载.docx
《矩阵实习作业Word格式文档下载.docx》由会员分享，可在线阅读，更多相关《矩阵实习作业Word格式文档下载.docx（17页珍藏版）》请在冰点文库上搜索。
49
1/3*a;
end
a
结果：
a=4.178866707295604e-024
2、前两项为
MATLAB程序:
functionwyc12
A=zeros（1,50）;
A
（1）=1;
A
（2）=1/3;
48
A（i+2）=（10/3）.*A（i+1）-A（i）;
A（50）
a=2.060436086367324e+006
经过程序发现第二种迭代过程结果偏离正常预期结果，通过观察发现当迭代到第18项左右时，由于出现了小数减去小数的状况进而造成了误差，进而可以得出这种迭代方法不可行。
二、利用迭代格式
及Aitken加速后的新迭代格式求方程
在
内的根
首先通过普通的迭代格式求求方程
内的根的MATLAB程序：
functionwyc2n
c=2*10^（-6）;
x=1;
x1=0;
whilec>
10^（-6）
x1=sqrt（10/（x+4））;
ifabs（x1-x）<
c=abs（x1-x）;
end
x=x1;
x
程序结果为：
x=1.36522998713958
Aitken加速后的新迭代格式求方程
functionwyc2a
f=inline（'
sqrt（10/（x+4））'
'
x'
）;
x1=f（f（x））-（（f（f（x））-f（x））^2）/（f（f（x））-2*f（x）+x）;
ifabs（x1-x）<
c=abs（x1-x）;
end
x=x1;
end
程序结果为：
x=1.36523001341410
通过对比两种迭代方法可以发现Aitken加速后的新迭代格式再不出现较大误差的情况下比普通的迭代速度更快，更加的高效。
三、解线性方程组
1．分别Jacobi迭代法和Gauss-Seidel迭代法求解线性方程组
迭代法计算停止的条件为：
．
Jacobi迭代法MATLAB程序：
functionwyc31j
D=diag（[6,5,8,7]）;
LU=[0,-2,-1,2;
-2,0,0,2;
2,0,0,-5;
-1,-3,-2,0;
];
B=inv（D）*LU;
b=[4;
7;
-1;
0;
f=inv（D）*b;
x=[0;
x1=[0;
x1=B*x+f;
fori=1:
3
ifabs（x1（i）-x（i））<
c=abs（x1（i）-x（i））;
end
x=x1;
x=[0.05205578320797
1.15094276230654
0.24463133345617
-0.57058744756179];
Gauss-Seidel迭代法MATLAB程序
functionwyc31gs
L=[0,0,0,0;
-2,0,0,0;
2,0,0,0;
-1,-3,-2,0;
U=[0,-2,-1,2;
0,0,0,2;
0,0,0,-5;
0,0,0,0;
B=inv（D-L）*U;
f=inv（D-L）*b;
ifabs（x1（i）-x（i））<
c=abs（x1（i）-x（i））;
a=[0.05204946628111
1.150********986
0.24463241176981
-0.57059206335719
Gauss-Seidel迭代法对已经计算出的信息加以充分利用，可以对Jacobi迭代法进行加速提高
2.用Gauss列主元消去法、QR方法求解如下方程组：
用Gauss列主元消去法求解方程MATLAB程序：
functionwyc32g
A1=[2,2,1,2;
4,1,3,-1;
-4,-2,0,1;
2,3,2,3;
b=[1;
2;
1;
0];
[n,m]=size（A1）;
x=zeros（n,1）;
A=[A1,b];
fork=1:
n-1
forp=k+1:
n
ifabs（A（p,k））>
abs（A（k,k））
forj=k:
n+1
t=A（k,j）;
A（k,j）=A（p,j）;
A（p,j）=t;
end
fori=k+1:
l=-A（i,k）/A（k,k）;
forj=k+1:
A（i,j）=A（i,j）+l*A（k,j）;
x（n）=A（n,n+1）/A（n,n）;
fori=n-1:
-1:
1
s=0;
forj=i+1:
s=s+A（i,j）*x（j）;
x（i）=（A（i,n+1）-s）/A（i,i）;
x=[1.54166666666667
-2.75000000000000
.0833********
1.66666666666667]
用QR方法求解方程MATLAB程序：
functionwyc32qr
A=[2,2,1,2;
2,3,2,3;
[n,m]=size（A）;
Q=zeros（m,n）;
R=zeros（m,n）;
[Q,R]=qr（A）;
b2=inv（Q）*b;
A2=[R,b2];
x（n）=A2（n,n+1）/A2（n,n）;
s=0;
s=s+A2（i,j）*x（j）;
x（i）=（A2（i,n+1）-s）/A2（i,i）;
x=[1.54166666666666
1.66666666666666]
四、已知一组数据点
，编写一程序求解三次样条插值函数
满足
并针对下面一组具体实验数据
0.25
0.3
0.39
0.45
0.53
0.5000
0.5477
0.6245
0.6708
0.7280
求解，其中边界条件为
.
functionwyc4
x=[0.25;
0.3;
0.39;
0.45;
0.53;
y=[0.5000;
0.5477;
0.6245;
0.6708;
0.7280;
s21=0;
s25=0;
[n,n1]=size（x）;
lamb=zeros（n,1）;
mu=zeros（n,1）;
h=zeros（n,1）;
g=zeros（n,1）;
m=zeros（n,1）;
s=zeros（4,1）;
symsp
h（i）=x（i+1）-x（i）;
forj=2:
lamb（j）=h（j）/（h（j-1）+h（j））;
mu（j）=h（j-1）/（h（j-1）+h（j））;
lamb（n）=1;
mu
（1）=1;
fork=2:
g（k）=3*（mu（k）*（y（k+1）-y（k））/h（k）+lamb（k）*（y（k）-y（k-1））/h（k-1））;
A=zeros（n,n）;
A（i,i）=2;
A（i,i+1）=mu（i）;
A（i+1,i）=lamb（i+1）;
A（n,n）=2;
g
（1）=3*（y
（2）-y
（1））/h
（1）-h
（1）*s21/2;
g（n）=3*（y（n）-y（n-1））/h（n-1）+h（n-1）*s25/2;
m=inv（A）*g;
4
（h（k）+2.*（p-x（k）））.*（（p-x（k+1））.^2）.*vpa（y（k）./（h（k）.^3））+（h（k）-2.*（p-x（k+1）））.*（（p-x（k））.^2）.*vpa（y（k+1）./（h（k）.^3））+（p-x（k））.*（（p-x（k+1））.^2）.*vpa（m（k）./（h（k）.^2））+（p-x（k+1））.*（（p-x（k））.^2）.*vpa（m（k+1）./（h（k）.^2））
当
s（p）=4000.*（-9/20+2*p）*（p-3/10）^2+4381.6*（13/20-2*p）*（p-1/4）^2+387.865164*（p-1/4）*（p-3/10）^2+369.069670*（p-3/10）*（p-1/4）^2
s（p）=751.303155*（-51/100+2*p）*（p-39/100）^2+856.652949*（87/100-2*p）*（p-3/10）^2+113.910392*（p-3/10）*（p-39/100）^2+98.670540*（p-39/100）*（p-3/10）^2
s（p）=2891.203704*（-18/25+2*p）*（p-9/20）^2+3105.555556*（24/25-2*p）*（p-39/100）^2+222.008716*（p-39/100）*（p-9/20）^2+206.2349887*（p-9/20）*（p-39/100）^2
s（p）=1310.15625*（-41/50+2*p）*（p-53/100）^2+1421.875*（57/50-2*p）*（p-9/20）^2+116.007181*（p-9/20）*（p-53/100）^2+109.574534*（p-53/100）*（p-9/20）^2
经过验证发现s（p）在各个区间段符合要求。
五、编写程序构造区间
上的以等分结点为插值结点的Newton插值公式，假设结点数为
（包括两个端点），给定相应的函数值，插值区间和等分的份数，该程序能快速计算出相应的插值公式。
以
，
为例计算其对应的插值公式，分别取不同的
值并画出原函数的图像以及插值函数的图像，观察当
增大时的逼近效果.
functionwyc5
n=4;
x=sym（'
real'
A=zeros（n+1,n）;
t=-1:
0.02:
z=1./（1+25.*（t.^2））;
1/（1+25*（y^2））'
y'
forj=1:
ifi==1
A（1,j）=2*（j-1）/（n-1）-1;
ifi==2
A（2,j）=f（A（1,j））;
if（i>
2）&
（i<
n+2）
ifj<
n+3-i
A（i,j）=（A（i-1,j+1）-A（i-1,j））/（（i-2）*2/（n-1））;
A
P3=ones（100）*A（2,1）;
q=zeros（n-1,1）;
100
x=0.02*k-1;
fori=1:
q（i）=A（i+2,1）;
i
q（i）=q（i）*（x-A（1,i））;
P3（k）=P3（k）+q（i）;
plot（P3）
当n=4时当n=6时
当n=8时当n=10时
当n=20时当n=30时
当n=50时
当n=100时
通过n取不同的值后观察发现当n的取值足够大时
区间内p（x）趋近于0，当x趋近于-1时p（x）趋近于无穷。
Newton插值公式不随插值点的增加而逼近，因此插值公式不适合对函数进行逼近。
升级会员 