From 9ce8b82bd3e7fd542237338f05dcb82f47ddb57f Mon Sep 17 00:00:00 2001 From: Richard West Date: Fri, 2 Aug 2024 10:03:56 -0400 Subject: [PATCH 01/14] Tweaking Chemkin species list parser for readability. This took me a while to figure out what was happening. So I changed variable name, removed a redundant if block, and added a couple of comments. --- rmgpy/chemkin.pyx | 14 ++++++-------- 1 file changed, 6 insertions(+), 8 deletions(-) diff --git a/rmgpy/chemkin.pyx b/rmgpy/chemkin.pyx index 068410dacfa..2a6c17cfea9 100644 --- a/rmgpy/chemkin.pyx +++ b/rmgpy/chemkin.pyx @@ -1209,18 +1209,16 @@ def read_species_block(f, species_dict, species_aliases, species_list): if token_upper == 'END': break - site_token = token.split('/')[0] - if site_token.upper() == 'SDEN': + species_name = token.split('/')[0] # CHO*/2/ indicates an adsorbate CHO* taking 2 surface sites + if species_name.upper() == 'SDEN': # SDEN/4.1e-9/ indicates surface site density continue # TODO actually read in the site density processed_tokens.append(token) - if token in species_dict: - species = species_dict[token] - elif site_token in species_dict: - species = species_dict[site_token] + if species_name in species_dict: + species = species_dict[species_name] else: - species = Species(label=token) - species_dict[token] = species + species = Species(label=species_name) + species_dict[species_name] = species species_list.append(species) From 9c04fd37732ae681479c7a1cd11eb59b3fb15277 Mon Sep 17 00:00:00 2001 From: Richard West Date: Fri, 2 Aug 2024 10:06:17 -0400 Subject: [PATCH 02/14] Cleaning up chemkinTest methods for multidentate adsorbate writing. It was fragile, depending on the precise ordering of lines in the test file. --- test/rmgpy/chemkinTest.py | 21 ++++++++------------- 1 file changed, 8 insertions(+), 13 deletions(-) diff --git a/test/rmgpy/chemkinTest.py b/test/rmgpy/chemkinTest.py index 09099b865c9..64b02b2d8ed 100644 --- a/test/rmgpy/chemkinTest.py +++ b/test/rmgpy/chemkinTest.py @@ -805,9 +805,9 @@ def test_write_bidentate_species(self): chemkin_path = os.path.join(folder, "surface", "chem-surface.inp") dictionary_path = os.path.join(folder, "surface", "species_dictionary.txt") chemkin_save_path = os.path.join(folder, "surface", "chem-surface-test.inp") - species, reactions = load_chemkin_file(chemkin_path, dictionary_path) - - surface_atom_count = species[3].molecule[0].get_num_atoms("X") + species, reactions = load_chemkin_file(chemkin_path, dictionary_path, use_chemkin_names=True) + chox3 = next(iter(s for s in species if s.label=="CHOX3")) + surface_atom_count = chox3.molecule[0].get_num_atoms("X") assert surface_atom_count == 3 save_chemkin_surface_file( chemkin_save_path, @@ -817,17 +817,12 @@ def test_write_bidentate_species(self): check_for_duplicates=False, ) - bidentate_test = " CH2OX2(52)/2/ \n" - tridentate_test = " CHOX3(61)/3/ \n" + bidentate_test = "CH2OX2(52)/2/" + tridentate_test = "CHOX3(61)/3/" with open(chemkin_save_path, "r") as f: - for i, line in enumerate(f): - if i == 3: - bidentate_read = line - if i == 4: - tridentate_read = line - - assert bidentate_test.strip() == bidentate_read.strip() - assert tridentate_test.strip() == tridentate_read.strip() + lines = [line.strip() for line in f] + assert bidentate_test in lines + assert tridentate_test in lines os.remove(chemkin_save_path) From 68ea1f4023bbafb845127027b55267d60b2a3b50 Mon Sep 17 00:00:00 2001 From: Richard West Date: Fri, 2 Aug 2024 00:54:12 -0400 Subject: [PATCH 03/14] Cleaning up and extending a chemkin test mechanism. To enable proper conversion in Cantera (for a future test) we need some gas phase parts to the chemkin mechanism. I added two extra types of surface reaction, an adsorption sticking coefficient, and an associative desorption. I also tidied up the species names and formatting. --- test/rmgpy/chemkinTest.py | 4 +- .../chemkin/chemkin_py/surface/chem-gas.inp | 54 +++++++++++++++++ .../chemkin_py/surface/chem-surface.inp | 50 ++++++++-------- .../chemkin_py/surface/species_dictionary.txt | 58 +++++++++++++++++-- 4 files changed, 135 insertions(+), 31 deletions(-) create mode 100644 test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-gas.inp diff --git a/test/rmgpy/chemkinTest.py b/test/rmgpy/chemkinTest.py index 64b02b2d8ed..69692ea4bbf 100644 --- a/test/rmgpy/chemkinTest.py +++ b/test/rmgpy/chemkinTest.py @@ -817,8 +817,8 @@ def test_write_bidentate_species(self): check_for_duplicates=False, ) - bidentate_test = "CH2OX2(52)/2/" - tridentate_test = "CHOX3(61)/3/" + bidentate_test = "CH2OX2/2/" + tridentate_test = "CHOX3/3/" with open(chemkin_save_path, "r") as f: lines = [line.strip() for line in f] assert bidentate_test in lines diff --git a/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-gas.inp b/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-gas.inp new file mode 100644 index 00000000000..0815daee86f --- /dev/null +++ b/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-gas.inp @@ -0,0 +1,54 @@ +ELEMENTS H C O N Ar END + +SPECIES + Ar + N2 + ethane + CH3 + CH4 + O2 +END + + +THERM ALL + 300.000 1000.000 5000.000 + +Ar Ar1 G200.000 6000.000 1000.00 1 + 2.50000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 0.00000000E+00 2 +-7.45375000E+02 4.37967000E+00 2.50000000E+00 0.00000000E+00 0.00000000E+00 3 + 0.00000000E+00 0.00000000E+00-7.45375000E+02 4.37967000E+00 4 + +N2 N 2 G200.000 6000.000 1000.00 1 + 2.95258000E+00 1.39690000E-03-4.92632000E-07 7.86010000E-11-4.60755000E-15 2 +-9.23949000E+02 5.87189000E+00 3.53101000E+00-1.23661000E-04-5.02999000E-07 3 + 2.43531000E-09-1.40881000E-12-1.04698000E+03 2.96747000E+00 4 + +ethane C 2 H 6 G100.000 5000.000 954.51 1 + 4.58979542E+00 1.41508364E-02-4.75965787E-06 8.60302924E-10-6.21723861E-14 2 +-1.27217507E+04-3.61718979E+00 3.78034578E+00-3.24276131E-03 5.52385395E-05 3 +-6.38587729E-08 2.28639990E-11-1.16203414E+04 5.21029728E+00 4 + +CH3 C 1 H 3 G100.000 5000.000 1337.62 1 + 3.54143640E+00 4.76790043E-03-1.82150130E-06 3.28880372E-10-2.22548587E-14 2 + 1.62239681E+04 1.66047100E+00 3.91546905E+00 1.84152818E-03 3.48746163E-06 3 +-3.32752224E-09 8.49972445E-13 1.62856393E+04 3.51736167E-01 4 + +CH4 C 1H 4 G 100.000 5000.000 1084.12 1 + 9.08256809E-01 1.14541005E-02-4.57174656E-06 8.29193594E-10-5.66316470E-14 2 +-9.71997053E+03 1.39931449E+01 4.20541679E+00-5.35559144E-03 2.51123865E-05 3 +-2.13763581E-08 5.97526898E-12-1.01619434E+04-9.21284857E-01 4 + +O2 O 2 G 100.000 5000.000 1074.56 1 + 3.15382774E+00 1.67803224E-03-7.69967755E-07 1.51273954E-10-1.08781177E-14 2 +-1.04082031E+03 6.16751905E+00 3.53732118E+00-1.21570202E-03 5.31615358E-06 3 +-4.89440364E-09 1.45843807E-12-1.03858843E+03 4.68368633E+00 4 + +END + + +REACTIONS KCAL/MOLE MOLES + +CH3 +CH3 = ethane 8.260e+17 -1.400 1.000 + +END + diff --git a/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-surface.inp b/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-surface.inp index 29dfedd651f..992e28bb709 100644 --- a/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-surface.inp +++ b/test/rmgpy/test_data/chemkin/chemkin_py/surface/chem-surface.inp @@ -1,61 +1,63 @@ ELEMENTS H - D /2.014/ - T /3.016/ C - CI /13.003/ O - OI /18.000/ - N - Ne Ar - He - Si - S - F - Cl - Br - I X /195.083/ END -SPECIES - X(1) - H*(10) - CH2OX2(52) - CHOX3(61) +SITE SDEN/2.4830E-09/ ! mol/cm^2 + X + HX + OX + CH4X + CH2OX2/2/ + CHOX3/3/ END THERM ALL 300.000 1000.000 5000.000 -X(1) X 1 G 100.000 5000.000 1554.80 1 +X X 1 G 100.000 5000.000 1554.80 1 1.60299900E-01-2.52235409E-04 1.14181275E-07-1.21471653E-11 3.85790025E-16 2 -7.08100885E+01-9.09527530E-01 7.10139498E-03-4.25619522E-05 8.98533016E-08 3 -7.80193649E-11 2.32465471E-14-8.76101712E-01-3.11211229E-02 4 -H*(10) H 1X 1 G 100.000 5000.000 952.91 1 +HX H 1X 1 G 100.000 5000.000 952.91 1 2.80339655E+00-5.41047017E-04 4.99507978E-07-7.54963647E-11 3.06772366E-15 2 -2.34636021E+03-1.59436787E+01-3.80965452E-01 5.47228709E-03 2.60912778E-06 3 -9.64961980E-09 4.63946753E-12-1.40561079E+03 1.01725550E+00 4 -CH2OX2(52) C 1H 2O 1X 2G 100.000 5000.000 777.56 1 +CH2OX2 C 1H 2O 1X 2G 100.000 5000.000 777.56 1 3.18259452E+00 1.04657053E-02-5.00464545E-06 9.79545329E-10-6.90993489E-14 2 -1.82783093E+04-1.50760119E+01-1.44645549E+00 3.42803139E-02-5.09483795E-05 3 4.03732306E-08-1.27356452E-11-1.75584787E+04 6.09156343E+00 4 -CHOX3(61) C 1H 1O 1X 3G 100.000 5000.000 689.33 1 +CHOX3 C 1H 1O 1X 3G 100.000 5000.000 689.33 1 3.47343833E+00 7.03348332E-03-3.51779073E-06 6.98047945E-10-4.93382204E-14 2 -1.53390119E+04-1.79870016E+01-1.32716780E+00 3.17587333E-02-5.05065871E-05 3 3.95520303E-08-1.17505741E-11-1.46027731E+04 3.92679635E+00 4 -END +OX O 1X 1 G 298.000 2000.000 1000.00 1 + 2.90244691E+00-3.38584457E-04 6.43372619E-07-3.66326660E-10 6.90093884E-14 2 +-1.70497471E+04-1.52559728E+01-2.94475701E-01 1.44162624E-02-2.61322704E-05 3 + 2.19005957E-08-6.98019420E-12-1.64619234E+04-1.99445623E-01 4 +CH4X C 1H 4X 1 G 298.000 2000.000 1000.00 1 + 9.54139378E+00-1.04025134E-02 1.83777400E-05-9.66765130E-09 1.71211379E-12 2 +-1.34475614E+04-3.55638434E+01 4.85496247E+00-5.54134984E-03 3.01198105E-05 3 +-2.99225917E-08 1.00502514E-11-1.17096278E+04-9.25620913E+00 4 +END REACTIONS KCAL/MOLE MOLES -X(1)+X(1)+CH2OX2(52)=H*(10)+CHOX3(61) 7.420000e+21 0.000 3.058 +X + X + CH2OX2 = HX + CHOX3 7.420000e+21 0.000 3.058 + +X + CH4 <=> CH4X 8.000e-03 0.000 0.000 + STICK + +OX + OX <=> X + X + O2 3.700000e+21 0.000 66.611 END diff --git a/test/rmgpy/test_data/chemkin/chemkin_py/surface/species_dictionary.txt b/test/rmgpy/test_data/chemkin/chemkin_py/surface/species_dictionary.txt index f6acc210007..c03ebd75413 100644 --- a/test/rmgpy/test_data/chemkin/chemkin_py/surface/species_dictionary.txt +++ b/test/rmgpy/test_data/chemkin/chemkin_py/surface/species_dictionary.txt @@ -1,11 +1,11 @@ -X(1) +X 1 X u0 p0 c0 -H*(10) +HX 1 H u0 p0 c0 {2,S} 2 X u0 p0 c0 {1,S} -CHOX3(61) +CHOX3 1 O u0 p2 c0 {2,S} {5,S} 2 C u0 p0 c0 {1,S} {3,S} {4,S} {6,S} 3 H u0 p0 c0 {2,S} @@ -13,10 +13,58 @@ CHOX3(61) 5 X u0 p0 c0 {1,S} 6 X u0 p0 c0 {2,S} -CH2OX2(52) +CH2OX2 1 O u0 p2 c0 {2,S} {6,S} 2 C u0 p0 c0 {1,S} {3,S} {4,S} {5,S} 3 H u0 p0 c0 {2,S} 4 H u0 p0 c0 {2,S} 5 X u0 p0 c0 {2,S} -6 X u0 p0 c0 {1,S} \ No newline at end of file +6 X u0 p0 c0 {1,S} + +Ar +1 Ar u0 p4 c0 + +N2 +1 N u0 p1 c0 {2,T} +2 N u0 p1 c0 {1,T} + +ethane +1 C u0 p0 c0 {2,S} {3,S} {4,S} {5,S} +2 C u0 p0 c0 {1,S} {6,S} {7,S} {8,S} +3 H u0 p0 c0 {1,S} +4 H u0 p0 c0 {1,S} +5 H u0 p0 c0 {1,S} +6 H u0 p0 c0 {2,S} +7 H u0 p0 c0 {2,S} +8 H u0 p0 c0 {2,S} + +CH3 +multiplicity 2 +1 C u1 p0 c0 {2,S} {3,S} {4,S} +2 H u0 p0 c0 {1,S} +3 H u0 p0 c0 {1,S} +4 H u0 p0 c0 {1,S} + +O2 +multiplicity 3 +1 O u1 p2 c0 {2,S} +2 O u1 p2 c0 {1,S} + +CH4 +1 C u0 p0 c0 {2,S} {3,S} {4,S} {5,S} +2 H u0 p0 c0 {1,S} +3 H u0 p0 c0 {1,S} +4 H u0 p0 c0 {1,S} +5 H u0 p0 c0 {1,S} + +CH4X +1 C u0 p0 c0 {2,S} {3,S} {4,S} {5,S} +2 H u0 p0 c0 {1,S} +3 H u0 p0 c0 {1,S} +4 H u0 p0 c0 {1,S} +5 H u0 p0 c0 {1,S} +6 X u0 p0 c0 + +OX +1 O u0 p2 c0 {2,D} +2 X u0 p0 c0 {1,D} \ No newline at end of file From b3906490e68dcb32bb2fdd944bb93c689bdef921 Mon Sep 17 00:00:00 2001 From: Nora Khalil Date: Mon, 13 Jun 2022 12:15:01 -0400 Subject: [PATCH 04/14] Updating methods that create in-memory Cantera objects of things. to_cantera for Reaction and Species and some surface kinetics types squashed commit from Richard: Should be able to use the to_cantera_kinetics methods. Not sure about all these nested elifs. I think some are confused. Co-authored-by: Nora Khalil Co-authored-by: Richard West --- rmgpy/kinetics/surface.pyx | 54 ++++++++++++++++++++++++++++++++++++++ rmgpy/reaction.py | 9 +++++++ rmgpy/species.py | 12 +++++++-- 3 files changed, 73 insertions(+), 2 deletions(-) diff --git a/rmgpy/kinetics/surface.pyx b/rmgpy/kinetics/surface.pyx index 679f7bb3389..cdabf620ecf 100644 --- a/rmgpy/kinetics/surface.pyx +++ b/rmgpy/kinetics/surface.pyx @@ -285,6 +285,34 @@ cdef class StickingCoefficient(KineticsModel): return string + def to_cantera_kinetics(self): + """ + Converts the Arrhenius object to a cantera Arrhenius object + + Arrhenius(A,b,E) where A is in units of m^3/kmol/s, b is dimensionless, and E is in J/kmol + """ + + import cantera as ct + + A = self._A.value_si + b = self._n.value_si + Ea = self._Ea.value_si * 1000 # convert from J/mol to J/kmol + + return ct.StickingArrheniusRate(A, b, Ea) + + + def set_cantera_kinetics(self, ct_reaction, species_list): + """ + Passes in a cantera Reaction() object and sets its + rate using to_cantera_kinetics(). + """ + import cantera as ct + + # Set the rate parameter to a cantera Arrhenius() object + ct_reaction.rate = self.to_cantera_kinetics() + + + ################################################################################ cdef class StickingCoefficientBEP(KineticsModel): """ @@ -589,6 +617,32 @@ cdef class SurfaceArrhenius(Arrhenius): ) + def to_cantera_kinetics(self): + """ + Converts the Arrhenius object to a cantera Arrhenius object + + Arrhenius(A,b,E) where A is in units of m^3/kmol/s, b is dimensionless, and E is in J/kmol + """ + + import cantera as ct + + A = self._A.value_si + b = self._n.value_si + Ea = self._Ea.value_si * 1000 # convert from J/mol to J/kmol + + return ct.Arrhenius(A, b, Ea) + + + def set_cantera_kinetics(self, ct_reaction, species_list): + """ + Passes in a cantera Reaction() object and sets its + rate using to_cantera_kinetics(). + """ + import cantera as ct + + # Set the rate parameter to a cantera Arrhenius() object + ct_reaction.rate = self.to_cantera_kinetics() + ################################################################################ cdef class SurfaceArrheniusBEP(ArrheniusEP): diff --git a/rmgpy/reaction.py b/rmgpy/reaction.py index 19136c60faf..64b34f457f5 100644 --- a/rmgpy/reaction.py +++ b/rmgpy/reaction.py @@ -268,6 +268,7 @@ def to_chemkin(self, species_list=None, kinetics=True): else: return rmgpy.chemkin.write_reaction_string(self) + def to_cantera(self, species_list=None, use_chemkin_identifier=False): """ Converts the RMG Reaction object to a Cantera Reaction object @@ -362,6 +363,14 @@ def to_cantera(self, species_list=None, use_chemkin_identifier=False): reactants=ct_reactants, products=ct_products, rate=rate ) + elif isinstance(self.kinetics, SurfaceArrhenius): + rate = self.kinetics.to_cantera_kinetics() + ct_reaction = ct.InterfaceReaction(equation=str(self), rate=rate) + + elif isinstance(self.kinetics, StickingCoefficient): + rate = self.kinetics.to_cantera_kinetics() + ct_reaction = ct.Reaction(equation=str(self), rate=rate) + elif isinstance(self.kinetics, Lindemann): high_rate = self.kinetics.arrheniusHigh.to_cantera_kinetics(arrhenius_class=True) low_rate = self.kinetics.arrheniusLow.to_cantera_kinetics(arrhenius_class=True) diff --git a/rmgpy/species.py b/rmgpy/species.py index 2969dac82b9..e2a9ff25814 100644 --- a/rmgpy/species.py +++ b/rmgpy/species.py @@ -451,9 +451,17 @@ def to_cantera(self, use_chemkin_identifier=False): else: element_dict[symbol] += 1 if use_chemkin_identifier: - ct_species = ct.Species(self.to_chemkin(), element_dict) + label = self.to_chemkin() else: - ct_species = ct.Species(self.label, element_dict) + label = self.label + + if self.contains_surface_site() and element_dict["X"] > 1: + # for multidentate adsorbates, 'size' is the same as 'sites'? for some reason, cantera won't take the input 'size,' so will need to use 'sites' + ct_species = ct.Species(label, element_dict, size=element_dict["X"]) + # hopefully this will be fixed soon, so that ct.Species can take a 'sites' parameter or that cantera can read input files with 'size' specified + else: + ct_species = ct.Species(label, element_dict) + if self.thermo: try: ct_species.thermo = self.thermo.to_cantera() From 46628d77a144fc93271dbf6e202878d509578f69 Mon Sep 17 00:00:00 2001 From: Richard West Date: Thu, 6 Jun 2024 10:02:06 -0400 Subject: [PATCH 05/14] Fixing docstrings in surface kinetics. --- rmgpy/kinetics/surface.pyx | 13 +++++++------ 1 file changed, 7 insertions(+), 6 deletions(-) diff --git a/rmgpy/kinetics/surface.pyx b/rmgpy/kinetics/surface.pyx index cdabf620ecf..5491c5ccca1 100644 --- a/rmgpy/kinetics/surface.pyx +++ b/rmgpy/kinetics/surface.pyx @@ -50,7 +50,7 @@ cdef class StickingCoefficient(KineticsModel): ======================= ============================================================= Attribute Description ======================= ============================================================= - `A` The preexponential factor + `A` The preexponential factor. Unitless (sticking probability) `T0` The reference temperature `n` The temperature exponent `Ea` The activation energy @@ -287,9 +287,9 @@ cdef class StickingCoefficient(KineticsModel): def to_cantera_kinetics(self): """ - Converts the Arrhenius object to a cantera Arrhenius object + Converts the Arrhenius object to a cantera StickingArrheniusRate object - Arrhenius(A,b,E) where A is in units of m^3/kmol/s, b is dimensionless, and E is in J/kmol + StickingArrheniusRate(A,b,E) where A is dimensionless, b is dimensionless, and E is in J/kmol """ import cantera as ct @@ -380,7 +380,7 @@ cdef class StickingCoefficientBEP(KineticsModel): self.Pmin, self.Pmax, self.coverage_dependence, self.comment)) property A: - """The preexponential factor.""" + """The preexponential factor. Dimensionless (sticking probability)""" def __get__(self): return self._A def __set__(self, value): @@ -423,7 +423,8 @@ cdef class StickingCoefficientBEP(KineticsModel): cpdef double get_sticking_coefficient(self, double T, double dHrxn=0.0) except -1: """ Return the sticking coefficient (dimensionless) at - temperature `T` in K and enthalpy of reaction `dHrxn` in J/mol. + temperature `T` in K and enthalpy of reaction `dHrxn` in J/mol. + Never exceeds 1.0. """ cdef double A, n, Ea, stickingCoefficient Ea = self.get_activation_energy(dHrxn) @@ -552,7 +553,7 @@ cdef class SurfaceArrhenius(Arrhenius): property A: """The preexponential factor. - This is the only thing different from a normal Arrhenius class.""" + This (and the coverage dependence) is the only thing different from a normal Arrhenius class.""" def __get__(self): return self._A def __set__(self, value): From 9fa619541753b3b53b2456cb032f9338594d2922 Mon Sep 17 00:00:00 2001 From: Sevy Harris Date: Wed, 8 Nov 2023 13:14:50 -0500 Subject: [PATCH 06/14] Added to_cantera_kinetics for SurfaceArrhenius and StickingCoefficient Nick add cm^5/(mol^2*s) dimensionality Co-authored-by: Sevy Harris Co-authored-by: Nicholas Tietje --- rmgpy/kinetics/surface.pyx | 86 ++++++++++++++++++++++++-------------- rmgpy/reaction.py | 8 +++- 2 files changed, 62 insertions(+), 32 deletions(-) diff --git a/rmgpy/kinetics/surface.pyx b/rmgpy/kinetics/surface.pyx index 5491c5ccca1..7787a82e479 100644 --- a/rmgpy/kinetics/surface.pyx +++ b/rmgpy/kinetics/surface.pyx @@ -287,31 +287,31 @@ cdef class StickingCoefficient(KineticsModel): def to_cantera_kinetics(self): """ - Converts the Arrhenius object to a cantera StickingArrheniusRate object - - StickingArrheniusRate(A,b,E) where A is dimensionless, b is dimensionless, and E is in J/kmol + Converts the RMG StickingCoefficient object to a cantera StickingArrheniusRate + + StickingArrheniusRate(A,b,E) where A is dimensionless, b is dimensionless, and E is in J/kmol """ - import cantera as ct - + import rmgpy.quantity + if type(self._A) != rmgpy.quantity.ScalarQuantity: + raise TypeError("A factor must be a dimensionless ScalarQuantity") A = self._A.value_si b = self._n.value_si - Ea = self._Ea.value_si * 1000 # convert from J/mol to J/kmol + E = self._Ea.value_si * 1000 # convert from J/mol to J/kmol - return ct.StickingArrheniusRate(A, b, Ea) + return ct.StickingArrheniusRate(A, b, E) def set_cantera_kinetics(self, ct_reaction, species_list): - """ - Passes in a cantera Reaction() object and sets its - rate using to_cantera_kinetics(). - """ - import cantera as ct - - # Set the rate parameter to a cantera Arrhenius() object - ct_reaction.rate = self.to_cantera_kinetics() + """ + Passes in a cantera Reaction() object and sets its + rate to a Cantera ArrheniusRate object. + """ + import cantera as ct + assert isinstance(ct_reaction.rate, ct.StickingArrheniusRate), "Must have a Cantera StickingArrheniusRate attribute" - + # Set the rate parameter to a cantera Arrhenius object + ct_reaction.rate = self.to_cantera_kinetics() ################################################################################ cdef class StickingCoefficientBEP(KineticsModel): @@ -620,29 +620,53 @@ cdef class SurfaceArrhenius(Arrhenius): def to_cantera_kinetics(self): """ - Converts the Arrhenius object to a cantera Arrhenius object + Converts the RMG SurfaceArrhenius object to a cantera InterfaceArrheniusRate - Arrhenius(A,b,E) where A is in units of m^3/kmol/s, b is dimensionless, and E is in J/kmol + InterfaceArrheniusRate(A,b,E) where A is in units like m^2/kmol/s (depending on dimensionality) + b is dimensionless, and E is in J/kmol """ - import cantera as ct - A = self._A.value_si - b = self._n.value_si - Ea = self._Ea.value_si * 1000 # convert from J/mol to J/kmol + rate_units_conversion = {'1/s': 1, + 's^-1': 1, + 'm^2/(mol*s)': 1000, + 'm^4/(mol^2*s)': 1000000, + 'cm^2/(mol*s)': 1000, + 'cm^4/(mol^2*s)': 1000000, + 'm^2/(molecule*s)': 1000, + 'm^4/(molecule^2*s)': 1000000, + 'cm^2/(molecule*s)': 1000, + 'cm^4/(molecule^2*s)': 1000000, + 'cm^5/(mol^2*s)': 1000000, + 'm^5/(mol^2*s)': 1000000, + } + + if self._T0.value_si != 1: + A = self._A.value_si / (self._T0.value_si) ** self._n.value_si + else: + A = self._A.value_si - return ct.Arrhenius(A, b, Ea) + try: + A *= rate_units_conversion[self._A.units] # convert from /mol to /kmol + except KeyError: + raise ValueError('Arrhenius A-factor units {0} not found among accepted units for converting to ' + 'Cantera Arrhenius object.'.format(self._A.units)) + b = self._n.value_si + E = self._Ea.value_si * 1000 # convert from J/mol to J/kmol + return ct.InterfaceArrheniusRate(A, b, E) def set_cantera_kinetics(self, ct_reaction, species_list): - """ - Passes in a cantera Reaction() object and sets its - rate using to_cantera_kinetics(). - """ - import cantera as ct - - # Set the rate parameter to a cantera Arrhenius() object - ct_reaction.rate = self.to_cantera_kinetics() + """ + Takes in a cantera Reaction object and sets its + rate to a cantera InterfaceArrheniusRate object. + """ + import cantera as ct + if not isinstance(ct_reaction.rate, ct.InterfaceArrheniusRate): + raise TypeError("ct_reaction.rate must be an InterfaceArrheniusRate") + + # Set the rate parameter to a cantera Arrhenius object + ct_reaction.rate = self.to_cantera_kinetics() ################################################################################ diff --git a/rmgpy/reaction.py b/rmgpy/reaction.py index 64b34f457f5..48bfd8a63a0 100644 --- a/rmgpy/reaction.py +++ b/rmgpy/reaction.py @@ -312,7 +312,10 @@ def to_cantera(self, species_list=None, use_chemkin_identifier=False): if self.kinetics: if isinstance(self.kinetics, Arrhenius): # Create an Elementary Reaction - ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.ArrheniusRate()) + if isinstance(self.kinetics, SurfaceArrhenius): # SurfaceArrhenius inherits from Arrhenius + ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.InterfaceArrheniusRate()) + else: + ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.ArrheniusRate()) elif isinstance(self.kinetics, MultiArrhenius): # Return a list of elementary reactions which are duplicates ct_reaction = [ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.ArrheniusRate()) @@ -388,6 +391,9 @@ def to_cantera(self, species_list=None, use_chemkin_identifier=False): reactants=ct_reactants, products=ct_products, rate=rate ) + elif isinstance(self.kinetics, StickingCoefficient): + ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.StickingArrheniusRate()) + else: raise NotImplementedError('Unable to set cantera kinetics for {0}'.format(self.kinetics)) From fed3a1e5bd7bc7c2580050f02ebcc0a27f569c52 Mon Sep 17 00:00:00 2001 From: Richard West Date: Thu, 1 Aug 2024 21:18:11 -0400 Subject: [PATCH 07/14] Adding to_cantera_kinetics for ArrheniusBM class. Cantera now has a BlowersMasel kinetics class. This requires Cantera > 2.6 As described in https://cantera.org/dev/python/kinetics.html#cantera.BlowersMaselRate the API may change suddenly on some versions. --- rmgpy/kinetics/arrhenius.pyx | 46 +++++++++++++++++++++++++++++++++--- 1 file changed, 43 insertions(+), 3 deletions(-) diff --git a/rmgpy/kinetics/arrhenius.pyx b/rmgpy/kinetics/arrhenius.pyx index f4b69179ed9..70508761e39 100644 --- a/rmgpy/kinetics/arrhenius.pyx +++ b/rmgpy/kinetics/arrhenius.pyx @@ -722,12 +722,52 @@ cdef class ArrheniusBM(KineticsModel): """ self._A.value_si *= factor + def to_cantera_kinetics(self): + """ + Converts the RMG ArrheniusBM object to a cantera BlowersMaselRate. + + BlowersMaselRate(A, b, Ea, W) where A is in units of m^3/kmol/s, + b is dimensionless, and Ea and W are in J/kmol + """ + import cantera as ct + + rate_units_conversion = {'1/s': 1, + 's^-1': 1, + 'm^3/(mol*s)': 1000, + 'm^6/(mol^2*s)': 1000000, + 'cm^3/(mol*s)': 1000, + 'cm^6/(mol^2*s)': 1000000, + 'm^3/(molecule*s)': 1000, + 'm^6/(molecule^2*s)': 1000000, + 'cm^3/(molecule*s)': 1000, + 'cm^6/(molecule^2*s)': 1000000, + } + + if self._T0.value_si != 1: + A = self._A.value_si / (self._T0.value_si) ** self._n.value_si + else: + A = self._A.value_si + + try: + A *= rate_units_conversion[self._A.units] # convert from /mol to /kmol + except KeyError: + raise ValueError(f'ArrheniusBM A-factor units {self._A.units} not found among accepted ' + 'units for converting to Cantera BlowersMaselRate object.') + + b = self._n.value_si + Ea = self._E0.value_si * 1000 # convert from J/mol to J/kmol + w = self._w0.value_si * 1000 # convert from J/mol to J/kmol + + return ct.BlowersMaselRate(A, b, Ea, w) + def set_cantera_kinetics(self, ct_reaction, species_list): """ - Sets a cantera Reaction() object with the modified Arrhenius object - converted to an Arrhenius form. + Accepts a cantera Reaction object and sets its rate to a Cantera BlowersMaselRate object. """ - raise NotImplementedError('set_cantera_kinetics() is not implemented for ArrheniusBM class kinetics.') + import cantera as ct + if not isinstance(ct_reaction.rate, ct.BlowersMaselRate): + raise TypeError("ct_reaction must have a cantera BlowersMaselRate as the rate attribute") + ct_reaction.rate = self.to_cantera_kinetics() ################################################################################ From 8f3b2f4c19675d975014a6ec44faabba4b56b5f3 Mon Sep 17 00:00:00 2001 From: Richard West Date: Fri, 2 Aug 2024 00:57:06 -0400 Subject: [PATCH 08/14] Testing conversion to Cantera for adsorbed species and reactions. This is in the canteramodelTest. We load a gas/surface mechanism saved in chemkin format (updated in prior commit) into both RMG and into Cantera. We then convert the RMG objects into Cantera objects and then compare them all - species and reactions. This mirrors the way the gas phase parts were tested. Hopefully increases test coverage of the new to_cantera_kinetics methods. --- test/rmgpy/tools/canteramodelTest.py | 45 +++++++++++++++++++++++++++- 1 file changed, 44 insertions(+), 1 deletion(-) diff --git a/test/rmgpy/tools/canteramodelTest.py b/test/rmgpy/tools/canteramodelTest.py index c2e1b1d57cf..81b7d6844e0 100644 --- a/test/rmgpy/tools/canteramodelTest.py +++ b/test/rmgpy/tools/canteramodelTest.py @@ -104,6 +104,27 @@ def setup_class(self): self.ctSpecies = job.model.species() self.ctReactions = job.model.reactions() + # Now load surface species and kinetics + folder = os.path.join(os.path.dirname(os.path.dirname(__file__)), "test_data", "chemkin", "chemkin_py") + chemkin_path = os.path.join(folder, "surface", "chem-gas.inp") + chemkin_surface_path = os.path.join(folder, "surface", "chem-surface.inp") + dictionary_path = os.path.join(folder, "surface", "species_dictionary.txt") + species, reactions = load_chemkin_file(chemkin_surface_path, dictionary_path) + self.rmg_surface_ct_species = [spec.to_cantera(use_chemkin_identifier=True) for spec in species] + self.rmg_surface_ct_reactions = [] + for rxn in reactions: + converted_reactions = rxn.to_cantera(species, use_chemkin_identifier=True) + if isinstance(converted_reactions, list): + self.rmg_surface_ct_reactions.extend(converted_reactions) + else: + self.rmg_surface_ct_reactions.append(converted_reactions) + job = Cantera() + job.surface = True + job.load_chemkin_model(chemkin_path, surface_file=chemkin_surface_path, quiet=True) + self.ct_surface_species = job.surface.species() + self.ct_surface_reactions = job.surface.reactions() + + def test_species_conversion(self): """ Test that species objects convert properly @@ -115,9 +136,31 @@ def test_species_conversion(self): def test_reaction_conversion(self): """ - Test that species objects convert properly + Test that reaction objects convert properly """ from rmgpy.tools.canteramodel import check_equivalent_cantera_reaction for i in range(len(self.ctReactions)): assert check_equivalent_cantera_reaction(self.ctReactions[i], self.rmg_ctReactions[i]) + + def test_surface_species_conversion(self): + """ + Test that surface species objects convert properly + """ + from rmgpy.tools.canteramodel import check_equivalent_cantera_species + + for i in range(len(self.ct_surface_species)): + #print("Chemkin-to-Cantera:", self.ct_surfaceSpecies[i].input_data) + #print("Chemkin-to-RMG-to-Cantera:", self.rmg_surface_ctSpecies[i].input_data) + assert check_equivalent_cantera_species(self.ct_surface_species[i], self.rmg_surface_ct_species[i]) + + def test_surface_reaction_conversion(self): + """ + Test that surface reaction objects convert properly + """ + from rmgpy.tools.canteramodel import check_equivalent_cantera_reaction + + for i in range(len(self.ct_surface_reactions)): + #print("Chemkin-to-Cantera:", self.ct_surfaceReactions[i].input_data) + #print("Chemkin-to-RMG-to-Cantera:", self.rmg_surface_ctReactions[i].input_data) + assert check_equivalent_cantera_reaction(self.ct_surface_reactions[i], self.rmg_surface_ct_reactions[i]) \ No newline at end of file From b84e568c417ffd0fe07ace00f11c570adc09a78e Mon Sep 17 00:00:00 2001 From: Sevy Harris Date: Wed, 27 Nov 2024 14:13:06 -0500 Subject: [PATCH 09/14] Fixup rmgpy/kinetics/arrhenius.pyx T0 isn't an attribute of ArrheniusBM (just Arrhenius), so this should be removed. Co-authored-by: Richard West --- rmgpy/kinetics/arrhenius.pyx | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/rmgpy/kinetics/arrhenius.pyx b/rmgpy/kinetics/arrhenius.pyx index 70508761e39..ea9e462d4f8 100644 --- a/rmgpy/kinetics/arrhenius.pyx +++ b/rmgpy/kinetics/arrhenius.pyx @@ -743,10 +743,7 @@ cdef class ArrheniusBM(KineticsModel): 'cm^6/(molecule^2*s)': 1000000, } - if self._T0.value_si != 1: - A = self._A.value_si / (self._T0.value_si) ** self._n.value_si - else: - A = self._A.value_si + A = self._A.value_si try: A *= rate_units_conversion[self._A.units] # convert from /mol to /kmol From cac63d1cdcf4c3f4140a4f1fa8d22c8ee28003d6 Mon Sep 17 00:00:00 2001 From: Sevy Harris Date: Wed, 27 Nov 2024 14:15:47 -0500 Subject: [PATCH 10/14] Fix rmgpy/reaction.py Since SurfaceArrhenius is a subclass of Arrhenius this elif will never be true. SurfaceArrhenius is dealt with on line 292 above. Co-authored-by: Richard West --- rmgpy/reaction.py | 3 --- 1 file changed, 3 deletions(-) diff --git a/rmgpy/reaction.py b/rmgpy/reaction.py index 48bfd8a63a0..12c2dfeb57f 100644 --- a/rmgpy/reaction.py +++ b/rmgpy/reaction.py @@ -366,9 +366,6 @@ def to_cantera(self, species_list=None, use_chemkin_identifier=False): reactants=ct_reactants, products=ct_products, rate=rate ) - elif isinstance(self.kinetics, SurfaceArrhenius): - rate = self.kinetics.to_cantera_kinetics() - ct_reaction = ct.InterfaceReaction(equation=str(self), rate=rate) elif isinstance(self.kinetics, StickingCoefficient): rate = self.kinetics.to_cantera_kinetics() From d8b93147b8e3ce0e92a23d58dc01eae81d3e7ef1 Mon Sep 17 00:00:00 2001 From: Richard West Date: Wed, 27 Nov 2024 15:32:34 -0500 Subject: [PATCH 11/14] In Reaction.to_cantera, don't make kinetics twice. The method calls the appropriate set_cantera_kinetics() method at the end, so we don't want to make them before then, because it's wasteful to do it twice. We just need to create a cantera kinetics object of the correct type. This was already done for several types of kinetics, but not all. --- rmgpy/reaction.py | 52 +++++++++++++++++++---------------------------- 1 file changed, 21 insertions(+), 31 deletions(-) diff --git a/rmgpy/reaction.py b/rmgpy/reaction.py index 12c2dfeb57f..810f8a6d75e 100644 --- a/rmgpy/reaction.py +++ b/rmgpy/reaction.py @@ -310,6 +310,9 @@ def to_cantera(self, species_list=None, use_chemkin_identifier=False): ct_collider[self.specific_collider.to_chemkin() if use_chemkin_identifier else self.specific_collider.label] = 1 if self.kinetics: + # Create the Cantera reaction object, + # with the correct type of kinetics object + # but don't actually set its kinetics (we do that at the end) if isinstance(self.kinetics, Arrhenius): # Create an Elementary Reaction if isinstance(self.kinetics, SurfaceArrhenius): # SurfaceArrhenius inherits from Arrhenius @@ -338,61 +341,47 @@ def to_cantera(self, species_list=None, use_chemkin_identifier=False): ct_reaction = ct.ThreeBodyReaction(reactants=ct_reactants, products=ct_products) elif isinstance(self.kinetics, Troe): - high_rate = self.kinetics.arrheniusHigh.to_cantera_kinetics(arrhenius_class=True) - low_rate = self.kinetics.arrheniusLow.to_cantera_kinetics(arrhenius_class=True) - A = self.kinetics.alpha - T3 = self.kinetics.T3.value_si - T1 = self.kinetics.T1.value_si - - if self.kinetics.T2 is None: - rate = ct.TroeRate( - high=high_rate, low=low_rate, falloff_coeffs=[A, T3, T1] - ) - else: - T2 = self.kinetics.T2.value_si - rate = ct.TroeRate( - high=high_rate, low=low_rate, falloff_coeffs=[A, T3, T1, T2] - ) - if ct_collider is not None: ct_reaction = ct.FalloffReaction( reactants=ct_reactants, products=ct_products, tbody=ct_collider, - rate=rate, + rate=ct.TroeRate() ) else: ct_reaction = ct.FalloffReaction( - reactants=ct_reactants, products=ct_products, rate=rate + reactants=ct_reactants, + products=ct_products, + rate=ct.TroeRate() ) - - elif isinstance(self.kinetics, StickingCoefficient): - rate = self.kinetics.to_cantera_kinetics() - ct_reaction = ct.Reaction(equation=str(self), rate=rate) - elif isinstance(self.kinetics, Lindemann): - high_rate = self.kinetics.arrheniusHigh.to_cantera_kinetics(arrhenius_class=True) - low_rate = self.kinetics.arrheniusLow.to_cantera_kinetics(arrhenius_class=True) - falloff = [] - rate = ct.LindemannRate(low_rate, high_rate, falloff) if ct_collider is not None: ct_reaction = ct.FalloffReaction( reactants=ct_reactants, products=ct_products, tbody=ct_collider, - rate=rate, + rate=ct.LindemannRate() ) else: ct_reaction = ct.FalloffReaction( - reactants=ct_reactants, products=ct_products, rate=rate + reactants=ct_reactants, + products=ct_products, + rate=ct.LindemannRate() ) + elif isinstance(self.kinetics, SurfaceArrhenius): + ct_reaction = ct.InterfaceReaction(reactants=ct_reactants, + products=ct_products, + rate=ct.InterfaceArrheniusRate()) + elif isinstance(self.kinetics, StickingCoefficient): - ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.StickingArrheniusRate()) + ct_reaction = ct.Reaction(reactants=ct_reactants, + products=ct_products, + rate=ct.StickingArrheniusRate()) else: - raise NotImplementedError('Unable to set cantera kinetics for {0}'.format(self.kinetics)) + raise NotImplementedError(f"Unable to set cantera kinetics for {self.kinetics}") # Set reversibility, duplicate, and ID attributes if isinstance(ct_reaction, list): @@ -407,6 +396,7 @@ def to_cantera(self, species_list=None, use_chemkin_identifier=False): ct_reaction.duplicate = self.duplicate ct_reaction.ID = str(self.index) + # Now we set the kinetics. self.kinetics.set_cantera_kinetics(ct_reaction, species_list) return ct_reaction From 100d74004e9cd52aef12cee158a2aac44f70ccc8 Mon Sep 17 00:00:00 2001 From: Richard West Date: Wed, 27 Nov 2024 15:33:15 -0500 Subject: [PATCH 12/14] Reduce whitespace in kinetics/falloff.pyx --- rmgpy/kinetics/falloff.pyx | 10 ---------- 1 file changed, 10 deletions(-) diff --git a/rmgpy/kinetics/falloff.pyx b/rmgpy/kinetics/falloff.pyx index 082ffc048d1..12290667be7 100644 --- a/rmgpy/kinetics/falloff.pyx +++ b/rmgpy/kinetics/falloff.pyx @@ -220,15 +220,12 @@ cdef class Lindemann(PDepKineticsModel): self.arrheniusLow.change_rate(factor) self.arrheniusHigh.change_rate(factor) - - def set_cantera_kinetics(self, ct_reaction, species_list): """ Sets the efficiencies and kinetics for a cantera reaction. """ import cantera as ct assert isinstance(ct_reaction.rate, ct.LindemannRate), "Must have a Cantera LindemannRate attribute" - ct_reaction.efficiencies = PDepKineticsModel.get_cantera_efficiencies(self, species_list) ct_reaction.rate = self.to_cantera_kinetics() @@ -400,11 +397,8 @@ cdef class Troe(PDepKineticsModel): for a cantera FalloffReaction. """ import cantera as ct - assert isinstance(ct_reaction.rate, ct.TroeRate), "Must have a Cantera TroeRate attribute" - ct_reaction.efficiencies = PDepKineticsModel.get_cantera_efficiencies(self, species_list) - ct_reaction.rate = self.to_cantera_kinetics() def to_cantera_kinetics(self): @@ -424,7 +418,3 @@ cdef class Troe(PDepKineticsModel): high = self.arrheniusHigh.to_cantera_kinetics(arrhenius_class=True) low = self.arrheniusLow.to_cantera_kinetics(arrhenius_class=True) return ct.TroeRate(high=high, low=low, falloff_coeffs=falloff) - - - - \ No newline at end of file From a3a67a110a422dc89022a4180958a48593e1da52 Mon Sep 17 00:00:00 2001 From: Sevy Harris Date: Wed, 27 Nov 2024 15:48:23 -0500 Subject: [PATCH 13/14] Fix comments in species.py --- rmgpy/species.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/rmgpy/species.py b/rmgpy/species.py index e2a9ff25814..c17bdd182a4 100644 --- a/rmgpy/species.py +++ b/rmgpy/species.py @@ -456,9 +456,11 @@ def to_cantera(self, use_chemkin_identifier=False): label = self.label if self.contains_surface_site() and element_dict["X"] > 1: - # for multidentate adsorbates, 'size' is the same as 'sites'? for some reason, cantera won't take the input 'size,' so will need to use 'sites' + # for multidentate adsorbates, 'size' is the same as 'sites'? + # for some reason,cantera won't take the input 'sites' so will need to use 'size' ct_species = ct.Species(label, element_dict, size=element_dict["X"]) - # hopefully this will be fixed soon, so that ct.Species can take a 'sites' parameter or that cantera can read input files with 'size' specified + # hopefully this will be fixed soon, so that ct.Species can take a 'sites' parameter + # or that cantera can read input files with 'size' specified else: ct_species = ct.Species(label, element_dict) From ffa0ca1b9be94d967cecc322e4569aadd9dbb5ec Mon Sep 17 00:00:00 2001 From: Richard West Date: Wed, 27 Nov 2024 16:09:36 -0500 Subject: [PATCH 14/14] Simple refactor in Reaction.to_cantera() By putting the "if not self.kinetics check first" and raising the exception early, instead of having that at the end in the else block, the rest of the code doesn't need to be inside the "if" block. --- rmgpy/reaction.py | 173 ++++++++++++++++++++++++---------------------- 1 file changed, 89 insertions(+), 84 deletions(-) diff --git a/rmgpy/reaction.py b/rmgpy/reaction.py index 810f8a6d75e..b8ab3951c8c 100644 --- a/rmgpy/reaction.py +++ b/rmgpy/reaction.py @@ -309,100 +309,105 @@ def to_cantera(self, species_list=None, use_chemkin_identifier=False): if self.specific_collider: # add a specific collider if exists ct_collider[self.specific_collider.to_chemkin() if use_chemkin_identifier else self.specific_collider.label] = 1 - if self.kinetics: - # Create the Cantera reaction object, - # with the correct type of kinetics object - # but don't actually set its kinetics (we do that at the end) - if isinstance(self.kinetics, Arrhenius): - # Create an Elementary Reaction - if isinstance(self.kinetics, SurfaceArrhenius): # SurfaceArrhenius inherits from Arrhenius - ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.InterfaceArrheniusRate()) - else: - ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.ArrheniusRate()) - elif isinstance(self.kinetics, MultiArrhenius): - # Return a list of elementary reactions which are duplicates - ct_reaction = [ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.ArrheniusRate()) - for arr in self.kinetics.arrhenius] + if not self.kinetics: + raise Exception('Cantera reaction cannot be created because there was no kinetics.') - elif isinstance(self.kinetics, PDepArrhenius): - ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.PlogRate()) + # Create the Cantera reaction object, + # with the correct type of kinetics object + # but don't actually set its kinetics (we do that at the end) + if isinstance(self.kinetics, Arrhenius): + # Create an Elementary Reaction + if isinstance(self.kinetics, SurfaceArrhenius): # SurfaceArrhenius inherits from Arrhenius + ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.InterfaceArrheniusRate()) + else: + ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.ArrheniusRate()) + elif isinstance(self.kinetics, MultiArrhenius): + # Return a list of elementary reactions which are duplicates + ct_reaction = [ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.ArrheniusRate()) + for arr in self.kinetics.arrhenius] - elif isinstance(self.kinetics, MultiPDepArrhenius): - ct_reaction = [ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.PlogRate()) - for arr in self.kinetics.arrhenius] + elif isinstance(self.kinetics, PDepArrhenius): + ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.PlogRate()) - elif isinstance(self.kinetics, Chebyshev): - ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.ChebyshevRate()) + elif isinstance(self.kinetics, MultiPDepArrhenius): + ct_reaction = [ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.PlogRate()) + for arr in self.kinetics.arrhenius] - elif isinstance(self.kinetics, ThirdBody): - if ct_collider is not None: - ct_reaction = ct.ThreeBodyReaction(reactants=ct_reactants, products=ct_products, third_body=ct_collider) - else: - ct_reaction = ct.ThreeBodyReaction(reactants=ct_reactants, products=ct_products) - - elif isinstance(self.kinetics, Troe): - if ct_collider is not None: - ct_reaction = ct.FalloffReaction( - reactants=ct_reactants, - products=ct_products, - tbody=ct_collider, - rate=ct.TroeRate() - ) - else: - ct_reaction = ct.FalloffReaction( - reactants=ct_reactants, - products=ct_products, - rate=ct.TroeRate() - ) - - elif isinstance(self.kinetics, Lindemann): - if ct_collider is not None: - ct_reaction = ct.FalloffReaction( - reactants=ct_reactants, - products=ct_products, - tbody=ct_collider, - rate=ct.LindemannRate() - ) - else: - ct_reaction = ct.FalloffReaction( - reactants=ct_reactants, - products=ct_products, - rate=ct.LindemannRate() - ) - - elif isinstance(self.kinetics, SurfaceArrhenius): - ct_reaction = ct.InterfaceReaction(reactants=ct_reactants, - products=ct_products, - rate=ct.InterfaceArrheniusRate()) - - elif isinstance(self.kinetics, StickingCoefficient): - ct_reaction = ct.Reaction(reactants=ct_reactants, - products=ct_products, - rate=ct.StickingArrheniusRate()) + elif isinstance(self.kinetics, Chebyshev): + ct_reaction = ct.Reaction(reactants=ct_reactants, products=ct_products, rate=ct.ChebyshevRate()) + elif isinstance(self.kinetics, ThirdBody): + if ct_collider is not None: + ct_reaction = ct.ThreeBodyReaction(reactants=ct_reactants, products=ct_products, third_body=ct_collider) else: - raise NotImplementedError(f"Unable to set cantera kinetics for {self.kinetics}") - - # Set reversibility, duplicate, and ID attributes - if isinstance(ct_reaction, list): - for rxn in ct_reaction: - rxn.reversible = self.reversible - # Set the duplicate flag to true since this reaction comes from multiarrhenius or multipdeparrhenius - rxn.duplicate = True - # Set the ID flag to the original rmg index - rxn.ID = str(self.index) + ct_reaction = ct.ThreeBodyReaction(reactants=ct_reactants, products=ct_products) + + elif isinstance(self.kinetics, Troe): + if ct_collider is not None: + ct_reaction = ct.FalloffReaction( + reactants=ct_reactants, + products=ct_products, + tbody=ct_collider, + rate=ct.TroeRate() + ) else: - ct_reaction.reversible = self.reversible - ct_reaction.duplicate = self.duplicate - ct_reaction.ID = str(self.index) - - # Now we set the kinetics. - self.kinetics.set_cantera_kinetics(ct_reaction, species_list) + ct_reaction = ct.FalloffReaction( + reactants=ct_reactants, + products=ct_products, + rate=ct.TroeRate() + ) + + elif isinstance(self.kinetics, Lindemann): + if ct_collider is not None: + ct_reaction = ct.FalloffReaction( + reactants=ct_reactants, + products=ct_products, + tbody=ct_collider, + rate=ct.LindemannRate() + ) + else: + ct_reaction = ct.FalloffReaction( + reactants=ct_reactants, + products=ct_products, + rate=ct.LindemannRate() + ) + + elif isinstance(self.kinetics, SurfaceArrhenius): + ct_reaction = ct.InterfaceReaction( + reactants=ct_reactants, + products=ct_products, + rate=ct.InterfaceArrheniusRate() + ) - return ct_reaction + elif isinstance(self.kinetics, StickingCoefficient): + ct_reaction = ct.Reaction( + reactants=ct_reactants, + products=ct_products, + rate=ct.StickingArrheniusRate() + ) else: - raise Exception('Cantera reaction cannot be created because there was no kinetics.') + raise NotImplementedError(f"Unable to set cantera kinetics for {self.kinetics}") + + # Set reversibility, duplicate, and ID attributes + if isinstance(ct_reaction, list): + for rxn in ct_reaction: + rxn.reversible = self.reversible + # Set the duplicate flag to true since this reaction comes from multiarrhenius or multipdeparrhenius + rxn.duplicate = True + # Set the ID flag to the original rmg index + rxn.ID = str(self.index) + else: + ct_reaction.reversible = self.reversible + ct_reaction.duplicate = self.duplicate + ct_reaction.ID = str(self.index) + + # Now we set the kinetics. + self.kinetics.set_cantera_kinetics(ct_reaction, species_list) + + return ct_reaction + + def get_url(self): """