ODE Problems

Defines an ordinary differential equation (ODE) problem. Documentation Page: https://diffeq.sciml.ai/stable/types/ode_types/

Mathematical Specification of an ODE Problem

To define an ODE Problem, you simply need to give the function $f$ and the initial condition $u_0$ which define an ODE:

\[M \frac{du}{dt} = f(u,p,t)\]

There are two different ways of specifying f:

• f(du,u,p,t): in-place. Memory-efficient when avoiding allocations. Best option for most cases unless mutation is not allowed.
• f(u,p,t): returning du. Less memory-efficient way, particularly suitable when mutation is not allowed (e.g. with certain automatic differentiation packages such as Zygote).

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.

For the mass matrix $M$, see the documentation of ODEFunction.

Problem Type

Constructors

ODEProblem can be constructed by first building an ODEFunction or by simply passing the ODE right-hand side to the constructor. The constructors are:

• ODEProblem(f::ODEFunction,u0,tspan,p=NullParameters();kwargs...)
• ODEProblem{isinplace}(f,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.

For specifying Jacobians and mass matrices, see the ODEFunction documentation.

Fields

• f: The function in the ODE.
• u0: The initial condition.
• tspan: The timespan for the problem.
• p: The parameters.
• kwargs: The keyword arguments passed onto the solves.

Example Problems

Example problems can be found in DiffEqProblemLibrary.jl.

To use a sample problem, such as prob_ode_linear, you can do something like:

#] add DiffEqProblemLibrary
using DiffEqProblemLibrary.ODEProblemLibrary
# load problems
ODEProblemLibrary.importodeproblems()
prob = ODEProblemLibrary.prob_ode_linear
sol = solve(prob)

ODEFunction{iip,F,TMM,Ta,Tt,TJ,JVP,VJP,JP,SP,TW,TWt,TPJ,S,S2,O,TCV} <: AbstractODEFunction{iip}

A representation of an ODE function f, defined by:

\[M \frac{du}{dt} = f(u,p,t)\]

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

ODEFunction{iip,recompile}(f;
                           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,
                           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!(du,u,p,t) or du = f(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 matrix M represented in the ODE function. Can be used to determine that the equation is actually a differential-algebraic equation (DAE) if M is singular. Note that in this case special solvers are required, see the DAE solver page for more details: https://diffeq.sciml.ai/stable/solvers/dae_solve/. Must be an AbstractArray or an AbstractSciMLOperator.
• 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 $\frac{\partial f(u,p,t)}{\partial t}$
• jac(J,u,p,t) or J=jac(u,p,t): returns $\frac{df}{du}$
• jvp(Jv,v,u,p,t) or Jv=jvp(v,u,p,t): returns the directional derivative$\frac{df}{du} v$
• vjp(Jv,v,u,p,t) or Jv=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 sized Tridiagonal matrix can be used as the prototype and integrators will specialize on this structure where possible. Non-structured sparsity patterns should use a SparseMatrixCSC with a correct sparsity pattern for the Jacobian. The default is nothing, which means a dense Jacobian.
• paramjac(pJ,u,p,t): returns the parameter Jacobian $\frac{df}{dp}$.
• syms: the symbol names for the elements of the equation. This should match u0 in size. For example, if u0 = [0.0,1.0] and syms = [:x, :y], this will apply a canonical naming to the values, allowing sol[:x] in the solution and automatically naming values in plots.
• indepsym: the canonical naming for the independent variable. Defaults to nothing, which internally uses t as the representation in any plots.
• colorvec: a color vector according to the SparseDiffTools.jl definition for the sparsity pattern of the jac_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 to nothing, 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

iip is the the optional boolean for determining whether a given function is written to be used in-place or out-of-place. In-place functions are f!(du,u,p,t) where the return is ignored and the result is expected to be mutated into the value of du. Out-of-place functions are du=f(u,p,t).

Normally this is determined automatically by looking at the method table for f and seeing the maximum number of arguments in available dispatches. For this reason, the constructor ODEFunction(f) generally works (but is type-unstable). However, for type-stability or to enforce correctness, this option is passed via ODEFunction{true}(f).

recompile: Controlling Compilation and Specialization

The recompile parameter controls whether the ODEFunction will fully specialize on the typeof(f). This causes recompilation of the solver for each new f function, but gives the maximum compiler information and runtime speed. By default recompile = true. If recompile = false, the ODEFunction uses Any type parameters for each of the functions, allowing for the reuse of compilation caches but adding a dynamic dispatch at the f call sites, potentially leading to runtime regressions.

Overriding the true default is done by passing a second type parameter after iip, for example ODEFunction{true,false}(f) is an in-place function with no recompilation specialization.

Fields

The fields of the ODEFunction type directly match the names of the inputs.

More Details on Jacobians

The following example creates an inplace ODEFunction whose jacobian is a Diagonal:

using LinearAlgebra
f = (du,u,p,t) -> du .= t .* u
jac = (J,u,p,t) -> (J[1,1] = t; J[2,2] = t; J)
jp = Diagonal(zeros(2))
fun = ODEFunction(f; jac=jac, jac_prototype=jp)

Note that the integrators will always make a deep copy of fun.jac_prototype, so there's no worry of aliasing.

In general the jacobian prototype can be anything that has mul! defined, in particular sparse matrices or custom lazy types that support mul!. A special case is when the jac_prototype is a AbstractDiffEqLinearOperator, in which case you do not need to supply jac as it is automatically set to update_coefficients!. Refer to the AbstractSciMLOperators documentation for more information on setting up time/parameter dependent operators.

Examples

Declaring Explicit Jacobians for ODEs

The most standard case, declaring a function for a Jacobian is done by overloading the function f(du,u,p,t) with an in-place updating function for the Jacobian: f_jac(J,u,p,t) where the value type is used for dispatch. For example, take the LotkaVolterra model:

function f(du,u,p,t)
  du[1] = 2.0 * u[1] - 1.2 * u[1]*u[2]
  du[2] = -3 * u[2] + u[1]*u[2]
end

To declare the Jacobian we simply add the dispatch:

function f_jac(J,u,p,t)
  J[1,1] = 2.0 - 1.2 * u[2]
  J[1,2] = -1.2 * u[1]
  J[2,1] = 1 * u[2]
  J[2,2] = -3 + u[1]
  nothing
end

Then we can supply the Jacobian with our ODE as:

ff = ODEFunction(f;jac=f_jac)

and use this in an ODEProblem:

prob = ODEProblem(ff,ones(2),(0.0,10.0))

Symbolically Generating the Functions

See the modelingtoolkitize function from ModelingToolkit.jl for automatically symbolically generating the Jacobian and more from the numerically-defined functions.

struct ODESolution{T, N, uType, uType2, DType, tType, rateType, P, A, IType, DE} <: SciMLBase.AbstractODESolution{T, N, uType}

Representation of the solution to an ordinary differential equation defined by an ODEProblem.

DESolution Interface

For more information on interacting with DESolution types, check out the Solution Handling page of the DifferentialEquations.jl documentation.

https://diffeq.sciml.ai/stable/basics/solution/

Fields

• u: the representation of the ODE solution. Given as an array of solutions, where u[i] corresponds to the solution at time t[i]. It is recommended in most cases one does not access sol.u directly and instead use the array interface described in the Solution Handling page of the DifferentialEquations.jl documentation.
• t: the time points corresponding to the saved values of the ODE solution.
• prob: the original ODEProblem that was solved.
• alg: the algorithm type used by the solver.
• destats: 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 (sol.retcode === :Success), whether it terminated due to a user-defined callback (sol.retcode === :Terminated), or whether it exited due to an error. For more details, see the return code section of the DifferentialEquations.jl documentation.

Linear ODE

\[\frac{du}{dt} = αu\]

with initial condition $u_0=\frac{1}{2}$, $α=1.01$, and solution

\[u(t) = u_0e^{αt}\]

with Float64s. The parameter is $α$

4x2 version of the Linear ODE

\[\frac{du}{dt} = αu\]

with initial condition $u_0$ as all uniformly distributed random numbers, $α=1.01$, and solution

\[u(t) = u_0e^{αt}\]

with Float64s

Linear ODE

\[\frac{du}{dt} = αu\]

with initial condition $u_0=\frac{1}{2}$, $α=1.01$, and solution

\[u(t) = u_0e^{αt}\]

with BigFloats

4x2 version of the Linear ODE

\[\frac{du}{dt} = αu\]

with initial condition $u_0$ as all uniformly distributed random numbers, $α=1.01$, and solution

\[u(t) = u_0e^{αt}\]

with BigFloats

100x100 version of the Linear ODE

\[\frac{du}{dt} = αu\]

with initial condition $u_0$ as all uniformly distributed random numbers, $α=1.01$, and solution

\[u(t) = u_0e^{αt}\]

with Float64s

4x2 version of the Linear ODE

\[\frac{du}{dt} = αu\]

with initial condition $u_0$ as all uniformly distributed random numbers, $α=1.01$, and solution

\[u(t) = u_0e^{αt}\]

on Float64. Purposefully not in-place as a test.

Fitzhugh-Nagumo (Non-stiff)

\[\frac{dv}{dt} = v - \frac{v^3}{3} - w + I_{est}\]

\[τ \frac{dw}{dt} = v + a -bw\]

with initial condition $v=w=1$

The ThreeBody problem as written by Hairer: (Non-stiff)

\[\frac{dy₁}{dt} = y₁ + 2\frac{dy₂}{dt} - \bar{μ}\frac{y₁+μ}{D₁} - μ\frac{y₁-\bar{μ}}{D₂}\]

\[\frac{dy₂}{dt} = y₂ - 2\frac{dy₁}{dt} - \bar{μ}\frac{y₂}{D₁} - μ\frac{y₂}{D₂}\]

\[D₁ = ((y₁+μ)^2 + y₂^2)^{3/2}\]

\[D₂ = ((y₁-\bar{μ})^2+y₂^2)^{3/2}\]

\[μ = 0.012277471\]

\[\bar{μ} =1-μ\]

From Hairer Norsett Wanner Solving Ordinary Differential Equations I - Nonstiff Problems Page 129

Usually solved on $t₀ = 0.0$ and $T = 17.0652165601579625588917206249$ Periodic with that setup.

Pleiades Problem (Non-stiff)

\[\frac{d^2xᵢ}{dt^2} = \sum_{j≠i} mⱼ(xⱼ-xᵢ)/rᵢⱼ\]

\[\frac{d^2yᵢ}{dt^2} = \sum_{j≠i} mⱼ(yⱼ-yᵢ)/rᵢⱼ\]

where

\[rᵢⱼ = ((xᵢ-xⱼ)^2 + (yᵢ-yⱼ)^2)^{3/2}\]

and initial conditions are

\[x₁(0) = 3\]

\[x₂(0) = 3\]

\[x₃(0) = -1\]

\[x₄(0) = -3\]

\[x₅(0) = 2\]

\[x₆(0) = -2\]

\[x₇(0) = 2\]

\[y₁(0) = 3\]

\[y₂(0) = -3\]

\[y₃(0) = 2\]

\[y₄(0) = 0\]

\[y₅(0) = 0\]

\[y₆(0) = -4\]

\[y₇(0) = 4\]

and with $\frac{dxᵢ(0)}{dt}=\frac{dyᵢ(0)}{dt}=0$ except for

\[\frac{dx₆(0)}{dt} = 1.75\]

\[\frac{dx₇(0)}{dt} = -1.5\]

\[\frac{dy₄(0)}{dt} = -1.25\]

\[\frac{dy₅(0)}{dt} = 1\]

From Hairer Norsett Wanner Solving Ordinary Differential Equations I - Nonstiff Problems Page 244

Usually solved from 0 to 3.

Van der Pol Equations

\[\frac{dx}{dt} = y\]

\[\frac{dy}{dt} = μ((1-x^2)y -x)\]

with $μ=1.0$ and $u_0=[0,\sqrt{3}]$

Non-stiff parameters.

Van der Pol Equations

\[\frac{dx}{dt} = y\]

\[\frac{dy}{dt} = μ((1-x^2)y -x)\]

with $μ=10^6$ and $u_0=[0,\sqrt{3}]$

Stiff parameters.

The Robertson biochemical reactions: (Stiff)

\[\frac{dy₁}{dt} = -k₁y₁+k₃y₂y₃\]

\[\frac{dy₂}{dt} = k₁y₁-k₂y₂^2-k₃y₂y₃\]

\[\frac{dy₃}{dt} = k₂y₂^2\]

where $k₁=0.04$, $k₂=3\times10^7$, $k₃=10^4$. For details, see:

Hairer Norsett Wanner Solving Ordinary Differential Equations I - Nonstiff Problems Page 129

Usually solved on $[0,1e11]$

Rigid Body Equations (Non-stiff)

\[\frac{dy₁}{dt} = I₁y₂y₃\]

\[\frac{dy₂}{dt} = I₂y₁y₃\]

\[\frac{dy₃}{dt} = I₃y₁y₂\]

with $I₁=-2$, $I₂=1.25$, and $I₃=-1/2$.

The initial condition is $y=[1.0;0.0;0.9]$.

From Solving Differential Equations in R by Karline Soetaert

or Hairer Norsett Wanner Solving Ordinary Differential Equations I - Nonstiff Problems Page 244

Usually solved from 0 to 20.

Hires Problem (Stiff)

It is in the form of

\[\frac{dy}{dt} = f(y)\]

with

\[ y(0)=y_0, \quad y \in ℝ^8, \quad 0 ≤ t ≤ 321.8122\]

where $f$ is defined by

$f(y) = \begin{pmatrix} −1.71y_1 & +0.43y_2 & +8.32y_3 & +0.0007y_4 & \\ 1.71y_1 & −8.75y_2 & & & \\ −10.03y_3 & +0.43y_4 & +0.035y_5 & & \\ 8.32y_2 & +1.71y_3 & −1.12y_4 & & \\ −1.745y_5 & +0.43y_6 & +0.43y_7 & & \\ −280y_6y_8 & +0.69y_4 & +1.71y_5 & −0.43y_6 & +0.69y_7 \\ 280y_6y_8 & −1.81y_7 & & & \\ −280y_6y_8 & +1.81y_7 & & & \end{pmatrix}$

Reference: demohires.pdf Notebook: Hires.ipynb

Orego Problem (Stiff)

It is in the form of $\frac{dy}{dt}=f(y), \quad y(0)=y0,$ with

\[y \in ℝ^3, \quad 0 ≤ t ≤ 360\]

where $f$ is defined by

$f(y) = \begin{pmatrix} s(y_2 - y_1(1-qy_1-y_2)) \\ (y_3 - y_2(1+y_1))/s \\ w(y_1-y_3) \end{pmatrix}$

where $s=77.27$, $w=0.161$ and $q=8.375⋅10^{-6}$.

Reference: demoorego.pdf Notebook: Orego.ipynb

Pollution Problem (Stiff)

This IVP is a stiff system of 20 non-linear Ordinary Differential Equations. It is in the form of

\[\frac{dy}{dt}=f(y)\]

with

\[y(0)=y0, \quad y \in ℝ^20, \quad 0 ≤ t ≤ 60\]

where $f$ is defined by

$f(y) = \begin{pmatrix} -\sum_{j∈{1,10,14,23,24}} r_j + \sum_{j∈{2,3,9,11,12,22,25}} r_j \\ -r_2 - r_3 - r_9 - r_12 + r_1 + r_{21} \\ -r_{15} + r_1 + r_{17} + r_{19} + r_{22} \\ -r_2 - r_{16} - r_{17} - r_{23} + r_{15} \\ -r_3 + 2r_4 + r_6 + r_7 + r_{13} + r_{20} \\ -r_6 - r_8 - r_{14} - r_{20} + r_3 + 2r_{18} \\ -r_4 - r_5 - r_6 + r_{13} \\ r_4 + r_5 + r_6 + r_7 \\ -r_7 - r_8 \\ -r_{12} + r_7 + r_9 \\ -r_9 - r_{10} + r_8 + r_{11} \\ r_9 \\ -r_{11} + r_{10} \\ -r_{13} + r_{12} \\ r_{14} \\ -r_{18} - r_{19} + r_{16} \\ -r_{20} \\ r_{20} \\ -r{21} - r_{22} - r_{24} + r_{23} + r_{25} \\ -r_{25} + r_{24} \end{pmatrix}$

with the initial condition of

\[y0 = (0, 0.2, 0, 0.04, 0, 0, 0.1, 0.3, 0.01, 0, 0, 0, 0 ,0, 0, 0, 0.007, 0, 0, 0)^T\]

Analytical Jacobian is included.

Reference: pollu.pdf Notebook: Pollution.ipynb

Nonlinear system of reactions with an analytical solution

\[\frac{dy_1}{dt} = -y_1\]

\[\frac{dy_2}{dt} = y_1 - y_2^2\]

\[\frac{dy_3}{dt} = y_2^2\]

with initial condition $y=[1;0;0]$ on a time span of $t \in (0,20)$

From

Liu, L. C., Tian, B., Xue, Y. S., Wang, M., & Liu, W. J. (2012). Analytic solution for a nonlinear chemistry system of ordinary differential equations. Nonlinear Dynamics, 68(1-2), 17-21.

The analytical solution is implemented, allowing easy testing of ODE solvers.

1D Brusselator

\[\frac{\partial u}{\partial t} = A + u^2v - (B+1)u + \alpha\frac{\partial^2 u}{\partial x^2}\]

\[\frac{\partial v}{\partial t} = Bu - u^2v + \alpha\frac{\partial^2 u}{\partial x^2}\]

and the initial conditions are

\[u(x,0) = 1+\sin(2π x)\]

\[v(x,0) = 3\]

with the boundary condition

\[u(0,t) = u(1,t) = 1\]

\[v(0,t) = v(1,t) = 3\]

From Hairer Norsett Wanner Solving Ordinary Differential Equations II - Stiff and Differential-Algebraic Problems Page 6

2D Brusselator

\[\frac{\partial u}{\partial t} = 1 + u^2v - 4.4u + \alpha(\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2}) + f(x, y, t)\]

\[\frac{\partial v}{\partial t} = 3.4u - u^2v + \alpha(\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2})\]

where

$f(x, y, t) = \begin{cases} 5 & \quad \text{if } (x-0.3)^2+(y-0.6)^2 ≤ 0.1^2 \text{ and } t ≥ 1.1 \\ 0 & \quad \text{else} \end{cases}$

and the initial conditions are

\[u(x, y, 0) = 22\cdot y(1-y)^{3/2}\]

\[v(x, y, 0) = 27\cdot x(1-x)^{3/2}\]

with the periodic boundary condition

\[u(x+1,y,t) = u(x,y,t)\]

\[u(x,y+1,t) = u(x,y,t)\]

From Hairer Norsett Wanner Solving Ordinary Differential Equations II - Stiff and Differential-Algebraic Problems Page 152

Filament PDE Discretization

Notebook: Filament.ipynb

In this problem is a real-world biological model from a paper entitled Magnetic dipole with a flexible tail as a self-propelling microdevice. It is a system of PDEs representing a Kirchhoff model of an elastic rod, where the equations of motion are given by the Rouse approximation with free boundary conditions.

Thomas' cyclically symmetric attractor equations

\[\begin{align} \frac{dx(t)}{dt} =& - b x\left( t \right) + \sin\left( y\left( t \right) \right) \\ \frac{dy(t)}{dt} =& - b y\left( t \right) + \sin\left( z\left( t \right) \right) \\ \frac{dz(t)}{dt} =& - b z\left( t \right) + \sin\left( x\left( t \right) \right) \end{align} \]

Reference

Wikipedia

Lorenz equations

\[\begin{align} \frac{dx(t)}{dt} =& \sigma \left( - x\left( t \right) + y\left( t \right) \right) \\ \frac{dy(t)}{dt} =& - y\left( t \right) + \left( \rho - z\left( t \right) \right) x\left( t \right) \\ \frac{dz(t)}{dt} =& x\left( t \right) y\left( t \right) - \beta z\left( t \right) \end{align} \]

Reference

Wikipedia

Aizawa equations

\[\begin{align} \frac{dx(t)}{dt} =& \left( - b + z\left( t \right) \right) x\left( t \right) - d y\left( t \right) \\ \frac{dy(t)}{dt} =& d x\left( t \right) + \left( - b + z\left( t \right) \right) y\left( t \right) \\ \frac{dz(t)}{dt} =& c - \frac{1}{3} \left( z\left( t \right) \right)^{3} + a z\left( t \right) - \left( 1 + e z\left( t \right) \right) \left( \left( x\left( t \right) \right)^{2} + \left( y\left( t \right) \right)^{2} \right) + \left( x\left( t \right) \right)^{3} f z\left( t \right) \end{align} \]

Reference

Dadras equations

\[\begin{align} \frac{dx(t)}{dt} =& - a x\left( t \right) + b y\left( t \right) z\left( t \right) + y\left( t \right) \\ \frac{dy(t)}{dt} =& c y\left( t \right) - x\left( t \right) z\left( t \right) + z\left( t \right) \\ \frac{dz(t)}{dt} =& d x\left( t \right) y\left( t \right) - e z\left( t \right) \end{align} \]

Reference

chen equations

\[\begin{align} \frac{dx(t)}{dt} =& a \left( - x\left( t \right) + y\left( t \right) \right) \\ \frac{dy(t)}{dt} =& c y\left( t \right) + \left( c - a \right) x\left( t \right) - x\left( t \right) z\left( t \right) \\ \frac{dz(t)}{dt} =& x\left( t \right) y\left( t \right) - b z\left( t \right) \end{align} \]

Reference

rossler equations

\[\begin{align} \frac{dx(t)}{dt} =& - y\left( t \right) - z\left( t \right) \\ \frac{dy(t)}{dt} =& a y\left( t \right) + x\left( t \right) \\ \frac{dz(t)}{dt} =& b + \left( - c + x\left( t \right) \right) z\left( t \right) \end{align} \]

Reference Wikipedia

rabinovich_fabrikant equations

\[\begin{align} \frac{dx(t)}{dt} =& b x\left( t \right) + \left( -1 + \left( x\left( t \right) \right)^{2} + z\left( t \right) \right) y\left( t \right) \\ \frac{dy(t)}{dt} =& \left( 1 - \left( x\left( t \right) \right)^{2} + 3 z\left( t \right) \right) x\left( t \right) + b y\left( t \right) \\ \frac{dz(t)}{dt} =& - 2 \left( a + x\left( t \right) y\left( t \right) \right) z\left( t \right) \end{align} \]

Reference

sprott equations

\[\begin{align} \frac{dx(t)}{dt} =& x\left( t \right) z\left( t \right) + a x\left( t \right) y\left( t \right) + y\left( t \right) \\ \frac{dy(t)}{dt} =& 1 + y\left( t \right) z\left( t \right) - \left( x\left( t \right) \right)^{2} b \\ \frac{dz(t)}{dt} =& - \left( x\left( t \right) \right)^{2} - \left( y\left( t \right) \right)^{2} + x\left( t \right) \end{align} \]

Reference

hindmarsh_rose equations

\[\begin{align} \frac{dx(t)}{dt} =& i - z\left( t \right) + \left( x\left( t \right) \right)^{2} b - \left( x\left( t \right) \right)^{3} a + y\left( t \right) \\ \frac{dy(t)}{dt} =& c - y\left( t \right) - \left( x\left( t \right) \right)^{2} d \\ \frac{dz(t)}{dt} =& r \left( s \left( - xr + x\left( t \right) \right) - z\left( t \right) \right) \end{align} \]

Reference