21 #ifndef MIO_GLCT_SECIR_MODEL_H
22 #define MIO_GLCT_SECIR_MODEL_H
50 template <
typename FP,
size_t NumExposed,
size_t NumInfectedNoSymptoms,
size_t NumInfectedSymptoms,
51 size_t NumInfectedSevere,
size_t NumInfectedCritical>
55 LctInfectionState<FP, InfectionState, 1, NumExposed, NumInfectedNoSymptoms, NumInfectedSymptoms,
56 NumInfectedSevere, NumInfectedCritical, 1, 1>,
57 mio::Populations<FP, LctInfectionState<FP, InfectionState, 1, NumExposed, NumInfectedNoSymptoms,
58 NumInfectedSymptoms, NumInfectedSevere, NumInfectedCritical, 1, 1>>,
64 NumInfectedSevere, NumInfectedCritical, 1,
87 if ((Eigen::Index)LctState::template get_num_subcompartments<InfectionState::Exposed>() !=
89 log_error(
"Constraint check: Dimension of the parameters does not match the number of subcompartments for "
94 if ((Eigen::Index)LctState::template get_num_subcompartments<InfectionState::InfectedNoSymptoms>() !=
97 "Constraint check: Dimension of the parameters does not match the number of subcompartments for the "
98 "InfectedNoSymptoms compartment.");
101 if ((Eigen::Index)LctState::template get_num_subcompartments<InfectionState::InfectedSymptoms>() !=
104 "Constraint check: Dimension of the parameters does not match the number of subcompartments for the "
105 "InfectedSymptoms compartment.");
108 if ((Eigen::Index)LctState::template get_num_subcompartments<InfectionState::InfectedSevere>() !=
110 log_error(
"Constraint check: Dimension of the parameters does not match the number of subcompartments for "
111 "the InfectedSevere "
115 if ((Eigen::Index)LctState::template get_num_subcompartments<InfectionState::InfectedCritical>() !=
118 "Constraint check: Dimension of the parameters does not match the number of subcompartments for the "
119 "InfectedCritical compartment.");
123 return (params.check_constraints() || this->populations.check_constraints());
138 void get_derivatives(Eigen::Ref<
const Eigen::VectorX<FP>> pop, Eigen::Ref<
const Eigen::VectorX<FP>> y, FP t,
139 Eigen::Ref<Eigen::VectorX<FP>> dydt)
const override
146 auto total_population = pop.sum() - pop[LctState::template get_first_index<InfectionState::Dead>()];
149 FP InfectedNoSymptoms_sum =
150 pop.segment(LctState::template get_first_index<InfectionState::InfectedNoSymptoms>(),
151 LctState::template get_num_subcompartments<InfectionState::InfectedNoSymptoms>())
154 FP InfectedSymptoms_sum =
155 pop.segment(LctState::template get_first_index<InfectionState::InfectedSymptoms>(),
156 LctState::template get_num_subcompartments<InfectionState::InfectedSymptoms>())
161 1 + params.template get<Seasonality<FP>>() *
162 sin(std::numbers::pi_v<ScalarType> * ((params.template
get<
StartDay<FP>>() + t) / 182.5 + 0.5));
164 -y[0] / total_population * season_val * params.template get<TransmissionProbabilityOnContact<FP>>() *
166 (params.template get<RelativeTransmissionNoSymptoms<FP>>() * InfectedNoSymptoms_sum +
170 dydt.segment(LctState::template get_first_index<InfectionState::Exposed>(),
171 LctState::template get_num_subcompartments<InfectionState::Exposed>()) -=
173 dydt.segment(LctState::template get_first_index<InfectionState::Exposed>(),
174 LctState::template get_num_subcompartments<InfectionState::Exposed>()) +=
176 y.segment(LctState::template get_first_index<InfectionState::Exposed>(),
177 LctState::template get_num_subcompartments<InfectionState::Exposed>());
181 dydt.segment(LctState::template get_first_index<InfectionState::InfectedNoSymptoms>(),
182 LctState::template get_num_subcompartments<InfectionState::InfectedNoSymptoms>()) =
184 Eigen::VectorX<FP>::Ones(LctState::template get_num_subcompartments<InfectionState::Exposed>()))
186 y.segment(LctState::template get_first_index<InfectionState::Exposed>(),
187 LctState::template get_num_subcompartments<InfectionState::Exposed>()) *
190 size_t dimensionInfectedNoSymptomsToInfectedSymptoms =
191 params.template get<TransitionMatrixInfectedNoSymptomsToInfectedSymptoms<FP>>().rows();
192 dydt.segment(LctState::template get_first_index<InfectionState::InfectedNoSymptoms>(),
193 dimensionInfectedNoSymptomsToInfectedSymptoms) +=
194 params.template get<TransitionMatrixInfectedNoSymptomsToInfectedSymptoms<FP>>().transpose() *
195 y.segment(LctState::template get_first_index<InfectionState::InfectedNoSymptoms>(),
196 dimensionInfectedNoSymptomsToInfectedSymptoms);
198 size_t dimensionInfectedNoSymptomsToRecovered =
199 params.template get<TransitionMatrixInfectedNoSymptomsToRecovered<FP>>().rows();
200 dydt.segment(LctState::template get_first_index<InfectionState::InfectedNoSymptoms>() +
201 dimensionInfectedNoSymptomsToInfectedSymptoms,
202 dimensionInfectedNoSymptomsToRecovered) +=
203 params.template get<TransitionMatrixInfectedNoSymptomsToRecovered<FP>>().transpose() *
204 y.segment(LctState::template get_first_index<InfectionState::InfectedNoSymptoms>() +
205 dimensionInfectedNoSymptomsToInfectedSymptoms,
206 dimensionInfectedNoSymptomsToRecovered);
208 dydt[LctState::template get_first_index<InfectionState::Recovered>()] +=
209 -(params.template get<TransitionMatrixInfectedNoSymptomsToRecovered<FP>>() *
210 Eigen::VectorX<FP>::Ones(dimensionInfectedNoSymptomsToRecovered))
212 y.segment(LctState::template get_first_index<InfectionState::InfectedNoSymptoms>() +
213 dimensionInfectedNoSymptomsToInfectedSymptoms,
214 dimensionInfectedNoSymptomsToRecovered);
218 dydt.segment(LctState::template get_first_index<InfectionState::InfectedSymptoms>(),
219 LctState::template get_num_subcompartments<InfectionState::InfectedSymptoms>()) =
221 (params.template get<TransitionMatrixInfectedNoSymptomsToInfectedSymptoms<FP>>() *
222 Eigen::VectorX<FP>::Ones(dimensionInfectedNoSymptomsToInfectedSymptoms))
224 y.segment(LctState::template get_first_index<InfectionState::InfectedNoSymptoms>(),
225 dimensionInfectedNoSymptomsToInfectedSymptoms);
227 size_t dimensionInfectedSymptomsToInfectedSevere =
228 params.template get<TransitionMatrixInfectedSymptomsToInfectedSevere<FP>>().rows();
229 dydt.segment(LctState::template get_first_index<InfectionState::InfectedSymptoms>(),
230 dimensionInfectedSymptomsToInfectedSevere) +=
231 params.template get<TransitionMatrixInfectedSymptomsToInfectedSevere<FP>>().transpose() *
232 y.segment(LctState::template get_first_index<InfectionState::InfectedSymptoms>(),
233 dimensionInfectedSymptomsToInfectedSevere);
235 size_t dimensionInfectedSymptomsToRecovered =
236 params.template get<TransitionMatrixInfectedSymptomsToRecovered<FP>>().rows();
237 dydt.segment(LctState::template get_first_index<InfectionState::InfectedSymptoms>() +
238 dimensionInfectedSymptomsToInfectedSevere,
239 dimensionInfectedSymptomsToRecovered) +=
240 params.template get<TransitionMatrixInfectedSymptomsToRecovered<FP>>().transpose() *
241 y.segment(LctState::template get_first_index<InfectionState::InfectedSymptoms>() +
242 dimensionInfectedSymptomsToInfectedSevere,
243 dimensionInfectedSymptomsToRecovered);
245 dydt[LctState::template get_first_index<InfectionState::Recovered>()] +=
246 -(params.template get<TransitionMatrixInfectedSymptomsToRecovered<FP>>() *
247 Eigen::VectorX<FP>::Ones(dimensionInfectedSymptomsToRecovered))
249 y.segment(LctState::template get_first_index<InfectionState::InfectedSymptoms>() +
250 dimensionInfectedSymptomsToInfectedSevere,
251 dimensionInfectedSymptomsToRecovered);
255 dydt.segment(LctState::template get_first_index<InfectionState::InfectedSevere>(),
256 LctState::template get_num_subcompartments<InfectionState::InfectedSevere>()) =
258 (params.template get<TransitionMatrixInfectedSymptomsToInfectedSevere<FP>>() *
259 Eigen::VectorX<FP>::Ones(dimensionInfectedSymptomsToInfectedSevere))
261 y.segment(LctState::template get_first_index<InfectionState::InfectedSymptoms>(),
262 dimensionInfectedSymptomsToInfectedSevere);
264 size_t dimensionInfectedSevereToInfectedCritical =
265 params.template get<TransitionMatrixInfectedSevereToInfectedCritical<FP>>().rows();
266 dydt.segment(LctState::template get_first_index<InfectionState::InfectedSevere>(),
267 dimensionInfectedSevereToInfectedCritical) +=
268 params.template get<TransitionMatrixInfectedSevereToInfectedCritical<FP>>().transpose() *
269 y.segment(LctState::template get_first_index<InfectionState::InfectedSevere>(),
270 dimensionInfectedSevereToInfectedCritical);
272 size_t dimensionInfectedSevereToDead = params.template get<TransitionMatrixInfectedSevereToDead<FP>>().rows();
273 dydt.segment(LctState::template get_first_index<InfectionState::InfectedSevere>() +
274 dimensionInfectedSevereToInfectedCritical,
275 dimensionInfectedSevereToDead) +=
276 params.template get<TransitionMatrixInfectedSevereToDead<FP>>().transpose() *
277 y.segment(LctState::template get_first_index<InfectionState::InfectedSevere>() +
278 dimensionInfectedSevereToInfectedCritical,
279 dimensionInfectedSevereToDead);
281 size_t dimensionInfectedSevereToRecovered =
282 params.template get<TransitionMatrixInfectedSevereToRecovered<FP>>().rows();
283 dydt.segment(LctState::template get_first_index<InfectionState::InfectedSevere>() +
284 dimensionInfectedSevereToInfectedCritical + dimensionInfectedSevereToDead,
285 dimensionInfectedSevereToRecovered) +=
286 params.template get<TransitionMatrixInfectedSevereToRecovered<FP>>().transpose() *
287 y.segment(LctState::template get_first_index<InfectionState::InfectedSevere>() +
288 dimensionInfectedSevereToInfectedCritical + dimensionInfectedSevereToDead,
289 dimensionInfectedSevereToRecovered);
291 dydt[LctState::template get_first_index<InfectionState::Recovered>()] +=
292 -(params.template get<TransitionMatrixInfectedSevereToRecovered<FP>>() *
293 Eigen::VectorX<FP>::Ones(dimensionInfectedSevereToRecovered))
295 y.segment(LctState::template get_first_index<InfectionState::InfectedSevere>() +
296 dimensionInfectedSevereToInfectedCritical + dimensionInfectedSevereToDead,
297 dimensionInfectedSevereToRecovered);
301 dydt.segment(LctState::template get_first_index<InfectionState::InfectedCritical>(),
302 LctState::template get_num_subcompartments<InfectionState::InfectedCritical>()) =
304 (params.template get<TransitionMatrixInfectedSevereToInfectedCritical<FP>>() *
305 Eigen::VectorX<FP>::Ones(dimensionInfectedSevereToInfectedCritical))
307 y.segment(LctState::template get_first_index<InfectionState::InfectedSevere>(),
308 dimensionInfectedSevereToInfectedCritical);
310 size_t dimensionInfectedCriticalToDead =
311 params.template get<TransitionMatrixInfectedCriticalToDead<FP>>().rows();
312 dydt.segment(LctState::template get_first_index<InfectionState::InfectedCritical>(),
313 dimensionInfectedCriticalToDead) +=
314 params.template get<TransitionMatrixInfectedCriticalToDead<FP>>().transpose() *
315 y.segment(LctState::template get_first_index<InfectionState::InfectedCritical>(),
316 dimensionInfectedCriticalToDead);
318 size_t dimensionInfectedCriticalToRecovered =
319 params.template get<TransitionMatrixInfectedCriticalToRecovered<FP>>().rows();
320 dydt.segment(LctState::template get_first_index<InfectionState::InfectedCritical>() +
321 dimensionInfectedCriticalToDead,
322 dimensionInfectedCriticalToRecovered) +=
323 params.template get<TransitionMatrixInfectedCriticalToRecovered<FP>>().transpose() *
324 y.segment(LctState::template get_first_index<InfectionState::InfectedCritical>() +
325 dimensionInfectedCriticalToDead,
326 dimensionInfectedCriticalToRecovered);
328 dydt[LctState::template get_first_index<InfectionState::Recovered>()] +=
329 -(params.template get<TransitionMatrixInfectedCriticalToRecovered<FP>>() *
330 Eigen::VectorX<FP>::Ones(dimensionInfectedCriticalToRecovered))
332 y.segment(LctState::template get_first_index<InfectionState::InfectedCritical>() +
333 dimensionInfectedCriticalToDead,
334 dimensionInfectedCriticalToRecovered);
337 dydt[LctState::template get_first_index<InfectionState::Dead>()] +=
338 -(params.template get<TransitionMatrixInfectedSevereToDead<FP>>() *
339 Eigen::VectorX<FP>::Ones(dimensionInfectedSevereToDead))
341 y.segment(LctState::template get_first_index<InfectionState::InfectedSevere>() +
342 dimensionInfectedSevereToInfectedCritical,
343 dimensionInfectedSevereToDead);
344 dydt[LctState::template get_first_index<InfectionState::Dead>()] +=
345 -(params.template get<TransitionMatrixInfectedCriticalToDead<FP>>() *
346 Eigen::VectorX<FP>::Ones(dimensionInfectedCriticalToDead))
348 y.segment(LctState::template get_first_index<InfectionState::InfectedCritical>(),
349 dimensionInfectedCriticalToDead);
368 log_error(
"Result does not match InfectionStates of the model.");
370 Eigen::VectorX<FP> error_output = Eigen::VectorX<FP>::Constant((Eigen::Index)
InfectionState::Count, -1);
372 return compartments_ts;
379 return compartments_ts;
An Index with more than one template parameter combines several Index objects.
Definition: index.h:191
Provides the functionality to be able to work with subcompartments in an LCT model.
Definition: lct_infection_state.h:42
static constexpr size_t Count
Definition: lct_infection_state.h:114
static Eigen::VectorX< FP > calculate_compartments(const Eigen::VectorX< FP > &subcompartments)
Cumulates a vector with the number of individuals in each subcompartment (with subcompartments accord...
Definition: lct_infection_state.h:93
Typesafe wrapper for a floating-point simulation time value (in days).
Definition: damping.h:78
stores vectors of values at time points (or some other abstract variable) the value at each time poin...
Definition: time_series.h:58
Eigen::Index get_num_elements() const
number of elements of vector at each time point
Definition: time_series.h:205
Eigen::Index get_num_time_points() const
number of time points in the series
Definition: time_series.h:197
FP & get_time(Eigen::Index i)
time of time point at index i
Definition: time_series.h:272
Eigen::Ref< Vector > add_time_point()
add one uninitialized time point
Definition: time_series.h:221
Class that defines an GLCT-SECIR model.
Definition: glct_secir/model.h:60
TimeSeries< FP > calculate_compartments(const TimeSeries< FP > &subcompartments_ts) const
Cumulates a simulation result with subcompartments to produce a result that divides the population on...
Definition: glct_secir/model.h:364
Pop Populations
Definition: compartmental_model.h:62
Params ParameterSet
Definition: compartmental_model.h:63
void get_derivatives(Eigen::Ref< const Eigen::VectorX< FP >> pop, Eigen::Ref< const Eigen::VectorX< FP >> y, FP t, Eigen::Ref< Eigen::VectorX< FP >> dydt) const override
Evaluates the right-hand-side f of the GLCT dydt = f(y, t).
Definition: glct_secir/model.h:138
bool check_constraints() const
Checks that the model satisfies all constraints (e.g.
Definition: glct_secir/model.h:82
Model()
Default constructor.
Definition: glct_secir/model.h:71
Parameters of an GLCT-SECIR model.
Definition: glct_secir/parameters.h:498
ad::internal::unary_intermediate< AD_TAPE_REAL, ad::internal::active_type< AD_TAPE_REAL, DATA_HANDLER_1 >, ad::operations::ad_sin< AD_TAPE_REAL > > sin(const ad::internal::active_type< AD_TAPE_REAL, DATA_HANDLER_1 > &x1)
Definition: ad.hpp:913
InfectionState
The InfectionState enum describes the basic categories for the infection state of persons.
Definition: glct_secir/infection_state.h:31
A collection of classes to simplify handling of matrix shapes in meta programming.
Definition: models/abm/analyze_result.h:30
auto i
Definition: io.h:810
void log_error(spdlog::string_view_t fmt, const Args &... args)
Definition: logging.h:114
constexpr std::tuple_element< I, std::tuple< Index< CategoryTags >... > >::type & get(Index< CategoryTags... > &i) noexcept
Retrieves the Index (by reference) at the Ith position of a MultiIndex.
Definition: index.h:304
CompartmentalModel is a template for a compartmental model for an array of initial populations and a ...
Definition: compartmental_model.h:59
Pop Populations
Definition: compartmental_model.h:62
Params ParameterSet
Definition: compartmental_model.h:63
mio::CompartmentalModel< FP, LctInfectionState< FP, InfectionState, 1, NumExposed, NumInfectedNoSymptoms, NumInfectedSymptoms, NumInfectedSevere, NumInfectedCritical, 1, 1 >, mio::Populations< FP, LctInfectionState< FP, InfectionState, 1, NumExposed, NumInfectedNoSymptoms, NumInfectedSymptoms, NumInfectedSevere, NumInfectedCritical, 1, 1 > >, Parameters< FP > >::parameters ParameterSet parameters
Definition: compartmental_model.h:157
The risk of infection from symptomatic cases in the GLCT-SECIR model.
Definition: glct_secir/parameters.h:433
The start day in the GLCT-SECIR model.
Definition: glct_secir/parameters.h:452
Vector with the probability to start in any of the subcompartments of the Exposed compartment.
Definition: glct_secir/parameters.h:42
Vector with the probability to start in any of the subcompartments of the InfectedCritical compartmen...
Definition: glct_secir/parameters.h:318
Vector with the probability to start in any of the subcompartments of the InfectedNoSymptoms compartm...
Definition: glct_secir/parameters.h:83
Vector with the probability to start in any of the subcompartments of the InfectedSevere compartment.
Definition: glct_secir/parameters.h:223
Vector with the probability to start in any of the subcompartments of the InfectedSymptoms compartmen...
Definition: glct_secir/parameters.h:153
Transition matrix of the Exposed compartment.
Definition: glct_secir/parameters.h:62