22Creating Input Files
33********************
44
5- There are three parts to a MEASURE input file, giving the species, path
6- reaction, and algorithm parameters. The species section must come before the
7- reaction section. Before discussing each of these sections, a brief word on the
8- input file syntax will be provided.
5+ There are four parts to a MEASURE input file, giving the species, important
6+ unimolecular and bimolecular configurations, path reactions, and algorithm
7+ parameters. The species section must come before the reaction section. Before
8+ discussing each of these sections, a brief word on the general input file
9+ syntax will be given.
910
1011Syntax
1112======
@@ -19,20 +20,20 @@ specified as text strings, numbers, or objects. Text strings must be wrapped in
1920either single or double quotes. There are two ways to specify numeric
2021quantities:
2122
22- * Dimensionless quantities or values already in SI units can be specified as
23- a simple literal, such as ``0.99 ``.
23+ * Dimensionless quantities can be specified as a simple literal, such as
24+ ``0.99 ``.
2425
25- * Quantities not in SI units are specified using a tuple of the form
26+ * Dimensional quantities are specified using a tuple of the form
2627 ``(value, units) ``, where ``value `` is either a single number or a list of
27- numbers, and ``units `` is a string, e.g. ``(-10.0,'kcal/mol') `` and
28+ numbers, and ``units `` is a string, e.g. ``(-10.0,'kcal/mol') `` or
2829 ``([1.0, 2.0, 3.0],"m") ``.
2930
3031.. note ::
3132
3233 MEASURE uses the ``quantities `` package to convert your numeric parameters
3334 into SI units, and therefore inherits all of the idiosyncracies from that
3435 package. In particular, the ``quantities `` package does *not *
35- follow the SI convention that all units after the circumflex are in the
36+ follow the SI convention that all units after the forward slash are in the
3637 denominator. For example, ``J/mol*K `` would be interpreted as ``(J/mol)*K ``
3738 rather than ``J/(mol*K) ``. Thus we recommend using parentheses where
3839 necessary to make your intentions explicit.
@@ -71,49 +72,104 @@ bimolecular reactant channels. When specifying the ``states`` parameter, use a
7172=========================== ====================================================
7273Parameter Description
7374=========================== ====================================================
74- ``rotationalConstants `` A list of the external rotational constants
75- ``symmetry `` The external symmetry number
76- ``frequencies `` A list of the vibrational frequencies
75+ ``rotations `` Parameters describing the external rotational motion, as a `` RigidRotor() `` object
76+ ``vibrations `` Parameters describing the internal vibrational motion, as a `` HarmonicOscillator() `` object
77+ ``torsions `` Parameters describing the internal torsional motion, as a list of `` HinderedRotor() `` objects
7778``frequencyScaleFactor `` The frequency scale factor to use (1.0 if not specified)
78- ``hinderedRotors `` A list of simple Pitzer hindered rotors, specified as 3-tuples (``frequency ``, ``barrier ``, ``symmetry ``)
7979``spinMultiplicity `` The ground-state spin multiplicity (degeneracy)
8080=========================== ====================================================
8181
82+ The ``RigidRotor() ``, ``HarmonicOscillator() ``, and ``HinderedRotor() ``
83+ constructors match the corresponding classes in the :mod: `rmgpy.statmech `
84+ module. The parameters for each are also summarized below:
85+
86+ =========================== ====================================================
87+ Parameter Description
88+ =========================== ====================================================
89+ ``RigidRotor() ``
90+ --------------------------- ----------------------------------------------------
91+ `linear ` ``True `` if the associated molecule is linear, ``False `` if nonlinear
92+ `inertia ` A list of the moment(s) of inertia of the molecule (1 if linear, 3 if nonlinear)
93+ `symmetry ` The total external rotational symmetry number
94+ --------------------------- ----------------------------------------------------
95+ ``HarmonicOscillator() ``
96+ --------------------------- ----------------------------------------------------
97+ `frequencies ` The set of vibrational frequencies
98+ --------------------------- ----------------------------------------------------
99+ ``HinderedRotor() ``
100+ --------------------------- ----------------------------------------------------
101+ `inertia ` The reduced moment of inertia of the hindered rotor
102+ `symmetry ` The symmetry number for the hindered rotation
103+ `barrier ` The barrier height of the cosine potential
104+ `fourier ` The :math: `2 \times C` array of Fourier coefficients for the Fourier series potential
105+ =========================== ====================================================
106+
107+ For each ``HinderedRotor() ``, you need only specify one of the barrier height
108+ or Fourier series coefficients.
109+
82110If ``states `` is specified and ``thermo `` is not, then the thermodynamic
83111parameters will be automatically computed. This is recommended unless you have
84112thermodynamic data that you believe to be more accurate than the molecular
85113degrees of freedom data. You can use any of the thermodynamics models in the
86- ``chempy .thermo `` module (from the ChemPy package) ; see that package for more
87- information on the available models and their syntax.
114+ ``rmgpy .thermo `` module; see that package for more information on the available
115+ models and their syntax.
88116
89117The following is an example of a typical species item, based on the acetylperoxy
90118radical :math: `\ce {CH3 C(=O)OO.}`::
91119
92120 species(
93- label=" acetylperoxy" ,
94- SMILES=" CC(=O)O[O]" ,
95- E0=(-144.766,"kJ /mol" ),
121+ label=' acetylperoxy' ,
122+ SMILES=' CC(=O)O[O]' ,
123+ E0=(-34.6,'kcal /mol' ),
96124 states=States(
97- rotationalConstants=([54.2978, 104.836, 156.049], "amu/angstrom^2"),
98- symmetry=1,
99- frequencies=([321.607, 503.468, 539.885, 547.148, 731.506, 979.187, 1043.98, 1126.42, 1188.62, 1399.43, 1458.2, 1463.42, 1881.7, 3055.28, 3115.45, 3155.14], "cm^-1"),
100- frequencyScaleFactor=1.0,
101- hinderedRotors=[
102- ((7.38359,"amu*angstrom^2"), (25.5921,"kJ/mol"), 1),
103- ((2.94725,"amu*angstrom^2"), (5.11105,"kJ/mol"), 3),
125+ rotations=RigidRotor(
126+ linear=False,
127+ inertia=([54.2978, 104.8364, 156.0495],"amu*angstrom^2"),
128+ symmetry=1,
129+ ),
130+ vibrations=HarmonicOscillator(
131+ frequencies=([321.607, 503.468, 539.885, 547.148, 731.506, 979.187, 1043.981, 1126.416, 1188.619, 1399.432, 1458.200, 1463.423, 1881.701, 3055.285, 3115.447, 3155.144], 'cm^-1'),
132+ ),
133+ torsions=[
134+ HinderedRotor(inertia=(7.38359,"amu*angstrom^2"), barrier=(6.11665,"kcal/mol"), symmetry=1),
135+ HinderedRotor(inertia=(2.94725,"amu*angstrom^2"), barrier=(1.22157,"kcal/mol"), symmetry=3),
104136 ],
137+ frequencyScaleFactor=0.99,
105138 spinMultiplicity=2,
106139 ),
107- thermo=ThermoGAModel(
108- Tdata=([303.231, 385.553, 467.875, 550.197, 632.519, 714.841, 797.163, 879.485, 961.807, 1044.13, 1126.45, 1208.77, 1291.09, 1373.42, 1455.74, 1538.06, 1620.38, 1702.7, 1785.03, 1867.35], "K"),
109- Cpdata=([80.917, 92.4545, 102.856, 111.821, 119.406, 125.81, 131.243, 135.887, 139.882, 143.341, 146.352, 148.984, 151.294, 153.33, 155.129, 156.725, 158.145, 159.413, 160.548, 161.567], "J/(mol*K)"),
110- H298=(-127.614, "kJ/mol"),
111- S298=(314.403, "J/(mol*K)"),
112- ),
113- lennardJones=LennardJones(sigma=(5.09e-10,"m"), epsilon=(6.53048e-21,"J")),
114- molecularWeight=(75.0434,"g/mol"),
140+ lennardJones=LennardJones(sigma=(5.09,'angstrom'), epsilon=(473,'K')),
115141 )
116142
143+ Molecular Configurations
144+ ========================
145+
146+ MEASURE is largely able to determine the molecular configurations that define
147+ the potential energy surface for your reaction network simply by inspecting the
148+ path reactions. However, you must indicate which unimolecular and bimolecular
149+ configurations you wish to include in the master equation formulation; all
150+ others will be treated as irreversible sinks.
151+
152+ * For a unimolecular configuration, use the ``isomer() `` method, passing as the
153+ only parameter a string containing the label of the species to treat as a
154+ unimolecular isomer.
155+
156+ * For a bimolecular configuration, use the ``reactants() `` method, passing as
157+ the two parameters a pair of strings containing the labels of the species to
158+ treat as a bimolecular reactant channel. (A reactant channel is allowed to be
159+ both a source and a sink, i.e. both association and dissociation pathways are
160+ kept.)
161+
162+ For example, the following input specifies acetylperoxy as a unimolecular
163+ isomer and acetyl + oxygen as a bimolecular reactant channel.
164+
165+ isomer('acetylperoxy')
166+
167+ reactants('acetyl', 'oxygen')
168+
169+ You do not need to specify the product channels (infinite sinks) in this
170+ manner, as any configuration not marked as an isomer or reactant channel will
171+ be treated as a product channel.
172+
117173Path Reaction Parameters
118174========================
119175
@@ -127,17 +183,9 @@ Parameter Required? Description
127183``reactants `` All reactions A list of strings indicating the labels of the reactant species
128184``products `` All reactions A list of strings indicating the labels of the product species
129185``transitionState `` All reactions Information about the transition state, using a ``TransitionState() `` block; see below
130- ``reversible `` ``True `` (default) if the reaction is reversible, ``False `` if not
131186``kinetics `` The high pressure-limit kinetics for the reaction
132187====================== ==================== ====================================
133188
134- The direction and reversibility of the reactions determine how bimolecular
135- configurations are divided into reactants. If the bimolecular configuration is
136- listed as the reactants or ``reversible `` is ``True ``, the configuration is
137- set as a reactant channel. If the bimolecular configuration is listed as the
138- products and ``reversible `` is ``False ``, the configuration is set as a product
139- channel.
140-
141189The type of information specified with each path reaction determines how the
142190microcanonical rate coefficient is computed:
143191
@@ -164,7 +212,6 @@ The following is an example of a typical reaction item, based on the reaction
164212 reaction(
165213 reactants=['acetylperoxy'],
166214 products=['ketene', 'hydroperoxyl'],
167- reversible=False,
168215 kinetics=Arrhenius(
169216 A=(2.62e9,'s^-1'),
170217 n=1.24,
@@ -173,11 +220,15 @@ The following is an example of a typical reaction item, based on the reaction
173220 transitionState=TransitionState(
174221 E0=(0.6,'kcal/mol'),
175222 states=States(
176- rotationalConstants=([55.4256, 136.1886, 188.2442], 'amu*angstrom^2'),
177- symmetry=1,
178- frequencies=([59.306, 205.421, 354.483, 468.861, 482.875, 545.574, 657.825, 891.898, 1023.947, 1085.617, 1257.494, 1316.937, 1378.552, 1688.566, 2175.346, 3079.822, 3154.325], 'cm^-1'),
223+ rotations=RigidRotor(
224+ linear=False,
225+ inertia=([55.4256, 136.1886, 188.2442],"amu*angstrom^2"),
226+ symmetry=1,
227+ ),
228+ vibrations=HarmonicOscillator(
229+ frequencies=([59.306, 205.421, 354.483, 468.861, 482.875, 545.574, 657.825, 891.898, 1023.947, 1085.617, 1257.494, 1316.937, 1378.552, 1688.566, 2175.346, 3079.822, 3154.325], 'cm^-1'),
230+ ),
179231 frequencyScaleFactor=0.99,
180- hinderedRotors=[],
181232 spinMultiplicity=2,
182233 ),
183234 frequency=(-1048.9950,'cm^-1'),
@@ -193,7 +244,7 @@ Collision Model
193244The collision model to use when constructing the master equation is specified
194245using a ``collisionModel() `` block, where you must specify both the ``type ``
195246of model, the value of any required ``parameters `` as a list of quantities,
196- and a list of tuples specifying the labels of the species in the bath gas and
247+ and a dictionary specifying the labels of the species in the bath gas and
197248their relative mole fractions. The collision models available are:
198249
199250* ``single exponential down `` - Specify ``alpha0 ``, ``T0 `` and ``n `` for the
@@ -267,11 +318,11 @@ Use a ``method()`` block to specify the approximate method to use when computing
267318:math: `k(T,P)` values from the full master equation. There are currently three
268319methods available: ``modified strong collision `` (MSC), ``reservoir state ``
269320(RS), and ``chemically-significant eigenvalues `` (CSE). The details of each of
270- these methods is provided in the :doc: `../theory/index `. In short: MSC is the
271- fastest but least accurate, RS usually provides a good balance of speed and
272- accuracy, while CSE is the most accurate but is slow and not robust. Currently
273- we recommend using MSC during initial explorations, then RS when more accurate
274- numbers are needed.
321+ these methods is provided in the :doc: `../../ theory/measure/ index `. In short:
322+ MSC is the fastest but least accurate, RS usually provides a good balance of
323+ speed and accuracy (except at very high temperatures) , while CSE is the most
324+ accurate but is slow and not robust. Currently we recommend using MSC during
325+ initial explorations, then RS when more accurate numbers are needed.
275326
276327An example ``method() `` block is given below::
277328
0 commit comments