Skip to content

Include root-finding equations in ODESystem #1337

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

Merged
merged 24 commits into from
Nov 11, 2021
Merged

Conversation

baggepinnen
Copy link
Contributor

@baggepinnen baggepinnen commented Nov 8, 2021

This PR adds a vector of equations to an ODESystem that can be used for root-finding when solving an ODEProblem.
The intended usage is to specify such equations when an ODESystem is created rather than when the ODEProblem is created, for example, when a component in a component library contains discountinous dynamics etc.

The provided tests pass locally, but there are a number of points I would like to get feedback on

  1. The root_eqs field was added to ODESystem only, how should all other AbstractSystems be treated?
  2. I wasn't able to use the inplace update for VectorContinuousCallback here, can the RuntimeGeneratedFunction somehow be generated in an inplace form as well?

Fixes #38

@YingboMa
Copy link
Member

YingboMa commented Nov 8, 2021

RuntimeGeneratedFunction handles arbitrary functions without closures in-place or not. Try to read build_function's implementation in Symbolics to see how in-place functions are expressed.

@baggepinnen baggepinnen changed the title WIP: Include root-finding equations in ODESystem Include root-finding equations in ODESystem Nov 9, 2021
@baggepinnen
Copy link
Contributor Author

I changed the implementation to use build_function instead of gen_nlsolve and it appears to work well. The remaining problem is related to other system types, e.g., NonlinearSystem does not have any root_eqs. The code in MTK assumes that all system types have the same fields in many places, for instance here https://github.com/SciML/ModelingToolkit.jl/pull/1337/files#diff-287e017b6c265948a34862d11b63bade810cb38f11173644350028d4ba836661R989
One can either introduce the field root_eqs for other system types, or add the keyword argument to their constructors, check that the value of the argument is empty/nothing. I took the latter route in ce89be9

@baggepinnen baggepinnen requested a review from YingboMa November 9, 2021 12:23
@baggepinnen
Copy link
Contributor Author

After offline discussions, it was decided that we should add affect!s to the interface as well. The commits above does this, but requires the user to write the affect! function using the standard integrator interface, which is inconvenient.

The remaining step is to allow symbolic affect expressions that are automatically compiled to an affect! function operator on the integrator interface.

@baggepinnen
Copy link
Contributor Author

baggepinnen commented Nov 10, 2021

One can now model events including effects like so

## Test bouncing ball with equation affect
@variables t x(t)=1 v(t)=0
D = Differential(t)

root_eqs = [x ~ 0]
affect   = [v ~ -v]

@named ball = ODESystem([
    D(x) ~ v
    D(v) ~ -9.8
], t, events = root_eqs => affect)

ball = structural_simplify(ball)

tspan = (0.0,5.0)
prob = ODEProblem(ball, Pair[], tspan)
sol = solve(prob,Tsit5())
@test 0 <= minimum(sol[x]) <= 1e-10 # the ball never went through the floor but got very close
plot(sol)

some restrictions apply on the affect equations, they are detailed in the docs, https://github.com/SciML/ModelingToolkit.jl/pull/1337/files#diff-3aa3a02fd9ea8d991b959e2b999e5e717502e8df1d1d57d4e49199fb1cfd1c5dR262

@ChrisRackauckas
Copy link
Member

While we're at it, to finish #38 we should have a way of designating an event condition in the discrete sense as well.

@isaacsas
Copy link
Member

Could the symbolic conditions and affects be stored in the current continuous and discrete callbacks provided by DiffEqBase? It might be nice to reuse those structures since users are already familiar with them.

@baggepinnen
Copy link
Contributor Author

baggepinnen commented Nov 10, 2021

Could the symbolic conditions and affects be stored in the current continuous and discrete callbacks provided by DiffEqBase? It might be nice to reuse those structures since users are already familiar with them.

