數(shù)值分析實(shí)驗(yàn)報(bào)告

1、.數(shù)值分析實(shí)驗(yàn)報(bào)告 目錄第二章11.【第一題】1第三章32.【第二題】3第四章6第六章14第七章32第七章習(xí)題43第八章46第二章1. 【第一題】在區(qū)間-1,1上分別取n=10、20用兩組等距節(jié)點(diǎn)對(duì)龍格函數(shù)做多項(xiàng)式插值及三次樣條插值，對(duì)每個(gè)n值，分別畫(huà)出插值函數(shù)及的圖形。解：先編制Lagrange插值函數(shù)的m文件：function y1=lagrange(x0,y0,x1)n=length(x0);syms x;for k=1:n l(k)=x/x; for p=1:n if p=k l(k)=l(k)*(x-x0(p)/(x0(k)-x0(p); end endendz=0;for k=1:

2、n z=z+l(k)*y0(k);endy1=subs(z,x,x1);然后對(duì)n=10的情形做Lagrange插值并畫(huà)圖：x0=-1:.2:1;y0=1./(1+25*x0.2);x=-1:.001:1;y1=lagrange(x0,y0,x);y=1./(1+25*x.2);plot(x,y,x,y1,'-.',x0,y0,'p');legend('Runge Function','插值函數(shù)','插值節(jié)點(diǎn)');title('n=10時(shí)的Lagrange插值的龍格現(xiàn)象');xlabel('x

3、');ylabel('y');再對(duì)n=20的情形做Lagrange插值并畫(huà)圖：x0=-1:.1:1;y0=1./(1+25*x0.2);x=-1:.001:1;y1=lagrange(x0,y0,x);y=1./(1+25*x.2);plot(x,y,x,y1,'-.',x0,y0,'p');legend('Runge Function','插值函數(shù)','插值節(jié)點(diǎn)');title('n=20時(shí)的Lagrange插值的龍格現(xiàn)象');xlabel('x');yla

4、bel('y');第三章2. 【第二題】下列數(shù)據(jù)節(jié)點(diǎn)的插值：01491625364964012345678可得到平方根函數(shù)的近似。x=0 1 4 9 16 25 36 49 64;y=0 1 2 3 4 5 6 7 8;x0=0:.5:64;y0=lagrange(x,y,x0);y1=sqrt(x0);plot(x0,y1,x0,y0,'-.',x,y,'p')legend('開(kāi)方函數(shù)','插值函數(shù)','插值節(jié)點(diǎn)'); x=0 1 4 9 16 25 36 49 64;y=0 1 2 3 4 5

5、6 7 8;x0=0:.5:64;y0=spline(x,y,x0);y1=sqrt(x0);plot(x0,y1,x0,y0,'-.',x,y,'p')legend('開(kāi)方函數(shù)','三次條插值函數(shù)','插值節(jié)點(diǎn)');x=0.0 .1 .2 .3 .5 .8 1.0;y=1.0 .41 .50 .61 .91 2.02 2.46;>> p3=polyfit(x,y,3);>> p4=polyfit(x,y,4);>> x0=0:.001:1;y3=polyval(p3,x0);y

6、4=polyval(p4,x0);plot(x,y,'p',x0,y3,'-.',x0,y4,'-');legend('擬合節(jié)點(diǎn)','三次擬合','四次擬合');xlabel('x');ylabel('y');title('多項(xiàng)式擬合');grid on第四章1.【第一題】用不同的數(shù)值方法計(jì)算積分.在MATLAB中先求精確解：a=int('sqrt(x)*log(x)',0,1)a =-4/9(1)取不同的步長(zhǎng)h，分別用復(fù)合梯形與復(fù)合辛

7、普森求積公式計(jì)算積分，給出誤差中關(guān)于h的函數(shù)，并與精確積分值比較兩個(gè)公式的精度，是否存在一個(gè)最小的h，使得精度不能再改善？采用復(fù)合梯形公式for k=50:50:8000x=linspace(0,1,k+1);x=x(2:k+1);y=x.(1/2).*log(x);b=trapz(x,y);a=a;k,b;endformat longa(:,2)+4/9ans = ；第 32 頁(yè) 共 32 頁(yè)略對(duì)于復(fù)合辛普生公式，先編制復(fù)合辛普生公式的M文件：function q=combinesimpson(F,a,b,n)%復(fù)合Simpson求積公式%F被積函數(shù)%a,b積分區(qū)間%n區(qū)間份數(shù)x=linsp

8、ace(a,b,n+1);q=0;for k=1:n q=q+F(x(k)+4*F(x(k)+x(k+1)/2)+F(x(k+1);endq=q*(b-a)/n/6;用辛普森公式有：F=inline('x.(1/2).*log(x)');a=;for k=100:100:5000b=combinesimpson(F,1/k,1,k);a=a;k,b;end>> a(:,2)+4/9ans = 略（2）用龍貝格求積計(jì)算完成問(wèn)題（1）.function q,n,T=Romberg(f,x,a,b,tol)%Romberg求積公式%f被積函數(shù)%x積分變量%a,b積分區(qū)間%

9、tol積分精度%n迭代次數(shù)T=;k=0;h=b-a;T(0+1,0+1)=h/2*(subs(f,x,a)+subs(f,x,b);while k=0|abs(T(0+1,k+1)-T(0+1,k-1+1)>tol k=k+1;h=h/2;s=0; for jj=1:2:2k-1 s=s+subs(f,x,a+jj*h); end T(k+1,0+1)=1/2*T(k-1+1,0+1)+h*s; for jj=1:k T(k-jj+1,jj+1)=4jj/(4jj-1)*T(k-jj+1+1,jj-1+1)-1/(4jj-1)*T(k-jj+1,jj-1+1); endendn=k;q=

10、T(0+1,k+1);for jj=2:n+1 T(jj:n+1,jj)=T(1:n-jj+2,jj); T(1:jj-1,jj)=zeros(jj-1,1);End123（3）用自適應(yīng)辛普生積分，使其精度達(dá)到.>> F=inline('sqrt(x).*log(x)');>> s=quad(F,0,1,1e-4)s = -0.44402.【第二題】計(jì)算二重積分.（1）若區(qū)域，試分別用復(fù)合辛普森公式及Gauss求積公式（取n=4）求積分。對(duì)于多重積分的復(fù)合辛普生公式，先編制多重積分復(fù)合辛普生公式的M文件：function q=combinesimpson

11、2(F,x0,a,b,n)%復(fù)合Simpson多元求積公式%F被積函數(shù)%x0被積函數(shù)自變量%a,b積分區(qū)間%n區(qū)間份數(shù)x=linspace(a,b,n+1);q=0;for k=1:n q=q+subs(F,x0,x(k)+4*subs(F,x0,(x(k)+x(k+1)/2)+subs(F,x0,x(k+1);endq=q*(b-a)/n/6;用辛普森公式計(jì)算多重積分有：>> syms x y;F=exp(-x.*y);>> s=combinesimpson2(combinesimpson2(F,'x',0,1,4),y,0,1,4)s = 0.796

12、6對(duì)于多重積分的Gauss公式，先編制Gauss公式的M文件，其中要用到函數(shù)GAUSSPOINTS，使用格式及作用如下：GAUSSPOINTS points for Gaussian quadrature nodes, weights = gausspoints(nquad) computes the nodes and weights for Gaussian quadrature with nquad (<=5) quadrature points on the standard interval -1,1編制的M文件為：function q=gaussquad(F,x0,a,b,n

13、)%Gauss求積公式%F被積函數(shù)%x0被積函數(shù)自變量%a,b積分區(qū)間%n節(jié)點(diǎn)個(gè)數(shù)syms t;F=subs(F,x0,(b-a)/2*t+(a+b)/2);x,A=gausspoints(n);q=(b-a)/2*sum(A.*subs(F,t,x);用Gauss求積公式計(jì)算多重積分有：>> clear>> syms x y;F=exp(-x.*y);>> s=gaussquad(gaussquad(F,x,0,1,4),y,0,1,4)s =0.7966（2）若區(qū)域，用復(fù)合辛普森公式（取n=4）求此積分.>> syms x y;F=exp(-

14、x.*y);>> s=combinesimpson2(combinesimpson2(F,y,0,sqrt(1-x2),4),x,0,1,4)s = 0.6701第五章1.function A,b,det=liezhuyuan(A,b)det=1;m=size(A);n=m(1);u=zeros(1,n);for k=1:n-1 ik=k; for jj=k+1:n if abs(A(jj,k)>abs(A(ik,k) ik=jj; end end if A(ik,k)=0 x=0;det=0;disp('無(wú)解');return; end if ik=k fo

15、r jj=k:n u(jj)=A(k,jj);A(k,jj)=A(ik,jj);A(ik,jj)=u(jj); end v=b(k);b(k)=b(ik);b(ik)=v;det=-det; end for ii=k+1:n A(ii,k)=A(ii,k)/A(k,k); for jj=k+1:n A(ii,jj)=A(ii,jj)-A(ii,k)*A(k,jj); end b(ii)=b(ii)-A(ii,k)*b(k); endendif A(n,n)=0 x=0;det=0;disp('無(wú)解');return;endb(n)=b(n)/A(n,n);for ii=n-1:

16、-1:1 s=0; for jj=ii+1:n s=s+A(ii,jj)*b(jj); end b(ii)=(b(ii)-s)/A(ii,ii);enddet=A(n,n)*det;function A,b,det=liezhuyuan(A,b)det=1;m=size(A);n=m(1);u=zeros(1,n);for k=1:n-1 ik=k; for jj=k+1:n if abs(A(jj,k)>abs(A(ik,k) ik=jj; end end if A(ik,k)=0 x=0;det=0;disp('無(wú)解');return; end if ik=k for

17、 jj=k:n u(jj)=A(k,jj);A(k,jj)=A(ik,jj);A(ik,jj)=u(jj); end v=b(k);b(k)=b(ik);b(ik)=v;det=-det; end for ii=k+1:n A(ii,k)=A(ii,k)/A(k,k); for jj=k+1:n A(ii,jj)=A(ii,jj)-A(ii,k)*A(k,jj); end b(ii)=b(ii)-A(ii,k)*b(k); endendif A(n,n)=0 x=0;det=0;disp('無(wú)解');return;endb(n)=b(n)/A(n,n);for ii=n-1:-

18、1:1 s=0; for jj=ii+1:n s=s+A(ii,jj)*b(jj); end b(ii)=(b(ii)-s)/A(ii,ii);enddet=A(n,n)*det;function A,b,Ip=liezhuyuanLU(A,b)m=size(A);n=m(1);s=zeros(1,n);tt=s;Ip=zeros(n,1);if det(A)=0 disp('無(wú)解');return;endfor r=1:n for ii=r:n A(ii,r)=A(ii,r)-A(ii,r:n)*A(r:n,r); end u=max(abs(A(ii,r) if ir=r

19、tt=A(r,:);A(r,:)=A(ir,:);A(ir,:)=tt; end A(r,r)=s(r); if r=n for ii=r+1:n A(ii,r)=A(ii,r)/A(r,r); u=0; for k=1:r-1 u=u+A(r,k)*A(k,ii); end A(r,ii)=A(r,ii)-u; end endendfor ii=1:n-1 u=Ip(ii); if ii=u v=b(u);b(u)=b(ii);b(ii)=v; endendfor ii=2:n u=0; for k=1:ii-1 u=u+A(ii,k)*b(k); end b(ii)=b(ii)-u;end

20、b(n)=b(n)/A(n,n);for ii=n-1:-1:1 u=0; for k=ii+1:n u=u+A(ii,k)*b(k); end b(ii)=(b(ii)-u)/A(ii,ii);end第六章1.【第一題】給出線性方程組，其中系數(shù)矩陣為Hilbert矩陣，假設(shè)，.若取n=6,8,10，分別用Jacobi迭代法及SOR迭代（=1,1.25,1.5）求解，比較計(jì)算結(jié)果。解：首先編制Jacobi迭代求線性方程解的m文件jacobiiteract.m：function x=jacobiiteract(A,b,x,n,flag)%Jacobi迭代求線性方程Ax=b的解%x-迭代初值%n-

21、迭代次數(shù)%flag-控制過(guò)程是否顯示L=0 x' norm(A*x-b);h=size(A);p=h(2);for k=1:n for s=1:p w=0; for t=1:p if t=s w=w+A(s,t)*x(t); end end z(s)=(b(s)-w)/A(s,s); end x=z' L=L;k x' norm(A*x-b);endif flag Lend順便編制Gauss Seidel迭代求線性方程解的m文件gauss_seidel_iteract.m：function x=gauss_seidel_iteract(A,b,x,n,flag)%Gau

22、ss Seidel迭代求線性方程Ax=b的解%x-迭代初值%n-迭代次數(shù)%flag-控制過(guò)程是否顯示L=0 x' norm(A*x-b);h=size(A);p=h(2);for k=1:n for s=1:p w=0; for t=1:s-1 w=w+A(s,t)*x(t); end for t=s+1:p w=w+A(s,t)*x(t); end x(s)=(b(s)-w)/A(s,s); end L=L;k x' norm(A*x-b);endif flag Lend再編制超松弛迭代求線性方程解的m文件SOR.m：function x=SOR(A,b,x,omega,n,

23、flag)%超松弛迭代求線性方程Ax=b的解%x-迭代初值%omega-松弛因子%n-迭代次數(shù)%flag-控制過(guò)程是否顯示L=0 x' norm(A*x-b);h=size(A);p=h(2);for k=1:n for s=1:p w=0; for t=1:p w=w+A(s,t)*x(t); end x(s)=x(s)+omega*(b(s)-w)/A(s,s); end L=L;k x' norm(A*x-b);endif flag Lend編制文件xiti6_1.m：for k=6:2:10xstar=ones(k,1);H=hilb(k);b=H*xstar;x=0.

24、8*ones(k,1);kdisp('Jacobi迭代結(jié)果');J=jacobiiteract(H,b,x,10,1)for v=1:.25:1.5disp('松弛因子'); vdisp('Gauss Seigel迭代結(jié)果');G=SOR(H,b,x,v,10,1)endend執(zhí)行>> xiti6_1結(jié)果略第七章1.【第一題】求下列方程的實(shí)根：（1）；（2）.要求（1）設(shè)計(jì)一種不動(dòng)點(diǎn)迭代法，要使迭代序列收斂，然后用Steffensen加速迭代，計(jì)算到為止。>> x=-10:.1:10;y1=x.2-3*x+2;y2=exp

25、(x);>> plot(x,y1,x,y2)>> x=-10:.1:6;y1=x.2-3*x+2;y2=exp(x);plot(x,y1,x,y2)>> x=-2:.05:2;y1=x.2-3*x+2;y2=exp(x);plot(x,y1,x,y2)>> x=0.255:.00005:0.26;y1=x.2-3*x+2;y2=exp(x);plot(x,y1,x,y2);grid on>> x0=0.257:0.00001:0.258;syms x;>> f1=sqrt(3*x-2+exp(x);>> f2=

26、(x2+2*x-exp(x)/3;>> f3=log(x2-3*x+2);>> f4=(exp(x)-2)/(x-3);>> f5=(exp(x)-2)/x+3;>> y1=subs(diff(f1),x,x0);>> y2=subs(diff(f2),x,x0);>> y3=subs(diff(f3),x,x0);>> y4=subs(diff(f4),x,x0);>> y5=subs(diff(f5),x,x0);>> subplot(3,1,1);plot(x0,y1)>&g

27、t; subplot(3,1,2);plot(x0,y2)>> subplot(3,1,3);plot(x0,y3)>> subplot(2,1,1);plot(x0,y4)>> subplot(2,1,2);plot(x0,y5)這表明f1 f3 f5迭代不收斂，而f4的收斂性要稍比f(wàn)2好，所以采用f4迭代。根據(jù)作圖，迭代初值設(shè)為0.2575。function xx,n,L=fixedpoint_iteract(f,x,x0,tol)%不動(dòng)點(diǎn)迭代求解非線性方程%f-方程字符串%x-方程變量%x0-迭代初值%tol-迭代精度%xx-迭代結(jié)果%n-迭代次數(shù)%L

28、-迭代過(guò)程y=subs(f,x,x0);n=1;L=n,y,y-x0;while abs(y-x0)>tol x0=y; y=subs(f,x,x0);n=n+1; U=n,y,y-x0;L=L;U;endxx=y;>> xx n L=fixedpoint_iteract(f4,x,0.2575,1e-8);>> xxxx = 0.25753028364411>> nn = 10>> LL = 1.00000000000000 0.25754172799259 0.00004172799259 2.00000000000000 0.2575

29、2596199544 -0.00001576599715 3.00000000000000 0.25753191898097 0.00000595698553 4.00000000000000 0.25752966823006 -0.00000225075090 5.00000000000000 0.25753051864311 0.00000085041304 6.00000000000000 0.25753019732752 -0.00000032131559 7.00000000000000 0.25753031873177 0.00000012140425 8.000000000000

30、00 0.25753027286101 -0.00000004587076 9.00000000000000 0.25753029019259 0.00000001733157 10.00000000000000 0.25753028364411 -0.00000000654847采用Steffensen加速迭代：function xx,n,L=steffensen_iteract(f,x,x0,tol)%steffensen迭代求解非線性方程%f-方程字符串%x-方程變量%x0-迭代初值%tol-迭代精度%xx-迭代結(jié)果%n-迭代次數(shù)%L-迭代過(guò)程F=x-(f-x)2/(subs(f,x,f

31、)-2*f+x);y=subs(F,x,x0);n=1;L=n,y,y-x0;while abs(y-x0)>tol x0=y; y=subs(F,x,x0);n=n+1; U=n,y,y-x0;L=L;U;endxx=y;>> xx n L=steffensen_iteract(f4,x,0.2575,1e-8);>> xxxx = 0.25753028543986>> nn = 2>> LL = 1.00000000000000 0.25753028534589 0.00003028534589 2.00000000000000 0.2

32、5753028543986 0.00000000009397（2）用牛頓迭代，同樣計(jì)算到。輸出迭代初值及各迭代值和迭代次數(shù)k，比較方法的優(yōu)劣。用Newton迭代法function xx,n,L=newton_iteract(f,x,x0,tol)%牛頓迭代求解非線性方程%f-方程字符串%x-方程變量%x0-迭代初值%tol-迭代精度%xx-迭代結(jié)果%n-迭代次數(shù)%L-迭代過(guò)程F=x-f/diff(f,x,1);y=subs(F,x,x0);n=1;L=n,y,y-x0;while abs(y-x0)>tol x0=y; y=subs(F,x,x0);n=n+1; U=n,y,y-x0;L

33、=L;U;endxx=y;>> f=x2-3*x+2-exp(x);>> xx n L=newton_iteract(f,x,0.2575,1e-8);>> xxxx = 0.25753028543986>> nn = 2>> LL = 1.00000000000000 0.25753028535414 0.00003028535414 2.00000000000000 0.25753028543986 0.00000000008572（1）>> x=-10:.1:10;y=x.3+2*x.2+10*x-20; plot(

34、x,y),grid on>> x=-4:.05:4;y=x.3+2*x.2+10*x-20; plot(x,y),grid on>> x=1:.005:2;y=x.3+2*x.2+10*x-20; plot(x,y),grid on>> x=1.3:.0001:1.4;y=x.3+2*x.2+10*x-20; plot(x,y),grid on>> x=1.36:.00005:1.38;y=x.3+2*x.2+10*x-20; plot(x,y),grid on>> x=1.365:.00001:1.370;y=x.3+2*x.2+1

35、0*x-20; plot(x,y),grid on>> x0=1.3685:0.000001:1.3690;syms x;>> f1=20/(x2+2*x+10);y1=subs(diff(f1),x,x0);>> plot(x0,y1)>> xx n L=fixedpoint_iteract(f1,x,1.3688,1e-8);>> xxxx = 1.36880810541671>> nn = 10>> LL =略>> f=x3+2*x2+10*x-20;>> xx n L=newto

36、n_iteract(f,x,1.3688,1e-8);>> xxxx = 1.36880810782137>> nn = 2>> LL = 1.00000000000000 1.36880810784040 0.00000810784040 2.00000000000000 1.36880810782137 -0.000000000019032.>> r=1 2 3 4 5 6 7 8 9 10;p=poly(r)p =略>> F=p;F(2)=F(2)+1e-6;>> x1=roots(F)x1 =略>> F

37、=p;F(2)=F(2)+1e-8;>> x2=roots(F)x2 =略>> F=p;F(2)=F(2)+1e-10;>> x3=roots(F)x3 =略第七章習(xí)題1. function xx,n,L=bicut(f,x,x1,x2,tol)%二分法迭代求解非線性方程%f-方程字符串%x-方程變量%x1,x2-迭代區(qū)間%tol-迭代精度%xx-迭代結(jié)果%n-迭代次數(shù)%L-迭代過(guò)程if subs(f,x,x1)*subs(f,x,x2)>=0 disp('迭代區(qū)間選取不合適，請(qǐng)重新選取');return;endn=0;L=;tol=m

38、in(abs(x2-x1)/2,tol);while abs(x2-x1)>tol xx=(x1+x2)/2;n=n+1; if subs(f,x,x1)*subs(f,x,xx)>0 x1=xx; elseif subs(f,x,xx)=0 return; else x2=xx; end U=n,xx,subs(f,x,xx);L=L;U;endsyms x;f=x2-x-1;xx n L=bicut(f,x,0,2,.05);xx = 1.59375000000000n = 6L = 1.00000000000000 1.00000000000000 -1.0000000000

39、0000 2.00000000000000 1.50000000000000 -0.25000000000000 3.00000000000000 1.75000000000000 0.31250000000000 4.00000000000000 1.62500000000000 0.01562500000000 5.00000000000000 1.56250000000000 -0.12109375000000 6.00000000000000 1.59375000000000 -0.053710937500002.（1）>> syms x;>> f=1+1/x2

40、;xx n L=fixedpoint_iteract(f,x,1.5,1e-5)略（2）>> f=(1+x2)(1/3);xx n L=fixedpoint_iteract(f,x,1.5,1e-5)略（3）function xx,n,L=fixedpoint_iteract(f,x,x0,tol)%不動(dòng)點(diǎn)迭代求解非線性方程%f-方程字符串%x-方程變量%x0-迭代初值%tol-迭代精度%xx-迭代結(jié)果%n-迭代次數(shù)%L-迭代過(guò)程if abs(subs(diff(f,x,1),x,x0)>=1 disp('迭代發(fā)散'); xx=x0; n=0; L=; ret

41、urn;endy=subs(f,x,x0);n=1;L=n,y,y-x0;uu=x0;while abs(y-x0)>tol x0=y; y=subs(f,x,x0); if abs(y-uu)>5|n>1e4 disp('迭代發(fā)散');xx=y;return; end n=n+1; U=n,y,y-x0;L=L;U;endxx=y;>> f=1/sqrt(x-1);xx n L=fixedpoint_iteract(f,x,1.5,1e-5)迭代發(fā)散第八章2【第二題】function V,lambda,L=mifa(A,n)%冪法求矩陣A的主特征

42、值lambda與其對(duì)應(yīng)的特征向量V%n迭代次數(shù)%L迭代過(guò)程k=size(A);v=ones(k(1),1);u=v;mu=max(v);L=0,u',mu;for k=1:n v=A*u;mu=max(v);u=v/mu;L=L;k,u',mu;endV=u;lambda=mu;>> A=1 1 .5;1 1 .25;.5 .25 2;>> V lambda L=mifa(A,20)略>> A=5 4 1 1;4 5 1 1;1 1 4 2;1 1 2 4;format long>> V lambda L=mifa(A,30)略f

43、unction V,lambda,L=mifa_rayleigh(A,n)%冪法求矩陣A的主特征值lambda與其對(duì)應(yīng)的特征向量V,并用Raleigh商加速%n迭代次數(shù)%L迭代過(guò)程k=size(A);v=ones(k(1),1);u=v;mu=max(v);L=0,u',mu;for k=1:n v=A*u;mu=(u'*A*u)/(u'*u); u=v/mu; L=L;k,u',mu;endV=u;lambda=mu;>> V lambda L=mifa_rayleigh(A,30)略function V,lambda,L,n=mifa1(A,to

44、l)%冪法求矩陣A的主特征值lambda與其對(duì)應(yīng)的特征向量V%n迭代次數(shù)%L迭代過(guò)程%tol特征值誤差限k=size(A);v=ones(k(1),1);u=v;mu0=max(v);L=0,u',mu0;v=A*u;mu=max(v);u=v/mu;L=L;1,u',mu;n=1;while abs(mu-mu0)>tol n=n+1; mu0=mu; v=A*u;mu=max(v);u=v/mu;L=L;n,u',mu;endV=u;lambda=mu;>> V lambda L n=mifa1(A,1e-8);>> nn =28fun

45、ction V,lambda,L,n=mifa_rayleigh1(A,tol)%冪法求矩陣A的主特征值lambda與其對(duì)應(yīng)的特征向量V,并用Raleigh商加速%n迭代次數(shù)%L迭代過(guò)程%tol特征值誤差限k=size(A);v=ones(k(1),1);u=v;mu0=max(v);L=0,u',mu0;v=A*u;mu=(u'*A*u)/(u'*u);u=v/mu;L=L;1,u',mu;n=1;while abs(mu-mu0)>tol n=n+1; mu0=mu; v=A*u;mu=(u'*A*u)/(u'*u); u=v/mu; L=L;n,u',mu;endV=u;lambda=mu;>> V lambda L n=mifa_rayleigh1(A,1e-8);>> nn = 15第九章function T,X=ode_euler(f,t