catmap.solvers package¶
Submodules¶
catmap.solvers.integrated_rate_control_solver module¶
-
class
catmap.solvers.integrated_rate_control_solver.
IntegratedRateControlSolver
(reaction_model=None)[source]¶ Bases:
catmap.solvers.solver_base.SolverBase
Class for estimating rates based on the degree of rate control screening method {citation after published}
-
__getattr__
(attr)¶ Return the value of the reaction model instance if its there. Otherwise return the instances own value (or none if the instance does not have the attribute defined and the attribute is not private)
-
__module__
= 'catmap.solvers.integrated_rate_control_solver'¶
-
__setattr__
(attr, val)¶ Set attribute for the instance as well as the reaction_model instance
-
catmap.solvers.mean_field_solver module¶
-
class
catmap.solvers.mean_field_solver.
MeanFieldSolver
(reaction_model=None)[source]¶ Bases:
catmap.solvers.solver_base.SolverBase
Class for handling mean-field type kinetic models. Can be sub-classed to get functionality for steady-state solutions, sabatier solutions, etc.
-
__getattr__
(attr)¶ Return the value of the reaction model instance if its there. Otherwise return the instances own value (or none if the instance does not have the attribute defined and the attribute is not private)
-
__module__
= 'catmap.solvers.mean_field_solver'¶
-
__setattr__
(attr, val)¶ Set attribute for the instance as well as the reaction_model instance
-
get_apparent_activation_energy
(rxn_parameters, epsilon=1e-10)[source]¶ returns apparent Arrhenius activation energies (in units of R) for production/consumption of each gas phase species. Calculated as E_app = T^2(dlnr_+/dT)=(T^2/r_+)(dr_+/dT), where r+ is the TOF :param rxn_parameters: reaction paramenters, see solver-base :param epsilon: degree of pertubation in temperature :type epsilon: float, optional
-
get_directional_rates
(rxn_parameters)[source]¶ get the exchange current density of a certain map entry on all elementary rxns for a given direction
type: direction: str (‘kf’ or ‘kr’)
-
get_elem_ec
(rxn_num, rxn_parameters, direction)[source]¶ return the ec on a certain coverage_map entry of a certain elementary step based on rxn_number and direction given
type: rxn_num: float direction: str (‘kf’ or ‘kr’)
-
get_empty_site_cvgs
()[source]¶ take the coverages at a certain coverage_map entry and return the dict of all the empty-sites coverages i.e. dict[site_name] = coverage
type: coverages: list
-
get_interacting_energies
(rxn_parameters)[source]¶ return the integral energy under high coverage with interactions :param rxn_parameters: reaction parameters, see solver-base
-
get_rate_control
(rxn_parameters)[source]¶ return list of degree of rate control for each reaction Ref: Stegelmann et al., DOI: 10.1021/ja9000097 :param rxn_parameters: reaction parameters, see solver-base
-
get_rxn_order
(rxn_parameters, epsilon=1e-10)[source]¶ return the reaction orders for the reactants :param rxn_parameters: reaction parameters, see solver-base :param epsilon: degree of perturbation in pressure :type epsilon: float, optional
-
get_rxn_rates
(coverages, rate_constants)[source]¶ returns list of reaction rate for each elementary reaction based on reaction constants & coverage .. todo:: coverages, rate_constants
-
get_selectivity
(rxn_parameters, weights=None)[source]¶ return list of selectivity of each reaction :param rxn_parameters: reaction parameters, see solver-base :param weights: weights for each species. Defaults to 1 for all species
-
get_selectivity_control
(rxn_parameters)[source]¶ return the list of degree of selectivity control for each rxn :param rxn_parameters: reaction parameters, see solver-base
-
get_turnover_frequency
(rxn_parameters, rates=None, verify_coverages=True)[source]¶ return list of turnover frequencies of all the gas-phase species :param rates: list of rates of each rxn :type rates: list :param verify_coverages: verify that the species has a certain value for the coverage :type verify_coverages: bool, optional :param rxn_parameters: reaction parameters, see solver-base
-
jacobian_equations
(adsorbate_interactions=True)[source]¶ Composes analytical expressions for the Jacobian matrix. Assumes:
kf is defined as a list of forward rate-constants kr is defined as a list of reverse rate-constants theta is defined as a list of coverages p is defined as a list of pressures
If the rate constants depend on coverage, use adsorbate_interactions = True. Assumes:
kB is defined as Boltzmann’s constant T is defined as the temperature dEf is defined as a list of lists where dEf[i][j] is the
derivative of forward activation free energy i wrt coverage j- dEr is defined as a list of lists where dEr[i][j] is the
- derivative of reverse activation free energy i wrt coverage j
Param: adsorbate_interactions: tell if need to include interactions Type: adsorbate_interactions: bool, optional
-
rate_equation_term
(species_list, rate_constant_string, d_wrt=None)[source]¶ Function to compose a term in the rate equation - e.g. kf[1]*theta[0]*p[0] :param species_list: list of species in rate equations :type species_list: list
-
rate_equations
()[source]¶ Compose analytical expressions for the reaction rates and change of surface species wrt time (dc/dt). Assumes:
kf is defined as a list of forward rate-constants kr is defined as a list of reverse rate-constants theta is defined as a list of coverages p is defined as a list of pressures
-
reaction_energy_equations
(adsorbate_interactions=True)[source]¶ Composes a list of analytical expressions which give the reaction and activation energies for elementary steps. Note that while this is useful primarily for models with adsorbate-interactions (otherwise these energetics can easily be obtained by the reaction model itself), they are technically valid for all mean-field models. Assumes:
- Gf is a list of formation energies ordered as
- adsorbate_names+transition_state_names
If model includes adsorbate interactions then use adsorbate_interactions = True to include dEa/dtheta in the output. Assumes:
- dGs is a matrix/array of derivatives of free energies wrt coverages
- such that dGs[:,i] is a vector of derivatives of the free energy of species i wrt each coverage ordered as adsorbate_names
Param: adsorbate_interaction: specify whether or not to include interactions Type: adsorbate_interaction: bool, optional
-
catmap.solvers.solver_base module¶
-
class
catmap.solvers.solver_base.
NewtonRoot
(f, x0, matrix, mpfloat, Axb_solver, **kwargs)[source]¶ Hacked from MDNewton in mpmath/calculus/optimization.py in order to allow for constraints on the solution.
Find the root of a vector function numerically using Newton’s method.
f is a vector function representing a nonlinear equation system.
x0 is the starting point close to the root.
J is a function returning the Jacobian matrix for a point.
Supports overdetermined systems.
Use the ‘norm’ keyword to specify which norm to use. Defaults to max-norm. The function to calculate the Jacobian matrix can be given using the keyword ‘J’. Otherwise it will be calculated numerically.
Please note that this method converges only locally. Especially for high- dimensional systems it is not trivial to find a good starting point being close enough to the root.
It is recommended to use a faster, low-precision solver from SciPy [1] or OpenOpt [2] to get an initial guess. Afterwards you can use this method for root-polishing to any precision.
[1] http://scipy.org
-
__module__
= 'catmap.solvers.solver_base'¶
-
maxsteps
= 10¶
-
-
class
catmap.solvers.solver_base.
SolverBase
(reaction_model=None)[source]¶ Bases:
catmap.ReactionModelWrapper
-
__getattr__
(attr)¶ Return the value of the reaction model instance if its there. Otherwise return the instances own value (or none if the instance does not have the attribute defined and the attribute is not private)
-
__init__
(reaction_model=None)[source]¶ Class for `solving’ for equilibrium coverages and rates as a function of reaction parameters. This class acts as a base class to be inherited by other solver classes, but is not functional on its own.
rxn_parameters: list of necessary parameters to solve the kinetic system. This will usually be populated by the scaler.
A functional derived solver class must also contain the methods:
- get_coverage(): a function which returns coverages for each
- adsorbate as a list [cvg_ads1,cvg_ads2,…]
- get_rate(): a function which returns steady-state reaction
- rates for each elementary step as a list [rate_rxn1,rate_rxn2,…]
- get_residual(): a function for computing the norm of the residual. This
- is the condition which will be minimized to reach steady-state.
compile(): a function to set-up/compile the solver.
-
__module__
= 'catmap.solvers.solver_base'¶
-
__setattr__
(attr, val)¶ Set attribute for the instance as well as the reaction_model instance
-
catmap.solvers.steady_state_solver module¶
-
class
catmap.solvers.steady_state_solver.
SteadyStateSolver
(reaction_model=None)[source]¶ Bases:
catmap.solvers.mean_field_solver.MeanFieldSolver
-
__getattr__
(attr)¶ Return the value of the reaction model instance if its there. Otherwise return the instances own value (or none if the instance does not have the attribute defined and the attribute is not private)
-
__module__
= 'catmap.solvers.steady_state_solver'¶
-
__setattr__
(attr, val)¶ Set attribute for the instance as well as the reaction_model instance
-
bisect_interaction_strength
(rxn_parameters, valid_strength, valid_coverages, target_strength, max_bisections, findrootArgs={})[source]¶ TODO:
-
get_apparent_activation_energy
(rxn_parameters, epsilon=1e-10)¶ returns apparent Arrhenius activation energies (in units of R) for production/consumption of each gas phase species. Calculated as E_app = T^2(dlnr_+/dT)=(T^2/r_+)(dr_+/dT), where r+ is the TOF :param rxn_parameters: reaction paramenters, see solver-base :param epsilon: degree of pertubation in temperature :type epsilon: float, optional
-
get_coverage
(rxn_parameters, c0=None, findrootArgs={})[source]¶ Return coverages for given reaction parameters and coverage constraints.
Parameters: - rxn_parameters ([float]) – Sequence of rxn_parameters
- c0 (TODO) – Coverage constraints.
- findrootArgs – deprecated
-
get_directional_rates
(rxn_parameters)¶ get the exchange current density of a certain map entry on all elementary rxns for a given direction
type: direction: str (‘kf’ or ‘kr’)
-
get_elem_ec
(rxn_num, rxn_parameters, direction)¶ return the ec on a certain coverage_map entry of a certain elementary step based on rxn_number and direction given
type: rxn_num: float direction: str (‘kf’ or ‘kr’)
-
get_empty_site_cvgs
()¶ take the coverages at a certain coverage_map entry and return the dict of all the empty-sites coverages i.e. dict[site_name] = coverage
type: coverages: list
-
get_ideal_coverages
(rxn_parameters, c0=None, refresh_rate_constants=True, findrootArgs={})[source]¶ Return
TODO:
-
get_initial_coverage
(rxn_parameters)[source]¶ Return coverages based on probabilties according to the Boltzmann distribution and the adsorption energies for a given sequence of rxn_parameters.
Parameters: rxn_parameters – Sequence of reaction parameters
-
get_interacting_coverages
(rxn_parameters, c0=None, interaction_strength=1.0, findrootArgs={})[source]¶ TODO:
-
get_interacting_energies
(rxn_parameters)¶ return the integral energy under high coverage with interactions :param rxn_parameters: reaction parameters, see solver-base
-
get_rate
(rxn_parameters, coverages=None, verify_coverages=True, **coverage_kwargs)¶
-
get_rate_constants
(rxn_parameters, coverages)[source]¶ Return rate constants for given sequence of reaction parameters and coverages.
Parameters:
-
get_rate_control
(rxn_parameters)¶ return list of degree of rate control for each reaction Ref: Stegelmann et al., DOI: 10.1021/ja9000097 :param rxn_parameters: reaction parameters, see solver-base
-
get_rxn_order
(rxn_parameters, epsilon=1e-10)¶ return the reaction orders for the reactants :param rxn_parameters: reaction parameters, see solver-base :param epsilon: degree of perturbation in pressure :type epsilon: float, optional
-
get_rxn_rates
(coverages, rate_constants)¶ returns list of reaction rate for each elementary reaction based on reaction constants & coverage .. todo:: coverages, rate_constants
-
get_selectivity
(rxn_parameters, weights=None)¶ return list of selectivity of each reaction :param rxn_parameters: reaction parameters, see solver-base :param weights: weights for each species. Defaults to 1 for all species
-
get_selectivity_control
(rxn_parameters)¶ return the list of degree of selectivity control for each rxn :param rxn_parameters: reaction parameters, see solver-base
-
get_steady_state_coverage
(rxn_parameters, steady_state_fn, jacobian_fn, c0=None, findrootArgs={})[source]¶ Return steady-state coverages using catmap.solvers.solver_base.NewtonRoot .
Parameters: - rxn_parameters ([float]) – Sequence of reaction parameters.
- steady_state_fn (TODO) – TODO
- jacobian_fn (TODO) – TODO
- c0 (TODO) – Coverage constraints
- findrootArgs – deprecated
-
get_turnover_frequency
(rxn_parameters, rates=None, verify_coverages=True)¶ return list of turnover frequencies of all the gas-phase species :param rates: list of rates of each rxn :type rates: list :param verify_coverages: verify that the species has a certain value for the coverage :type verify_coverages: bool, optional :param rxn_parameters: reaction parameters, see solver-base
-
jacobian_equations
(adsorbate_interactions=True)¶ Composes analytical expressions for the Jacobian matrix. Assumes:
kf is defined as a list of forward rate-constants kr is defined as a list of reverse rate-constants theta is defined as a list of coverages p is defined as a list of pressures
If the rate constants depend on coverage, use adsorbate_interactions = True. Assumes:
kB is defined as Boltzmann’s constant T is defined as the temperature dEf is defined as a list of lists where dEf[i][j] is the
derivative of forward activation free energy i wrt coverage j- dEr is defined as a list of lists where dEr[i][j] is the
- derivative of reverse activation free energy i wrt coverage j
Param: adsorbate_interactions: tell if need to include interactions Type: adsorbate_interactions: bool, optional
-
optimize_analytical_function
(func_name, func_string, insertion_line, indention_level, *test_args)[source]¶ Replace some common multiplication terms to speed up functions.
-
rate_equation_term
(species_list, rate_constant_string, d_wrt=None)¶ Function to compose a term in the rate equation - e.g. kf[1]*theta[0]*p[0] :param species_list: list of species in rate equations :type species_list: list
-
rate_equations
()¶ Compose analytical expressions for the reaction rates and change of surface species wrt time (dc/dt). Assumes:
kf is defined as a list of forward rate-constants kr is defined as a list of reverse rate-constants theta is defined as a list of coverages p is defined as a list of pressures
-
reaction_energy_equations
(adsorbate_interactions=True)¶ Composes a list of analytical expressions which give the reaction and activation energies for elementary steps. Note that while this is useful primarily for models with adsorbate-interactions (otherwise these energetics can easily be obtained by the reaction model itself), they are technically valid for all mean-field models. Assumes:
- Gf is a list of formation energies ordered as
- adsorbate_names+transition_state_names
If model includes adsorbate interactions then use adsorbate_interactions = True to include dEa/dtheta in the output. Assumes:
- dGs is a matrix/array of derivatives of free energies wrt coverages
- such that dGs[:,i] is a vector of derivatives of the free energy of species i wrt each coverage ordered as adsorbate_names
Param: adsorbate_interaction: specify whether or not to include interactions Type: adsorbate_interaction: bool, optional
-
site_string_list
()¶ Function to compose an analytic expression for the coverage of empty sites
-
substitutions_dict
()¶ Dictionary of substitutions needed for static compiled functions
-
summary_text
()¶ Stub for producing solver summary.
-