-
Notifications
You must be signed in to change notification settings - Fork 0
/
SantiagoHenao_S11C1.cpp
114 lines (93 loc) · 2.21 KB
/
SantiagoHenao_S11C1.cpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
/**
*
* Modelo SIR solucionado con el método de Euler
*
* ./a.out t_0 t_f step N I_0 beta sigma (R_0=0)
*
* t_0 double Tiempo inicial
* t_f double Tiempo final
* step double Resolución temporal
* N int Número total de personas. "Constante" en este modelo
* I_0 int Personas infectadas en t=0
* beta double Tasa de propagación
* sigma double Tasa de remoción (personas que se vuelven inmunes)
* (R_0) int Personas infectadas en t=0. Valor por defecto como 0.
*
**/
#include <iostream>
#include <fstream>
#include <math.h>
#include <vector>
#include <algorithm>
#include <iterator>
using namespace std;
// vectores para guardar los datos
vector<double> S,I,R,t;
int t_i=0;
// constantes del modelo
double beta=0;
double sigma=0;
/*
* Pasos individuales del método de Euler
*/
double Euler_step( vector<double>X ,double (*func)(int) , int time_index , double dt)
{
return X[time_index]+dt*func(time_index);
}
/*
* Función de la ecuación de S
*/
double dR(int time_index)
{
return sigma*I[time_index];
}
/*
* Función de la ecuación de S
*/
double dS(int time_index)
{
return -beta*I[time_index]*S[time_index];
}
/*
* Función de la ecuación de S
*/
double dI(int time_index)
{
return -dS(time_index)-dR(time_index);
}
int main(int argc, char *argv[])
{
// Leer condiciones inciales dadas por parámetro
double t_0=atof(argv[1]);
double t_f=atof(argv[2]);
double step=atof(argv[3]);
int N = atoi(argv[4]);
I.push_back(atof(argv[5]));
beta+=atof(argv[6]);
sigma +=atof(argv[7]);
if(argc==9)
{
R.push_back(atof(argv[8]));
}
else
{
R.push_back(0.0);
}
// Calcular S inicial
S.push_back(N-I[0]-R[0]);
t.push_back(t_0);
while(t[t.size()-1]<=t_f)
{
R.push_back(Euler_step(R,dR,t_i,step));
S.push_back(Euler_step(S,dS,t_i,step));
I.push_back(Euler_step(I,dI,t_i,step));
t_i++;
t.push_back(t[t.size()-1]+step);
}
ofstream result("MIR.dat");
for(int i=0; i < S.size(); i++)
{
result << t[i] << "\t" << S[i] << "\t" << R[i]<< "\t" << I[i] << endl;
}
result.close();
}