微分方程數值解法的matlab程序

歐拉法

clear clc syms x y f=2*x; a=0; b=10; y0=0; h=0.1; n=(b-a)/h+1; xx=zeros(n,1); yy=zeros(n,1); for i=1:nxx(i)=a+h*(i-1); end yy=zeros(n,1); yy(1)=y0; for i=2:nyy(i)=yy(i-1)+h*subs(f,{x,y},{xx(i-1),yy(i-1)}); end ff1=[xx,yy];

改進歐拉法

clear clc syms x y f=3*x^2; a=0; b=10; y0=0; h=0.1; n=(b-a)/h+1; xx=zeros(n,1); yy=zeros(n,1); for i=1:nxx(i)=a+h*(i-1); end yy=zeros(n,1); yy(1)=y0; for i=2:nyy(i)=yy(i-1)+h*subs(f,{x,y},{xx(i-1),yy(i-1)});p=zeros(4,1);p(1)=yy(i);for j=2:4p(j)=yy(i-1)+(h/2)*(subs(f,{x,y},{xx(i-1),yy(i-1)})+subs(f,{x,y},{xx(i),p(j-1)}));endyy(i)=p(4); end ff1=[xx,yy];

龍格庫塔法

clear clc syms x y f=3*x^2; a=0; b=10; y0=0; h=0.1; n=(b-a)/h+1; xx=zeros(n,1); yy=zeros(n,1); for i=1:nxx(i)=a+h*(i-1); end yy=zeros(n,1); yy(1)=y0; for i=2:nk1=subs(f,{x,y},{xx(i-1),yy(i-1)});k2=subs(f,{x,y},{xx(i-1)+h/2,yy(i-1)+h/2*k1});k3=subs(f,{x,y},{xx(i-1)+h,yy(i-1)-h*k1+2*h*k2});yy(i)=yy(i-1)+h/6*(k1+4*k2+k3); end ff1=[xx,yy];

隱式龍格庫塔法

%隱式比顯式有更高的精度，但是因為求解非線性或者線性方程組導致計算量過大。 %下面編寫三階隱式 clear clc syms x y f=3*x^2; a=0; b=10; y0=0; h=0.1; n=(b-a)/h+1; xx=zeros(n,1); yy=zeros(n,1); for i=1:nxx(i)=a+h*(i-1); end yy=zeros(n,1); yy(1)=y0; c=[(5-15^0.5)/10,0.5,(5+15^0.5)/10]; b=[5/18,4/9,5/18]; aa=[5/36,(10-3*15^0.5)/45,(25-6*15^0.5)/180;(10+3*15^0.5)/72,2/9,(10-3*15^0.5)/72;(25+6*15^0.5)/180,(10+3*15^0.5)/45,5/36]; for i=2:nsyms k1 k2 k3k1=subs(f,{x,y},{xx(i-1)+c(1)*h,yy(i-1)+h*(aa(1,1)*k1+aa(1,2)*k2+aa(1,3)*k3)});k2=subs(f,{x,y},{xx(i-1)+c(2)*h,yy(i-1)+h*(aa(2,1)*k1+aa(2,2)*k2+aa(2,3)*k3)});k3=subs(f,{x,y},{xx(i-1)+c(3)*h,yy(i-1)+h*(aa(3,1)*k1+aa(3,2)*k2+aa(3,3)*k3)});yy(i)=yy(i-1)+h/6*(k1+4*k2+k3); end ff1=[xx,yy];

半隱式龍格庫塔法

%半隱式比隱式計算量更小，并且精度相當。 %下面編寫二階半隱式 clear clc syms x y f=6*x; A=diff(f,y,1); a=0; b=10; y0=0; h=0.1; n=(b-a)/h+1; xx=zeros(n,1); yy=zeros(n,1); for i=1:nxx(i)=a+h*(i-1); end yy=zeros(n,1); yy(1)=y0; c=[(5-15^0.5)/10,0.5,(5+15^0.5)/10]; b=[5/18,4/9,5/18]; aa=[5/36,(10-3*15^0.5)/45,(25-6*15^0.5)/180;(10+3*15^0.5)/72,2/9,(10-3*15^0.5)/72;(25+6*15^0.5)/180,(10+3*15^0.5)/45,5/36]; omig1=(6+6^0.5)/6; omig2=(6-6^0.5)/6; c2=0.17378667; b1=-0.41315432; b2=1-b1; for i=2:nk1=(1-omig1*subs(A,{x,y},{xx(i-1),yy(i-1)}))^(-1)*subs(f,{x,y},{xx(i-1),yy(i-1)});k2=(1-omig2*subs(A,{x,y},{xx(i-1)+c2*h,yy(i-1)+h*c2*k1}))^(-1)*subs(f,{x,y},{xx(i-1)+c2*h,yy(i-1)+h*c2*k1});yy(i)=yy(i-1)+h*(b1*k1+b2*k2); end ff1=[xx,yy];

總結

以上是生活随笔為你收集整理的微分方程数值解法的matlab程序的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。