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

House

" {7 ~3 J& n( i  U
) ]9 N8 ~3 ^* t: }  {' O, uMATLAB源程序代码分享：MATLAB实现四阶龙格库塔法求解常微分方程组
4 F. o. J% q: I: [$ T2 y) {" jfunction F=f(t,Y)
% 定义待求的微分方程组
5 D0 K; Y$ V1 G0 |/ cx=Y(1);
y=Y(2);
f1=-x^2*(5-y);
f2=y*(4-3*x);
# T' @2 C- Y# r" b( E+ b* \F=[f1;f2];
end

% J6 u- r1 j% b. i5 t' r; }0 i) Z
8 n" E" h$ _% E1 X
0 W! p) Q% C% ?2 G2 l
* w8 s" x  l% u# W" b5 y%% 定义计算的步长, 设置变量的初始值
clear;clc;close all
* r( ^3 {! W4 {; ^+ g0 |8 d, F3 JDelta=0.001;     % 步长
t=0: Delta :8;     % 定义自变量 t
n=length(t);
7 Y& ~; O. i; Z" oY(:,1)=[0.5;3];  % 变量 x y 的初始值
( Y. Q- a8 q" G$ v& @( F
%% 自定义龙格库塔法, 求解微分方程组
for k=1: n-1
    z1=f(t(k),Y(:,k));
! b5 E. G4 s" s/ f    z2=f(t(k)+Delta/2,Y(:,k)+z1*Delta/2);
    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;
( u6 A2 u& i) R* d' Q0 ?$ q: tend
x=Y(1,: );
y=Y(2,: );

: W7 P# n1 k0 P8 y%% 绘制 x y 的求解结果
figure
set(gcf,'units','normalized','position',[0.15 0.2 0.7 0.6]);  % 设置 figure 窗口的位置和尺寸
subplot(1,2,1)
4 ~* ~! B5 ], s  J. Wplot(t,x,'b')
xlabel('t')
ylabel('x')

subplot(1,2,2)
) ?% X) k3 [( T; |* ]; ]plot(t,y,'r')
/ R( R7 v# e/ Txlabel('t')
5 W+ `0 y( L% lylabel('y')

$ U6 L2 S" U5 Q/ `

helendcany

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