DAE Problems
#
SciMLBase.DAEProblem
— Type
Defines an implicit ordinary differential equation (ODE) or differential-algebraic equation (DAE) problem. Documentation Page: https://docs.sciml.ai/DiffEqDocs/stable/types/dae_types/
Mathematical Specification of an DAE Problem
To define a DAE 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(resid,du,u,p,t)
). 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.
Problem Type
Constructors
-
DAEProblem(f::DAEFunction,du0,u0,tspan,p=NullParameters();kwargs...)
-
DAEProblem{isinplace,specialize}(f,du0,u0,tspan,p=NullParameters();kwargs...)
: Defines the DAE with the specified functions.isinplace
optionally sets whether the function is inplace or not. This is determined automatically, but not inferred.specialize
optionally controls the specialization level. See the specialization levels section of the documentation for more details. The default isAutoSpecialize
.
For more details on the in-place and specialization controls, see the ODEFunction documentation.
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.
For specifying Jacobians and mass matrices, see the DiffEqFunctions page.
Fields
-
f
: The function in the ODE. -
du0
: The initial condition for the derivative. -
u0
: The initial condition. -
tspan
: The timespan for the problem. -
differential_vars
: A logical array which declares which variables are the differential (non-algebraic) vars (i.e.du'
is in the equations for this variable). Defaults to nothing. Some solvers may require this be set if an initial condition needs to be determined. -
p
: The parameters for the problem. Defaults toNullParameters
-
kwargs
: The keyword arguments passed onto the solves.
Example Problems
Examples problems can be found in DiffEqProblemLibrary.jl.
To use a sample problem, such as prob_dae_resrob
, you can do something like:
#] add DAEProblemLibrary
using DAEProblemLibrary
prob = DAEProblemLibrary.prob_dae_resrob
sol = solve(prob,IDA())
#
SciMLBase.DAEFunction
— Type
DAEFunction{iip,F,Ta,Tt,TJ,JVP,VJP,JP,SP,TW,TWt,TPJ,S,S2,S3,O,TCV} <: AbstractDAEFunction{iip,specialize}
A representation of an implicit DAE function f
, defined by:
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
DAEFunction{iip,specialize}(f;
analytic = __has_analytic(f) ? f.analytic : 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,
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 function f
itself is required. This function should be given as f!(out,du,u,p,t)
or out = f(du,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:
-
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. -
jac(J,du,u,p,gamma,t)
orJ=jac(du,u,p,gamma,t)
: returns the implicit DAE Jacobian defined as -
jvp(Jv,v,du,u,p,gamma,t)
orJv=jvp(v,du,u,p,gamma,t)
: returns the directional derivative$\frac{df}{du} v$ -
vjp(Jv,v,du,u,p,gamma,t)
orJv=vjp(v,du,u,p,gamma,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. -
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 DAEFunction type directly match the names of the inputs.
Examples
Declaring Explicit Jacobians for DAEs
For fully implicit ODEs (DAEProblem
s), a slightly different Jacobian function is necessary. For the DAE
The Jacobian should be given in the form gamma*dG/d(du) + dG/du
where gamma
is given by the solver. This means that the signature is:
f(J,du,u,p,gamma,t)
For example, for the equation
function testjac(res,du,u,p,t)
res[1] = du[1] - 2.0 * u[1] + 1.2 * u[1]*u[2]
res[2] = du[2] -3 * u[2] - u[1]*u[2]
end
we would define the Jacobian as:
function testjac(J,du,u,p,gamma,t)
J[1,1] = gamma - 2.0 + 1.2 * u[2]
J[1,2] = 1.2 * u[1]
J[2,1] = - 1 * u[2]
J[2,2] = gamma - 3 - u[1]
nothing
end
Symbolically Generating the Functions
See the modelingtoolkitize
function from ModelingToolkit.jl for automatically symbolically generating the Jacobian and more from the numerically-defined functions.
Solution Type
#
SciMLBase.DAESolution
— Type
struct DAESolution{T, N, uType, duType, uType2, DType, tType, P, A, ID, S} <: SciMLBase.AbstractDAESolution{T, N, uType}
Representation of the solution to an differential-algebraic equation defined by an DAEProblem.
DESolution Interface
For more information on interacting with DESolution
types, check out the Solution Handling page of the DifferentialEquations.jl documentation.
Fields
-
u
: the representation of the DAE solution. Given as an array of solutions, whereu[i]
corresponds to the solution at timet[i]
. It is recommended in most cases one does not accesssol.u
directly and instead use the array interface described in the Solution Handling page of the DifferentialEquations.jl documentation. -
du
: the representation of the derivatives of the DAE solution. -
t
: the time points corresponding to the saved values of the DAE solution. -
prob
: the original DAEProblem that was solved. -
alg
: the algorithm type used by the solver. -
stats
: statistics of the solver, such as the number of function evaluations required, number of Jacobians computed, and more. -
retcode
: the return code from the solver. Used to determine whether the solver solved successfully, whether it terminated early due to a user-defined callback, or whether it exited due to an error. For more details, see the return code documentation.
Example Problems
Examples problems can be found in DiffEqProblemLibrary.jl.
To use a sample problem, such as prob_dae_resrob
, you can do something like:
using DiffEqProblemLibrary.DAEProblemLibrary, Sundials
prob = DAEProblemLibrary.prob_dae_resrob
sol = solve(prob, IDA())
#
DAEProblemLibrary.prob_dae_resrob
— Constant
The Robertson biochemical reactions in DAE form
where , , . For details, see: Hairer Norsett Wanner Solving Ordinary Differential Equations I - Nonstiff Problems Page 129 Usually solved on ]