Please explain this code about differential equations.

Here is the matlab code

theta=0.10895;
YF=0.0667;
alpha= 0.29;
beta= 0.68;
gamma1=450;
gamma2=11.25;
X0=0; Y0=0.0667; Z0=0;
f=@(t,y)[-y(1)/theta+(1+alpha)*gamma1*(1-y(1))*y(2)^2+beta*gamma1*(1-y(1))*y(3)^2;...
    (YF-y(2))/theta+(1-alpha)*gamma1*(1-y(1))*y(2)^2-gamma2*y(2);...
    -y(3)/theta+beta*gamma1*(1-y(1))*y(3)^2+2*alpha*gamma1*(1-y(1))*y(2)^2-gamma2*y(3)/beta];
[T,Y]=ode45(f,[100 120],[X0,Y0,Z0]);
plot(T,Y(:,1),'-',T,Y(:,3),'-.',T,Y(:,3),'.');

Can you please tell me where the underlined numbers (below) come from in f=@(t,y) ??????????

-y(1)

(1-y(1))*y(2)

(1-y(1))*y(3)^2

(YF-y(2))

(1-y(1))*y(2)^2-gamma2*y(2)

-y(3)/theta+beta

(1-y(1))*y(3)^2+2

y(2)^2-gamma2*y(3)/beta

Answer 1

[Star Strider]

They refer to the function values returned by ode45 (in this code).

y(1) is X

y(2) is Y

y(3) is Z

That can easily be inferred by comparing the code to the symbolic differential equation system.

Answer 2

[Sam Chak]

@Mr.DDWW

The code is annotated now. Hopefully sufficient for you to understand. By the way, I've fixed the plot line because you plotted Z(t) twice.

% Parameters
theta  = 0.10895;
YF     = 0.0667;
alpha  = 0.29;
beta   = 0.68;
gamma1 = 450;
gamma2 = 11.25;

% Initial values
X0 = 0; 
Y0 = 0.0667; 
Z0 = 0;

% A system of differential equations
f = @(t,y)[-y(1)/theta+(1+alpha)*gamma1*(1-y(1))*y(2)^2+beta*gamma1*(1-y(1))*y(3)^2;...              % this is dX/dt
           (YF-y(2))/theta+(1-alpha)*gamma1*(1-y(1))*y(2)^2-gamma2*y(2);...                          % this is dY/dt
           -y(3)/theta+beta*gamma1*(1-y(1))*y(3)^2+2*alpha*gamma1*(1-y(1))*y(2)^2-gamma2*y(3)/beta]; % this is dZ/dt

% Solving the system f from time 100 s to 120 s with the initial values using the ode45 solver
[T, Y] = ode45(f, [100 120], [X0, Y0, Z0]);

% plotting the solutions for X(t), Y(t), Z(t)
plot(T, Y(:,1), '-', T, Y(:,2), '-.', T, Y(:,3), '.');

% additional stuff to display label, title, and legends
xlabel('Time, t [sec]')
title('Time responses of the System')
legend({'$X(t)$', '$Y(t)$', '$Z(t)$'}, 'Interpreter', 'latex', 'location', 'best')

Result: