Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

simplify steadystate logic #1485

Closed
wants to merge 17 commits into from
Closed

Conversation

FFroehlich
Copy link
Member

No description provided.

@codecov
Copy link

codecov bot commented Apr 3, 2021

Codecov Report

Merging #1485 (99d0cfb) into develop (3d61ac6) will increase coverage by 0.12%.
The diff coverage is 97.95%.

❗ Current head 99d0cfb differs from pull request most recent head 2c4c380. Consider uploading reports for the commit 2c4c380 to get more accurate results
Impacted file tree graph

@@             Coverage Diff             @@
##           develop    #1485      +/-   ##
===========================================
+ Coverage    79.43%   79.55%   +0.12%     
===========================================
  Files           66       66              
  Lines        10370    10374       +4     
===========================================
+ Hits          8237     8253      +16     
+ Misses        2133     2121      -12     
Flag Coverage Δ
cpp 75.94% <97.95%> (+0.18%) ⬆️
petab 67.85% <ø> (ø)
python 68.15% <ø> (ø)

Flags with carried forward coverage won't be shown. Click here to find out more.

Impacted Files Coverage Δ
python/amici/numpy.py 82.08% <ø> (ø)
src/steadystateproblem.cpp 88.08% <97.91%> (+2.05%) ⬆️
include/amici/steadystateproblem.h 100.00% <100.00%> (ø)
src/newton_solver.cpp 76.75% <0.00%> (-2.17%) ⬇️
src/solver_cvodes.cpp 70.41% <0.00%> (+0.20%) ⬆️
src/spline.cpp 71.57% <0.00%> (+6.31%) ⬆️
src/exception.cpp 91.30% <0.00%> (+8.69%) ⬆️

@paulstapor
Copy link
Contributor

