Matlab numerical simulation for ode system by changing the parameter value in range graph not running error coming
Afficher commentaires plus anciens
function su
options = odeset('RelTol',1e-6,'Stats','on');
%initial conditions
Xo = [0.005; 0.0007; 0.001; 0.0001; 0.001; 0; 0; 0.007];
tspan = [0,120];
tic
[t,X] = ode45(@TestFunction,tspan,Xo,options);
toc
%figure
%plot(t, X(:,1), 'b')
plot(t, X(:,2), 'b')
%plot(t, X(:,3), 'g')
%plot(t, X(:,4), 'm')
%plot(t, X(:,5), 'k')
%plot(t, X(:,6), 'c')
%plot(t, X(:,7), 'y')
%plot(t, X(:,8), 'y')
hold on
% legend('x1','x2')
% ylabel('x - Population')
% xlabel('t - Time')
hold on
return
function [dx_dt]= TestFunction(~,x)
% Parameters
s = 0.038;
alpha = 0.02;
gamma = 0.10;
r = 0.03;
dH = 0.0083;
dX = 0.0125;
rho = 0.07;
K = 500;
beta = 0.0005;
theta = 0.03;
eta = 0.015;
muH = 0.015;
muI = 0.08;
UP = 0.20;
deltaP = 0.033;
UL = 0.50;
deltaL = 0.05;
wP = 10000;
wL = 20000;
wV = 1;
LC = 0.05;
dx_dt(1) = s - (alpha * x(3) * x(1)) / (1 + gamma * x(1)) - dH * x(1) + r * x(2);
dx_dt(2) = (alpha * x(3) * x(1)) / (1 + gamma * x(1)) - (r+dX) * x(2);
dx_dt(3) = rho * x(3) * (1 - (x(3) + x(4)) / K) - beta * x(3) * x(5) - muH * x(3) - LC * x(8) * x(3);
dx_dt(4) = beta * x(3) * x(5) - muI * x(4) - LC * x(8) * x(4);
dx_dt(5) = theta * x(4) - eta * x(5) - beta * x(3) * x(5);
dx_dt(6) = UP - deltaP * x(6);
dx_dt(7) = UL - deltaL * x(7);
dx_dt(8) = wP * x(6) + wL * x(7) + wV * x(5);
% dx_dt(1) = (s) - ((alpha1 * (x(3)+x(4)) * x(1)) / (1 + gamma1 * x(1))) - ((alpha2 * x(5) * x(1)) / (1 + gamma2 * x(1)));
%dx_dt(2) = ((alpha1 * (x(3)+x(4)) * x(1)) / (1 + gamma1 * x(1)) + (alpha2 * x(5) * x(1)) / (1 + gamma2 * x(1))) - ((k1 * x(3) + k2 * x(5)) * x(2)) - (r * x(2));
%dx_dt(3) = (gamma1 * x(3) * (1 - ((x(3)+x(4)) / kR))) - (muR * x(3)) - (beta * x(3) * x(4));
%dx_dt(4) = (beta * x(3) * x(4)) - ((muR + gamma * x(6)) * x(4));
%dx_dt(5) = (gamma2 * x(5) * (1 - x(5) / kW)) - (muW * x(5)) - (theta * (x(3)+x(4)) * x(5));
% dx_dt(6) = (delta_v) + (rho * x(4)) - (eta * x(6));
dx_dt = dx_dt';
return
hold on In this code error coming. This is the error.
Operation terminated by user during ode45
In summ (line 7)
[t,X] = ode45(@TestFunction,tspan,Xo,options);
Réponse acceptée
Your ODE system is stiff. Use "ode15s" instead of "ode45".
The reason for the error you receive seems to be that you interrupted the computation (maybe because it took too long).
su()
function su
options = odeset('RelTol',1e-6,'Stats','on');
%initial conditions
Xo = [0.005; 0.0007; 0.001; 0.0001; 0.001; 0; 0; 0.007];
tspan = [0,120];
tic
[t,X] = ode15s(@TestFunction,tspan,Xo,options);
toc
%figure
%plot(t, X(:,1), 'b')
plot(t, X(:,2), 'b')
%plot(t, X(:,3), 'g')
%plot(t, X(:,4), 'm')
%plot(t, X(:,5), 'k')
%plot(t, X(:,6), 'c')
%plot(t, X(:,7), 'y')
%plot(t, X(:,8), 'y')
%hold on
% legend('x1','x2')
% ylabel('x - Population')
% xlabel('t - Time')
%hold on
end
function [dx_dt]= TestFunction(~,x)
% Parameters
s = 0.038;
alpha = 0.02;
gamma = 0.10;
r = 0.03;
dH = 0.0083;
dX = 0.0125;
rho = 0.07;
K = 500;
beta = 0.0005;
theta = 0.03;
eta = 0.015;
muH = 0.015;
muI = 0.08;
UP = 0.20;
deltaP = 0.033;
UL = 0.50;
deltaL = 0.05;
wP = 10000;
wL = 20000;
wV = 1;
LC = 0.05;
dx_dt(1) = s - (alpha * x(3) * x(1)) / (1 + gamma * x(1)) - dH * x(1) + r * x(2);
dx_dt(2) = (alpha * x(3) * x(1)) / (1 + gamma * x(1)) - (r+dX) * x(2);
dx_dt(3) = rho * x(3) * (1 - (x(3) + x(4)) / K) - beta * x(3) * x(5) - muH * x(3) - LC * x(8) * x(3);
dx_dt(4) = beta * x(3) * x(5) - muI * x(4) - LC * x(8) * x(4);
dx_dt(5) = theta * x(4) - eta * x(5) - beta * x(3) * x(5);
dx_dt(6) = UP - deltaP * x(6);
dx_dt(7) = UL - deltaL * x(7);
dx_dt(8) = wP * x(6) + wL * x(7) + wV * x(5);
% dx_dt(1) = (s) - ((alpha1 * (x(3)+x(4)) * x(1)) / (1 + gamma1 * x(1))) - ((alpha2 * x(5) * x(1)) / (1 + gamma2 * x(1)));
%dx_dt(2) = ((alpha1 * (x(3)+x(4)) * x(1)) / (1 + gamma1 * x(1)) + (alpha2 * x(5) * x(1)) / (1 + gamma2 * x(1))) - ((k1 * x(3) + k2 * x(5)) * x(2)) - (r * x(2));
%dx_dt(3) = (gamma1 * x(3) * (1 - ((x(3)+x(4)) / kR))) - (muR * x(3)) - (beta * x(3) * x(4));
%dx_dt(4) = (beta * x(3) * x(4)) - ((muR + gamma * x(6)) * x(4));
%dx_dt(5) = (gamma2 * x(5) * (1 - x(5) / kW)) - (muW * x(5)) - (theta * (x(3)+x(4)) * x(5));
% dx_dt(6) = (delta_v) + (rho * x(4)) - (eta * x(6));
dx_dt = dx_dt';
end
1 commentaire
For your records, state no.8 grows rapidly to become very large (at the magnitude of 
), which may contribute to the system becoming increasingly stiff over time.
), which may contribute to the system becoming increasingly stiff over time. options = odeset('RelTol', 1e-6);
Xo = [0.005; 0.0007; 0.001; 0.0001; 0.001; 0; 0; 0.007];
tspan = [0, 120];
[t, X] = ode15s(@TestFunction, tspan, Xo, options);
figure
plot(t, X(:,8), 'b')
ylabel('x - Population')
xlabel('t - Time')
dx_dt = TestFunction(t(end), X(end,:)')
function [dx_dt]= TestFunction(t, x)
% Parameters
s = 0.038;
alpha = 0.02;
gamma = 0.10;
r = 0.03;
dH = 0.0083;
dX = 0.0125;
rho = 0.07;
K = 500;
beta = 0.0005;
theta = 0.03;
eta = 0.015;
muH = 0.015;
muI = 0.08;
UP = 0.20;
deltaP = 0.033;
UL = 0.50;
deltaL = 0.05;
wP = 10000;
wL = 20000;
wV = 1;
LC = 0.05;
dx_dt(1)= s - (alpha * x(3) * x(1)) / (1 + gamma * x(1)) - dH * x(1) + r * x(2);
dx_dt(2)= (alpha * x(3) * x(1)) / (1 + gamma * x(1)) - (r+dX) * x(2);
dx_dt(3)= rho * x(3) * (1 - (x(3) + x(4)) / K) - beta * x(3) * x(5) - muH * x(3) - LC * x(8) * x(3);
dx_dt(4)= beta * x(3) * x(5) - muI * x(4) - LC * x(8) * x(4);
dx_dt(5)= theta * x(4) - eta * x(5) - beta * x(3) * x(5);
dx_dt(6)= UP - deltaP * x(6);
dx_dt(7)= UL - deltaL * x(7);
dx_dt(8)= wP * x(6) + wL * x(7) + wV * x(5);
dx_dt = dx_dt';
end
Plus de réponses (0)
Catégories
En savoir plus sur Ordinary Differential Equations dans Centre d'aide et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
