Runge-Kutta 4th order method for two variable

12 views (last 30 days)
clc
clear
syms x1 x2
a=3;
b=5;
c=7;
d=2;
t(1)=0;
X1(1)=1;
X2(1)=-1;
TG=5;
n=100;
h=(TG-t(1))/n;
fx1=inline(x2);
fx2=inline(-a*x1+b*x1^2);
for i=1:n
K1_1=h*fx1(X2(i));
K2_1=h*fx1(X2(i)+K1_1/2);
K3_1=h*fx1(X2(i)+K2_1/2);
K4_1=h*fx1(X2(i)+K3_1);
X1(i+1)=X1(i)+(1/6)*(K1_1+2*K2_1+2*K3_1+K4_1);
K1=h*fx2(X1(i));
K2=h*fx2(X1(i)+K1_1/2);
K3=h*fx2(X1(i)+K2_1/2);
K4=h*fx2(X1(i)+K3_1);
X2(i+1)=X2(i)+(1/6)*(K1+2*K2+2*K3+K4);
t(i+1)=t(i)+h;
disp(' ')
disp(['t(',num2str(i),'): ',num2str(t(i+1))])
disp(['X1(',num2str(i),'): ',num2str(X1(i+1))])
disp(['X2(',num2str(i),'): ',num2str(X2(i+1))])
end
disp(' ')
subplot(1,3,1)
plot(t,'r--*')
title('t')
subplot(1,3,2)
plot(X1,'-.o')
title('X1')
subplot(1,3,3)
plot(X2,'-.*')
title('X2')

Answers (0)

Tags

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!