Library

PDMPProblem(F, R, Delta, nu, xc0, xd0, p, tspan)

Create a PDMP problem to be simulated. To define a PDMP Problem, you first need to give the function F and the initial condition xc0 which define an ODE: dxc/dt = F(xc(t),xd(t),p,t). Jumps are defined as a Jump process which changes states at some rate R. Note, that in between jumps, xd is constant but xc is allowed to evolve.

Arguments

• F: inplace function F(du, u, p, t) representing the vector field
• R: the function to compute the transition rates. It should be specified in-place as R(rate::AbstractVector, xc, xd, p, t, issum::Bool) where it mutates rate. Note that a boolean issum is provided and the behavior of R should be as follows
  1. if issum == true, we only require R to return the total rate, e.g. sum(rate). We use this formalism because sometimes you can compute the sum without mutating rate.
  2. if issum == false, R must populate rate with the updated rates

We then need to provide the way the jumps affect the state variable. There are two possible ways here: 1. either give a transition matrix nu: it will only affect the discrete component xd and leave xc unaffected. 2. give a function to implement jumps Delta(xc, xd, parms, t, ind_reaction::Int64) where you can mutate xc,xd or parms. The argument ind_reaction is the index of the reaction at which the jump occurs. See examples/pdmp_example_eva.jl for an example.

• Delta [Optional]: the function to effect the jumps
• nu [Optional]: the transition matrix
• xc0: the initial condition of the continuous part
• xd0: the initial condition of the discrete part
• p: the parameters to be provided to the functions F, R, Delta
• tspan: The timespan for the problem.

Constructors

• PDMPProblem(F,R,Delta,nu,xc0,xd0,p,tspan)
• PDMPProblem(F,R,nu,xc0,xd0,p,tspan) when ones does not want to provide the function Delta
• PDMPProblem(F,R,Delta,reaction_number::Int64,xc0,xd0,p,tspan) when ones does not want to provide the transition matrix nu. The length reaction_number of the rate vector must then be provided.

We also provide a wrapper to JumpProcesses.jl. This is quite similar to how a JumpProblem would be created.

• PDMPProblem(prob, jumps...) where prob can be an ODEProblem. For an example, please consider example/examplediffeqjumpwrapper.jl.

solve(problem::PDMPProblem, algo; verbose = false, n_jumps = Inf64, save_positions = (false, true), reltol = 1e-7, abstol = 1e-9, save_rate = false, finalizer = finalize_dummy, kwargs...)

Simulate the PDMP problem using the CHV algorithm.

Arguments

• problem::PDMPProblem
• alg can be CHV(ode) (for the CHV algorithm), Rejection(ode) for the Rejection algorithm and RejectionExact() for the rejection algorithm in case the flow in between jumps is known analytically. In this latter case, prob.F is used for the specification of the Flow. The ODE solver ode can be any solver of DifferentialEquations.jl like Tsit5() for example or anyone of the list [:cvode, :lsoda, :adams, :BDF, :euler]. Indeed, the package implement an iterator interface which does not work yet with ode = LSODA(). In order to have access to the ODE solver LSODA(), one should use ode = :lsoda.
• verbose display information during simulation
• n_jumps maximum number of jumps to be computed
• save_positions which jump position to record, pre-jump (savepositions[1] = true) and/or post-jump (savepositions[2] = true).
• reltol: relative tolerance used in the ODE solver
• abstol: absolute tolerance used in the ODE solver
• ind_save_c: which indices of xc should be saved
• ind_save_d: which indices of xd should be saved
• save_rate = true: requires the solver to save the total rate. Can be useful when estimating the rate bounds in order to use the Rejection algorithm as a second try.
• X_extended = zeros(Tc, 1 + 1): (advanced use) options used to provide the shape of the extended array in the CHV algorithm. Can be useful in order to use StaticArrays.jl for example.
• finalizer = finalize_dummy: allows the user to pass a function finalizer(rate, xc, xd, p, t) which is called after each jump. Can be used to overload / add saving / plotting mechanisms.
• kwargs keyword arguments passed to the ODE solver (from DifferentialEquations.jl)

Same as the solve for CHV(::DiffEqBase.DEAlgorithm) but for CHV(::Symbol). This is an old implementation of the CHV algorithm which can be used with :lsoda. For all other solvers, use the the new solver.