《應(yīng)用計(jì)算方法教程》matlab作業(yè)二

1、作業(yè)六6-1 試驗(yàn)?zāi)康?計(jì)算特征值，實(shí)現(xiàn)算法試驗(yàn)內(nèi)容：隨機(jī)產(chǎn)生一個(gè)10階整數(shù)矩陣，各數(shù)均在-5和5之間。（1） 用MATLAB函數(shù)“eig”求矩陣全部特征值。（2） 用冪法求A的主特征值及對(duì)應(yīng)的特征向量。（3） 用基本QR算法求全部特征值（可用MATLAB函數(shù)“qr”實(shí)現(xiàn)矩陣的QR分解）。原理冪法：設(shè)矩陣A的特征值為并設(shè)A有完全的特征向量系(它們線性無(wú)關(guān))，則對(duì)任意一個(gè)非零向量所構(gòu)造的向量序列有，其中表示向量的第j個(gè)分量。為避免逐次迭代向量不為零的分量變得很大（ 時(shí)）或很?。?時(shí)），將每一步的按其模最大的元素進(jìn)行歸一化。具體過(guò)程如下：選擇初始向量，令，當(dāng)充分大時(shí)。QR法求全部特征值： 由于此題

2、的矩陣是10階的，上述算法計(jì)算時(shí)間過(guò)長(zhǎng)，考慮采用改進(jìn)算法移位加速。迭代格式如下： 計(jì)算右下角的二階矩陣的特征值，當(dāng)為實(shí)數(shù)時(shí)，選為中最接近的。程序A=-5+round(10*rand(10);V,D=eig(A)lamda u=lab6_2_power(A,1;1;1;1;1;1;1;1;1;1,10(-5),1000)d=lab6_3_qr2(A,10(-5)function lamda u=lab6_2_power(a,v,eps,N)lamda=0;err=1;k=1;while(keps) u=a*v; m j=max(abs(u); dc=abs(lamda-m); u=u/m; dv

3、=norm(u-v); err=max(dc,dv); v=u; lamda=m; k=k+1;endfunction D=lab6_3_qr2(A,eps)n,n=size(A);m=n;D=zeros(n,1);B=A;while(m1) while(abs(B(m,m-1)=eps*(abs(B(m-1,m-1)+abs(B(m,m) S=eig(B(m-1:m,m-1:m); j,k=min(abs(B(m,m)-S(1),abs(B(m,m)-S(2); Q,U=qr(B-S(k)*eye(m); B=U*Q+S(k)*eye(m); end A(1:m,1:m)=B; m=m-1;

4、 B=A(1:m,1:m);endD=diag(A);界面(1)(2)(3)作業(yè)七7-1 試驗(yàn)?zāi)康模菏煜ご鷶?shù)插值試驗(yàn)內(nèi)容：已知在f(x)在7個(gè)點(diǎn)的函數(shù)值如下表所示，分別使用拉格朗日插值法和牛頓插值法求f(0.596)與f(0.906)的近似值。0.40.50.60.70.80.91.011.751.962.192.442.713.00原理拉格朗日插值多項(xiàng)式： 牛頓插值多項(xiàng)式： 其中。程序function y1=lab7_1_Lagrange(x,y,x1)y1=0;m n=size(x);n=n-1;for k=1:n+1 t=1; for i=1:n+1 if(i=k) t=t*(x1-x(

5、i)/(x(k)-x(i); end end y1=y1+t*y(k);endfunction y1=lab7_2_Newton(x,y,x1)m n=size(x);n=n-1;for j=1:n for i=n+1:-1:j+1 y(i)=(y(i)-y(i-1)/(x(i)-x(i-j); endendy1=y(n+1);for j=n:-1:1 y1=y(j)+(x1-x(j)*y1;end界面作業(yè)八8-1 試驗(yàn)?zāi)康模菏煜ぷ钚《朔〝M合多項(xiàng)式試驗(yàn)內(nèi)容：給定數(shù)據(jù)點(diǎn)(,)，X0.400.550.650.800.901.05f(x)0.410750.578150.696750.888111.

6、026521.25386用3次最小二乘多項(xiàng)式擬合數(shù)據(jù)，并求平方誤差。原理要作三次最小二乘擬合，令，計(jì)算法方程，其中，。其平方誤差為。 程序x=0.4 0.55 0.65 0.80 0.90 1.05;f=0.41075 0.57815 0.69675 0.88811 1.02652 1.25386;G=zeros(4,4);for j=1:4 for k=1:4 for i=1:6 G(j,k)=G(j,k)+x(i)(j+k-2); end endendd=zeros(4,1);for k=1:4 for i=1:6 d(k,1)=d(k,1)+f(i)*x(i)(k-1); endenda

7、=Gds=0;for i=1:6 s=s+(f(i)-(a(1)+a(2)*x(i)+a(3)*x(i)2+a(4)*x(i)3)2;ends界面作業(yè)九9-1 試驗(yàn)?zāi)康模菏煜?shù)值積分公式，掌握數(shù)值計(jì)算定積分的方法試驗(yàn)內(nèi)容：采用不同方法數(shù)值計(jì)算積分 編寫復(fù)合梯形公式和復(fù)合Simpson公式通用子程序，分別采用4，8，16，32，64等分區(qū)間計(jì)算。原理復(fù)合梯形公式：將區(qū)間a,b作n等分，結(jié)點(diǎn)，復(fù)合Simpson公式：將區(qū)間a,b作2n等分，記， 程序function y=lab9_f(x)y=(log(1+x)/x;function y=lab9_1_fTrapezoid(a,b,eps,n)f0

8、=0;h=(b-a)/n;for i=0:n-1 f0=f0+lab9_f(a+i*(b-a)/n)+lab9_f(a+(i+1)*(b-a)/n);endy=f0*h/2;function y=lab9_2_fSimpson(a,b,eps,n)f0=0;h=(b-a)/n;k=n/2;for i=0:k-1 f0=f0+lab9_f(a+2*i*(b-a)/n)+4*lab9_f(a+(2*i+1)*(b-a)/n)+lab9_f(a+(2*i+2)*(b-a)/n);endy=f0*h/3;界面作業(yè)十10-1 試驗(yàn)?zāi)康模簩W(xué)會(huì)用Euler法、改進(jìn)Euler法、經(jīng)典的4階Runge-Kutt

9、a法求解常微分方程初值問(wèn)題。試驗(yàn)內(nèi)容：分別用1) Euler法（步長(zhǎng)h=0.025）2) 改進(jìn)Euler法（步長(zhǎng)h=0.05）3) 4階Runge-Kutta（步長(zhǎng)h=0.1）求解下面的初值問(wèn)題： 比較公共節(jié)點(diǎn)解的誤差。精確解為。原理Euler法：令，（）。改進(jìn)Euler法（梯形公式）： 4階Runge-Kutta：程序y0=-2;h=0.025;n=2/h;y(1)=y0-h*(y0+1)*(y0+3);for i=2:n y(i)=y(i-1)-h*(y(i-1)+1)*(y(i-1)+3);endyfor j=1:n y1(j)=-3+2/(1+exp(-4*j/n);endy1x=0.

10、025:0.025:2;plot(x,y1,r,x,y,b)y0=-2;h=0.05;n=2/h;k1=-(y0+1)*(y0+3);k2=-(y0+h*k1+1)*(y0+h*k1+3);y21(1)=y0+h/2*(k1+k2);for i=2:n k1=-(y21(i-1)+1)*(y21(i-1)+3); k2=-(y21(i-1)+h*k1+1)*(y21(i-1)+h*k1+3); y21(i)=y21(i-1)+h/2*(k1+k2);endy21for j=1:n y22(j)=-3+2/(1+exp(-4*j/n);endy22x=0.05:0.05:2;plot(x,y21,r,x,y22,b)y0=-2;h=0.1;n=2/h;k1=-h*(y0+1)*(y0+3);k2=-h*(y0+k1/2+1)*(y0+k1/2+3);k3=-h*(y0+k2/2+1)*(y0+k2/2+3);k4=-h*(y0+k3+1)*(y0+k3+3);y31(1)=y0+(k1+2*k2+2*k3+k4)/6;for i=2:n k1=-h*(y31(i-1)+1)*(y31(i-1)+3); k2=-h*(y31(i-1)+k1/2+1)*(y31(i-1)+k1/2+3); k3=-h*(y31(i-1)+k2/2+1)*(y31(i