Solve system of differential equations using Euler forward and ode45

I have a basic SIR model which is described by three differential equations:

and I want to solve these using Euler forward and ode45. I have never worked with these types of problems before but from research I have found that for Euler forward I should use the equations:

But all the examples I have seen using this are given initial values for the different parameters, which I haven't. So I don't even know where to start. I would appreciate if someone could push me in the right direction of how to solve this. Thank you!

Answer 1

[Ameer Hamza]

When using Euler forward, you start with

S(1) = s0; % initial value of S
I(1) = i0; % initial value of I
R(1) = r0; % initial value of R

and then run for loop using the equations in your question

N = % total number of time steps
dT = % length of each time step

for i = 2:N
    S(i) = S(i-1) + % equation for Sn in question
    I(i) = I(i-1) + % equation for In in question
    R(i) = R(i-1) + % equation for Rn in question
end

t = 0:dT:(N-1)*dT;

plot(t, S, t, I, t, R)

See these example of how Euler method can be written in MATLAB

https://www.mathworks.com/matlabcentral/fileexchange/45664-euler-method

https://www.mathworks.com/matlabcentral/fileexchange/72522-euler-method

[Katara So]

So should the for loop equations should be as below?

for i = 2:N
    S(i) = S(i-1) + S(N)*(1-dT(beta*I(N))/N)
    I(i) = I(i-1) + I(N)*(1-dT(beta*S(N)-(gamma*N))/N)
    R(i) = R(i-1) + R(N)+dT*gamma*I(N)
end

But then again I don't have any values for the parameters so how will I be able to solve these?

[Ameer Hamza]

No, the question mentions two relations for Sn, In, and Rn. If you want to use the second one, then write it like this

for i = 2:N
    S(i) = S(i-1)*(1-dT(beta*I(i-1))/N)
    I(i) = I(i-1)*(1-dT(beta*S(i-1)-(gamma*N))/N)
    R(i) = R(i-1)+dT*gamma*I(i-1)
end

You can start with any arbitrary value of parameters.

[Katara So]

Thank you for taking the time to help me.

I have another question. Say I was given:

tMax = 20;
timeSpan = [0 tMax];

S0 = 1.0;
I0 = 0.0;
R0 = 0.0;
dt = 0.01;

time_vector = 0:dt:tMax;
nIterations = length(time_vector);

tau = 0.6;
h = 0.5;
rho = 0.8;
r = 0.2;

And asked to calculate the differential equations I posted in my original question using Euler forward and ode45. Here the equations that I found for Sn+1 etc. are applicable since I don't have a parameter gamma. So I approached the problem in another way however I don't seem to get the correct plot. I would truly appreciate it if you could take a look at my code and see if what I have done is correct or perhaps spot where I have gone wrong. Thank you!

tMax = 20;
timeSpan = [0 tMax];
dt = 0.01;
tau = 0.6;
h = 0.5;
rho = 0.8;
r = 0.2;
beta = (h*exp(-h*tau))/(1-exp(-h*tau));

time_vector = 0:dt:tMax;
nIters = length(time_vector);
t = 0:dt:tMax;

% Initial conditions
S(1) = 1.0;
I(1) = 0.0;
R(1) = 0.0;

n=0;
for i = time_vector
    n=n+1;
    S(n+1) = S(n)-dt*h*S(n)+dt*rho*I(n)+dt*beta*R(n);
    I(n+1) = I(n)+dt*h*S(n)-dt*rho*I(n)-dt*r*I(n);
    R(n+1) = R(n)+dt*r*I(n)-dt*beta*R(n);
end 
hold on
plot (t,S(1:end-1),'-r');
plot (t,I(1:end-1),'-b');
plot (t,R(1:end-1),'-g');
 
 
xlabel ('t');
ylabel ('population');
title ('Spread of Malaria');
legend ('S','I','R');

and the plot is

[Ameer Hamza]

Can you explain what is wrong with this plot? This seems to follow the expected trajectories of the SIR model.