1003 Add ode seir metapop bindings

cpp/examples/ode_seir_metapop.cpp

@@ -68,12 +68,10 @@ int main()
     for (size_t index = 0; index < (size_t)interpolated_result.get_num_time_points(); index++) {
         printf("\n %f", interpolated_result.get_times()[index]);
         for (size_t i = 0; i < (size_t)model.parameters.get_num_regions(); i++) {
-
-            printf("\t %.5f ",
-                   interpolated_result.get_value(index)[(size_t)model.parameters.get_num_regions() *
-                                                            ((size_t)mio::oseirmetapop::InfectionState::Count - 1) +
-                                                        i] /
-                       model.populations.get_group_total(mio::regions::Region(i)) * 100);
+            auto infected_idx = model.populations.get_flat_index(
+                {mio::regions::Region(i), mio::AgeGroup(0), mio::oseirmetapop::InfectionState::Infected});
+            printf("\t %.5f ", interpolated_result.get_value(index)[infected_idx] /
+                                   model.populations.get_group_total(mio::regions::Region(i)) * 100);
         }
     }
     printf("\n");

cpp/memilio/math/stepper_wrapper.h

@@ -53,8 +53,9 @@ struct step_adjuster : public boost::numeric::odeint::default_step_adjuster<Valu
  * @brief This is an adaptive IntegratorCore. It creates and manages an instance of a
  * boost::numeric::odeint::controlled_runge_kutta integrator, wrapped as mio::IntegratorCore.
  */
