Skip to content

Logging for isoutofdomain and add interface for domain_checks - #4114

Open
Shreyas-Ekanathan wants to merge 3 commits into
SciML:masterfrom
Shreyas-Ekanathan:domain_logging
Open

Shreyas-Ekanathan wants to merge 3 commits into
SciML:masterfrom
Shreyas-Ekanathan:domain_logging

Conversation

@Shreyas-Ekanathan

Copy link
Copy Markdown
Contributor

This PR does three things:

  1. Creates an interfance called domain_checks that behaves similar to variable bounds in MTK directly in OrdinaryDiffEq, for reporting of variable domains. This is enforced on every RHS call, as opposed to isoutofdomain, which is enforced after each step.
  2. Adds logging for isoutofdomain by building an AST like structure that can be traced to determine which premises failed. This has very little overhead (if any) relative to the current scheme by using a macro, but needs some new user input (specifically isoutofdomain needs to be provided using @isoutofdomain).
  3. Adds a hook for MTK's variable bounds to directly be inserted into the solve via the isoutofdomain kwarg. A PR (coming shortly) to MTK will make use of this hook.

Some code examples:

domain_checks and logging:

f(u, p, t) = -u
prob = ODEProblem(f, [1.0], (0.0, 1.0))
sol = solve(
    prob, Tsit5(); dtmin = 1e-10,
    domain_checks = [(u, p, t) -> u[1] >= 1.0 => "u must not decrease"]
)

┌ Warning: Verbosity toggle: dt_min_unstable
│  dt(1.0e-10) <= dtmin(1.0e-10) at t=0.0. Aborting. There is either an error in your model specification or the true solution is unstable.
│ 
│ Diagnostics:
│ 
│ State Analysis:
│   All 1 state variables are non-finite (NaN/Inf)
│ 
│ Domain Check Analysis: Gates inputs to function evaluations:1 domain_checks predicate(s) failed:
│   [1] u must not decrease
│ 
│ Error Analysis:
│   step error estimate EEst = NaN (a step is accepted when EEst <= 1)
└ @ SciMLBase ~/Documents/GitHub/SciMLBase.jl/src/integrator_interface.jl:1007

isoutofdomain and logging:

f(u, p, t) = -u
prob = ODEProblem(f, [1.0], (0.0, 1.0))
integrator = init(
    prob, Tsit5(); dtmin = 1e-10,
    isoutofdomain = @isoutofdomain (u, p, t) -> u[1] > 0.9 || u[1] < -5.0
)
solve!(integrator)

┌ Warning: Verbosity toggle: dt_min_unstable 
│  dt(1.0e-10) <= dtmin(1.0e-10) at t=0.0. Aborting. There is either an error in your model specification or the true solution is unstable.
│ 
│ Diagnostics:
│ 
│ Is Out of Domain: Gates results of integrator steps:
│   isoutofdomain predicate returned true for the proposed state, causing failure of the step
│   predicate: u[1] > 0.9 || u[1] < -5.0
│     [1] u[1] > 0.9  => true  (u[1] = 1)
│     [2] u[1] < -5.0 => false  (u[1] = 1)
│ 
│ Error Analysis:
│   step error estimate EEst = 1.065e-23 (a step is accepted when EEst <= 1)
└ @ SciMLBase ~/Documents/GitHub/SciMLBase.jl/src/integrator_interface.jl:1007

Comment thread lib/DiffEqBase/src/solve.jl Outdated
Comment on lines +555 to +556
* `domain_checks`: Specifies a vector of predicates, each `(u,p,t) -> Bool` or
`(u,p,t) -> Pair{Bool,<:AbstractString}` (to attach a custom message), checked

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

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

This can have function wrapping issues? Why not require the user supplies the struct instead?

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

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

that would work but has heavy user-side implications right? realistically a user doesn't want to supply the struct.

an alternative that is cleanly function wrapped would be to wrap the function in these domain checks prior to function wrapping itself happening, which is then stable throughout. This is a fairly complicated change though, I've scoped it out and pushed it so we can choose if we like this. Let me know your thoughts, I'm open to changes here since this is a big design choice.

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