Skip to content

Commit 4886fcf

Browse files
lenaploetzkeHenrZu
andauthored
947 example for ide secir model not working as expected (#948)
- Changed the initial transitions so that when using a smaller step size, the resulting number of people in a compartment is no longer negative.. - small corrections of the comments Co-authored-by: Henrik Zunker <69154294+HenrZu@users.noreply.github.com>
1 parent ffa248d commit 4886fcf

File tree

1 file changed

+15
-13
lines changed

1 file changed

+15
-13
lines changed

cpp/examples/ide_secir.cpp

Lines changed: 15 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -20,13 +20,13 @@
2020

2121
#include "ide_secir/model.h"
2222
#include "ide_secir/infection_state.h"
23-
#include "ide_secir/parameters.h"
2423
#include "ide_secir/simulation.h"
2524
#include "memilio/config.h"
2625
#include "memilio/math/eigen.h"
2726
#include "memilio/utils/time_series.h"
2827
#include "memilio/epidemiology/uncertain_matrix.h"
29-
#include <iostream>
28+
#include "memilio/epidemiology/state_age_function.h"
29+
#include "memilio/data/analyze_result.h"
3030

3131
int main()
3232
{
@@ -35,14 +35,14 @@ int main()
3535
ScalarType tmax = 10;
3636
ScalarType N = 10000;
3737
ScalarType deaths = 13.10462213;
38-
ScalarType dt = 1;
38+
ScalarType dt = 1e-2;
3939

4040
int num_transitions = (int)mio::isecir::InfectionTransition::Count;
4141

42-
// create TimeSeries with num_transitions elements where transitions needed for simulation will be stored
42+
// Create TimeSeries with num_transitions elements where transitions needed for simulation will be stored.
4343
mio::TimeSeries<ScalarType> init(num_transitions);
4444

45-
// add time points for initialization of transitions
45+
// Add time points for initialization of transitions.
4646
Vec vec_init(num_transitions);
4747
vec_init[(int)mio::isecir::InfectionTransition::SusceptibleToExposed] = 25.0;
4848
vec_init[(int)mio::isecir::InfectionTransition::ExposedToInfectedNoSymptoms] = 15.0;
@@ -54,21 +54,22 @@ int main()
5454
vec_init[(int)mio::isecir::InfectionTransition::InfectedSevereToRecovered] = 1.0;
5555
vec_init[(int)mio::isecir::InfectionTransition::InfectedCriticalToDead] = 1.0;
5656
vec_init[(int)mio::isecir::InfectionTransition::InfectedCriticalToRecovered] = 1.0;
57-
// add initial time point to time series
57+
vec_init = vec_init * dt;
58+
// Add initial time point to time series.
5859
init.add_time_point(-10, vec_init);
59-
// add further time points until time 0
60-
while (init.get_last_time() < 0) {
61-
vec_init *= 1.01;
60+
// Add further time points until time 0.
61+
while (init.get_last_time() < -dt / 2) {
6262
init.add_time_point(init.get_last_time() + dt, vec_init);
6363
}
6464

6565
// Initialize model.
6666
mio::isecir::Model model(std::move(init), N, deaths);
6767

68+
// Uncomment one of these lines to use a different method to initialize the model using the TimeSeries init.
6869
// model.m_populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Susceptible] = 1000;
6970
// model.m_populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Recovered] = 0;
7071

71-
// Set working parameters
72+
// Set working parameters.
7273
mio::SmootherCosine smoothcos(2.0);
7374
mio::StateAgeFunctionWrapper delaydistribution(smoothcos);
7475
std::vector<mio::StateAgeFunctionWrapper> vec_delaydistrib(num_transitions, delaydistribution);
@@ -97,7 +98,8 @@ int main()
9798
mio::isecir::Simulation sim(model, 0, dt);
9899
sim.advance(tmax);
99100

100-
sim.get_result().print_table({"S", "E", "C", "I", "H", "U", "R", "D "}, 16, 8);
101-
sim.get_transitions().print_table({"S->E", "E->C", "C->I", "C->R", "I->H", "I->R", "H->U", "H->R", "U->D", "U->R"},
102-
16, 8);
101+
auto interpolated_results = mio::interpolate_simulation_result(sim.get_result(), dt / 2);
102+
interpolated_results.print_table({"S", "E", "C", "I", "H", "U", "R", "D "}, 16, 8);
103+
// Uncomment this line to print the transitions.
104+
// sim.get_transitions().print_table({"S->E", "E->C", "C->I", "C->R", "I->H", "I->R", "H->U", "H->R", "U->D", "U->R"}, 16, 8);
103105
}

0 commit comments

Comments
 (0)