That one was on my to-do-list anyway... Would have been good to discuss this prior to you putting effort in it...
The logic here isn't simple, but that's not (let's say, maybe not only) because it's not perfectly programmed, but rather because the problem as such is intricate. You probably saw it from the fact that you want to simplify the code, but actually added more lines than you removed...

So, give me some background please, as the PR has no description:

@FFroehlich
Copy link
Member Author

FFroehlich commented Apr 6, 2021

That one was on my to-do-list anyway... Would have been good to discuss this prior to you putting effort in it...
The logic here isn't simple, but that's not (let's say, maybe not only) because it's not perfectly programmed, but rather because the problem as such is intricate. You probably saw it from the fact that you want to simplify the code, but actually added more lines than you removed...

The big issue is that there was a single function that served 3 different purposes (check whether simulation needs to be stored, whether simulation needs to be run with sensitivities, and whether sensitivities are to be computed using the implicit function theorem). This is in violation of the single-responsibility principle (https://en.wikipedia.org/wiki/Single-responsibility_principle), giving rise to code execution that is so complicated that eve you got confused about it #1461 .

So, give me some background please, as the PR has no description:

  • What is precisely you aim of the code? What does the current code not do? What should the new code do (better)?
  1. the simulationFSA sensitivity mode is intended to be used when the Jacobian is singular. In that case, sensitivity computation via the implicit function theorem should never be used.
  2. Having a mode where you use forward sensitivities only when simulation is run is not practically relevant, at least I cannot tell any situations when the user would want to do that. Primarily, since simulation has to be run until convergence of sensitivities. Moreover, it is questionable whether forward sensitivities will be accurate when newtons method was used prior to simulation, since the sensitivity initialization are not updated by newtons method. While you may get lucky and they may converge eventually, but whether this will work in multistable systems is unclear.
  3. You cannot run adjoint sensitivity analysis and then use the same solver for postequilibration using sensitivity analysis. I added a check that detects that situation and raises an issue
  4. Many of the test that were supposed to test combinations of forward and adjoint methods simply used newtons method to compute sensitivity and hence didn't really test what they were supposed to test. I added checks for preeq and posteq status flags that verifies that the right method were run.
  5. Convergence of steadystate sensitivities were based on norm of x_, not sx_, which honestly just doesn't make too much sense.

It partially addresses the item Add flag "simulationOnly" to steadyStateSensitivityMethod, which does not even try to launch the Newton solver.

Copy link
Contributor

@paulstapor paulstapor left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

That's somewhat orthogonal to what I have planned to do with the steady state solver, therefore I'm not completely happy about it
If I'm not mistaken, the SteadyStateContext flag isn't used any more by now. In this case, please remove it also from include/amici/defines.h

Comment on lines 45 to 52

if (solver.getSensitivityMethodPreequilibration()
== SensitivityMethod::adjoint &&
model.getSteadyStateSensitivityMode() ==
SteadyStateSensitivityMode::simulationFSA)
throw AmiException("Preequilibration using adjoint sensitivities "
"is not compatible with simulationFSA "
"sensitivity mode.");
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I understand that setting SteadyStateSensitivityMode::simulationFSA and allowing to use adjoint preequilibration is counter-intuitive. However, it would make much sense to allow the user to switch off the Newton solver and using adjoint preequilibration... Would be happy if there was a way to encode this...

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I understand that setting SteadyStateSensitivityMode::simulationFSA and allowing to use adjoint preequilibration is counter-intuitive. However, it would make much sense to allow the user to switch off the Newton solver and using adjoint preequilibration... Would be happy if there was a way to encode this...

How about setting the allowed newton steps to 0?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ideally, there would be the options newtonOnly, simulationFSA, and simulationASA and something that first tries Newton and uses simulation only if Newton fails... But no need to implement this now.

Something else: I'd advocate to move the steadyStateSensitivitiyMode from model to solver. The value may be model dpeendent, but it encodes the behavior of the solver. Therefore, I'd be happy if we could move this at some point. What's your opinion on this?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ideally, there would be the options newtonOnly, simulationFSA, and simulationASA and something that first tries Newton and uses simulation only if Newton fails... But no need to implement this now.

I think most of these settings can already be achieved with existing options.

Something else: I'd advocate to move the steadyStateSensitivitiyMode from model to solver. The value may be model dpeendent, but it encodes the behavior of the solver. Therefore, I'd be happy if we could move this at some point. What's your opinion on this?

Yes, I agree, I found it confusing that this is in model, not solver.

Comment on lines 69 to 72
if (solver->getSensitivityOrder() >= SensitivityOrder::first &&
model->getSteadyStateSensitivityMode() ==
SteadyStateSensitivityMode::simulationFSA &&
solver->getSensitivityMethod() != SensitivityMethod::forward) {
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Again not quite getting why you want to exclude this, but okay if you feel like it's really necessary...
Anyhow: All these quantities should have been set prior to the creation of the SteadyStateProblem object and could be caught in the constructor (together with the additional check you implemented there), ideally by using a method like check_steady_state_solver_options or similar...

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This check may have become obsolete after adding the respective checks to the constructor.

@@ -462,7 +472,7 @@ bool SteadystateProblem::checkConvergence(const Solver *solver,
sx_ = solver->getStateSensitivity(t_);
model->fsxdot(t_, x_, dx_, ip, sx_[ip], dx_, xdot_);
wrms_ = getWrmsNorm(
x_, xdot_, solver->getAbsoluteToleranceSteadyStateSensi(),
sx_[ip], xdot_, solver->getAbsoluteToleranceSteadyStateSensi(),
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

👍

include/amici/steadystateproblem.h Outdated Show resolved Hide resolved
@FFroehlich
Copy link
Member Author

That's somewhat orthogonal to what I have planned to do with the steady state solver, therefore I'm not completely happy about it

If you plan stuff for a year or so and never follow through, you kinda have to live with the possibility that someone starts working on it beforehand.

If I'm not mistaken, the SteadyStateContext flag isn't used any more by now. In this case, please remove it also from include/amici/defines.h

👍

@paulstapor
Copy link
Contributor

That's somewhat orthogonal to what I have planned to do with the steady state solver, therefore I'm not completely happy about it

If you plan stuff for a year or so and never follow through, you kinda have to live with the possibility that someone starts working on it beforehand.

No problem with somebody else working on it at all. But I created a milestone and organized issues for this for the precise reason that one can discuss the way how its implemented prior to it being done...

@FFroehlich
Copy link
Member Author

That's somewhat orthogonal to what I have planned to do with the steady state solver, therefore I'm not completely happy about it

@paulstapor I added a mixed mode which restores the possibility to use combinations of FSA with newton methods

@sonarcloud
Copy link

sonarcloud bot commented Apr 11, 2021

@FFroehlich
Copy link
Member Author

okay at this point I simply give up and consider the simulationFSA functionality broken.

for setting1, setting2 in itertools.product(settings, settings):
never ran because once a generator is used, its empty.

@FFroehlich FFroehlich closed this Apr 12, 2021
@dweindl dweindl deleted the simplify_steadystate_logic branch April 8, 2022 21:15
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants