Skip to content

Commit

Permalink
Merge pull request #2700 for better to_cantera methods.
Browse files Browse the repository at this point in the history 
Add and improve to_cantera methods for conversion to in-memory Cantera objects.
  • Loading branch information
@rwest
rwest authored Dec 3, 2024
2 parents ac71689 + ffa0ca1 commit 2a1f910
Show file tree
Hide file tree
Showing 11 changed files with 416 additions and 153 deletions.
14 changes: 6 additions & 8 deletions rmgpy/chemkin.pyx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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)


Expand Down
43 changes: 40 additions & 3 deletions rmgpy/kinetics/arrhenius.pyx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -722,12 +722,49 @@ 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,
}

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()

################################################################################

Expand Down
10 changes: 0 additions & 10 deletions rmgpy/kinetics/falloff.pyx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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()

Expand Down Expand Up @@ -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):
Expand All @@ -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)




87 changes: 83 additions & 4 deletions rmgpy/kinetics/surface.pyx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -285,6 +285,34 @@ cdef class StickingCoefficient(KineticsModel):

return string

def to_cantera_kinetics(self):
"""
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
E = self._Ea.value_si * 1000 # convert from J/mol to J/kmol

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 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):
"""
Expand Down Expand Up @@ -352,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):
Expand Down Expand Up @@ -395,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)
Expand Down Expand Up @@ -524,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):
Expand Down Expand Up @@ -589,6 +618,56 @@ cdef class SurfaceArrhenius(Arrhenius):
)


def to_cantera_kinetics(self):
"""
Converts the RMG SurfaceArrhenius object to a cantera InterfaceArrheniusRate
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

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

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):
"""
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()

################################################################################

cdef class SurfaceArrheniusBEP(ArrheniusEP):
Expand Down
Loading

0 comments on commit 2a1f910

Please sign in to comment.