Solvers

solve(bp::T, alg::A; kwargs...) where {T<:ConstantRateBranchingProblem, A<:Union{SciMLBase.AbstractSciMLAlgorithm,Nothing}}

Solve a branching stochastic process with constant branching rate defined by the ConstantRateBranchingProblem bp. The positional argument alg and optional keyword arguments kwargs... are passed to the solver used to sample trajectories of the underlying SDE problem.

When the keyword argument reduction is provided (e.g. reduction=sum), the on-the-fly solver solve_and_reduce is used instead of solve_and_split, and a ReducedBranchingProcessSolution is returned directly without storing any BranchingProcessNodes. The keyword arguments transform (default: identity) and output_dt (default: 0.01) control the transform applied to particle values and the time-grid spacing of the reduced solution, respectively. Any remaining keyword arguments (including dt) are passed to the underlying SDE/jump process solver. See solve_and_reduce for details on supported reduction functions.

When reduction is not provided (the default), returns a BranchingProcessSolution containing the problem definition and the resulting tree structure. If bp.ndim is 1, 2, or 3, spatial grid positions are automatically assigned to all BranchingProcessNodes using tissue_growth! after the tree is built.

See also: ConstantRateBranchingProblem, solve_and_split, solve_and_reduce, tissue_growth!, common solver options

solve_and_split(bp::ConstantRateBranchingProblem, alg=nothing; kwargs...)

Recursively solve a branching stochastic process defined by the ConstantRateBranchingProblem bp. The positional argument alg and optional keyword arguments kwargs... are passed to the solver used to sample the trajectory of each particle.

The timespan of bp.prob defines the total time interval for the branching process. A lifetime for the first particle is sampled from the provided lifetime distribution. If the lifetime is larger than the total time interval, the problem is solved until the end of the original interval and a solution node is returned without children. If the lifetime is smaller than the total time interval, the problem is solved until the sampled lifetime, and a solution node is returned with recursively solved children for the remaining time interval.

Returns a BranchingProcessNode representing the tree structure.

See also: ConstantRateBranchingProblem, remake, sample_lifetime, common solver options

solve_and_split(prob::T, lifetime::L, nchild::O, alg::A; kwargs...) where {T<:SciMLBase.AbstractDEProblem, L<:UnivariateDistribution, O<:Union{Integer,DiscreteUnivariateDistribution}, A<:Union{SciMLBase.AbstractSciMLAlgorithm,Nothing}}

Recursively solve a branching stochastic process where the single-particle dynamics is defined by the SDE problem prob, the lifetime distribution is lifetime, and the number of children nchild of each particle is either a non-negative integer or a discrete distribution from which the number of children is sampled. The positional argument alg and optional keyword arguments kwargs... are passed to the solver used to sample the trajectory of each particle.

The timespan of the problem prob defines the total time interval for the branching process. A lifetime for the first particle is sampled from the provided lifetime distribution. If the lifetime is larger than the total time interval, the problem is solved until the end of the original interval and a solution node is returned without children. If the lifetime is smaller than the total time interval, the problem is solved until the sampled lifetime, and a solution node is returned with recursively solved children for the remaining time interval.

Returns a BranchingProcessNode representing the tree structure.

See also: SDE problems, sample_lifetime, common solver options

solve_and_reduce(bp::ConstantRateBranchingProblem, alg=nothing; reduction=sum, transform=identity, output_dt=0.01, kwargs...)

Recursively solve a branching process and construct a reduced time series on the fly, without storing any BranchingProcessNodes.

Unlike reduce_tree, which requires first sampling and storing a complete tree (whose size grows exponentially with the length of the time span), this function builds the reduced time series incrementally during the recursive solve. It exploits two key properties of branching processes:

1. A particle cannot contribute to the reduced time series outside its own time span.
2. Given their initial value (inherited from their mother particle), particles evolve independently of each other.

