Dynamical, Hamiltonian and 2nd Order ODE Problems
#
SciMLBase.DynamicalODEProblem
— Type
Defines a dynamical ordinary differential equation (ODE) problem. Documentation Page: https://docs.sciml.ai/DiffEqDocs/stable/types/dynamical_types/
Dynamical ordinary differential equations, such as those arising from the definition of a Hamiltonian system or a second order ODE, have a special structure that can be utilized in the solution of the differential equation. On this page, we describe how to define second order differential equations for their efficient numerical solution.
Mathematical Specification of a Dynamical ODE Problem
These algorithms require a Partitioned ODE of the form:
This is a Partitioned ODE partitioned into two groups, so the functions should be specified as f1(dv,v,u,p,t)
and f2(du,v,u,p,t)
(in the inplace form), where f1
is independent of v
(unless specified by the solver), and f2
is independent of u
and t
. This includes discretizations arising from SecondOrderODEProblem
s where the velocity is not used in the acceleration function, and Hamiltonians where the potential is (or can be) time-dependent, but the kinetic energy is only dependent on v
.
Note that some methods assume that the integral of f2
is a quadratic form. That means that f2=v'*M*v
, i.e. , giving du = v
. This is equivalent to saying that the kinetic energy is related to . The methods which require this assumption will lose accuracy if this assumption is violated. Methods listed make note of this requirement with "Requires quadratic kinetic energy".
Constructor
DynamicalODEProblem(f::DynamicalODEFunction,v0,u0,tspan,p=NullParameters();kwargs...)
DynamicalODEProblem{isinplace}(f1,f2,v0,u0,tspan,p=NullParameters();kwargs...)
Defines the ODE with the specified functions. isinplace
optionally sets whether the function is inplace or not. This is determined automatically, but not inferred.
Parameters are optional, and if not given, then a NullParameters()
singleton will be used which will throw nice errors if you try to index non-existent parameters. Any extra keyword arguments are passed on to the solvers. For example, if you set a callback
in the problem, then that callback
will be added in every solve call.
Fields
-
f1
andf2
: The functions in the ODE. -
v0
andu0
: The initial conditions. -
tspan
: The timespan for the problem. -
p
: The parameters for the problem. Defaults toNullParameters
-
kwargs
: The keyword arguments passed onto the solves.
#
SciMLBase.SecondOrderODEProblem
— Type
Defines a second order ordinary differential equation (ODE) problem. Documentation Page: https://docs.sciml.ai/DiffEqDocs/stable/types/dynamical_types/
Mathematical Specification of a 2nd Order ODE Problem
To define a 2nd Order ODE Problem, you simply need to give the function and the initial condition which define an ODE:
f
should be specified as f(du,u,p,t)
(or in-place as f(ddu,du,u,p,t)
), and u₀
should be an AbstractArray (or number) whose geometry matches the desired geometry of u
. Note that we are not limited to numbers or vectors for u₀
; one is allowed to provide u₀
as arbitrary matrices / higher dimension tensors as well.
From this form, a dynamical ODE:
is generated.
Constructors
SecondOrderODEProblem{isinplace}(f,du0,u0,tspan,callback=CallbackSet())
Defines the ODE with the specified functions.
Fields
-
f
: The function for the second derivative. -
du0
: The initial derivative. -
u0
: The initial condition. -
tspan
: The timespan for the problem. -
callback
: A callback to be applied to every solver which uses the problem. Defaults to nothing.
#
SciMLBase.DynamicalODEFunction
— Type
DynamicalODEFunction{iip,F1,F2,TMM,Ta,Tt,TJ,JVP,VJP,JP,SP,TW,TWt,TPJ,S,S2,S3,O,TCV} <: AbstractODEFunction{iip,specialize}
A representation of an ODE function f
, defined by:
as a partitioned ODE:
and all of its related functions, such as the Jacobian of f
, its gradient with respect to time, and more. For all cases, u0
is the initial condition, p
are the parameters, and t
is the independent variable.
Constructor
DynamicalODEFunction{iip,specialize}(f1,f2;
mass_matrix = __has_mass_matrix(f) ? f.mass_matrix : I,
analytic = __has_analytic(f) ? f.analytic : nothing,
tgrad= __has_tgrad(f) ? f.tgrad : nothing,
jac = __has_jac(f) ? f.jac : nothing,
jvp = __has_jvp(f) ? f.jvp : nothing,
vjp = __has_vjp(f) ? f.vjp : nothing,
jac_prototype = __has_jac_prototype(f) ? f.jac_prototype : nothing,
sparsity = __has_sparsity(f) ? f.sparsity : jac_prototype,
paramjac = __has_paramjac(f) ? f.paramjac : nothing,
syms = __has_syms(f) ? f.syms : nothing,
indepsym= __has_indepsym(f) ? f.indepsym : nothing,
paramsyms = __has_paramsyms(f) ? f.paramsyms : nothing,
colorvec = __has_colorvec(f) ? f.colorvec : nothing,
sys = __has_sys(f) ? f.sys : nothing)
Note that only the functions f_i
themselves are required. These functions should be given as f_i!(du,u,p,t)
or du = f_i(u,p,t)
. See the section on iip
for more details on in-place vs out-of-place handling.
All of the remaining functions are optional for improving or accelerating the usage of f
. These include:
-
mass_matrix
: the mass matrixM_i
represented in the ODE function. Can be used to determine that the equation is actually a differential-algebraic equation (DAE) ifM
is singular. Note that in this case special solvers are required, see the DAE solver page for more details: https://docs.sciml.ai/DiffEqDocs/stable/solvers/dae*solve/. Must be an AbstractArray or an AbstractSciMLOperator. Should be given as a tuple of mass matrices, i.e.(M*1, M_2)
for the mass matrices of equations 1 and 2 respectively. -
analytic(u0,p,t)
: used to pass an analytical solution function for the analytical solution of the ODE. Generally only used for testing and development of the solvers. -
tgrad(dT,u,p,t)
or dT=tgrad(u,p,t): returns -
jac(J,u,p,t)
orJ=jac(u,p,t)
: returns -
jvp(Jv,v,u,p,t)
orJv=jvp(v,u,p,t)
: returns the directional derivative$\frac{df}{du} v$ -
vjp(Jv,v,u,p,t)
orJv=vjp(v,u,p,t)
: returns the adjoint derivative$\frac{df}{du}^\ast v$ -
jac_prototype
: a prototype matrix matching the type that matches the Jacobian. For example, if the Jacobian is tridiagonal, then an appropriately sizedTridiagonal
matrix can be used as the prototype and integrators will specialize on this structure where possible. Non-structured sparsity patterns should use aSparseMatrixCSC
with a correct sparsity pattern for the Jacobian. The default isnothing
, which means a dense Jacobian. -
paramjac(pJ,u,p,t)
: returns the parameter Jacobian . -
syms
: the symbol names for the elements of the equation. This should matchu0
in size. For example, ifu0 = [0.0,1.0]
andsyms = [:x, :y]
, this will apply a canonical naming to the values, allowingsol[:x]
in the solution and automatically naming values in plots. -
indepsym
: the canonical naming for the independent variable. Defaults to nothing, which internally usest
as the representation in any plots. -
paramsyms
: the symbol names for the parameters of the equation. This should matchp
in size. For example, ifp = [0.0, 1.0]
andparamsyms = [:a, :b]
, this will apply a canonical naming to the values, allowingsol[:a]
in the solution. -
colorvec
: a color vector according to the SparseDiffTools.jl definition for the sparsity pattern of thejac_prototype
. This specializes the Jacobian construction when using finite differences and automatic differentiation to be computed in an accelerated manner based on the sparsity pattern. Defaults tonothing
, which means a color vector will be internally computed on demand when required. The cost of this operation is highly dependent on the sparsity pattern.
iip: In-Place vs Out-Of-Place
For more details on this argument, see the ODEFunction documentation.
specialize: Controlling Compilation and Specialization
For more details on this argument, see the ODEFunction documentation.
Fields
The fields of the DynamicalODEFunction type directly match the names of the inputs.
Solution Type
Dynamical ODE solutions return an ODESolution
. For more information, see the ODE problem definition page for the ODESolution
docstring.
Hamiltonian Problems
HamiltonianProblem
s are provided by DiffEqPhysics.jl and provide an easy way to define equations of motion from the corresponding Hamiltonian. To define a HamiltonianProblem
one only needs to specify the Hamiltonian:
and autodifferentiation (via ForwardDiff.jl) will create the appropriate equations.