MATLAB源程序代码分享：MATLAB实现四阶龙格库塔法求解常微分方程组

House

  K( k9 e, L. E! J& H! y
MATLAB源程序代码分享：MATLAB实现四阶龙格库塔法求解常微分方程组
0 h5 b3 i% P3 p. k- R& \function F=f(t,Y)
, j/ e$ z( W- N3 D0 n% 定义待求的微分方程组
3 H& v" m3 z, l( k1 y8 V2 j$ A5 ~x=Y(1);
4 e& H5 S$ Z! {* Y8 `/ R2 S5 iy=Y(2);
f1=-x^2*(5-y);
f2=y*(4-3*x);
F=[f1;f2];
end
9 m1 Y% F) \( E% q7 Q

3 ^  R) k# K  w" ]" J  |2 }
$ @* j' q9 F: m4 ]$ Q: X
%% 定义计算的步长, 设置变量的初始值
, b8 Y7 H) d  J" W! Dclear;clc;close all
Delta=0.001;     % 步长
, o! Z9 x/ B2 @* v$ c8 q- h  Pt=0: Delta :8;     % 定义自变量 t
n=length(t);
2 V6 S; D- j" e  J# gY(:,1)=[0.5;3];  % 变量 x y 的初始值
" u! d* y0 S6 J7 Q" o( z0 k
. v/ e" d9 e3 A& P9 r! b%% 自定义龙格库塔法, 求解微分方程组
6 D2 U$ i8 x: c) B4 E/ ]1 o5 dfor k=1: n-1
    z1=f(t(k),Y(:,k));
    z2=f(t(k)+Delta/2,Y(:,k)+z1*Delta/2);
6 v( p5 p  z' Y0 m4 ?' L    z3=f(t(k)+Delta/2,Y(:,k)+z2*Delta/2);
    z4=f(t(k)+Delta,Y(:,k)+z3*Delta);
    Y(:,k+1)=Y(:,k)+Delta*(z1+2*z2+2*z3+z4)/6;
end
% I) r" y8 d4 e( E9 ~x=Y(1,: );
y=Y(2,: );
( B# W5 F1 [+ U+ B# i
+ D( u2 U# l, q# t7 |: D%% 绘制 x y 的求解结果
- g" L5 k3 U: v/ r) h$ D3 \" kfigure
0 A, O3 h  L. V& o% wset(gcf,'units','normalized','position',[0.15 0.2 0.7 0.6]);  % 设置 figure 窗口的位置和尺寸
  I0 o# N* b7 ^  M- @0 b- Osubplot(1,2,1)
plot(t,x,'b')
xlabel('t')
ylabel('x')

subplot(1,2,2)
plot(t,y,'r')
- N! M* `# O# Q8 }" d3 x: l# \xlabel('t')
ylabel('y')
+ @) |, Z" ]# h% F& K
- u  U7 E4 w" T

helendcany

MATLAB实现四阶龙格库塔法求解常微分方程组