Integrators¶
Contains integrator methods that will be used to integrate ODEs or SDEs.
-
class
epipack.integrators.
IntegrationMixin
[source]¶ Bases:
object
A helper MixIn that enables the base class to set initial conditions and integrate numerical ODEs.
Expects the base class to have the following methods and attributes:
get_numerical_dydt()
get_compartment_id()
compartments
N_comp
-
integrate
(time_points, return_compartments=None, *args, **kwargs)[source]¶ Returns values of the given compartments at the demanded time points (as a dictionary).
- Parameters
time_points (np.ndarray) -- An array of time points at which the compartment values should be evaluated and returned.
return_compartments (list, default = None) -- A list of compartments for which the result should be returned. If
return_compartments
is None, all compartments will be returned.integrator (str, default = 'dopri5') -- Which method to use for integration. Currently supported are
'euler'
and'dopri5'
. If'euler'
is chosen, \(\delta t\) will be determined by the difference of consecutive entries intime_points
.adopt_final_state (bool, default = False) -- Whether or not to adopt the final state of the integration as new initial conditions.
diffusion_constants (numpy.ndarray) -- Scalar and constant diffusion coefficients as prefactors for each compartment's Wiener process (has to be of same length as y0)
- Returns
result -- Dictionary mapping compartments to time series.
- Return type
-
integrate_and_return_by_index
(time_points, return_compartments=None, integrator='dopri5', adopt_final_state=False, diffusion_constants=None)[source]¶ Returns values of the given compartments at the demanded time points (as a numpy.ndarray of shape
(return_compartments), len(time_points)
.- Parameters
time_points (np.ndarray) -- An array of time points at which the compartment values should be evaluated and returned.
return_compartments (list, default = None) -- A list of compartments for which the result should be returned. If
return_compartments
is None, all compartments will be returned.integrator (str, default = 'dopri5') -- Which method to use for integration. Currently supported are
'euler'
and'dopri5'
. If'euler'
is chosen, \(\delta t\) will be determined by the difference of consecutive entries intime_points
.adopt_final_state (bool, default = False) -- Whether or not to adopt the final state of the integration
diffusion_constants (numpy.ndarray) -- Scalar and constant diffusion coefficients as prefactors for each compartment's Wiener process (has to be of same length as y0)
- Returns
result -- First axis accesses system states, second axis refers to the system state corresponding to the time points in
time_points
.- Return type
numpy.ndarray
-
integrate_and_return_by_index_until
(t0, stop_condition, return_compartments=None, adopt_final_state=False, solve_ivp_kwargs={})[source]¶ Integrate the system until a stop condition is reached, using
scipy.integrate.solve_ivp
. Returns values of the given compartments at the demanded time points (as a numpy.ndarray of shape(return_compartments), len(time_points)
.- Parameters
t0 (float) -- Initial time
stop_condition (function) --
A function
stop_condition(t,y)
that flips from positive to negative values when a certain condition is reached. For instance, when an element reaches a certain threshold value:stop_condition = lambda t, y: 0.3 - y[2]
return_compartments (list, default = None) -- A list of compartments for which the result should be returned. If
return_compartments
is None, all compartments will be returned.adopt_final_state (bool, default = False) -- Whether or not to adopt the final state of the integration
solve_ivp_kwargs (dict) -- Keyword arguments that will be passed to
scipy.integrate.solve_ivp
.
- Returns
t (float) -- Times when the system state was recorded The final time corresponds to the time when the stop condition was reached.
y (numpy.ndarray) -- First axis accesses system states, second axis refers to the system state corresponding to the time points in t.
-
integrate_until
(t0, stop_condition, return_compartments=None, adopt_final_state=False, solve_ivp_kwargs={})[source]¶ Integrate the system until a stop condition is reached, using
scipy.integrate.solve_ivp
.- Parameters
t0 (float) -- Initial time
stop_condition (function) --
A function
stop_condition(t,y)
that flips from positive to negative values when a certain condition is reached. For instance, when an element reaches a certain threshold value:stop_condition = lambda t, y: 0.3 - y[2]
return_compartments (list, default = None) -- A list of compartments for which the result should be returned. If
return_compartments
is None, all compartments will be returned.adopt_final_state (bool, default = False) -- Whether or not to adopt the final state of the integration
solve_ivp_kwargs (dict) -- Keyword arguments that will be passed to
scipy.integrate.solve_ivp
.
- Returns
t (float) -- Times when the system state was recorded The final time corresponds to the time when the stop condition was reached.
result (dict) -- Dictionary mapping compartments to time series.
-
set_initial_conditions
(initial_conditions, initial_time=0.0, allow_nonzero_column_sums=False)[source]¶ Set the initial conditions for integration
- Parameters
initial_conditions (dict) -- A dictionary that maps compartments to a fraction of the population. Compartments that are not set in this dictionary are assumed to have an initial condition of zero.
allow_nonzero_column_sums (bool, default = False) -- If True, an error is raised when the initial conditions do not sum to the population size.
-
epipack.integrators.
integrate_SDE
(dydt, t, y0, diffusion_constants, *args)[source]¶ Integrate an SDE system with Euler's method.
- Parameters
dydt (function) -- A function returning the momenta of the deterministic ODE system.
t (numpy.ndarray of float) -- Array of time points at which the functions should be evaluated. Time steps must be equally spaced.
y0 (numpy.ndarray) -- Initial conditions
diffusion_constants (numpy.ndarray) --
Scalar and constant diffusion coefficients as prefactors for each compartment's Wiener process (has to be of same length as y0)
corresponds to \(D_i\) in
\[dY_i = f_i(\mathbf Y,t) dt + D_i dW_i\]*args (
list
) -- List of parameters that will be passed to the momenta function.
-
epipack.integrators.
integrate_dopri5
(dydt, t, y0, *args)[source]¶ Integrate an ODE system with the Runge-Kutte Dormand Prince method (with step-size control).
- Parameters
dydt (function) -- A function returning the momenta of the ODE system.
t (numpy.ndarray of float) -- Array of time points at which the functions should be evaluated.
y0 (numpy.ndarray) -- Initial conditions
*args (
list
) -- List of parameters that will be passed to the momenta function.
-
epipack.integrators.
integrate_euler
(dydt, t, y0, *args)[source]¶ Integrate an ODE system with Euler's method.
- Parameters
dydt (function) -- A function returning the momenta of the ODE system.
t (numpy.ndarray of float) -- Array of time points at which the functions should be evaluated.
y0 (numpy.ndarray) -- Initial conditions
*args (
list
) -- List of parameters that will be passed to the momenta function.
-
epipack.integrators.
integrate_ivp_until
(dydt, t0, y0, stop_condition, solve_ivp_kwargs={})[source]¶ Integrate an ODE system with the
scipy.integrate.solve_ivp
method.- Parameters
dydt (function) -- A function returning the momenta of the ODE system.
t0 (float) -- Initial time
y0 (numpy.ndarray) -- Initial conditions
stop_condition (function) --
A function
stop_condition(t,y)
that flips from positive to negative values when a certain condition is reached. For instance, when an element reaches a certain threshold value:stop_condition = lambda t, y: 0.3 - y[2]
solve_ivp_kwargs (dict) -- Keyword arguments that will be passed to
scipy.integrate.solve_ivp
.
- Returns
t (float) -- Times when the system state was recorded The final time corresponds to the time when the stop condition was reached.
y (numpy.ndarray) -- First axis accesses system states, second axis refers to the system state corresponding to the time points in
t
.
-
epipack.integrators.
time_leap_ivp
(t0, y0, get_event_rates, rand=None)[source]¶ Compute a time leap for time-varying rate functions based on Gillespie's SSA using Newton's method on a quad-integrator.
- Parameters
- Returns
new_t -- A new time.
- Return type