For every particle that is solved, its values are added to the reduced time series at the relevant time steps using an incremental binary operation derived from reduction. This approach works for reductions that form a monoid under element-wise application: sum, prod, maximum, and minimum (or the equivalent strings "sum", "prod", "max", "min").

Arguments

• bp: A ConstantRateBranchingProblem defining the single-particle dynamics, lifetime distribution, and offspring distribution.
• alg: The solver algorithm passed to the underlying SDE/jump process solver (optional).

Keyword Arguments

• reduction: The reduction function. Must be one of sum, prod, maximum, minimum (Julia built-ins), or the equivalent strings "sum", "prod", "max", "min". Default: sum.
• transform: A function applied to each particle's values before they are combined. Default: identity.
• output_dt: Spacing of the output time grid for the reduced solution. Default: 0.01.
• Additional keyword arguments (e.g. dt, saveat, reltol, abstol) are passed to the underlying SDE/jump process solver.

Returns

A ReducedBranchingProcessSolution containing the reduced time series.

Examples

using BranchingProcesses, StochasticDiffEq

f(u, p, t) = 0.0
g(u, p, t) = 0.5
prob = SDEProblem(f, g, 1.0, (0.0, 5.0))
bp = ConstantRateBranchingProblem(prob, 1.0, 2)

# Via solve (convenient wrapper)
sol = solve(bp; reduction=sum)
sol = solve(bp; reduction=sum, output_dt=0.05, transform=log)

# Directly
sol = solve_and_reduce(bp; reduction=sum, output_dt=0.05)
sol = solve_and_reduce(bp; reduction=maximum, output_dt=0.1)

See also: reduce_tree, solve_and_split, ConstantRateBranchingProblem

SciMLBase.remake(bp::ConstantRateBranchingProblem; prob=missing, lifetime=missing, nchild=missing, ndim=missing, kwargs...)

Remake a ConstantRateBranchingProblem with modified fields. The fields prob, lifetime, nchild, and ndim can be replaced directly. Any additional keyword arguments (e.g. u0, tspan, p) that are not fields of ConstantRateBranchingProblem are passed to SciMLBase.remake(bp.prob; kwargs...) to modify the single-particle dynamics problem. This works for both SDEProblem and JumpProblem inner problems.

Examples

# Change the number of children
new_bp = remake(bp, nchild=3)

# Change initial condition of the inner problem (shortcut syntax)
new_bp = remake(bp, u0=1.0)

# Change both the lifetime and initial condition
new_bp = remake(bp, lifetime=Exponential(0.5), u0=1.0)

# Enable 2D spatial growth
new_bp = remake(bp, ndim=2)

sample_lifetime(lifetime::T) where T <: UnivariateDistribution

Sample the lifetime of a particle from the provided lifetime distribution. This function draws a random sample from the given distribution.

For backward compatibility, when a ConstantRateBranchingProblem is constructed with a positive real number λ as the second argument, it is interpreted as a branch rate and converted to an exponential distribution with parameter 1/λ.

sample_offspring(nchild::O) where O<:Union{Integer,DiscreteUnivariateDistribution}

Sample the number of offspring for a particle. If the input nchild is an integer, that integer is returned. If nchild is a discrete distribution, a sample from that distribution is returned.

get_timespan(prob::P) where P<:SciMLBase.AbstractDEProblem

Get the timespan of the problem prob. Works for SDEProblems and JumpProblems. Throws an error for NoiseProblems.

get_timespan(tree::BranchingProcessSolution)

Get the timespan of a BranchingProcessSolution tree, defined as the time from the initial time of the root to the latest final time of any of the leaves.

remake_initial_condition(prob::P, tspan, u0=nothing) where P<:SciMLBase.AbstractDEProblem

Remake the problem prob with a new timespan tspan and, optionally, a new initial condition u0. Works for SDEProblems and JumpProblems. Throws an error for NoiseProblems.