PositiveIntegrators.jl API

ConservativePDSProblem(P, u0, tspan, p = NullParameters();
                       p_prototype = nothing,
                       analytic = nothing,
                       std_rhs = nothing)

A structure describing a conservative system of ordinary differential equation in form of a production-destruction system (PDS). P denotes the function defining the production matrix $P$. The diagonal of $P$ contains production terms without destruction counterparts. u0 is the vector of initial conditions and tspan the time span (t_initial, t_final) of the problem. The optional argument p can be used to pass additional parameters to the function P.

The function P can be given either in the out-of-place form with signature production_terms = P(u, p, t) or the in-place form P(production_terms, u, p, t).

Keyword arguments:

• p_prototype: If P is given in in-place form, p_prototype or copies thereof are used to store evaluations of P. If p_prototype is not specified explicitly and P is in-place, then p_prototype will be internally set to zeros(eltype(u0), (length(u0), length(u0))).
• analytic: The analytic solution of a PDS must be given in the form f(u0,p,t). Specifying the analytic solution can be useful for plotting and convergence tests.
• std_rhs: The standard ODE right-hand side evaluation function callable as du = std_rhs(u, p, t) for the out-of-place form and as std_rhs(du, u, p, t) for the in-place form. Solvers that do not rely on the production-destruction representation of the ODE, will use this function instead to compute the solution. If not specified, a default implementation calling P is used.

References

• Hans Burchard, Eric Deleersnijder, and Andreas Meister. "A high-order conservative Patankar-type discretisation for stiff systems of production-destruction equations." Applied Numerical Mathematics 47.1 (2003): 1-30. DOI: 10.1016/S0168-9274(03)00101-6

PDSProblem(P, D, u0, tspan, p = NullParameters();
           p_prototype = nothing,
           analytic = nothing,
           std_rhs = nothing)

A structure describing a system of ordinary differential equations in form of a production-destruction system (PDS). P denotes the function defining the production matrix $P$. The diagonal of $P$ contains production terms without destruction counterparts. D is the function defining the vector of destruction terms $D$ without production counterparts. u0 is the vector of initial conditions and tspan the time span (t_initial, t_final) of the problem. The optional argument p can be used to pass additional parameters to the functions P and D.

The functions P and D can be used either in the out-of-place form with signature production_terms = P(u, p, t) or the in-place form P(production_terms, u, p, t).

Keyword arguments:

• p_prototype: If P is given in in-place form, p_prototype or copies thereof are used to store evaluations of P. If p_prototype is not specified explicitly and P is in-place, then p_prototype will be internally set to zeros(eltype(u0), (length(u0), length(u0))).
• analytic: The analytic solution of a PDS must be given in the form f(u0,p,t). Specifying the analytic solution can be useful for plotting and convergence tests.
• std_rhs: The standard ODE right-hand side evaluation function callable as du = std_rhs(u, p, t) for the out-of-place form and as std_rhs(du, u, p, t) for the in-place form. Solvers that do not rely on the production-destruction representation of the ODE, will use this function instead to compute the solution. If not specified, a default implementation calling P and D is used.

References

• Hans Burchard, Eric Deleersnijder, and Andreas Meister. "A high-order conservative Patankar-type discretisation for stiff systems of production-destruction equations." Applied Numerical Mathematics 47.1 (2003): 1-30. DOI: 10.1016/S0168-9274(03)00101-6

prob_pds_bertolazzi

Positive and conservative autonomous nonlinear PDS

\[\begin{aligned} \mathbf{u}'=\begin{pmatrix}2 &-1 &-1\\-1 &2 &-1\\-1& -1& 2\end{pmatrix}\begin{pmatrix}5u_2u_3/(10^{-2} + (u_2u_3)^2) + u_2u_3/(10^{-16} + u_2u_3(10^{-8} + u_2u_3))\\ 10u_1u_3^2\\ 0.1(u_3 - u_2 - 2.5)^2u_1u_2\end{pmatrix} \end{aligned}\]

with initial value $\mathbf{u}_0 = (0.0, 1.0, 2.0)^T$ and time domain $(0.0, 1.0)$. There is one independent linear invariant, e.g. $u_1+u_2+u_3 = 3.0$.

References

• Enrico Bertolazzi. "Positive and conservative schemes for mass action kinetics." Computers and Mathematics with Applications 32 (1996): 29-43. DOI: 10.1016/0898-1221(96)00142-3

prob_pds_brusselator

Positive and conservative autonomous nonlinear PDS

\[\begin{aligned} u_1' &= -u_1,\\ u_2' &= -u_2u_5,\\ u_3' &= u_2u_5,\\ u_4' &= u_5,\\ u_5' &= u_1 - u_2u_5 + u_5^2u_6 - u_5,\\ u_6' &= u_2u_5 - u_5^2u_6, \end{aligned}\]

with initial value $\mathbf{u}_0 = (10.0, 10.0, 0.0, 0.0, 0.1, 0.1)^T$ and time domain $(0.0, 20.0)$. There are two independent linear invariants, e.g. $u_1+u_4+u_5+u_6 = 10.2$ and $u_2+u_3 = 10.0$.

References

• Luca Bonaventura, and Alessandro Della Rocca. "Unconditionally Strong Stability Preserving Extensions of the TR-BDF2 Method." Journal of Scientific Computing 70 (2017): 859 - 895. DOI: 10.1007/s10915-016-0267-9

prob_pds_linmod

Positive and conservative autonomous linear PDS

\[\begin{aligned} u_1' &= u_2 - 5u_1,\\ u_2' &= 5u_1 - u_2, \end{aligned}\]

with initial value $\mathbf{u}_0 = (0.9, 0.1)^T$ and time domain $(0.0, 2.0)$. There is one independent linear invariant, e.g. $u_1+u_2 = 1$.

References

• Hans Burchard, Eric Deleersnijder, and Andreas Meister. "A high-order conservative Patankar-type discretisation for stiff systems of production-destruction equations." Applied Numerical Mathematics 47.1 (2003): 1-30. DOI: 10.1016/S0168-9274(03)00101-6

prob_pds_linmod_inplace

Same as prob_pds_linmod but with in-place computation.

prob_pds_minmapk

Positive and nonconservative autonomous nonlinear PDS

\[\begin{aligned} u_1' &= k_6u_6-k_7u_1-k_1u_1u_2 +k_2u_4,\\ u_2' &= k_5u_3-k_1u_2u_2,\\ u_3' &= k_2u_4-k_3u_1u_3+k_4u_5-k_5u_3,\\ u_4' &= k_1u_1u_2-k_2u_4,\\ u_5' &= k_3u_1u_3-k_4u_5,\\ u_6' &= k_7u_1-k_6u_6, \end{aligned}\]

with constants

\[\begin{aligned} k_1 &=\frac{100}{3}, & k_2 &=\frac{1}{3}, & k_3 &=50,\\ k_4 &=0.5, & k_5 &=\frac{10}{3} , & k_6 &= 0.1,\\ k_7 &= 0.1. \end{aligned}\]

The initial value is $\mathbf{u}_0 = (0.1, 0.175, 0.15, 1.15, 0.81, 0.5)^T$ and the time domain $(0, 200)$. There are two independent linear invariants, e.g. $u_1+u_4+u_6=1.75$ and $u_2+u_3+u_4+u_5 =2.285$.

References

• Sergio Blanes, Arieh Iserles, and Shev Macnamara. "Positivity preserving methods for ordinary differential equations." ESAIM: Mathematical Modelling and Numerical Analysis 56 (2022): 1843–1870. DOI: 10.1051/m2an/2022042
• Otto Hadač, František Muzika, Vladislav Nevoral, Michal Přibyl, and Igor Schreiber "Minimal oscillating subnetwork in the Huang-Ferrell model of the MAPK cascade." PLoS ONE 12 (2017): e0178457. DOI: 10.1371/journal.pone.0178457

prob_pds_nonlinmod

Positive and conservative autonomous nonlinear PDS

\[\begin{aligned} u_1' &= -\frac{u_1u_2}{u_1 + 1.0},\\ u_2' &= \frac{u_1u_2}{u_1 + 1.0} - 0.3u_2,\\ u_3' &= 0.3 u_2, \end{aligned}\]

with initial value $\mathbf{u}_0 = (9.98, 0.01, 0.01)^T$ and time domain $(0.0, 30.0)$. There is one independent linear invariant, e.g. $u_1+u_2+u_3 = 10.0$.

References

• Hans Burchard, Eric Deleersnijder, and Andreas Meister. "A high-order conservative Patankar-type discretisation for stiff systems of production-destruction equations." Applied Numerical Mathematics 47.1 (2003): 1-30. DOI: 10.1016/S0168-9274(03)00101-6

prob_pds_npzd

Positive and conservative autonomous nonlinear PDS

\[\begin{aligned} u_1' &= 0.01u_2 + 0.01u_3 + 0.003u_4 - \frac{u_1u_2}{0.01 + u_1},\\ u_2' &= \frac{u_1u_2}{0.01 + u_1}- 0.01u_2 - 0.5( 1 - e^{-1.21u_2^2})u_3 - 0.05u_2,\\ u_3' &= 0.5(1 - e^{-1.21u_2^2})u_3 - 0.01u_3 - 0.02u_3,\\ u_4' &= 0.05u_2 + 0.02u_3 - 0.003u_4 \end{aligned}\]

with initial value $\mathbf{u}_0 = (8.0, 2.0, 1.0, 4.0)^T$ and time domain $(0.0, 10.0)$. There is one independent linear invariant, e.g. $u_1+u_2+u_3+u_4 = 15.0$.

References

• Hans Burchard, Eric Deleersnijder, and Andreas Meister. "Application of modified Patankar schemes to stiff biogeochemical models for the water column." Ocean Dynamics 55 (2005): 326-337. DOI: 10.1007/s10236-005-0001-x

prob_pds_robertson

Positive and conservative autonomous nonlinear PDS

\[\begin{aligned} u_1' &= -0.04u_1+10^4 u_2u_3,\\ u_2' &= 0.04u_1-10^4 u_2u_3-3⋅10^7 u_2^2,\\ u_3' &= 3⋅10^7 u_2^2, \end{aligned}\]

with initial value $\mathbf{u}_0 = (1.0, 0.0, 0.0)^T$ and time domain $(0.0, 10^{11})$. There is one independent linear invariant, e.g. $u_1+u_2+u_3 = 1.0$.

References

• Ernst Hairer, Gerd Wanner. "Solving Ordinary Differential Equations II - Stiff and Differential-Algebraic Problems." 2nd Edition, Springer (2002): Section IV.1.

prob_pds_sir

Positive and conservative autonomous nonlinear PDS

\[\begin{aligned} u_1' &= -2u_1u_2,\\ u_2' &= 2u_1u_2 - u_2,\\ u_3' &= u_2, \end{aligned}\]

with initial value $\mathbf{u}_0 = (0.99, 0.005, 0.005)^T$ and time domain $(0.0, 20.0)$. There is one independent linear invariant, e.g. $u_1+u_2+u_3 = 1.0$.

References

• Ronald E. Mickens, and Talitha M. Washington. "NSFD discretizations of interacting population models satisfying conservation laws." Computers and Mathematics with Applications 66 (2013): 2307-2316. DOI: 10.1016/j.camwa.2013.06.011

prob_pds_stratreac

Positive and nonconservative autonomous nonlinear PDS

\[\begin{aligned} u_1' &= r_5 - r_6 - r_7,\\ u_2' &= 2r_1 - r_2 + r_3 - r_4 + r_6 - r_9 + r_{10} - r_{11},\\ u_3' &= r_2 - r_3 - r_4 - r_5 - r_7 - r_8,\\ u_4' &= -r_1 -r_2 + r_3 + 2r_4+r_5+2r_7+r_8+r_9,\\ u_5' &= -r_8+r_9+r_{10}-r_{11},\\ u_6' &= r_8-r_9-r_{10}+r_{11}, \end{aligned}\]

with reaction rates

\[\begin{aligned} r_1 &=2.643⋅ 10^{-10}σ^3 u_4, & r_2 &=8.018⋅10^{-17}u_2 u_4 , & r_3 &=6.12⋅10^{-4}σ u_3,\\ r_4 &=1.567⋅10^{-15}u_3 u_2 , & r_5 &= 1.07⋅ 10^{-3}σ^2u_3, & r_6 &= 7.11⋅10^{-11}⋅ 8.12⋅10^6 u_1,\\ r_7 &= 1.2⋅10^{-10}u_1 u_3, & r_8 &= 6.062⋅10^{-15}u_3 u_5, & r_9 &= 1.069⋅10^{-11}u_6 u_2,\\ r_{10} &= 1.289⋅10^{-2}σ u_6, & r_{11} &= 10^{-8}u_5 u_2, \end{aligned}\]

where

\[\begin{aligned} T &= t/3600 \mod 24,\quad T_r=4.5,\quad T_s = 19.5,\\ σ(T) &= \begin{cases}1, & T_r≤ T≤ T_s,\\0, & \text{otherwise}.\end{cases} \end{aligned}\]

The initial value is $\mathbf{u}_0 = (9.906⋅10^1, 6.624⋅10^8, 5.326⋅10^{11}, 1.697⋅10^{16}, 4⋅10^6, 1.093⋅10^9)^T$ and the time domain $(4.32⋅ 10^{4}, 3.024⋅10^5)$. There are two independent linear invariants, e.g. $u_1+u_2+3u_3+2u_4+u_5+2u_6=(1,1,3,2,1,2)\cdot\mathbf{u}_0$ and $u_5+u_6 = 1.097⋅10^9$.

References

• Stephan Nüsslein, Hendrik Ranocha, and David I. Ketcheson. "Positivity-preserving adaptive Runge-Kutta methods." Communications in Applied Mathematics and Computer Science 16 (2021): 155-179. DOI: 10.2140/camcos.2021.16.155

MPE([linsolve = ..., small_constant = ...])

The first-order modified Patankar-Euler algorithm for production-destruction systems. This one-step, one-stage method is first-order accurate, unconditionally positivity-preserving, and linearly implicit.

The scheme was introduced by Burchard et al. for conservative production-destruction systems. For nonconservative production–destruction systems we use the straight forward extension

$u_i^{n+1} = u_i^n + Δt \sum_{j, j≠i} \biggl(p_{ij}^n \frac{u_j^{n+1}}{u_j^n}-d_{ij}^n \frac{u_i^{n+1}}{u_i^n}\biggr) + {\Delta}t p_{ii}^n - Δt d_{ii}^n\frac{u_i^{n+1}}{u_i^n}$,

where $p_{ij}^n = p_{ij}(t^n,\mathbf u^n)$ and $d_{ij}^n = d_{ij}(t^n,\mathbf u^n)$.

The modified Patankar-Euler method requires the special structure of a PDSProblem or a ConservativePDSProblem.

You can optionally choose the linear solver to be used by passing an algorithm from LinearSolve.jl as keyword argument linsolve. You can also choose the parameter small_constant which is added to all Patankar-weight denominators to avoid divisions by zero. You can pass a value explicitly, otherwise small_constant is set to floatmin of the floating point type used.

The current implementation only supports fixed time steps.

References

• Hans Burchard, Eric Deleersnijder, and Andreas Meister. "A high-order conservative Patankar-type discretisation for stiff systems of production-destruction equations." Applied Numerical Mathematics 47.1 (2003): 1-30. DOI: 10.1016/S0168-9274(03)00101-6

MPRK22(α; [linsolve = ..., small_constant = ...])

A family of second-order modified Patankar-Runge-Kutta algorithms for production-destruction systems. Each member of this family is an adaptive, one-step, two-stage method which is second-order accurate, unconditionally positivity-preserving, and linearly implicit. In this implementation the stage-values are conservative as well. The parameter α is described by Kopecz and Meister (2018) and studied by Izgin, Kopecz and Meister (2022) as well as Torlo, Öffner and Ranocha (2022).

This method supports adaptive time stepping, using the Patankar-weight denominators $σ_i$, see Kopecz and Meister (2018), as first order approximations to estimate the error.

The scheme was introduced by Kopecz and Meister for conservative production-destruction systems. For nonconservative production–destruction systems we use a straight forward extension analogous to MPE.

This modified Patankar-Runge-Kutta method requires the special structure of a PDSProblem or a ConservativePDSProblem.

You can optionally choose the linear solver to be used by passing an algorithm from LinearSolve.jl as keyword argument linsolve. You can also choose the parameter small_constant which is added to all Patankar-weight denominators to avoid divisions by zero. You can pass a value explicitly, otherwise small_constant is set to floatmin of the floating point type used.

References

• Hans Burchard, Eric Deleersnijder, and Andreas Meister. "A high-order conservative Patankar-type discretisation for stiff systems of production-destruction equations." Applied Numerical Mathematics 47.1 (2003): 1-30. DOI: 10.1016/S0168-9274(03)00101-6
• Stefan Kopecz and Andreas Meister. "On order conditions for modified Patankar-Runge-Kutta schemes." Applied Numerical Mathematics 123 (2018): 159-179. DOI: 10.1016/j.apnum.2017.09.004
• Thomas Izgin, Stefan Kopecz, and Andreas Meister. "On Lyapunov stability of positive and conservative time integrators and application to second order modified Patankar-Runge-Kutta schemes." ESAIM: Mathematical Modelling and Numerical Analysis 56.3 (2022): 1053-1080. DOI: 10.1051/m2an/2022031
• Davide Torlo, Philipp Öffner, and Hendrik Ranocha. "Issues with positivity-preserving Patankar-type schemes." Applied Numerical Mathematics 182 (2022): 117-147. DOI: 10.1016/j.apnum.2022.07.014

SSPMPRK22(α, β; [linsolve = ..., small_constant = ...])

A family of second-order modified Patankar-Runge-Kutta algorithms for production-destruction systems. Each member of this family is an adaptive, one-step, two-stage method which is second-order accurate, unconditionally positivity-preserving, and linearly implicit. The parameters α and β are described by Huang and Shu (2019) and studied by Huang, Izgin, Kopecz, Meister and Shu (2023). The difference to MPRK22 is that this method is based on the SSP formulation of an explicit second-order Runge-Kutta method. This family of schemes contains the MPRK22 family, where MPRK22(α) = SSMPRK22(0, α) applies.

This method supports adaptive time stepping, using the first order approximations $(σ_i - u_i^n) / τ + u_i^n$ with $τ=1+(α_{21}β_{10}^2)/(β_{20}+β_{21})$, see (2.7) in Huang and Shu (2019), to estimate the error.

The scheme was introduced by Huang and Shu for conservative production-destruction systems. For nonconservative production–destruction systems we use the straight forward extension analogous to MPE.

This modified Patankar-Runge-Kutta method requires the special structure of a PDSProblem or a ConservativePDSProblem.

You can optionally choose the linear solver to be used by passing an algorithm from LinearSolve.jl as keyword argument linsolve. You can also choose the parameter small_constant which is added to all Patankar-weight denominators to avoid divisions by zero. You can pass a value explicitly, otherwise small_constant is set to floatmin of the floating point type used.

References

• Juntao Huang and Chi-Wang Shu. "Positivity-Preserving Time Discretizations for Production–Destruction Equations with Applications to Non-equilibrium Flows." Journal of Scientific Computing 78 (2019): 1811–1839 DOI: 10.1007/s10915-018-0852-1
• Juntao Huang, Thomas Izgin, Stefan Kopecz, Andreas Meister and Chi-Wang Shu. "On the stability of strong-stability-preserving modified Patankar-Runge-Kutta schemes." ESAIM: Mathematical Modelling and Numerical Analysis 57 (2023):1063–1086 DOI: 10.1051/m2an/2023005

MPRK43I(α, β; [linsolve = ..., small_constant = ...])

A family of third-order modified Patankar-Runge-Kutta schemes for production-destruction systems, which is based on the two-parameter family of third order explicit Runge–Kutta schemes. Each member of this family is an adaptive, one-step method with four-stages which is third-order accurate, unconditionally positivity-preserving, conservative and linearly implicit. In this implementation the stage-values are conservative as well. The parameters α and β must be chosen such that the Runge–Kutta coefficients are nonnegative, see Kopecz and Meister (2018) for details.

These methods support adaptive time stepping, using the Patankar-weight denominators $σ_i$, see Kopecz and Meister (2018), as second order approximations to estimate the error.

The scheme was introduced by Kopecz and Meister for conservative production-destruction systems. For nonconservative production–destruction systems we use the straight forward extension analogous to MPE.

These modified Patankar-Runge-Kutta methods require the special structure of a PDSProblem or a ConservativePDSProblem.

You can optionally choose the linear solver to be used by passing an algorithm from LinearSolve.jl as keyword argument linsolve. You can also choose the parameter small_constant which is added to all Patankar-weight denominators to avoid divisions by zero. You can pass a value explicitly, otherwise small_constant is set to floatmin of the floating point type used.

References

• Stefan Kopecz and Andreas Meister. "Unconditionally positive and conservative third order modified Patankar–Runge–Kutta discretizations of production–destruction systems." BIT Numerical Mathematics 58 (2018): 691–728. DOI: 10.1007/s10543-018-0705-1

MPRK43II(γ; [linsolve = ..., small_constant = ...])

A family of third-order modified Patankar-Runge-Kutta schemes for production-destruction systems, which is based on the one-parameter family of third order explicit Runge–Kutta schemes with non-negative Runge–Kutta coefficients. Each member of this family is an adaptive, one-step method with four stages which is third-order accurate, unconditionally positivity-preserving, conservative and linearly implicit. In this implementation the stage-values are conservative as well. The parameter γ must satisfy 3/8 ≤ γ ≤ 3/4. Further details are given in Kopecz and Meister (2018).

This method supports adaptive time stepping, using the Patankar-weight denominators $σ_i$, see Kopecz and Meister (2018), as second order approximations to estimate the error.

The scheme was introduced by Kopecz and Meister for conservative production-destruction systems. For nonconservative production–destruction systems we use the straight forward extension analogous to MPE.

These modified Patankar-Runge-Kutta methods require the special structure of a PDSProblem or a ConservativePDSProblem.

You can optionally choose the linear solver to be used by passing an algorithm from LinearSolve.jl as keyword argument linsolve. You can also choose the parameter small_constant which is added to all Patankar-weight denominators to avoid divisions by zero. To display the default value for data type type evaluate MPRK43II(gamma).small_constant_function(type), where type can be, e.g., Float64.

References

• Stefan Kopecz and Andreas Meister. "Unconditionally positive and conservative third order modified Patankar–Runge–Kutta discretizations of production–destruction systems." BIT Numerical Mathematics 58 (2018): 691–728. DOI: 10.1007/s10543-018-0705-1

SSPMPRK43([linsolve = ..., small_constant = ...])

A third-order modified Patankar-Runge-Kutta algorithm for production-destruction systems. This scheme is a one-step, four-stage method which is third-order accurate, unconditionally positivity-preserving, and linearly implicit. The scheme is described by Huang, Zhao and Shu (2019) and studied by Huang, Izgin, Kopecz, Meister and Shu (2023). The difference to MPRK43I or MPRK43II is that this method is based on the SSP formulation of an explicit third-order Runge-Kutta method.

The scheme was introduced by Huang, Zhao and Shu for conservative production-destruction systems. For nonconservative production–destruction systems we use the straight forward extension analogous to MPE.

This modified Patankar-Runge-Kutta method requires the special structure of a PDSProblem or a ConservativePDSProblem.

You can optionally choose the linear solver to be used by passing an algorithm from LinearSolve.jl as keyword argument linsolve. You can also choose the parameter small_constant which is added to all Patankar-weight denominators to avoid divisions by zero. To display the default value for data type type evaluate SSPMPRK43. small_constant_function(type), where type can be, e.g., Float64.

The current implementation only supports fixed time steps.

References

• Juntao Huang, Weifeng Zhao and Chi-Wang Shu. "A Third-Order Unconditionally Positivity-Preserving Scheme for Production–Destruction Equations with Applications to Non-equilibrium Flows." Journal of Scientific Computing 79 (2019): 1015–1056 DOI: 10.1007/s10915-018-0881-9
• Juntao Huang, Thomas Izgin, Stefan Kopecz, Andreas Meister and Chi-Wang Shu. "On the stability of strong-stability-preserving modified Patankar-Runge-Kutta schemes." ESAIM: Mathematical Modelling and Numerical Analysis 57 (2023):1063–1086 DOI: 10.1051/m2an/2023005

isnegative(sol::ODESolution)

Returns true if sol contains negative elements.

See also isnonnegative.

isnonnegative(u)

Negation of isnegative.