MATLAB中已知系统微分方程及初始值用欧拉法和龙格库塔法解一阶微分方程

function Euler

%欧拉法和龙格库塔算法解一阶常微分方程源代码

%例子dy/dx=-y+x+1

f=inline('-y+x+1','x','y'); %微分方程的右边项

dx=05; %x方向步长

xleft=0; %区域的左边界

xright=10; %区域的右边界

xx=xleft:dx:xright; %一系列离散的点

n=length(xx); %点的个数

y0=1;

%%(1)欧拉法

Euler=y0;

for i=2:n

Euler(i)=Euler(i-1)+dxf(xx(i-1),Euler(i-1));

end

%%(2)龙格库塔法

RK=y0;

for i=2:n

k1=f(xx(i-1),RK(i-1));

k2=f(xx(i-1)+dx/2,RK(i-1)+k1dx/2);

k3=f(xx(i-1)+dx/2,RK(i-1)+k2dx/2);

k4=f(xx(i-1)+dx,RK(i-1)+k3dx);

RK(i)=RK(i-1)+dx(k1+2k2+2k3+k4)/6;

end

%%Euler和Rk法结果比较

plot(xx,Euler,xx,RK)

hold on

%精确解用作图

syms x

rightsolve=dsolve('Dy=-y+x+1','y(0)=1','x');%求出解析解

rightdata=subs(rightsolve,xx);%将xx代入解析解,得到解析解对应的数值

plot(xx,rightdata,'r')

legend('Euler','Runge-Kutta','analytic')

function dy=test(t,y)

dy=[-12cos(y(2))-120cos(2082pi/360-3y(2));(12sin(y(2))+120sin(2082pi/360-3y(2)))/y(1);];

[t,y]=ode45('test',[001,1],[1,1])

plot(t,y(:,1),t,y(:,2));

x=3000sin(702pi/360)-y(1)sin(y(2));

z=3000cos(702pi/360)-12t-y(1)cos(y(2));

plot(t,x,t,z)

最后的这图像不好。你可再看看。方法是没问题的。

你可到MATLAB中文论坛去请教关于MATLAB的问题，那儿高手多的很。

function [Y] = RK45(t,X,f,h)

K1=f(t,X);

K2=f(t+h/2,X+h/2K1);

K3=f(t+h/2,X+h/2K2);

K4=f(t+h,X+hK3);

Y=X+h/6(K1+2K2+2K3+K4);

end

以上是4阶龙格库塔法的代码:

自己写函数，存为fm

function dxdt = f (t,x)

dxdt(1)=exp(x(1)sin(t))+x(2);

dxdt(2)=exp(x(2)cos(t))+x(1); % x(1)是你的f，x(2)是你的g

dxdt=dxdt(:);

end

自己给出t0,x0,h的值（初始时间，初值，步长）

如果求t0到t1的轨迹的话:给个例子如下

t0=0;t1=5;h=002;x0=[-1;-1];

T=t0:h:t1;X=zeros(length(x0),length(T));X(:,1)=x0;

for j=1:length(T)-1

X(:,j+1)=RK45(T(j),X(:,j),@(t,x) f(t,x),h);

end

plot(T,X(1,:));

hold on;

plot(T,X(2,:),'r');

具体参数自己设置

用m4rkutta函数求解，见附件

输入以下程序

df=@(x,y)3y/(1+x);

[x,y]=m4rkutta(df,[0 1],1,02)

得到的结果为：

x =

0    02000    04000    06000    08000    10000

y =

10000    17275    27430    40942    58292    79960

再matlab命令窗口输入 doc ode45 可以查看龙格库塔算法的详细解释和用法

[T,Y] = ode45(@vdp1000,[0 3000],[2 0]); 这是龙格库塔4阶算法的示例，[2,0]为初值。

以上就是关于MATLAB中已知系统微分方程及初始值用欧拉法和龙格库塔法解一阶微分方程全部的内容，包括:MATLAB中已知系统微分方程及初始值用欧拉法和龙格库塔法解一阶微分方程、matlab用四阶龙格库塔法解微分方程组；、matlab龙格库塔法解微分方程等相关内容解答，如果想了解更多相关内容，可以关注我们，你们的支持是我们更新的动力！