Time-dependent immunity and drug resistant SIR model
6 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
I have to solve this SIR model with Euler's forward. But my teacher said that it is wrong with my graph. Can anyone help me ?
Here is my code
%initial values
S(1)= 1;
I_1(1)= 0;
I_2 (1) = 0;
R(1)= 0;
tMax= 50;
dt = 0.01;
time_vector = 0:dt:tMax;
nIterations = length(time_vector);
t= 0:dt:tMax;
v= 0.0024;
tau = 0.6;
u= 0.8;
h= 0.5;
rho = 0.8;
R = 0.2;
r= re(R,v,t);
b = BetaL(h,tau);
n=0 ;
for i = time_vector
n = n+1;
S(n+1)= S(n)-dt*h*S(n)+dt*rho*I_1(n)+dt*b.*R(n);
I_1(n+1) =I_1(n)+ dt*(1-u)*S(n) - dt*rho*I_1(n) - dt*u*h*I_1(n);
I_2(n+1) =I_2(n)+ dt*u*h*S(n) + dt*u*h*I_1(n) - dt*r(n)*I_2(n);
R(n+1)= R(n) + dt*r(n)*I_2(n) - dt*b*R(n);
end
hold on
plot (t,S(1:end-1),'-r');
plot (t,I_1(1:end-1),'-b');
plot (t,I_2(1:end-1),'.b');
plot (t,R(1:end-1),'-g');
xlabel ('t');
ylabel ('population');
title ('Drug resistant model');
legend ('S','I_1','I_2','R');
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/245222/image.png)
1 Commento
Aquatris
il 28 Ott 2019
What is BetaL function?
In the line with "S(n+1) = ..." you do not use the b(h)
In the line "I2(n+1) = ..." you are missing (1-exp(-vt^2).
In the line "R(n+1) = ..." you are missing (1-exp(-vt^2)
Risposta accettata
ME
il 28 Ott 2019
Modificato: ME
il 28 Ott 2019
Assuming you have your BetaL and re functions setup correctly then there are the errors listed in the comment from Aquatris but you are also missing a multiple "h" in the I_1(n+1) equation.
Your equations at the bottom of your question also state that
=0.2 but you appear to have set R=0.2 instead and then defined a function r(n) instead of using
.
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/245236/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/245236/image.png)
If fixing those items doesn't help then I think you might need to give some more information as to what the BetaL and re functions are supposed to be.
Più risposte (0)
Vedere anche
Categorie
Scopri di più su 2-D and 3-D Plots in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!