Chemical networks

This module provides routines to read, compare and manipulate chemical networks. The data structure used to represent a chemical network are also used to read in the chemical network used with a model (i.e. the Reactions.out see ReactionNetworkPout).

The new chemical networks can than also be written to files in various formats.

Usage example

This examples reads a network from a ProDiMo Reaction.in file, renames a species in it and writes the network back to a new file in the same format.

import prodimopy.chemistry.network as pcnet

# Create a new network and load the rates from the Reactions.in (local directory)
reacin=pcnet.ReactionNetworkPin()
reacin.load_reactions("Reactions.in")

# now do the renaming
reacin.renameSpecies("C2H6","CH3CH3")

# write it again in a new file, be careful if you directly want to overwrite your old file.
reacin.write_reactions("Reactions.in.renamed")

Source documentation

Created on 15 Feb 2021

@author: Christian Rab

class prodimopy.chemistry.network.Reaction[source]

Bases: object

General representation of one chemical reaction.

id: int: The Reaction Id as is.

idU: int: The Umist reaction id … it is unclear what this actually is.

type: str: The Reaction type identifier as used in ProDiMo

gasphase: bool: Gas phase reaction or not.

reactants: list[str]: The list of reaction reactants (Species names)

products: list[str]: The list of reaction products (Species names)

coeffs: list[ndarray[tuple[int, ...], dtype[float64]]]: The list of ([alpha,beta,gamma]) coefficients per temperature range.

temps: list[ndarray[tuple[int, ...], dtype[float64]]]: The list of temperature ranges given by (min,max).

rate: float | None: The real rate for a given set of physical conditions. Only exists if the network is loaded from an actual model.

method: str | None: And indicator on how the rated was derived. Usually UMIST method identifier (M: measured; C: calculated; E: estimated; L: literature).

accuracy: str | None: A string describing the quality of the reactions coefficients.

comment: str: Any relevant comment for this particular reaction.

arrhenius_rate()[source]

Calculate the Arrhenius rate for the given (from the data) temperature range of the reaction.

Returns:: A tuple containing two lists of numpy arrays. The first list contains the temperature points for each temperature range, the second list contains the corresponding Arrhenius rates for each temperature range.
Return type:: tuple[list[NDArray[np.float64]], list[NDArray[np.float64]]]

compare(reaction, specsonly=False)[source]

Compares the Reaction to the passed reaction.

The routines evaluates for a selection of the properties from class:`Reaction individually if it is equal or not. See below for the return value.^

Parameters:

reaction (Reaction) – Another reaction to compare with.
specsonly (bool) – If True, only the reactants and products are compared. Default: False In that case the returned dictionary only contains the keys reactants and products.

Returns:

A dictionary with boolean values with the following keys. The keys have the same names as the properties of the Reaction class.

type is the reaction type the sames
gasphase is the gasphase flag the same
reactants does have the reaction the same reactants (order does not matter)
products same as reactants but for the products
coeffs are the coefficients (alpha,beta,gamma) the same
temps are the temperature ranges the sames

Return type:

dict[str, bool]

..todo:

- the routine is quite slow if speconly=False
- deal with comparison of the network if actual rates only exist in one network
- make the tolerances for comparison of the rate a parameter somehow

class prodimopy.chemistry.network.ReactionNetwork(name, filename=None, filename_dictionary=None)[source]

Bases: ABC

General representation of a chemical reaction network.

Todo

make it an abstract class that requires the child classes to implement the load_reactions routine.
properly define what is an ice species, e.g. add routines such as species_is_ice which can be overwritten.
properly define non-species, also allow to overwrite that in subclass e.g. routine species_is_nonspecies

Initialize a new chemical reaction network.

Parameters:

name (str) – A name for the reaction network.
filename (str | None) – The name/path of the file containing the Reaction network. Can be None if the network is created from scratch. Or if the data is read later on. If it is not None the network is loaded from the file in the constructor. Default: None
filename_dictionary (str | None) – The name/path of the dictionary file containing the mapping of the names used in this reaction network, to the new names. Can be None if no mapping is used. If it is not None the mapping is applied to the network in the constructor. after reading the network. Default: None

name: str: A name for the reaction network

filename: str | None: The name/path of the file containing the reaction data.

filename_dictionary: str | None: The name/path of the dictionary file containing the mapping of the names used in this reaction network.

reactions: list[Reaction]: The list of reactions included in the network.

non_species: list[str]: A list of Species names that are considered as non species (not really a chemical species) but might be included in the network.

property species: list[str]: The list of individual species in the network. The list is construction from the reactions data.

property elements: list[str]

The list of unique elements in the network, sorted by atomic number

Todo

use parse_stoichiometry() to be consistent in the parsing of the species names.

duplicates()[source]

Returns possible duplicate reactions in the network.

Currently a two reactions are a potential duplicate if they have the same reactants and products.

Returns:: Returns a list of tuples containing the duplicate Reactions or an empty list if there are no duplicates
Return type:: list[tuple[Reaction, Reaction]]

rename_species(old, new, variants='')[source]

Renames species with name old with species name new in the whole network.

Parameters:

old (str) – The species name of the species that should be renamed.
new (str) – The new name of the species.
variants (tuple[str] | None) – Postfix to the species names that should also be considered for renaming. For example if variants=[“”,”+”,”++”,”#”] and the old name is SP and the new name is SPNEW the routine would rename SP,SP+,SP# with SPNEW,SPNEW+,SPNEW# , respectively. If None -> variants=(“”), which means only identical names are considered,. Default: (“”)
log (boolean) – log detailed changes on the screen. Default: False

Returns:

The total number of changed Reactions.

Return type:

int

map_species_names(filename_dictionary)[source]

Maps the species names in the reaction network to new names based on a provided dictionary file. Note that only the variants (“”,”#”) are considered (see rename_species()).

Parameters:: filename_dictionary (str) – The path to the dictionary file containing the mapping of the names used in this reaction network, to the new names. E.g. Here KIDA to UMIST. see also read_species_dictionary().

compareSpecies(reactionNetwork)[source]

Compares the species list of two networks.

Parameters:

reactionNetwork (ReactionNetwork) – Another reaction network

Returns:

(tuple)
* boolean … equal or not
* list … list of species only in old network
* list … list of species only in new network

find_species(searchstr, checkelements=False, ignore_non_species=False)[source]

Find all species of the network that contain the given search string.

Parameters:

searchstr (str | list[str]) – The string or list of strings to search for in the species names.
checkelements (bool) – If True the search string is interpreted as an element symbol and the routine checks additionally the stoichiometry for that element in the species. Default: False
ignore_non_species (bool) – If True the routine ignores species that are in the non_species list of the network. Default: False

Returns:

A list of species names that contain the search string.

Return type:

list[str]

find_isomers(species=None)[source]

Find all species that have isomers (i.e. the same stoichiometry) in the network. If a species name is passed, the routine only returns the isomer group for this particular species. The passed species does not have to be in the network. If the species is the only isomer, an empty tuple is returned. If the species is not in the network, it is not included in the returned isomer group.

By default the routine searches for all species in the network for isomers, and returns a list of tuples containing those isomer groups.

Ice species ‘#’ are not considered as isomers, and are ignored in the search.

Parameters:

network – The reaction network to analyze.
species (str | None) – The specific species name to find isomers for, by default None -> search for all isomer groups in the network

Returns:

A list of tuples where each tuple contains isomer species names (sorted) or just a tuple (if only one species is passed).

Return type:

list[tuple[str]] or tuple[str]

find(species, reac_type=None, formation=True, destruction=True)[source]

Find all reactions that involve any of the given species.

Parameters:

species (str | list[str]) – The species name(s) to search for in the reactions.
reac_type (str | None) – If provided, only reactions of this type will be returned.
formation (bool) – If True, include formation reactions.
destruction (bool) – If True, include destruction reactions.

Returns:

A list of reactions that involve the given species. Please note those are references to the original Reactions in the network, so if you update those they are also updated in the network. If you want to avoid that you need to copy the reactions.

Return type:

list[Reaction]

remove(id)[source]

Remove a reactions or a list of reactions from the network.

Parameters:: id (int | list[int] | Reaction | list[Reaction]) – The reaction(s) to delete.
Returns:: The deleted reactions.
Return type:: list[Reaction]

compare(reactionNetwork, printresults=False, eqfunc=None, chgfunc=None)[source]

Compares the network to the given reaction network. Currently simply prints out the results to stdout.

Parameters:

reactionNetwork (ReactionNetwork) – The reaction network to compare with.
printresults (bool) – Print results to stdout.
eqfunc (Callable[[dict[str, bool]], bool] | None) – Pass a function that decides if two reactions are equal. The function must take as an argument the outcome of Reaction.compare().
chgfunc (Callable[[dict[str, bool]], bool] | None) – Pass a function that decides if two reactions are the same reaction, but some other quantities differ (coefficients, type etc). The function must take as an argument the outcome of Reaction.compare().

Returns:

tuple containing:

bool: are the networks equal or not
list of changed reactions (both, and the equal dictionary)
list of reactions only in this network
list of reactions only in the passed network

Return type:

tuple[bool, list[tuple[Reaction, Reaction, dict[str, bool]]], list[Reaction], list[Reaction]]

abstractmethod load_reactions(filename=None)[source]

In this routine the reading of the Reaction network needs to be implemented.

Fills in the reactions field in this class.

Parameters:: filename (str) – The name/path of the file to read

abstractmethod load_rates(filename=None)[source]

In this routine the reading of rate coefficients need to be implemented.

Parameters:: filename (str) – The name/path of the file to read

write_reactions(filename, fmt='csv', coeffsfmts=('{:8.5}', '{:8.6}', '{:12.4f}'), header=None)[source]

Writes the network to a file. The default format is the csv format as used in ProDiMo (i.e. Reactions.in.csv). Multiple temperature fits are only supported in the csv format.

Warning

Not well tested … so be careful.

Parameters:

filename (str) – The name/path of the output file. And existing one will be overwritten.
fmt (str) –
- csv … write the network in a csv format, as used by ProDiMo.
- oldprodimo … write the network in the old ProDiMo format, as used in Reactions.in. This format is not recommended any more.
header (str) – A string that is written as a header to the file. Only used for fmt=’csv’. The string will be written as is, but using header.strip(), it is recommended to use triple-quoted strings. Default: None

class prodimopy.chemistry.network.ReactionNetworkPout(name='ReactionNetworkPout', modeldir=None)[source]

Bases: ReactionNetwork

Implementation of ReactionNetwork for the Reactions.out of a ProDiMo model.

Parameters:

modeldir (str) – If set the init routine tries to load the Reactions.out and the rates log file (at the moment) from the given modeldir directory. But continues if it does not work.
name (str)

load_reactions(filename='Reactions.out')[source]

Reads a reaction network in the format of the ProDiMo Reactions.out.

Fills in the reactions field in this class.

Parameters:: filename (str) – The name/path of the file to read.

load_rates(filename=None)[source]

Load the real rates for a given set of physical conditions. Can be used for comparison. Has to be run after load_reactions.

Parameters:: filename (str | None)

class prodimopy.chemistry.network.ReactionNetworkPin(name, filename=None, filename_dictionary=None, modeldir=None)[source]

Bases: ReactionNetwork

Implementation of ReactionNetwork for the Reactions.in and Reactions.in.csv formats of ProDiMo. As the csv format is the same as for the UMIST data files in ProDiMo it can/should also be used for those.

Parameters:

name (str) – a name for the network.
filename (str | None) – The name/path of the file to read.
filename_dictionary (str | None) – The name/path of the dictionary file to read.
modeldir (str | None) – If set, the init routine tries to load the Reactions.in from within modeldir.

load_reactions(filename='Reactions.in', fmt=None, threeReactants=True)[source]

Reads a reaction network in the format of the ProDiMo Reactions.in or Reactions.in.csv format.

Fills reactions field in this class.

Todo

Include reading of T-dependent rates for the csv format.
more sophisticated guessing of the file format.

Parameters:

filename (str) – The name/path of the file to read.
fmt (str) – Format of the file. Currently either in for old Reactions.in or csv for the UMIST csv format. If None the routine tries to guess the format (primitive at the moment).
threeReactants (bool)

updateAdsorptionDesorptionReactions(species, comment='')[source]

Adds or updates existing adsorption and desorption reactions for the given species. Only add and update operations are supported.

For e.g. CO# it would create the following reactions:

` 1281,IC,CO,dust,,CO#,dust,,,0.00E+00,+0.00000,0.0000,L,10,41000,B, 1282,DT,CO#,dust,,CO,dust,,,0.00E+00,+0.00000,899.0000,L,10,41000,B,Martin 899/1150GH06 1283,DC,CO#,CRP,,CO,,,,0.00E+00,+0.00000,899.0000,L,10,41000,B,GH06 1284,DP,CO#,PHOTON,,CO,,,,0.00E+00,+0.00000,899.0000,L,10,41000,B,GH06 ` :param species: A list of species names or a dictionary with species names as keys and their corresponding binding energies (in K) as values. :type species: list[str,float] or dict[str,float]

Returns:: The list of the created and added reactions.
Return type:: list[Reaction]
Parameters:: species (list[tuple[str, float]] | dict[str, float])

load_rates(filename=None)[source]

In this routine the reading of rate coefficients need to be implemented.

Parameters:: filename (str) – The name/path of the file to read

class prodimopy.chemistry.network.ReactionNetworkUMIST2022(name, filename=None, filename_dictionary=None)[source]

Bases: ReactionNetwork

Implementation of ReactionNetwork for the UMIST 2022 database format.

Initialize a new chemical reaction network.

Parameters:

name (str) – A name for the reaction network.
filename (str | None) – The name/path of the file containing the Reaction network. Can be None if the network is created from scratch. Or if the data is read later on. If it is not None the network is loaded from the file in the constructor. Default: None
filename_dictionary (str | None) – The name/path of the dictionary file containing the mapping of the names used in this reaction network, to the new names. Can be None if no mapping is used. If it is not None the mapping is applied to the network in the constructor. after reading the network. Default: None

load_reactions(filename='rate22_final.rates.csv')[source]

Reads a UMIST2022 reaction network in the format of the UMIST 2022 csv format

Fills reactions field in this class.

Parameters:: filename (str) – The name/path of the file to read.

load_rates(filename=None)[source]

In this routine the reading of rate coefficients need to be implemented.

