:py:mod:`BioSANS2020.propagation.stochastic.mystiffcle` ======================================================= .. py:module:: BioSANS2020.propagation.stochastic.mystiffcle .. autoapi-nested-parse:: This module is the mystiffcle module This can propagate non-stiff to moderately stiff stochastic simulation using the chemical langevine equation. Here two versions are provided 1) Tau-adaptive CLE 2) Fix-inreval CLE The following are the list of function for this module. 1. cle_model 2. cle_calculate 3. cle2_calculate Module Contents --------------- Functions ~~~~~~~~~ .. autoapisummary:: BioSANS2020.propagation.stochastic.mystiffcle.cle_model BioSANS2020.propagation.stochastic.mystiffcle.cle_calculate BioSANS2020.propagation.stochastic.mystiffcle.cle2_calculate .. py:function:: cle_model(sp_comp, ks_dict, conc, r_dict, p_dict, stch_var, dtime, del_coef, reg=False) This functions prepare the CLE model for integration Args: sp_comp (dict): dictionary of appearance or position of species or component in the reaction tag of BioSANS topology file. For example; #REACTIONS A => B, -kf1 # negative means to be estimated B => C, kf2 The value of sp_comp is sp_comp = {'A': {0}, 'B': {0, 1}, 'C': {1}} A appears in first reaction with index 0 B appears in first and second reaction with index 0, 1 C appears in second reaction with index 1 ks_dict (dict): dictionary of rate constant that appear in each reactions. For example; #REACTIONS A => B , 0.3 # first reaction B <=> C, 0.1, 0.2 # second reaction The value of ks_dict is ks_dict = { 0 : [0.3], # first reaction 1 : [0.1, 0.2] # second reaction } conc (dict): dictionary of initial concentration. For example; {'A': 100.0, 'B': -1.0, 'C': 0.0} negative means unknown or for estimation r_dict (dict): dictionary of reactant stoichiometry. For example r_dict = { 0: {'A': 1}, # first reaction, coefficient of A is 1 1: {'B': 1} # second reaction, coefficient of B is 1 } p_dict (dict): dictionary of product stoichiometry. For example p_dict = { 0: {'B': 1}, # first reaction, coefficient of B is 1 1: {'C': 1} # second reaction, coefficient of C is 1 } stch_var (numpy.ndarray): stoichiometric matrix. For example v_stoich = np.array([ [ -1, 0 ] # species A [ 1, -1 ] # species B [ 0, 1 ] # species C #1st rxn 2nd rxn ]) dtime (float): step-size del_coef (float): step-size factor or modifier reg (bool, optional): If True, the model is for fix-interval CLE . Defaults to False. Returns: np.ndarray: f_d = fofx * dtime + gofx * sqdt .. py:function:: cle_calculate(tvar, sp_comp, ks_dict, sconc, r_dict, p_dict, stch_var, del_coef=10, rand_seed=1, implicit=False, rfile='') This functions performs the tau-adaptive CLE integration Args: tvar (list): time stamp of simulation sp_comp (dict): dictionary of appearance or position of species or component in the reaction tag of BioSANS topology file. For example; #REACTIONS A => B, -kf1 # negative means to be estimated B => C, kf2 The value of sp_comp is sp_comp = {'A': {0}, 'B': {0, 1}, 'C': {1}} A appears in first reaction with index 0 B appears in first and second reaction with index 0, 1 C appears in second reaction with index 1 ks_dict (dict): dictionary of rate constant that appear in each reactions. For example; #REACTIONS A => B , 0.3 # first reaction B <=> C, 0.1, 0.2 # second reaction The value of ks_dict is ks_dict = { 0 : [0.3], # first reaction 1 : [0.1, 0.2] # second reaction } sconc (dict): dictionary of initial concentration. For example; {'A': 100.0, 'B': -1.0, 'C': 0.0} negative means unknown or for estimation r_dict (dict): dictionary of reactant stoichiometry. For example r_dict = { 0: {'A': 1}, # first reaction, coefficient of A is 1 1: {'B': 1} # second reaction, coefficient of B is 1 } p_dict (dict): dictionary of product stoichiometry. For example p_dict = { 0: {'B': 1}, # first reaction, coefficient of B is 1 1: {'C': 1} # second reaction, coefficient of C is 1 } stch_var (numpy.ndarray): stoichiometric matrix. For example v_stoich = np.array([ [ -1, 0 ] # species A [ 1, -1 ] # species B [ 0, 1 ] # species C #1st rxn 2nd rxn ]) del_coef (float): step-size factor or modifier rand_seed (float): random seed value picked at random for each trajectory. They have been sampled from the calling program. implicit (bool, optional): True means report in time intervals similar to the input time intervals even if actual step is more or less. Defaults to False. rfile (string, optional): name of topology file where some parameters or components are negative indicating they have to be estimated. Defaults to None. Returns: tuple: (time, trajectories) .. py:function:: cle2_calculate(tvar, sp_comp, ks_dict, sconc, r_dict, p_dict, stch_var, del_coef=1, rand_seed=1, rfile='') This functions performs the fix-interval CLE integration Args: tvar (list): time stamp of simulation sp_comp (dict): dictionary of appearance or position of species or component in the reaction tag of BioSANS topology file. For example; #REACTIONS A => B, -kf1 # negative means to be estimated B => C, kf2 The value of sp_comp is sp_comp = {'A': {0}, 'B': {0, 1}, 'C': {1}} A appears in first reaction with index 0 B appears in first and second reaction with index 0, 1 C appears in second reaction with index 1 ks_dict (dict): dictionary of rate constant that appear in each reactions. For example; #REACTIONS A => B , 0.3 # first reaction B <=> C, 0.1, 0.2 # second reaction The value of ks_dict is ks_dict = { 0 : [0.3], # first reaction 1 : [0.1, 0.2] # second reaction } sconc (dict): dictionary of initial concentration. For example; {'A': 100.0, 'B': -1.0, 'C': 0.0} negative means unknown or for estimation r_dict (dict): dictionary of reactant stoichiometry. For example r_dict = { 0: {'A': 1}, # first reaction, coefficient of A is 1 1: {'B': 1} # second reaction, coefficient of B is 1 } p_dict (dict): dictionary of product stoichiometry. For example p_dict = { 0: {'B': 1}, # first reaction, coefficient of B is 1 1: {'C': 1} # second reaction, coefficient of C is 1 } stch_var (numpy.ndarray): stoichiometric matrix. For example v_stoich = np.array([ [ -1, 0 ] # species A [ 1, -1 ] # species B [ 0, 1 ] # species C #1st rxn 2nd rxn ]) del_coef (float): step-size factor or modifier rand_seed (float): random seed value picked at random for each trajectory. They have been sampled from the calling program. rfile (string, optional): name of topology file where some parameters or components are negative indicating they have to be estimated. Defaults to None. Returns: tuple: (time, trajectories)