% 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('-------------------------------------------------------------------')
M=zeros(1001,2,6);
N=140;
for i=1:6
z=zeros(1,2);
z(2)=(i-1)*20+20;   % I
z(1)=N-z(2); % S

a=0.01; % EXPLORAR         a = 0.015; 0.02; 0.03; 0.04; 0.05 } -> S > b/a 
b=0.60; % EXPLORAR LUEGO   a = 0.009; 0.008; 0.007; 0.0065; 0.006; 0.005 } -> S < b/a
tf=100;

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);
M(:,:,i)=z;
end
A=0:N;
B=N*ones(1,length(A))-A;
C=(b/a)*ones(1,N);
D=1:N;
figure(1)
hold on
for i=1:6
plot(M(:,1,i),M(:,2,i),'k')
end
plot(C,D,'--k',A,B,'k')
axis([0 150 0 150])
ylabel('I')
xlabel('S')

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