Internal API
Index Handlers
Index Handler Interface
User-defined index handlers should inherit from AbstractIndexHandler and implement the following methods:
For matrix conversions they should additionally implement:
FiniteStateProjection.singleindices — Functionsingleindices(idxhandler::AbstractIndexHandler, arr)Returns all indices I in arr. Defaults to CartesianIndices, but can be overloaded for arbitrary index handlers.
FiniteStateProjection.pairedindices — Functionpairedindices(idxhandler::AbstractIndexHandler, arr, shift::CartesianIndex)Returns all pairs of indices (I .- shift, I) in arr.
FiniteStateProjection.getsubstitutions — Functiongetsubstitutions(idxhandler::AbstractIndexHandler, rs::ReactionSystem; state_sym::Symbol)Returns a dict of the form S_i => f_i(state_sym), where each f_i is an expression for the abundance of species S_i in terms of the state variable state_sym.
FiniteStateProjection.build_rhs_header — Functionbuild_rhs_header(sys::FSPSystem)Return initialisation code for the RHS function, unpacking the parameters p supplied by DifferentialEquations. The default implementation just unpacks parameters from p.
See also: unpackparams, build_rhs
Base.LinearIndices — TypeLinearIndices(idxhandler::AbstractIndexHandler, arr)Returns an object lind which converts indices returned from singleindices and pairedindices to linear indices compatible with vec via lind[idx_cart] = idx_lin. The indices are related via
arr[idx_cart] == vec(idxhandler, arr)[idx_lin]See also: vec
Built-in implementations
FiniteStateProjection.getsubstitutions — Methodgetsubstitutions(sys::FSPSystem{DefaultIndexHandler}; state_sym::Symbol)::DictDefines the abundance of species $S_i$ to be state_sym[i] - offset.
Function Building
FiniteStateProjection.build_rhs — Functionbuild_rhs(sys::FSPSystem)Builds the function f(du,u,p,t) that defines the right-hand side of the CME for use with DifferentialEquations.jl.
FiniteStateProjection.unpackparams — Functionunpackparams(sys::FSPSystem, psym::Symbol)Returns code unpacking the parameters of the system from the symbol psym in the form (p1, p2, ...) = psym. This should be called in all overloads of build_rhs_header. It is assumed that the variable psym represents an AbstractVector.
See also: build_rhs_header, build_rhs
FiniteStateProjection.build_ratefuncs — Functionbuild_ratefuncs(rs, ih; state_sym::Symbol, combinatoric_ratelaw::Bool)::VectorReturn the rate functions converted to Julia expressions in the state variable state_sym. Abundances of the species are computed using getsubstitutions.
See also: getsubstitutions, build_rhs
FiniteStateProjection.build_rhs_firstpass — Functionbuild_rhs_firstpass(sys::FSPSystem)Return code for the first pass of the RHS function, for the time-dependent FSP. Goes through all reactions and computes the negative part of the CME (probability flowing out of states). This is a simple array traversal and can be done in one go for all reactions.
See also: build_rhs
FiniteStateProjection.build_rhs_secondpass — Functionbuild_rhs_secondpass(sys::FSPSystem)Return code for the second pass of the RHS function. Goes through all reactions and computes the positive part of the CME (probability flowing into states). This requires accessing du and u at different locations depending on the net stoichiometries. In order to reduce random memory access reactions are processed one by one.
See also: build_rhs
Steady-State Functions
FiniteStateProjection.build_rhs_ss — Functionbuild_rhs_ss(sys::FSPSystem)Builds the function f(du,u,p,t) that defines the right-hand side of the CME for use with SteadyStateProblems.
FiniteStateProjection.build_rhs_singlepass_ss — Functionbuild_rhs_singlepass_ss(sys::FSPSystem)Return code for the RHS function in a single pass, for use with steady state problems. Transitions out of the truncated state space are ignored.
See also: build_rhs_ss