From 08724c9ac3e338853d111333552041081e3a0b30 Mon Sep 17 00:00:00 2001 From: Ray Speth Date: Sat, 6 Jul 2013 21:44:10 +0000 Subject: [PATCH] [ck2cti] Translate "chemically activated" reactions These are reactions with an additional set of Arrhenius parameters labeled 'HIGH' in the input file. They may also have a set of falloff function parameters. --- interfaces/python/ck2cti.py | 94 +++++++++++++++++++-- test/data/chemically-activated-reaction.inp | 38 +++++++++ test/python/testConvert.py | 9 ++ 3 files changed, 132 insertions(+), 9 deletions(-) create mode 100644 test/data/chemically-activated-reaction.inp diff --git a/interfaces/python/ck2cti.py b/interfaces/python/ck2cti.py index d2e3629fd..2e0e33a6b 100755 --- a/interfaces/python/ck2cti.py +++ b/interfaces/python/ck2cti.py @@ -627,6 +627,62 @@ class Falloff(ThirdBody): return '\n'.join(lines) +class ChemicallyActivated(ThirdBody): + """ + A kinetic model of a phenomenological rate coefficient k(T, P) using the + expression + + .. math:: k(T,P) = k_0(T) \\left[ \\frac{1}{1 + P_\\mathrm{r}} \\right] F + + where + + .. math:: + + P_\\mathrm{r} &= \\frac{k_0(T)}{k_\\infty(T)} [\\ce{M}] + + k_0(T) &= A_0 T^{n_0} \\exp \\left( - \\frac{E_0}{RT} \\right) + + k_\\infty(T) &= A_\\infty T^{n_\\infty} \\exp \\left( - \\frac{E_\\infty}{RT} \\right) + + and :math:`[\\ce{M}] \\approx P/RT` is the concentration of the bath gas. + The Arrhenius expressions :math:`k_0(T)` and :math:`k_\\infty(T)` + represent the low-pressure and high-pressure limit kinetics, respectively. + The former is necessarily one reaction order higher than the latter. The + allowable parameterizations for the function *F* are the same as for the + `Falloff` class. A collision efficiency can be used to further correct the + value of :math:`k(T,P)`. + + The attributes are: + + =============== ======================= ==================================== + Attribute Type Description + =============== ======================= ==================================== + `arrheniusLow` :class:`Arrhenius` The Arrhenius kinetics at the low-pressure limit + `arrheniusHigh` :class:`Arrhenius` The Arrhenius kinetics at the high-pressure limit + `efficiencies` ``dict`` A mapping of species to collider efficiencies + `F` Falloff function parameterization + =============== ======================= ==================================== + """ + def __init__(self, arrheniusLow=None, F=None, **kwargs): + ThirdBody.__init__(self, **kwargs) + self.arrheniusLow = arrheniusLow + self.F = F + + def to_cti(self, reactantstr, arrow, productstr, indent=0): + rxnstr = reactantstr + arrow + productstr + prefix = ' '*(indent + 30) + lines = ['chemically_activated_reaction({0!r},'.format(rxnstr)] + lines.append(prefix + 'kLow={0},'.format(self.arrheniusLow.rateStr())) + lines.append(prefix + 'kHigh={0},'.format(self.arrheniusHigh.rateStr())) + if self.efficiencies: + lines.append(prefix + 'efficiencies={0!r},'.format(self.efficiencyString())) + if self.F: + lines.append(prefix + 'falloff={0},'.format(self.F.to_cti())) + + lines[-1] = lines[-1][:-1] + ')' + return '\n'.join(lines) + + class Troe(object): """ For the Troe model the parameter :math:`F` is computed via @@ -1052,10 +1108,9 @@ class Parser(object): klow_units = self.getRateConstantUnits(length_dim + 3, 'cm', quantity_dim + 1, quantity_units) - # The rest of the first line contains the high-P limit Arrhenius parameters (if available) - #tokens = lines[0][52:].split() + # The rest of the first line contains Arrhenius parameters tokens = lines[0].split()[1:] - arrheniusHigh = Arrhenius( + arrhenius = Arrhenius( A=(A,kunits), n=n, Ea=(Ea, energy_units), @@ -1064,11 +1119,12 @@ class Parser(object): ) if len(lines) == 1: - # If there's only one line then we know to use the high-P limit kinetics as-is - reaction.kinetics = arrheniusHigh + # If there's only one line then we know to use the kinetics as-is + reaction.kinetics = arrhenius else: # There's more kinetics information to be read arrheniusLow = None + arrheniusHigh = None falloff = None chebyshev = None pdepArrhenius = None @@ -1084,7 +1140,7 @@ class Parser(object): reaction.duplicate = True elif 'low' in line.lower(): - # Low-pressure-limit Arrhenius parameters + # Low-pressure-limit Arrhenius parameters for "falloff" reaction tokens = tokens[1].split() arrheniusLow = Arrhenius( A=(float(tokens[0].strip()),klow_units), @@ -1094,6 +1150,20 @@ class Parser(object): parser=self ) + elif 'high' in line.lower(): + # High-pressure-limit Arrhenius parameters for "chemically + # activated" reaction + tokens = tokens[1].split() + arrheniusHigh = Arrhenius( + A=(float(tokens[0].strip()),kunits), + n=float(tokens[1].strip()), + Ea=(float(tokens[2].strip()),energy_units), + T0=(1,"K"), + parser=self + ) + # Need to fix units on the base reaction: + arrhenius.A = (arrhenius.A[0], klow_units) + elif 'rev' in line.lower(): reaction.reversible = False @@ -1215,17 +1285,23 @@ class Parser(object): parser=self ) elif arrheniusLow is not None: - reaction.kinetics = Falloff(arrheniusHigh=arrheniusHigh, + reaction.kinetics = Falloff(arrheniusHigh=arrhenius, arrheniusLow=arrheniusLow, F=falloff, parser=self, efficiencies=efficiencies) + elif arrheniusHigh is not None: + reaction.kinetics = ChemicallyActivated(arrheniusHigh=arrheniusHigh, + arrheniusLow=arrhenius, + F=falloff, + parser=self, + efficiencies=efficiencies) elif thirdBody: - reaction.kinetics = ThirdBody(arrheniusHigh=arrheniusHigh, + reaction.kinetics = ThirdBody(arrheniusHigh=arrhenius, parser=self, efficiencies=efficiencies) else: - reaction.kinetics = arrheniusHigh + reaction.kinetics = arrhenius return reaction, revReaction diff --git a/test/data/chemically-activated-reaction.inp b/test/data/chemically-activated-reaction.inp new file mode 100644 index 000000000..c45ab18a3 --- /dev/null +++ b/test/data/chemically-activated-reaction.inp @@ -0,0 +1,38 @@ +elements +c h n o ar +!to modify +end + +species +ch3 oh ch2o h2 n2 +end + +thermo all + 300.000 1000.000 5000.000 +ch3 iu0702c 1 h 3 0 0 g 200.000 6000.00 1000.00 1 + 0.29781206e+01 0.57978520e-02-0.19755800e-05 0.30729790e-09-0.17917416e-13 2 + 0.16509513e+05 0.47224799e+01 0.36571797e+01 0.21265979e-02 0.54583883e-05 3 +-0.66181003e-08 0.24657074e-11 0.16422716e+05 0.16735354e+01 0.17643935e+05 4 +oh iu3/03o 1 h 1 0 0g 200.000 6000.00 1000.00 1 + 2.83853033e+00 1.10741289e-03-2.94000209e-07 4.20698729e-11-2.42289890e-15 2 + 3.69780808e+03 5.84494652e+00 3.99198424e+00-2.40106655e-03 4.61664033e-06 3 +-3.87916306e-09 1.36319502e-12 3.36889836e+03-1.03998477e-01 4.48613328e+03 4 +ch2o g 8/88h 2 c 1 o 1 0 g 200.000 6000.00 1000.00 1 + 3.16952665e+00 6.19320560e-03-2.25056366e-06 3.65975660e-10-2.20149458e-14 2 +-1.44922756e+04 6.04207898e+00 4.79372312e+00-9.90833322e-03 3.73219990e-05 3 +-3.79285237e-08 1.31772641e-11-1.43227879e+04 6.02798058e-01-1.30729291e+04 4 +h2 tpis78h 2 0 0 0g 200.000 6000.00 1000.00 1 + 2.93286575e+00 8.26608026e-04-1.46402364e-07 1.54100414e-11-6.88804800e-16 2 +-8.13065581e+02-1.02432865e+00 2.34433112e+00 7.98052075e-03-1.94781510e-05 3 + 2.01572094e-08-7.37611761e-12-9.17935173e+02 6.83010238e-01 0.00000000e+00 4 +n2 g 8/02n 2 0 0 0g 200.000 6000.00 1000.00 1 + 2.95257637e+00 1.39690040e-03-4.92631603e-07 7.86010195e-11-4.60755204e-15 2 +-9.23948688e+02 5.87188762e+00 3.53100528e+00-1.23660988e-04-5.02999433e-07 3 + 2.43530612e-09-1.40881235e-12-1.04697628e+03 2.96747038e+00 0.00000000e+00 4 +end + +reactions +ch3+oh(+m)<=>ch2o+h2(+m) 282320.078 1.46878 -3270.56495 +high / 5.88E-14 6.721 -3022.227 / +troe / 1.671E+00 434.782 2934.21 3.919E+03 / +end diff --git a/test/python/testConvert.py b/test/python/testConvert.py index 4039f3f62..28a50bdf4 100644 --- a/test/python/testConvert.py +++ b/test/python/testConvert.py @@ -118,6 +118,15 @@ class chemkinConverterTest(utilities.CanteraTest): 'sri-falloff.cti') self.checkKinetics(ref, gas, [300, 800, 1450, 2800], [5e3, 1e5, 2e6]) + def test_chemically_activated(self): + name = 'chemically-activated-reaction' + convertMech('../data/{0}.inp'.format(name), + outName='{0}.cti'.format(name), quiet=True) + + ref, gas = self.checkConversion('../data/{0}.xml'.format(name), + '{0}.cti'.format(name)) + self.checkKinetics(ref, gas, [300, 800, 1450, 2800], [5e3, 1e5, 2e6, 1e7]) + def test_explicit_third_bodies(self): convertMech('../data/explicit-third-bodies.inp', thermoFile='../data/dummy-thermo.dat',