% Sistema de 2 Ecuaciones Diferenciales
% MODELO SIR "SIMPLIFICADO" (Sin ecuación de "Recuperados")

clc
close all
clear all

disp('Sistema de 2 Ecuaciones Diferenciales')
disp('MODELO SIR "SIMPLIFICADO" (Sin ecuación de "Recuperados")')
 
disp('...................................................................')
 
disp('Ecuaciones: ')
 
%disp('S -> dS/dt = -a*S*I') 
%disp('I -> dI/dt = a*S*I - b*I')

disp('S -> dS/dt = -a*S*I') 
disp('I -> dI/dt = a*S*I - b*I')

% S: número de Susceptibles.
% I: número de Infectados.
% a, b: parámetros.

disp('...................................................................')
 
%disp('Valores por defecto:')
%disp('Valor inicial de S [z(1)]: 100) 
%disp('Valor inicial de I [z(2)]: 2)
%disp('Valor de a: 0.01; ')
%disp('Valor de b: 0.63; ')
%disp('tiempo final, tf: 100')
 
%disp('-------------------------------------------------------------------')

z=zeros(1,2);

z(1)=100; % S
z(2)=2;   % I
a=0.01; % EXPLORAR         a = 0.015; 0.02; 0.03; 0.04; 0.05 } -> S > b/a 
b=0.63; % EXPLORAR LUEGO   a = 0.009; 0.008; 0.007; 0.0065; 0.006; 0.005 } -> S < b/a
tf=40;

S0=z(1)
I0=z(2)
format loose
b_dividido_a=b/a

% z es una Matriz de 2 Columnas
% z(:,1) (todos los elementos de la 1era. Columna) representará los sucesivos valores de la función variable "S"
% z(:,2) (todos los elementos de la 2da.  Columna) representará los sucesivos valores de la función variable "I"

% Se definen las funciones "f(t,S,I)" y "g(t,S,I)" -> Ecuaciones Diferenciales "dS/dt" y "dI/dt"

fg=@(t,z) [(-a*z(1)*z(2));(a*z(1)*z(2)-b*z(2))]; % z(1) es "S", z(2) es "I"

tspan=[0:0.1:tf];
[t,z]=ode45(fg,tspan,z);

%**************************************************************************

% "S" vs. "t"
S=z(:,1);
figure(1)
hold on
plot(t,S,'r','LineWidth',2)
hold off
grid on
%axis([-10 110 -10 110])
xlabel('t (días)')
ylabel('S (Susceptibles)');
title('S vs. t')

%**************************************************************************

% "I" vs. "t"
I=z(:,2);
figure(2)
hold on
plot(t,I,'b','LineWidth',2)
hold off
grid on
%axis([-10 tf+10 -10 70])
xlabel('t (días)')
ylabel('I (Infectados)');
title('I vs. t')

%**************************************************************************

% "S" e "I" vs. "t"
figure(3)
plot(t,S,'r',t,I,'b','LineWidth',2)
grid on
%axis([-10 tf+10 -10 S0+I0+10])
xlabel('t (días)')
ylabel('S e I (Susceptibles e Infectados)');
legend('S','I');
title('S e I vs. t')

%**************************************************************************

% I vs. S  (Diagrama de Fases)
figure(4)
hold on
plot(S,I,'g','LineWidth',2)
hold off
grid on
%axis([-10 110 -10 80])
xlabel('S (Susceptibles)');
ylabel('I (Infectados)');
title('I vs. S (Diagrama de Fases)')

%**************************************************************************
%%
% R vs.t
R=S0+I0-S-I;
figure(5)
hold on
plot(t,S,'r',t,I,'b',t,R,'k','LineWidth',2)
grid on
%axis([-10 tf+10 -10 S0+I0+10])
xlabel('t (días)')
ylabel('S e I (Susceptibles e Infectados)');
legend('S','I','R');
%title('S e I vs. t')