Source code for gillespy2.solvers.cpp.tau_leaping_c_solver

# GillesPy2 is a modeling toolkit for biochemical simulation.
# Copyright (C) 2019-2023 GillesPy2 developers.

# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.

# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# GNU General Public License for more details.

# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <>.
import numpy as np

from gillespy2.solvers.cpp.c_decoder import IterativeSimDecoder
from gillespy2.solvers.utilities import solverutils as cutils
from gillespy2.core import GillesPySolver, Model
from gillespy2.core.gillespyError import *
from gillespy2.core import Results 

from .c_solver import CSolver, SimulationReturnCode

from import SanitizedModel

[docs]class TauLeapingCSolver(GillesPySolver, CSolver): """ A Tau Leaping solver for GillesPy2 models. This solver uses an algorithm that calculates multiple reactions in a single step over a given tau step size. The change in propensities over this step are bounded by relative change in state, yielding greatly improved run-time performance with very little trade-off in accuracy. """ name = "TauLeapingCSolver" target = "tau_leap" def __init__(self, model = None, output_directory = None, delete_directory = True, resume=None, variable = False, constant_tau_stepsize=None): self.constant_tau_stepsize = constant_tau_stepsize super().__init__(model=model, output_directory=output_directory, delete_directory=delete_directory, resume=resume, variable=variable) def _build(self, model, simulation_name, variable, debug = False, custom_definitions=None): sanitized_model = SanitizedModel(model, variable=variable) # determine if a constant stepsize has been requested if self.constant_tau_stepsize is not None: sanitized_model.options['GPY_CONSTANT_TAU_STEPSIZE'] = str(float(self.constant_tau_stepsize)) else: sanitized_model.options['GPY_CONSTANT_TAU_STEPSIZE'] = '0' return super()._build(sanitized_model, simulation_name, variable, debug)
[docs] @classmethod def get_solver_settings(cls): """ Returns a list of arguments supported by :returns: Tuple of strings, denoting all keyword argument for this solvers run() method. :rtype: tuple """ return ('model', 't', 'number_of_trajectories', 'timeout', 'increment', 'seed', 'debug', 'profile')
[docs] def run(self=None, model: Model = None, t: int = None, number_of_trajectories: int = 1, timeout: int = 0, increment: int = None, seed: int = None, debug: bool = False, profile: bool = False, variables={}, resume=None, live_output: str = None, live_output_options: dict = {}, tau_tol=0.03, constant_tau_stepsize=None, **kwargs): """ :param model: The model on which the solver will operate. (Deprecated) :type model: gillespy2.Model :param t: End time of simulation. :type t: int :param number_of_trajectories: Number of trajectories to simulate. By default number_of_trajectories = 1. :type number_of_trajectories: int :param timeout: If set, if simulation takes longer than timeout, will exit. :type timeout: int :param increment: Time step increment for plotting. :type increment: float :param seed: The random seed for the simulation. Optional, defaults to None. :type seed: int :param variables: Dictionary of species and their data that will override existing species data. :type variables: dict :param resume: Result of a previously run simulation, to be resumed. :type resume: gillespy2.Results :param live_output: The type of output to be displayed by solver. Can be "progress", "text", or "graph". :type live_output: str :param live_output_options: dictionary contains options for live_output. By default {"interval":1}. "interval" specifies seconds between displaying. "clear_output" specifies if display should be refreshed with each display. :type live_output_options: dict :param tau_tol: Tolerance level for Tau leaping algorithm. Larger tolerance values will result in larger tau steps. Default value is 0.03. :type tau_tol: float :param constant_tau_stepsize: If set, overrides the automatic stepsize selection and uses the given value as the stepsize on each step. :returns: A result object containing the results of the simulation :rtype: gillespy2.Results """ if self is None: # Post deprecation block # raise SimulationError("TauLeapingCSolver must be instantiated to run the simulation") # Pre deprecation block log.warning( """ `` is deprecated. You should use ` Future releases of GillesPy2 may not support this feature. """ ) self = TauLeapingCSolver(model, resume=resume) if model is not None: log.warning('model = gillespy2.model is deprecated. Future releases ' 'of GillesPy2 may not support this feature.') if self.model is None: if model is None: raise SimulationError("A model is required to run the simulation.") self._set_model(model=model) self.constant_tau_stepsize = constant_tau_stepsize self.model.compile_prep() self.validate_model(self.model, model) self.validate_sbml_features(model=self.model) self.validate_tspan(increment=increment, t=t) if increment is None: increment = self.model.tspan[-1] - self.model.tspan[-2] if t is None: t = self.model.tspan[-1] # Validate parameters prior to running the model. self._validate_type(variables, dict, "'variables' argument must be a dictionary.") self._validate_variables_in_set(variables, self.species + self.parameters) self._validate_resume(t, resume) self._validate_kwargs(**kwargs) if resume is not None: t = abs(t - int(resume["time"][-1])) number_timesteps = int(round(t / increment + 1)) args = { "trajectories": number_of_trajectories, "timesteps": number_timesteps, "tau_tol": tau_tol, "end": t, "interval": str(number_timesteps), } if self.variable: populations = cutils.update_species_init_values(self.model.listOfSpecies, self.species, variables, resume) parameter_values = cutils.change_param_values(self.model.listOfParameters, self.parameters, self.model.volume, variables) args.update({ "init_pop": populations, "parameters": parameter_values }) seed = self._validate_seed(seed) if seed is not None: args.update({ "seed": seed }) if live_output is not None: live_output_options['type'] = live_output display_args = { "model": self.model, "number_of_trajectories": number_of_trajectories, "timeline": np.linspace(0, t, number_timesteps), "live_output_options": live_output_options, "resume": bool(resume) } else: display_args = None args = self._make_args(args) decoder = IterativeSimDecoder.create_default(number_of_trajectories, number_timesteps, len(self.model.listOfSpecies)) sim_exec = self._build(self.model,, self.variable, False) sim_status = self._run(sim_exec, args, decoder, timeout, display_args) if sim_status == SimulationReturnCode.FAILED: raise ExecutionError("Error encountered while running simulation C++ file:\n" f"Return code: {int(sim_status)}.\n") trajectories, time_stopped = decoder.get_output() simulation_data = self._format_output(trajectories) if sim_status == SimulationReturnCode.PAUSED: simulation_data = self._make_resume_data(time_stopped, simulation_data, t) if resume is not None: simulation_data = self._update_resume_data(resume, simulation_data, time_stopped) self.result = simulation_data self.rc = int(sim_status) return Results.build_from_solver_results(self, live_output_options)