Parameters:: filename (str) – The name/path of the file to read

class prodimopy.chemistry.network.ReactionNetworkKIDA(name, filename=None, filename_dictionary=None)[source]

Bases: ReactionNetwork

Implementation of ReactionNetwork for the KIDA database format.

Parameters:

name – a name for the network.
filename (str | None) – The name/path of the file to read.
filename_dictionary (str | None) – The name/path of the dictionary file to read.

load_reactions(filename='KIDA_uva_2024.dat')[source]

Reads a KIDA reaction network in the format of the KIDA database.

Fills reactions list in this object.

Parameters:: filename (str) – The name/path of the file to read.

write_reactions(filename='Reactions.in.new', fmt='KIDA')[source]

Writes a KIDA reaction network in the format of the KIDA database.

Parameters:

filename (str) – The name/path of the file to write.
fmt (str) – Format of the file. Currently only KIDA is supported.

load_rates(filename=None)[source]

In this routine the reading of rate coefficients need to be implemented.

Parameters:: filename (str) – The name/path of the file to read

prodimopy.chemistry.network.fit_double_arrhenius_ranges(reactions, min_points=3)[source]

Fit two single Arrhenius laws with non-overlapping, contiguous temperature ranges to the total (summed) rate of a reaction with a double/triple … Arrhenius fit.

The split point between the two ranges is chosen to minimise the total log-space least-squares residual across both segments.

Returns a new Reaction with the coefficents for the two temperature ranges, The returned Reaction is a new object identical to the input (first Reaction in reactions list) except that its coeffs and temps are replaced by the two fitted ranges (one entry each), and its comment is extended with “Newly fitted with 2 temperature fit”.

Routine was partly written by Claude.

Parameters:

reactions (list[Reaction]) – A list of Reaction objects that have the same reactants and products, but multie Arrhenius fits for the same temperature range (only one per Reaction). It is assumend that the actual total rate of the Reaction is the sum of those Rreactions rates. And this is what will be fitted.
min_points (int) – Minimum number of positive-rate points required in each segment to perform the fit. If the total number of positive-rate points is less than 2 * min_points, the function returns None and does not perform the fit.

Return type:

A new Reaction object with the fitted coefficients and temperature ranges, or None if the fit cannot be performed (too few positive-rate points).

prodimopy.chemistry.network.parse_stoichiometry(species)[source]

Parse the species names and calculates the soichiometry and charge of the species.

Several special cases are handled:

ice suffix ‘#’ is stripped before parsing.

exc (e.g. H2exc) is stripped before parsing

l-, c-, t- isomer prefixes are stripped before parsing

Z and PAH are treated as elements

