diff --git a/.github/pull_request_template.md b/.github/pull_request_template.md index e1c24f0b2a..a53c48c1b2 100644 --- a/.github/pull_request_template.md +++ b/.github/pull_request_template.md @@ -22,7 +22,8 @@ Please check our [git workflow](https://memilio.readthedocs.io/en/latest/develop - [ ] New code adheres to [coding guidelines](https://memilio.readthedocs.io/en/latest/development.html#coding-guidelines) - [ ] No large data files have been added (files should in sum not exceed 100 KB, avoid PDFs, Word docs, etc.) - [ ] Tests are added for new functionality and a local test run was successful (with and without OpenMP) -- [ ] Appropriate **documentation** for new functionality has been added (Doxygen in the code and explanations in the online documentation) +- [ ] Appropriate **documentation within the code** (Doxygen) for new functionality has been added in the code +- [ ] Appropriate **external documentation** (ReadTheDocs) for new functionality has been added to the online documentation - [ ] Proper attention to licenses, especially no new third-party software with conflicting license has been added - [ ] (For ABM development) Checked [benchmark results](https://memilio.readthedocs.io/en/latest/development.html#agent-based-model-development) and ran and posted a local test above from before and after development to ensure performance is monitored. diff --git a/README.md b/README.md index 5333ab438a..f009c4b64a 100644 --- a/README.md +++ b/README.md @@ -9,17 +9,17 @@ MEmilio implements various models for infectious disease dynamics, from simple c If you use MEmilio, please cite our work -- Bicker J, Kerkmann D, Korf S, Plötzke L, Schmieding R, Wendler A, Zunker H et al. (2025) *MEmilio - a High Performance Modular Epidemics Simulation Software*. Available at `https://github.com/SciCompMod/memilio` and `https://elib.dlr.de/213614/`. +- Bicker J, Kerkmann D, Korf S, Plötzke L, Schmieding R, Wendler A, Zunker H et al. (2025) *MEmilio - a High Performance Modular Epidemics Simulation Software*. Available at https://github.com/SciCompMod/memilio and https://elib.dlr.de/213614/ . and, in particular, for -- Ordinary differential equation-based (ODE) and Graph-ODE models: Zunker H, Schmieding R, Kerkmann D, Schengen A, Diexer S, et al. (2024). *Novel travel time aware metapopulation models and multi-layer waning immunity for late-phase epidemic and endemic scenarios*. *PLOS Computational Biology* 20(12): e1012630. `https://doi.org/10.1371/journal.pcbi.1012630` -- Integro-differential equation-based (IDE) models: Wendler AC, Plötzke L, Tritzschak H, Kühn MJ. (2024). *A nonstandard numerical scheme for a novel SECIR integro differential equation-based model with nonexponentially distributed stay times*. Submitted for publication. `https://arxiv.org/abs/2408.12228` -- Agent-based models (ABMs): Kerkmann D, Korf S, Nguyen K, Abele D, Schengen A, et al. (2025). *Agent-based modeling for realistic reproduction of human mobility and contact behavior to evaluate test and isolation strategies in epidemic infectious disease spread*. *Computers in Biology and Medicine* 193: 110269. `DOI:10.1016/j.compbiomed.2025.110269 `_ -- Hybrid agent-metapopulation-based models: Bicker J, Schmieding R, Meyer-Hermann M, Kühn MJ. (2025). *Hybrid metapopulation agent-based epidemiological models for efficient insight on the individual scale: A contribution to green computing*. *Infectious Disease Modelling* 10(2): 571-590. `https://doi.org/10.1016/j.idm.2024.12.015` -- Graph Neural Networks: Schmidt A, Zunker H, Heinlein A, Kühn MJ. (2024). *Towards Graph Neural Network Surrogates Leveraging Mechanistic Expert Knowledge for Pandemic Response*. arXiv. `https://arxiv.org/abs/2411.06500` -- ODE-based models with Linear Chain Trick: Plötzke L, Wendler A, Schmieding R, Kühn MJ. (2024). *Revisiting the Linear Chain Trick in epidemiological models: Implications of underlying assumptions for numerical solutions*. Submitted for publication. `https://doi.org/10.48550/arXiv.2412.09140` -- Behavior-based ODE models: Zunker H, Dönges P, Lenz P, Contreras S, Kühn MJ. (2025). *Risk-mediated dynamic regulation of effective contacts de-synchronizes outbreaks in metapopulation epidemic models*. Chaos, Solitons & Fractals. `https://doi.org/10.1016/j.chaos.2025.116782` +- Ordinary differential equation-based (ODE) and Graph-ODE models: Zunker H, Schmieding R, Kerkmann D, Schengen A, Diexer S, et al. (2024). *Novel travel time aware metapopulation models and multi-layer waning immunity for late-phase epidemic and endemic scenarios*. *PLOS Computational Biology* 20(12): e1012630. https://doi.org/10.1371/journal.pcbi.1012630 +- Integro-differential equation-based (IDE) models: Wendler A, Plötzke L, Tritzschak H, Kühn MJ. (2026). *A nonstandard numerical scheme for a novel SECIR integro differential equation-based model with nonexponentially distributed stay times*. *Applied Mathematics and Computation* 509: 129636. https://doi.org/10.1016/j.amc.2025.129636 +- Agent-based models (ABMs): Kerkmann D, Korf S, Nguyen K, Abele D, Schengen A, et al. (2025). *Agent-based modeling for realistic reproduction of human mobility and contact behavior to evaluate test and isolation strategies in epidemic infectious disease spread*. *Computers in Biology and Medicine* 193: 110269. https://doi.org/10.1016/j.compbiomed.2025.110269 +- Hybrid agent-metapopulation-based models: Bicker J, Schmieding R, Meyer-Hermann M, Kühn MJ. (2025). *Hybrid metapopulation agent-based epidemiological models for efficient insight on the individual scale: A contribution to green computing*. *Infectious Disease Modelling* 10(2): 571-590. https://doi.org/10.1016/j.idm.2024.12.015 +- Graph Neural Networks: Schmidt A, Zunker H, Heinlein A, Kühn MJ. (2025). *Graph Neural Network Surrogates to leverage Mechanistic Expert Knowledge towards Reliable and Immediate Pandemic Response*. Submitted for publication. https://doi.org/10.48550/arXiv.2411.06500 +- ODE-based models with Linear Chain Trick: Plötzke L, Wendler A, Schmieding R, Kühn MJ. (2024). *Revisiting the Linear Chain Trick in epidemiological models: Implications of underlying assumptions for numerical solutions*. Submitted for publication. https://doi.org/10.48550/arXiv.2412.09140 +- Behavior-based ODE models: Zunker H, Dönges P, Lenz P, Contreras S, Kühn MJ. (2025). *Risk-mediated dynamic regulation of effective contacts de-synchronizes outbreaks in metapopulation epidemic models*. Chaos, Solitons & Fractals. https://doi.org/10.1016/j.chaos.2025.116782 **Getting started** @@ -28,6 +28,12 @@ and, in particular, for The documentation for MEmilio can be found at https://memilio.readthedocs.io/en/latest/index.html +**Publication simulations** + +Simulations used for publications, along with their specific plotting and processing scripts, +are available in a separate repository: +https://github.com/SciCompMod/memilio-simulations + **Development** The coding guidelines and git workflow description can be found in the documentation at diff --git a/cpp/examples/smm.cpp b/cpp/examples/smm.cpp index cb63adf39b..1e7f5d285e 100644 --- a/cpp/examples/smm.cpp +++ b/cpp/examples/smm.cpp @@ -19,6 +19,7 @@ */ #include "smm/simulation.h" +#include "smm/model.h" #include "smm/parameters.h" #include "memilio/data/analyze_result.h" #include "memilio/epidemiology/adoption_rate.h" diff --git a/docs/source/_static/custom.css b/docs/source/_static/custom.css index fc86d8fd05..d6198044dd 100644 --- a/docs/source/_static/custom.css +++ b/docs/source/_static/custom.css @@ -25,11 +25,13 @@ dl.class > dt:first-of-type { .rst-content dl .sig { width: -webkit-fill-available; } + .rst-content .viewcode-link { display: inline-flex; float: inline-end; margin-right: 1.5em; } + .rst-content .headerlink { position: absolute; right: 0.5em; @@ -39,3 +41,48 @@ dl.class > dt:first-of-type { .wy-table-responsive table th { white-space: normal; } + +/* Team page styles */ +.developer-photo { + border-radius: 8px; + margin-right: 20px; + margin-bottom: 10px; + /* box-shadow: 0 2px 4px rgba(0, 0, 0, 0.1); */ +} + +/* Clear floats after each developer section */ +.section h3 { + clear: both; +} + +/* Add some spacing between developer sections */ +.section>.section { + margin-bottom: 40px; + padding-bottom: 20px; + border-bottom: 1px solid #e1e4e5; +} + +/* Style the links section */ +.section ul { + list-style-type: none; + padding-left: 0; +} + +.section ul li { + margin-bottom: 5px; +} + +.section ul li a { + color: #2980b9; + text-decoration: none; + font-weight: 500; +} + +.section ul li a:hover { + text-decoration: underline; +} + +.section p { + text-align: justify; + line-height: 1.6; +} diff --git a/docs/source/_static/funder/bmbf.jpg b/docs/source/_static/funder/bmbf.jpg new file mode 100644 index 0000000000..e91c985687 Binary files /dev/null and b/docs/source/_static/funder/bmbf.jpg differ diff --git a/docs/source/_static/funder/bmdv.png b/docs/source/_static/funder/bmdv.png new file mode 100644 index 0000000000..42dbba3894 Binary files /dev/null and b/docs/source/_static/funder/bmdv.png differ diff --git a/docs/source/_static/funder/hdslee.png b/docs/source/_static/funder/hdslee.png new file mode 100644 index 0000000000..08b012c3c3 Binary files /dev/null and b/docs/source/_static/funder/hdslee.png differ diff --git a/docs/source/_static/funder/helmholtz.jpg b/docs/source/_static/funder/helmholtz.jpg new file mode 100644 index 0000000000..60bd2e120b Binary files /dev/null and b/docs/source/_static/funder/helmholtz.jpg differ diff --git a/docs/source/_static/funder/mfund.png b/docs/source/_static/funder/mfund.png new file mode 100644 index 0000000000..e0913145a3 Binary files /dev/null and b/docs/source/_static/funder/mfund.png differ diff --git a/docs/source/_static/team/max_mustermann.jpg b/docs/source/_static/team/max_mustermann.jpg new file mode 100644 index 0000000000..3c56a98195 Binary files /dev/null and b/docs/source/_static/team/max_mustermann.jpg differ diff --git a/docs/source/citation.rst b/docs/source/citation.rst index a8e0a26779..674774dffe 100644 --- a/docs/source/citation.rst +++ b/docs/source/citation.rst @@ -3,14 +3,15 @@ Citing MEmilio If you use MEmilio, please cite our work -- Kühn, Martin Joachim et al. (2025). *MEmilio - a High Performance Modular Epidemics Simulation Software (2022)*. Available at `https://github.com/SciCompMod/memilio `_ and `https://elib.dlr.de/213614/ `_. +- Bicker J, Kerkmann D, Korf S, Ploetzke L, Schmieding R, Wendler A, Zunker H et al. (2025). *MEmilio - a High Performance Modular Epidemics Simulation Software*. Available at `https://github.com/SciCompMod/memilio `_ and `https://elib.dlr.de/213614/ `_. and, in particular, for - **Ordinary differential equation-based (ODE) and Graph-ODE models**: Zunker H, Schmieding R, Kerkmann D, Schengen A, Diexer S, et al. (2024). *Novel travel time aware metapopulation models and multi-layer waning immunity for late-phase epidemic and endemic scenarios*. *PLOS Computational Biology* 20(12): e1012630. `DOI:10.1371/journal.pcbi.1012630 `_ -- **Integro-differential equation-based (IDE) models**: Wendler AC, Plötzke L, Tritzschak H, Kühn MJ. (2024). *A nonstandard numerical scheme for a novel SECIR integro differential equation-based model with nonexponentially distributed stay times*. Submitted for publication. `arXiv:2408.12228 `_ +- **Integro-differential equation-based (IDE) models**: Wendler A, Plötzke L, Tritzschak H, Kühn MJ. (2026). *A nonstandard numerical scheme for a novel SECIR integro differential equation-based model with nonexponentially distributed stay times*. *Applied Mathematics and Computation* 509: 129636. `DOI:10.1016/j.amc.2025.129636 `_ - **Agent-based models (ABMs)**: Kerkmann D, Korf S, Nguyen K, Abele D, Schengen A, et al. (2025). *Agent-based modeling for realistic reproduction of human mobility and contact behavior to evaluate test and isolation strategies in epidemic infectious disease spread*. *Computers in Biology and Medicine* 193: 110269. `DOI:10.1016/j.compbiomed.2025.110269 `_ - **Hybrid agent-metapopulation-based models**: Bicker J, Schmieding R, Meyer-Hermann M, Kühn MJ. (2025). *Hybrid metapopulation agent-based epidemiological models for efficient insight on the individual scale: A contribution to green computing*. *Infectious Disease Modelling* 10(2): 571-590. `DOI:10.1016/j.idm.2024.12.015 `_ -- **Graph Neural Networks**: Schmidt A, Zunker H, Heinlein A, Kühn MJ. (2024).*Towards Graph Neural Network Surrogates Leveraging Mechanistic Expert Knowledge for Pandemic Response*. arXiv. `arXiv:2411.06500 `_ -- **ODE-based models with Linear Chain Trick**: Plötzke L, Wendler A, Schmieding R, Kühn MJ. (2024). *Revisiting the Linear Chain Trick in epidemiological models: Implications of underlying assumptions for numerical solutions*. Submitted for publication. `DOI:10.48550/arXiv.2412.09140 `_ -- **Behavior-based ODE models**: Zunker H, Dönges P, Lenz P, Contreras S, Kühn MJ. (2025). *Risk-mediated dynamic regulation of effective contacts de-synchronizes outbreaks in metapopulation epidemic models*. arXiv. `arXiv:2502.14428 `_ +- **Graph Neural Networks**: Schmidt A, Zunker H, Heinlein A, Kühn MJ. (2025). *Graph Neural Network Surrogates to leverage Mechanistic Expert Knowledge towards Reliable and Immediate Pandemic Response*. Submitted for publication. `arXiv:2411.06500 `_ +- **ODE-based models with Linear Chain Trick**: Plötzke L, Wendler A, Schmieding R, Kühn MJ. (2024). *Revisiting the Linear Chain Trick in epidemiological models: Implications of underlying assumptions for numerical solutions*. Accepted for publication. `DOI:10.48550/arXiv.2412.09140 `_ +- **Behavior-based ODE models**: Zunker H, Dönges P, Lenz P, Contreras S, Kühn MJ. (2025). *Risk-mediated dynamic regulation of effective contacts de-synchronizes outbreaks in metapopulation epidemic models*. *Chaos, Solitons & Fractals* 199:116782. `DOI:10.1016/j.chaos.2025.116782 `_ + diff --git a/docs/source/cpp/aggregated_models.rst b/docs/source/cpp/aggregated_models.rst index abfecc855a..5ab44e9b07 100644 --- a/docs/source/cpp/aggregated_models.rst +++ b/docs/source/cpp/aggregated_models.rst @@ -6,8 +6,8 @@ There are different equation-based models implemented in MEmilio that consider a .. toctree:: :maxdepth: 1 - Ordinary differential equations - Linear Chain Trick - Generalized Linear Chain Trick - Stochastic-differential equations - Integro-differential equations \ No newline at end of file + Ordinary differential equations (ODE) + Linear Chain Trick (LCT) + Generalized Linear Chain Trick (GLCT) + Integro-differential equations (IDE) + Stochastic-differential equations (SDE) \ No newline at end of file diff --git a/docs/source/cpp/cpp_makros_preprocessing.rst b/docs/source/cpp/cpp_makros_preprocessing.rst new file mode 100644 index 0000000000..15c34247c0 --- /dev/null +++ b/docs/source/cpp/cpp_makros_preprocessing.rst @@ -0,0 +1,104 @@ +Configuration macros +------------------- + +MEmilio provides several compile-time configuration options through CMake that control which features and dependencies are available during compilation. These configurations are defined as preprocessor macros and can be used to conditionally compile code sections based on available libraries and enabled features. +The following macros are automatically defined by CMake during the build process based on the availability of dependencies and build options: + +.. list-table:: + :header-rows: 1 + :widths: 30 70 + + * - Macro Name + - Description + * - :code:`MEMILIO_HAS_HDF5` + - Defined when HDF5 library is available. Enables HDF5 file I/O functionality for reading and writing results. + * - :code:`MEMILIO_HAS_JSONCPP` + - Defined when JsonCpp library is available. Enables JSON file I/O functionality for parameter import and configuration files. + * - :code:`MEMILIO_ENABLE_MPI` + - Defined when MPI (Message Passing Interface) support is enabled. Allows distributed computing and parallel simulations across multiple processes or compute nodes. + * - :code:`MEMILIO_ENABLE_OPENMP` + - Defined when OpenMP support is enabled. Enables shared-memory parallelization for multi-threaded execution within a single process. + * - :code:`MEMILIO_ENABLE_PROFILING` + - Defined when profiling support is enabled. Activates performance monitoring and timing instrumentation throughout the codebase. + +Usage in Code +^^^^^^^^^^^^^^ + +These macros are used with preprocessor conditionals to enable or disable specific functionality. You should use these macros to ensure that your code can compile and run correctly depending on the available libraries and features. +When using these macros, you should always consider fallbacks or alternative implementations when a feature is not available. This ensures that your code remains robust and can handle cases where optional dependencies are not present. +Here are common usage patterns: + +**Conditional Compilation for Optional Dependencies** + +.. code-block:: cpp + + #ifdef MEMILIO_HAS_JSONCPP + #include "memilio/io/epi_data.h" + #include "memilio/io/result_io.h" + + // JSON-based parameter I/O functions + template + IOResult read_divi_data(const std::string& path, + const std::vector& vregion, + Date date, + std::vector& vnum_icu) { + } + #endif // MEMILIO_HAS_JSONCPP + +**HDF5-Specific Functionality** + +.. code-block:: cpp + + #ifdef MEMILIO_HAS_HDF5 + template + IOResult export_input_data_county_timeseries( + std::vector models, + const std::string& results_dir, + /* ... other parameters ... */) { + // HDF5-based time series export + } + #else + template + IOResult export_input_data_county_timeseries( + std::vector models, + const std::string& results_dir, + /* ... other parameters ... */) { + return failure(StatusCode::UnknownError, + "HDF5 not available"); + } + #endif // MEMILIO_HAS_HDF5 + +**Parallel Computing** + +.. code-block:: cpp + + #ifdef MEMILIO_ENABLE_OPENMP + #include + + void parallel_simulation() { + #pragma omp parallel for + for (int i = 0; i < num_regions; ++i) { + // Parallel execution of regional simulations + } + } + #endif // MEMILIO_ENABLE_OPENMP + +**Check Multiple Features Simultaneously** + +.. code-block:: cpp + + // Check if both JSON and HDF5 are available + #if defined(MEMILIO_HAS_JSONCPP) && defined(MEMILIO_HAS_HDF5) + IOResult read_and_export_data() { + // Implementation using both JSON input and HDF5 output + } + #elif defined(MEMILIO_HAS_JSONCPP) + IOResult read_and_export_data() { + // JSON-only implementation + } + #else + IOResult read_and_export_data() { + return failure(StatusCode::UnknownError, + "Neither JSON nor HDF5 available"); + } + #endif diff --git a/docs/source/cpp/data_types.rst b/docs/source/cpp/data_types.rst new file mode 100644 index 0000000000..006397a460 --- /dev/null +++ b/docs/source/cpp/data_types.rst @@ -0,0 +1,51 @@ +Common data types +----------------- + +The following list explains the nonstandard data types that are used throughout MEmilio. + +.. list-table:: + :header-rows: 1 + :widths: 20 60 + + * - Data type name + - Description + * - :code:`FP` + - A floating point type. Usually :code:`double` is used, but for instance in the optimization using optimal control :code:`FP` is equal to :code:`Ipopt::Number`, see :doc:`models/oseair` and `examples/ode_seair_optimization.cpp `_. + * - :code:`UncertainValue` + - This data type describes a value sampled from a given distribution. The value is initialized with a given :code:`FP` and can be (re)sampled with the :code:`draw_sample()` function. + * - :code:`AgeGroup` + - A typesafe ``size_t``, i.e. an integer that cannot be confused with other integer types, so operations like assignment, addition etc. only work with other :code:`AgeGroup`\s. Derived from :code:`mio::Index`. + * - :code:`Region` + - A typesafe ``size_t``, derived from :code:`mio::Index`. + * - :code:`CustomIndexArray` + - A contiguous array whose values can be accessed by a multi-index, e.g. a :code:`mio::Index` with one or more categories. This datatype is, for example, used in the parameter :code:`mio::abm::TimeExposedToNoSymptoms` making it dependent on :code:`mio::abm::VirusVariant` and :code:`mio::AgeGroup`. Its values can then be set for a specific :code:`virus_variant` and :code:`age_group` using :code:`model.parameters.template get()[{virus_variant, age_group}]`. + * - :code:`Populations` + - Is a :code:`mio::CustomIndexArray` with :code:`mio::UncertainValue` as values. Adds some convenient functions like :code:`get_group_total`. + * - :code:`TimeSeries` + - Stores vectors of values at time points. Each time point has a vector of values of the same size with operations like adding time points, retrieving values, exporting to CSV, etc. It's also used for storing and analyzing simulation results over time. + * - :code:`Graph` + - A generic graph structure that represents a network of nodes connected by edges. Each node and edge can have associated properties. The graph is used to model geographical regions connected by mobility patterns (e.g., commuting), where each node is represented by its own epidemiological model. + * - :code:`Node` + - Represents a node in a graph with a unique ID and associated properties. + * - :code:`Edge` + - Represents a directed connection between two nodes in a graph with associated properties. + * - :code:`EdgeBase`, :code:`InEdgeBase`, :code:`OutEdgeBase` + - Base classes for Edge that define start and end node indices for connections in the graph. + * - :code:`SimulationNode` + - Represents a simulation in one node of a mobility graph. Contains a simulation model of any type and keeps track of the last state and time point. + * - :code:`MobilityCoefficients` + - Time-dependent mobility coefficients used to model how populations move between nodes in a graph. + * - :code:`MobilityCoefficientGroup` + - A collection of time-dependent mobility coefficients that differentiate between various sources of mobility. + * - :code:`MobilityParameters` + - Parameters that influence mobility between nodes, including coefficients and dynamic nonpharmaceutical interventions (NPIs). + * - :code:`MobilityEdge` + - Represents mobility between two nodes in a graph. Handles the movement of populations between nodes, tracks mobile populations, and applies mobility returns according to epidemiological models. + * - :code:`ContactMatrix` + - Time-dependent contact frequencies between groups, derived from ``DampingMatrixExpression``. Models how the contact rates between different age groups change over time due to interventions. + * - :code:`ContactMatrixGroup` + - A collection of contact matrices that represent different contexts (e.g., home, school, work) whose sum is the total number of contacts, derived from ``DampingMatrixExpressionGroup``. + * - :code:`DampingMatrixExpression` + - Represents a coefficient-wise matrix expression :math:`B - D \odot (B - M)`, where :math:`B` is a baseline matrix, :math:`M` is a minimum matrix, :math:`D` is a time-dependent complex damping factor, and :math:`\odot` is element wise multiplication. Used as the base for time-dependent contact matrices. + * - :code:`DampingMatrixExpressionGroup` + - Represents a collection of ``DampingMatrixExpression``\s that are summed up. Used for representing multiple sources of contacts or mobility. diff --git a/docs/source/cpp/development.rst b/docs/source/cpp/development.rst index ff19249dae..b1db5347f1 100644 --- a/docs/source/cpp/development.rst +++ b/docs/source/cpp/development.rst @@ -1,96 +1,11 @@ Development -================================= +=========== -.. _performance-monitoring-cpp: +Here we collect some information that is relevant for the development processes of MEmilio. -Performance monitoring ---------------------------------- +.. toctree:: + :maxdepth: 2 -LIKWID -~~~~~~ - -To measure the performance of the code, we use LIKWID. We recommend measuring the performance on an HPC system with LIKWID installed or refer to ``_ on how to get LIKWID. -Run ``likwid-perfctr ./bin/my_example`` to measure the performance of the entire example. To measure the performance of a specific part, use the LIKWID marker API: - -.. code-block:: cpp - - #include - - // ... - - LIKWID_MARKER_INIT; - - LIKWID_MARKER_START("my_marker"); - // code to be measured - LIKWID_MARKER_STOP("my_marker"); - - LIKWID_MARKER_CLOSE; - -Set the CMake variable ``MEMILIO_USE_LIKWID=ON`` to enable LIKWID support and run ``likwid-perfctr -m ./bin/my_example``. -For more details see the LIKWID documentation, available `here `_. - -Agent-based model benchmarks -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -There is a suite of benchmarks for the ABM that can be used to check performance. The suite contains setups of different sizes. If you added a new feature (i.e., you didn't just fix a bug in an existing feature), make sure the feature is actually used by the benchmark. Add it to the benchmark if necessary, then run the benchmark to see if the cost for the new feature is acceptable and as expected. - -Most new features will add some overhead, but this needs to be limited and in proportion to the added value of the feature so runtime doesn't grow out of control. Optional features that can be disabled should only incur minimal overhead. If you did not add any new feature, just run the benchmark before and after your changes to make sure there are no performance regressions. This process will hopefully be automated soon by running benchmarks in the CI. - -Build the benchmarks by defining the CMake variable ``MEMILIO_BUILD_BENCHMARKS=ON`` in the build. Make sure to use a **Release** build to test performance. - -.. code-block:: bash - - cmake .. -DMEMILIO_BUILD_BENCHMARKS=ON -DCMAKE_BUILD_TYPE=Release - cmake --build . - -Run the benchmark executable: - -.. code-block:: bash - - ./build/bin/abm_benchmark - -Each benchmark is run for a number of iterations and the average time is reported. - -.. code-block:: text - - Benchmark Time CPU Iterations - --------------------------------------------------------------------------- - abm_benchmark/abm_benchmark_50k 7583 ms 7583 ms 1 - abm_benchmark/abm_benchmark_100k 18216 ms 18214 ms 1 - abm_benchmark/abm_benchmark_200k 41492 ms 41489 ms 1 - -You may get a warning: - -.. code-block:: text - - ***WARNING*** CPU scaling is enabled, the benchmark real time measurements may be noisy and will incur extra overhead. - -If possible, disable CPU scaling to improve the consistency of results. See the Google Benchmark documentation here: -https://google.github.io/benchmark/reducing_variance.html - -Also, try to reduce other system load during the benchmark run. - -If it is not possible to disable frequency scaling, increase the runtime of the benchmark using the commands below. Constant CPU frequency is necessary to get the most reliable results and to measure small differences. - -**REMINDER:** Don't forget to re-enable CPU scaling after you ran the benchmarks to save energy. Rebooting may restore the settings as well. - -The benchmark executable has a number of command line arguments that customize execution. Use ``--help`` to list them all. - -Two important options for consistency and stability: - -- ``--benchmark_min_time=``: Iterate each benchmark so that the total runtime is at least ``T`` seconds. - Default is 1 second, which may not be enough. - Try 60 seconds for better stability (you may need to experiment). - -- ``--benchmark_repetitions=``: Repeat each benchmark ``N`` times and report mean, median, and variance. - (Repetitions are **not** iterations; a benchmark can be repeated 10 times with 5 iterations each. Each repetition runs for at least the minimum time.) - -``benchmark_repetitions`` is useful to check timing consistency, as it reports variance. -However, it can be expensive because long-running benchmarks are repeated. -``benchmark_min_time`` increases iterations only for fast-running benchmarks, which tend to be less stable. - -**Suggested workflow:** - -1. Run with 5–10 repetitions to check variance. -2. Increase ``benchmark_min_time`` until variance is acceptable. -3. Continue benchmarking with 1 repetition and the adjusted minimum time. \ No newline at end of file + performance_monitoring + data_types + cpp_makros_preprocessing diff --git a/docs/source/cpp/diffusive_abm.rst b/docs/source/cpp/diffusive_abm.rst index f099ea05f9..eec7238630 100644 --- a/docs/source/cpp/diffusive_abm.rst +++ b/docs/source/cpp/diffusive_abm.rst @@ -10,6 +10,9 @@ with :math:`X` the vector of all agents' positions and :math:`Z` the vector of a The Model class needs an Implementation class as template argument which provides the domain agents move and interact in. A quadwell potential and a singlewell potential given in the classes **QuadWell** and **SingleWell** respectively are implemented, but any other suitable potential can be used as implementation. +In the following, we present more details of the diffusive agent-based model, including code examples. +An overview of nonstandard but often used data types can be found under :doc:`data_types`. + Infection states ---------------- @@ -35,7 +38,7 @@ Using the infection states Susceptible (S), Exposed (E), Carrier (C), Infected ( Infection state transitions --------------------------- -The infection state transitions are explicitly given by the adoption rates and are therefore subject to user input. Adoption rates always depend on their source infection state. If an adoption event requires interaction of agents (e.g. disease transmission), the corresponding rate depends not only on the source infection state, but also on other infection states, the **Influence**s. An adoption rate that only depends on the source infection state, e.g. recovery or worsening of disease symptoms, is called `first-order` adoption rate and an adoption rate that has influences is called `second-order` adoption rate. Adoption rates are region-dependent; therefore it is possible to have different rates in two regions for the same state transition which can be useful when having e.g. region-dependent interventions or contact behavior. +The infection state transitions are explicitly given by the adoption rates and are therefore subject to user input. Adoption rates always depend on their source infection state. If an adoption event requires interaction of agents (e.g. disease transmission), the corresponding rate depends not only on the source infection state, but also on other infection states, the **Influence**\s. An adoption rate that only depends on the source infection state, e.g. recovery or worsening of disease symptoms, is called `first-order` adoption rate and an adoption rate that has influences is called `second-order` adoption rate. Adoption rates are region-dependent; therefore it is possible to have different rates in two regions for the same state transition which can be useful when having e.g. region-dependent interventions or contact behavior. Using the infection states from above and only one region, there are five first-order and one second-order adoption rate that can be set via: @@ -122,7 +125,7 @@ There are no non-pharmaceutical interventions (NPIs) explicitly implemented in t .. code-block:: cpp //Reduce the transmission risk by 10% - model.get_adoption_rates().at({mio::mpm::Region(0), Status::S, Status::E}).factor *= 0.9; + model.get_adoption_rates().at({mio::regions::Region(0), Status::S, Status::E}).factor *= 0.9; Simulation ----------- diff --git a/docs/source/cpp/glct.rst b/docs/source/cpp/glct.rst index e0d1639e19..39ef4abf2c 100644 --- a/docs/source/cpp/glct.rst +++ b/docs/source/cpp/glct.rst @@ -1,18 +1,17 @@ Generalized Linear Chain Trick model ==================================== -MEmilio implements a SECIR-type model utilizing the Generalized Linear Chain Trick (GLCT). This is a generalization of -the LCT model and allows for phase-type distributed stay times in the compartments. Phase-type distributions are dense -in the field of all positive-valued distributions. Therefore, for any positive-valued distribution, a phase-type -distribution of arbitrary precision can be identified. Note that the resulting system can still be described by ODEs. +MEmilio implements a SECIR-type model utilizing the Generalized Linear Chain Trick (GLCT). This is a generalization of the LCT model. In contrast to simpler ODE models that assume (possibly unrealistic) exponentially distributed stay times, the GLCT allows for more realistic, phase-type distributed stay times in the compartments through the use of subcompartments. Phase-type distributions are dense in the field of all positive-valued distributions. Therefore, for any positive-valued distribution, a phase-type +distribution of arbitrary precision can be identified. Note that the resulting system can still be described by an ODE-system. -In the following, we present the general structure of the GLCT model. The particular model documentation with examples -is linked at the bottom of this page. +In the following, we present the general structure of the GLCT model. The particular model documentation with examples is linked at the bottom of this page. + +An overview of nonstandard but often used data types can be found under :doc:`data_types`. Infection states ---------------- -The model contains a list of **InfectionState**\s that define particular features of the subpopulations in the particular state. +The model contains a list of **InfectionState**s that define particular features of the subpopulations in the particular state. .. code-block:: RST @@ -27,7 +26,7 @@ To make use of the GLCT, we additionally need to define the numbers of subcompar `Number of subcompartments of State2` `...` -Note that the model is implemented as **CompartmentalModel**. +The model is implemented as **CompartmentalModel**. Parameters @@ -36,38 +35,30 @@ Parameters The parameters of the model are defined as structs and are combined in a class ``ParameterSet``. We use different types of parameters to represent epidemiological parameters such as the mean stay times in a compartment or the contact rates between different age groups. Most model parameters are constants that describe -pathogen-specific characteristics (possibly resolved by sociodemographic groups) and are represented by a vector with a -value for each sociodemographic group. To model different contact rates between different sociodemographic groups, we -use a parameter denoted **ContactPatterns** of type **UncertainContactMatrix**. The **UncertainContactMatrix** contains -a set of contact matrices of arbitrary length and which can represent the different contact locations in the model like -schools, workplaces, or homes. The matrices can be loaded or stored in the particular example. +pathogen-specific characteristics (possibly resolved by sociodemographic groups) and are represented by a vector with +a value for each sociodemographic group. To model different contact rates between different sociodemographic groups, +we use a parameter denoted **ContactPatterns** of type **UncertainContactMatrix**. The **UncertainContactMatrix** contains an +arbitrary large set of contact matrices which can e.g. represent the different contact locations in the model +like schools, workplaces, or homes. The matrices can be loaded or stored in the particular example. In the **ContactPatterns**, each matrix element stores baseline contact rates :math:`c_{i,j}` between sociodemographic -group :math:`i` to group :math:`j`. The dimension of the matrix is automatically defined by the model initiation and it is reduced -to one value if no stratifcation is used. The values can be adjusted during the simulation, e.g., through implementing -nonpharmaceutical interventions, see the section on :ref:`Nonpharmaceutical Interventions`. -Parameters can get accessed via ``model.parameters.get>()`` and set via either -``model.parameters.get>() = value`` or ``model.parameters.set>(value)``. +group :math:`i` to group :math:`j`. The dimension of the matrix is automatically defined by the model initialization +and it is reduced to one value if no stratification is used. The values can be adjusted during the simulation, e.g., +through implementing nonpharmaceutical interventions, see the section on :ref:`Nonpharmaceutical Interventions GLCT`. +Parameters can be accessed via ``model.parameters.get>()`` and set via either +``model.parameters.set>(value)`` or ``model.parameters.get>() = value``. Initial conditions ------------------ -The initial conditions of the model are represented by a class **LctPopulations** that gives the number of individuals in -each sociodemographic group and each subcompartment of each **InfectionState**. For more details, see -:doc:`Model Creation `. Before the simulation, set the initial conditions for each **InfectionState** and -sociodemographic group. +The initial conditions of the model are represented by a class **LctPopulations** that gives the number of individuals in each sociodemographic group and each subcompartment of each **InfectionState**. For more details, see :doc:`Model Creation `. Before the simulation, the initial conditions for each **InfectionState** and sociodemographic group must be set. -.. _Nonpharmaceutical Interventions: +.. _Nonpharmaceutical Interventions GLCT: Nonpharmaceutical interventions ------------------------------- -Contact rates can be adjusted during the simulation to model nonpharmaceutical interventions (NPIs) such as lockdowns, -school closures, or social distancing. This is done by adding **Damping**\s to the **ContactPatterns** of the model. A -**Damping** is defined by a time point at which the intervention starts and a matrix of the same size as the -**ContactMatrix**. While in many cases, the reduction matrix is given by a constant matrix with factor :math:`r`, also -group-specific reductions are possible through setting particular rows or columns differently. With a constant -reduction factor :math:`r`, the reduced contact rate is :math:`(1-r) * c_{i,j}`. +Contact rates can be adjusted during the simulation to model nonpharmaceutical interventions (NPIs) such as lockdowns, school closures, or social distancing. This is done by adding **Damping**\s to the **ContactPatterns** of the model. A **Damping** is defined by a time point at which the intervention starts and a matrix of the same size as the **ContactMatrix**. While in many cases, the reduction matrix is given by a constant matrix with factor :math:`r`, also group-specific reductions are possible through setting particular rows or columns differently. With a constant reduction factor :math:`r`, the reduced contact rate is :math:`(1-r) * c_{i,j}`. .. dropdown:: :fa:`gears` Expert's settings @@ -86,47 +77,28 @@ reduction factor :math:`r`, the reduced contact rate is :math:`(1-r) * c_{i,j}`. Simulation ---------- -Once the model is setup, run a simple simulation from time ``t0`` to ``tmax`` with an initial step size ``dt`` using the -``mio::simulation()`` function. This will run a simulation of type **Simulation** that saves the sizes of each -subcompartment over time. -You can run a simulation using either fixed or adaptive solution schemes with an absolute or relative tolerance. By -default, the simulation uses an adaptive solution scheme of the boost library and an absolute tolerance of 1e-10 and a -relative tolerance of 1e-5. For more details on the possible integration schemes, see . +Once the model is setup, one can run a simple simulation from time ``t0`` to ``tmax`` with initial step size ``dt`` using the ``mio::simulate()`` function. This will run a simulation of type **Simulation** that saves the sizes of each subcompartment over time. +You can run a simulation using either fixed or adaptive integration schemes with an absolute or relative tolerance. By default, the simulation uses an adaptive solution scheme of the boost library with an absolute tolerance of 1e-10 and a relative tolerance of 1e-5. For more details on the possible integration schemes, see . Output ------ -The output of the **Simulation** is a **TimeSeries** ``result`` containing the sizes of each subcompartment at each time point. +The output of the **Simulation** is a ``mio::TimeSeries`` containing the sizes of each subcompartment at each time point. To obtain a result with respect to the compartments, the subcompartments can be accumulated via the function ``calculate_compartments()``. A simple table can be printed using the ``print_table()`` function of the -``TimeSeries`` class. The compartment sizes can be printed with ``model.calculate_compartments(result).print_table()``. -As adaptive step size methods are used by default, the output will not be available on equidistant time points like `dt` -or days. To obtain outputs on days or user-defined time points, you can interpolate the results to days or -any other series of times points with ``mio::interpolate_simulation_result()``. +``mio::TimeSeries`` class. The compartment sizes can be printed with ``model.calculate_compartments(result).print_table()``. +As adaptive step size methods are used by default, the output will not be available on equidistant time points like `dt` or days. To obtain outputs on days or user-defined time points, you can interpolate the results to days or any other series of times points with ``mio::interpolate_simulation_result()``. Visualization ------------- -To visualize the results of a simulation, you can use the Python package :doc:`memilio_plot <../python/memilio_plot>` -and its documentation. - -In the following, we give detailed explanations of the GLCT-SECIR model. - -Introduction -------------- - -This model is based on the Generalized Linear Chain Trick (GLCT). - -The GLCT provides the option to use phase-type distributed stay times in the compartments through the use of subcompartments. The Generalized Linear Chain Trick is an extension of the Linear Chain Trick (as the name already suggests). Phase-type distributions are dense in the field of all positive-valued distributions. Therefore, for any positive-valued distribution, a phase-type distribution of arbitrary precision can be identified. -The normal ODE models have (possibly unrealistic) exponentially distributed stay times. -The GLCT model can still be described by an ordinary differential equation system. - +To visualize the results of a simulation, you can use the Python package :doc:`m-plot <../python/m-plot>` and its documentation. List of models ------------------------------ +-------------- .. toctree:: :titlesonly: diff --git a/docs/source/cpp/graph_abm.rst b/docs/source/cpp/graph_abm.rst index 99edb6d657..1d30a044d2 100644 --- a/docs/source/cpp/graph_abm.rst +++ b/docs/source/cpp/graph_abm.rst @@ -4,7 +4,18 @@ Graph-based agent-based model Introduction ------------- -This model realizes multiple instances of the mobility-based agent-based model ``abm::Model`` (see :doc:`mobility_based_abm` for documentation) as nodes in a directed graph. One local model represents a geographical region. The regions are connected by the graph edges. Mobility within one node and via the graph edges follows the same mobility rules that can be handed as argument to ``mio::GraphABModel``. Therefore this graph-based agent-based model can be reduced to a single mobility-based agent-based if simulation time steps within the whole graph, i.e. the step size of each node and the step size of the edge exchange, are equal. +This model realizes multiple instances of the mobility-based agent-based model (ABM) ``abm::Model`` (see :doc:`mobility_based_abm` for the documentation) +as nodes in a directed graph. One local model represents a geographical region. The regions are connected by the graph edges. Mobility within one node +and via the graph edges follows the same mobility rules that can be handed as argument to ``mio::GraphABModel``. +Therefore this graph-based agent-based (graph-ABM) model can be reduced to a single mobility-based agent-based model if +simulation time steps within the whole graph, i.e. the step size of each node and the step size of the edge exchange, are equal. +Preimplemented mobility rules can be found in ``_. +The motivation behind the graph-ABM is to have multiple ABMs run independently from each other in parallel for different regions and only synchronize, +i.e. exchange agents via edges, in fixed time intervals. The synchronization time steps should be bigger than the internal +ABM time steps to reduce communication between nodes as much as possible. + +An overview of nonstandard but often used data types can be found under :doc:`data_types`. + Simulation ----------- @@ -122,7 +133,9 @@ Assigning infection states and locations to persons in all models can be done vi adult2.set_assigned_location(mio::abm::LocationType::Work, work, model2.get_id()); adult3.set_assigned_location(mio::abm::LocationType::Work, work, model2.get_id()); -For initializing the graph nodes and edges a ``mio::Graph`` is created which gets ``mio::ABMSimulationNode`` and ``mio::ABMMobilityEdge`` as templates. Additionally, every node needs a ``mio::History`` object to log its results during the simulation. See :doc:`io` for information on how to use ``mio::History``. Below, ``mio::abm::LogInfectionState`` is used as logger. +For initializing the graph nodes and edges, a ``mio::Graph`` is created which gets ``mio::ABMSimulationNode`` and ``mio::ABMMobilityEdge`` as templates. +Additionally, every node needs a ``mio::History`` object to log its results during the simulation. See :ref:`history` for information on how to use ``mio::History``. +Below, ``mio::abm::LogInfectionState`` is used as logger. .. code-block:: cpp @@ -135,7 +148,8 @@ For initializing the graph nodes and edges a ``mio::Graph`` is created which get graph.add_edge(model1.get_id(), model2.get_id()); graph.add_edge(model2.get_id(), model1.get_id()); -To simulate the model from `start_date` to `end_date` with given graph step size `exchange_time_span`, a GraphSimulation has to be created. The step size is used to regularly exchange agents via the graph edges. Advancing the simulation until `end_date` is done as follows: +To simulate the model from `start_date` to `end_date` with given graph step size `exchange_time_span`, a GraphSimulation has to be created. +The step size is used to regularly exchange agents via the graph edges. Advancing the simulation until `end_date` is done as follows: .. code-block:: cpp diff --git a/docs/source/cpp/graph_metapop.rst b/docs/source/cpp/graph_metapop.rst index eaa57fa1c3..a7df8d2275 100644 --- a/docs/source/cpp/graph_metapop.rst +++ b/docs/source/cpp/graph_metapop.rst @@ -4,7 +4,16 @@ Graph-based metapopulation model Introduction ------------- -This model realizes a metapopulation model using a directed graph with instances of an ODE-based model as graph nodes. Any ODE-based model implemented in MEmilio can be used in the graph nodes, see :doc:`ode`, representing one geographical region. Currently, only ODE-based models are supported as nodes, but support for other models will be added soon. Mobility between graph nodes is modelled via graph edges. There are two different types of mobility edges implemented in MEmilio: ``mio::MobilityEdge`` and ``mio::MobilityEdgeStochastic``. The first one implements a deterministic instant mobility approach where subpopulations are exchanged between graph nodes instantaneously in fixed time intervals. The latter implements mobility as stochastic jumps between graph nodes using a temporal Gillespie algorithm for simulation. +This model realizes a metapopulation model using a directed graph with instances of an ODE-based model as graph nodes. +Any ODE-based model implemented in MEmilio can be used in the graph nodes, see :doc:`ode`, representing one geographical region. +Currently, only ODE-based models are supported as nodes, but support for other models will be added soon. +Mobility between graph nodes is modelled via graph edges. There are two different types of mobility edges implemented +in MEmilio: ``mio::MobilityEdge`` and ``mio::MobilityEdgeStochastic``. The first one implements a deterministic instant +mobility approach where subpopulations are exchanged between graph nodes instantaneously in fixed time intervals. +The latter implements mobility as stochastic jumps between graph nodes using a temporal Gillespie algorithm for simulation. + +An overview of nonstandard but often used data types can be found under :doc:`data_types`. + Simulation ----------- @@ -15,7 +24,7 @@ During the simulation, graph nodes are advanced independently from each other ac The instant mobility approach was introduced by Kühn et al. (2021) and mainly represents daily instant commuting activities between different regions. Using realistic commuting data, such as the absolute number of commuters between regions, we can model realistic exchanges of individuals in the metapopulation context. -In this approach, the commuters are exchanged twice daily in an instantaneous manner, allowing the integration of age-specific mobility restrictions. The implementation also allows the restriction of mobility activities for certain infection states, such as symptomatic infected individuals. An important feature implemented in Kühn et al. (2022) is the testing of commuters, which enables the detection and potential isolation of asymptomatic or symptomatic infected individuals before they can infect people in other regions. +In this approach, the commuters are exchanged twice daily in an instantaneous manner, allowing the integration of age-specific mobility restrictions. The implementation also allows for the restriction of mobility activities for certain infection states, such as symptomatic infected individuals. An important feature implemented in Kühn et al. (2022) is the testing of commuters, which enables the detection and potential isolation of asymptomatic or symptomatic infected individuals before they can infect people in other regions. The instant mobility approach assumes instant exchange between regions without considering travel times. The length of stay is fixed to half a day. Since all regions exchange commuters simultaneously, this scheme offers great properties for parallelization and has fewer constraints on the integration time step compared to more detailed mobility approaches. @@ -29,9 +38,9 @@ For further details, please refer to: The detailed mobility scheme was introduced in Zunker et al. (2024) and further develops the instant approach. The detailed mobility model addresses some of the constraints of the instant approach by adding realistic travel and stay times. -The core idea of this detailed approach is that each region has a second, local model, which can be parameterized independently from the existing model, to represent mobility within the region. This allows for more detailed modeling of population movement. +The core idea of this detailed approach is that each region has a second local model, which can be parameterized independently from the existing model, to represent mobility within the region. This allows for more detailed modeling of population movement. -To have a realistic representation of the exchange process, we calculate the centroids of each region, which is represented by a polygon. When modeling travel between two non-adjacent regions (nodes k and l), the model creates a path connecting their centroids, which may pass through other regions. The predefined travel time between regions is then distributed across all regions along this path. +To have a realistic representation of the exchange process, we calculate the centroids of each region, which is represented by a polygon. When modeling travel between two non-adjacent regions (nodes :math:`k` and :math:`l`), the model creates a path connecting their centroids, which may pass through other regions. The predefined travel time between regions is then distributed across all regions along this path. Commuters move along these paths and pass through various mobility models during their trip. This allows for potential interactions with other commuters during transit. The detailed mobility approach enables modeling of scenarios that cannot be addressed with the instant approach, such as the impact of different mask type usage in public transport. However, this increased realism comes at the cost of greater computational effort and requires more detailed input data. For further details, please refer to: @@ -40,7 +49,7 @@ For further details, please refer to: **Stochastic mobility approach** -In the stochastic mobility approach, transitions of individuals between regions, i.e. graph nodes, are modeled as stochastic jumps. The frequency of individuals transitioning from region i to region j is determined by a rate :math:`\lambda_{ij}` which can be interpreted as the (expected) fraction of the population in region i that commutes to region j within one day. In contrast to the instant and detailed mobility approach, the time points and the number of transitions between two regions are stochastic. +In the stochastic mobility approach, transitions of individuals between regions, i.e. graph nodes, are modeled as stochastic jumps. The frequency of individuals transitioning from region :math:`i` to region :math:`j` is determined by a rate :math:`\lambda_{ij}` which can be interpreted as the (expected) fraction of the population in region :math:`i` that commutes to region :math:`j` within one day. In contrast to the instant and detailed mobility approach, the time points and the number of transitions between two regions are stochastic. Behavior-based ODE models @@ -57,19 +66,18 @@ The extension and inclusion of regional and global information is based on the f How to: Set up a graph and run a graph simulation ------------------------------------------------- -The graph simulation couples multiple simulation instances representing different geographic regions via a graph-based approach. In this example, a compartment model based on ODEs (ODE-SECIR) is used for each region. The simulation proceeds by advancing each region independently and then exchanging portions of the population between regions along the graph edges. +The graph simulation couples multiple simulation instances representing different geographic regions via a graph-based approach. In this example, a compartmental model based on ODEs (ODE-SECIR) is used for each region. The simulation proceeds by advancing each region independently and then exchanging portions of the population between regions along the graph edges. The following steps detail how to configure and execute a graph simulation: -1. **Initialize the local compartment model:** +1. **Initialize the local compartmental model:** - First, set up the compartment model by initializing the populations and parameters. In this example, we use a single age group and start with 10,000 susceptible individuals. Additionally, the necessary epidemiological parameters (e.g. time periods, transmission probabilities) are set. + First, set up the compartmental model by initializing the parameters that are equal in every region. In this example, we use a single age group and set the necessary epidemiological parameters (e.g. time periods, transmission probabilities). .. code-block:: cpp const size_t num_groups = 1; mio::osecir::Model model(num_groups); - model.populations[{mio::AgeGroup(0), mio::osecir::InfectionState::Susceptible}] = 10000; model.parameters.set(0); model.parameters.set>(0.2); @@ -94,7 +102,7 @@ The following steps detail how to configure and execute a graph simulation: 2. **Create simulation groups and adjust contact patterns:** - To represent different geographic regions, clone the base model into separate model groups. In this example, two model groups are created. The first group is modified by applying a contact damping to simulate contact restrictions. + To represent different geographic regions, clone the base model into separate model groups. In this example, two model groups are created. The first group is modified by applying a contact damping to simulate contact restrictions. Additionally, set the populations in all models. .. code-block:: cpp @@ -111,6 +119,8 @@ The following steps detail how to configure and execute a graph simulation: model_group1.populations[{mio::AgeGroup(0), mio::osecir::InfectionState::Susceptible}] = 9990; model_group1.populations[{mio::AgeGroup(0), mio::osecir::InfectionState::Exposed}] = 100; + model_group2.populations[{mio::AgeGroup(0), mio::osecir::InfectionState::Susceptible}] = 10000; + 3. **Define compartments to save from edges:** To extract the mobility results, define the compartments to save from the edges. In this example, the compartments for infected individuals with and without symptoms are saved for each region. @@ -147,7 +157,7 @@ The following steps detail how to configure and execute a graph simulation: g.add_edge(1, 0, Eigen::VectorXd::Constant((size_t)mio::osecir::InfectionState::Count, 0.1), indices_save_edges); -For the stochastic mobility, ``mio::MobilityEdgeStochastic`` has to be used as edge type for the graph. The rates or mibility coefficients can be set as follows: + For the stochastic mobility, ``mio::MobilityEdgeStochastic`` has to be used as edge type for the graph. The rates or mobility coefficients can be set as follows: .. code-block:: cpp @@ -172,7 +182,7 @@ For the stochastic mobility, ``mio::MobilityEdgeStochastic`` has to be used as e .. code-block:: cpp const auto tmax = 30.; - const auto dt = 0.5; // time step or Mobility (daily mobility occurs every second step) + const auto dt = 0.5; // time step for Mobility (daily mobility occurs every second step) auto sim = mio::make_mobility_sim(t0, dt, std::move(g)); sim.advance(tmax); @@ -185,6 +195,4 @@ For the stochastic mobility, ``mio::MobilityEdgeStochastic`` has to be used as e auto& edge_1_0 = sim.get_graph().edges()[1]; auto& results = edge_1_0.property.get_mobility_results(); results.print_table({"Commuter INS", "Commuter ISy", "Commuter Total"}); - - return 0; - } + diff --git a/docs/source/cpp/hybrid_models.rst b/docs/source/cpp/hybrid_models.rst index 2778b7df3c..c957a00b61 100644 --- a/docs/source/cpp/hybrid_models.rst +++ b/docs/source/cpp/hybrid_models.rst @@ -1,7 +1,7 @@ Hybrid models ============== -Hybrid models combine two modeling approaches to balance the trade-off that arises when choosing an appropriate model for a given application. The motivation behind combining two modeling approaches is that multiple models have different advantages and disadvantages, making it difficult to identify a single best option. In MEmilio a temporal hybridization is implemented where the concept is to switch between two models according to a predefined condition during the course of the simulation. +Hybrid models combine two modeling approaches to balance the trade-off that arises when choosing an appropriate model for a given application. The motivation behind combining two modeling approaches is that multiple models have different advantages and disadvantages, making it difficult to identify a single best option. In MEmilio, a temporal hybridization is implemented where the concept is to switch between two models according to a predefined condition during the course of the simulation. .. toctree:: :maxdepth: 1 diff --git a/docs/source/cpp/ide.rst b/docs/source/cpp/ide.rst index b5dc11ea9b..e09b0d06d0 100644 --- a/docs/source/cpp/ide.rst +++ b/docs/source/cpp/ide.rst @@ -9,7 +9,10 @@ In the following, we present the general structure of the IDE-SECIR model. The IDE-SEIR model follows a different structure, an introduction and a detailed example are provided in a separate section below. -Infection States +An overview of nonstandard but often used data types can be found under :doc:`data_types`. + + +Infection states ---------------- Every model contains a list of **InfectionState**\s that define particular features of the subpopulations in the particular state. @@ -21,7 +24,7 @@ Every model contains a list of **InfectionState**\s that define particular featu `...` -Infection State Transitions +Infection state transitions --------------------------- Additionally, we define **InfectionTransition**\s that define the possible transitions between the **InfectionState**\s. @@ -34,7 +37,7 @@ When solving the model, we obtain results for the **InfectionTransition**\s as w `...` -Sociodemographic Stratification +Sociodemographic stratification ------------------------------- For the IDE-SECIR model, the population can also be stratified by one sociodemographic dimension. This dimension is denoted @@ -45,17 +48,17 @@ Parameters The parameters of the model are defined as structs and are combined in a class ``ParameterSet``. We use different types of parameters to represent epidemiological parameters such as the distributions of the stay times in a -compartment or the contact rates between different age groups. Most model parameters are constants that describe +compartment or the contact rates between different sociodemographic groups. Most model parameters are constants that describe pathogen-specific characteristics (possibly resolved by sociodemographic groups) and are represented by a vector with a -value for each sociodemographic group. To model different contact rates between different sociodemographic groups, we -use a parameter denoted **ContactPatterns** of type **UncertainContactMatrix**. The **UncertainContactMatrix** contains -a set of contact matrices of arbitrary length and which can represent the different contact locations in the model like +value for each group. To model different contact rates between different sociodemographic groups, we +use a parameter denoted **ContactPatterns** of type **UncertainContactMatrix**. The **UncertainContactMatrix** contains an +arbitrary large set of contact matrices which can represent the different contact locations in the model like schools, workplaces, or homes. The matrices can be loaded or stored in the particular example. In the **ContactPatterns**, each matrix element stores baseline contact rates :math:`c_{i,j}` between sociodemographic -group :math:`i` to group :math:`j`. The dimension of the matrix is automatically defined by the model initiation and it is reduced -to one value if no stratifcation is used. The values can be adjusted during the simulation, e.g., through implementing -nonpharmaceutical interventions, see the section on :ref:`Nonpharmaceutical Interventions`. +group :math:`i` and group :math:`j`. The dimension of the matrix is automatically defined by the model initiation and it is reduced +to one value if no stratification is used. The values can be adjusted during the simulation, e.g., through implementing +nonpharmaceutical interventions, see the section on :ref:`Nonpharmaceutical Interventions IDE`. An important feature of our IDE-based model is that we can choose the transition distributions in a flexible way. The default distribution is a smoother cosine function as it provides good testing qualities. For more realistic simulations, @@ -65,11 +68,10 @@ e.g. ``SmootherCosine``. Any class that is derived from ``StateAgeFunction`` can ``StateAgeFunctionWrapper`` object that is then passed to the model. Parameters can get accessed via ``model.parameters.get>()`` and set via either -``model.parameters.get>() = value``. - +``model.parameters.get>() = value`` or ``model.parameters.set>(value)``. -Initial Conditions +Initial conditions ------------------ The initial conditions consist of a **TimeSeries** that contains the flows per time step for some time interval before @@ -79,8 +81,8 @@ the documentation of **StateAgeFunction**. Note that the last time point of the start time of the simulation. -.. _Nonpharmaceutical Interventions: -Nonpharmaceutical Interventions +.. _Nonpharmaceutical Interventions IDE: +Nonpharmaceutical interventions ------------------------------- Contact rates can be adjusted during the simulation to model nonpharmaceutical interventions (NPIs) such as lockdowns, @@ -126,80 +128,30 @@ the results to days or any other series of times points with ``mio::interpolate_ Visualization ------------- -To visualize the results of a simulation, you can use the Python package :doc:`memilio_plot <../python/memilio_plot>` +To visualize the results of a simulation, you can use the Python package :doc:`m-plot <../python/m-plot>` and its documentation. -List of models ------------------------ - -.. toctree:: - :titlesonly: - - models/iseir - models/isecir - - IDE-SEIR model ---------------- - -Introduction -~~~~~~~~~~~~~ -The four compartments - -- `Susceptible` (:math:`S`), may become exposed at any time -- `Exposed` (:math:`E`), becomes infected after some time -- `Infected` (:math:`I`), will recover after some time -- `Recovered` (:math:`R`) - -are used to simulate the spread of the disease. - -Simulation -~~~~~~~~~~~ - -The simulation runs in discrete time steps using a trapezoidal rule. The model and the numerical scheme is based on the paper `"Modeling infectious diseases using integro-differential equations: Optimal -control strategies for policy decisions and Applications in COVID-19" by Keimer and Pflug, 2020 `_. +-------------- -For a detailed description and application of the model, see: +The IDE-SEIR roughly follows the same structure but differs in a few points: -Plötzke L (2021) Modellierung epidemischer Infektionskrankheiten auf der Basis von gewöhnlichen und Integro-Differentialgleichungen. Bachelor thesis, University of Cologne. https://elib.dlr.de/143504/ +- We do not consider infection state transitions but only infection states. +- The initial conditions are defined by the compartment sizes at the simulation start :math:`t_0`. -How to: Set up and run a simulation of the IDE-SEIR model -~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ +Note that the IDE-SEIR model is solved with a different numerical scheme than the IDE-SECIR model. For further information, +see :doc:`IDE-SEIR <./models/iseir>`. +List of models +-------------- -Here, we set the contact matrix used in the simulation. One can define multiple matrices for different locations. The size of each of these matrices is defined by the number of age groups. -Below, we use only one contact matrix for one location. As we only consider one age group in our example, we set the corresponding contact rate to :math:`10`. - -.. code-block:: cpp - - mio::ContactMatrixGroup contact_matrix = mio::ContactMatrixGroup(1, 1); - contact_matrix[0] = mio::ContactMatrix(Eigen::MatrixXd::Constant(1, 1, 10.)); - -To simulate the implementation of nonpharmaceutical interventions, we add dampings to the contact rate. Here, we apply a damping of :math:`0.7` after :math:`10` days, meaning that the contact rate is reduced to :math:`30%` of the initial value. - -.. code-block:: cpp - - contact_matrix[0].add_damping(0.7, mio::SimulationTime(10.)); - model.parameters.get>() = mio::UncertainContactMatrix(contact_matrix); - -After defining :math:`t_{\max}`, we can simulate, which means that we calculate the value for the compartment :math:`S`. - -.. code-block:: cpp - - int tmax = 15; - model.simulate(tmax); - -The values of the remaining compartments :math:`E`, :math:`I` and :math:`R` are calculated using the parameters ``LatencyTime`` and ``InfectiousTime`` and obtain a time series containing the values of all compartments. - -.. code-block:: cpp - - auto result = model.calculate_EIR(); - -Finally, we can print our results. +.. toctree:: + :titlesonly: + + models/isecir + models/iseir -.. code-block:: cpp - result.print_table({"S", "E", "I", "R"}); diff --git a/docs/source/cpp/ide_creation.rst b/docs/source/cpp/ide_creation.rst index ae4166f024..53314a59b6 100644 --- a/docs/source/cpp/ide_creation.rst +++ b/docs/source/cpp/ide_creation.rst @@ -1,4 +1,4 @@ -IDE Models +IDE models ========== Currently, the implemented solver for integro-differential equations is tailored to the SECIR-type model. From a theoretical point of view, this type of solver can also be applied to other epidemic IDE-based models. Implementing this into our software is work in progress. \ No newline at end of file diff --git a/docs/source/cpp/installation.rst b/docs/source/cpp/installation.rst index 97a9eb4699..9e0bfebea6 100644 --- a/docs/source/cpp/installation.rst +++ b/docs/source/cpp/installation.rst @@ -1,4 +1,246 @@ +Build instructions +================== + +The MEmilio core library (MEmilio C++) is written in C++ and uses `CMake `_ as build system. Before +installing MEmilio C++, make sure a C++20 compiler, CMake and a build tool (like GNU Make or Ninja) is installed on your +device. The following guide will make use of the command line, but you can use graphical build tools from an IDE as +well. + +Quick start +----------- + +These are minimal build instructions. More detailed steps and explanations are given below. + +.. code:: bash + + git clone https://github.com/SciCompMod/memilio # download the project + cd memilio # go into the project directory + cmake -S cpp -B cpp/build # *configure* the project, creating a "build" directory under cpp/ + cmake --build cpp/build -j 2 # *build* all targets from the project with 2 threads + +After the build process is done, you can run files from "cpp/build/bin", for example our test suite: + +.. code:: bash + + ./cpp/build/bin/memilio-test + +This will run several tests and should write out ``[ PASSED ]`` in the end. + +Requirements +------------ +MEmilio C++ is regularly tested with the following compilers (list will be extended over time): + +- GCC, versions 11 and 13 +- Clang, version 14 and 17 +- MSVC, version 19.43 (Visual Studio 2022) + +MEmilio C++ is regularly tested on GitHub runners using Ubuntu 22.04 and 24.04 and Windows Server 2022 and 2025. It is +expected to run on any comparable Linux or Windows system. It is currently not tested on macOS. + +The following table lists the dependencies that are used. Most of them are required, but some are optional. The library +can be used without them but with slightly reduced features. CMake will warn about them during configuration. Most of +them are bundled with this library and do not need to be installed manually. Bundled libraries are either included with +this project or loaded from the web on demand. For each dependency, there is a CMake option to use an installed version +instead. Version compatibility needs to be ensured by the user, the version we currently use is included in the table. + +.. list-table:: + :header-rows: 1 + + * - Library + - Version + - Required + - Bundled + - Notes + * - spdlog + - 1.15.0 + - Yes + - Yes (git repo) + - https://github.com/gabime/spdlog + * - Eigen + - 3.4.0 + - Yes + - Yes (git repo) + - http://gitlab.com/libeigen/eigen + * - Boost + - 1.84.0 + - Yes + - Yes (git repo) + - https://github.com/boostorg/boost + * - JsonCpp + - 1.9.6 + - No + - Yes (git repo) + - https://github.com/open-source-parsers/jsoncpp + * - HDF5 + - 1.12.0 + - No + - No + - https://www.hdfgroup.org/, package libhdf5-dev on apt (Ubuntu) + * - GoogleTest + - 1.10 + - For Tests only + - Yes (git repo) + - https://github.com/google/googletest + * - LibSBML + - 5.20.2 + - No + - No + - https://sbml.org/software/libsbml/ (For SBML integration only) + +See the `thirdparty directory `_ for more details. + +Step-by-step instructions +------------------------- + +Download +~~~~~~~~ + +Start by downloading the newest version (or a specific release) of MEmilio from our +`github repository `_ to a directory of your choice, or clone it directly using +git by first opening a terminal in that directory and then running + +.. code:: bash + + git clone https://github.com/SciCompMod/memilio + +.. dropdown:: :fa:`gears` Note for developers + + If you need to push changes to the main repo, register an ssh key with GitHub and clone from + *git@github.com:SciCompMod/memilio.git* instead. You can also change to the ssh address later using + + .. code:: bash + + git remote set-url origin git@github.com:SciCompMod/memilio.git + +This will create a new directory called "memilio". Change into this directory. + +.. code:: bash + + cd memilio + +Configuration +~~~~~~~~~~~~~ + +Before we can *build* anything, we need to *configure* the project first. If you want to use its default options, +simply run + +.. code:: bash + + cmake -S cpp -B cpp/build + +Additional options can be specified by appending one or more ``-D