Bimolecular Master Equations for a Single and Multiple Potential Wells

Mar 27, 2018 - The analytic solutions, that is, populations, are derived for the K-adiabatic and K-active bimolecular master equations, separately, fo...
0 downloads 8 Views 913KB Size
Article Cite This: J. Phys. Chem. A XXXX, XXX, XXX−XXX

pubs.acs.org/JPCA

Bimolecular Master Equations for a Single and Multiple Potential Wells with Analytic Solutions Nima Ghaderi* Department of Physics, Beckman Institute, and Noyes Laboratory of Chemical Physics, California Institute of Technology, 1200 East California Boulevard, Pasadena, California 91125, United States ABSTRACT: The analytic solutions, that is, populations, are derived for the K-adiabatic and Kactive bimolecular master equations, separately, for a single and multiple potential wells and reaction channels, where K is the component of the total angular momentum J along the axis of least moment of inertia of the recombination products at a given energy E. The analytic approach provides the functional dependence of the population of molecules on its K-active or K-adiabatic dissociation, association rate constants and the intermolecular energy transfer, where the approach may complement the usual numerical approaches for reactions of interest. Our previous work, Part I, considered the solutions for a single potential well, whereby an assumption utilized there is presently obviated in the derivation of the exact solutions and farther discussed. At the high-pressure limit, the K-adiabatic and K-active bimolecular master equations may each reduce, respectively, to the K-adiabatic and K-active bimolecular Rice-Ramsperger-KasselMarcus theory (high-pressure limit expressions) for bimolecular recombination rate constant, for a single potential well, and augmented by isomerization terms when multiple potential wells are present. In the low-pressure limit, the expression for population above the dissociation limit, associated with a single potential well, becomes equivalent to the usual presumed detailed balance between the association and dissociation rate constants, where the multiple well case is also considered. When the collision frequency of energy transfer, ZLJ, between the chemical intermediate and bath gas is sufficiently less than the dissociation rate constant kd(E′J′K′) for postcollision (E′J′K), then the solution for population, g(EJK)+, above the critical energy further simplifies such that depending on ZLJ, the dissociation and association rate constant kr(EJK), as g(EJK)+ = kr(EJK)A·BC/ [ZLJ+kd(EJK)], where A and BC are the reactants, for example, relevant for O3 formation from O + O2 + Ar up to ∼100 bar; otherwise, additional contributions from postcollision are present and especially relevant at high pressures. In the aforementioned regime ZLJ < kd(E′J′K) the physical connection of recombination rate constants, krec based on either utilizing population from the master equation approach or a collision based bimolecular RRKM theory is traced and elucidated analytically that the rate constants are equal. Recombination rate constants, krec, based on the population, are also given and considered for an adiabatic or active K. For example, for O3 formation in Ar bath gas, the K-adiabatic-based krec for O3 yields 4.0 × 10−34 cm6 molecule−2 s−1 at T = 300 K and 1 bar, in agreement with the experimental value, where the contribution from the population of metastable ozone is discussed and the adiabaticity of K highlighted. A facile numerical implementation of the formalism for g(EJK) and krec for O3 is noted. The application of the expressions to ozone recombination as a function of pressure and temperature is given elsewhere, beyond the selection considered here, for studying the physical features, including the contributions from the K and intermolecular energy transfer to the krec, and the puzzles reported from experiments for this reaction. thermal bimolecular sampling),14 and its survival probability for the important energy and angular momentum region of interest may be such that multiple collisions with a bath gas may take place, for example, by 1000 ps at E = 0.2kT (above the dissociation limit) and J = 15, the survival probability of Ο*3 is ∼50%,14 where then ∼10 collisions between Ο*3 and Ar may take place at T = 300 K and 1 bar, before either dissociating back to the reactants via its entrance or exchange channels or the intermediate become stabilized (e.g., Ο*3 +M → O3 + M, with a bath gas M).14 We treat this case using an analytic approach. The most general case includes the possibility of the complex isomerizing to another complex, a stabilization, or a dissociation either back to reactants or to one of several sets of bimolecular product channels.5 In this paper, we also consider this general

