949 Correct handling of dt_min in integrators

[reneSchm]

now the integrator makes a time step of at least dt_min, instead of potentially doing nothing except decreasing dt without limit. this breaks fsal steppers (in part. rk-dopri5, only used in tests and benchmarks), so they were disabled for now.

Changes and Information

Please briefly list the changes (main added features, changed items, or corrected bugs) made:

• made changes to ControlledStepperWrapper such that the step method
  • correctly uses dt_min
  • returns a dt in [dt_min, dt_max]
  • always makes an integration step of at least size dt_min, even if the step got rejected due to not reaching tolerances
• updated RKIntegratorCore to behave similarly.

If need be, add additional information and what the reviewer should look out for in particular:

• Please check that the requirements listed in the documentation of IntegratorCore::step are understandable and sufficient.
• Also consider whether the compare-with tests are correctly updated, or if the testing method might need to change.

Merge Request - Guideline Checklist

Please check our git workflow. Use the draft feature if the Pull Request is not yet ready to review.

Checks by code author

• Every addressed issue is linked (use the "Closes #ISSUE" keyword below)
• New code adheres to 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 Markdown files if necessary)
• Proper attention to licenses, especially no new third-party software with conflicting license has been added
• (For ABM development) Checked benchmark results and ran and posted a local test above from before and after development to ensure performance is monitored.

Checks by code reviewer(s)

• Corresponding issue(s) is/are linked and addressed
• Code is clean of development artifacts (no deactivated or commented code lines, no debugging printouts, etc.)
• Appropriate unit tests have been added, CI passes, code coverage and performance is acceptable (did not decrease)
• No large data files added in the whole history of commits(files should in sum not exceed 100 KB, avoid PDFs, Word docs, etc.)
• On merge, add 2-5 lines with the changes (main added features, changed items, or corrected bugs) to the merge-commit-message. This can be taken from the briefly-list-the-changes above (best case) or the separate commit messages (worst case).

Closes #949

[reneSchm]

I fixed the issue using another version of the try_step method. However, this breaks the runge_kutta_dopri5 stepper due to a bug in boost. Since we currently only use it in testing and benchmarks, I decided to prevent it (and other FSAL steppers) to be used in the wrapper.

@mknaranja Do you think its okay to disable FSAL steppers for now? Otherwise, the ControlledStepperWrapper would have to make a manual integration step whenever the normal try_step loop fails. This is not a terrible solution, but we potentially make a repeated integration step whenever dt==dt_min.

The test TestOdeSecir.compareAgeResWithPreviousRunCashKarp currently fails since it uses an integration result with a time step < dt_min as reference. @MaxBetzDLR is it safe to just update the csv with a new integration result?

[reneSchm]

Currently, both ControlledStepperWrapper and RKIntegratorCore are allowed to make steps of size > dt_max, only if the function argument dt > dt_max.
We could either keep this behaviour, I don't think it will matter much or have a negative impact. But it may be unexpected, so we could warn/exit on a bad input; or allow any input value, but restrict it to [dt_min, dt_max] internally.
I would prefer restricting dt internally, as we already enforce the lower bound of dt_min.

[codecov]

Codecov Report

Attention: Patch coverage is 97.67442% with 1 lines in your changes are missing coverage. Please review.

Project coverage is 96.23%. Comparing base (058b3cc) to head (31f4ae7).
Report is 3 commits behind head on main.

Files	Patch %	Lines
cpp/memilio/math/integrator.cpp	0.00%	1 Missing ⚠️

Additional details and impacted files

@@           Coverage Diff           @@
##             main     #960   +/-   ##
=======================================
  Coverage   96.22%   96.23%           
=======================================
  Files         124      122    -2     
  Lines        9699     9684   -15     
=======================================
- Hits         9333     9319   -14     
+ Misses        366      365    -1     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

[mknaranja]

Currently, both ControlledStepperWrapper and RKIntegratorCore are allowed to make steps of size > dt_max, only if the function argument dt > dt_max. We could either keep this behaviour, I don't think it will matter much or have a negative impact. But it may be unexpected, so we could warn/exit on a bad input; or allow any input value, but restrict it to [dt_min, dt_max] internally. I would prefer restricting dt internally, as we already enforce the lower bound of dt_min.

@reneSchm This is not good design we did. We should not allow that, adapt it and inform/warn on it.

[mknaranja]

Hi @reneSchm
Thank you very much. I had looked at almost everything. I have mostly commenting and documentation suggestions. If someone could have a more detailed look on the new step method, that would be great though. If that's done, we can merge it.

[HenrZu]

I looked at the step functions and they look correct and also well commented:)
One comment. Should we also adjust the default values for the minimum step sizes? (also in stepper_wrapper.h)

memilio/cpp/memilio/math/adapt_rk.h

Line 108 in afd4982

, m_dt_min(std::numeric_limits<double>::min())

[reneSchm]

One comment. Should we also adjust the default values for the minimum step sizes? (also in stepper_wrapper.h)

I don't think we should set a different dt_min, which scale would we even choose?

Maybe we could set a dt_min in the simulation, but I would prefer to keep it model/executable specific.

[reneSchm]

So I added a small warning when the given dt was out of bounds, and strangely enough it might have improved performance for some steppers?

• RKIntegratorCore aka adapt_rk and ControlledStepperWrapper<runge_kutta_fehlberg78> aka rkf78 are consistently about 50ns faster. Which is nice, but confusing. Apart from the couple lines for the logging, the assembly looked like only some addresses got reordered, at least I didn't see any significant changes. My best guess is that I got lucky and some reordering caused better cache hits, but I don't know.
• ControlledStepperWrapper<runge_kutta_cash_carp54> aka rk_ck54 seems to be roughly 10ns slower, but that is within the margin of error for these measurements. This is as expected, and I think well acceptable.

Before logging (dae153c):

simulate SecirModel adapt_rk            1008 ns         1007 ns       696182
simulate SecirModel boost rk_ck54        746 ns          746 ns       938309
simulate SecirModel boost rkf78         2196 ns         2195 ns       317635

With logging (b93c0e6):

simulate SecirModel adapt_rk             935 ns          934 ns       749251
simulate SecirModel boost rk_ck54        759 ns          759 ns       921252
simulate SecirModel boost rkf78         2141 ns         2141 ns       323717