Skip to content

Commit afd4982

Browse files
936 implement the class infectionstate of the lct model more efficiently (#941)
- LctInfectionState is now a class template and has only static members (variables and functions). More computation is now done at compile time. - LctInfectionState can be used by various LCT models (not just SECIR models). Co-authored-by: reneSchm <[email protected]>
1 parent c484b89 commit afd4982

File tree

11 files changed

+476
-625
lines changed

11 files changed

+476
-625
lines changed

cpp/examples/lct_secir.cpp

+69-26
Original file line numberDiff line numberDiff line change
@@ -25,42 +25,85 @@
2525
#include "memilio/utils/time_series.h"
2626
#include "memilio/epidemiology/uncertain_matrix.h"
2727
#include "memilio/math/eigen.h"
28+
#include "memilio/utils/logging.h"
29+
2830
#include <vector>
2931

3032
int main()
3133
{
32-
/** Simple example to demonstrate how to run a simulation using an LCT SECIR model.
33-
Parameters, initial values and subcompartments are not meant to represent a realistic scenario. */
34+
// Simple example to demonstrate how to run a simulation using an LCT SECIR model.
35+
// Parameters, initial values and the number of subcompartments are not meant to represent a realistic scenario.
3436

35-
// Set vector that specifies the number of subcompartments.
36-
std::vector<int> num_subcompartments((int)mio::lsecir::InfectionStateBase::Count, 1);
37-
num_subcompartments[(int)mio::lsecir::InfectionStateBase::Exposed] = 2;
38-
num_subcompartments[(int)mio::lsecir::InfectionStateBase::InfectedNoSymptoms] = 3;
39-
num_subcompartments[(int)mio::lsecir::InfectionStateBase::InfectedCritical] = 5;
40-
mio::lsecir::InfectionState infection_state(num_subcompartments);
37+
using Model = mio::lsecir::Model<2, 3, 1, 1, 5>;
38+
using LctState = Model::LctState;
4139

4240
ScalarType tmax = 20;
4341

44-
// Define initial distribution of the population in the subcompartments.
45-
Eigen::VectorXd init(infection_state.get_count());
46-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::Susceptible)] = 750;
47-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::Exposed)] = 30;
48-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::Exposed) + 1] = 20;
49-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::InfectedNoSymptoms)] = 20;
50-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::InfectedNoSymptoms) + 1] = 10;
51-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::InfectedNoSymptoms) + 2] = 10;
52-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::InfectedSymptoms)] = 50;
53-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::InfectedSevere)] = 50;
54-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::InfectedCritical)] = 10;
55-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::InfectedCritical) + 1] = 10;
56-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::InfectedCritical) + 2] = 5;
57-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::InfectedCritical) + 3] = 3;
58-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::InfectedCritical) + 4] = 2;
59-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::Recovered)] = 20;
60-
init[infection_state.get_firstindex(mio::lsecir::InfectionStateBase::Dead)] = 10;
42+
// Define the initial value vector init with the distribution of the population into subcompartments.
43+
// This method of defining the vector using a vector of vectors is a bit of overhead, but should remind you how
44+
// the entries of the initial value vector relate to the defined template parameters of the model or the number of subcompartments.
45+
// It is also possible to define the initial value vector directly.
46+
std::vector<std::vector<ScalarType>> initial_populations = {{750}, {30, 20}, {20, 10, 10}, {50},
47+
{50}, {10, 10, 5, 3, 2}, {20}, {10}};
48+
49+
// Assert that initial_populations has the right shape.
50+
if (initial_populations.size() != (int)LctState::InfectionState::Count) {
51+
mio::log_error("The number of vectors in initial_populations does not match the number of InfectionStates.");
52+
return 1;
53+
}
54+
if ((initial_populations[(int)LctState::InfectionState::Susceptible].size() !=
55+
LctState::get_num_subcompartments<LctState::InfectionState::Susceptible>()) ||
56+
(initial_populations[(int)LctState::InfectionState::Exposed].size() !=
57+
LctState::get_num_subcompartments<LctState::InfectionState::Exposed>()) ||
58+
(initial_populations[(int)LctState::InfectionState::InfectedNoSymptoms].size() !=
59+
LctState::get_num_subcompartments<LctState::InfectionState::InfectedNoSymptoms>()) ||
60+
(initial_populations[(int)LctState::InfectionState::InfectedSymptoms].size() !=
61+
LctState::get_num_subcompartments<LctState::InfectionState::InfectedSymptoms>()) ||
62+
(initial_populations[(int)LctState::InfectionState::InfectedSevere].size() !=
63+
LctState::get_num_subcompartments<LctState::InfectionState::InfectedSevere>()) ||
64+
(initial_populations[(int)LctState::InfectionState::InfectedCritical].size() !=
65+
LctState::get_num_subcompartments<LctState::InfectionState::InfectedCritical>()) ||
66+
(initial_populations[(int)LctState::InfectionState::Recovered].size() !=
67+
LctState::get_num_subcompartments<LctState::InfectionState::Recovered>()) ||
68+
(initial_populations[(int)LctState::InfectionState::Dead].size() !=
69+
LctState::get_num_subcompartments<LctState::InfectionState::Dead>())) {
70+
mio::log_error("The length of at least one vector in initial_populations does not match the related number of "
71+
"subcompartments.");
72+
return 1;
73+
}
74+
75+
// Transfer the initial values in initial_populations to the vector init.
76+
Eigen::VectorXd init = Eigen::VectorXd::Zero(LctState::Count);
77+
init[(int)LctState::get_first_index<LctState::InfectionState::Susceptible>()] =
78+
initial_populations[(int)LctState::InfectionState::Susceptible][0];
79+
for (unsigned int i = 0; i < LctState::get_num_subcompartments<LctState::InfectionState::Exposed>(); i++) {
80+
init[(int)LctState::get_first_index<LctState::InfectionState::Exposed>() + i] =
81+
initial_populations[(int)LctState::InfectionState::Exposed][i];
82+
}
83+
for (unsigned int i = 0; i < LctState::get_num_subcompartments<LctState::InfectionState::InfectedNoSymptoms>();
84+
i++) {
85+
init[(int)LctState::get_first_index<LctState::InfectionState::InfectedNoSymptoms>() + i] =
86+
initial_populations[(int)LctState::InfectionState::InfectedNoSymptoms][i];
87+
}
88+
for (unsigned int i = 0; i < LctState::get_num_subcompartments<LctState::InfectionState::InfectedSymptoms>(); i++) {
89+
init[(int)LctState::get_first_index<LctState::InfectionState::InfectedSymptoms>() + i] =
90+
initial_populations[(int)LctState::InfectionState::InfectedSymptoms][i];
91+
}
92+
for (unsigned int i = 0; i < LctState::get_num_subcompartments<LctState::InfectionState::InfectedSevere>(); i++) {
93+
init[(int)LctState::get_first_index<LctState::InfectionState::InfectedSevere>() + i] =
94+
initial_populations[(int)LctState::InfectionState::InfectedSevere][i];
95+
}
96+
for (unsigned int i = 0; i < LctState::get_num_subcompartments<LctState::InfectionState::InfectedCritical>(); i++) {
97+
init[(int)LctState::get_first_index<LctState::InfectionState::InfectedCritical>() + i] =
98+
initial_populations[(int)LctState::InfectionState::InfectedCritical][i];
99+
}
100+
init[(int)LctState::get_first_index<LctState::InfectionState::Recovered>()] =
101+
initial_populations[(int)LctState::InfectionState::Recovered][0];
102+
init[(int)LctState::get_first_index<LctState::InfectionState::Dead>()] =
103+
initial_populations[(int)LctState::InfectionState::Dead][0];
61104