Emits a UserWarning for any characters that cannot be matched to a known element.

Parameters:: species (str) – The species name to parse.
Returns:: A dictionary containing the stoichiometry and charge of the species.
Return type:: dict[str, int]

Examples

>>> parse_stoichiometry('CH3ON')
({'H': 3, 'C': 1, 'N': 1, 'O': 1,'charge': 0})
>>> parse_stoichiometry('H3O+')
({'H': 3, 'O': 1, 'charge': 1})
>>> parse_stoichiometry('H-')
({'H': 1, 'charge': -1})
>>> parse_stoichiometry('H2O#')
({'H': 2, 'O': 1, 'charge': 0})

prodimopy.chemistry.network.readSpeciesin(filename='Species.in')[source]

Reads the list of species and their abundances from a ProDiMo Species.in file.

Parameters:: filename (str, optional) – The input file name, by default “Species.in”
Returns:: A list of species names and their abundances.
Return type:: list[tuple[str, float]]

prodimopy.chemistry.network.updateSpeciesin(filename, newSpecies)[source]

Updates the list of species and their abundances in a ProDiMo Species.in file. Supports only add an overwrite operations.

Parameters:

filename (str) – The input file name.
newSpecies (list[tuple[str, float]] | dict[str, float]) – A list or dictionary of new species names and their abundances.

prodimopy.chemistry.network.writesSpeciesin(network, filename='Species.in', excludeElems=[], excludeSpecies=[])[source]

Writes the list of species in the chemical network to a file in the ProDiMo format. All species have zero (very low) initial concentrations.

Parameters:

network (ReactionNetwork) – A Chemical network
filename – The output file name, by default “Species.in”

prodimopy.chemistry.network.readElementsin(filename='Elements.in')[source]

Reads the list of elements and their abundances and masses from a ProDiMo Elements.in file.

Parameters:: filename (str, optional) – The input file name, by default “Elements.in”
Returns:: A list of element names and their abundances and masses.
Return type:: list[tuple[str, float, float]]

prodimopy.chemistry.network.updateElementsin(filename, newElements)[source]

Updates the list of elements in a ProDiMo Elements.in file with new elements and their abundances and masses.

Parameters:

filename (str) – The input file name, by default “Elements.in”
newElements (list[tuple[str, float, float]], optional) – A list of new elements to add to the file, by default []

prodimopy.chemistry.network.writesElementsin(network, filename='Elements.in', excludeElems=[])[source]

Writes the list of elements in the chemical network to a file in the ProDiMo format.

Parameters:

network (ReactionNetwork) – A Chemical network
filename (str, optional) – The output file name, by default “Elements.in”

prodimopy.chemistry.network.read_AdsorptionEnergies(filename='AdsorptionEnergies.in')[source]

Reads the adsorption energies data from a file and returns a dictionary mapping species to their adsorption energy values. If the value given for a species in the file is not a valid float, the value will be set to None, the species is still included in the dictionary.

Parameters:: filename (str, optional) – The input file name, by default “AdsorptionEnergies.in”
Returns:: A dictionary mapping species names to their adsorption energy values (in K). If a value is not a valid float, it will be set to None.
Return type:: dict[str, float | None]

prodimopy.chemistry.network.read_heat_of_formation(filename='Heat_of_formation.dat')[source]

Reads the heat of formation data from a file and returns a dictionary mapping species to their heat of formation values. If the value given for a species in the file is not a valid float, the value will be set to None, the species is still included in the dictionary.

Parameters:

filename – The input file name, by default “Heat_of_formation.dat”.

Returns:

dict[str, float | None] – A dictionary mapping species names to their heat of formation values (in K). If a value is not a valid float, it will be set to None.
.. todo:: –
- Allow for different format for the input files, e.g. without number of lines.

Return type:

dict[str, float | None]

prodimopy.chemistry.network.read_species_dictionary(filename)[source]

Reads dictionary data from a file and returns a dictionary, mapping species to their mapped names.

The format has to be a two column file, with the first column being the species name and the second column being the mapped name, separated by whitespace characters. After this two columns, separated by whitespace, any comment is also allowed, but will not be read.

The first line of the file will always be ignored, as it is assumed to be a header line. A # as the first character of a line indicates a comment line, which will also be ignored.

Parameters:: filename (str)
Return type:: dict[str, str]