-template <typename FP, template <class State, class Value, class Deriv, class Time, class Algebra, class Operations,
-                                 class Resizer> class ControlledStepper>
+template <typename FP,
+          template <class State, class Value, class Deriv, class Time, class Algebra, class Operations, class Resizer>
+          class ControlledStepper>
 class ControlledStepperWrapper : public mio::OdeIntegratorCore<FP>
 {
     using Algebra      = boost::numeric::odeint::vector_space_algebra;
@@ -204,8 +205,9 @@ class ControlledStepperWrapper : public mio::OdeIntegratorCore<FP>
  * @brief This is a non-adaptive IntegratorCore. It creates and manages an instance of an explicit stepper from
  * boost::numeric::odeint, wrapped as mio::IntegratorCore.
  */
-template <typename FP, template <class State, class Value, class Deriv, class Time, class Algebra, class Operations,
-                                 class Resizer> class ExplicitStepper>
+template <typename FP,
+          template <class State, class Value, class Deriv, class Time, class Algebra, class Operations, class Resizer>
+          class ExplicitStepper>
 class ExplicitStepperWrapper : public mio::OdeIntegratorCore<FP>
 {
 public:

cpp/models/ode_seir_metapop/infection_state.h

@@ -17,7 +17,6 @@
 * See the License for the specific language governing permissions and
 * limitations under the License.
 */
-
 #ifndef SEIRMETAPOP_INFECTIONSTATE_H
 #define SEIRMETAPOP_INFECTIONSTATE_H
 

pycode/examples/simulation/ode_seir_metapop.py

@@ -0,0 +1,83 @@
+#############################################################################
+# Copyright (C) 2020-2026 MEmilio
+#
+# Authors: Carlotta Gerstein
+#
+# Contact: Martin J. Kuehn <Martin.Kuehn@DLR.de>
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+#     http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+#############################################################################
+
+import numpy as np
+
+import memilio.simulation as mio
+import memilio.simulation.oseir_metapop as oseir_metapop
+from memilio.simulation import AgeGroup, LogLevel, set_log_level
+set_log_level(LogLevel.Off)
+
+t0 = 0
+tmax = 10
+dt = 0.1
+
+model = oseir_metapop.Model(3, 1)
+
+for i in range(model.parameters.num_regions.get()):
+    # Set infection state stay times (in days)
+    model.populations[mio.Region(i), AgeGroup(
+        0), oseir_metapop.InfectionState.Susceptible] = 10000
+
+model.populations[mio.Region(0), AgeGroup(
+    0), oseir_metapop.InfectionState.Infected].value += 100
+model.populations[mio.Region(0), AgeGroup(
+    0), oseir_metapop.InfectionState.Susceptible].value -= 100
+
+mobility_data_commuter = np.array(
+    [[0.4, 0.3, 0.3], [0.2, 0.7, 0.1], [0.4, 0.1, 0.5]])
+model.set_commuting_strengths(mobility_data_commuter)
+
+# Set contact frequency
+model.parameters.ContactPatterns.cont_freq_mat[0].baseline = np.ones(
+    (3, 3)) * 2.7
+
+model.parameters.TimeExposed[mio.Region(0), AgeGroup(0)] = 3.
+model.parameters.TimeExposed[mio.Region(1), AgeGroup(0)] = 4.
+model.parameters.TimeExposed[mio.Region(2), AgeGroup(0)] = 5.
+model.parameters.TimeInfected[mio.Region(0), AgeGroup(0)] = 7.
+model.parameters.TimeInfected[mio.Region(1), AgeGroup(0)] = 8.
+model.parameters.TimeInfected[mio.Region(2), AgeGroup(0)] = 9.
+model.parameters.TransmissionProbabilityOnContact[mio.Region(
+    0), AgeGroup(0)] = 0.07333
+model.parameters.TransmissionProbabilityOnContact[mio.Region(
+    1), AgeGroup(0)] = 0.07333
+model.parameters.TransmissionProbabilityOnContact[mio.Region(
+    2), AgeGroup(0)] = 0.07333
+
+result = oseir_metapop.simulate(t0, tmax, dt, model)
+interpolated_result = oseir_metapop.interpolate_simulation_result(result)
+
+print("Infected individuals per Region over time [%]:")
+print("".join(f"\t Region {i}" for i in range(
+    model.parameters.num_regions.get())))
+
+for t_idx in range(interpolated_result.get_num_time_points()):
+    row = f" {interpolated_result.get_times()[t_idx]:.6f}"
+    values = interpolated_result.get_value(t_idx)
+
+    for i in range(model.parameters.num_regions.get()):
+        population = model.populations.get_group_total_Region(mio.Region(i))
+        infected_idx = model.populations.get_flat_index(
+            (mio.Region(i), AgeGroup(0), oseir_metapop.InfectionState.Infected))
+        percentage = values[infected_idx] / population * 100
+        row += f"\t {percentage:.5f} "
+
+    print(row)

pycode/memilio-simulation/CMakeLists.txt

@@ -110,6 +110,11 @@ add_pymio_module(_simulation_oseir
     SOURCES memilio/simulation/bindings/models/oseir.cpp
 )
 
+add_pymio_module(_simulation_oseir_metapop
+    LINKED_LIBRARIES memilio ode_seir_metapop
+    SOURCES memilio/simulation/bindings/models/oseir_metapop.cpp
+)
+
 add_pymio_module(_simulation_oseirdb
     LINKED_LIBRARIES memilio ode_seirdb
     SOURCES memilio/simulation/bindings/models/oseirdb.cpp

pycode/memilio-simulation/memilio/simulation/__init__.py

@@ -39,6 +39,9 @@ def __getattr__(attr):
     elif attr == "oseir":
         import memilio.simulation.oseir as oseir
         return oseir
+    elif attr == "oseir_metapop":
+        import memilio.simulation.oseir_metapop as oseir_metapop
+        return oseir_metapop
     elif attr == "oseirdb":
         import memilio.simulation.oseirdb as oseirdb
         return oseirdb

pycode/memilio-simulation/memilio/simulation/bindings/epidemiology/populations.h

@@ -24,6 +24,7 @@
 #include "utils/custom_index_array.h"
 #include "memilio/utils/custom_index_array.h"
 #include "memilio/epidemiology/populations.h"
+#include "geography/regions.h"
 
 #include "pybind11/pybind11.h"
 

pycode/memilio-simulation/memilio/simulation/bindings/geography/regions.h

@@ -0,0 +1,38 @@
+/* 
+* Copyright (C) 2020-2026 MEmilio
+*
+* Authors: Carlotta Gerstein
+*
+* Contact: Martin J. Kuehn <Martin.Kuehn@DLR.de>
+*
+* Licensed under the Apache License, Version 2.0 (the "License");
+* you may not use this file except in compliance with the License.
+* You may obtain a copy of the License at
+*
+*     http://www.apache.org/licenses/LICENSE-2.0
+*
+* Unless required by applicable law or agreed to in writing, software
+* distributed under the License is distributed on an "AS IS" BASIS,
+* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+* See the License for the specific language governing permissions and
+* limitations under the License.
+*/
+#ifndef PYMIO_REGION_H
+#define PYMIO_REGION_H
+
+#include "memilio/geography/regions.h"
+
+#include "pybind11/pybind11.h"
+
+namespace pymio
+{
+
+template <>
+inline std::string pretty_name<mio::regions::Region>()
+{
+    return "Region";
+}
+
+} // namespace pymio
+
+#endif //PYMIO_REGION_H

pycode/memilio-simulation/memilio/simulation/bindings/models/oseir_metapop.cpp

@@ -0,0 +1,108 @@
+/* 
+* Copyright (C) 2020-2025 MEmilio
+*
+* Authors: Henrik Zunker
+*
+* Contact: Martin J. Kuehn <Martin.Kuehn@DLR.de>
+*
+* Licensed under the Apache License, Version 2.0 (the "License");
+* you may not use this file except in compliance with the License.
+* You may obtain a copy of the License at
+*
+*     http://www.apache.org/licenses/LICENSE-2.0
+*
+* Unless required by applicable law or agreed to in writing, software
+* distributed under the License is distributed on an "AS IS" BASIS,
+* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+* See the License for the specific language governing permissions and
+* limitations under the License.
+*/
+
+// Includes from pymio
+#include "pybind_util.h"
+#include "compartments/simulation.h"
+#include "compartments/flow_simulation.h"
+#include "compartments/compartmental_model.h"
+#include "data/analyze_result.h"
+#include "epidemiology/age_group.h"
+#include "epidemiology/populations.h"
+#include "utils/custom_index_array.h"
+#include "utils/parameter_set.h"
+#include "utils/index.h"
+
+// Includes from MEmilio
+#include "models/ode_seir_metapop/model.h"
+#include "models/ode_seir_metapop/infection_state.h"
+#include "models/ode_seir_metapop/parameters.h"
+#include "memilio/compartments/simulation.h"
+#include "memilio/data/analyze_result.h"
+
+#include "pybind11/pybind11.h"
+#include "pybind11/eigen.h"
+
+namespace py = pybind11;
+
+namespace pymio
+{
+// specialization of pretty_name
+template <>
+inline std::string pretty_name<mio::oseirmetapop::InfectionState>()
+{
+    return "InfectionState";
+}
+
+} // namespace pymio
+
+PYBIND11_MODULE(_simulation_oseir_metapop, m)
+{
+    pymio::bind_interpolate_result_methods(m);
+
+    pymio::iterable_enum<mio::oseirmetapop::InfectionState>(m, "InfectionState", py::module_local{})
+        .value("Susceptible", mio::oseirmetapop::InfectionState::Susceptible)
+        .value("Exposed", mio::oseirmetapop::InfectionState::Exposed)
+        .value("Infected", mio::oseirmetapop::InfectionState::Infected)
+        .value("Recovered", mio::oseirmetapop::InfectionState::Recovered);
+
+    pymio::bind_ParameterSet<mio::oseirmetapop::ParametersBase<double>, pymio::EnablePickling::Required>(
+        m, "ParametersBase");
+
+    pymio::bind_CustomIndexArray<mio::UncertainValue<double>, mio::oseirmetapop::Region, mio::AgeGroup>(
+        m, "RegionAgeGroupArray");
+    using RegionAgeGroupPopulations = mio::Populations<double, mio::oseirmetapop::Region, mio::AgeGroup>;
+    pymio::bind_Population(m, "RegionAgeGroupPopulations", mio::Tag<RegionAgeGroupPopulations>{});
+
+    pymio::bind_class<mio::oseirmetapop::Parameters<double>, pymio::EnablePickling::Never,
+                      mio::oseirmetapop::ParametersBase<double>>(m, "Parameters")
+        .def(py::init<mio::oseirmetapop::Region, mio::AgeGroup>(), py::arg("num_regions"), py::arg("num_agegroups"))
+        .def("check_constraints", &mio::oseirmetapop::Parameters<double>::check_constraints)
+        .def("apply_constraints", &mio::oseirmetapop::Parameters<double>::apply_constraints)
+        .def_property_readonly("num_regions", &mio::oseirmetapop::Parameters<double>::get_num_regions)
+        .def_property_readonly("num_agegroups", &mio::oseirmetapop::Parameters<double>::get_num_agegroups);
+
+    using Populations =
+        mio::Populations<double, mio::oseirmetapop::Region, mio::AgeGroup, mio::oseirmetapop::InfectionState>;
+
+    pymio::bind_Population(m, "Populations", mio::Tag<mio::oseirmetapop::Model<double>::Populations>{});
+    pymio::bind_CompartmentalModel<mio::oseirmetapop::InfectionState, Populations,
+                                   mio::oseirmetapop::Parameters<double>, pymio::EnablePickling::Never>(m, "ModelBase");
+
+    using Model = mio::oseirmetapop::Model<double>;
+    pymio::bind_class<Model, pymio::EnablePickling::Never,
+                      mio::CompartmentalModel<double, mio::oseirmetapop::InfectionState, Populations,
+                                              mio::oseirmetapop::Parameters<double>>>(m, "Model")
+        .def(py::init<int, int>(), py::arg("num_regions"), py::arg("num_agegroups"))
+        .def("set_commuting_strengths", py::overload_cast<const Eigen::MatrixXd&>(&Model::set_commuting_strengths),
+             py::arg("commuting_strengths"))
+        .def("set_commuting_strengths_identity", py::overload_cast<>(&Model::set_commuting_strengths));
+
+    pymio::bind_Simulation<mio::Simulation<double, Model>>(m, "Simulation");
+    pymio::bind_Flow_Simulation<mio::FlowSimulation<double, Model>>(m, "FlowSimulation");
+
+    m.def("simulate", &mio::simulate<double, Model>, "Simulates an ODE SEIR metapopulation model from t0 to tmax.",
+          py::arg("t0"), py::arg("tmax"), py::arg("dt"), py::arg("model"), py::arg("integrator") = py::none());
+    m.def("simulate_flows", &mio::simulate_flows<double, Model>,
+          "Simulates an ODE SEIR metapopulation model with flows from t0 to tmax.", py::arg("t0"), py::arg("tmax"),
+          py::arg("dt"), py::arg("model"), py::arg("integrator") = py::none());
+
+    m.attr("__version__") = "dev";
+}

pycode/memilio-simulation/memilio/simulation/bindings/pybind_util.h

@@ -255,7 +255,7 @@ auto bind_Range(pybind11::module_& m, const std::string& class_name)
         .def(
             "__iter__",
             [](Range& self) {
-                return self;
+                return Iterator{{self.begin(), self.end()}};
             },
             pybind11::keep_alive<1, 0>{}) //keep alive the Range as long as there is an iterator
         .def(

pycode/memilio-simulation/memilio/simulation/bindings/simulation.cpp

@@ -59,6 +59,9 @@ PYBIND11_MODULE(_simulation, m)
     pymio::bind_CustomIndexArray<mio::UncertainValue<double>, mio::AgeGroup>(m, "AgeGroupArray");
     pymio::bind_class<mio::AgeGroup, pymio::EnablePickling::Required, mio::Index<mio::AgeGroup>>(m, "AgeGroup")
         .def(py::init<size_t>());
+    pymio::bind_Index<mio::regions::Region>(m, "IndexRegion");
+    pymio::bind_class<mio::regions::Region, pymio::EnablePickling::Never, mio::Index<mio::regions::Region>>(m, "Region")
+        .def(py::init<size_t>());
 
     pymio::bind_CustomIndexArray<double, mio::AgeGroup, mio::SimulationDay>(m, "AgeGroupSimulationDayArray");
     pymio::bind_class<mio::SimulationDay, pymio::EnablePickling::Never, mio::Index<mio::SimulationDay>>(m,

pycode/memilio-simulation/memilio/simulation/bindings/utils/index.h

@@ -51,6 +51,7 @@ void bind_Index(pybind11::module_& m, std::string const& name)
     c.def(pybind11::init<size_t>(), pybind11::arg("value"));
     c.def(pybind11::self == pybind11::self);
     c.def(pybind11::self != pybind11::self);
+    c.def("get", &mio::Index<Tag>::get);
 
     bind_Index_members_if_enum(c);
 }

pycode/memilio-simulation/memilio/simulation/oseir_metapop.py

@@ -0,0 +1,25 @@
+#############################################################################
+# Copyright (C) 2020-2025 MEmilio
+#
+# Authors: Henrik Zunker
+#
+# Contact: Martin J. Kuehn <Martin.Kuehn@DLR.de>
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+#     http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+#############################################################################
+
+"""
+Python bindings for the MEmilio ODE SEIR metapopulation model.
+"""
+
+from memilio.simulation._simulation_oseir_metapop import *