Why do the ode45 solver is not running in my code?

1 view (last 30 days)
clear; clc;
a = 1; %alpha
p = 7.5; %rho
b = 4; %beta
g = 16; %gamma
u = 1/g; %mue
q = 1; %taking q(m) = 1 or m* = 1
k = 1500;
s = .5; %sigma
eta = 1;
c1 = 150; %C*
% V(1) = V, V(2) = F, V(3) = C,V(4) = M
f = @(t,V) [a*V(1)-p*V(1)*V(2);b*V(3)-g*p*V(1)*V(2)-a*V(2);-u*(V(3)-c1)+q*k*V(1)*V(2);s*V(1)-eta*V(4)];
Y0 = [10000;9850;c1;0];
[time,sol] = ode45(f,[0 8],Y0);
plot(time,sol(:,1),'b')
hold on
plot(time,sol(:,2),'r')
hold on
plot(time,sol(:,3),'g')
hold on
plot(time,sol(:,4),'m')

Accepted Answer

Stephan
Stephan on 14 May 2021
Your problem appears to be stiff - use ode15s instead:
clear; clc;
a = 1; %alpha
p = 7.5; %rho
b = 4; %beta
g = 16; %gamma
u = 1/g; %mue
q = 1; %taking q(m) = 1 or m* = 1
k = 1500;
s = .5; %sigma
eta = 1;
c1 = 150; %C*
% V(1) = V, V(2) = F, V(3) = C,V(4) = M
f = @(t,V) [a*V(1)-p*V(1)*V(2);b*V(3)-g*p*V(1)*V(2)-a*V(2);-u*(V(3)-c1)+q*k*V(1)*V(2);s*V(1)-eta*V(4)];
Y0 = [10000;9850;c1;0];
[time,sol] = ode15s(f,[0 8],Y0);
plot(time,sol(:,1),'b')
hold on
plot(time,sol(:,2),'r')
hold on
plot(time,sol(:,3),'g')
hold on
plot(time,sol(:,4),'m')
  1 Comment
Aman Gupta
Aman Gupta on 14 May 2021
Edited: Aman Gupta on 14 May 2021
Thanks.
Any idea why the first plot is not getting plotted.

Sign in to comment.

More Answers (0)

Products

Community Treasure Hunt

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

Start Hunting!