62105
// Initialize model.
63-
mio::lsecir::Model model(std::move(init), infection_state);
106+
Model model(std::move(init));
64107

65108
// Set Parameters.
66109
model.parameters.get<mio::lsecir::TimeExposed>() = 3.2;

cpp/memilio/CMakeLists.txt

+1
Original file line numberDiff line numberDiff line change
@@ -16,6 +16,7 @@ add_library(memilio
1616
epidemiology/damping_sampling.cpp
1717
epidemiology/dynamic_npis.h
1818
epidemiology/dynamic_npis.cpp
19+
epidemiology/lct_infection_state.h
1920
geography/regions.h
2021
geography/regions.cpp
2122
epidemiology/simulation_day.h
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,86 @@
1+
/*
2+
* Copyright (C) 2020-2024 MEmilio
3+
*
4+
* Authors: Lena Ploetzke
5+
*
6+
* Contact: Martin J. Kuehn <[email protected]>
7+
*
8+
* Licensed under the Apache License, Version 2.0 (the "License");
9+
* you may not use this file except in compliance with the License.
10+
* You may obtain a copy of the License at
11+
*
12+
* http://www.apache.org/licenses/LICENSE-2.0
13+
*
14+
* Unless required by applicable law or agreed to in writing, software
15+
* distributed under the License is distributed on an "AS IS" BASIS,
16+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
17+
* See the License for the specific language governing permissions and
18+
* limitations under the License.
19+
*/
20+
#ifndef MIO_EPI_LCT_INFECTION_STATE_H
21+
#define MIO_EPI_LCT_INFECTION_STATE_H
22+
23+
#include <array>
24+
25+
namespace mio
26+
{
27+
/**
28+
* @brief Provides the functionality to be able to work with subcompartments in an LCT model.
29+
*
30+
* @tparam InfectionStates An enum class that defines the basic infection states.
31+
* @tparam Ns Number of subcompartments for each infection state defined in InfectionState.
32+
* The number of given template arguments must be equal to the entry Count from InfectionState.
33+
*/
34+
template <class InfectionStates, unsigned int... Ns>
35+
class LctInfectionState
36+
{
37+
public:
38+
using InfectionState = InfectionStates;
39+
static_assert((unsigned int)InfectionState::Count == sizeof...(Ns),
40+
"The number of integers provided as template parameters must be "
41+
"the same as the entry Count of InfectionState.");
42+
43+
static_assert(((Ns > 0) && ...), "The number of subcompartments must be at least 1.");
44+
45+
/**
46+
* @brief Gets the number of subcompartments in an infection state.
47+
*
48+
* @tparam State: Infection state for which the number of subcompartments should be returned.
49+
* @return Number of subcompartments for State.
50+
*/
51+
template <InfectionState State>
52+
static constexpr unsigned int get_num_subcompartments()
53+
{
54+
static_assert(State < InfectionState::Count, "State must be a a valid InfectionState.");
55+
return m_subcompartment_numbers[(int)State];
56+
}
57+
58+
/**
59+
* @brief Gets the index of the first subcompartment of an infection state.
60+
*
61+
* In a simulation, the number of individuals in the subcompartments are stored in vectors.
62+
* Accordingly, the index of the first subcompartment of State in such a vector is returned.
63+
* @tparam State: Infection state for which the index should be returned.
64+
* @return Index of the first subcompartment for a vector with one entry per subcompartment.
65+
*/
66+
template <InfectionState State>
67+
static constexpr unsigned int get_first_index()
68+
{
69+
static_assert(State < InfectionState::Count, "State must be a a valid InfectionState.");
70+
unsigned int index = 0;
71+
for (int i = 0; i < (int)(State); i++) {
72+
index = index + m_subcompartment_numbers[i];
73+
}
74+
return index;
75+
}
76+
77+
static constexpr unsigned int Count{(... + Ns)};
78+
79+
private:
80+
static constexpr const std::array<unsigned int, sizeof...(Ns)> m_subcompartment_numbers{
81+
Ns...}; ///< Vector which defines the number of subcompartments for each infection state of InfectionState.
82+
};
83+
84+
} // namespace mio
85+
86+
#endif

cpp/models/lct_secir/CMakeLists.txt

+1-2
Original file line numberDiff line numberDiff line change
@@ -2,8 +2,7 @@ add_library(lct_secir
22
infection_state.h
33
model.h
44
model.cpp
5-
simulation.h
6-
simulation.cpp
5+
simulation.h
76
parameters.h
87
)
98
target_link_libraries(lct_secir PUBLIC memilio)

cpp/models/lct_secir/README.md

+5-5
Original file line numberDiff line numberDiff line change
@@ -35,11 +35,11 @@ Below is an overview of the model architecture and its compartments.
3535
| $\xi_{I_{Sy}}$ | `RiskOfInfectionFromSymptomatic` | Proportion of infected people with symptomps who are not isolated. |
3636
| $N$ | `m_N0` | Total population. |
3737
| $D$ | `D` | Number of death people. |
38-
| $n_E$ | Defined via `InfectionState` | Number of subcompartments of the Exposed compartment. |
39-
| $n_{NS}$ | Defined via `InfectionState` | Number of subcompartments of the InfectedNoSymptoms compartment. |
40-
| $n_{Sy}$ | Defined via `InfectionState` | Number of subcompartments of the InfectedSymptoms compartment. |
41-
| $n_{Sev}$ | Defined via `InfectionState` | Number of subcompartments of the InfectedSevere compartment.|
42-
| $n_{Cr}$ | Defined via `InfectionState` | Number of subcompartments of the InfectedCritical compartment. |
38+
| $n_E$ | `NumExposed` | Number of subcompartments of the Exposed compartment. |
39+
| $n_{NS}$ | `NumInfectedNoSymptoms` | Number of subcompartments of the InfectedNoSymptoms compartment. |
40+
| $n_{Sy}$ | `NumInfectedSymptoms` | Number of subcompartments of the InfectedSymptoms compartment. |
41+
| $n_{Sev}$ |`NumInfectedSevere` | Number of subcompartments of the InfectedSevere compartment.|
42+
| $n_{Cr}$ | `NumInfectedCritical` | Number of subcompartments of the InfectedCritical compartment. |
4343
| $T_E$ | `TimeExposed` | Average time in days an individual stays in the Exposed compartment. |
4444
| $T_{I_{NS}}$ | `TimeInfectedNoSymptoms` | Average time in days an individual stays in the InfectedNoSymptoms compartment. |
4545
| $T_{I_{Sy}}$ | `TimeInfectedSymptoms` | Average time in days an individual stays in the InfectedSymptoms compartment. |

0 commit comments

Comments
 (0)