diff --git a/cpp/memilio/data/analyze_result.h b/cpp/memilio/data/analyze_result.h index 430cc0e28d..1cd9f37039 100644 --- a/cpp/memilio/data/analyze_result.h +++ b/cpp/memilio/data/analyze_result.h @@ -20,12 +20,12 @@ #ifndef MEMILIO_DATA_ANALYZE_RESULT_H #define MEMILIO_DATA_ANALYZE_RESULT_H +#include "memilio/config.h" #include "memilio/utils/time_series.h" #include "memilio/mobility/metapopulation_mobility_instant.h" #include "memilio/math/interpolation.h" #include "memilio/io/io.h" -#include #include namespace mio @@ -42,7 +42,8 @@ namespace mio * @return interpolated time series */ template -TimeSeries interpolate_simulation_result(const TimeSeries& simulation_result, const FP abs_tol = 1e-14); +TimeSeries interpolate_simulation_result(const TimeSeries& simulation_result, + const FP abs_tol = FP{100.} * Limits::zero_tolerance()); /** * @brief interpolate time series with freely chosen time points that lie in between the time points of the given time series up to a given tolerance. @@ -185,9 +186,9 @@ TimeSeries interpolate_simulation_result(const TimeSeries& simulation_re const auto t0 = simulation_result.get_time(0); const auto t_max = simulation_result.get_last_time(); // add another day if the first time point is equal to day_0 up to absolute tolerance tol - const auto day0 = (t0 - abs_tol < ceil(t0) - 1) ? floor(t0) : ceil(t0); + const auto day0 = (t0 - floor(t0) < abs_tol) ? floor(t0) : ceil(t0); // add another day if the last time point is equal to day_max up to absolute tolerance tol - const auto day_max = (t_max + abs_tol > floor(t_max) + 1) ? ceil(t_max) : floor(t_max); + const auto day_max = (ceil(t_max) - t_max < abs_tol) ? ceil(t_max) : floor(t_max); // create interpolation_times vector with all days between day0 and day_max std::vector tps(static_cast(day_max) - static_cast(day0) + 1); diff --git a/cpp/memilio/math/integrator.h b/cpp/memilio/math/integrator.h index 1d9aa83332..a5f6ad8a48 100755 --- a/cpp/memilio/math/integrator.h +++ b/cpp/memilio/math/integrator.h @@ -188,11 +188,12 @@ class SystemIntegrator using std::fabs; using std::max; using std::min; - const FP t0 = results.get_last_time(); + const FP t0 = results.get_last_time(); + const FP eps = Limits::zero_tolerance(); assert(tmax > t0); assert(dt > 0); - const size_t num_steps = + const auto num_steps = static_cast(ceil((tmax - t0) / dt)); // estimated number of time steps (if equidistant) results.reserve(results.get_num_time_points() + num_steps); @@ -204,7 +205,7 @@ class SystemIntegrator FP dt_min_restore = m_core->get_dt_min(); // used to restore dt_min, if it was decreased to reach tmax FP t = t0; - for (size_t i = results.get_num_time_points() - 1; fabs((tmax - t) / (tmax - t0)) > 1e-10; ++i) { + for (size_t i = results.get_num_time_points() - 1; fabs((tmax - t) / (tmax - t0)) > eps; ++i) { // We don't make time steps too small as the error estimator of an adaptive integrator //may not be able to handle it. this is very conservative and maybe unnecessary, //but also unlikely to happen. may need to be reevaluated. @@ -225,7 +226,7 @@ class SystemIntegrator results.get_last_time() = t; // if dt has been changed by step, register the current m_core as adaptive. - m_is_adaptive |= !floating_point_equal(dt, dt_copy, Limits::zero_tolerance()); + m_is_adaptive |= !floating_point_equal(dt, dt_copy, eps); } m_core->get_dt_min() = dt_min_restore; // restore dt_min // if dt was decreased to reach tmax in the last time iteration, @@ -236,7 +237,7 @@ class SystemIntegrator if (!step_okay) { log_warning("Adaptive step sizing failed. Forcing an integration step of size dt_min."); } - else if (fabs((tmax - t) / (tmax - t0)) > 1e-14) { + else if (fabs((tmax - t) / (tmax - t0)) > eps) { log_warning("Last time step too small. Could not reach tmax exactly."); } else { diff --git a/pycode/memilio-simulation/memilio/simulation/bindings/data/analyze_result.h b/pycode/memilio-simulation/memilio/simulation/bindings/data/analyze_result.h new file mode 100644 index 0000000000..ac5acb6ed7 --- /dev/null +++ b/pycode/memilio-simulation/memilio/simulation/bindings/data/analyze_result.h @@ -0,0 +1,60 @@ +/* +* Copyright (C) 2020-2026 MEmilio +* +* Authors: Rene Schmieding +* +* Contact: Martin J. Kuehn +* +* 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_DATA_ANALYZE_RESULT_H +#define PYMIO_DATA_ANALYZE_RESULT_H + +#include "memilio/data/analyze_result.h" +#include "pybind_util.h" + +#include "pybind11/pybind11.h" + +namespace py = pybind11; + +namespace pymio +{ + +/** + * @brief Bind functions interpolate_simulation_result and interpolate_ensemble_results. + */ +void bind_interpolate_result_methods(py::module_& m) +{ + + m.def( + "interpolate_simulation_result", + [](const mio::TimeSeries& ts) { + return mio::interpolate_simulation_result(ts); + }, + py::arg("ts")); + m.def("interpolate_simulation_result", + static_cast (*)(const mio::TimeSeries&, const double)>( + &mio::interpolate_simulation_result), + py::arg("ts"), py::arg("abs_tol")); + + m.def("interpolate_simulation_result", + static_cast (*)(const mio::TimeSeries&, const std::vector&)>( + &mio::interpolate_simulation_result), + py::arg("ts"), py::arg("interpolation_times")); + + m.def("interpolate_ensemble_results", &mio::interpolate_ensemble_results>); +} + +} // namespace pymio + +#endif // PYMIO_DATA_ANALYZE_RESULT_H diff --git a/pycode/memilio-simulation/memilio/simulation/bindings/models/omseirs4.cpp b/pycode/memilio-simulation/memilio/simulation/bindings/models/omseirs4.cpp index 13af660254..0cd77f8223 100644 --- a/pycode/memilio-simulation/memilio/simulation/bindings/models/omseirs4.cpp +++ b/pycode/memilio-simulation/memilio/simulation/bindings/models/omseirs4.cpp @@ -25,6 +25,7 @@ #include "epidemiology/populations.h" #include "utils/parameter_set.h" #include "utils/index.h" +#include "data/analyze_result.h" // Includes from MEmilio #include "ode_mseirs4/model.h" @@ -51,17 +52,7 @@ inline std::string pretty_name() PYBIND11_MODULE(_simulation_omseirs4, m) { // interpolation helpers - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const double)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("abs_tol") = 1e-14); - - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const std::vector&)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("interpolation_times")); - - m.def("interpolate_ensemble_results", &mio::interpolate_ensemble_results>); + pymio::bind_interpolate_result_methods(m); // InfectionState enum pymio::iterable_enum(m, "InfectionState") diff --git a/pycode/memilio-simulation/memilio/simulation/bindings/models/osecir.cpp b/pycode/memilio-simulation/memilio/simulation/bindings/models/osecir.cpp index 799477ba9d..46facf0d03 100644 --- a/pycode/memilio-simulation/memilio/simulation/bindings/models/osecir.cpp +++ b/pycode/memilio-simulation/memilio/simulation/bindings/models/osecir.cpp @@ -33,6 +33,7 @@ #include "mobility/metapopulation_mobility_instant.h" #include "io/mobility_io.h" #include "io/result_io.h" +#include "data/analyze_result.h" //Includes from MEmilio #include "ode_secir/model.h" @@ -87,17 +88,7 @@ PYBIND11_MAKE_OPAQUE(std::vector); PYBIND11_MODULE(_simulation_osecir, m) { // https://github.com/pybind/pybind11/issues/1153 - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const double)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("abs_tol") = 1e-14); - - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const std::vector&)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("interpolation_times")); - - m.def("interpolate_ensemble_results", &mio::interpolate_ensemble_results>); + pymio::bind_interpolate_result_methods(m); m.def("ensemble_mean", &mio::ensemble_mean); m.def("ensemble_percentile", &mio::ensemble_percentile); diff --git a/pycode/memilio-simulation/memilio/simulation/bindings/models/osecirvvs.cpp b/pycode/memilio-simulation/memilio/simulation/bindings/models/osecirvvs.cpp index 3ef274963e..288c95fd95 100755 --- a/pycode/memilio-simulation/memilio/simulation/bindings/models/osecirvvs.cpp +++ b/pycode/memilio-simulation/memilio/simulation/bindings/models/osecirvvs.cpp @@ -31,6 +31,7 @@ #include "epidemiology/populations.h" #include "io/mobility_io.h" #include "io/result_io.h" +#include "data/analyze_result.h" //Includes from MEmilio #include "ode_secirvvs/model.h" @@ -79,15 +80,7 @@ PYBIND11_MAKE_OPAQUE(std::vector); PYBIND11_MODULE(_simulation_osecirvvs, m) { - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const double)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("abs_tol") = 1e-14); - - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const std::vector&)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("interpolation_times")); + pymio::bind_interpolate_result_methods(m); pymio::iterable_enum(m, "InfectionState") .value("SusceptibleNaive", mio::osecirvvs::InfectionState::SusceptibleNaive) diff --git a/pycode/memilio-simulation/memilio/simulation/bindings/models/oseir.cpp b/pycode/memilio-simulation/memilio/simulation/bindings/models/oseir.cpp index d62c8c9493..42bd8edc4c 100755 --- a/pycode/memilio-simulation/memilio/simulation/bindings/models/oseir.cpp +++ b/pycode/memilio-simulation/memilio/simulation/bindings/models/oseir.cpp @@ -28,6 +28,7 @@ #include "compartments/compartmental_model.h" #include "epidemiology/age_group.h" #include "epidemiology/populations.h" +#include "data/analyze_result.h" //Includes from MEmilio #include "ode_seir/model.h" @@ -51,17 +52,7 @@ inline std::string pretty_name() PYBIND11_MODULE(_simulation_oseir, m) { - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const double)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("abs_tol") = 1e-14); - - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const std::vector&)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("interpolation_times")); - - m.def("interpolate_ensemble_results", &mio::interpolate_ensemble_results>); + pymio::bind_interpolate_result_methods(m); pymio::iterable_enum(m, "InfectionState") .value("Susceptible", mio::oseir::InfectionState::Susceptible) diff --git a/pycode/memilio-simulation/memilio/simulation/bindings/models/oseirdb.cpp b/pycode/memilio-simulation/memilio/simulation/bindings/models/oseirdb.cpp index 1a28476729..68c1a991a1 100644 --- a/pycode/memilio-simulation/memilio/simulation/bindings/models/oseirdb.cpp +++ b/pycode/memilio-simulation/memilio/simulation/bindings/models/oseirdb.cpp @@ -28,6 +28,7 @@ #include "compartments/compartmental_model.h" #include "epidemiology/age_group.h" #include "epidemiology/populations.h" +#include "data/analyze_result.h" //Includes from MEmilio #include "ode_seirdb/model.h" @@ -51,17 +52,7 @@ inline std::string pretty_name() PYBIND11_MODULE(_simulation_oseirdb, m) { - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const double)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("abs_tol") = 1e-14); - - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const std::vector&)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("interpolation_times")); - - m.def("interpolate_ensemble_results", &mio::interpolate_ensemble_results>); + pymio::bind_interpolate_result_methods(m); pymio::iterable_enum(m, "InfectionState") .value("Susceptible", mio::oseirdb::InfectionState::Susceptible) diff --git a/pycode/memilio-simulation/memilio/simulation/bindings/models/osir.cpp b/pycode/memilio-simulation/memilio/simulation/bindings/models/osir.cpp index 4f386e4eb0..6fe0cc20c7 100644 --- a/pycode/memilio-simulation/memilio/simulation/bindings/models/osir.cpp +++ b/pycode/memilio-simulation/memilio/simulation/bindings/models/osir.cpp @@ -28,6 +28,7 @@ #include "compartments/compartmental_model.h" #include "epidemiology/age_group.h" #include "epidemiology/populations.h" +#include "data/analyze_result.h" //Includes from MEmilio #include "ode_sir/model.h" @@ -51,17 +52,7 @@ inline std::string pretty_name() PYBIND11_MODULE(_simulation_osir, m) { - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const double)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("abs_tol") = 1e-14); - - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const std::vector&)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("interpolation_times")); - - m.def("interpolate_ensemble_results", &mio::interpolate_ensemble_results>); + pymio::bind_interpolate_result_methods(m); pymio::iterable_enum(m, "InfectionState") .value("Susceptible", mio::osir::InfectionState::Susceptible) diff --git a/pycode/memilio-simulation/memilio/simulation/bindings/models/ssir.cpp b/pycode/memilio-simulation/memilio/simulation/bindings/models/ssir.cpp index 9503153e46..56d510ccd0 100644 --- a/pycode/memilio-simulation/memilio/simulation/bindings/models/ssir.cpp +++ b/pycode/memilio-simulation/memilio/simulation/bindings/models/ssir.cpp @@ -25,6 +25,7 @@ #include "utils/custom_index_array.h" #include "utils/parameter_set.h" #include "epidemiology/populations.h" +#include "data/analyze_result.h" //Includes from MEmilio #include "sde_sir/model.h" @@ -52,17 +53,7 @@ inline std::string pretty_name() PYBIND11_MODULE(_simulation_ssir, m) { - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const double)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("abs_tol") = 1e-14); - - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const std::vector&)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("interpolation_times")); - - m.def("interpolate_ensemble_results", &mio::interpolate_ensemble_results>); + pymio::bind_interpolate_result_methods(m); pymio::iterable_enum(m, "InfectionState") .value("Susceptible", mio::ssir::InfectionState::Susceptible) diff --git a/pycode/memilio-simulation/memilio/simulation/bindings/models/ssirs.cpp b/pycode/memilio-simulation/memilio/simulation/bindings/models/ssirs.cpp index 77904208ef..e85d0372ed 100644 --- a/pycode/memilio-simulation/memilio/simulation/bindings/models/ssirs.cpp +++ b/pycode/memilio-simulation/memilio/simulation/bindings/models/ssirs.cpp @@ -25,6 +25,7 @@ #include "utils/custom_index_array.h" #include "utils/parameter_set.h" #include "epidemiology/populations.h" +#include "data/analyze_result.h" //Includes from MEmilio #include "sde_sirs/model.h" @@ -52,17 +53,7 @@ inline std::string pretty_name() PYBIND11_MODULE(_simulation_ssirs, m) { - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const double)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("abs_tol") = 1e-14); - - m.def("interpolate_simulation_result", - static_cast (*)(const mio::TimeSeries&, const std::vector&)>( - &mio::interpolate_simulation_result), - py::arg("ts"), py::arg("interpolation_times")); - - m.def("interpolate_ensemble_results", &mio::interpolate_ensemble_results>); + pymio::bind_interpolate_result_methods(m); pymio::iterable_enum(m, "InfectionState") .value("Susceptible", mio::ssirs::InfectionState::Susceptible)