System of 4 ODEs solver

14 views (last 30 days)
Duong Quang Chi
Duong Quang Chi on 6 May 2024
Commented: Duong Quang Chi on 6 May 2024
I'm solving a system of 4 ODEs but it wasn't working. It showed "Unable to find explicit solution". What should I do next?
Here is the ODE system:
syms V(t)
syms X(t)
syms S(t)
syms P(t)
Q = 30;
nuy_max = 0.2;
Ks = 0.5;
nuy = nuy_max * S / (Ks + S);
nuy_g = 0.05;
K = 0;
Y = 0.3;
S_F = 300;
ode1 = diff(V) == Q;
ode2 = diff(X) == (nuy - Q/V)*X;
ode3 = diff(P*V) == nuy_g * X * V - K*P*V;
ode4 = diff(S*V) == S_F * Q - nuy * X * V / Y;
odes = [ode1; ode2; ode3; ode4];
cond1= V(0) == 200;
cond2 = S(0) == 1;
cond3 = X(0) == 20;
cond4 = P(0) == 0.1;
conds = [cond1; cond2 ; cond3; cond4];
Str = dsolve(odes, conds);
V.Sol(t) = Str.V
X.Sol(t) = Str.X
S.Sol(t) = Str.S
P.Sol(t) = Str.P
The result is
Please help me. Thanks in advance.

Accepted Answer

Sam Chak
Sam Chak on 6 May 2024
Edited: Sam Chak on 6 May 2024
If 'dsolve' cannot find explicit solutions for , , , and , we can use 'ode45' to numerically solve the ODEs and plot the results.
syms X(t) V(t) P(t) S(t)
Q = 30;
nuy_max = 0.2;
Ks = 0.5;
nuy = nuy_max*S/(Ks + S);
nuy_g = 0.05;
K = 0;
Y = 0.3;
S_F = 300;
eqn1 = diff(V) == Q;
eqn2 = diff(X) == (nuy - Q/V)*X;
eqn3 = diff(P*V) == nuy_g*X*V - K*P*V;
eqn4 = diff(S*V) == S_F*Q - nuy*X*V/Y;
eqns = [eqn1 eqn2 eqn3 eqn4];
[V, S] = odeToVectorField(eqns)
V = 
S = 
odefun = matlabFunction(V, 'vars', {'t','Y'});
tspan = [0 20];
X0 = 20;
V0 = 200;
P0 = 0.1;
S0 = 1;
y0 = [X0; V0; P0; S0];
[t, y] = ode45(odefun, tspan, y0);
tb = tiledlayout(2, 2, 'TileSpacing', 'Compact');
nexttile
plot(t, y(:,1)), grid on, title('X(t)')
nexttile
plot(t, y(:,2)), grid on, title('V(t)')
nexttile
plot(t, y(:,3)), grid on, title('P(t)')
nexttile
plot(t, y(:,4)), grid on, title('S(t)')
xlabel(tb, 'Time')
  1 Comment
Duong Quang Chi
Duong Quang Chi on 6 May 2024
Omg it works :)) Though the explicit solutions cannot be found, these demonstration might be sufficient for my projects. Thanks a lot, I thought I was doomed.

Sign in to comment.

More Answers (0)

Products


Release

R2019a

Community Treasure Hunt

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

Start Hunting!