Solve the following simultaneous equations:

I am trying to make a code and plot and Draw the diagram of X and Y vs. t for 0<t<5

here is my code. I think there is somthoing worng

theta=0.10895;
YF=0.0667;
X=0;
alpha= 0.29;
beta= 0.68;
y1=450;
y2=11.25;
t=0; Y=0.0667; Z=0;
f1=@(t,y)[-X/theta+(1+alpha)*y1(1-X)*Y^2+beta*y1(1-X)*Z^2;(YF-Y/theta)+(1+alpha)*y1(1-X)*Y^2-y2*Y;-Z/theta+beta*y1(1-X)*Z^2+2*alpha*y1(1-X)*Y^2-((y2*Z/beta));];

[T,u]=ode45(f1,[100 120],[0 0 0.0667]);
plot(T,u(:,1),'-',T,u(:,3),'-.',T,u(:,3),'.');

[Jan]

Please mention why you think, that there is a problem. This is a very useful information.

If you get error message, post them here.

[Mr.DDWW]

% I tried to run the code but it is not working. here is the code
clc;clf;clear all;close all; 
syms X(t) Y(t) Z(t) YF alpha beta gamma1  gamma2 T theta

% The given data
YF=0.0667;
alpha= 0.29;
beta= 0.68;
gamma1=450;
gamma2=11.25;
theta=0.10895;


Ode1 = diff(X) == -X/theta + (1+alpha)*gamma1*(1-X)*Y^2 + beta*gamma1*(1-X)*Z^2;
Ode2 = diff(Y) == (YF-Y)/theta+(1-alpha)*gamma1*(1-X)*Y^2-gamma2*Y;
Ode3 = diff(Z) == -Z/theta+beta*gamma1*(1-X)*Z^2+2*alpha*gamma1*(1-X)*Y^2-(gamma2*Z)/beta;

Odes = [Ode1;Ode2; Ode3];
[VF, Subs] = OdeToVectorField ([Ode1, Ode2,Ode3]);

M = matlabFunction (VF, 'Vars',{'T','Y'});

% The  given condtions
[t,y]=ode45(M, [0 5],[0.0667 0 0]);
figure
plot(t,y(:,1:2))
gride on
title('X and Y')
legend (string(subs))

[t,y]=ode45(M, [100 120],[0.0667 0 0]);
figure
plot(t,y(:,1:2))
gride on
title('X and Y')
legend (string(subs))

[t,y]=ode45(M, [100 120],[0.0667 0 0]);
figure
plot(t,y(:,1),t,y(:,2))
gride on
title('X vs Y')
legend (string(subs))

Answer 1

[Jan]

f1 = @(t,y) [-X/theta+(1+alpha)*y1(1-X)*Y^2+beta*y1(1-X)*Z^2; ...
    ... %                         ^                ^   * are missing
             (YF-Y/theta)+(1+alpha)*y1(1-X)*Y^2-y2*Y; ...
    ... %    ^^^^^^^^^^^^             ^  and a missing * again
             -Z/theta+beta*y1(1-X)*Z^2+2*alpha*y1(1-X)*Y^2-((y2*Z/beta));];
    ... %                    ^                   ^

y1(1-X) would be indexing in the vector y1. You mean y1 * (1 - X) at all 5 locations.

At the 2nd marked location the closing parenthesis is at the wrong position:

(YF - Y) / theta

Using some spaces would increase the readability and support the debugging:

f1 = @(t,y) ...
    [-X / theta + (1+alpha) * y1 * (1-X) * Y^2 + beta * y1 * (1-X) * Z^2; ...
     (YF-Y) / theta + (1+alpha) * y1 * (1-X) * Y^2 - y2 * Y; ...
     -Z / theta + beta * y1 * (1-X) * Z^2 + 2 * alpha * y1 * (1-X) * Y^2 - y2 * Z / beta];

Answer 2

[Torsten]

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),'.');

[Mr.DDWW]

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];

can you please tell me how did you get the underline numbers ??????????

[Jan]

In the image of the original formula the vriables X, Y, Z are used. In the function f=@(t,y) the variables are provided as a vector:

y = [X, Y, Z]

Then X is y(1), Y is y(2) and Z is y(3).

[Torsten]

@Mr.DDWW

The MATLAB solvers expect the unknowns be defined as a vector (here: y), not as a collection of scalar variables (here: X, Y and Z). So one has to decide which of your solution variables (X, Y and Z) to placed at which position of this vector. I decided to take X = y(1), Y = y(2), Z = y(3).