1. INTRODUCTION Bimolecular reactions are fundamental physical processes in chemical reactions1−3 and may occur in various ways,1−7 in which statistical theories have been developed for their treatments.1,2,8−14 In perhaps the simplest case in the course of a reaction, a simple barrier may separate a set of bimolecular reactants from a set of bimolecular products, for example, A + BC → AB + C.1 In this case, the rate coefficient is a function of temperature but not of pressure, since the intermediate complex may not form with a requisite lifetime to interact with a bath gas, and hence the bath gas would not assist in the stabilization of any intermediate or facilitate decay into the product channels. Another case, a common occurrence, is for one or more potential energy wells to lie along the reaction path, and so the rate constant can exhibit a pressure dependence.4,5 In the presence of a potential energy well, the intermediate complex formed from the reactant may form a long-lived complex (e.g., an average half-time of ∼50 ps for an excited ozone, Ο3*, from a © XXXX American Chemical Society

Received: September 16, 2017 Revised: December 15, 2017

A

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

equations are derived, using an algebraic method, for a single potential energy well and separately for a multiple potential energy wells, including isomerization events and dissociation into bimolecular product channels. The solutions, like before,30 may accommodate various intermolecular energy transfer models such as a single exponential,1 double exponential,1,36 Gaussian,1,37,38 stepladder,39−41 and near singularity,42 as the functional form of each is noted in the Appendix A of Part I30 and a selection described in the Appendix B of the present work. The case when the energy transfer function is explicitly factorizable into pre- and postcollision energies, for example, for a single exponential model, is also noted, and derivation of solution of the master equation is given in Appendix E. In the regime where the collision frequency that leads to an intermolecular energy transfer for the excited intermediate is sufficiently less than the dissociation rate constant of the intermediate, the solution describing the population above the dissociation limit further simplifies as discussed in the main text, yet still retains information on collision, that is, reduces to a modified equilibrium constant with a facile numerical implementation. In unimolecular and bimolecular reactions, the experimental observables amenable to measurement are the rate constants, quantum state distribution of products, evolution of population, and correlations of vector and scalar quantities involved in the reaction.1,15,23,24 The magnitude of population of molecules near the dissociation limit may be investigated experimentally;43−45 for example, with a pump−dump−probe experimental setup, the HFCO molecules have been excited to a dissociative level by stimulated emission pumping, and then the decay of the excited population to HF + CO was monitored by measuring the fluorescence intensity of HFCO molecules by the probe laser as a function of time delay between the dump and probe pulses and so measuring the temporal evolution of the population.44 Herein, we solve for the steady-state populations of a bimolecular reaction for below and above the dissociation limit separately. Thence, the predictions for the magnitude of these populations in both energy regions can be investigated by experiments. Expressions for the recombination rate constants based on the analytic-based populations are also given in the present work. We now pause to comment on the validity of the steady-state assumption15 in the study of chemical reactions and also note some examples of it in chemical reactions, since the analytic solutions of the bimolecular master equations are derived in the present work for this regime. The steady-state hypothesis,15 that is, d[intermediate]/dt = 0, may be applicable to reactions, for

case as may apply to reactions of interest and derive analytic solutions (populations) for the steady-state (SS) regime,15,16 where after an initial induction period the population has relaxed toward a steady-state,16 for a single and multiple potential wells K-adiabatic and K-active bimolecular master equations, where K is the component of the total angular momentum J along the axis of least moment of inertia of the recombination products. The validity of the steady-state regime is described in more detail later in this section. Some reactions of interest with multiple potential wells studied by master equations include OH + CH2O, OH + CH2CH2, OH + CH2NH, and H + SO2.17,18 The number of possible intermediate isomers increases rapidly with the total number of atoms, (e.g., a five-atom system like HCO + NO can have an energy diagram more complicated than two isomers).19 Reactions with both a barrier and multiple potential wells have been noted to be particularly interesting,5 as they may exhibit a multimodal lifetime distribution and possible manifestations of nonstatistical behavior.5 The analytic solution for the multiple potential wells presented in the present work permits the treatment of any given barrier for the determination of populations, as the number of states at the transition state (TS)20 and the density of states in each potential well may be counted (e.g., using a Beyer−Swinehart algorithm,21 with any anharmonic considerations)22,23 and then utilized in RiceRamsperger-Kassel-Marcus (RRKM) theory1,9,10,15,23,24 to determine the rate constants of dissociation and association processes when determining the populations expressed by the analytic-based formulas derived from the K-adiabatic and Kactive bimolecular master equations. Any of the nonstatistical effects, such as recrossings of the transition state and phase space inaccessibility,1,14,23 may be included as corrections to the number of states at TS and density of states, respectively, when present. The dissociation and association rate constants may also be determined by other methods of interest as selected by the investigator and be used in the analytic solutions. When K is a dynamically slow variable and shares energy less freely with the remaining relevant degrees of freedom, then it is referred to as K-adiabatic, in contrast to the K-active case, whereby the Coriolis coupling may mix the possible (2J + 1) K levels for the J quantum number.9,10,14,22,25,26 The ozone molecule is an example of such a system, where its K has been investigated and also shown to exhibit adiabaticity,14,27−29 for example, the role of K in the intermolecular energy transfer with a bath gas,27 K-diffusion in the metastable ozone,28,29 and its role in the recombination rate constant krec for Ο + Ο2 → Ο3.14 Expressions for various alternative combinations for the treatment of angular momentum in RRKM theory,1,9,10,15,23,24 including K-active and K-adiabatic, have been delineated by Hase26 and explored elsewhere,10 and they may be considered for a reaction of interest. The analytic-based K-adiabatic and Kactive populations, considered in this work, may elucidate the population’s interwoven dependencies of association rate constants, dissociation rate constants, collision rates, and intermolecular energy transfers with the bath gases for reactions of interest, at a given pressure and temperature. In a previous work, for brevity referred to as Part I,30 the Kadiabatic and K-active analytic solutions for a single potential well for the bimolecular master equations were acquired using properties of linear inhomogeneous integral equations31 as developed by Fredholm32 and Hilbert−Schmidt theory,33,34 upon utilizing an assumption described in footnote 35. In the present work, this assumption35 is omitted, and analytic solutions for the K-adiabatic and K-active bimolecular master

k1

k3

example, A + BC ⇄ ABC*, ABC* + M → ABC + M, in k2

particular, if the concentration of the intermediate ABC* is much less than the sum of the concentrations of the reactant and products at all times, where M is the bath gas. This will be true if k2 + k3[M] ≫ k1[BC]. For example, in considering ozone formation, where Ο*3 is the metastable intermediate, using the rate constants from the literature based on experiment46 has k2 + k3[Ar] ≈ 1 × 1010 sec−1 and k1[O2] ≈ 1 × 106 sec−1, and therefore the inequality is satisfied for concentrations of O2 and Ar used in experiment,47 and the SS hypothesis is rendered applicable. The time necessary for the steady state of the intermediate to be reached, that is, the elapsed induction period, can be estimated by integrating the equation d[ABC*]/dt = k1[A][BC] − k2[ABC*] − k3[ABC*][M] over time with the B

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

population is also analyzed. The high-pressure limiting forms of the K-adiabatic and K-active master equations and their reduction to the high-pressure limiting form of bimolecular RRKM theory9 for a single potential well, and supplanted by isomerization processes for the case of the multiple potential wells, are given in Section 4, with derivations provided in Appendix H. Expressions for recombination rate constants based on the analytic-based populations are given in Section 5. Select results regarding ozone formation are given in Section 6. Relevant features of the population are discussed in Section 7. The summary and concluding remarks are given in Section 8.

assumption that [A] is constant during the induction period. Combining the latter result with setting d[ABC*]/dt = 0 yields the condition that steady state is reached when (k2 + k3[M])t ≫ 1.48 For ozone formation, the latter condition is satisfied when t is ∼1 × 10−8 sec or greater. When the time scale of the experiment, texp, is greater than the time to reach the steady state of the intermediate, t, that is, texp > t ≫ 1/(k2 + k3[M]), then the SS approximation spans the requisite range of time of the experiment after the induction period; for example, in a radiopulse experiment for studying ozone formation, the time scale of the experiment is ∼1000 times greater than the noted t and so applicable over the near entire time range of ozone formation, that is, after the induction period, in the experiment.48 Methane oxidation at low temperatures provides another example, whereby the intermediate species, CH3O, H, O, and OH are produced in small concentrations relative to the other participating species, and so the SS approximation may be applied to the rate equation describing each intermediate in the reaction mechanism.49 A similar prescription of a steady-state analysis may be extended to reactions possessing multiple potential wells with isomerization processes. Time-dependent master equation simulations for multiple wells may inform on the time profile of the intermediates and final product, for example, studied on cyclopropene system and propargyl production, achieving steady populations by 0.1 ms.50 In contrast, shock-tube experiments represent a case where the reaction is fast, and the time resolution of the observation in experiment may also be short, and in the latter regime the contribution of the transient to reaction rate may not be negligible compared to the steady-state component,51,52 with additional discussion, including the relative length of time, noted in ref 51 and the references noted there. The analytic results on population for a single potential well and recombination rate constants in the steady state may then be utilized, for example, to study O3 formation, that is, O + O2 + Ar → O3 + Ar, as a function of temperature and pressure,41,53−59 to understand the key issues, including the anomalous trends in the experiments,53 for example, “swoop-up” or the unexpected increase of recombination rate constants as pressure increases toward 1000 bar at temperature T = 300 K, and the rate constants at low temperature (T = 163 K) and low pressure (1 bar), along with the role of the energy transfer (ET) mechanism, that is, a Lindemann form,1,15 with the reaction steps noted in the next section, where the reaction steps are for general reaction. The application of the expressions and detailed results on the formation of O3 as a function of temperature and pressure will be presented elsewhere, while select numerical results are shared for elucidation on population and its contribution to recombination of O3 in the present work and thence also highlighting the utility of the expressions. The paper is organized as follows. In Section 2, the expressions for the K-adiabatic and K-active continuum bimolecular master equations with a single and multiple potential wells are reviewed30 with additional expounding, and the expressions also serve to be referred upon when deriving their analytic solutions. In Section 3, the analytic solutions (i.e., populations) for the K-adiabatic and K-active bimolecular master equations for a single and multiple potential wells are derived, and proofs for each are given in Appendices C and D, respectively, and augment the results from Part I.30,35 The case when the energy transfer function is explicitly factorizable into pre- and postcollision energies is also noted, and derivation of solution is given in Appendix E. The low-pressure limit of

2. THE BIMOLECULAR MASTER EQUATIONS 2.1. Single Chemical Intermediate. For the convenience of the reader, we recap the origin and the steps leading to the Kadiabatic and K-active bimolecular master equations from Part I30 for a single potential well in this section and for the multiple potential wells in the next, Section 2.2, to which also assists in referring to the requisite equations in this section, purposeful, when in the course of deriving the expressions for their solutions in Section 3. Additional expounding, for example, microscopic reversibility for activating and deactivating collisions with the bath gas (Appendix A) and collision-related rotational states is noted. The treatment of the bimolecular master equations and connection to previous works,5,41 appearing later in refs 60 and 61, is also noted. The (EJK) resolved kinetic scheme for a bimolecular reaction and energy transfer steps, as a Lindemann type,1,15 may be described by A + BC → ABC(EJK )

k r(EJK )

ABC(EJK ) + M → ABC(E′J ′K ′) + M

(1)

Z(EJK , E′J ′K ′) (2)

ABC(E′J ′K ′) + M → ABC(EJK ) + M

Z(E′J ′K ′, EJK ) (3)

ABC(E′J ′K ′) → A + BC

kden(E′J ′K ′)

(4)

ABC(E′J ′K ′) → AB + C

kdex(E′J ′K ′)

(5)

where A, B, and C are atomic or molecular species; for example, for simplicity BC may be a diatomic. The en and ex denote the entrance and the exchange channels, respectively. The Z(EJK,E′J′K′) in reaction step 2 is defined as the number of collisions per unit time, per unit energy, with energy transfer between vibrationally and rotationally excited intermediate and a bath gas. The Z(E′J′K′,EJK) in reaction step 3 has a similar definition. The frequency of activating collision is related to deactivation by microscopic reversibility (Appendix A). In a usual procession, a Z(EJK,E′J′K′) may comprise the normalized probability of energy transferred with the bath gas, multiplied by the total collision Lennard-Jones frequency, ZLJ = ωHSΩ, where ωHS is the hard-sphere collision frequency, and Ω, unitless, is a realistic molecule−bath interaction determined from a collision integral, for example, Lennard-Jones, which corrects the ωHS.1 Various models for intermolecular energy transfer1 for the reaction steps 2 and 3 are noted in the Appendix B. Reaction steps 1, 4, and 5 may take place only when E > E*(JK), where E*(JK) is the critical energy for reaction. Let g(EJK) be the population density of molecular products ABC (EJK) formed from a process A + BC→ABC, for a given E, J, and K.60 The K is the component of the total angular momentum J along the body-fixed axis with the least moment of C

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A inertia of the recombination products. The K in g(EJK) is assumed adiabatic when the Coriolis coupling does not mix the (2J + 1) K levels for a given J and is treated as such when resolved in the presentation throughout the paper, where otherwise the K is considered active in g(EJ) and so then absent in the presentation, as having been averaged over. For an ABC molecule that initially forms from a bimolecular reaction, for example, reaction step 1, its initial energy is, of course, greater than E*(JK), prior to any collision with the bath gas. The equation for the distribution function g(EJK) satisfies the bimolecular master equation30

EJK , Q JK =

JK

and gJ = 2J + 1, for a statistical distribution of J and K for any given E of ABC. However, prior to collisions with the bath gas, the rotational angular momentum of BC and the orbital angular momentum do not always add statistically when forming ABC,64,65 and thus any nonstatistical contribution to the rotational state of ABC that may remain after collisions would then serve to correct P(JK) in eq 10 and P(J′K′) in eq 11. We write, similarly: Z(EJK , E′J ′K ′) = P(J ′K ′)Z(E , E′)

dg (EJK , t ) = k r(EJK )A(t ) ·BC(t ) − g (EJK , t ) dt J′K′

+

Z(EJK , E′J ′K ′)dE′

Eo

∑∫



Z(E′J ′K ′, EJK )g (E′J ′K ′, t )dE′

Eo

J′K′

− kd(EJK )g (EJK , t )

(6)

for E > E*(JK), where the lower limit in the integrals for E′ corresponds to the minimum energy of the intermediate, where it may be set to Eo = −∞. The A(t)·BC(t) denotes the concentration of the reactants A and BC as a function of time.61 The association rate constant41 is given by k r(EJK ) = N *(EJK )exp( − E /kT )/hQ r

0 = k r(EJK )A·BC − g (EJK ) ∑ +

J′K′

− g (EJK , t ) ∑ J′K′

(8)

(12)

then eq 12 further simplifies to k r(EJK )A·BC = g (EJK )

∑∫

∫E





Z(E , E′)dE′ − P(JK )

o

Z(E′, E)g (E′J ′K ′)dE′

Eo

J′K′

+ kd(EJK )g (EJK )

(13)

A K-active g(EJ) can be derived following a similar procedure for the K-adiabatic case, except now the K degree of freedom is not explicitly resolved, or rather assumed to have appropriately been averaged over to yield the K-active master equation, given by1,30



Z(EJK , E′J ′K ′)dE′ (9)

whereupon the dissociation and association reaction channels are now closed, and so each is absent in eq 9. A recurrence relation for evolving g(EJK,t) as a function of time, for both energy regions E > E*(JK) and E ≤ E*(JK), can be considered.30 We next introduce a statistical assumption for the rotational state of ABC after a collision1,63 Z(E′J ′K ′, EJK ) = P(JK )Z(E′, E)

Z(E′, E)P(JK )g (E′J ′K ′)dE′ − kd(EJK )

Eo

J′K′

Z(E′J ′K ′, EJK )g (E′J ′K ′, t )dE′

o

Z(E , E′)P(J ′K ′)dE′

o

∑ P(J′K ′) = 1

Eo

∫E





The magnitude of the concentration of A·BC in eq 12 is assumed to be sufficiently in excess relative to that of the intermediate in the SS regime and so also appears as timeindependent in the forthcoming steps leading to the solution of the master equation. Next, factorizing P in the second and third term on the right-hand side in eq 12 and noting the normalization condition



∑∫

J′K′

∫E

g (EJK )

(7)

for the entrance and exchange channels en and ex, respectively, with units of (time)−1 per (EJK). The N*(EJK) is the number of states at the TS,20 and ρ(EJK) is the density of states of the metastable product. Mixed models of K,26 for example, K-active and K-adiabatic for N*(EJ) and ρ(EJK), respectively, may also be considered for the determination of ken,ex and kr and d subsequently in the analytic expressions for the populations in Section 3, as may be suited for a chemical reaction of interest. Alternatively, ken,ex and kr may be determined by other methods d of interest. For a process where E ≤ E*(JK), we have instead dg (EJK , t ) = dt

∑∫ J′K′

with partition functions Qr for the reactant pairs, with an overall unit of (concentration)−1 (time)−1 per (EJK). The RRKM theory1,9,10 dissociation rate constant is given by62 kden,ex(EJK ) = N *(EJK )/hρ(EJK )

(11)

In the present work, P(JK) and P(J′K′) are assumed normalized and can be represented by other distributions, where it would not affect the progression of the forthcoming derivations besides assuming an equilibrium distribution, and its application invites a proper selection for a given reaction of interest. For the steady-state approximation,15,16,42 dg(EJK,t)/dt = 0, where its validity was discussed in Section 1, then eq 6 becomes



∑∫

∑ gJ × exp( − EJK /kT )



∫E Z(E , E′)dE′ − P(J) ∞ ∑ ∫ Z(E′, E)g(E′J′)dE′ + kd(EJ )g(EJ ) E

k r(EJ )A·BC = g (EJ )

o

J′

o

(14)

which is the analogue of eq 13. The analytic solutions g(EJK) of the bimolecular master equation, eq 13, for the energy regions of E > E*(JK) and E ≤ E*(JK), and for its K-active counterpart, eq 14, are derived in Section 3. The recurrence relations for finding the energy

(10)

where P(JK) is the microcanonical equilibrium distribution of J and K at a given E, that is, gJ × exp(−EJK/kT)/QJK with rotational energy D

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A dependence of g(EJK) in eq 13 and for the K-active g(EJ) in eq 14 can also be determined.30 2.2. Multiple Chemical Intermediates. In this section, we consider chemical species that can be identified with local minima, potential wells, on the potential energy hypersurface. For a physical situation with multiple intermediates and reaction channels, the K-adiabatic master equation or its K-active counterpart can be written for each potential well, and the equations then can be coupled via the chemical reaction terms. Each reaction channel is either associated with another well or with fragmentation products. A master equation describing the time evolution of a population for N distribution functions gi(EJK,t), i = 1, ···, N, for N chemical configurations of A·BC is given by30 dgi(EJK , t ) dt

where all the terms in eq 16 are now the K-active counterparts of those defined earlier for eq 15. The E-resolved version of eq 16 is found in refs 1, 5, and 66 and in those cited therein. In the usual treatment, the term kri appearing in the master equation is rewritten as kdi and the equilibrium constant Keqi for A + BC ⇆ i using detailed balance condition.1 In the present treatment, the kri may be retained and determined by using eq 7 or other methods of interest, for each chemical intermediate. 61 Recurrence relations for treating the time evolution of the Kactive or K-adiabatic populations can also be considered.30 For the steady-state approximation,15,16,42 dgi(EJK,t)/dt = 0, eq 15 becomes k ri(EJK )A·BC = gi(EJK )

∑∫

Eo

J′K′

+

Zi(E′, E)gi(E′J ′K ′)dE′

q≠i

∑ kqi(EJK )gq(EJK ) + k d (EJK )gi(EJK ) i

q≠i

q≠i N

+

np

∑ k p(EJK )gi(EJK ) i

p=1

∑ kqi(EJK )gq(EJK , t ) − k d (EJK )gi

(17)

i

A K-active gi(EJ) can be derived following a similar procedure for the K-adiabatic case, except now the K-degree of freedom is not explicitly resolved, or rather it is assumed to have appropriately been averaged over to yield the K-active master equation, given by

np

∑ k p(EJK )gi(EJK , t ) p=1

i

(15)

where A(t)·BC(t) denotes the time-dependent bimolecular reactants A and BC, with an association rate constant kri(EJK) into well i. The kiq(EJK) is the unimolecular rate coefficient for isomerization of the complex from well i to the complex in well q, and kqi(EJK) is the rate constant for the reverse process for the complex from well q to that in well i. The kdi(EJK) is the dissociation rate constant for the ith chemical intermediate, for well i, dissociating to the original reactants (A and BC); kpi(EJK) is the analogous dissociation rate constant from well i to a set of bimolecular product p, for np product sets. Here, for brevity of notation, the Eo denotes E(i) o , the ground-state energy of isomer i, and may be noted in subsequent appearances in the case of multiple potential wells. The K-active version of eq 15 is given by5,30

k ri(EJ )A·BC = gi(EJ )

∑∫

∑∫ ∑∫

q≠i

∑ k p(EJ )gi(EJ ) i

(18)

where all the terms in eq 18 are now the K-active counterpart of those defined earlier for eq 17.

3. ANALYTIC SOLUTIONS FOR THE K-ADIABATIC AND K-ACTIVE BIMOLECULAR MASTER EQUATIONS FOR A SINGLE AND MULTIPLE POTENTIAL WELLS 3.1. Single Potential Well. 3.1.A. K-Adiabatic Population g(EJK) for a Single Potential Well. We next wish to acquire the analytic solution, that is, population, of the bimolecular master equation, where the single potential well case is considered in this Section 3.1. The adopted strategy is to solve eq 13 for g(EJK) by recasting eq 13 form of Fredholm equation and then find its solution using an algebraic method. Recasting the Kadiabatic master equation, eq 13, as a Fredholm of the second kind equation,31−34,67

N q≠i np

∑ k p(EJ )gi(EJ , t ) p=1

N

∑ kiq(EJ )gi(EJ ) − ∑ kqi(EJ )gq(EJ )

p=1

∑ kiq(EJ )gi(EJ , t )+∑ kqi(EJ )gq(EJ , t )

− k di(EJ )gi(EJ , t ) −

Zi(E′, E)gi(E′J ′)dE′

+ k di(EJ )gi(EJ ) +

Zi(E′J ′, EJ )gi(E′J ′, t )dE′

q≠i



np

Zi(EJ , E′J ′)dE′+

N

Zi(E , E′)dE′−Pi(J )

o

q≠i



Eo

J′



N

+



Eo

J′

∫E

Eo

J′

= k ri(EJ )A(t ) ·BC(t ) − gi(EJ , t )



o

N

∑ kiq(EJK )gi(EJK , t )

(EJK , t ) −

dt



∑ kiq(EJK )gi(EJK )−

+

Zi(E′J ′K ′, EJK )gi(E′J ′K ′, t )dE′

q≠i

dgi(EJ , t )

∫E

N

Zi(EJK , E′J ′K ′)dE′+

N



Zi(E , E′)dE′

o

J′K′



Eo ∞

J′K′



−Pi(JK ) ∑

= k ri(EJK )A(t ) ·BC(t ) − gi(EJK , t )

∑∫

∫E

i

y(x) − λ

(16) E

∫a

b

k(x , t ) × y(t )dt = f (x) DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

is governed by the conditions of E ≥ E′ or E′ > E, as each region is considered for Z(E′,E) of interest. On the basis of eqs 19−25, we obtain the following results regarding the conditions on the solution g(EJK) for each energy region: (1) If 1 ≠ λ̂, (i.e., E > E*(JK)), then for an arbitrary right-hand side of eq 19 there exists a unique solution for eq 19. Upon inserting  , eq 23, into eq 22, yields the form69

a ≤ x ≤ b, yields g (EJK ) − λ(EJK ) ∑ J′K′

= f (EJK )

∫E



Z(E′, E)g (E′J ′K ′)dE′

o

Eo ≤ E ≤ ∞

(19)

where f(EJK) = 0 for E ≤ E*(JK), and upon comparing the terms between eq 13 and eq 19 assists to identify k (EJK )A·BC f (EJK ) = r Z LJ + kd(EJK )

g (EJK )+ = f (EJK ) + (20)

∫E



o

⎤−1 Z(E , E′)dE′ + kd(EJK )⎥ ⎦

(21)

In the main text, we present the derivation for the solutions of eq 19, for the energy regions of E > E*(JK) and E ≤ E*(JK), with proofs given in Appendix C, and make a note when Z(E′,E) in eq 19 is factorizable into Z(E)Z(E′) at the appropriate junctions. In the Appendix E, we present a derivation when Z(E′,E) is factorizable, where then Z(E) would be brought outside the integral in eq 19 and get absorbed in λ(EJK), and then derive the solution for the single exponential model. The contribution of population to the recombination rate constant, based on whether Z(E′,E) is factorized or not in g(EJK)+, is found to be similar in magnitude and further discussed in ref 68. We seek a solution of the K-adiabatic master equation, eq 19, in the form g (EJK ) = f (EJK ) + λ(EJK )Â

g (EJK )− = λ(EJK )− Ĉ

(22)

(23)

Z(E′, E)f (E′J ′K ′)dE′

E

(24)

and λ̂ =



∑∫* J′K′

E

(29)

where λ(EJK)− ≠ 0. The introduction of a delta function in cE ensures the solution, eq 27, satisfies the master equation, while information is retained about E′ via c′E due to collisional energy transfer process, where c′E has been determined in previous works for various energy transfer models16,42 and further delineated in ref 42. Upon noting eqs 23−29, then the following boundary condition holds between the g(EJK)+ and g(EJK)− domains:



J′K′

(28)

Ĉ = P(JK )geq (0)c Ee−E / kT /λ(EJK )−

∑∫*

(27)

where now f(EJK) = 0 in the right-hand side of the master eq 19, since kr(EJK) = 0 in the numerator of f(EJK), as the associative and dissociative channels are absent for E ≤ E*(JK). Eq 27 is a solution of the master equation, eq 19, as confirmed upon its substitution into the master equation, for E ≤ E*(JK) (proof provided in Appendix C), where Ĉ in eq 27 is a constant to be determined. To determine Ĉ , we utilize the boundary condition,30,42 g(EJK)− = geq(EJK), when E → −D (D is the dissociation energy of ABC molecule measured from the bottom of its potential well), where geq(EJK) denotes the equilibrium probability density that ABC has energy E, defined by geq(EJK) = ρ(EJK)e−E/kT/Q′, where ρ(EJK) is the density of states of the molecule, and Q′ is the partition functions of ABC in the centerof-mass system of coordinates. Corrections to the equilibrium density have also been derived for various energy transfer models and are noted here as c′E.16,42 For example, the contribution to population from a collisional energy transfer process with a single exponential model is c′E = (1 − γ′/γ eE/kT).16,42 Writing geq(EJK) as P(JK)geq(0)e−E/kT and equating the latter with g(EJK)− in eq 27, so to be P(JK)geq(0)cEe−E/kT = λ(EJK)−Ĉ , and including the cE = δ(E − E′)c′E, yields for Ĉ

where the identifications are made for f1 =

E ≤ E*(JK )

λ(EJK )− = P(JK )

f1 [1 − λ ]̂

(26)

with

where g(EJK) is later below partitioned into g(EJK)+ and g(EJK)− for the energy regions of E > E*(JK) and E ≤ E*(JK), respectively, and the conditions noted for each and the matching boundary condition are given for the populations. The form of the solution with two superposed terms emerges, because, upon the consideration for the region of E > E*(JK), the first term f(EJK) in eq 22 assists in canceling the same on the right-hand side of the master equation in eq 19, while the second term in eq 22 works in tandem with its first term inside the integrand of eq 19 in satisfying the master equation by rendering the equality put forth by eq 19, as the proofs inform in the Appendices C and E. The  in eq 22 is a coefficient determined by substituting g(EJK), eq 22, into the bimolecular master eq 19, (steps are delineated in Appendix C) to yield  =

E > E*(JK )

and is a solution of the bimolecular master equation, eq 19, where its proof is given in Appendix C, with the properties g(EJK)+ = g(EJK) for E > E*(JK), and g(EJK)+ = 0, E ≤ E*(JK). (2) If 1 = λ̂ and f1 = 0, (i.e., E ≤ E*(JK)), where the dissociation and recombination events are absent in the master equation but only the energy transfer steps are available (e.g., reaction steps 2 and 3), then a solution of eq 19 can be represented in the form

and ⎡ λ(EJK ) = P(JK )⎢ ⎣

λ(EJK )f1 1 − λ̂

Z(E′, E)λ(E′J ′K ′)dE′ (25)

g (EJK )+ f (EJK ) + λ(EJK )Â

The integrals in eq 23 are supposed to exist, and indeed do, on physical grounds, for λ̂ and f1, where 1 ≠ λ̂ for E > E*(JK) and a separate solution found for E ≤ E*(JK). The counterpart of eqs 23−25 for a factorizable Z(E′,E) is given in Appendix E, as eqs E7−E9. The domain of E′ in the integral of eqs 24 and 25

=

g (EJK )− λ(EJK )− Ĉ

(30)

Further exposition on the boundary condition is discussed in Appendix F. Cases (1) and (2) are further treated in the present paper for their interesting physical relevance. F

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A Upon considering case (1), the E > E*(JK) region, and inserting the expressions for f(EJK), λ(EJK), f1, and λ̂, eqs 20, 21, 24, and 25, respectively, into g(EJK)+, eq 26, then yields70 ∞

P(JK )∑J ′ K ′ ∫ * f (E′J ′K ′)Z(E′, E)dE′ k (EJK )A·BC E g (EJK )+ = r +[Z LJ + kd(EJK )]−1 ∞ Z LJ + kd(EJK ) 1 − ∑J ′ K ′ ∫ * Z(E′, E)P(J ′K ′)[Z LJ + kd(E′J ′K ′)]−1 dE′ E

g (EJK )− = P(JK )geq (0)c Ee−E / kT

The transition probability Z(E,E′) in ⎡ λ(EJK ) = P(JK ) × ⎢ ⎣

∫E



o

⎤−1 Z(E , E′)dE′ + kd(EJK )⎥ ⎦

E



Z(E , E′)dE′ = Z LJ

where the lower integration limit on the first integral is extended to −∞ for mathematical convenience and the unity multiplied by the ZLJ. Upon considering both the up and down energy transition, then

∫E

(32)

to be the equilibrium distribution of population multiplied by a correction factor cE given elsewhere for a single exponential, double exponential, and near singularity models for intermolecular energy transfer.16,42 This solution is consistent with previous results using a Wiener−Hopf method,42 or the trial solution given by Troe,16 for the E-resolved case, when the sum over P(JK) is made to yield unity as it is normalized. The single exponential model for intermolecular energy transfer, which may be used in the expression for population, eq 31, or in the bimolecular master equations, is given by1,16

satisfies the completeness requirement15,71,72 in the form

∫−∞ Z(E , E′)dE′ + ∫E

(31)



Z(E , E′)dE′

⎧ Z e−(E − E′)/ γ ⎪ o Z (E ′ , E ) = ⎨ ⎪ Z e−(E′− E)/ γ′ ⎩ o

o

may be set to ZLJ, as done in eq 31, where Z(E,E′) is assumed normalized, multiplied by the Lennard-Jones frequency. When Z(E′,E) may be factored as Z(E)Z(E′) in the master equation eq 19, then Z(E′) would appear in lieu of Z(E′,E) in the expression of the numerator in eq 31, as its derivation is shown in Appendix E (cf. eq E14). On checking the upper limiting value of g(EJK)+, in the region of E > E*(JK) away from the reaction threshold as energy increases, we find as E → ∞, [ZLJ + kd(EJK)] > kr(EJK), since the collision frequency via ZLJ increases, and kr → 0 according to eq 7, so the first term and, similarly, the second term approach zero in eq 31, hence, g(EJK)+ → 0 as expected. We now consider case (2), E ≤ E*(JK), for g(EJK)−. Upon substituting for λ(EJK)− and Ĉ , eqs 28 and 29, respectively, into eq 27 for g(EJK)− yields

E′ ≤ E E′ ≥ E

(33)

where γ and γ′ are the deactivation and activation constants, which are related to each other by detailed balance, that is, γ′−1 = γ−1 + kT−1.42 The Zo contains the normalization constant, 1/(γ + γ′), for the transition probability of a single exponential model and the Lennard-Jones frequency. The functional forms of Z(E′,E) for other intermolecular energy transfer models, for example, double exponential,1,36 Gaussian,1,37,38 and stepladder,39−41 are given in the Appendix B. For the region of E > E*(JK), upon making the substitution for f(E′J′K′) and Z(E′,E) from eqs 20 and 33 into eq 31 yields ∞

( )

P(JK )Zo∑J ′ K ′ ∫ * e−|E − E′| / γ ′ k r(E′J ′K ′)A·BC[Z LJ + kd(E′J ′K ′)]−1 dE′ k r(EJK )A·BC E g (EJK )+ = + ( ) ∞ Z LJ + kd(EJK ) [Z LJ + kd(EJK )][1 − Zo∑J ′ K ′ ∫ * e−|E − E′| / γ ′ P(J ′K ′)[Z LJ + kd(E′J ′K ′)]−1 dE′] E

for E > E*(JK). The reaction scheme 1−5 may be augmented by an additional channel described by a kr(E′J′K′). However, when the concentration of the reactants A·BC is deemed independent of (E′J′K′), for example, when the concentration of reactants produced in reaction step 4 of the kinetic scheme, with energy E′, is negligible compared to the starting initial concentration of reactants in the state (EJK), and the rotational distributions independent of J′K′, then for simplicity A·BC may be factored outside of the integral sign in the second term of eq 34, contained in the f(E′J′K′) term. Likewise, the prime appearing in kr(E′J′K′) for E′ may be dropped, and as well the primes on J′K′ (again when making the assumption of negligible populations from the counterpart primed channels relative to a dominating reaction step 1), in the second term under the integral sign. The reporting of calculations on ozone formation in the present work retains the primes in kr(E′J′K′).

(34)

For the evaluation of the integral in the second term of eq 34, the conditions from eq 33 are observed for the energy regions, where γ(′)in eq 34 denotes either selecting γ for the region E′ ≤ E, or selecting γ′ for E′ > E, with the appropriate limits of integration for the selection of either region,73 and with kr and kd defined by eqs 7 and 8, respectively, or determined by using other methods of interest. Eq 34 is derived under the consideration of both the up and down intermolecular energy transfers. In the case when Z(E′,E) is factorizable, for example, for the single exponential model, then a Z(E′) = eE′/γ would ( )

appear instead of e−|E−E′|/γ ′ in the integral of the numerator eq 34 for E > E′, or Z(E′) = e−E′/γ′ when E′ > E, and the corresponding Z(E) would be factored outside the integral, while in the integral of the denominator in eq 31 or 34, namely, in its λ̂ term eq 23, a Z(E′) would appear in lieu of Z(E′,E) that G

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

g(EJ)+ = kr(EJ)A·BC/kd(EJ), similar to the K-adiabatic result. For the region energy E ≤ E*(J)

ends up canceling as informed by eq E9 in Appendix E, to yield eq E14 or E15 as exact solutions, as counterparts to eqs 31 and 34, respectively. In the low-pressure limit, when ZLJ = 0, then the second term of g(EJK)+ is zero, since Zo in its numerator is zero, as the latter receives contribution from ZLJ. As pressure increases, in the lowpressure regime, but ZLJ < kd(E′J′K′), for example, at 1 bar and T = 300 K, ZLJ = 7.4 × 109 sec−1 between Ar and O3 and kd(E′J′K′) = 1.7 × 1012 sec−1, for example, when E = 1kT with T = 300 K, J = 20, K = 3, and approximately so for other neighboring values of interest, then the second term of g(EJK)+, expressed by eq 31 or 34, is approximately a factor of 1 × 10−4 to 1 × 10−6 smaller than the first term of g(EJK)+, depending on the (EJK) values of interest. In the latter regime, we signify that ZLJ < kd(E′J′K′) sufficiently, and so then the second term of g(EJK)+ is negligible,74 as further discussed in ref 74 and in the discussion in Section 7.1, and hence eq 31 or 34 reduce to k (EJK )A·BC g (EJK )+ = r Z LJ + kd(EJK )

g (EJ )− = P(J )geq (0)c Ee−E / kT

yi (x) − λi

×

1 − ∑J ′ ∫



E*

∫E



o

Zi(E′, E)gi(E′J ′K ′)dE′

= fi (EJK )

(39)

for Eo ≤ E ≤ ∞, and upon comparing the terms between eqs 17 and 39 assists to identify N −1

fi (EJK ) = λi(EJK )Pi(JK ) [∑ kqi(EJK )gq(EJK ) q≠i

+ k ri(EJK )A·BC]

(40)

and

∫E



λi(EJK ) = Pi(JK )[

Zi(E , E′)dE′

o

np

N

+

∑ kiq(EJK ) + k d (EJK ) + ∑ k p(EJK )]−1 i

q≠i

p=1

i

(41)

We seek a solution of the K-adiabatic master equation, eq 39, in the form gi(EJK ) = fi (EJK ) + λi(EJK )Â i

Z(E′, E)f (E′J ′)dE′ −1

Z(E′, E)P(J ′)[Z LJ + kd(E′J ′)] dE′

(42)

where gi(EJK) is later below partitioned into gi(EJK)+ and gi(EJK)− for the energy regions of E > E*(JK) and E ≤ E*(JK), respectively, and the conditions noted for each, and the matching boundary condition are given for the populations. The  i in eq 42 is a coefficient determined by substituting gi(EJK), eq 42, into the master equation, eq 39, (steps are delineated in Appendix D) to yield

where in the derivation, the K degree of freedom is appropriately averaged over at the outset. Previous comments following eq 31 regarding the factorizability of Z(E′,E) also apply to eq 36. The domain of E′ in the integrals of eq 36 is again selected according to when E ≥ E′ or E′ > E is under consideration for the energy transfer model of interest. Similar to the K-adiabatic result, eq 36 further simplifies to k r(EJ )A·BC Z LJ + kd(EJ )

ki(x , t ) × yi (t )dt = fi (x)

J′K′

(36)

g (EJ )+ =

b

gi(EJK ) − λi(EJK ) ∑

k (EJ )A·BC g (EJ )+ = r +[Z LJ + kd(EJ )]−1 Z LJ + kd(EJ ) ∞

∫a

a ≤ x ≤ b, yields

Eq 35 refers to a physical region, where the effect of collisional energy transfer from g(E′J′K′)+ to g(EJK)+ is negligible, because the population in the states g(E′J′K′)+ dominantly dissociates, since kd(E′J′K′) > ZLJ sufficiently to preempt a viable contribution to repopulate g(EJK)+. Its application to ozone formation is given in the results in Section 6. In the low-pressure limit, as ZLJ → 0, the second term in the right-hand side of eqs 31 and 34 becomes zero, as noted earlier, but further remark, that it then leaves the first term of the expressions for g(EJK)+, eqs 31 and 34, and along with eq 35, to further simplify to g(EJK)+ = kr(EJK)A·BC/kd(EJK), where the ZLJ that would have appeared in the denominator of the latter expression is now zero and so rendered absent regarding its appearance. The latter expression for g(EJK)+ can be readily evaluated numerically, and note the expression also resembles the familiar detailed balance condition, Keq = kr/kd, where Keq is the equilibrium constant for the association reaction.75,76 3.1.B. K-Active Population g(EJ) for a Single Potential Well. The K-active counterpart of g(EJK)+, for the energy region of E > E*(J), is

E*

(38)

is likewise acquired. 3.2. Multiple Potential Wells. 3.2.A. K-Adiabatic Populations gi(EJK) for Multiple Potential Wells. We now consider the analytic solution for the multiple potential well case. We solve the K-adiabatic bimolecular master equation with multiple potential wells, eq 17, for the population of the ith well gi(EJK), by recasting eq 17 form of Fredholm equation and then obtain its solution using an algebraic method, a strategy similar to the single potential well case. Recasting the K-adiabatic master eq 17 as a coupled system of Fredholm of the second kind equations

(35)

P(J )∑J ′ ∫

E ≤ E*(J )

 i =

f1

i

1 − λî

(43)

where the identifications are made for

(37)



with a negligible second term in eq 36, when ZLJ < kd(E′J′) sufficiently, according to the physical condition discussed in the paragraph containing eq 35. In the low-pressure limit, as ZLJ → 0, then the expression for g(EJ)+, eqs 36 and 37, simplifies to

f1 = i

∑∫* J′K′

E

Zi(E′, E)fi (E′J ′K ′)dE′

(44)

and H

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A λî =



∑∫* J′K′

N

Zi(E′, E)λi(E′J ′K ′)dE′

E

∑ kqi(EJK ) × gq(EJK )

(45)

q≠i

The integrals in eq 43 are supposed to exist for λî and f1 , and

appearing in f i(E′J′K′)−, in eq 49, is further investigated in the Appendix G. In the case (b) that kpi channel is also closed for a given channel pi or for all values of pi, 1 through np, then correspondingly kpi = 0 in the terms appearing in eqs 47−52. The following boundary condition holds between the gi(EJK)+ and gi(EJK)− domains

i

indeed do on physical grounds, where kri ≠ 0 and kdi ≠ 0 in f1i, and kdi ≠ 0 in λî , for E > E*(JK). The counterpart of eqs 43−45 for a factorizable Zi(E′,E) in the master equation eq 39 would have Zi(E′) appear in lieu of Zi(E′,E), and the Zi(E) assigned to appear in the multiwell analogue of λ̅(E′J′K′), namely, as λ̅i(E′J′K′) = Zi(E)λi(E′J′K′) (cf. eq E2), similar to the singlewell case (Appendix E). The domain of E′ in the integrals of eqs 44 and 45 is informed by E ≥ E′ or E′>E, for the consideration of each case in Zi(E′,E). On the basis of eqs 39−45, we obtain the following results regarding the conditions on the solution gi(EJK): (1) If 1 ≠ λ̂i, (for E > E*(JK), then for an arbitrary right-hand side of eq 39 there exists a unique solution for eq 39. Upon substituting for  i into eq 42, then it can be written in the form gi(EJK )+ = fi (EJK ) + λi(EJK )

f1

gi(EJK )+ fi (EJK ) + λi(EJK )Â i

i

i



λî =

p=1

but the latter equality is made invalid by kdi(EJK) > 0 for this region E > E*(JK), and so renders 1 ≠ λî . The condition 1 = λ−î , for E ≤ E*(JK), implies that np

N

∑ kiq(EJK ) + ∑ k p(EJK ) = 0 i

q≠i

p=1

and so N

∑ kiq(EJK ) = 0 q≠i

and

(47)

np

∑ k p(EJK ) = 0 i

f 1−

p=1

i

1 − λî



that is, these channels are closed at the given E(JK), then when possible and physically elucidating for the reactions, as the potential energy surface may also inform, the E*(JK) may be instead selected such that at least one channel associated with either kiq or kpi is open, to render 1 ≠ λ−î , and the solutions from cases (1) and (2) would then follow, or the case be excluded. Similar to the single potential well case, when f i(EJK)+ > λi(EJK)+Â i sufficiently such that the contribution of the second term in eq 46 to gi(EJK)+ is negligible, then eq 46 reduces to

(48)

with f 1− =

i

q≠i

where Âi − =

(53)

np

N

as a solution of eq 39, where its proof is given in Appendix D. If the distribution functions for gq(EJK) are known, or when can be satisfactorily approximated by a valid manner, then they can be utilized in the f i(EJK) in the solution. (2) If 1 ≠ λî and E ≤ E*(JK) where (a) the dissociation kdi(EJK) and recombination kri(EJK) channels are absent in the master equation but only the energy transfer steps, the isomerization and the kpi channels are open, then the solution of the eq 39 still has the form of eq 46, but now the kdi(EJK) = 0 and kri(EJK) = 0 in the terms appearing in eq 46. So, the solution becomes E ≤ E*(JK )

fi (EJK )− + λi(EJK )− Â i−

∑ kiq(EJK ) + k d (EJK ) + ∑ k p(EJK ) = 0 (46)

gi(EJK )− = fi (EJK )− + λi(EJK )− Â i −

gi(EJK )−

upon noting 40, 41 and 43−52. Further exposition on the boundary condition is discussed in Appendix F. For the last conditions on the solution of gi(EJK), (3) the condition 1 = λî , for E > E*(JK), does not physically arise, since it implies that

E > E*(JK )

i

1 − λî

=

∑∫ J′K′

Zi(E′, E)fi (E′J ′K ′)− dE′

Eo

∑∫ J′K′

E*

(49)

E*

Eo

Zi(E′, E)λi(E′J ′K ′)− dE′ (50)

gi(EJK )+ = fi (EJK )+

N −1

fi (EJK )− = λi(EJK )− Pi(JK )

∑ kqi(EJK )gq(EJK ) q≠i

and similarly eq 47 would simplify to

(51)

gi(EJK )− = fi (EJK )−

and

∫E

λi(EJK )− = Pi(JK )[

Zi(E , E′)dE′ +

o

∑ kiq(EJK ) q≠i

np

+

∑ k p(EJK )]−1 p=1

i

(55)

for its consideration. In the low-pressure limit as ZLJi → 0, then the second term of the right-hand side of eq 46 is zero, since f1i = 0 and 1−λ̂i = 1, so to then yield

N



(54)

lim g (EJK )+ = fi (EJK )

Z LJ → 0 i

(52)

i

with proof is provided in Appendix D. The term

where now I

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

∫E



4.2. Multiple Potential Wells. At the high-pressure limit, as ZLJi → ∞, the K-adiabatic bimolecular master equation for multiple potential wells, eq 17, may reduce to

Zi(E , E′)dE′ = 0 0

in f i(EJK), eq 40, and in eq 54. Similarly, for the energy region of E ≤ E*(JK),

k r+i, ∞(EJK ) = Ni*(EJK )e−E / kT /hQ i np

lim g (EJK )− = fi (EJK )−

Z LJ → 0 i

+

i

∑ N *p (EJK )e−E /kT /hQ i p=1

where the collision term describing energy transfer is now absent in f i(EJK)−, eq 51 and eq 55. 3.2.B. K-Active Populations gi(EJ) for Multiple Potential Wells. The K-active counterpart of gi(EJK)+, eq 46, is gi(EJ )+ = fi (EJ ) + λi(EJ )

f1

N

+

N



gi(EJ )− = fi (EJ )− + λi(EJ )−

− 1 − λî

np

k r+i, ∞(EJ ) = Ni*(EJ )e−E / kT /hQ i +

E ≤ E*(J ) (57)

N

+

i

N

∑ kiq(EJ )gi(EJ )+ − ∑ kqi(EJ )gq(EJ )+ q≠i

q≠i

(61)

per reactant concentrations, where the K degree of freedom was averaged at the outset in the derivation. In effect, the K-adiabatic and K-active recombination rate constants are now described by terms resembling the high-pressure RRKM theory limiting forms (i.e., the first and second terms of eqs 60 and 61) and augmented by isomerization terms. The high pressure and temperature limits of the rate constants of eqs 58 and 60 are discussed in Section 7.2. For the energy region of E ≤ E*(JK), kr−i,∞(EJK) = 0, as noted earlier for the single potential well case, since the association channel is closed by definition, and similarly k−ri,∞(EJ) = 0, valid, of course, also at the high-pressure limit. However, if the bimolecular product channels are open, then the relationship between any isomerization processes coupled with the processes associated with the bimolecular product channels is given by

4. BIMOLECULAR MASTER EQUATIONS AT THE HIGH-PRESSURE LIMIT AND RECOMBINATION RATE CONSTANTS 4.1. Single Potential Well. At the high-pressure limit, as ZLJ → ∞, the K-adiabatic bimolecular master equation, eq 13, may reduce to yield a recombination rate constant (58)

per reactant concentrations, for the energy region of E > E*(JK), where Q is the partition function of the reactants in the centerof-mass system of coordinates. The derivation is given in Appendix H, with the tracing of the units discussed therein. The unit of k+r,∞(EJK) is (concentration)−1 (time)−1 per (EJK). It may be noted that eq 58 displays the high-pressure limiting form of bimolecular RRKM theory,24 as the latter consideration emerges and is noted from the initial derivation by Marcus,9 where presently considered for an (EJK) resolved representation. A similar analysis for the high-pressure limit expression for the K-active master equation yields k r,+∞(EJ ) = N *(EJ )e−E / kT /hQ

∑ N *p (EJ )e−E /kT /hQ i p=1

is acquired. The conditions on 1 − λî and 1 − λ−î are similar to the aforementioned K-adiabatic case. The K-active analogue of eqs 54 and 55 may be noted upon the condition f i(EJ) > λi(EJ)Ai and f i(EJ)− > λi(EJ)−Ai−, respectively. Also, the low-pressure limits of gi(EJ)+ and gi(EJ)− yield the K-active version of the results given earlier for the Kadiabatic case.

k r,+∞(EJK ) = N *(EJK )e−E / kT /hQ

(60)

per reactant concentrations, for the energy region of E > E*(JK). The derivation is given in the Appendix H and unit discussed therein. A similar analysis for the high-pressure limit expression for the K-active master equation yields

where in the derivation the K degree of freedom was averaged over at the outset, in each term. Likewise for the energy region E ≤ E*(J) i

∑ kqi(EJK )gq(EJK )+ q≠i

(56)

f 1−

∑ kiq(EJK )gi(EJK )+ q≠i

E > E*(J )

i

1 − λî

i

np

∑ k p(EJK )gi(EJK )− p=1

i

N

=

N

∑ kqi(EJK )gq(EJK )− − ∑ kiq(EJK )gi(EJK )− q≠i

q≠i

(62)

at the high-pressure limit, with the derivation given in Appendix H. In the case when the pi channels are closed, then kpi = 0, and eq 62 becomes N

(59)

per reactant concentrations, where the K degree of freedom was averaged at the outset in the derivation. The averaging over the (EJ) or (EJK) degrees of freedom then, of course, yields the thermal high-pressure limit kr,∞. For the energy region E ≤ E*(JK), k−r,∞(EJK) = 0, as the association channel is closed by definition, and similarly k−r,∞(EJ) = 0, at the high-pressure limit.

N

∑ kqi(EJK )gq(EJK )− =

∑ kiq(EJK )gi(EJK )−

q≠i

q≠i

(63)

In this scenario, detailed balance77 may hold between gi(EJK)− and gq(EJK)−.77,78 The potential energy surface may inform if any isomerization process is taking place in this energy region, for example, when E is greater than a low activation barrier between the two potential wells, while E ≤ E*(JK), and so to establish some population via activating collisions above J

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

then typically the ratio of the second to the first term of the population is to within 5% at 100 bar, (cf. at 1 bar, the ratio is in the range of ∼1 × 10−4 to 1 × 10−6 as noted in Section 2). So, in this pressure region near 100 bar, the population of metastable Ο3*(EJK) still receives a small contribution from the energy transfer events associated with Ο3*(E′J′K′), up to ∼5%, since the dissociation events of the latter dominates the energy transfer events associated with it. This is analyzed in Section 7.1. The maximum magnitude of the population for O*3 (EJK) takes place at E = 1kT above the dissociation limit, J = 20, K = 3, where g(EJK)+ = 1.0 × 10−27 cm3/wavenumber per (EJK) and per reactant concentrations for T = 300 K and 1 bar, where the units is in accord with ref 60 and the (2J + 1) weighting factor for N*(EJK) is included. Similar results are found as pressure increases to 100 bar. This (EJK) region with its higher energy relative to a lower E = 0.2kT region still makes a contribution to the krec of O3, that is ∼10% at 1 bar and up to 100 bar, after taking into account the exponential weighting e−E/γ due to the intermolecular energy transfer contributions for a single exponential model in the integral of eq 65, where results on the krec and γ are given in Section 6.2 for T = 300 K and 1 bar. Additional results on the temperature and pressure dependence of population will be given elsewhere. We note, as earlier in Section 1, the ozone molecule is an example of a system where its K has been shown to exhibit adiabaticity27 with diffusive aspects28,29 at the lower pressures, and in the range of T = 373 K to 130 K in one study,14 and so the K-adiabatic have been of initial focus here, where the K-active contributions to krec may also further be investigated as pressure increases towards 1000 bar. Lastly, we remark that the numerical implementation of eqs E15, 31, 34, and 35 is found to be facile.79 6.2. krec. The application of eq 65, that is, the expression for the thermal krec accommodating general pressure employing g(EJK)+, to ozone formation at low pressure, 1 bar, in the bath gas Ar yields krec = 4.0 × 10−34 cm6 molecule−2 s−1, in agreement with the same value from experiment,53 at T = 300 K and low pressure. The agreement is rendered upon the selection of a deactivation energy transfer γ = 67 cm−1 for a single exponential model1,16 of intermolecular energy transfer. The activation energy γ′ is determined from detailed balance, that is, γ′−1 = γ−1 + kT−1.42 The value of γ is an approximation and one that is found to also approximately agree from the values reported from classical trajectory59 and RRKM-based bimolecular studies14 for ozone formation in Ar. The numerical implementation of eq 65 is noted to be facile.79 The temperature and pressure dependence of krec for ozone formation, using eq 65, will be investigated elsewhere for their experimentally interesting aspects. 6.3. High-Pressure Limit k+r,∞(EJK) and Thermal kr,∞. We consider the application of the high-pressure limit of the recombination rate constant k+r,∞(EJK), expressed by eq 58, to ozone formation. As an example, the calculation yields k+r,∞(EJK) = 1.3 × 10−16 cm3 molecule−1 s−1 for O3 when E = 0.2kT above the dissociation limit, J = 10, K = 3 for a region of interest. Upon averaging the k+r,∞(EJK), over the (EJK) region then yields the thermal high-pressure rate constant, kr,∞ = 1.5 × 10−11cm3 molecule−1 s−1 at T = 300 K. The experimental kr,∞ for O3 in Ar recorded at the highest pressure, 1000 bar, is reported to be 4.6 × 10−12 cm3 molecule−1 s−1 at T = 300 K,53 and so has not yet reached the high-pressure limiting value from RRKM theory.

the energy barrier for a given potential well at thermal conditions, where increasing temperatures may further assist the process.

5. RECOMBINATION RATE CONSTANT AS A FUNCTION OF PRESSURE AND TEMPERATURE Hitherto, we have been considering the (EJK) or (EJ) resolved populations and microcanonical rate constants, where the latter appeared in the expression defining the K-adiabatic or K-active populations. We next consider averaging over these microcanonical degrees of freedom of the population so to acquire the thermal recombination rate constant krec for product formation from the bimolecular reaction as described by the reaction scheme 1−5. The motivation behind this quest is that the thermal rate constant is an observable that is typically measured in kinetic-based experiments. The thermal recombination rate constant based on the analytic-based population is given by k rec =

∑ ∫ ∫ Z(E , E′)g(EJK )dEdE′ (64)

JK

where the population g(EJK) may, for example, be represented by using eq 26 or eq E3 in Appendix E, upon relevance. The limits of integration for E is (E*,∞) and for E′ is (−D,E*) for a successful third-body collision, where D is the dissociation energy of the molecule measured from the bottom of the potential well, and E* is the critical energy for the reaction at a given (JK). The pressure and temperature dependence of krec is contained in g(EJK) and Z(E,E′). The unit of krec is volume2 molecule−2 time−1. This expression, eq 64, is the analogue of the krec based on RRKM theory of ref 14, where now g(EJK) obtained from the bimolecular master equation appears in lieu of ρ(EJK)e−E/kT/Qr in eq 16 there, where Qr is the partition function of the colliding reactants in the center of mass system of coordinates. The K-active analogue of eq 64 would have K appropriately averaged in its expression. Upon considering a single exponential model for Z(E,E′), that is, eq 33, into eq 64 then E′ may be analytically integrated over to yield k rec =

Zo ∑ (γ + γ ′) JK



∫E*

g (EJK )e−E / γ dE (65)

where those E′ less than E*(JK) contribute to the bimolecular formation of the stable product, and here Zo is the LennardJones collision number with units of volume molecule−1 time−1. The krec may, of course, also be considered for a K-active g(EJ). The g(EJK) contributes a volume·molecule−1, and so krec acquires the unit of volume2 molecule−2 time−1.

6. RESULTS ON OZONE FORMATION 6.1. g(EJK)+. An applicable regime of eq 35, that is, when ZLJ < kd(E′J′K′), is to the reaction of ozone formation, that is, O + O2 + M →O3 + M, up to ∼100 bar at T = 300 K, where the population via eq 35 contributes 96% to the thermal recombination rate constant of O3 via eq 65 at 100 bar in bath gas M = Ar, and the remaining 4% is contributed to the krec from the second term of the population from either eq 31 or 34. At the latter pressure, the ZLJ = 6.7 × 1011 sec−1 in Ar bath gas is less than kd(E′J′K′) = 1.7 × 1012 sec−1 of Ο3*(E′J′K′), by a factor of 2.5, for the (E′J′K′) region of interest, for example, E′ = 1kT above the dissociation limit at T = 300 K, J′ = 20, K′ = 3, where K

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

7. DISCUSSION 7.1. Analysis of g(EJK)+. It was noted earlier in Sections 3.1 and 6.1 that, when ZLJ < kd(E′J′K′) sufficiently such that the second term of the solution makes a negligible contribution to g(EJK)+ for general reaction of interest, then g(EJK)+ described by eq 31 or 34 may further be simplified to g(EJK)+ = kr(EJK)A· BC[ZLJ + kd(EJK)]−1, eq 35, by retaining the first term of either equation only. In the absence of the second term, the numerical implementation of g(EJK)+, of course, also simplifies. We next further examine the validity of this regime by illustrating it using the ozone formation reaction: for a region of interest regarding (E′J′K′), for example, E′ = 1kT with T = 300 K, J′ = 20, K′ = 3, we have kd(E′J′K′) = 1.7 × 1012 s−1 for Ο3*(E′J′K′) and ZLJ = 7.4 × 109 s−1 between O3 and Ar, at 1 bar, where for consistency the same (E′J′K′) values are considered as those given earlier in the result Section 6.1. If say 1000 molecules comprise the initial metastable population of Ο*3 (E′J′K′), then by 200 ps, the kd(E′J′K′) informs that 340 molecules or 34% of the initial population have dissociated to the reactants O + O2, while the ZLJ informs that ∼1 collision between Ο*3 (E′J′K′) and Ar has taken place during this time to contribute to g(EJK)+ of O3(EJK) from g(E′J′K′)+, which in turn the metastable complex may or may not lead to form a stable O3 molecule. Then, after 3 × 200 or 600 ps, all the population of ozone in the state (E′J′K′) has dissociated to O + O2, while during this same time (7.4 × 109 collisions sec−1) × (600 ps) or ∼4 collisions may have taken place to contribute to g(EJK)+ of O3, and so the maximum relative possible contribution to g(EJK)+ from g(E′J′K′)+ via collisions is 4/1000 or 0.4%, as the dissociation mechanism dominates. Thus, in this regime of ZLJ < kd(E′J′K′), the second term makes a negligible contribution to g(EJK)+, compared to the first term, and so eq 35 would be valid. In the case of ozone formation, this regime is operative at low pressures and temperatures ranging from T = 163 to 1000 K, where at T = 300 K, the significant contribution from (EJK)-based population of Ο3* to the thermal krec arises from regions about the excess energy above the dissociation limit (E = 1 kT, J = 20), (E = 0.5kT, J = 10), and (E = 0.2 kT, J = 5), for K = 3. The natural cutoff of the adiabatic K associated with the flux at the transition state is ∼5, upon noting conservation of E and J in the Hamiltonian of ozone with a hindered rotor80 and using a potential energy surface81 of O3. Eq 35 may also be arrived at, in the aforementioned regime, when Z(E′,E) is assumed factorizable in the derivation of g(EJK)+, and the neglect of the second term in eq E15 is again corroborated by numerical results regarding O3 formation. An alternative physical path for arriving at eq 35 is discussed next. One may foundationally arrive at this simpler expression for g(EJK)+, eq 35, via an alternative physical approach, in the following way: let the probability of survival of a metastable molecule, for example, Ο3*, in the presence of collisions be Pc(EJK,t). It is defined as the product of two independent probabilities: the survival probability in the absence of collisions P(EJK,t), multiplied by the probability Pnon(EJK,t) that the molecule has not undergone a collision satisfying −dPnon(EJK,t)/dt = Pnon(EJK,t)ZLJ, where ZLJ is the frequency of collisions with third bodies that depletes the energy of an intermediate, for example, Ο*3 , to below the critical energy E*(JK). The solution of the latter differential equation is, of course, Pnon(EJK,t) = e−ZLJt, and so Pc(EJK,t) = P(EJK,t) Pnon(EJK,t) = P(EJK,t)e−ZLJt. We next incorporate the expression for Pc(EJK,t) into the rate expression14 with (EJK) resolved as

k rec(EJK ) =

1 h

∫ ∫ P(EJK , t )W (EJK )

× Z(E′, E)e−Z LJt e−E / kT dE′dt /Q

where W(EJK) is a dimensionless quantity related to the flux and equivalent to N*(EJK),13,14,82 and upon letting P(EJK,t) = e−kd(EJK)t, where kd(EJK) = N*(EJK)/hρ(EJK) from RRKM theory1,15,23,24 and analytically integrating over the time t in the range (0,∞) then yields N*(EJK)e−E/kT/Q[ZLJ + kd(EJK)], which is identified as g(EJK)+ = kr(EJK)[ZLJ + kd(EJK)]−1, equivalent to eq 35, upon noting kr(EJK) = N*(EJK)e−E/kT/Q and taking account of the reactant concentrations, so to complete the derivation as aimed to show. In turn, the krec, eq 64, may accommodate this simpler g(EJK)+ for this regime, with a facile numerical implementation. We also conclude, based on the aforementioned analysis, that the RRKM-based bimolecular krec, eq 9 in ref 14, there described as a simple “bird’s eye view” of the pressure dependence, is now understood to return an equivalent result for krec as that based on the population from the bimolecular master equation, since its integrand is now identified to contain g(EJK)+, eq 35, when the regime ZLJ < kd(EJK) is sufficiently valid, for example, up to ∼100 bar for ozone formation at T = 300 K. Herein, we briefly comment on the rate constants found in the expressions of g(EJK)+ for the single well or multiple potential wells. The microcanonical rate constants k(EJK) or k(EJ), describing dissociation, association, or isomerization in the multiple well case, found in the analytic expression of the populations, via 47, 56, or 57, of course, pertain to an elementary reaction such as for an isomerization step between two potential wells in a multiple-well case, treated by the master equation in the present formalism. The contribution to the population in well i, from a specific well q, may precisely be ascertained by isolating the channel and discounting the contribution of other channel(s) in the analysis of the expressions. A phenomenological rate constant, as a function of thermodynamic variables,83 based on the population may then be defined for the connecting wells of interest, to which the wells may or may not be adjacent, 84 that is, nonlocal, as done by Miller and Klippenstein85 in describing a “well-mediated”84 reaction. For example, a numerical-based treatment of the master-equation formalism6 has utilized its population to yield phenomenological rate coefficients to account for isomerization reactions that “skip a well” in a well-mediated reaction, for example, allene with propene, in a cyclopropene-mediated reaction using an eigenvalue method for obtaining populations and thereafter acquiring the phenomenological rate coefficient(s).6,85 7.2. High Pressure and Temperature Limits of Rate Constants. At the high-pressure limit, a potential well that directly connects with the bimolecular reactant channel is expected to become stably populated at lower temperatures, for example, at T = 300 K, where a deeper well would accommodate a more stabilized product. This is so because eq 58, defining k+r,∞(EJK), or the K-active version eq 59, assumes that a reactive flux has passed the TS and is on its way to form a stable product in its associated well for a finite temperature of interest, moduli any recrossings of the TS, which then may receive a correction. However, as temperature increases toward infinity, deviations from the high-pressure limit k+r,∞ at finite temperature is expected by the way of the complex not stabilizing in its potential well. We commence to investigate the latter point further by using an analytic analysis of k+r,∞ and also draw upon numerical results reported in the literature. L

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

reaction yields the dissociative product, phenyl radical + H.66 The distinct physics that has been noted66 to be responsible for the aforementioned scenario is that the average intermolecular activation energy transfer to an intermediate complex from collisions with bath gas has become comparable to the average deactivation energy in this temperature range, for a net ⟨ΔE⟩ = 0, so to yield no stable product with energy less than the critical energy of any of the potential wells. Although, in the study66 the highest temperature in this regime is less than ∞, where the thermal rate coefficient for product formation was noted to be independent of pressure in the range studied, the outcome of the propargyl recombination reactions reflects the physical behavior of a bimolecular reaction expected at T → ∞. To expound, the high-pressure limiting expression of rate constant of a stable product from the master equation, eq 60 or eq 61 does not possess terms in its expression regarding energy transfer and assumes all of the reactive flux that pass the TS forms product, unless T → ∞, which then no stable product forms as the limiting behavior of the expression shows. Similar in outcome, no stable population of an isomer is obtained from the master equation at the higher T values, greater than 2000 K, where the ⟨ΔE⟩ = 0 is found at high but finite pressures of 100 atm, and so in effect rendering the energy-transfer process ineffective to cause stabilization of an isomeric complex. However, as temperature is lowered to the range of 500 K < T < 2000 K, a region of interest where combustion takes place, then various isomers of C6H6 are found to stabilize by the multiple collisions.66 The microcanonical rate constant for isomerization between isomers in the high-pressure limit displayed in the third and fourth terms in eqs 60 and 61 may be finite, rather than zero, and an isomerization event may take place above or below near the critical energy of reaction, with the likelihood of the event increasing as temperature increases, as illustrated by two scenarios: (a) Collisions may assist the reaction coordinate of the complex, for example, the internal rotation of a halogenated substituted ethane from a gauche to an anti staggered conformation with a rotational barrier height of 3−4 kcal/mol for a monosubstituted complex, with the initial energy less or greater than, but near, the critical energy in a potential well A, to then traverse the transition state and then reach well B to descend below the critical energy or be excess of it, and then vice versa from B to A. In this manner, technically, an isomerization event has taken place, to above or below E*(JK). A key physical ingredient in the admission of this process is the breaking of a perfect isotropic collisional interaction between each atom in a complex with the bath gas at any given time, and a sufficiently high temperature. Otherwise, the net force on a given atom by the bath gas at a given instant may vectorially cancel in the highpressure limit. Upon this admission, when a metastable complex initially forms at the TS connecting two potential wells, then a higher-symmetried structure of the activated complex at the TS may favor on average a relatively enhanced statistical branching into either well compared to one of a lower symmetry. (b) A low-symmetry structured complex populated at an energy near but less than the critical energy may have its energy increased and the reaction coordinate reach the TS by the assistance of collisions, where “nonstatistical” dynamical effect may take place at sufficiently high temperature, for example, collision-assisted roaming in the densely packed environment, where ABC becomes AB−C at the TS, evolves to BA−C as a function of time, and then forms the isomer BAC. Since collision frequency with the bath gas is infinite at the high-pressure limit, the

We consider the behavior of the recombination rate constant for a product formation at the high-pressure limit as temperature increases to infinity for reactions with a single well and multiple potential wells. For illustration, we consider the (EJK) case, whereof similar forthcoming considerations would also apply to the K-active case. For the single potential well case, as T → ∞ in eq 58, then the high-pressure limit k+r,∞(EJK) → 0, and so no stable complex with energy less than the critical energy E*(JK) is expected to form in the single potential well. This may be arrived at by noting that

lim e−E / kT = 1

T →∞

in the numerator of eq 58, while lim Q = ∞

T →∞

in the denominator of eq 58 upon using l’Hôpital’s rule86,87 on Q (delineated in ref 88), where Q is the reactant partition function, and so88 lim k r,+∞(EJK ) ∝ 1/∞ = 0

T →∞

So, the bimolecular reaction involving reactants A and BC would simply yield a fragmented product, for example, AB and C, rather than a stable complex ABC. For the multiple potential well case, lim k r+i, ∞(EJK ) → 0

T →∞

is also expected regarding a complex formation in the ith potential well, by the following analysis of k+ri,∞(EJK), eq 60. The first and second terms of eq 60, similar to those of eq 58, each approach 0 as T → ∞ using the same analysis used for the single-well case. The remaining terms for consideration in eq 60 are the third and the fourth terms, N

∑ kiq(EJK ) × gi(EJK )+ q≠i

and N

∑ kqi(EJK ) × gq(EJK )+ q≠i

respectively, describing isomerization processes. Upon assuming the gi and gq to have a functional form ρ(EJK)e−E/kT/Q, where Q is the partition function for the complex in its center of mass, then lim Q = ∞

T →∞

and so contributions from the isomerization terms approach 0.89 Thus, as T → ∞, a stable complex does not form in the potential well i, but rather a fragmented product is to be expected from the initial bimolecular reaction. The latter result is evidenced in the trends found in propargyl (C3H3) radical recombination, which is an important cyclization step in flames of aliphatic fuels.66 At T > 2000 K and high but finite pressure, 100 atm, the numerical results from a master equation analysis show that as radical C3H3 combines with radical C3H3, no stable complex forms in any of the isomeric potential wells, or expressed more precisely, 13 isomers exhibit a favorable energetic accessibility in the purview of study,66 with benzene, fulvene, and 2-ethynyl-1,3butadiene having the deepest wells, and instead the bimolecular M

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A rotation of AB to BA by π radians to form BAC may then be assisted by collisions at high temperature and the asymmetry of the complex, characterized by its asymmetric mass distribution or structure, for the execution of the π rotation of say the atom or diatom. Otherwise, from a statistical perspective of a strictly isotropic interaction between bath gas and the complex at each instant, the net torque on AB from clockwise and counterclockwise directions due to the bath gas would be zero, and isomerization would not take place from ABC to BAC, for example, if the identity of A and B is the same and there exists a perfect cancellation of forces on a given atom. As T → ∞, a dissociative product is to be expected as noted earlier, and since a stable population is absent in a given well, then concomitantly an isomerization event would then be absent. 7.3. N*(EJK) and Transition State(s). The number of states N* at the TS that appears in the terms kr and kd in eqs 34 and 35 and others may be determined using variational RRKM theory,1,9,10 that is, N(EJK) versus R, upon the appropriate consideration of a one- or a two-state TS for the reaction of interest. With the introduction of K as a state variable, providing a sharper distribution in the TS, it becomes increasingly important to introduce a two-state TS form of TS theory. For example, in ozone formation, two TS manifest when K is adiabatic and the total rotational energy is a substantial fraction of the total energy above the asymptote, for example, for K near 3 as region of interest and J = 15, E = 0.2kT above the asymptote. The inner TS is the relatively tighter one with the shorter reaction coordinate R nearer to the unimolecular intermediate, and the outer TS may be an orbiting one described by a longer reaction coordinate, for example, atom−diatom center-to-center distance of R = 8.5 Bohr radius for an outer TS corresponding to an O2···O roaming complex; that is, where a bond becomes greatly extended relative to its equilibrium length and prior to its shortening, the O2 moiety possibly reorient, and a tighter TS located at 5.2 Bohr radius for Ο*3 is determined from variational RRKM theory, informing on the position of the two TS. The dynamical roaming in O3 is confirmed by classical trajectory.14 During the process of Ο3* dissociation, the occurrence of roaming becomes more probable as the excess energy above the dissociation limit is lowered, for example, from 1kT to ∼0.2kT at T = 300 K while J increasing from 1 to 15, where the lifetime of the complex increases to hundreds of picoseconds, and the energy of the roaming complex is predominantly locked in the orbital motion of the extended oxygen atom at the constant J, for example, 10, in the absence of collisions with the bath gas. The unimolecular reaction of C4H8+· is another system whereby two TS have been implicated,90 and its inclusion yielded an improved dissociation rate constant in close agreement with experimental rate constants.90 Roaminginduced product branching by trajectory simulation has also been observed for the decomposition of H2CO, CH3CHO, CH3OOH, and CH3CCH.91 The analysis of reactive flux and potential energy surface has also informed on the separability of tight and roaming pathways in the decomposition of H2CO, CH3CHO, MgH2, and HNNOH.92 The two TS have been posited to be the more general case, rather than an exception,90 where the aforementioned examples of the reactions point in that direction. Thus, the resolution of K degree of freedom, and in the anticipation of any roaming at the outer TS contingent on the (EJK), then renders a two TS treatment of reactive flux relevant for consideration. For considering the N*(EJK) for two TS, it is

given by the expression N*(EJK)eff = N*EJK1 N*EJK2 /[N*EJK1 + N*EJK2 − NEJK *1 NEJK *2 /Nmax *1 and NEJK *2 represent the number of EJK ], where NEJK states at the two minima in a plot of N(EJK) versus the reaction coordinate R, and Nmax EJK represents the maximum in the number of states in between these two minima.12,14,90,93,94

8. SUMMARY AND CONCLUDING REMARKS The analytic solutions, populations, for the bimolecular master equations are derived for a single and multiple potential wells, separately, for the K-adiabatic and K-active cases. The solutions may accommodate various intermolecular energy-transfer models, for example, single-exponential, double-exponential, Gaussian, stepladder, and near singularity. Considerations were given when Z(E′,E) is factorizable, or not, in the derivation of population. On the latter point, the magnitude of the krec based on either way of deriving the population is found to be similar, for example, within 1% for the single-exponential model for ozone formation. The present method extends Part I30 by (1) utilizing an algebraic approach in acquiring the solutions where the implicit assumption35 regarding the (E′J′K′) independence of λ̂ and  in the second term of solution made in Part I30 is now obviated in acquiring the exact solutions and so favorable for utility, and (2) treating the multiple potential wells case. The utilization of these solutions may inform, by inspection of the expressions, on how the population of a given molecule is expected to be modulated as the association rate constant, dissociation rate constant, collision frequency with a given bath gas, and the type of intermolecular energy transfer model vary at a given pressure and temperature in the steady-state regime. The possession of the functional behavior of population may serve for acquiring the steady-state populations for a single or multiple potential well reactions, and the predictions on magnitude may be compared with results from experiments, for example, using pump-dump-probe method for energies of interest,43 to quantify population. A numerical mapping of g(EJK)+ in the (EJK) space and the krec(EJK) based on it elucidate the important regions contributing to the thermal krec. Elsewhere, the expressions are utilized to investigate the population and the recombination rate constant for ozone formation as a function of temperature and pressure for various intermolecular energy-transfer models. At the high-pressure limit, the K-adiabatic and K-active bimolecular master equations each reduce to the high-pressure limiting form of the K-adiabatic and K-active RRKM theory, respectively, for a single potential well and are augmented by terms of isomerization and bimolecular product channels for the case of multiple potential wells. When T → ∞, at the highpressure limit, the analysis shows that a stable molecule is not expected to form in its potential well, but rather the bimolecular reactants lead to dissociative products, while, at the low-pressure limit, the expression for population becomes equivalent to the familiar presumed detailed balance condition between the association and dissociation rate constants regarding a single potential well, where the multiple potential-well case contains additional terms when analyzed separately for the energy regions above and below the dissociation limit. When the collision frequency of the complex with the bath gas is ZLJ ≠ 0 but finite, and ZLJ < kd(E′J′K′) as discussed in the text, the expression for population may further simplify to its first term, reducing to eq 35, while still containing information on collision, in which for O3 covers pressures up to ∼100 bar at T = 300 K, where ∼96% contribution or more to krec comes from the first term of the population via eq 31 or 34 for a single-exponential model, with mostly similar contribution for the other models of N

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

−1 where EJ = J2/2I2, EK = K2(I−1 1 − I2 )/2 for I2 > I1, and the factor of 2 in eq A4 accounts for the degeneracy of K or K′, for its range in the plus or minus direction. The normalization of P canceled in eq A4, as due by contribution from each side of the equation.

intermolecular energy transfer. At higher pressures reaching 1000 bar, postcollision energies become more important in their contribution to krec and are explored elsewhere. The numerical implementation of the expression for the populations and the recombination rate constant based on them is facile, whereby RRKM theory, or other methods of interest, may be utilized for determining the rate constants appearing in the analytic expressions and quantum state counting by using for example a Beyer−Swinehart algorithm.21 Any nonstatistical effects may be included when applicable, such as for the recrossings of the TS by the reaction coordinate from reactant to product, and phase space inaccessibility for the density of states from a bimolecular process, for example, ∼20% and ∼10% for ozone, respectively.14 Resolving the K degree of freedom in a molecular system of interest, and in the flux across the transition state as may be the case, along with its numerical implementation in the analytic approach is also facile and simply entails the inclusion of K in the counting of states in the determination of the K-adiabatic population and recombination constants, and computationally efficient. The analytic-based populations, which include collisional contributions, may also serve for studying isotopic reactions of interest, such as a primary, that is, intermolecular or intramolecular, or a secondary kinetic isotope effect as a function of pressure and temperature, and the predictions of isotope effects on population may be compared with results from experiments. The relative magnitude of the populations for a pair of isotopomers, and contribution to the associated recombination rate constants, may also readily be predicted for energies less or greater than the critical energy of the reaction. An application to isotopic ozone formation, to study its pressure and temperature dependence of enrichment and mass-independent effects, will be presented elsewhere.





Intermolecular Energy Transfer Models for Z(E′,E)

In this section, additional selection of intermolecular energy transfer models, Z(E′,E), are given that may appear in the bimolecular master equation, to augment the single exponential model, eq 33, given in the text. For a double exponential model of energy transfer between a molecule and its bath gas we have1 Z(E′, E) = Zo[exp( −(E − E′)/γ ) + c exp( − (E − E′)/d)]

+ c exp( − (E′ − E)/d′)]

⎧ −(E − E′)2 / γ 2 ⎪ Zoe Z (E ′ , E ) = ⎨ ⎪ Z e−(E′− E)2 / γ′2 ⎩ o

Z(E′J ′K ′ → EJK ) = ωdP(J ′K ′ → JK |E)δ(E′ − E − ΔE) (A1)

and Z(EJK → E′J ′K ′) = ωaP(JK → J ′K ′|E)δ(E − E′ − ΔE)

Z(E′, E) = Zoσ

(A2)

Z(E′J ′K ′ → EJK )ρ(E′J ′K ′)e−E′/ kT =

JJ ′ KK ′

Z(EJK → E′J ′K ′)ρ(EJK )e−E / kT (A3)



For way of demonstration, upon introducing a strong collision for the transfer of rotational energy and inserting eqs A1 and A2 into eq A3 yields



ωdδ(E′ − E − ΔE)ρ(E′J ′K ′)e

(2J + 1)e

× 2e

=



ωaδ(E − E′ + ΔE)ρ(EJK )e

E′ ≥ E

(B3)

(B4)

APPENDIX C

We first proceed to derive the coefficient  , eq 23, which appears in g(EJK),eq 22, and in g(EJK)+, eq 26, in which the  is needed and used in the proof for the solution following it. Substituting g(EJK) into the master equation in the form of Fredholm integral equation yields

−E / kT

JJ ′ KK ′

× (2J ′ + 1)e−E J′/ kT 2e−E K′/ kT

E′ ≤ E

Proofs for g(EJK)+ and g(EJK)−, Eqs 26 and 27, as Solutions of the Bimolecular Master Equation with a Single Potential Well, Eq 19, and Derivation for  , Eq 23

−EJ / kT

JJ ′ KK ′ −EK / kT

(B2)

where Zo is a normalization constant, and σ is a constant when |E − E′| = ⟨ΔE⟩d; otherwise, Z(E′,E) = 0, and the activation transition is given by detailed balance.39−41 An activation and deactivation frequency, ωa and ωd, may be defined, respectively, according to ωd = ωae−E/kT, when the density of states is assumed constant. A complex stepladder model with a distribution of steps may also be considered.39

According to microscopic reversibility, we have

−E′/ kT

E′ ≥ E

where γ and γ′ are again the deactivation and activation constants, which are related to each other by detailed balance. The normalization constant in Zo may be determined numerically. For a stepladder model for intermolecular energy transfer,39−41 we have

We note the energy transfer processes

JJ ′ KK ′

(B1)

where Zo is a constant that contains a normalization constant multiplied by the Lennard-Jones frequency, and γ and γ′ are related by detailed balance and likewise for d and d′. The γ and γ′ are the deactivation and activation constants with units of energy, and similarly so for d and d′. For a Gaussian model of intermolecular energy transfer,37,38

Microscopic Reversibility in the (EJK) Regime



E′ ≤ E

Z(E′, E) = Zo[exp( −(E′ − E)/γ ′)

APPENDIX A



APPENDIX B

(A4) O

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

equation, eq 19, where g(EJK) = g(EJK)+ holds for this region. Upon inserting g(EJK)+, eq 26, into the K-adiabatic master equation, eq 19, yields

f (EJK ) + λ(EJK )Â − λ(EJK ) ∞

∑∫*

Z(E′, E)[f (E′J ′K ′) + λ(E′J ′K ′)Â ]dE′

E

J′K′

= f (EJK )

where the lower limit of the integration is set to E* for the commensurability of presentation in the text when considering f1 and λ̂, and without loss of generality, or the analysis may also proceed with the limit set to Eo, with additional comment following eq C4. Strictly speaking, the domain of E′ is determined by the cases E′ > E or E ≥ E′, and the procession may accommodate either case. Canceling f(EJK) from the right and left hand sides of eq C1, and distributing the

J′K′

J′K′

= f (EJK )



∑∫*

Z (E ′ , E ) × d E ′

term then yields

λ(EJK )Â − λ(EJK ) ∑ J′K′ ∞

∑∫*



∫E*

f1

Z(E′, E)f (E′J ′K ′)dE′−λ(EJK )



(C2)

Upon readily eliminating λ(EJK) in eq C2 and factoring  then yields

where

Z(E′ , E)λ(E′J ′K ′)dE′]

E

J′K′

∑∫*

Z(E′, E)f (E′J ′K ′)dE′

E

J′K′

(C3)

Solving for  in eq C3 then yields  =

∑J ′ K ′ ∫



E*

1 − ∑J ′ K ′ ∫

Z(E′, E)f (E′J ′K ′)dE′ ∞

E*

Z(E′, E)λ(E′J ′K ′)dE′

∑∫*

Z(E′, E)f (E′J ′K ′)dE′

E

J′K′

f1



1−λ

∑∫* ̂ J′K′

Z(E′, E)λ(E′J ′K ′)dE′

E

(C6)

f1 1 − λ̂

is independent of (E′J′K′) as informed by eq C4 and

eqs 23−25, and so was brought outside of the integral and the summation sign in the third term, whereupon insisting that E in Z(E′,E), and in the f1 and λ̂, retain its pre-collision definition. Otherwise, a derivation based on a factorizable Z(E′,E) = Z(E) Z(E′) is given in the Appendix E. Noting that the second term in eq C6 is f1, defined in eq 24, and then factoring the f1 and eliminating it in eq C6 then yields



=



=0



∑∫*



1 − λ̂

Z(E′, E)λ(E′J ′K ′)Â dE′ = 0

E

J′K′

Z (E ′ , E )

E

J′K′

E

 [1 −

(C5)

Upon canceling f(EJK)from both sides of eq C5, canceling λ(EJK), and distributing the

term then yields

×

− λ(EJK ) 1 − λ̂ ∞ ⎡ f ⎤ Z(E′, E)⎢f (E′J ′K ′) + λ(E′J ′K ′) 1 ⎥dE′ E* ⎣ 1 − λ̂ ⎦

∑∫



∑∫*

f1

f (EJK ) + λ(EJK )

(C1)

1 1 −1− ∑ ̂ 1−λ 1 − λ ̂ J′K′

(C4)

The numerator in eq C4 is the f1 defined in eq 24, and



∫E*

Z(E′, E)λ(E′J ′K ′)dE′ = 0 (C7)



∑∫* J′K′

Z(E′, E) × λ(E′J ′K ′)dE′

Also, upon noting that

E

is λ̂ defined in eq 25, where both substitutions into eq C4 then yield eq 23, Â =

λ̂ =



∑∫* J′K′

f1

Z(E′, E) × λ(E′J ′K ′) × dE′

E

from eq 25, then eq C7 becomes

[1 − λ ]̂

1 λ̂ =0 −1− 1 − λ̂ 1 − λ̂

The lower limit of the integrals in eqs C1−C4 could have been set to Eo and so be reflected in  . As noted earlier, more precisely, the domain of E′ is determined by the cases E′ > E or E ≥ E′, and the procession may accommodate either case, for each energy region in the integrals of eqs C1−C4. Thus, g(EJK), eq 22, is set up as a solution of the bimolecular master equation with its coefficient  , and next we prove it for E > E*(JK), where g(EJK) = g(EJK)+, where g(EJK)+ = 0 for E ≤ E*(JK), since  and f(EJK) are each zero in eq 22 or 26 for the latter energy region. A similar analysis may be used regarding g(EJ) as a solution of its K-active master equation, eq 14, or of the analogue of eq 19. We now prove, by substitution, that g(EJK)+, eq 26, for the energy region E > E*(JK), is a solution of the K-adiabatic master

(C8)

which upon rearrangement yields 1 − λ̂ =1 1 − λ̂

that then completes the proof. An analogous analysis may prove that the K-active g(EJ)+ is a solution of its master equation. We next prove that, for the energy region of E ≤ E*(JK) g(EJK)−, eq 27, is a solution of its K-adiabatic master equation, eq 19. Upon inserting g(EJK)− = λ(EJK)−Ĉ , eq 27, into the lefthand side of K-adiabatic master equation, eq 19, that is, the first term and inside the integral sign of the second term, yields P

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A λ(EJK )− Ĉ −

λ(EJK )− ∑ Z LJ J ′ K ′

∫E

E*

fi (EJK ) + λi(EJK )Â i − λi(EJK ) ∑

Z(E′, E)λ(E′J ′K ′)− Ĉ dE′

J′K′

0

Since E ≤ E*(JK), the λ(EJK) positioned outside the integral in the master equation, eq 19, becomes λ(EJK)−/ZLJ, since kd(EJK) = 0 in λ(EJK), eq 21, leaving λ(EJK)− and the ZLJ in the denominator in eq C9. Also, the f(EJK) that would have appeared on the right-hand side of the master equation is set to zero in eq C9, since the association and dissociation rate constants in the f(EJK) are each zero below the critical energy limit. We next substitute λ(EJK)−Ĉ = P(JK)geq(0)cEe−E/kT, into the master eq 9, where cE = δ(E − E′)c′E, and c′E is the contribution from the collisional energy transfer to the populations as previously determined,16,42 for example, for a single exponential model for a collision transition probability

λi(EJK )Â i − λi(EJK ) ∑ J′K′ ∞

∫E*

−λi(EJK ) ∑ J′K′

∫E



∫E*

Zi(E′, E)fi (E′J ′K ′)dE′

Zi(E′, E)λi(E′J ′K ′)Â i dE′ = 0 (D2)

Upon readily eliminating λi(EJK) in eq D2 and factoring  i then yields

was obtained using a Wiener−Hopf method,42 where γ′ and γ are the activation and deactivation constants, respectively, with units of energy, to yield P(JK )geq (0)c Ee−E / kT

(D1)

where the lower limit of the integration is set to E* for the convenience of the presentation rendered for f1i and λî in the text, or the analysis may also proceed with Eo, with additional comment following eq D4. Strictly speaking, the domain of E′ is determined by the cases E′ > E or E ≥ E′, and the procession may accommodate either case. Canceling f i(EJK) from the rightand left-hand sides of eq D1 and distributing the terms inside the integral sign then yields

⎞ ⎛ γ′ c′E = ⎜1 − e E / kT ⎟ γ ⎠ ⎝

P(JK ) ∑ − Z LJ J ′ K ′

Zi(E′, E)

× [fi (E′J ′K ′) + λi(E′J ′K ′)Â i ]dE′ = fi (EJK )

(C9)

=0



∫E*

⎡ Â i ⎢1 − ⎢⎣

E*



∑∫* E

J′K′

⎤ Zi(E′, E)λi(E′J ′K ′)dE′⎥ ⎥⎦



Z(E′, E)P(J ′K ′)

=

0

∑∫* J′K′

× geq (0)c′E δ(E − E′)e−E′/ kT dE′ = 0

E

Zi(E′, E)fi (E′J ′K ′)dE′

(D3)

Solving for  i in eq D3 then yields

In the latter equation P(JK) may be factored and cancel from the equation, and

 i =

∑J ′ K ′ ∫



E*

1 − ∑J ′ K ′ ∫

Zi(E′, E)fi (E′J ′K ′)dE′ ∞

E*

∑ P(J′K ′) = 1

Zi(E′, E)λi(E′J ′K ′)dE′

(D4)

The numerator in eq D4 is the f1i defined in eq 44, and

J′K′



∑∫*

according to normalization. Furthermore, integrating over E′ in the second term yields ZLJgeq(0)c′Ee−E/kT, which then the master equation simplifies to geq(0)c′Eδ(E−E′)e−E/kT−geq(0)c′Ee−E/kT = 0. Canceling terms in the latter equation leaves δ(E − E′) − 1 = 0, which is true when E = E′, as assumed at the outset, which then completes the proof that the solution satisfies the master equation. As noted above, although δ(E − E′) appeared in the integral as a device in acquiring the solution for this energy region, information about the collision energy transfer is retained in the solution as it appears in c′E, as calculated elsewhere for various types of energy transfer models.16,42 An analogous analysis may prove that the K-active g(EJ)− is a solution of its master equation.



J′K′

E

Zi(E′, E) × λi(E′J ′K ′) × dE′

is λî , eq 45, where substituting both into eq D4 then yields eq 43, Â i =

f1

i

[1 − λî ]

The lower limit of the integrals in eqs D1−D4 could have been set to Eo and so be reflected in  i. As noted earlier, more precisely, the domain of E′ is determined by the cases E′ > E or E ≥ E′, and the procession may accommodate either case, for each energy region in the integrals of eqs C1−C4. Thus, gi(EJK), eq 42, is set up as solution to the bimolecular master equation with its coefficient  i, and next we prove it for E > E*(JK). A similar analysis may be used regarding gi(EJ) as a solution for its K-active master equation. We prove, by substitution, that gi(EJK)+, eq 46, for the energy region E > E*(JK), is a solution of the K-adiabatic master equation, eq 39, where gi(EJK) = gi(EJK)+ holds for this region. Upon inserting gi(EJK)+, eq 46, into the K-adiabatic master equation, eq 39, yields

APPENDIX D

Proofs for gi(EJK)+ and gi(EJK)−, Eqs 46 and 47, as Solutions of the Master Equation with Multiple Potential Wells, Eq 39, and Derivation for  i, Eq 43

We retain the strategy employed for the single-well case and first proceed to derive the coefficient  i, eq 43, for the multiple-well case that appears in gi(EJK), eq 42, and in gi(EJK)+, eq 46, in which it is needed and used in the proof for the solution following it. Hence, we consider the bimolecular master equation in the form of the Fredholm integral equation, eq 39, and substitute gi(EJK) into its left-hand side and also into its right-hand side, inside the integral sign, to yield

fi (EJK ) + λi(EJK )Â i − λi(EJK ) ∑ J′K′



∫E*

Zi(E′, E)

× [fi (E′J ′K ′) + λi(E′J ′K ′)Â i ]dE′ = fi (EJK ) Q

(D5)

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

region, E ≥ E′ or E′ > E, for defining the domain of E′. The  i−, eq 48, is algebraically derived from eq D9 in a similar manner as delineated for  i, eq D4. Upon substituting for  i−, eq 48, into eq D9, canceling f i(EJK)− from both sides, canceling λi(EJK)−, and distributing the

Upon substituting  i, eq D4 or 43, into eq D5, canceling f i(EJK) from both sides of eq D5, canceling λi(EJK), and distributing the ∞

∑∫*

Zi(E′, E)

E

J′K′

f1

1 − λî −



∑∫* E

J′K′

f1

i

Zi(E′, E)fi (E′J ′K ′)dE′

∑∫* ̂ J′K′

E

where

i

1 − λî

f 1− i

Zi(E′, E)λi(E′J ′K ′)dE′

1 − λî (D6)

=0 f1





∑∫

Noting that

E*

i J′K′

f 1− i

− 1 − λî

dE ′ = 0

f1− i 1 − λî



is independent of (E′J′K′), as informed by

∫E

E*

Zi(E′, E)λi(E′J ′K ′)− dE′

o

(D11)

=0

Also, upon noting that

Also, upon noting that ∞

E

Zi(E′, E)λi(E′J ′K ′)−

Eo

1 1 − − 1 − −∑ ̂ 1 − λi 1 − λî J ′ K ′

(D7)

J′K′

Zi(E′, E)fi (E′J ′K ′)− dE′

eqs 49 and 50, it can be brought outside of the integral in the third term, and the second term in eq D10 is f1i− as defined in eq 49, and so then factoring the f1i− and eliminating it in eq D10 then yields



∑∫*

Eo

J′K′

∑∫

E*

(D10)

44, and 45, and so was brought outside both of the integral and the summation sign in the third term of eq D6, whereupon insisting that E in Z(E′,E), in the f1i and  i, retain its pre-collision definition. Otherwise, a derivation based on a factorizable Zi(E′,E) = Zi(E)Zi(E′) may be considered, similar to the singlewell case (Appendix E). Noting that the second term in eq D6 is f1i, defined in eq 44, and then factoring the f1i and eliminating it in eq D6 then yields ∞ 1 1 −1− Zi(E′, E)λi(E′J ′K ′)dE′ = 0 ∑ E* 1 − λ̂ 1 − λ̂

λî =



J′K′

is independent of (E′J′K′) as informed by eqs D4,

i

Zi(E′, E)

term then yields



1 − λi

Eo

J′K′



i

E*

∑∫

term, then yields

Zi(E′, E) × λi(E′J ′K ′) × dE′



∑∫

λî =

J′K′

from eq 45, then eq D7 becomes

E*

Eo

Zi(E′, E) × λi(E′J ′K ′)− × dE′

from eq 50, then eq D11 becomes

λî 1 =0 −1− ̂ 1 − λi 1 − λî

− λî 1 − − 1 − − = 0 1 − λî 1 − λî

(D8)

which upon rearrangement yields

(D12)

which upon rearrangement yields

1 − λî =1 1 − λ̂

− 1 − λî − = 1 1 − λ̂

i

and then completes the proof. An analogous analysis may prove that the K-active gi(EJ)+ is a solution of its bimolecular master equation. We next prove that, for the energy region of E ≤ E*(JK) gi(EJK)−, eq 47, is a solution of its K-adiabatic master equation, eq 39, where gi(EJK) = gi(EJK)− holds for this region, and the f i(EJK)− and λi(EJK)−, eqs 51 and 52, respectively, appear in lieu of their counterparts f i(EJK) and λi(EJK) in the master equation. Upon inserting gi(EJK)−, eq 47, into the K-adiabatic master equation, eq 39, onto its left-hand side and the righthand side inside the integral sign, yields fi (EJK )− + λi(EJK )− Â i− − λi(EJK )− ∑ J′K′

∫E

i

and then completes the proof. An analogous analysis may prove that the K-active gi(EJ)− is a solution of its bimolecular master equation, analogue of eq 39.



Analytic Solution g(EJK)+ of the Bimolecular Master Equation, Eq 19, When Z(E′,E) = Z(E)Z(E′) and Proof

We proceed to first derive the analytic solution of the bimolecular master equation, eq 19, for the case when Z(E′,E) is factorized as Z(E)Z(E′) for the energy region of E > E*(JK) and consider for it the single-exponential model. Thereafter, the proof of the solution is presented. The solution for the region of E ≤ E*(JK), eq 27, and its proof remain as before, presented in Appendix C, where an expounded comment is made at the end of the section. We write the bimolecular master equation, eq 19, as

E*

Zi(E′, E)

o

× [fi (E′J ′K ′)− + λi(E′J ′K ′)− Â i− ]dE′ = fi (EJK )−

APPENDIX E

(D9)

Note the upper limit of the integral in eq D9 is E* for this energy region. The procession may accommodate either energy R

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A g (EJK ) − λ ̅ (EJK ) ∑ J′K′

∫E



λ̂ =

Z(E′)g (E′J ′K ′)dE′ = f (EJK )

o

 =

f1 [1 − λ ]̂

For a single exponential model, when E > E′, then we have

λ ̅ (EJK ) = Z(E)λ(EJK )

(E2)

E

f1 =

appearing outside of the integral in eq E1, where for a singleexponential model, we have Z(E) = e−E/γ and Z(E′) = eE′/γ when

g (EJK )+ = f (EJK ) + λ ̅ (EJK )Â

λ̂ =

∫E

E

for a Z(E′) = eE′/γ. In the integrand of  , the eE′/γ end up canceling the e−E′/γ in λ̅(E′J′K′). We next prove, by substitution, that g(EJK)+, eq E3, for the energy region of E > E*(JK), is a solution of the K-adiabatic master equation, eq E1. Upon inserting eq E3 into the Kadiabatic master equation, eq E1, then yields f (EJK ) + λ ̅ (EJK )



Z(E′)

f1 1 − λ̂

− λ ̅ (EJK ) ∑ J′K′

(E4)

Canceling f(EJK) from both sides of the equation, and distributing the E



1 − λ̂

λ ̅ (EJK )Â − λ ̅ (EJK ) ∑ J′K′ ∞

∑∫*





∫E*

Z(E′)f (E′J ′K ′)dE′−λ ̅ (EJK )



where

(E5)

J′K′



1−λ



∑∫* ̂ J′K′

Z(E′)λ (̂ E′J ′K ′)dE′

E

f1

(E11)

1 − λ̂

is independent of (E′J′K′) as informed by eqs E7 to

1 1 −1− ∑ ̂ 1−λ 1 − λ ̂ J′K′

Z(E′)λ ̅ (E′J ′K ′)dE′]

E

J′K′

f1

Z(E′)f (E′J ′K ′)dE′

E

E9, and so was brought outside of the integral in the third term. Noting that the second term in eq E11 is f1, defined in eq E8, and then factoring the f1 and eliminating it in eq E11 then yields

Upon readily eliminating λ̅(EJK) in eq E5 and factoring  then yields

∑∫*

∑∫*

=0

Z(E′)λ ̅ (E′J ′K ′)Â dE′ = 0

E

J′K′



∫E*

Z(E′)λ ̅ (E′J ′K ′)dE′ = 0



∑∫* J′K′

(E12)

Z(E′)f (E′J ′K ′)dE′

E

Also, upon noting that

(E6)



Solving for  in eq E6 yields ∑J ′ K ′ ∫



E*

1 − ∑J ′ K ′ ∫

λ̅ =

Z(E′)λ ̅ (E′J ′K ′)dE′

λ̂ 1 −1− =0 ̂ 1−λ 1 − λ̂



J′K′

E

Z(E′) × λ ̅ (E′J ′K ′) × dE′

E

from eq E9, then eq E12 becomes (E7)

where we define

∑∫*

∑∫* J′K′

Z(E′)f (E′J ′K ′)dE′ ∞

E*

f1 =

(E10)

E

f1

Z(E′) × dE′

then yields

 =

Z(E′)

Upon canceling f(EJK)from both sides of eq E10, canceling ∞ λ̅(EJK), and distributing the ∑J ′ K ′ ∫ * Z(E′) term then yields



∑∫*

=



∫E*

⎡ f ⎤ ⎢f (E′J ′K ′) + λ ̅ (E′J ′K ′) 1 ⎥dE′ = f (EJK ) ⎣ 1 − λ̂ ⎦

o

× [f (E′J ′K ′) + λ ̅ (E′J ′K ′)Â ]dE′ = f (EJK )

E

∑ ∫ * e E′/γ′ × λ ̅ (E′J′K ′) × dE′ J′K′

for the energy region of E > E*(JK). We derive an expression for  in eq E3 by inserting eq E3 into eq E1 and displaying the steps. Inserting eq E3 into E1 yields

J′K′

E

and

(E3)

f (EJK ) + λ ̅ (EJK )Â − λ ̅ (EJK ) ∑

∑ ∫ * e E′/γ′ × f (E′J′K ′) × dE′ J′K′

E > E′, and Z(E) = eE/γ′ and Z(E′) = e−E′/γ′ when E′ > E. For convenience the Z(E′) is assigned to carry the normalization and the collision frequency. The procession may apply for either when E > E′ or E′ > E, separately, where each defines the domain of E′ in the integral of eq E1. The f(EJK) in eq E1 is defined by eq 20 in the main text. We again seek a solution of eq E1 in the form of

 [1 −

(E9)

and so eq E7 may compactly be written as

(E1)

where after factoring Z(E,E′) = Z(E)Z(E′), and combining the Z(E) with the λ(EJK) term, eq 21, outside adjacent to the integral in eq 19, we then have

×

Z(E′)λ ̅ (E′J ′K ′)dE′

E

J′K′

Eo ≤ E ≤ ∞

J′K′



∑∫*

(E13)

which upon rearrangement yields

Z(E′)f (E′J ′K ′)dE′

1 − λ̂ =1 1 − λ̂

(E8)

and S

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A and then completes the proof. An analogous analysis may prove that the K-active g(EJ)+ is a solution of its master equation. When E′ > E, then as noted earlier, for a single exponential

g (EJK )+ =

×

model we have Z(E) = eE/γ′ and Z(E′) = e−E′/γ′, and these would appear in lieu of their E > E′ counterpart in the derivation. The same proof of the solution applies. The g(EJK)+, eq E3, for a factorized Z(E′,E) = Z(E)Z(E′) then appears in full form as

k r(EJK )A·BC + [Z LJ + kd(EJK )]−1 Z LJ + kd(EJK )

Z(E)P(JK )∑J ′ K ′ ∫



E*

1 − ∑J ′ K ′ ∫



E*

f (E′J ′K ′)Z(E′)dE′

P(J ′K ′)[Z LJ + kd(E′J ′K ′)]−1 dE′ (E14)

For the single-exponential model, eq E1, with a factorized Z(E′,E), the g(EJK)+ defined by eq E3 or E14, in full form then appears as ∞

( )

( )

e∓E / γ ′ P(JK )Zo∑J ′ K ′ ∫ * e±E′/ γ ′ k r(E′J ′K ′)A·BC[Z LJ + kd(E′J ′K ′)]−1 dE′ k r(EJK )A·BC E g (EJK )+ = + ∞ Z LJ + kd(EJK ) [Z LJ + kd(EJK )][1 − Zo∑J ′ K ′ ∫ * P(J ′K ′)[Z LJ + kd(E′J ′K ′)]−1 dE′]

(E15)

E

( )

(′)

where the −E/γ and +E′/γ are used in e∓E/ γ ′ and e±E′/ γ , respectively, for the region E > E′, or use +E/γ′ and −E′/γ′ for E′ > E, in elucidating the usage of the notations γ(′), ∓, and ± in eq E15, upon the consideration of the domains of E′. Each domain also determines the limits of the integrals in eq E15. The solution for the region energy E ≤ E*(JK), as expressed by eq 27, and its proof remain as before, presented in Appendix C, since the factorization of Z(E′,E) does not alter the premise and the course of the approach partaken.



*



|gi(EJK )+ |
v′, |gi(zJK)+| has no singularity. So, gi(zJK)+ is an analytic function in the half plane of the z-plane for which Imz > v′. Now upon considering gi(EJK)−: we Fourier transform ∞

|gi(zJK )− | =

APPENDIX F

∫−∞ e2πizE × gi(EJK )−dE

Again letting z = u + iv, we have

Boundary Conditions for gi(EJK)+ and gi(EJK)− for a Single and Multiple Potential Wells and Analytic Continuation

E*

|gi(zJK )− | =

Herein, we proceed to elucidate the boundary conditions of populations for multiple potential wells whereupon gi(EJK)+ = gi(EJK)− for the K-adiabatic case, by using complex analysis to check for singularity and convergence of each function for each energy domain E > E*(JK) and E ≤ E*(JK) and analytic continuation across the boundary.95,96 A similar analysis may be applied to a K-active case. Thereafter, the boundary condition for a single potential well case is considered. We first consider gi(EJK)+ for the energy region of E > E*(JK). Upon applying a Fourier transform,

∫−∞ e2πiuE × e−2πvEgi(EJK )−dE E*



∫−∞ e−2πvE ×

gi(EJK )− dE

We seek a solution such that |gi(zJK )− | < x 2ie 2πv ″ E

as E → −∞ where v″ > 0, and hence gi(EJK)− → 0 as E → −∞, and thus we have



f (z ) =

∫−∞ e2πizE × f (E)dE

E*

|gi(zJK )− |
E*(JK )

If a solution exists for gi(zJK)+ such that

v′ − v < 0

(F3)

and

|gi(EJK )+ | < x1ie 2πv ′ E

gi(EJK )− = x 2ie 2πEv ″|Ei[−2π (E − E*)(v″ − v)]| /2π

as E → ∞ where v′ < 0, so |gi(EJK)+| → 0 as E → ∞, then we have

E ≤ E*(JK ) T

v″ − v > 0

(F4) DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A



The units of v, v′, and v″ is inverse energy. The Ei denotes an exponential integral95,97 in eqs F3 and F4, defined as ∞

Ei(z) =

∫−z

Article

APPENDIX G N

Analytic Expression for ∑q ≠ i kqi(EJK )gq(EJK ) in Eq 49

We proceed to derive an analytic expression for

e −t / t dt

N

∑ kqi(EJK ) × gq(EJK )

where the principal value of the integral is taken. The inverse transform was considered under the restriction of the conditions noted in eqs F3 and F4. At the boundary for a given (EJK), we have gi(EJK)+ = gi(EJK)−, upon using eqs F3 and F4, with the constants x1i and x2i selected to ensure the equality, and obeying the conditions v′ − v < 0 and v″ − v > 0 for the two energy regions of E. Inserting gi(EJK)+ and gi(EJK)− from eqs F3 and F4 into the boundary condition eq 53 then establishes the relationship between x1i, x2i, f i(EJK), f i(EJK)−, λi(EJK), λi(EJK)−, λî , and λ−î . The aforementioned analysis also applies to the K-active boundary condition, whereby now the K degree of freedom has appropriately been averaged over at the outset and would be absent in each step of the above presentation. A similar analysis may be performed for elucidating the boundary conditions for a single potential well, where now the subindex i would be dropped in each step of the analysis and in eqs F1−F4,30 without the subindex i appearing there, would be applicable. Then, inserting g(EJK)+ and g(EJK)− from the analogue of eqs F3 and F4 into the boundary condition for the single well, eq 30, then establishes the relationship between x1, x2, f(EJK), λ(EJK), λ1, and  .98

q≠i

found in eq 49, for the energy region of E ≤ E*(JK). We utilize the boundary condition gi(EJK)− = geqi(EJK) when E → −Di (Di is the dissociation energy of ABC molecule measured from the bottom of its potential well i), where geqi(EJK) denotes the equilibrium probability density that ABC has energy E, defined by geqi(EJK)− = ρi(EJK)e−E/kT/Q′i, where ρi(EJK) is the density of states of the molecule in well i, and Q′i is the partition function of the molecule in the center-of-mass system of coordinates in potential well i. Corrections to the equilibrium density have also been derived for various intermolecular energy transfer models16,42 and may be included as c′Ei for potential well i. Rewriting geqi(EJK)− as geqi(0)e−E/kT, where geqi(0) = ρi(EJK)/ Q′i, and substituting into the left-hand side of eq 47 yields geq (0)c′Ei e i

−E / kT

= fi (EJK )− + λi(EJK )−

f 1− i

− 1 − λî

(G1)

Upon substituting the definitions for f1i−, λ−î , and f i(EJK)−, eqs 49 50, and 51, respectively, into the right-hand side of eq G1 yields

⎡N ⎤ geq (0)c′Ei e−E / kT = λi(EJK )− Pi(JK )−1⎢∑ kqi(EJK )gq(EJK )⎥ + λi(EJK )− ⎢⎣ q ≠ i ⎥⎦ i E* ⎡ N ⎤ ∑J ′ K ′ ∫ Zi(E′, E)λi(E′J ′K ′)− Pi(J ′K ′)−1⎣∑q ≠ i kqi(E′J ′K ′)gq(E′J ′K ′)⎦dE′ Eo 1 − ∑J ′ K ′ ∫

E*

Eo

Zi(E′, E)λi(E′J ′K ′)− dE′

(G2)

Upon dividing eq G2 by λi(EJK)− and solving for N

∑ kqi(EJK ) × gq(EJK ) q≠i

yields N

∑ kqi(EJK )gq(EJK ) =

Pi(JK )geq (0)c′Ei e−E / kT i

λi(EJK )−

q≠i



Pi(JK )∑J ′ K ′ ∫

E*

Eo

N

Zi(E′, E)λi(E′J ′K ′)− Pi(J ′K ′)−1[∑q ≠ i kqi(E′J ′K ′)gq(E′J ′K ′)]dE′ 1 − ∑J ′ K ′ ∫

E*

Eo

Zi(E′, E)λi(E′J ′K ′)− dE′

(G3)

well q in the center-of-mass system of coordinates. As E′ → E*, anharmonic contributions to density of state may also be taken into account. Rewriting geqq(E′J′K′)− as geqq(0)e−E′/kT, where geqq(0) = ρq(E′J′K′)/Q′q, and including c′E′q as the correction to the equilibrium density for potential well q, then substituting into the gq(E′J′K′) of the right-hand side of eq G3 yields

Also applying the boundary condition, g q(E′J′K′) = geqq(E′J′K′) when E′ → −Dq (Dq is the energy from the bottom of potential well q to the lowest energy transition state associated with potential well q), then we have geqq(E′J′K′)− = ρq(E′J′K′)e−E′/kT/Q′q where ρq(E′J′K′) is the density of states of the isomer in well q, Q′q is the partition function of the isomer in

U

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A Pi(JK )geq (0)c′Ei e−E / kT

N

∑ kqi(EJK )gq(EJK ) =

i

λi(EJK )−

q≠i

Pi(JK )∑J ′ K ′ ∫

E*

Eo



⎤ ⎡ N Zi(E′, E)λi(E′J ′K ′)− Pi(J ′K ′)−1⎢∑q ≠ i kqi(E′J ′K ′)geq (0)c′E′q e−E′/ kT ⎥dE′ ⎦ ⎣ q 1 − ∑J ′ K ′ ∫

E*

Eo



∫E

APPENDIX H



∫E

∫E



Z(E , E′) × dE′

in the numerator cancels

∫E



Z′(E , E′) × dE′

o

in the denominator, leaving a ZLJ in the numerator and upon the latter’s reabsorption with the remaining part of the first term (i.e.,



Z(E′, E)g (E′J ′K ′)+ dE′ (H1)

and first consider the energy region of E > E*(JK), where the lower limit of integration is set to E* for this energy region, where more strictly, the domain of E′ is defined by E ≥ E′ or E′ > E. Next, we take the high-pressure limit of eq H1 as ZLJ → ∞. Noting the definition of f(EJK), eq 20, we have

Z LJ ∑ J′K′

P(JK ) ∑ J′K′

=0

The application of l’Hôpital’s rule (infinity-over-infinity case) to the second term of the right-hand side eq H1 yields g (EJK )+ =



E* ∞

Z′(E′, E)g (E′J ′K ′)+ dE′

∫E Z′(E , E′)dE′ o

(H2)

upon noting the ZLJ in

∫E



∫E



Z(E , E′)dE′ = Z LJ

o

Z′(E , E′) × dE′

that appears as the first term in λ(EJK), eq 21. Then inserting g(EJK)+, eq H2, into the first term of the right-hand side of the master equation, eq 13, yields a master equation k r+(EJK )A·BC =

∑J ′ K ′ ∫

E*

×

∫E

k r+(EJK ) = N *(EJK )e−E / kT /hQ

Z′(E′, E)g (E′J ′K ′)+ dE′

∫E Z′(E , E′)dE′

o

+ kd(EJK )g (EJK )+

∑∫* J′K′

E

(H5)

(H6)

per reactant concentrations, for the energy region of E > E*(JK), whereupon noting eq H5, and as noted earlier that the first and second terms of eq H3 canceled each other. The K-active version of eq H6 may also be obtained by averaging over the Kdegree of freedom in the derivation. The units of k+r (EJK), in eq H6, is (concentration)−1 (time)−1 per (EJK). The unit arises by noting that the translational partition function defined per unit volume contributes a per unit concentration and that h has units of (energy·time) for a given (JK). For the energy region of E ≤





Z(E , E′)dE′ −

(H4)

Thus, the K-adiabatic master equation, eq H3, at the highpressure limit, reduces to

o



Z(E′, E)g (E′J ′K ′)+ dE′

kd(EJK )g (EJK )+ = N *(EJK )e−E / kT /hQ

o





∫E*

Therefore, the first and second terms in the right-hand side of the master equation in eq H3, or in eq 13, cancel, which leaves the remaining third term in eq H3, or in eq 13, for further consideration. We next consider the remaining term, kd(EJK)g(EJK)+ in eq H3. Upon noting the definition, kd(EJK) = N*(EJK)/hρ(EJK) from RRKM theory1,9,10,15,23,24 and assuming g(EJK)+ = ρ(EJK) e−E/kT/Q, or otherwise introduce a correction factor to the density of states, where Q is the partition function of the reactants in the center-of-mass system of coordinates, then we have for the last term of eq H3

86,87,99

P(JK )∑J ′ K ′ ∫

Z(E′, E) × g (E′J ′K ′)+ dE ′)

E

serves to reduce the first term of the right-hand side of eq H3 to become

since lim [Z LJ + kd(EJK )]

Z′(E′, E) × g (E′J ′K ′)+ dE′

∑∫* J′K′

lim f (EJK ) = k r(EJK )[Z LJ + kd(EJK )]−1 A·BC = 0

Z LJ →∞



∫E*



=

Z LJ →∞

−1

Z′(E , E′) × dE′

o

o

g (EJK )+ = f (EJK ) + λ(EJK ) E



then in the first term of eq H3, the

H.1. Single Potential Well. We investigate the high-pressure limit K-adiabatic association rate constant of the bimolecular master equation. The results on the K-active case follow by analogy. We commence by rearranging eq 19 for g(EJK)+

∑∫*

(G4)

Z(E , E′) × dE′ = Z LJ

o

Derivation for the Bimolecular Master Equations Reducing to the High-Pressure Limit Form of RRKM Theory, Eqs 58 and 59, with Isomerization Processes Eqs 60 and 61

J′K′

Zi(E′, E)λi(E′J ′K ′)− dE′

Z(E′, E)g (E′J ′K ′)+ dE′ (H3)

Noting that V

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A E*(JK), k−r (EJK) = 0, since the association channel is closed including at the high-pressure limit. In Part I,30 the K-adiabatic and K-active bimolecular master equations each reduced exactly to the high-pressure limiting forms of RRKM theory at the high-pressure limit, upon using an analytic solution based on an assumption involving E′, delineated in ref 35 of the present work, in its derivation. When omitting the latter assumption in the present work, the master equations again reduce to the high-pressure limit form of RRKM theory for the association rate constant at the highpressure limit. H.2. Multiple Potential Well. We investigate the highpressure limit association rate constants of bimolecular master equations with multiple potential wells, eq 17, in which we provide a derivation that at the high-pressure limit, its Kadiabatic kri(EJK), may reduce, in form, to the rate constant based on the high-pressure limit form of RRKM theory and additional terms describing isomerization, and by extension a similar derivation is valid for the K-active case, where the K degree of freedom is averaged over. We first consider the energy region, E > E*(JK), and commence by demonstrating that the first two terms in the right-hand side of eq 17 cancel at the highpressure limit, and in turn the third term emerges as the bimolecular rate constant of the high-pressure limit of RRKM theory, with the remaining terms of isomerization associated with the multiple potential wells and the bimolecular products. We commence by rearranging the bimolecular master equation, eq 39, for gi(EJK)+

where the lower limit of integration is set to E* for this energy region. Specifically, the domain of E′ is defined by according to whether E′ > E or E ≥ E′. We take the high-pressure limit of eq H7, as ZLJi → ∞. Using the definition of f i(EJK), eq 40, we have lim f (EJK ) = λi(EJK )[∑ kqi(EJK ) × gq (EJK )

Z LJ →∞ i

+ k ri(EJK )A·BC] = 0

since lim λi(EJK ) = 0

Z LJ →∞ i

upon noting the ZLJi in

∫E

∫E

i



Z′i (E , E′) × dE′

o



Pi(JK )∑J ′ K ′ ∫

E* ∞

gi(EJK )+ =

Z′i (E′, E)gi(E′J ′K ′)+ dE′

∫E Z′i (E , E′)dE′ o

(H8)

Upon inserting gi(EJK)+, eq H8, into the first term of the right-hand side of the master equation, eq 17, yields a master equation

Zi(E′, E)gi(E′J ′K ′)+ dE′

J′K′ E*

(H7)

Pi(JK )∑J ′ K ′ ∫



E* ∞

=

Zi(E , E′)dE′ = Z LJ

that appears as the first term in λi(EJK), eq 41. The application of l’Hôpital’s rule86,87 (infinity-over-infinity case) to the second term of the right-hand side of eq H7 yields



k r+i (EJK )A·BC



o

gi(EJK )+ = fi (EJK ) + λi(EJK )

∑∫

q≠i

i

Z′i (E′, E)gi(E′J ′K ′)+ dE′

∫E Z′i (E , E′)dE′

∫E



Zi(E , E′)dE′−Pi(JK ) ∑

o

J′K′



∫E*

Zi(E′, E)gi(E′J ′K ′)+ dE′

o

N

np

N

∑ kiq(EJK )gi(EJK )+ −∑ kqi(EJK )gq(EJK )+ + k d (EJK )gi(EJK )+ + ∑ k p(EJK )gi(EJK )+

+

i

q≠i

q≠i

p=1

Z LJ ∑

Noting that

i

∫E



Zi(E , E′)dE′ = Z LJ

∫E

i

o



=

o



∫E*

Z′i (E′, E) × gi(E′J ′K ′)+ dE′

∑∫* J′K′

Zi(E′, E) × gi(E′J ′K ′)+ dE ′)

E

serves to reduce the first term of eq H9 to become



Pi(JK ) ∑

Zi(E , E′) × dE′

o

J′K′



∫E*

Zi(E′, E)gi(E′J ′K ′)+ dE′

(H10)

Therefore, the first and second terms in the right-hand side of the master equation in eq H9, or in eq 17, cancel, which leaves the remaining terms in eq H9 for further consideration. We next consider the fifth term, kdi(EJK)gi(EJK)+, in eq H9. Upon noting the definition, kdi(EJK) = Ni*(EJK)/hρi(EJK) from RRKM theory and assuming gi(EJK)+ = ρi(EJK)e−E/kT/Qi, or otherwise introduce a correction factor to the density of states, then we have for the fifth term of eq H9

in the numerator cancels

∫E

(H9)



Z′i (E , E′) × dE′

then in the first term of the right-hand side of eq H9, the

∫E

J′K′

i



Z′i (E , E′) × dE′

o

in the denominator, leaving a ZLJi in the numerator and upon the latter’s absorption with the remaining part of the first term (i.e., W

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A k di(EJK )gi(EJK )+ = Ni*(EJK )e−E / kT /hQ i

where the upper limit in the integration is set to E* for this energy region. Next, we take the high-pressure limit of eq H14, as ZLJi → ∞. Using the definition of f i(EJK)−, eq 51, we have

(H11)

The Qi is the partition function of the reactant in the centerof-mass system of coordinates, for the ith potential well. For the sixth term of eq H9,

lim f (EJK )− = λi(EJK )− [∑ kqi(EJK ) × gq(EJK )] = 0

Z LJ →∞ i

q≠i

i

np

∑ k p(EJK ) × gi(EJK )+

since

i

p=1

lim λi(EJK )− = 0

Z LJ →∞

we have kpi(EJK) = N*pi (EJK)/hρpi(EJK) and ρpi = ρi, since the product channels pi share the same density of states with well i. Also, again assuming the equilibrium definition for gi(EJK)+, or introducing a correction factor for it when applicable, then np

np

∑ k p(EJK ) × gi(EJK )+ = ∑ p=1

i

p=1

N *p (EJK ) i

hρi (EJK )

i

upon noting the ZLJi in

∫E



Zi(E , E′)dE′ = Z LJ

i

o

to yield np

∫E Z′i (E , E′)dE′ (H15)

(H12)

p=1

Upon inserting gi(EJK)−, eq H15, into the first term of the right-hand side of the master equation, eq 17, yields a master equation

for the sixth term of eq H9, which is the high-pressure limiting form of RRKM theory for the set of bimolecular product channels, pi associated with potential well i. Thus, the Kadiabatic bimolecular master equation, eq H9 or 17, at the highpressure limit, reduces to

0=

Pi(JK )∑J ′ K ′ ∫

i

×

∫E Zi′(E , E′)dE′

∑ kiq(EJK )gi(EJK )+

×

q≠i

(H13)

+

Eo

N

Zi(E′, E)gi(E′J ′K ′)− dE′+∑q ≠ i kiq(EJK ) N

∑q≠ i kqi(EJK )gq(EJK )−

np

∑p = 1 k p(EJK )gi(EJK )−

(H16)

i

np

∑ k p(EJK )gi(EJK )− p=1

i

N

=

N

∑ kqi(EJK )gq(EJK )−−∑ kiq(EJK )gi(EJK )− q≠i

q≠i

(H17)

at the high-pressure limit, which is eq 62.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected].

E*

Eo

E*

where k−ri (EJK) = 0 and k−di(EJK) = 0, which initially appeared in the left and right-hand side of the eq H16, respectively, since the association and dissociation channels are closed for the energy region of E ≤ E*(JK). Similar to the analysis of eq H9, again, the first two terms on the right-hand side of eq H16 cancel, and eq H16 further simplifies to

gi(EJK )− = fi (EJK )− + λi(EJK )− J′K′

∑∫

× gi(EJK )− −

per reactant pair concentrations, for the energy region of E > E*(JK), which is that reported as eq 6, whereupon noting eqs H11 and H12, and also as noted earlier that the first and second terms of eq H9 canceled each other. The consideration of units discussed earlier for the single potential well is applicable to eq H13 as well. Each term on the right-hand side of eq H13 has again units of (concentration)−1 (time)−1. The kiq(EJK) and kqi(EJK) in the third and fourth terms of eq H13 are the unimolecular rates of isomerization between wells with units of (time)−1. For the energy region of E ≤ E*(JK), k−ri (EJK) = 0, since the association channel is closed including at the high-pressure limit. However, if the bimolecular product channels are open, then we may seek the relationship between the isomerization processes with the processes associated with the bimolecular product channels, according to the delineation of the following steps, in so deriving eq 62: we rearrange eq 39 for gi(EJK)−

Zi(E′, E)gi(E′J ′K ′)− dE′

Zi(E , E′)dE′ − Pi(JK )

J′K′

N q≠i

∫E



o

N

∑ kqi(EJK )gq(EJK )+

Zi′(E′, E)gi(E′J ′K ′)dE′

o

∑ N *p (EJK ) p=1

E*

Eo ∞

np

∑∫

Z′i (E′, E)gi(E′J ′K ′)− dE′

o

i



E*

Eo ∞

∑ N *p (EJK )e−E /kT /hQ i

× e−E / kT /hQ i +

Z′i (E , E′)dE′

o

Pi(JK )∑J ′ K ′ ∫

gi(EJK )− =

= Ni*(EJK )e−E / kT /hQ i +



that appears as the first term in λi(EJK)−, eq 52. The application of l’Hôpital’s rule86,87 (infinity-over-infinity case) to the second term of eq H14 yields

ρi (EJK )

× e−E / kT /Q i

k r+i (EJK )

∫E

ORCID (H14)

Nima Ghaderi: 0000-0002-6257-1373 X

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A Notes

(18) Murakami, Y.; Onishi, S.; Fujii, M. Shock Tube Kinetic Study for the Reaction of H Atoms with SO2: Comparison between Experiments and Theory. J. Phys. Chem. A 2004, 108, 8141−8144. (19) Xu, Z. F.; Hsu, C.-H.; Lin, M. C. Ab initio Kinetics of the Reaction of HCO with NO: Abstraction versus Association/ Elimination Mechanism. J. Chem. Phys. 2005, 122, 234308. (20) Wigner, E. Calculation of the Rate of Elementary Association Reactions. J. Chem. Phys. 1937, 5, 720−725. (21) Beyer, T.; Swinehart, D. F. Algorithm 448: Number of MultiplyRestricted Partitions. Commun. ACM 1973, 16, 379. (22) Lourderaj, U.; Hase, W. L. Theoretical and Computational Studies of Non-RRKM Unimolecular Dynamics. J. Phys. Chem. A 2009, 113, 2236−2253. (23) Baer, T.; Hase, W. L. Unimolecular Reaction Dynamics; Oxford University Press: Oxford, U.K., 1996. (24) Holbrook, K. A.; Pilling, M. J.; Robertson, S. H. Unimolecular Reactions, 2nd ed.; John Wiley & Sons: Chichester, U.K., 1996. (25) Zhu, L.; Hase, W. L. Comparison of Models for Calculating the RRKM Unimolecular Rate Constant k(E, J). Chem. Phys. Lett. 1990, 175, 117−124. (26) Zhu, L.; Chen, W.; Hase, W. L.; Kaiser, E. W. Comparison of Models for Treating Angular Momentum in RRKM Calculations with Vibrator Transition States. Pressure and Temperature Dependence of Cl + C2H2 Association. J. Phys. Chem. 1993, 97, 311−322. (27) Ivanov, M. V.; Grebenshchikov, S.Yu.; Schinke, R. Intra- and Intermolecular Energy Transfer in Highly Excited Ozone Complexes. J. Chem. Phys. 2004, 120, 10015−10024. (28) Kryvohuz, M.; Marcus, R. A. Coriolis Coupling as a Source of Non-RRKM Effects in Triatomic Near-Symmetric Top Molecules: Diffusive Intramolecular Energy Exchange between Rotational and Vibrational Degrees of Freedom. J. Chem. Phys. 2010, 132, 224304. (29) Kryvohuz, M.; Marcus, R. A. Coriolis Coupling as a Source of Non-RRKM Effects in Ozone Molecule: Lifetime Statistics of Vibrationally Excited Ozone Molecules. J. Chem. Phys. 2010, 132, 224305. (30) Ghaderi, N. Bimolecular Recombination Reactions: K-adiabatic and K-active Forms of the Bimolecular Master Equations and Analytic Solutions. J. Chem. Phys. 2016, 144, 124114. (31) Courant, R.; Hilbert, D. Methods of Mathematical Physics; Interscience Publishers: New York, 1989; Vol. 1. (32) Fredholm, I. Sur Une Classe D’équations Fonctionnelles. Acta Math. 1903, 27, 365−390. (33) Hilbert, D. Grundzüge einer Allgemeinen Theorie der Linearen Integral-Gleichungen. Nachr. Ges. Wiss. Göttingen, Math.-Phys. Kl 1904, 49−91. (34) Schmidt, E. Zur Theorie der Linearen und Nichtlinearen Integralgleichungen. I. Teil: Entwicklung Willkürlicher Funktionen nach Systemen Vorgeschriebener. Math. Ann. 1907, 63, 433−476. (35) In Part I,30 the analytic solution for population g(EJK), eq 20 there, was derived assuming that λ and  depended on (EJK), but not on (E′J′K′), and so were brought outside the integral in the third term of eq C2 there, relevant for region E > E*(JK). A numerical analysis of this approximation on the population, for E > E*(JK), informs on three physical regimes responsible for any deviation from the exact magnitude of population, and the correction factor is then assigned for each regime, to correct its results if necessary, whereby one regime appears to stand out. The key physical ingredients in λ and  for consideration of analysis is ZLJ + kd(E′J′K′): (1) when ZLJ ≫ kd(E′J′K′), then the approximation of factoring ZLJ + kd(E′J′K′) out of the integral in eq 17 of ref 30 or to retain it inside the integral for acquiring the exact result for populations g(EJK), yields an error that approaches zero, (2) when ZLJ ≈ kd(E′J′K′), then retaining ZLJ + kd(E′J′K′) inside the integral versus factoring it out leads to ∼15% error in underestimating the populations, (3) when ZLJ ≪ kd(E′J′K′), then retaining ZLJ + kd(E′J′K′) inside the integral versus factoring it out leads to errors approaching 0% as ZLJ→0, since in this case the remaining term in eq H9, f(EJK), of ref 30 dominates the integral term, that is,

The author declares no competing financial interest.



ACKNOWLEDGMENTS It is a great pleasure to acknowledge Prof. R. A. Marcus for the valuable and constructive comments, and kindly noting a great teacher and a remarkable mind. The critically valuable comments of Dr. S. J. Klippenstein are appreciated very much. The very thoughtful and the insightful comments of the anonymous reviewer carrying depth and knowledge are greatly appreciated and duly noted. The various thoughtful comments and/or support of Profs. W. A. Goddard III, H. B. Gray, S. I. Chan, J. Troe, B. Widom, and T. F. Rosenbaum are kindly noted and appreciated very much. The author gratefully thanks Mrs. M. Vojdani (grandmother) for the inquiries, the stimulating discussions regarding bimolecular reactions, and the kind encouragement, as noted with much appreciation. It is a pleasure to gratefully acknowledge support from CalTech and kindly acknowledge the NSF.



REFERENCES

(1) Gilbert, R. G.; Smith, S. C. Theory of Unimolecular and Recombination Reactions; Blackwell Scientific: Oxford, U.K., 1990. (2) Nikitin, E. E. Theory of Elementary Atomic and Molecular Processes in Gases; Oxford University Press: New York, 1974. (3) Seinfeld, J. H.; Pandis, S. N. Atmospheric Chemistry and Physics; John Wiley & Sons: New York, 1998. (4) Barker, J. R.; Golden, D. M. Master Equation Analysis of PressureDependent Atmospheric Reactions. Chem. Rev. 2003, 103, 4577−4592. (5) Fernández-Ramos, A.; Miller, J. A.; Klippenstein, S. J.; Truhlar, D. G. Modeling the Kinetics of Bimolecular Reactions. Chem. Rev. 2006, 106, 4518−4584. (6) Miller, J. A.; Klippenstein, S. J. Master Equation Methods in Gas Phase Chemical Kinetics. J. Phys. Chem. A 2006, 110, 10528−10544. (7) Rissanen, M. P.; Eskola, A. J.; Nguyen, T. L.; Barker, J. R.; Liu, J.; Liu, J.; Halme, E.; Timonen, R. S. CH2NH2 + O2 and CH3CHNH2 + O2 Reaction Kinetics: Photoionization Mass Spectrometry Experiments and Master Equation Calculations. J. Phys. Chem. A 2014, 118, 2176− 2186. (8) Wigner, E. The Transition State Method. Trans. Faraday Soc. 1938, 34, 29−41. (9) Marcus, R. A. Unimolecular Dissociations and Free Radical Recombination Reactions. J. Chem. Phys. 1952, 20, 359−364. (10) Marcus, R. A. Recombination of Methyl Radicals and Atomic Cracking of Ethyl Radicals. J. Chem. Phys. 1952, 20, 364−368. (11) Light, J. C. Statistical Theory of Bimolecular Exchange Reactions. Discuss. Faraday Soc. 1967, 44, 14−29. (12) Miller, W. H. Unified Statistical Model for “Complex’’ and “Direct” Reaction Mechanisms. J. Chem. Phys. 1976, 65, 2216−2223. (13) Ghaderi, N.; Marcus, R. A. Bimolecular Recombination Reactions: Low Pressure Rates in Terms of Time-Dependent Survival Probabilities, Total J Phase Space Sampling of Trajectories, and Comparison with RRKM. J. Phys. Chem. B 2011, 115, 5625−5633. (14) Ghaderi, N.; Marcus, R. A. Bimolecular Recombination Reactions: K-Adiabatic and K-Active Forms of RRKM Theory, Nonstatistical Aspects, Low-Pressure Rates, and Time-Dependent Survival Probabilities with Application to Ozone. 2. J. Phys. Chem. A 2014, 118, 10166−10178. (15) Forst, W. Theory of Unimolecular Reactions; Academic Press: New York, 1973. (16) Troe, J. Theory of Thermal Unimolecular Reactions at Low Pressures. I. Solutions of the Master Equation. J. Chem. Phys. 1977, 66, 4745−4757. (17) Akbar Ali, M.; Barker, J. R. Comparison of Three Isoelectronic Multiple-Well Reaction Systems: OH + CH2O, OH + CH2CH2, and OH + CH2NH. J. Phys. Chem. A 2015, 119, 7578−7592. Y

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

which ref 15 and those cited therein, for example, ref 52, may be consulted for further information. Importantly, the solution consists of a transient term and a steady-state term. For a coarse time resolution, k′ is sufficiently large in the exponent, so the first term, that is, the transient term, decays quickly. In shock-tube experiments, the first term on the right-hand side of the above expression may not be negligible compared to the second term, when the time of observation in experiment is short relative to the transient process. As an example of time scales, considering a simpler unimolecular decomposition interacting with a bath gas, in cyclohexene decomposition, the pulsed laser flash absorption (PLFA) techniques with time resolutions ( λ∑ J′K′

∫E



Z(E′, E) × g (E′J ′K ′)dE′

o

The ratio of the latter inequality favors f(EJK) in eq 17, for example, can be a factor of ∼5 when ZLJ/kd(E′J′K′) ≈ 0.5 for a deactivating intermolecular energy of γ = 50 cm−1 (activation γ′ obtained from detailed balance),42 T = 300 K, and hence effectively masking the error associated with the integral by reducing its contribution to g(EJK) in eq 17, and increasingly so as ZLJ→0. In the above analysis, when ZLJ + kd(E′J′K′) was retained in the integral, it was assumed that kd(E′J′K′) changed by a factor of 3 over 100 cm−1 when integrated over E′, to take into account the sensitivity of the energy dependence of kd. The aforementioned results and comments also apply to the K-active case. In the present paper, the exact analytic solutions are acquired when Z(E′,E) is separable into Z(E′)Z(E) (Appendix E), and the assumption made in ref 30, as discussed above, is now absent. The proof for the solutions, populations, when Z(E′,E) is unfactorizable into Z(E′)Z(E) are given in Appendices C and D, for a single and multiple potential well cases, respectively, in the present work. (36) Brown, N. J.; Miller, J. A. Collisional Energy Transfer in the LowPressure-Limit Unimolecular Dissociation of HO2. J. Chem. Phys. 1984, 80, 5568−5580. (37) Tardy, D. C.; Rabinovitch, B. S. Intermolecular Vibrational Energy Transfer in Thermal Unimolecular Systems. Chem. Rev. 1977, 77, 369−408. (38) Oref, I.; Tardy, D. C. Energy Transfer in Highly Excited Large Polyatomic Molecules. Chem. Rev. 1990, 90, 1407−1445. (39) Setser, D. W.; Rabinovitch, B. S.; Simons, J. W. Collisional Transition Probabilities for Deactivation of Highly Vibrationally Excited Cyclopropane and Dimethylcyclopropane. J. Chem. Phys. 1964, 40, 1751−1767. (40) Rabinovitch, B. S.; Flowers, M. C. Chemical Activation. Q. Rev., Chem. Soc. 1964, 18, 122−167. (41) Gao, Y. Q.; Marcus, R. A. On the Theory of the Strange and Unconventional Isotopic Effects in Ozone Formation. J. Chem. Phys. 2002, 116, 137−154. (42) Zhu, Z.; Marcus, R. A. On Collisional Energy Transfer in Recombination and Dissociation Reactions: A Wiener−Hopf Problem and the Effect of a Near Elastic Peak. J. Chem. Phys. 2008, 129, 214106. (43) Silva, M.; Jongma, R.; Field, R. W.; Wodtke, A. M. The Dynamics of “Stretched Molecules”: Experimental Studies of Highly Vibrationally Excited Molecules With Stimulated Emission Pumping. Annu. Rev. Phys. Chem. 2001, 52, 811−852. and references cited therein. (44) Choi, Y. S.; Moore, C. B. State-Specific Unimolecular Reaction Dynamics of HFCO. I. Dissociation Rates. J. Chem. Phys. 1992, 97, 1010−1021. (45) Choi, Y. S.; Moore, C. B. State-Specific Unimolecular Dissociation Dynamics of HFCO. II. CO Rotational Distribution and Doppler Widths. J. Chem. Phys. 1995, 103, 9981−9988. (46) Klein, F. S.; Herron, J. T. Mass-Spectrometric Study of the Reactions of O Atoms with NO and NO2. J. Chem. Phys. 1964, 41, 1285−1290. (47) In the experiments48 for monitoring ozone formation, a range of concentration for Ar was considered, from 0.45 to 4.3 M. For O2, a range of pressure was considered for up to a total pressure of 100 atm.48 (48) Sauer, M. C., Jr.; Dorfman, L. M. Pulse Radiolysis of Gaseous Argon-Oxygen Solutions. Rate Constant for the Ozone Formation Reaction. J. Am. Chem. Soc. 1965, 87, 3801−3806. (49) Steinfeld, J. I.; Francisco, J. S.; Hase, W. L. Chemical Kinetics and Dynamics; Prentice Hall: New Jersey, 1989. (50) Frankcombe, T. J.; Smith, S. C.; Gates, K. E.; Robertson, S. H. A Master Equation Model for Bimolecular Reaction via Multi-Well Isomerizing Intermediates. Phys. Chem. Chem. Phys. 2000, 2, 793−803. (51) The kinetic scheme described by the Lindemann form as noted in the paragraph yields a solution −d[A]/dt = C2(k′ − k)e−k′t + k[A], where [A] is the concentration of reactant A, and C2, k, k′ are complicated coefficients that involve the individual rate constants, in Z

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A (65) Wang, H.; Hase, W. L. Kinetics of F− + CH3Cl SN2 Nucleophilic Substitution. J. Am. Chem. Soc. 1997, 119, 3093−3102. (66) Miller, J. A.; Klippenstein, S. J. The Recombination of Propargyl Radicals: Solving the Master Equation. J. Phys. Chem. A 2001, 105, 7254−7266, and references cited therein. (67) A Fredholm equation of the second kind has the form

y(x) − λ

∫a

2P(JK )Zo∑J ′ K ′ ∫

∞ −|E − E′| / γ (′)

E*

1 − Zo∑J ′ K ′ ∫

[Z LJ + kd(E′J ′K ′)]−1 dE′

e

∞ −|E − E′| / γ (′)

E*

e

P(J ′K ′)[Z LJ + kd(E′J ′K ′)]−1 dE′

of the second term, where the factor of 2 is noted. When ZLJ = 0, that is, the low-pressure limit, where ZLJ appears in Zo, then the latter expression yields zero due to the Zo appearing in its numerator, and what remains is the first term, that is, eq 35. As pressure increases, but ZLJ < kd(E′J′K′), for example, at low pressure, 1 bar, at T = 300 K with the values of ZLJ and kd(E′J′K′) from the region of (E′J′K′) interest noted in the main text for O3, then the second term reduces to P(JK) (0.002/0.99), where the normalized P(JK) is ∼0.007 for J and K of interest, that is, values of ∼20 and 3, respectively, so the second term is ∼1 × 10−4 to 1 × 10−6 relative to unity of the first term at 1 bar depending on the (EJK) values of interest. The exact numerical evaluation of the first and second term also yields the aforementioned result. So, in this regime eq 35 is a valid approximation of g(EJK)+. Then, in the latter regime of ZLJ relative to kd, the contribution of the second term of g(EJK)+ to the krec of the stable product is also negligible. (75) At the macroscopic level, when the standard practice is made that kd be equal to that for irreversible dissociation, then it has been noted76 accurate when Keqnm(1 − f ne) ≪ 1, where nm is the concentration of the excess reactant for the association reaction, and f ne is the fractional contribution of the association rate coefficient of the slowest relaxing eigenmode of the system.76 Detailed discussion is found in ref 76. (76) Miller, J. A.; Klippenstein, S. J. Some Observations Concerning Detailed Balance in Association/Dissociation Reactions. J. Phys. Chem. A 2004, 108, 8296−8306. (77) Miller, J. A.; Klippenstein, S. J.; Robertson, S. H.; Pilling, M. J.; Green, N. J. B. Detailed Balance in Multiple-Well Chemical Reactions. Phys. Chem. Chem. Phys. 2009, 11, 1128−1137. (78) When a nonequilibrium factor77 for an ith well, at a given E and J, deviates from unity, then it has been found to indicate that the associated microcanonical rate constant might not satisfy detailed balance, where the deviations may occur more readily for very shallow wells. Additional insightful discussions on detailed balance on multiple wells are found in ref 77. (79) The N*(EJK) and ρ(EJK) for O3 are each determined separately, upon drawing from ref 13, and an interpolating function is fitted to each in (EJK) space, which in turn is used in a function call for determining kr(EJK) and kd(EJK) in eq 35, where in the cases of eqs 31, 34, and E15, the g(EJK) dependence regarding their second terms was further noted. The evaluation of g(EJK)+, eq 35, utilized 7.8 × 10−5 seconds in real time with a CPU at 3.40 GHz and 15.6 GB of RAM memory, on a 64 bit linux type of operating system. The calculation of the second term of eqs 31, 34, or E15 for g(EJK)+ utilized time on the order of minutes encompassing the entire region of interest for (E′J′K′) and (EJK). The calculation of the thermal krec, eq 65, for O3 utilized 0.04 seconds, after acquiring g(EJK)+ with its second term present or not depending on the pressure region, as the overall magnitude of the population in (EJK) space does not practically affect the computation time of krec. A forthcoming investigation of the temperature and pressure effects on ozone formation will also issue forth the details of calculation encompassing those in this work. (80) For the orbital hindered part of the Hamiltonian,14 we have H = j2/2μBCr2 + l2/2μABCR2 + V(r,R,θ), apart from the kinetic energy term, where all terms are accounted for in the calculation. The j is the action variable for the angular momentum j/2π of the diatom O2, and l is the action variable for the orbital angular momentum l/2π of the two colliding partners, O and O2. The V is the three-body potential energy81 for O3 as a function of the Jacobi coordinates: the diatomic instantaneous bond length r, center-to-center distance R, and the angle θ between the vectors r and R. The reduced masses for the diatomic and the colliding pairs, that is, a triatomic, are μBC and μABC, respectively. (81) Varandas, A. J. C.; Murrell, J. N. Dynamics of the 18O + 16O2 (υ = 0)Exchange Reaction on a New Potential Energy Surface for GroundState Ozone. Chem. Phys. Lett. 1982, 88, 1−6. The reader may consult

b

k(x , t ) × y(t )dt = f (x)

where y(x) is the unknown function, f(x) and the kernel k(x, t) are known, and λ is a parameter.31−34 In the present work, the unknown function is g(EJK), that is, the population in the K-adiabatic case, and is acquired using algebraic method. The solution for the K-active case is also given. (68) The calculated thermal recombination rate constant, eq 65, for O3 formation, based on a population that utilizes the factoring of Z(E′,E), eq 15, into Z(E)Z(E′) versus not as in eq 34, yields values of krec, by eq 65, that are within 1% of each other at 1 bar, and within 2% at 1000 bar, at a given temperature, upon using the single exponential model1,16 for intermolecular energy transfer. (69) The g(EJK)+, eq 26, is a solution of the bimolecular master equation eq 19, in which is utilized an algebraic approach for its derivation, with its proof given in Appendix C, and for the separable Z(E′,E) case considered in Appendix E, eq E3. In Part I,30 an approximate solution similar in form, eq 24 there, was derived using property of inhomogeneous linear integral equations (with Fredholm and Hilbert−Schmidt theory),31−34 which used an assumption noted in footnote 35, whereupon the correction factor to the population there is known a priori for a given calculation, when needed, to acquire the near exact result, as discussed therein for it in footnote 35. The solution, eq E3, in the present paper omits these assumptions and consequences thereof and yields the exact solution of the bimolecular master equation. (70) The solution for the population, g(EJK)+, eq 31, has the same first term as the approximate solution from Part I,30 eq 29 there. The second term of the solution, eq 31, is augmented by terms in its denominator, in comparison to its counterpart from Part I.30,35 (71) Montroll, E. W.; Shuler, K. E. The Application of the Theory of Stochastic Processes to Chemical Kinetics. Adv. Chem. Phys. 1958, 1, 361−399. (72) Widom, B. Collision Theory of Chemical Reaction Rates. Adv. Chem. Phys. 2007, 5, 353−386. (73) For the evaluation of the integral in the second term of eq 34, for the energy region E′ ≤ E, the lower limit of integration is the least value of E′ such that E′ > E*(JK), so to validate the condition kd(E′J′K′) > 0, and its upper limit is E. For region E′ > E, the lower limit of the integral is a selected E, such that E > E*(JK), while its upper limit is ∞, or can be set, for practical purposes, to such a value that would ensure the convergence of the integral. (74) The comparison of the first with the second term of g(EJK)+, eq 31 or 34, informs that ZLJ + kd(EJK) appears in the denominator of each term and A·BC appears in the numerator of each term. Then, the constituent terms that remain in the comparison between the first and second term is the kr(EJK) of the first term with the remaining terms in the second term. The kr(E′J′K′) in the numerator of the second term may change by a factor of ∼2 in the integrand as E′ increases in the region of interest from 0.1 to 1 kT, with T = 300 K, in excess of the critical energy for J′ and K′ of interest for the primed channel, and so kr(EJK) of the first term approximately cancels the kr(E′J′K′) of the second term when the latter is brought outside the integral, but the factor ∼2 retained in the second term in the analysis to note the change. If the channel associated with kr(E′J′K′) is assumed negligible when applicable, then kr(E′J′K′) would be brought out of the integral and treated as kr(EJK) and the factor of 2 dropped. When the latter factor is retained, the comparison is now between unity from the first term and the AA

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A lim Q = ∞

the Supporting Information of ref 14 for a discussion on the treatment of the potential energy surface. (82) The W(EJK) = 2J∫ ···∫ Θ(EJK,Γ6)dΓ6, where Γ6 denotes the set of generalized initial rotational-orbital-vibrational action-angle variables (j,l,n,wj,wl,wn) at any given total E, J, and K for the atom−diatom system. The step function Θ(EJK,Γ6) is unity if an intermediate is formed from the colliding partners and zero otherwise. An initial integration over J, wJ, and wM yielded the 2J, unity, and unity, respectively, in the total J representation. A transformation from l to mj can be made using the relation l2 = j2 + J2 − 2mjJ, with a Jacobian of J/l, where the action mj is identified as K, in the approximation of using the line of centers of the atom and diatom as the body-fixed principal axis associated with the least moment of inertia. The W(EJK) is equivalent to

∫ ···∫ δ(K − mj) × Θ(E −

N *(EJK ) = min 2J R

T →∞

Delineation entails: factoring T in the numerator of Q yields T5/2e−hv/2kTQelecC where C denotes the remaining terms independent of T. The denominator of interest for Q is 1 − e−hv/kT, and so then lim Q =

T →∞

where the chain rule84 is employed in treating the numerator. Since Q appears in the denominator of eq 58 then

lim k r,+∞(EJK ) ∝ 1/∞ = 0

T →∞

as noted in the text. (89) We show the lim Q = ∞, where Q is the product of partition

H′)dΓ6J

T →∞

functions, Q = QrotQvib, for the complex in its center of mass. For example, for a triatomic,

where H′ is the total energy minus the PR2/2μ at the TS with dΓJ6 = dwjdwldjdldrdpr/h. Additional elucidations are found in refs 13 and 14. (83) Miller, J. A.; Klippenstein, S. J.; Robertson, S. H.; Pilling, M. J.; Shannon, R.; Zádor, J.; Jasper, A. W.; Goldsmith, C. F.; Burke, M. P. Comment on “When Rate Constants Are Not Enough. J. Phys. Chem. A 2016, 120, 306−312. (84) The term “well-skipping” or “well-mediated” applies to the phenomenological description of a reaction, even though the reactants may pass through an intermediate well, or an intermediate that may be long-lived above the dissociation limit, to then reach products, for example, C2H5 + O2 → C2H4 + HO2, even though the latter reaction has the reactant pass through C2H5O2.83 At the microscopic level, the C2H5O2 is not skipped en route to products. The phenomenological rate constant or rate coefficient is a function of macroscopic thermodynamic parameters, for example, temperature or pressure, distinct from a microscopic rate coefficient k(EJK) or k(EJ) being a function of microscopic molecular constants of motion. The former description of the rate constant may describe a reaction between configurations that are not adjacent to each other (well-skipping), for example, allene with propyne,85 and in general not interpreted as the probability per unit time that a reactant(s) will react to form products.83 (85) Miller, J. A.; Klippenstein, S. J. From the Multiple-Well Master Equation to Phenomenological Rate Coefficients: Reactions on a C3H4 Potential Energy Surface. J. Phys. Chem. A 2003, 107, 2680−2692. (86) Stein, S. K. Calculus and Analytic Geometry, 4th ed.; McGrawHill: New York, 1987. (87) The l’Hôpital’s rule (infinity-over-infinity case)86 states that if f and g are defined and differentiable for all x larger than some fixed number, then, if

3

Q vib =

where Θvn = hvn/k with the Boltzmann constant k, Planck constant h, and the frequencies vn(Hz) for each mode, n, and Qrot = π1/2(8π2I1kT/ h2)(8π2I3kT/h2)1/2/σ, where σ is a symmetry number, and moments of inertias I3 < I1 ≈ I2 for a prolate near-symmetric top, for example, ozone. The application of the l’Hôpital’s rule86,87 and chain rule86 to Q then yields

lim Q = ∞

T →∞

Since Q appears in the denominator of ρ(EJK)e−E/kT/Q, as the assumed functional form for g in the isomerization terms of eq 60, then

lim k r+i, ∞(EJK ) ∝ 1/∞ = 0

T →∞

as noted in the text. (90) Chesnavich, W. J.; Bass, L.; Su, T.; Bowers, M. T. Multiple Transition States in Unimolecular Reactions: A Transition State Switching Model. Application to the C4H8+· System. J. Chem. Phys. 1981, 74, 2228−2246. (91) Klippenstein, S. J.; Georgievskii, Y.; Harding, L. B. Statistical Theory for the Kinetics and Dynamics of Roaming Reactions. J. Phys. Chem. A 2011, 115, 14370−14381. (92) Harding, L. B.; Klippenstein, S. J.; Jasper, A. W. Separability of Tight and Roaming Pathways to Molecular Decomposition. J. Phys. Chem. A 2012, 116, 6967−6982. (93) Hirschfelder, J. O.; Wigner, E. Some Quantum-Mechanical Considerations in the Theory of Reactions Involving an Activation Energy. J. Chem. Phys. 1939, 7, 616−628. (94) Klippenstein, S. J.; Marcus, R. A. Application of Unimolecular Reaction Rate Theory for Highly Flexible Transition States to the Dissociation of CH2CO into CH2 and CO. J. Chem. Phys. 1989, 91, 2280−2292. (95) Marsden, J. E.; Hoffman, M. J. Basic Complex Analysis, 3rd ed.; W.H. Freeman and Company: New York, 1999. (96) Churchill, R. V. Complex Variables and Applications, 2nd ed.; McGraw-Hill: New York, 1960. (97) Press, W. H.; Teukolsky, S. A.; Vetterling, W. T.; Flannery, B. P. Numerical Recipes: The Art of Scientific Computing, 3rd ed.; Cambridge University Press: New York, 2007. (98) In Part I,30 the analysis for the boundary condition and analytic continuation were given for a single potential well. There, the solution for g utilized an assumption of the independence of λ and A from E′, as noted in footnote 35. Also, there, in ref 30, eq D1 should possess an E* in its exponent on the right-hand side, and eq D2 should possess

x →∞

and lim

x →∞

f ′(x) =L g ′(x)

then it follows that

lim

x →∞

∏ e−Θvn /2T /(1 − e−Θvn / T ) n=1

lim f (x) = ∞ , lim g (x) = ∞

x →∞

d d [T 5/2 × e−hv /2kT Q elecC ]/ [1 − e−hv /2kT ] = ∞ dT dT

f (x) =L g (x)

For further discussion please see ref 86. (88) The product of the partition functions for the atom−diatom reactants, Q = QrotQvibQtransQelec, is given by Q = (kT/σB)(e−hv/2kT/[1 − e−hv/kT])(2πμkT/h2)3/2 (15 + 9e−158.5/kT + 3e−226.5/kT), where the first, second, third, and fourth term correspond to Qrot, Qvib, Qtrans, (per volume), and Qelec, respectively. For illustration, Qelec for ozone is considered. The σ is the usual symmetry number, for example, 2 for O2, the B is the rotational constant, ν is the vibrational frequency, and μ is the reduced mass. The application of l’Hôpital’s rule (infinityoverinfinity)84,85 case to Q yields

e2π(v″−v)E*. The attainment of the subsequent results, eqs D3 and D4 in ref 30 assumed the presence of these factors. (99) The application of the l’Hôpital’s rule86,87 to eq H1 gives AB

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A lim g (EJK )+ =

Z LJ →∞

⎡ d ⎢ ∑ dZ LJ ⎢⎣ J ′ K ′ /



∫E*

d ×[ dZ LJ

⎤ Z(E′, E) × g (E′J ′K ′)dE′⎥ ⎥⎦



∫E*

Z(E , E′)dE′ + kd(EJK )]

where the numerator in the latter equation yields eq H2, and the denominator is equal to 1.



NOTE ADDED AFTER ASAP PUBLICATION This paper was published ASAP on March 27, 2018, with errors in two equations in Reference 88. The corrected paper was reposted on March 29, 2018.

AC

DOI: 10.1021/acs.jpca.7b09244 J. Phys. Chem. A XXXX, XXX, XXX−XXX