Such callbacks are generated internally and are available through the created ODEProblem as prob.kwargs[:callback].
https://github.com/SciML/ModelingToolkit.jl/pull/1337/files#diff-4ff6b77a6fa947c70cfd3bcaf373b11d61ad04f2f1a74022ea9250e386459f74R156-R180

@isaacsas
Copy link
Member

isaacsas commented Nov 10, 2021

I meant that users create and pass the symbolic callbacks using the same structures (but built with symbolic equations), and use callback as a kwarg to ODESystems (as opposed to the events kwarg used currently and it’s format). This would keep the symbolic workflow and terminology the same as the DifferentialEquations.jl workflow and terminology from a user’s perspective.

@ChrisRackauckas
Copy link
Member

@baggepinnen we should at least thing of the naming scheme so that way the next PR isn't immediately breaking. Is it continuous_events and discrete_events?

@baggepinnen
Copy link
Contributor Author

continuous_events and discrete_events sound reasonable. The internal structure EqAffectPair can be renamed SymbolicContinousCallback as well to align more with the standard diffeq interface.

@ChrisRackauckas
Copy link
Member

Sounds reasonable to me. Yeah, let's get that together and merge, and then follow up with the discrete.

@isaacsas
Copy link
Member

Could you add a test with multiple affects for a given condition? i.e. changing more than one variable?

@ChrisRackauckas
Copy link
Member

How is this working when a system with an event is used as a subsystem in another system? Where is it appropriately renamespacing the events?

@baggepinnen
Copy link
Contributor Author

There is a test that the event equations have been renamespaced here
https://github.com/SciML/ModelingToolkit.jl/pull/1337/files#diff-4ff6b77a6fa947c70cfd3bcaf373b11d61ad04f2f1a74022ea9250e386459f74R133

and there is a method to renamespace events here
https://github.com/SciML/ModelingToolkit.jl/pull/1337/files#diff-287e017b6c265948a34862d11b63bade810cb38f11173644350028d4ba836661R188

@ChrisRackauckas
Copy link
Member

Missed it thanks.

Comment on lines 159 to 161
callbacks = map(cbs) do eq_aff
generate_rootfinding_callback(eq_aff, sys, dvs, ps; kwargs...)
end
Copy link
Member

Choose a reason for hiding this comment

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

If you map here, won't it not fuse into a single VectorContinuousCallback?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

It only returns VectorContinuousCallback if there is an event with multiple equations. If there are several events, they are treated as separate callbacks. So this kind of event leads to VectorContinuousCallback
https://github.com/SciML/ModelingToolkit.jl/pull/1337/files#diff-4ff6b77a6fa947c70cfd3bcaf373b11d61ad04f2f1a74022ea9250e386459f74R194

while this kind
https://github.com/SciML/ModelingToolkit.jl/pull/1337/files#diff-4ff6b77a6fa947c70cfd3bcaf373b11d61ad04f2f1a74022ea9250e386459f74R132-R133
leads to a CallbackSet.

I kept them apart in order to make the following check
https://github.com/SciML/ModelingToolkit.jl/pull/1337/files#diff-7195a73d9d71f7a993fe71086614d700629ebc22306be7b602f60fbc054bcdb7R219-R220

I couldn't prove for myself that there couldn't be bad outcomes from having one state affected by multiple equations in the same event since the ordering of the affect equations would then matter etc.

Copy link
Member

Choose a reason for hiding this comment

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

How do you make a VectorContinuousCallback which have a different affect based on the condition that is fired? Why not fuse anyways since it's more efficient?

@ChrisRackauckas
Copy link
Member

@YingboMa you got anything else?

@YingboMa
Copy link
Member

LGTM

@YingboMa YingboMa merged commit 52bdc81 into SciML:master Nov 11, 2021
@baggepinnen baggepinnen deleted the rootfuns branch November 11, 2021 19:42
@tshort
Copy link

tshort commented Nov 12, 2021

Cool stuff. Any chance it could be expanded to cover affect_neg!?

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.

Support for events
5 participants