How to solve SIR model with using DTM (Differential Transform Method)

11 次查看（过去 30 天）
mohammad amin akhlaghi
回答： Athanasios Paraskevopoulos 2024-5-16
I am trying to use dtm for solving SIR model, Although my code is run but I think the DTM part is wrong. I need help for DTM part
here is my code
function sir_model_dtm()
% Parameters
alpha = 0.3; % Infection rate
beta = 0.1; % Recovery rate
N = 1000; % Total population
I0 = 1; % Initial number of infected individuals
R0 = 0; % Initial number of recovered individuals
S0 = N - I0 - R0; % Initial number of susceptible individuals
% Time parameters
T = 100; % Total simulation time
dt = 1; % Time step
% Initialize arrays to store results
S = zeros(1, T);
I = zeros(1, T);
R = zeros(1, T);
% Initial conditions
S(1) = S0;
I(1) = I0;
R(1) = R0;
% Simulate the SIR model using DTM
for t = 2:T
% Compute new values
S(t) = S(t-1) - alpha*S(t-1)*I(t-1)/N * dt;
I(t) = I(t-1) + (alpha*S(t-1)*I(t-1)/N - beta*I(t-1)) * dt;
R(t) = R(t-1) + beta*I(t-1) * dt;
end
% Plot the results
t = 0:dt:T-dt;
plot(t, S, 'b', t, I, 'r', t, R, 'g', 'LineWidth', 2);
legend('Susceptible', 'Infectious', 'Recovered');
xlabel('Time');
ylabel('Number of individuals');
title('SIR Model');
end

请先登录，再回答此问题。

回答（1 个）

Athanasios Paraskevopoulos
function sir_model_dtm()
% Parameters
alpha = 0.3; % Infection rate
beta = 0.1; % Recovery rate
N = 1000; % Total population
I0 = 1; % Initial number of infected individuals
R0 = 0; % Initial number of recovered individuals
S0 = N - I0 - R0; % Initial number of susceptible individuals
% Time parameters
T = 100; % Total simulation time
dt = 1; % Time step
% Initialize arrays to store results
S = zeros(1, T);
I = zeros(1, T);
R = zeros(1, T);
% Initial conditions
S(1) = S0;
I(1) = I0;
R(1) = R0;
% Simulate the SIR model using DTM
for t = 2:T
% Current values
S_curr = S(t-1);
I_curr = I(t-1);
R_curr = R(t-1);
% Compute intermediate values using DTM (predictor step)
S_inter = S_curr - alpha * S_curr * I_curr / N * dt;
I_inter = I_curr + (alpha * S_curr * I_curr / N - beta * I_curr) * dt;
R_inter = R_curr + beta * I_curr * dt;
% Compute new values (corrector step)
S(t) = S_curr - 0.5 * dt * (alpha * S_curr * I_curr / N + alpha * S_inter * I_inter / N);
I(t) = I_curr + 0.5 * dt * ((alpha * S_curr * I_curr / N - beta * I_curr) + (alpha * S_inter * I_inter / N - beta * I_inter));
R(t) = R_curr + 0.5 * dt * (beta * I_curr + beta * I_inter);
end
% Plot the results
t = 0:dt:T-dt;
plot(t, S, 'b', t, I, 'r', t, R, 'g', 'LineWidth', 2);
legend('Susceptible', 'Infectious', 'Recovered');
xlabel('Time');
ylabel('Number of individuals');
title('SIR Model');
end

类别

在 Help Center 和 File Exchange 中查找有关 Optimization Toolbox 的更多信息

标签

产品


版本

R2023b

Community Treasure Hunt

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

Start Hunting!

Translated by