diff --git a/cpp/examples/ide_initialization.cpp b/cpp/examples/ide_initialization.cpp index 619c7dc140..bf2508e6e5 100644 --- a/cpp/examples/ide_initialization.cpp +++ b/cpp/examples/ide_initialization.cpp @@ -84,7 +84,7 @@ int main(int argc, char** argv) init.add_time_point(init.get_last_time() + dt, Vec::Constant((int)mio::isecir::InfectionTransition::Count, 1. * dt)); } - model.m_transitions = init; + model.flows = init; } else { // Use the real data for initialization. @@ -100,8 +100,8 @@ int main(int argc, char** argv) sim.advance(2.); // Print results. - 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); + sim.get_flows().print_table({"S->E", "E->C", "C->I", "C->R", "I->H", "I->R", "H->U", "H->R", "U->D", "U->R"}, 16, + 8); return 0; } \ No newline at end of file diff --git a/cpp/models/ide_secir/model.cpp b/cpp/models/ide_secir/model.cpp index f826931d1c..f173b696d3 100644 --- a/cpp/models/ide_secir/model.cpp +++ b/cpp/models/ide_secir/model.cpp @@ -35,43 +35,43 @@ namespace mio namespace isecir { -Model::Model(TimeSeries&& init, CustomIndexArray N_init, - CustomIndexArray deaths, size_t num_agegroups, - CustomIndexArray total_confirmed_cases) +Model::Model(TimeSeries&& flows_init, CustomIndexArray N_init, + CustomIndexArray deaths_init, size_t num_agegroups, + CustomIndexArray total_confirmed_cases_init) : parameters{Parameters(AgeGroup(num_agegroups))} - , m_transitions{std::move(init)} - , m_populations{TimeSeries(Eigen::Index(InfectionState::Count) * num_agegroups)} - , m_total_confirmed_cases{total_confirmed_cases} + , flows{std::move(flows_init)} + , populations{TimeSeries(Eigen::Index(InfectionState::Count) * num_agegroups)} + , total_confirmed_cases{total_confirmed_cases_init} , m_N{N_init} , m_num_agegroups{num_agegroups} { - if (m_transitions.get_num_time_points() > 0) { - // Add first time point in m_populations according to last time point in m_transitions which is where we start + if (flows.get_num_time_points() > 0) { + // Add first time point in populations according to last time point in flows which is where we start // the simulation. - m_populations.add_time_point>( - m_transitions.get_last_time(), + populations.add_time_point>( + flows.get_last_time(), TimeSeries::Vector::Constant((size_t)InfectionState::Count * m_num_agegroups, 0.)); } else { - // Initialize m_populations with zero as the first point of time if no data is provided for the transitions. + // Initialize populations with zero as the first point of time if no data is provided for the flows. // This can happen for example in the case of initialization with real data. - m_populations.add_time_point>( + populations.add_time_point>( 0, TimeSeries::Vector::Constant((size_t)InfectionState::Count * m_num_agegroups, 0.)); } // Set deaths at simulation start time t0. for (AgeGroup group = AgeGroup(0); group < AgeGroup(m_num_agegroups); group++) { - int Di = get_state_flat_index(Eigen::Index(InfectionState::Dead), group); - m_populations[Eigen::Index(0)][Di] = deaths[group]; + int Di = get_state_flat_index(Eigen::Index(InfectionState::Dead), group); + populations[Eigen::Index(0)][Di] = deaths_init[group]; } } bool Model::check_constraints(ScalarType dt) const { - if (!((size_t)m_transitions.get_num_elements() == (size_t)InfectionTransition::Count * m_num_agegroups)) { - log_error("A variable given for model construction is not valid. Number of elements in transition vector " + if (!((size_t)flows.get_num_elements() == (size_t)InfectionTransition::Count * m_num_agegroups)) { + log_error("A variable given for model construction is not valid. Number of elements in vector of flows" "does not match the required number."); return true; } @@ -80,7 +80,7 @@ bool Model::check_constraints(ScalarType dt) const for (int i = 0; i < (int)InfectionState::Count; i++) { int index = get_state_flat_index(i, group); - if (m_populations[0][index] < 0) { + if (populations[0][index] < 0) { log_error("Initialization failed. Initial values for populations are less than zero."); return true; } @@ -88,32 +88,32 @@ bool Model::check_constraints(ScalarType dt) const } // It may be possible to run the simulation with fewer time points, but this number ensures that it is possible. - if (m_transitions.get_num_time_points() < (Eigen::Index)std::ceil(get_global_support_max(dt) / dt)) { - log_error("Initialization failed. Not enough time points for transitions given before start of " + if (flows.get_num_time_points() < (Eigen::Index)std::ceil(get_global_support_max(dt) / dt)) { + log_error("Initialization failed. Not enough time points for flows given before start of " "simulation."); return true; } for (AgeGroup group = AgeGroup(0); group < AgeGroup(m_num_agegroups); ++group) { - for (int i = 0; i < m_transitions.get_num_time_points(); i++) { + for (int i = 0; i < flows.get_num_time_points(); i++) { for (int j = 0; j < (int)InfectionTransition::Count; j++) { int index = get_transition_flat_index(j, group); - if (m_transitions[i][index] < 0) { - log_error("Initialization failed. One or more initial value for transitions is less than zero."); + if (flows[i][index] < 0) { + log_error("Initialization failed. One or more initial value for flows is less than zero."); return true; } } } } - if (m_transitions.get_last_time() != m_populations.get_last_time()) { - log_error("Last time point of TimeSeries for transitions does not match last time point of " + if (flows.get_last_time() != populations.get_last_time()) { + log_error("Last time point of TimeSeries for flows does not match last time point of " "TimeSeries for " "compartments. Both of these time points have to agree for a sensible simulation."); return true; } - if (m_populations.get_num_time_points() != 1) { + if (populations.get_num_time_points() != 1) { log_error("The TimeSeries for the compartments contains more than one time point. It is unclear how to " "initialize."); return true; @@ -179,11 +179,11 @@ void Model::compute_compartment_from_flows(ScalarType dt, Eigen::Index idx_Infec Eigen::Index calc_time_index = (Eigen::Index)std::ceil(calc_time / dt) - 1; - Eigen::Index num_time_points = m_transitions.get_num_time_points(); + Eigen::Index num_time_points = flows.get_num_time_points(); - // Index referring to m_transitions. + // Index referring to flows. int flow_index = get_transition_flat_index(idx_IncomingFlow, (size_t)group); - // Index referring to m_populations. + // Index referring to populations. int state_index = get_state_flat_index(idx_InfectionState, (size_t)group); for (Eigen::Index i = num_time_points - 1 - calc_time_index; i < num_time_points - 1; i++) { @@ -194,10 +194,10 @@ void Model::compute_compartment_from_flows(ScalarType dt, Eigen::Index idx_Infec parameters.get()[group][idx_TransitionDistribution1].eval(state_age) + (1 - parameters.get()[group][idx_TransitionDistribution1]) * parameters.get()[group][idx_TransitionDistribution2].eval(state_age)) * - m_transitions[i + 1][flow_index]; + flows[i + 1][flow_index]; } - m_populations.get_last_value()[state_index] = sum; + populations.get_last_value()[state_index] = sum; } void Model::initial_compute_compartments_infection(ScalarType dt) @@ -242,7 +242,7 @@ void Model::initial_compute_compartments(ScalarType dt) bool susceptibles_given = true; for (AgeGroup group = AgeGroup(0); group < AgeGroup(m_num_agegroups); ++group) { int Si = get_state_flat_index(Eigen::Index(InfectionState::Susceptible), group); - if (m_populations[Eigen::Index(0)][Si] < 1e-12) { + if (populations[Eigen::Index(0)][Si] < 1e-12) { susceptibles_given = false; break; } @@ -250,15 +250,15 @@ void Model::initial_compute_compartments(ScalarType dt) bool recovered_given = true; for (AgeGroup group = AgeGroup(0); group < AgeGroup(m_num_agegroups); ++group) { int Ri = get_state_flat_index(Eigen::Index(InfectionState::Recovered), group); - if (m_populations[Eigen::Index(0)][Ri] < 1e-12) { + if (populations[Eigen::Index(0)][Ri] < 1e-12) { recovered_given = false; break; } } // We check which Initialization method we want to use. - if (!(m_total_confirmed_cases == CustomIndexArray()) && - std::all_of(m_total_confirmed_cases.begin(), m_total_confirmed_cases.end(), [](ScalarType x) { + if (!(total_confirmed_cases == CustomIndexArray()) && + std::all_of(total_confirmed_cases.begin(), total_confirmed_cases.end(), [](ScalarType x) { return x > 1e-12; })) { m_initialization_method = 1; @@ -272,16 +272,16 @@ void Model::initial_compute_compartments(ScalarType dt) int Ri = get_state_flat_index(Eigen::Index(InfectionState::Recovered), group); int Di = get_state_flat_index(Eigen::Index(InfectionState::Dead), group); // The scheme of the ODE model for initialization is applied here. - m_populations[Eigen::Index(0)][Ri] = m_total_confirmed_cases[group] - m_populations[Eigen::Index(0)][ISyi] - - m_populations[Eigen::Index(0)][ISevi] - - m_populations[Eigen::Index(0)][ICri] - - m_populations[Eigen::Index(0)][Eigen::Index(InfectionState::Dead)]; - - m_populations[Eigen::Index(0)][Si] = - m_N[group] - m_populations[Eigen::Index(0)][Ei] - m_populations[Eigen::Index(0)][INSi] - - m_populations[Eigen::Index(0)][ISyi] - m_populations[Eigen::Index(0)][ISevi] - - m_populations[Eigen::Index(0)][ICri] - m_populations[Eigen::Index(0)][Ri] - - m_populations[Eigen::Index(0)][Di]; + populations[Eigen::Index(0)][Ri] = total_confirmed_cases[group] - populations[Eigen::Index(0)][ISyi] - + populations[Eigen::Index(0)][ISevi] - + populations[Eigen::Index(0)][ICri] - + populations[Eigen::Index(0)][Eigen::Index(InfectionState::Dead)]; + + populations[Eigen::Index(0)][Si] = m_N[group] - populations[Eigen::Index(0)][Ei] - + populations[Eigen::Index(0)][INSi] - populations[Eigen::Index(0)][ISyi] - + populations[Eigen::Index(0)][ISevi] - + populations[Eigen::Index(0)][ICri] - populations[Eigen::Index(0)][Ri] - + populations[Eigen::Index(0)][Di]; } } @@ -298,11 +298,11 @@ void Model::initial_compute_compartments(ScalarType dt) int ICri = get_state_flat_index(Eigen::Index(InfectionState::InfectedCritical), group); int Ri = get_state_flat_index(Eigen::Index(InfectionState::Recovered), group); int Di = get_state_flat_index(Eigen::Index(InfectionState::Dead), group); - m_populations[Eigen::Index(0)][Ri] = - m_N[group] - m_populations[Eigen::Index(0)][Si] - m_populations[Eigen::Index(0)][Ei] - - m_populations[Eigen::Index(0)][INSi] - m_populations[Eigen::Index(0)][ISyi] - - m_populations[Eigen::Index(0)][ISevi] - m_populations[Eigen::Index(0)][ICri] - - m_populations[Eigen::Index(0)][Di]; + populations[Eigen::Index(0)][Ri] = m_N[group] - populations[Eigen::Index(0)][Si] - + populations[Eigen::Index(0)][Ei] - populations[Eigen::Index(0)][INSi] - + populations[Eigen::Index(0)][ISyi] - + populations[Eigen::Index(0)][ISevi] - + populations[Eigen::Index(0)][ICri] - populations[Eigen::Index(0)][Di]; } } else if (recovered_given) { @@ -318,11 +318,11 @@ void Model::initial_compute_compartments(ScalarType dt) int ICri = get_state_flat_index(Eigen::Index(InfectionState::InfectedCritical), group); int Ri = get_state_flat_index(Eigen::Index(InfectionState::Recovered), group); int Di = get_state_flat_index(Eigen::Index(InfectionState::Dead), group); - m_populations[Eigen::Index(0)][Si] = - m_N[group] - m_populations[Eigen::Index(0)][Ei] - m_populations[Eigen::Index(0)][INSi] - - m_populations[Eigen::Index(0)][ISyi] - m_populations[Eigen::Index(0)][ISevi] - - m_populations[Eigen::Index(0)][ICri] - m_populations[Eigen::Index(0)][Ri] - - m_populations[Eigen::Index(0)][Di]; + populations[Eigen::Index(0)][Si] = m_N[group] - populations[Eigen::Index(0)][Ei] - + populations[Eigen::Index(0)][INSi] - populations[Eigen::Index(0)][ISyi] - + populations[Eigen::Index(0)][ISevi] - + populations[Eigen::Index(0)][ICri] - populations[Eigen::Index(0)][Ri] - + populations[Eigen::Index(0)][Di]; } } else { @@ -350,16 +350,15 @@ void Model::initial_compute_compartments(ScalarType dt) with a number less than zero, so that a log_error is thrown. However, this initialization method is consistent with the numerical solver of the model equations, so it may sometimes make sense to use this method. */ - m_populations[Eigen::Index(0)][Si] = - m_transitions.get_last_value()[StEi] / (dt * m_forceofinfection[group]); + populations[Eigen::Index(0)][Si] = flows.get_last_value()[StEi] / (dt * m_forceofinfection[group]); // Recovered; calculated as the difference between the total population and the sum of the other // compartment sizes. - m_populations[Eigen::Index(0)][Ri] = - m_N[group] - m_populations[Eigen::Index(0)][Si] - m_populations[Eigen::Index(0)][Ei] - - m_populations[Eigen::Index(0)][INSi] - m_populations[Eigen::Index(0)][ISyi] - - m_populations[Eigen::Index(0)][ISevi] - m_populations[Eigen::Index(0)][ICri] - - m_populations[Eigen::Index(0)][Di]; + populations[Eigen::Index(0)][Ri] = + m_N[group] - populations[Eigen::Index(0)][Si] - populations[Eigen::Index(0)][Ei] - + populations[Eigen::Index(0)][INSi] - populations[Eigen::Index(0)][ISyi] - + populations[Eigen::Index(0)][ISevi] - populations[Eigen::Index(0)][ICri] - + populations[Eigen::Index(0)][Di]; } } else { @@ -378,7 +377,7 @@ void Model::initial_compute_compartments(ScalarType dt) for (AgeGroup group = AgeGroup(0); group < AgeGroup(m_num_agegroups); ++group) { for (Eigen::Index i = 0; i < (Eigen::Index)InfectionState::Count; i++) { int idx = get_state_flat_index(i, group); - if (m_populations[0][idx] < 0) { + if (populations[0][idx] < 0) { m_initialization_method = -2; log_error( "Initialization failed. One or more initial values for populations are less than zero. It may " @@ -396,12 +395,11 @@ void Model::compute_susceptibles(ScalarType dt) { // We need to compute to compute the Susceptibles in every AgeGroup. for (AgeGroup group = AgeGroup(0); group < AgeGroup(m_num_agegroups); ++group) { - Eigen::Index num_time_points = m_populations.get_num_time_points(); + Eigen::Index num_time_points = populations.get_num_time_points(); // Using number of Susceptibles from previous time step and force of infection from previous time step: - // Compute current number of Susceptibles and store Susceptibles in m_populations. - int Si = get_state_flat_index(Eigen::Index(InfectionState::Susceptible), group); - m_populations.get_last_value()[Si] = - m_populations[num_time_points - 2][Si] / (1 + dt * m_forceofinfection[group]); + // Compute current number of Susceptibles and store Susceptibles in populations. + int Si = get_state_flat_index(Eigen::Index(InfectionState::Susceptible), group); + populations.get_last_value()[Si] = populations[num_time_points - 2][Si] / (1 + dt * m_forceofinfection[group]); } } @@ -424,17 +422,17 @@ void Model::compute_flow(Eigen::Index idx_InfectionTransitions, Eigen::Index idx for (Eigen::Index i = current_time_index - calc_time_index; i < current_time_index; i++) { // (current_time_index - i) is the index corresponding to time the individuals have already spent in this state. sum += m_transitiondistributions_derivative[group][idx_InfectionTransitions][current_time_index - i] * - m_transitions[i + 1][idx_IncomingFlow]; + flows[i + 1][idx_IncomingFlow]; } - m_transitions.get_value(current_time_index)[transition_idx] = + flows.get_value(current_time_index)[transition_idx] = (-dt) * parameters.get()[group][idx_InfectionTransitions] * sum; } void Model::compute_flow(Eigen::Index idx_InfectionTransitions, Eigen::Index idx_IncomingFlow, ScalarType dt, AgeGroup group) { - Eigen::Index current_time_index = m_transitions.get_num_time_points() - 1; + Eigen::Index current_time_index = flows.get_num_time_points() - 1; compute_flow(idx_InfectionTransitions, idx_IncomingFlow, dt, current_time_index, group); } @@ -446,7 +444,7 @@ void Model::flows_current_timestep(ScalarType dt) int Si = get_state_flat_index(Eigen::Index(InfectionState::Susceptible), group); // Calculate flow SusceptibleToExposed with force of infection from previous time step and Susceptibles from current time step. - m_transitions.get_last_value()[StEi] = dt * m_forceofinfection[group] * m_populations.get_last_value()[Si]; + flows.get_last_value()[StEi] = dt * m_forceofinfection[group] * populations.get_last_value()[Si]; // Calculate all other flows with compute_flow. // Flow from Exposed to InfectedNoSymptoms @@ -485,17 +483,17 @@ void Model::update_compartment_from_flow(InfectionState infectionState, { int state_idx = get_state_flat_index(Eigen::Index(infectionState), group); - Eigen::Index num_time_points = m_populations.get_num_time_points(); - ScalarType updated_compartment = m_populations[num_time_points - 2][state_idx]; + Eigen::Index num_time_points = populations.get_num_time_points(); + ScalarType updated_compartment = populations[num_time_points - 2][state_idx]; for (const InfectionTransition& inflow : IncomingFlows) { int inflow_idx = get_transition_flat_index(Eigen::Index(inflow), group); - updated_compartment += m_transitions.get_last_value()[inflow_idx]; + updated_compartment += flows.get_last_value()[inflow_idx]; } for (const InfectionTransition& outflow : OutgoingFlows) { int outflow_idx = get_transition_flat_index(Eigen::Index(outflow), group); - updated_compartment -= m_transitions.get_last_value()[outflow_idx]; + updated_compartment -= flows.get_last_value()[outflow_idx]; } - m_populations.get_last_value()[state_idx] = updated_compartment; + populations.get_last_value()[state_idx] = updated_compartment; } void Model::update_compartments() @@ -558,15 +556,15 @@ void Model::compute_forceofinfection(ScalarType dt, bool initialization) ScalarType current_time; if (initialization) { - // Determine m_forceofinfection at time t0-dt which is the penultimate timepoint in m_transitions. - num_time_points = m_transitions.get_num_time_points() - 1; - // Get time of penultimate timepoint in m_transitions. - current_time = m_transitions.get_time(num_time_points - 1); + // Determine m_forceofinfection at time t0-dt which is the penultimate timepoint in flows. + num_time_points = flows.get_num_time_points() - 1; + // Get time of penultimate timepoint in flows. + current_time = flows.get_time(num_time_points - 1); } else { - // Determine m_forceofinfection for current last time in m_transitions. - num_time_points = m_transitions.get_num_time_points(); - current_time = m_transitions.get_last_time(); + // Determine m_forceofinfection for current last time in flows. + num_time_points = flows.get_num_time_points(); + current_time = flows.get_last_time(); } //We compute the Season Value. ScalarType season_val = @@ -598,10 +596,10 @@ void Model::compute_forceofinfection(ScalarType dt, bool initialization) ScalarType deaths_j; if (initialization) { // Determine the number of individuals in Dead compartment at time t0-dt. - deaths_j = m_populations[Eigen::Index(0)][Dj] - m_transitions.get_last_value()[ICrtDj]; + deaths_j = populations[Eigen::Index(0)][Dj] - flows.get_last_value()[ICrtDj]; } else { - deaths_j = m_populations.get_last_value()[Dj]; + deaths_j = populations.get_last_value()[Dj]; } ScalarType sum = 0; @@ -610,15 +608,14 @@ void Model::compute_forceofinfection(ScalarType dt, bool initialization) for (Eigen::Index l = num_time_points - 1 - calc_time_index; l < num_time_points - 1; l++) { Eigen::Index state_age_index = num_time_points - 1 - l; ScalarType state_age = state_age_index * dt; - sum += season_val * parameters.get()[i].eval(state_age) * - parameters.get().get_cont_freq_mat().get_matrix_at(current_time)( - static_cast((size_t)i), static_cast((size_t)j)) * - (m_transitiondistributions_in_forceofinfection[j][0][num_time_points - l - 1] * - m_transitions[l + 1][EtINSj] * - parameters.get()[j].eval(state_age) + - m_transitiondistributions_in_forceofinfection[j][1][num_time_points - l - 1] * - m_transitions[l + 1][INStISyj] * - parameters.get()[j].eval(state_age)); + sum += + season_val * parameters.get()[i].eval(state_age) * + parameters.get().get_cont_freq_mat().get_matrix_at(current_time)( + static_cast((size_t)i), static_cast((size_t)j)) * + (m_transitiondistributions_in_forceofinfection[j][0][num_time_points - l - 1] * + flows[l + 1][EtINSj] * parameters.get()[j].eval(state_age) + + m_transitiondistributions_in_forceofinfection[j][1][num_time_points - l - 1] * + flows[l + 1][INStISyj] * parameters.get()[j].eval(state_age)); } const double divNj = (m_N[j] - deaths_j < Limits::zero_tolerance()) ? 0.0 : 1.0 / (m_N[j] - deaths_j); diff --git a/cpp/models/ide_secir/model.h b/cpp/models/ide_secir/model.h index 9e02bdd50c..6460b3c182 100644 --- a/cpp/models/ide_secir/model.h +++ b/cpp/models/ide_secir/model.h @@ -42,7 +42,7 @@ class Model /** * @brief Constructor to create an IDE-SECIR model. * - * @param[in, out] init TimeSeries with the initial values of the number of individuals, + * @param[in, out] flows_init TimeSeries with the initial values of the number of individuals, * which transit within one timestep dt from one compartment to another. * Possible transitions are specified in #InfectionTransition%s. * Considered time points should have the distance dt. The last time point determines the start time t0 of the @@ -50,14 +50,14 @@ class Model * The time history must reach a certain point in the past so that the simulation can be performed. * A warning is displayed if the condition is violated. * @param[in] N_init A vector, containing the populations of the considered region, for every AgeGroup. - * @param[in] deaths A vector, containing the total number of deaths at time t0, for every AgeGroup. + * @param[in] deaths_init A vector, containing the total number of deaths at time t0, for every AgeGroup. * @param[in] num_agegroups The number of AgeGroups. - * @param[in] total_confirmed_cases A vector, containing the total confirmed cases at time t0 can be set if it + * @param[in] total_confirmed_cases_init A vector, containing the total confirmed cases at time t0 can be set if it * should be used for initialization, for every AgeGroup. */ - Model(TimeSeries&& init, CustomIndexArray N_init, - CustomIndexArray deaths, size_t num_agegroups, - CustomIndexArray total_confirmed_cases = CustomIndexArray()); + Model(TimeSeries&& flows_init, CustomIndexArray N_init, + CustomIndexArray deaths_init, size_t num_agegroups, + CustomIndexArray total_confirmed_cases_init = CustomIndexArray()); // ---- Additional functionality such as constraint checking, setters and getters, etc. ---- /** @@ -67,10 +67,10 @@ class Model bool check_constraints(ScalarType dt) const; /** - * @brief Returns a flat index for the InfectionTransition TimeSeries. + * @brief Returns a flat index for the TimeSeries flows which contains values for the InfectionTransitions. * - * In the TimeSeries we store a vector for each time point. In this vector we store the different - * #InfectionTransition%s for every AgeGroup. + * In the TimeSeries we store a vector for each time point. In this vector we store values for the different + * InfectionTransitions for every AgeGroup. * This function is used to get the right index in this vector for a specific AgeGroup and InfectionTransition. * * @param[in] transition_idx Index determining which InfectionTransition we want to evaluate. @@ -83,10 +83,10 @@ class Model } /** - * @brief Returns a flat index for the InfectionState TimeSeries. + * @brief Returns a flat index for the TimeSeries populations which contains values for the InfectionStates. * - * In the TimeSeries we store a vector for each time point. In this vector we store the different InfectionStates - * for every AgeGroup. + * In the TimeSeries we store a vector for each time point. In this vector we store values for the + * different InfectionStates for every AgeGroup. * This function is used to get the right index in this vector for a specific AgeGroup and InfectionState. * * @param[in] state_idx Index at which InfectionState we want to evaluate. @@ -142,15 +142,15 @@ class Model // ---- Public parameters. ---- ParameterSet parameters{AgeGroup(m_num_agegroups)}; ///< ParameterSet of Model Parameters. - // Attention: m_populations and m_transitions do not necessarily have the same number of time points due to the + // Attention: populations and flows do not necessarily have the same number of time points due to the // initialization part. TimeSeries - m_transitions; ///< TimeSeries containing points of time and the corresponding number of transitions for every - // AgeGroup. - TimeSeries m_populations; ///< TimeSeries containing points of time and the corresponding number of + flows; ///< TimeSeries containing points of time and the corresponding number of individuals transitioning from + // one InfectionState to another as defined in InfectionTransitions for every AgeGroup. + TimeSeries populations; ///< TimeSeries containing points of time and the corresponding number of // people in defined #InfectionState%s for every AgeGroup. CustomIndexArray - m_total_confirmed_cases; ///< CustomIndexArray that contains the total number of confirmed cases at time t0 for every AgeGroup. + total_confirmed_cases; ///< CustomIndexArray that contains the total number of confirmed cases at time t0 for every AgeGroup. private: // ---- Functionality to calculate the sizes of the compartments for time t0. ---- @@ -165,7 +165,7 @@ class Model * @param[in] dt Time discretization step size. * @param[in] idx_InfectionState Specifies the considered #InfectionState * @param[in] group The AgeGroup for which we want to compute. - * @param[in] idx_IncomingFlow Specifies the index of the incoming flow to #InfectionState in m_transitions. + * @param[in] idx_IncomingFlow Specifies the index of the incoming flow to #InfectionState in flows. * @param[in] idx_TransitionDistribution1 Specifies the index of the first relevant TransitionDistribution, * related to a flow from the considered #InfectionState to any other #InfectionState. * This index is also used for related probability. @@ -185,7 +185,7 @@ class Model * * The values for the compartments Exposed, InfectedNoSymptoms, InfectedSymptoms, InfectedSevere and * InfectedCritical for time t_0 are calculated using the initial data in form of flows. - * Calculated values are stored in m_populations. + * Calculated values are stored in populations. * * @param[in] dt Time discretization step size. */ @@ -214,10 +214,10 @@ class Model // ---- Functionality for the iterations of a simulation. ---- /** - * @brief Computes number of Susceptible%s for the current last time in m_populations. + * @brief Computes number of Susceptible%s for the current last time in populations. * * Number is computed using previous number of Susceptible%s and the force of infection (also from previous timestep). - * Number is stored at the matching index in m_populations. + * Number is stored at the matching index in populations. * @param[in] dt Time discretization step size. */ void compute_susceptibles(ScalarType dt); @@ -240,10 +240,10 @@ class Model Eigen::Index current_time_index, AgeGroup group); /** - * @brief Computes size of a flow for the current last time value in m_transitions. + * @brief Computes size of a flow for the current last time value in flows. * * Computes size of one flow from #InfectionTransition, specified in idx_InfectionTransitions, for the current - * last time value in m_transitions. + * last time value in flows. * * @param[in] idx_InfectionTransitions Specifies the considered flow from #InfectionTransition. * @param[in] idx_IncomingFlow Index of the flow in #InfectionTransition, which goes to the considered starting @@ -256,9 +256,9 @@ class Model AgeGroup group); /** - * @brief Sets all required flows for the current last timestep in m_transitions. + * @brief Sets all required flows for the current last timestep in flows. * - * New values are stored in m_transitions. Most values are computed via the function compute_flow(). + * New values are stored in flows. Most values are computed via the function compute_flow(). * * @param[in] dt Time step. */ @@ -267,7 +267,7 @@ class Model /** * @brief Updates the values of one compartment using flows. * - * New value is stored in m_populations. The value is calculated using the compartment size in the previous + * New value is stored in populations. The value is calculated using the compartment size in the previous * time step and the related flows of the current time step. * Therefore the flows of the current time step should be calculated before using this function. */ @@ -278,7 +278,7 @@ class Model /** * @brief Updates the values of all compartments except Susceptible at initialization. * - * New values are stored in m_populations. The values are calculated using the compartment size in the previous + * New values are stored in populations. The values are calculated using the compartment size in the previous * time step and the related flows of the current time step. * Therefore the flows of the current time step should be calculated before using this function. * @@ -286,14 +286,14 @@ class Model void update_compartments(); /** - * @brief Computes force of infection for the current last time in m_transitions. + * @brief Computes force of infection for the current last time in flows. * * Computed value is stored in m_forceofinfection. * * @param[in] dt Time discretization step size. * @param[in] initialization If true we are in the case of the initialization of the model. * For this we need forceofinfection at time point t0-dt and not at the current last time - * (given by m_transitions) as in the other time steps. + * (given by flows) as in the other time steps. */ void compute_forceofinfection(ScalarType dt, bool initialization = false); diff --git a/cpp/models/ide_secir/parameters_io.cpp b/cpp/models/ide_secir/parameters_io.cpp index 39fcd07dc3..beab230825 100644 --- a/cpp/models/ide_secir/parameters_io.cpp +++ b/cpp/models/ide_secir/parameters_io.cpp @@ -71,14 +71,14 @@ IOResult set_initial_flows(Model& model, ScalarType dt, std::string const& // Get the number of AgeGroups. const size_t num_age_groups = ConfirmedCasesDataEntry::age_group_names.size(); - // m_transitions should be empty at the beginning. + // flows should be empty at the beginning. - if (model.m_transitions.get_num_time_points() > 0) { - model.m_transitions = TimeSeries(Eigen::Index(InfectionTransition::Count) * num_age_groups); + if (model.flows.get_num_time_points() > 0) { + model.flows = TimeSeries(Eigen::Index(InfectionTransition::Count) * num_age_groups); } - if (model.m_populations.get_time(0) != 0) { - model.m_populations.remove_last_time_point(); - model.m_populations.add_time_point( + if (model.populations.get_time(0) != 0) { + model.populations.remove_last_time_point(); + model.populations.add_time_point( 0, TimeSeries::Vector::Constant((int)InfectionState::Count * num_age_groups, 0)); } @@ -101,8 +101,8 @@ IOResult set_initial_flows(Model& model, ScalarType dt, std::string const& for (AgeGroup group = AgeGroup(0); group < AgeGroup(num_age_groups); group++) { // Set the Dead compartment to 0 so that RKI data can be added correctly. - int Di = model.get_state_flat_index(Eigen::Index(InfectionState::Dead), group); - model.m_populations[0][Di] = 0; + int Di = model.get_state_flat_index(Eigen::Index(InfectionState::Dead), group); + model.populations[0][Di] = 0; mean_ExposedToInfectedNoSymptoms[group] = model.parameters @@ -138,17 +138,17 @@ IOResult set_initial_flows(Model& model, ScalarType dt, std::string const& // Create TimeSeries with zeros. The index of time zero is start_shift. for (Eigen::Index i = -start_shift; i <= last_time_index_needed; i++) { // Add time point. - model.m_transitions.add_time_point( + model.flows.add_time_point( i * dt, TimeSeries::Vector::Constant((size_t)InfectionTransition::Count * num_age_groups, 0.)); } // Setting the entries in m_total_confirmed_cases to zero before overwriting it with the RKI data. - model.m_total_confirmed_cases = CustomIndexArray(AgeGroup(num_age_groups), 0.); - //--- Calculate the flow InfectedNoSymptomsToInfectedSymptoms using the RKI data and store in the m_transitions object.--- + model.total_confirmed_cases = CustomIndexArray(AgeGroup(num_age_groups), 0.); + //--- Calculate the flow InfectedNoSymptomsToInfectedSymptoms using the RKI data and store in the flows object.--- ScalarType min_offset_needed = std::ceil( - model.m_transitions.get_time(0) - + model.flows.get_time(0) - 1); // Need -1 if first time point is integer and just the floor value if not, therefore use ceil and -1 - ScalarType max_offset_needed = std::ceil(model.m_transitions.get_last_time()); + ScalarType max_offset_needed = std::ceil(model.flows.get_last_time()); bool min_offset_needed_avail = false; bool max_offset_needed_avail = false; @@ -169,32 +169,31 @@ IOResult set_initial_flows(Model& model, ScalarType dt, std::string const& max_offset_needed_avail = true; } // Smallest index for which the entry is needed. - idx_needed_first = - Eigen::Index(std::max(std::floor((offset - model.m_transitions.get_time(0) - 1) / dt), 0.)); + idx_needed_first = Eigen::Index(std::max(std::floor((offset - model.flows.get_time(0) - 1) / dt), 0.)); // Biggest index for which the entry is needed. - idx_needed_last = Eigen::Index(std::min(std::ceil((offset - model.m_transitions.get_time(0) + 1) / dt), - double(model.m_transitions.get_num_time_points() - 1))); + idx_needed_last = Eigen::Index(std::min(std::ceil((offset - model.flows.get_time(0) + 1) / dt), + double(model.flows.get_num_time_points() - 1))); int INStISyi = model.get_transition_flat_index( Eigen::Index(InfectionTransition::InfectedNoSymptomsToInfectedSymptoms), group); for (Eigen::Index i = idx_needed_first; i <= idx_needed_last; i++) { - time_idx = model.m_transitions.get_time(i); + time_idx = model.flows.get_time(i); if (offset == int(std::floor(time_idx))) { - model.m_transitions[i][INStISyi] += + model.flows[i][INStISyi] += (1 - (time_idx - std::floor(time_idx))) * scale_confirmed_cases * entry.num_confirmed; } if (offset == int(std::ceil(time_idx))) { - model.m_transitions[i][INStISyi] += + model.flows[i][INStISyi] += (time_idx - std::floor(time_idx)) * scale_confirmed_cases * entry.num_confirmed; } if (offset == int(std::floor(time_idx - dt))) { - model.m_transitions[i][INStISyi] -= + model.flows[i][INStISyi] -= (1 - (time_idx - dt - std::floor(time_idx - dt))) * scale_confirmed_cases * entry.num_confirmed; } if (offset == int(std::ceil(time_idx - dt))) { - model.m_transitions[i][INStISyi] -= + model.flows[i][INStISyi] -= (time_idx - dt - std::floor(time_idx - dt)) * scale_confirmed_cases * entry.num_confirmed; } } @@ -205,7 +204,7 @@ IOResult set_initial_flows(Model& model, ScalarType dt, std::string const& if (offset == std::floor(-mean_InfectedSymptomsToInfectedSevere[group] - mean_InfectedSevereToInfectedCritical[group] - mean_InfectedCriticalToDead[group])) { - model.m_populations[0][Di] += + model.populations[0][Di] += (1 - (-mean_InfectedSymptomsToInfectedSevere[group] - mean_InfectedSevereToInfectedCritical[group] - mean_InfectedCriticalToDead[group] - @@ -216,7 +215,7 @@ IOResult set_initial_flows(Model& model, ScalarType dt, std::string const& if (offset == std::ceil(-mean_InfectedSymptomsToInfectedSevere[group] - mean_InfectedSevereToInfectedCritical[group] - mean_InfectedCriticalToDead[group])) { - model.m_populations[0][Di] += + model.populations[0][Di] += (-mean_InfectedSymptomsToInfectedSevere[group] - mean_InfectedSevereToInfectedCritical[group] - mean_InfectedCriticalToDead[group] - std::floor(-mean_InfectedSymptomsToInfectedSevere[group] - @@ -225,7 +224,7 @@ IOResult set_initial_flows(Model& model, ScalarType dt, std::string const& } if (offset == 0) { - model.m_total_confirmed_cases[group] = scale_confirmed_cases * entry.num_confirmed; + model.total_confirmed_cases[group] = scale_confirmed_cases * entry.num_confirmed; } } } @@ -297,10 +296,10 @@ IOResult set_initial_flows(Model& model, ScalarType dt, std::string const& int INStISyi = model.get_transition_flat_index( Eigen::Index(InfectionTransition::InfectedNoSymptomsToInfectedSymptoms), group); for (Eigen::Index i = -2 * global_support_max_index; i <= 0; i++) { - model.m_transitions[i + start_shift][EtINSi] = + model.flows[i + start_shift][EtINSi] = (1 / model.parameters.get()[group][Eigen::Index( InfectionTransition::InfectedNoSymptomsToInfectedSymptoms)]) * - model.m_transitions[i + start_shift + index_shift_mean][INStISyi]; + model.flows[i + start_shift + index_shift_mean][INStISyi]; } // Compute flow SusceptibleToExposed for -global_support_max, ..., 0. @@ -311,10 +310,10 @@ IOResult set_initial_flows(Model& model, ScalarType dt, std::string const& int StEi = model.get_transition_flat_index(Eigen::Index(InfectionTransition::SusceptibleToExposed), group); for (Eigen::Index i = -global_support_max_index; i <= 0; i++) { - model.m_transitions[i + start_shift][StEi] = + model.flows[i + start_shift][StEi] = (1 / model.parameters.get()[group][Eigen::Index( InfectionTransition::InfectedNoSymptomsToInfectedSymptoms)]) * - model.m_transitions[i + start_shift + index_shift_mean][INStISyi]; + model.flows[i + start_shift + index_shift_mean][INStISyi]; } // InfectedNoSymptomsToRecovered for -global_support_max, ..., 0. @@ -328,10 +327,10 @@ IOResult set_initial_flows(Model& model, ScalarType dt, std::string const& } // At the end of the calculation, delete all time points that are not required for the simulation. - auto transition_copy(model.m_transitions); - model.m_transitions = TimeSeries(Eigen::Index(InfectionTransition::Count) * num_age_groups); + auto flows_copy(model.flows); + model.flows = TimeSeries(Eigen::Index(InfectionTransition::Count) * num_age_groups); for (Eigen::Index i = -global_support_max_index; i <= 0; i++) { - model.m_transitions.add_time_point(i * dt, transition_copy.get_value(i + start_shift)); + model.flows.add_time_point(i * dt, flows_copy.get_value(i + start_shift)); } return mio::success(); diff --git a/cpp/models/ide_secir/simulation.cpp b/cpp/models/ide_secir/simulation.cpp index 3743feb34b..190f90c86e 100644 --- a/cpp/models/ide_secir/simulation.cpp +++ b/cpp/models/ide_secir/simulation.cpp @@ -30,18 +30,18 @@ namespace isecir void Simulation::advance(ScalarType tmax) { - mio::log_info("Simulating IDE-SECIR from t0 = {} until tmax = {} with dt = {} .", - m_model->m_transitions.get_last_time(), tmax, m_dt); + mio::log_info("Simulating IDE-SECIR from t0 = {} until tmax = {} with dt = {}.", m_model->flows.get_last_time(), + tmax, m_dt); m_model->set_transitiondistributions_support_max(m_dt); m_model->set_transitiondistributions_derivative(m_dt); m_model->set_transitiondistributions_in_forceofinfection(m_dt); m_model->initial_compute_compartments(m_dt); // For every time step: - while (m_model->m_transitions.get_last_time() < tmax - m_dt / 2) { + while (m_model->flows.get_last_time() < tmax - m_dt / 2) { - m_model->m_transitions.add_time_point(m_model->m_transitions.get_last_time() + m_dt); - m_model->m_populations.add_time_point(m_model->m_populations.get_last_time() + m_dt); + m_model->flows.add_time_point(m_model->flows.get_last_time() + m_dt); + m_model->populations.add_time_point(m_model->populations.get_last_time() + m_dt); // compute Susceptibles: m_model->compute_susceptibles(m_dt); diff --git a/cpp/models/ide_secir/simulation.h b/cpp/models/ide_secir/simulation.h index 6c070894ef..a90713ebc4 100644 --- a/cpp/models/ide_secir/simulation.h +++ b/cpp/models/ide_secir/simulation.h @@ -62,7 +62,7 @@ class Simulation */ TimeSeries get_result() { - return m_model->m_populations; + return m_model->populations; } /** @@ -72,17 +72,17 @@ class Simulation */ const TimeSeries& get_result() const { - return m_model->m_populations; + return m_model->populations; } /** - * @brief Get the transitions between the different #InfectionState%s. + * @brief Get the computed flows between the different #InfectionState%s. * - * @return TimeSeries with stored transitions calculated in the simulation. + * @return TimeSeries with stored flows calculated in the simulation. */ - TimeSeries const& get_transitions() + TimeSeries const& get_flows() { - return m_model->m_transitions; + return m_model->flows; } /** diff --git a/cpp/tests/test_ide_parameters_io.cpp b/cpp/tests/test_ide_parameters_io.cpp index dcad2031f9..4afe24f019 100644 --- a/cpp/tests/test_ide_parameters_io.cpp +++ b/cpp/tests/test_ide_parameters_io.cpp @@ -96,16 +96,16 @@ TEST(TestIDEParametersIo, RKIcompareWithPreviousRun) std::vector deaths = { 1, 2.471428571455, 26.34999999999, 603.621428571465, 3972.41428571431, 7668.84999999998}; - std::vector total_confirmed_cases = {10269.2857142857, 29615.8571428571, 185321.571428571, - 215386.428571429, 77163.5714285714, 35588.4285714286}; + std::vector total_confirmed_cases_test = {10269.2857142857, 29615.8571428571, 185321.571428571, + 215386.428571429, 77163.5714285714, 35588.4285714286}; for (mio::AgeGroup group = mio::AgeGroup(0); group < mio::AgeGroup(num_agegroups); ++group) { int Di = model.get_state_flat_index((Eigen::Index)mio::isecir::InfectionState::Dead, group); - EXPECT_NEAR(model.m_populations.get_value(0)[Di], deaths[size_t(group)], 1e-4); - EXPECT_NEAR(model.m_total_confirmed_cases[group], total_confirmed_cases[size_t(group)], 1e-4); + EXPECT_NEAR(model.populations.get_value(0)[Di], deaths[size_t(group)], 1e-4); + EXPECT_NEAR(model.total_confirmed_cases[group], total_confirmed_cases_test[size_t(group)], 1e-4); } - // Compare transitions at last time point with results from a previous run that are given here. + // Compare flows at last time point with results from a previous run that are given here. Eigen::VectorX compare(num_transitions * num_agegroups); compare << 336.428571428600, 328.285714285701, 162.000000000000, 163.071428571425, 80.130989648839, 79.803571428575, 39.476374533415, 39.476374533415, 19.550404043081, 19.550404043081, 1105.714285714297, 1069.857142857200, @@ -119,10 +119,10 @@ TEST(TestIDEParametersIo, RKIcompareWithPreviousRun) 80.130989648839, 79.803571428575, 39.476374533415, 39.476374533415, 19.550404043081, 19.550404043081; mio::isecir::Simulation sim(model, dt); - ASSERT_EQ(compare.size(), model.m_transitions.get_last_value().size()); + ASSERT_EQ(compare.size(), model.flows.get_last_value().size()); for (int j = 0; j < compare.size(); j++) { - ASSERT_NEAR(compare[j], model.m_transitions.get_last_value()[j], 1e-7); + ASSERT_NEAR(compare[j], model.flows.get_last_value()[j], 1e-7); } } @@ -172,8 +172,8 @@ TEST(TestIDEParametersIo, ParametersIoRKIFailure) for (size_t group = 0; group < num_agegroups; ++group) { int INStISy = model.get_transition_flat_index( (Eigen::Index)mio::isecir::InfectionTransition::InfectedNoSymptomsToInfectedSymptoms, group); - for (Eigen::Index i = 0; i < model.m_transitions.get_num_time_points() - 2; i++) { - EXPECT_EQ(0., model.m_transitions.get_value(i)[INStISy]); + for (Eigen::Index i = 0; i < model.flows.get_num_time_points() - 2; i++) { + EXPECT_EQ(0., model.flows.get_value(i)[INStISy]); } } diff --git a/cpp/tests/test_ide_secir.cpp b/cpp/tests/test_ide_secir.cpp index 6eeb39e5d2..a6df999cae 100755 --- a/cpp/tests/test_ide_secir.cpp +++ b/cpp/tests/test_ide_secir.cpp @@ -146,23 +146,23 @@ TEST_F(ModelTestIdeSecir, compareWithPreviousRunTransitions) mio::isecir::Simulation sim(*model, dt); sim.advance(5); - auto transitions = sim.get_transitions(); + auto flows = sim.get_flows(); size_t iter_0 = 0; - while (transitions.get_time(iter_0) < compare[0][0]) { + while (flows.get_time(iter_0) < compare[0][0]) { iter_0++; } for (size_t i = 0; i < compare.size(); i++) { - ASSERT_EQ(compare[i].size(), static_cast(transitions.get_num_elements()) + 1) << "at row " << i; - ASSERT_NEAR(transitions.get_time(i + iter_0), compare[i][0], 1e-7) << "at row " << i; + ASSERT_EQ(compare[i].size(), static_cast(flows.get_num_elements()) + 1) << "at row " << i; + ASSERT_NEAR(flows.get_time(i + iter_0), compare[i][0], 1e-7) << "at row " << i; for (size_t j = 1; j < compare[i].size(); j++) { - ASSERT_NEAR(transitions.get_value(i + iter_0)[j - 1], compare[i][j], 1e-7) << " at row " << i; + ASSERT_NEAR(flows.get_value(i + iter_0)[j - 1], compare[i][j], 1e-7) << " at row " << i; } } } -// Check that the start time of the simulation is determined by the given time points for the transitions. +// Check that the start time of the simulation is determined by the given time points for the flows. TEST(IdeSecir, checkStartTime) { using Vec = mio::TimeSeries::Vector; @@ -178,11 +178,11 @@ TEST(IdeSecir, checkStartTime) int num_transitions = (int)mio::isecir::InfectionTransition::Count; - // Create TimeSeries with num_transitions * num_agegroups elements where transitions needed for simulation + // Create TimeSeries with num_transitions * num_agegroups elements where flows needed for simulation // will be stored. mio::TimeSeries init(num_transitions * num_agegroups); - // Define transitions that will be used for initialization. + // Define flows that will be used for initialization. Vec vec_init = Vec::Constant(num_transitions * num_agegroups, 0.); vec_init[(int)mio::isecir::InfectionTransition::SusceptibleToExposed] = 1.0; vec_init[(int)mio::isecir::InfectionTransition::InfectedNoSymptomsToInfectedSymptoms] = 8.0; @@ -200,8 +200,8 @@ TEST(IdeSecir, checkStartTime) mio::isecir::Simulation sim(model, dt); // Check that the last time point of transitions is equal to t0. - mio::TimeSeries transitions = sim.get_transitions(); - EXPECT_NEAR(t0, transitions.get_last_time(), 1e-8); + mio::TimeSeries flows = sim.get_flows(); + EXPECT_NEAR(t0, flows.get_last_time(), 1e-8); // Carry out simulation and check that first time point of resulting compartments is equal to t0. sim.advance(tmax); @@ -269,18 +269,18 @@ TEST(IdeSecir, checkSimulationFunctions) // Carry out simulation. mio::isecir::Simulation sim(model, dt); sim.advance(tmax); - mio::TimeSeries secihurd_simulated = sim.get_result(); - mio::TimeSeries transitions_simulated = sim.get_transitions(); + mio::TimeSeries secihurd_simulated = sim.get_result(); + mio::TimeSeries flows_simulated = sim.get_flows(); - // Define vectors for compartments and transitions with values from example + // Define vectors for compartments and flows with values from example // (calculated by hand, see internal Overleaf document). // TODO: Add link to material when published. Vec secihurd0((int)mio::isecir::InfectionState::Count); Vec secihurd1((int)mio::isecir::InfectionState::Count); - Vec transitions1(num_transitions); + Vec flows1(num_transitions); secihurd0 << 4995, 0.5, 0, 4, 0, 0, 4990.5, 10; secihurd1 << 4994.00020016, 0.49989992, 0.49994996, 0.12498749, 1.03124687, 0.25781172, 4993.45699802, 10.12890586; - transitions1 << 0.99979984, 0.99989992, 0.24997498, 0.24997498, 2.06249374, 2.06249374, 0.51562344, 0.51562344, + flows1 << 0.99979984, 0.99989992, 0.24997498, 0.24997498, 2.06249374, 2.06249374, 0.51562344, 0.51562344, 0.12890586, 0.12890586; // Compare SECIHURD compartments at times 0 and 1. @@ -289,9 +289,9 @@ TEST(IdeSecir, checkSimulationFunctions) EXPECT_NEAR(secihurd_simulated[1][i], secihurd1[i], 1e-8); } - // Compare transitions at time 1. + // Compare flows at time 1. for (Eigen::Index i = 0; i < num_transitions; i++) { - EXPECT_NEAR(transitions_simulated[transitions_simulated.get_num_time_points() - 1][i], transitions1[i], 1e-8); + EXPECT_NEAR(flows_simulated[flows_simulated.get_num_time_points() - 1][i], flows1[i], 1e-8); } } @@ -311,7 +311,7 @@ TEST(IdeSecir, checkInitializations) int num_transitions = (int)mio::isecir::InfectionTransition::Count; - // Create TimeSeries with num_transitions elements where transitions needed for simulation will be stored. + // Create TimeSeries with num_transitions elements where flows needed for simulation will be stored. mio::TimeSeries init(num_transitions * num_agegroups); // Add initial time point to time series. init.add_time_point(-10, Vec::Constant(num_transitions, 3.0)); @@ -347,8 +347,8 @@ TEST(IdeSecir, checkInitializations) mio::CustomIndexArray(mio::AgeGroup(num_agegroups), 0.)); model2.parameters.get() = mio::UncertainContactMatrix(contact_matrix); - model2.m_populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Susceptible] = 5000; - model2.m_populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Recovered] = 0; + model2.populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Susceptible] = 5000; + model2.populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Recovered] = 0; // Carry out simulation. mio::isecir::Simulation sim2(model2, dt); @@ -363,8 +363,8 @@ TEST(IdeSecir, checkInitializations) mio::CustomIndexArray(mio::AgeGroup(num_agegroups), 0.)); model3.parameters.get() = mio::UncertainContactMatrix(contact_matrix); - model3.m_populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Susceptible] = 0; - model3.m_populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Recovered] = 1000; + model3.populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Susceptible] = 0; + model3.populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Recovered] = 1000; // Carry out simulation. mio::isecir::Simulation sim3(model3, dt); @@ -394,8 +394,8 @@ TEST(IdeSecir, checkInitializations) mio::CustomIndexArray(mio::AgeGroup(num_agegroups), 0.)); model5.parameters.get() = mio::UncertainContactMatrix(contact_matrix); - model5.m_populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Susceptible] = 0; - model5.m_populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Recovered] = 0; + model5.populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Susceptible] = 0; + model5.populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Recovered] = 0; // Carry out simulation. mio::isecir::Simulation sim5(model5, dt); @@ -412,8 +412,8 @@ TEST(IdeSecir, checkInitializations) mio::CustomIndexArray(mio::AgeGroup(num_agegroups), 0.)); model6.parameters.get() = mio::UncertainContactMatrix(contact_matrix); - model6.m_populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Susceptible] = 0; - model6.m_populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Recovered] = 0; + model6.populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Susceptible] = 0; + model6.populations.get_last_value()[(Eigen::Index)mio::isecir::InfectionState::Recovered] = 0; // Carry out simulation. mio::isecir::Simulation sim6(model6, dt); @@ -467,7 +467,7 @@ TEST(IdeSecir, testModelConstraints) // --- Test with negative number of deaths. // Create TimeSeries with num_transitions elements. mio::TimeSeries init(num_transitions); - // Add time points for initialization of transitions. + // Add time points for initialization of flows. Vec vec_init = Vec::Constant(num_transitions, 0.); vec_init[(int)mio::isecir::InfectionTransition::ExposedToInfectedNoSymptoms] = 10.0; init.add_time_point(-3, vec_init); @@ -479,14 +479,14 @@ TEST(IdeSecir, testModelConstraints) // Initialize a model. mio::isecir::Model model_negative_deaths(std::move(init), N, deaths, num_agegroups); - // Return true for negative entry in m_populations. + // Return true for negative entry in populations. constraint_check = model_negative_deaths.check_constraints(dt); EXPECT_TRUE(constraint_check); // --- Test with too few time points. // Create TimeSeries with num_transitions elements. mio::TimeSeries init_few_timepoints(num_transitions); - // Add time points for initialization of transitions. + // Add time points for initialization of flows. init_few_timepoints.add_time_point(-3, vec_init); while (init_few_timepoints.get_last_time() < 0) { init_few_timepoints.add_time_point(init_few_timepoints.get_last_time() + dt, vec_init); @@ -501,40 +501,40 @@ TEST(IdeSecir, testModelConstraints) std::vector vec_delaydistrib(num_transitions, delaydistribution); model.parameters.set(vec_delaydistrib); - // Return true for not enough time points given for the initial transitions. + // Return true for not enough time points given for the initial flows. constraint_check = model.check_constraints(dt); EXPECT_TRUE(constraint_check); - // --- Test with negative transitions. + // --- Test with negative flows. // Create TimeSeries with num_transitions elements. - mio::TimeSeries init_negative_transitions(num_transitions); - // Add time points for initialization of transitions. - init_negative_transitions.add_time_point(-3, vec_init); - while (init_negative_transitions.get_last_time() < 0) { - init_negative_transitions.add_time_point(init_negative_transitions.get_last_time() + dt, (-1) * vec_init); + mio::TimeSeries init_negative_flows(num_transitions); + // Add time points for initialization of flows. + init_negative_flows.add_time_point(-3, vec_init); + while (init_negative_flows.get_last_time() < 0) { + init_negative_flows.add_time_point(init_negative_flows.get_last_time() + dt, (-1) * vec_init); } // Initialize a model. - mio::isecir::Model model_negative_transitions(std::move(init_negative_transitions), N, deaths, num_agegroups); + mio::isecir::Model model_negative_flows(std::move(init_negative_flows), N, deaths, num_agegroups); - // Return true for negative entries in the initial transitions. - constraint_check = model_negative_transitions.check_constraints(dt); + // Return true for negative entries in the initial flows. + constraint_check = model_negative_flows.check_constraints(dt); EXPECT_TRUE(constraint_check); - // --- Test with last time point of transitions not matching last time point of populations. + // --- Test with last time point of flows not matching last time point of populations. // Create TimeSeries with num_transitions elements. mio::TimeSeries init_different_last_time(num_transitions); - // Add enough time points for initialization of transitions but with different last time point - // than before so that it does not match last time point of m_populations (that was set in + // Add enough time points for initialization of flows but with different last time point + // than before so that it does not match last time point of populations (that was set in // when constructing model above). init_different_last_time.add_time_point(-4, vec_init); while (init_different_last_time.get_last_time() < 1) { init_different_last_time.add_time_point(init_different_last_time.get_last_time() + dt, vec_init); } - model.m_transitions = init_different_last_time; + model.flows = init_different_last_time; - // Return true for not last time points of compartments and transitions not matching. + // Return true for not last time points of compartments and flows not matching. constraint_check = model.check_constraints(dt); EXPECT_TRUE(constraint_check); @@ -548,7 +548,7 @@ TEST(IdeSecir, testModelConstraints) populations_many_timepoints.add_time_point(populations_many_timepoints.get_last_time() + dt, vec_populations); } - model.m_populations = populations_many_timepoints; + model.populations = populations_many_timepoints; // Return true for too many time points given for populations. constraint_check = model.check_constraints(dt); @@ -560,7 +560,7 @@ TEST(IdeSecir, testModelConstraints) // Add one time point. correct_populations.add_time_point(1, vec_populations); - model.m_populations = correct_populations; + model.populations = correct_populations; constraint_check = model.check_constraints(dt); EXPECT_FALSE(constraint_check); @@ -719,10 +719,10 @@ TEST(IdeSecir, checkProportionRecoveredDeath) int num_transitions = (int)mio::isecir::InfectionTransition::Count; - // Create TimeSeries with num_transitions elements where transitions needed for simulation will be stored. + // Create TimeSeries with num_transitions elements where flows needed for simulation will be stored. mio::TimeSeries init(num_transitions * num_agegroups); - // Add time points for initialization for transitions. + // Add time points for initialization for flows. Vec vec_init = Vec::Constant(num_transitions * num_agegroups, 0.); vec_init[(int)mio::isecir::InfectionTransition::ExposedToInfectedNoSymptoms] = 10.0; vec_init[(int)mio::isecir::InfectionTransition::InfectedSymptomsToInfectedSevere] = 10.0; @@ -815,10 +815,10 @@ TEST(IdeSecir, compareEquilibria) int num_transitions = (int)mio::isecir::InfectionTransition::Count; - // Create TimeSeries with num_transitions elements where transitions needed for simulation will be stored. + // Create TimeSeries with num_transitions elements where flows needed for simulation will be stored. mio::TimeSeries init(num_transitions); - // Add time points for initialization for transitions. + // Add time points for initialization for flows. Vec vec_init = Vec::Constant(num_transitions * num_agegroups, 0.); vec_init[(int)mio::isecir::InfectionTransition::ExposedToInfectedNoSymptoms] = 10.0; vec_init[(int)mio::isecir::InfectionTransition::InfectedSymptomsToInfectedSevere] = 10.0; diff --git a/cpp/tests/test_ide_secir_ageres.cpp b/cpp/tests/test_ide_secir_ageres.cpp index e7048f4181..8c72ba3e0d 100644 --- a/cpp/tests/test_ide_secir_ageres.cpp +++ b/cpp/tests/test_ide_secir_ageres.cpp @@ -28,10 +28,10 @@ TEST(TestIdeAgeres, compareWithPreviousRun) int num_transitions = (int)mio::isecir::InfectionTransition::Count; int num_compartments = (int)mio::isecir::InfectionState::Count; - // Create TimeSeries with num_transitions * num_agegroups elements where transitions needed for simulation will be stored. + // Create TimeSeries with num_transitions * num_agegroups elements where flows needed for simulation will be stored. mio::TimeSeries init(num_transitions * num_agegroups); - // Define transitions that will be used for initialization. + // Define flows that will be used for initialization. Vec vec_init = Vec::Constant(num_transitions * num_agegroups, 1.); // First AgeGroup. vec_init[(int)mio::isecir::InfectionTransition::SusceptibleToExposed] = 20.0; @@ -140,19 +140,19 @@ TEST(TestIdeAgeres, compareWithPreviousRun) ASSERT_NEAR(compartments.get_last_value()[j], compare_compartments[j], 1e-7); } - // Compare transitions at last time point with results from a previous run that are given here. + // Compare flows at last time point with results from a previous run that are given here. - mio::TimeSeries transitions = sim.get_transitions(); - Eigen::VectorX compare_transitions(num_transitions * num_agegroups); - compare_transitions << 31.5370062111, 30.6497959470, 14.1231866958, 14.7543908776, 6.6982921386, 6.6982921386, + mio::TimeSeries flows = sim.get_flows(); + Eigen::VectorX compare_flows(num_transitions * num_agegroups); + compare_flows << 31.5370062111, 30.6497959470, 14.1231866958, 14.7543908776, 6.6982921386, 6.6982921386, 3.1606794140, 3.1606794140, 1.4742153411, 1.4742153411, 31.5370062111, 29.5087817552, 14.7543908776, 14.1231866958, 6.3213588280, 6.3213588280, 2.9484306823, 2.9484306823, 1.3533839877, 1.3533839877, 39.4212577639, 30.1092463410, 14.4459531059, 14.4459531059, 6.5114345346, 6.5114345346, 3.0573621415, 3.0573621415, 1.4155355888, 1.4155355888; - ASSERT_EQ(compare_transitions.size(), static_cast(transitions.get_last_value().size())); + ASSERT_EQ(compare_flows.size(), static_cast(flows.get_last_value().size())); - for (int j = 0; j < compare_transitions.size(); j++) { - ASSERT_NEAR(transitions.get_last_value()[j], compare_transitions[j], 1e-7); + for (int j = 0; j < compare_flows.size(); j++) { + ASSERT_NEAR(flows.get_last_value()[j], compare_flows[j], 1e-7); } }