矩阵与数值分析作业.docx
- 文档编号:8002686
- 上传时间:2023-01-27
- 格式:DOCX
- 页数:12
- 大小:489.91KB
矩阵与数值分析作业.docx
《矩阵与数值分析作业.docx》由会员分享,可在线阅读,更多相关《矩阵与数值分析作业.docx(12页珍藏版)》请在冰豆网上搜索。
矩阵与数值分析作业
2011级工科硕士研究生
《矩阵与数值分析》课程数值实验题目
一、对于数列
,有如下两种生成方式
1、首项为
,递推公式为
;
2、前两项为
,递推公式为
;
给出利用上述两种递推公式生成的序列的第50项。
解:
matlab编程
1、
c1
(1)=1;
fort1=2:
50
c1(t1)=c1(t1-1)/3;
end
c1(50)
运算结果:
ans=4.1788e-024
2、
c2
(1)=1;c2
(2)=1/3;
fort2=3:
50
c2(t2)=c2(t2-1)*10/3-c2(t2-2);
end
c2(50)
运算结果:
ans=-4.9661e+006
二、利用迭代格式
及Aitken加速后的新迭代格式求方程
在
内的根
解:
matlab编程
k=1;x
(1)=1.3;f=1;
while(abs(f)>1.0e-6)
x(k+1)=sqrt(10/(x(k)+4));
k=k+1;
f=x(k)^3+4*x(k)^2-10;
end
x,f
clear
k=1;x
(1)=1.3;f=1;
while(abs(f)>1.0e-6)
y(k)=sqrt(10/(x(k)+4));
z(k)=sqrt(10/(y(k)+4));
x(k+1)=z(k)-(z(k)-y(k))^2/(z(k)-2*y(k)+x(k));
k=k+1;
f=x(k)^3+4*x(k)^2-10;
end
x,f
运算结果:
x=1.30001.37351.36411.36531.36521.36511.36511.3651
f=5.8603e-007
x=1.30001.36511.3651
f=2.5046e-012
分析:
原迭代格式运算8次,达到精度要求;而是用Aitken加速后,只用3次迭代,即Aitken加速可以明显加快迭代速度。
三、解线性方程组
1.分别Jacobi迭代法和Gauss-Seidel迭代法求解线性方程组
迭代法计算停止的条件为:
.
解:
matlab编程
1
Xj=[0;0;0;0];Xg=[0;0;0;0];
A=[621-2;
250-2;
-2085;
1327];
b=[4;7;-1;0];
D=diag(diag(A));
L=-tril(A)+D;
U=-triu(A)+D;
Bj=inv(D)*(L+U);
f=inv(D)*b;
k=1;flag=1;
Xj(:
k+1)=Bj*Xj(:
k)+f;
whileflag==1
ifabs(Xj(:
k+1)-Xj(:
k))<1e-6
flag=0;
else
k=k+1;
Xj(:
k+1)=Bj*Xj(:
k)+f;
end
end
Xj
Bg=inv(D-L)*U;
fg=inv(D-L)*b;
k=1;flag=1;
Xg(:
k+1)=Bg*Xg(:
k)+fg;
whileflag==1
ifabs(Xg(:
k+1)-Xg(:
k))<1e-6
flag=0;
else
k=k+1;
Xg(:
k+1)=Bg*Xg(:
k)+fg;
end
end
Xg
运算结果:
Xj=Columns1through11
00.66680.22080.06220.08380.05670.05370.05350.05220.05220.0521
01.40001.13331.04781.16351.14421.14781.15151.15051.15091.1509
0-0.12500.04160.34240.22120.24330.24780.24320.24480.24470.2446
00-0.6593-0.5292-0.5557-0.5738-0.5681-0.5704-0.5706-0.5705-0.5706
Columns12through16
0.05210.05210.05210.05200.0520
1.15091.15091.15091.15091.1509
0.24450.24450.24450.24450.2445
-0.5707-0.5707-0.5707-0.5707-0.5707
Xg=
00.66670.08430.05510.05230.05210.05210.05200.0520
01.13321.12911.14901.15081.15081.15091.15091.1509
00.04170.26660.24640.24480.24460.24460.24460.2446
0-0.5928-0.5721-0.5707-0.5706-0.5706-0.5706-0.5706-0.5706
2.用Gauss列主元消去法、QR方法求解如下方程组:
2、
A=[2212;
413-1;
-4-201;
2323];
b=[1;2;1;0];
Ab=[Ab];
forj=1:
3
k=j;
max(j)=Ab(j,j);
fori=(j+1):
4
ifmax(j) k=i; max(j)=Ab(i,j); end end ifk~=j t=Ab(j,: ); Ab(j,: )=Ab(k,: ); Ab(k,: )=t; end fori=(j+1): 4 Ab(i,: )=Ab(i,: )-Ab(j,: )*Ab(i,j)/Ab(j,j); end end fori=4: -1: 1 Abt=Ab(i,5); forj=4: -1: (i+1) Abt=Abt-x(j)*Ab(i,j); end x(i)=Abt/Ab(i,i); end Xg=x' [Q,R]=qr(A); Xqr=R\Q'*b 运算结果: Xg= 1.5417 -2.7500 0.0832 1.6667 Xqr= 1.5417 -2.7500 0.0832 1.6667 分析: G-S迭代法与Jacob迭代法相比结果相同,但减少很多迭代步数。 四、已知一组数据点 ,编写一程序求解三次样条插值函数 满足 并针对下面一组具体实验数据 0.25 0.3 0.39 0.45 0.53 0.5000 0.5477 0.6245 0.6708 0.7280 求解,其中边界条件为 . 解: n=5; X=[0.250.30.390.450.53]; Y=[0.50000.54770.62450.67080.7280]; f00=0; fnn=0; fork=1: (n-1) h(k)=X(k+1)-X(k); end g (1)=3*(Y (2)-Y (1))/h (1)-f00*h (1)/2; g(n)=3*(Y(n)-Y(n-1))/h(n-1)+fnn*h(n-1)/2; fork=2: (n-1) v(k)=h(k)/(h(k)+h(k-1)); u(k)=h(k-1)/(h(k)+h(k-1)); g(k)=3*u(k)*(Y(k+1)-Y(k))/h(k)+3*v(k)*(Y(k)-Y(k-1))/h(k-1); end u (1)=1; v (1)=[];v(n-1)=1; A=2*eye(n)+diag(u,1)+diag(v,-1); m=A\g'; symsxSx; fork=1: (n-1) fprintf('%g Sx=(h(k)+2*(x-X(k)))*(x-X(k+1))^2*Y(k)/h(k)^3+... (h(k)-2*(x-X(k+1)))*(x-X(k))^2*Y(k+1)/h(k)^3+... (x-X(k))*(x-X(k+1))^2*m(k)/h(k)^2+... (x-X(k+1))*(x-X(k))^2*m(k+1)/h(k)^2; Sx=vpa(expand(Sx),4) End 运算结果: 0.25 Sx= 4.699*x^2-0.2051*x+0.3555-6.265*x^3 0.3 Sx= -2.633*x^2+1.995*x+0.1355+1.881*x^3 0.39 Sx= 0.1064*x^2+0.9261*x+0.2744-0.4600*x^3 0.45 Sx= -3.409*x^2+2.508*x+0.3710e-1+2.144*x^3 分析: 三次样条函数具有良好的光滑性和几何光顺性,同时收敛性和稳定性也较好。 五、编写程序构造区间 上的以等分结点为插值结点的Newton插值公式,假设结点数为 (包括两个端点),给定相应的函数值,插值区间和等分的份数,该程序能快速计算出相应的插值公式。 以 , 为例计算其对应的插值公式,分别取不同的 值并画出原函数的图像以及插值函数的图像,观察当 增大时的逼近效果. 解: matlab编程 functionLagrange clc;clear;closeall; fori=1: 3; ifi==1 N=5; elseifi==2 N=8; else N=10; end f=inline('1/(1+25*x^2)'); x1=zeros(1,N+1); a=-1; b=1; fori=1: N+1 x1(i)=a+(i-1)*(b-a)/N; y1(i)=feval(f,x1(i)); end symsx ff=0; fori=1: N+1 f=1; forj=1: i-1 f=f*(x-x1(j))/(x1(i)-x1(j)); end forj=i+1: N+1 f=f*(x-x1(j))/(x1(i)-x1(j)); end ff=f*y1(i)+ff; f=1; end ff=collect(ff,x); ff=vpa(ff,4); y=ff; p=ezplot(y,[a,b]);grid YLIM([-0.10.6]); ifN==5 set(p,'Color','black'); set(p,'LineStyle','--'); lagrange_5=y elseifN==8 set(p,'Color','r'); set(p,'LineStyle','--'); lagrange_8=y else set(p,'Color','b'); set(p,'LineStyle','--') lagrange_10=y end holdon; xlabel('x');ylabel('y'); title('y=p(x)');holdon; end Lag_Cheb(); x=-1: 0.01: 1; y=1./(1+25*x.^2); acu=plot(x,y);gridon;holdon set(acu,'Color','m'); set(acu,'LineStyle','-'); legend('N=5','N=8','N=10','Cheb,N=10','¾«È·Öµ'); functionLag_Cheb() f=inline('1/(1+25*x^2)'); N=10; x1=zeros(1,N+1); a=-1; b=1; fori=1: N+1 x1(i)=cos((2*i-1)*pi/(2*N)); y1(i)=feval(f,x1(i)); end symsx ff=0; fori=1: N+1 f=1; forj=1: i-1 f=f*(x-x1(j))/(x1(i)-x1(j)); end forj=i+1: N+1 f=f*(x-x1(j))/(x1(i)-x1(j)); end ff=f*y1(i)+ff; f=1; end ff=collect(ff,x); ff=vpa(ff,4); yy=ff; ff=collect(ff,x); yy=ff; lagrange_chebshev_10=yy cheb=ezplot(yy,[a,b]);gridon YLIM([-0.10.6]); set(cheb,'Color','g'); set(cheb,'LineStyle',': '); 结果输出 lagrange_5=0.5673+1.202*x^4-1.731*x^2 lagrange_8=1.+53.69*x^8-102.8*x^6+61.37*x^4-13.20*x^2 lagrange_10=1.-220.9*x^10+494.9*x^8-381.4*x^6+123.4*x^4-16.86*x^2 lagrange_chebshev_10=0.7413-5.359*x^10-0.4e-2*x^9+18.96*x^8-0.1321*x^7- 25.78*x^6+0.10*x^5+16.81*x^4+0.5e-2*x^3-5.336*x^2+0.5288e-3*x 分析: (1)拉格朗日插值法随着插值节点增多而更加接近真实曲线,但是还是会在区间 [-1,-0.4]及[0.4,1]上,由于高次振荡产生Runge现象。 (2)分段插值法则优于前者,因为分段插值在每个小区间上都可以最大限度的满足边界条件,但光滑度较低。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 矩阵 数值 分析 作业