资源描述：
《解常微分方程组(system》由会员上传分享，免费在线阅读，更多相关内容在应用文档-天天文库。
1、第六章解常微分方程組(SystemsofOrdinaryDiff.Equations)本章先解一階聯立常微分方程式,其數值方法與解一階常微分方程式類似.然後再利用解一階聯立常微分方程式的方法,解高階的常微分方程式(組).在本章中包含Matlab的m-file1.taylor4sys.m2.rk4sys.m3.rk4sysnew.m;fxsys.m4.adamsmsys.m將須要的m-file之檔案夾加入搜尋路徑中path('d:umerical',path)註:如果你有安裝MatlabNotebook要執行下列inputcells(綠色指令敘述)之前必須先執行上面的cell–[path
2、(…)]藍色的內容是Matlab[outputcells]1.taylor4sys.m–利用泰勒級數方法取到第4階X(t+h)=X(t)+hX'(t)+h2/2X''(t)+h3/3!X'''(t)+h4/4!X(4)(t)來估算聯立常微分方程式X'(t)=F(X,t)的數值解顯示taylor4sys.m的內容typetaylor4sys.mfunctionrs=taylor4sys(X0,a,b,m)%toreturntheapproximationvaluesofX(t)%byTaylorseriesmethodoforder4,X0istheinitial%tisin[a,b]n=l
3、ength(X0);rs=zeros(n,m);d1=zeros(n,1);d2=d1;d3=d1;d4=d1;t=a;X=X0;h=(b-a)/m;fork=1:m%foryourF(t,X)computethederivativesofX',X'',X'''andX4attd1(1)=X(1)-X(2)+t*(2-t*(1+t));%d1isX'd1(2)=X(1)+X(2)+t^2*(-4+t);d2(1)=d1(1)-d1(2)+2-t*(2+3*t);%d2isX''d2(2)=d1(1)+d1(2)+t*(-8+3*t);d3(1)=d2(1)-d2(2)-2-6*t;%d3i
4、sX'''d3(2)=d2(1)+d2(2)-8+6*t;d4(1)=d3(1)-d3(2)-6;%d4isX4d4(2)=d3(1)+d3(2)+6;%computeX(t+h)fori=1:nX(i)=X(i)+h*(d1(i)+h*(d2(i)+h*(d3(i)+h*d4(i)/4)/3)/2);rs(i,k)=X(i);end%forit=t+h;end%forkchapter6page9例題1:解聯立常微分方程式x'(t)=x(t)-y(t)+2t-t2–t3;y'(t)=x(t)+y(t)-4t2+t3初值x(0)=1,y(0)=0,0≦t≦1a=0;b=1;m=10;X0=[
5、1;0];Rs=taylor4sys(X0,a,b,m);Rs=[X0Rs]%增加初值的資料Rs=Columns1through71.00001.10971.23711.37961.53411.69691.863900.10930.23470.37190.51690.66540.8129Columns8through112.03022.19052.33892.46870.95431.08451.19771.2874比較上面的結果與真正解x(t)=etcos(t)+t2;y(t)=etsin(t)–t3;並畫出圖形t=0:0.1:1;sol1=exp(t).*cos(t)+t.^2sol2=
6、exp(t).*sin(t)-t.^3tt=0:0.05:1;xt=exp(tt).*cos(tt)+tt.^2;yt=exp(tt).*sin(tt)-tt.^3;plot(Rs(1,:),Rs(2,:),'r*',xt,yt)sol1=Columns1through71.00001.10961.23711.37961.53411.69691.8639Columns8through112.03022.19052.33892.4687sol2=Columns1through700.10930.23470.37190.51690.66540.8128Columns8through110.95
7、431.08451.19771.2874chapter6page9依據上面的數據我們發現利用taylor4sys得到的結果,可正確到小數點第四位,此外從圖形上顯示結果亦算蠻好的.2.rk4sys.m--利用Runge-Kutta方法X(t+h)=X(t)+h/6(K1+2K2+2K3+K4),K1=F(t,X),K2=F(t+h/2,X+h/2K1),K2=F(t+h/2,X+h/2K2),K4=F(t+h,X+hK3)
