Ab Initio Potential Energy Surfaces and Quantum Dynamics for

Mar 26, 2018 - There has been great progress in the development of potential energy surfaces (PESs) and quantum dynamics calculations in the gas phase...
3 downloads 0 Views 2MB Size
Perspective Cite This: J. Chem. Theory Comput. 2018, 14, 2289−2303

pubs.acs.org/JCTC

Ab Initio Potential Energy Surfaces and Quantum Dynamics for Polyatomic Bimolecular Reactions Bina Fu* and Dong H. Zhang* State Key Laboratory of Molecular Reaction Dynamics and Center for Theoretical and Computational Chemistry, Dalian Institute of Chemical Physics, Chinese Academy of Sciences, Dalian 116023, China ABSTRACT: There has been great progress in the development of potential energy surfaces (PESs) and quantum dynamics calculations in the gas phase. The establishment of a fitting procedure for highly accurate PESs and new developments in quantum reactive scattering on reliable PESs allow accurate characterization of reaction dynamics beyond triatomic systems. This review will give the recent development in our group in constructing ab initio PESs based on neural networks and the timedependent wave packet calculations for bimolecular reactions beyond three atoms. Bimolecular reactions of current interest to the community, namely, OH + H2, H + H2O, OH + CO, H + CH4, and Cl + CH4, are focused on. Quantum mechanical characterization of these reactions uncovers interesting dynamical phenomena with an unprecedented level of sophistication and has greatly advanced our understanding of polyatomic reaction dynamics.

I. INTRODUCTION Chemical reactivity is essentially a process of the breakage of an old chemical bond and the formation of a new chemical bond. For a bimolecular reaction, the collision of two reactants leads to new product molecules with the rearrangement of atoms. Characterizing chemical reactions in the gas phase at the most fundamental level has been an important and fundamental goal of modern physical chemistry. Reaction dynamics can interpret reactive collisional processes with quantum state resolution and provide essential mechanistic information on chemical reactions. Tremendous progress has been made for theoretical and experimental research in chemical reaction dynamics in the recent past, moving from the very thorough studies of some triatomic reactions to truly polyatomic reactions beyond three atoms.1−32 The theoretical and computational approach to polyatomic reaction dynamics has been challenged in two major ways. The first is the potential energy surface (PES), and the second is the dynamics, for which the accurate quantum reaction dynamics is essential. A single PES governing the motion of the nuclei is defined as a function of nuclear coordinates in the Born−Oppenheimer (BO) approximation. In the past few decades, development of electronic structure theory and fitting techniques of ab initio electronic structure data has made the construction of accurate global PESs possible for multidimensional molecular systems beyond the chemical accuracy.33−43 Besides, the quantum dynamical methodologies and computational algorithms have been greatly advanced along with the fast growth of numerical computing power. This has made full quantum dynamics calculations for polyatomic reactions based on first principles truly possible.23−29 The construction of an accurate PES needs a large number of high-level ab initio points over a large configuration space and a faithful analytical representation of the discrete ab initio data. With the advances of computational capacity, high-level ab initio theories, such as the coupled-cluster singles, doubles, and perturbative triples (CCSD(T))44 and multireference configurational interaction (MRCI) methods,45 together with large © 2018 American Chemical Society

basis sets, can basically give reliable determination of electronic energy. However, generating a highly accurate analytical representation of the discrete ab initio data for multidimensional systems poses a big challenge. Splines and nonlinear fitting based on many-body expansions, functional forms, and so on46,47 are efficient in representing the PESs for systems with no more than three atoms but do not scale well for larger reaction systems. For reaction systems with more than three atoms, some general fitting approaches have been developed to represent high-dimensional PESs.33−42 The interpolating moving least-squares (IMLS)33 and modified Shepard interpolation methods34 are examples of interpolation approaches to represent the PESs for polyatomic reactions. This group once constructed full-dimensional PESs for the CH5 and OH3 reaction systems using the modified Shepard interpolation scheme.23,48 However, the evaluation of these PESs is extremely slow due to the nature of the highdimensional interpolation method. Perhaps the first significant step in constructing accurate high-dimensional PESs was taken by Bowman and co-workers, who developed a method based on permutation-invariant polynomials (PIPs),35,36 which has been applied to a large number of permutationally invariant polyatomic PESs. The PIP method uses primary and secondary invariants as the fitting basis to generate permutation invariant polynomials while the coefficients are obtained by linear leastsquares fitting. Based on the PIP approach, we constructed accurate global PESs for the O(3P) + C2H4 and O(1D) + CH4 reactions.49−51 Due to the complexity of these PESs, the fitting is challenging to describe accurately all of the stationary points and reaction channels. In addition, the QCT calculations on these PESs produced dynamical results in good accord with the experimental results. Another mathematically equivalent method based on the symmetrized monomials was also developed by Bowman and co-workers.52 More recent efforts Received: January 3, 2018 Published: March 26, 2018 2289

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation

packet (TDWP) approach captures the quantum reaction dynamics information by solving the time-dependent Schrödinger equation, which is numerically much more efficient than the time-independent methods.5,6,25,60,61 The TDWP approach makes large-scale reactive scattering calculations possible, and thus it has been widely used in multidimensional quantum dynamics studies on polyatomic reactions.23−30 The dynamical information, such as the initial state-selected total reaction probabilities and cross-sections, can be obtained from the applications of the initial state-specific wave packet (ISSWP) approach. Since the ISSWP method resolves an initial value problem, the computational effort as well as the gained dynamical information corresponds to the action on a specific initial state of interest. One can rely on the reactant Jacobi coordinates for the accurate and efficient wave packet propagation either with the time or the Chebyshev evolution operator,6,61−66 which is similar to that for nonreactive problems. To avoid boundary reflections, the wave function is absorbed at the edges of the grid which is accomplished by simply multiplying the wave function by an exponentially decaying function of a coordinate near the boundary at the end of each propagation step. The ISSWP calculation can now reach far beyond atom−diatom reactions, and it has been extended to larger polyatomic reactions, such as the reactions of the X + YCZ3 type.24,26−29 The state-to-state dynamical information is much more difficult to capture, mainly owing to the lack of optimal coordinates for simultaneous description of the reactant channel, strongly interacting region, and product channel. The “reactant coordinates based (RCB)” and “product coordinates based (PCB)” approaches were used to investigate the state-to-state dynamics for atom−diatom reactions, where the reagent to product coordinate transformation is carried out at one time step around the transition state region.67−73 The reactant−product decoupling (RPD) approach, originally proposed by Peng and Zhang,74 has been successfully applied in recent state-to-state calculations beyond three-atom reactions.23,75−78 This method first propagates the wave packet in the reactant channel and transforms the no-return part of the reacted wave packet continuously in time from reactant to product coordinates with the help of absorbing potentials, and further propagates the wave packet into the product asymptote. In addition, we can perform the wave function absorption after multiple propagation steps,75,76 instead of at each propagation step, which makes the RPD approach more appealing for stateto-state dynamics studies of polyatomic chemical reactions. The RPD approach resolves, to a great extent, the problem of the choice of coordinates in state-to-state TDWP calculations, in particular for direct activated reactions. The fully converged differential cross-sections (DCSs) for four-atom reactions with all six degrees of freedom treated exactly were calculated efficiently by the TDWP method based on the RPD scheme in this group.23 The calculated DCS for the HD + OH → H + HOD reaction agrees quite well with a crossed-molecular beam experiment,23 representing a milestone in quantum reactive scattering calculations because of the first complete state-tostate treatment for reactions of beyond three atoms. Manthe and co-workers recently proposed a transition state wave packet (TSWP) approach to state-to-state reaction dynamics,79,80 which is based on the quantum transition state theory of Miller.81 The TSWP method has been employed with both the multiconfiguration time-dependent Hartree (MCTDH)28,79 and traditional grid/basis framework82−84 in

have been made to construct the highly accurate multidimensional PESs based on artificial neural networks (NNs).53 A key advantage of the NN approach is its more faithful representation of ab initio points. A systematic procedure based on the NN fitting was recently proposed in this group to construct accurate PESs, as recently demonstrated for OH3, HOCO, and CH5 systems with very small root-mean-square error (RMSE).38,39,54 These NN PESs did not enforce the permutational symmetry of identical atoms explicitly, which can produce accurate results in quantum dynamics calculations but may introduce some errors in QCT calculations. An exchange scheme to solve the problem was proposed, so that the conservation of the total energy can be substantially improved with the exchange scheme in the QCT calculations.38 Guo and co-workers were recently inspired by PIP and proposed the PIP-NN fitting method.40,41 The PIP-NN PESs use a set of PIPs instead of pairwise distances as input vector in the neural network, which can automatically guarantee the permutational symmetry. The input PIPs include all the polynomials truncated by a given degree. However, the number of polynomials increases considerably with the degree bound. It is impractical to include all the invariant polynomials with the highest degree bound of secondary invariants in the input. Although the number of polynomials can be reduced with a lower degree bound, a large number of secondary polynomials have to be discarded. Very recently, a more flexible NN approach using the fundamental invariants (FIs) as the input vector was proposed in this group to construct the PESs with permutational sysmmetry.42 FI-NN is mathematically precise, and it can approximate any permutation invariant PESs to arbitrary accuracy. Different from the primary and secondary polynomials in PIP and the polynomials in the PIP-NN method, FI minimizes the size of input invariants which can generate all the invariant polynomials. Therefore, FI-NN can efficiently reduce the evaluation time of potential energy compared to the corresponding PIP-NN PESs, in particular for larger molecular systems with more identical atoms. The FIs for all possible molecular systems of up to five atoms were provided. The PESs for OH3 and CH4 were constructed with FI-NN with very small RMSE, whose accuracy was further confirmed by fulldimensional quantum dynamics scattering and bound state calculations.42 Significant progress has also been made for multidimensional quantum reactive scattering calculations. Quantum reactive scattering aims to extract all observable quantities of chemical reactions by solving the time-independent or time-dependent Schrödinger equation for nuclear motion. It not only provides highly accurate theoretical results which can be directly compared in detail with the experimental observations but also gives an in depth understanding of reaction dynamics. In addition, quantum scattering calculations on chemical reactions also serve as quantum mechanical benchmarks for testing approximate theories which can more readily be applied to more complicated reactions. In the past decades, state-to-state dynamics problems have been resolved for many direct atom− diatom reactions by using hyperspherical coordinate based time-independent quantum dynamics methods,55,56 which was often accomplished through the application of the ABC program.57 However, due to the dramatic rise of the computational costs with the increase of basis functions, it is difficult to extend the time-independent approach to larger polyatomic reactive systems.58,59 The time-dependent wave 2290

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation

The threshold of the jth neuron of the lth layer is determined by the biases blj, and the transfer functions f1 and f 2 are taken as hyperbolic tangent functions. During the NN fitting procedure, we have to perform many tests by using a different number of neurons for the two hidden layers based on a specific set of data points, in order to get the smallest fitting error, mainly because the structure of the NN intensely influences the quality of fit. For a specific NN structure, the weights and biases in eqs 1−3 can be updated through proper NN training using the Levenberg−Marquardt algorithm.95 The root-mean-square error

investigating the state-to-state dynamics for atom−diatom reactions and some reactions of beyond three atoms. Besides, a new RCB method for computing S-matrix elements and DCSs for tetratomic reactions was proposed by Guo and coworkers.85−87 It involves the interpolation of the timedependent wave packets, using a collocation method at selected time intervals on the product grid that naturally defines the product asymptotic states. This Perspective will be short reviews of primary research progress on PESs and quantum dynamics of polyatomic bimolecular reactions in the gas phase in this group in recent years. Readers can refer to ref 30 for a more complete review on the field of quantum scattering calculations for polyatomic bimolecular reactions published in the past 10 years. In sections II and III, we will give briefly the neural network based fitting methods and TDWP methodologies, as well as some detailed illustrations of PESs and quantum reaction dynamics for polyatomic reactions developed in this group. We conclude in section IV with a brief summary of this review.

RMSE =

I



i=1



∑ (w1j ,ixi)⎟⎟,

j = 1, 2, ..., J (1)

and the output of kth neuron in the second hidden layer is yk2



=

f 2 ⎜⎜bk2 ⎝

J

+



(wk2, j

j=1



y1j )⎟⎟ , ⎠

k = 1, 2, ..., K (2)

Therefore, the final result is given by K

y = b13 +

∑ (w1,3k yk2 ) k=1

n

∑ (Efit − Eab initio)2 i=1

(4)

is used to measure the fitting error. Further, we used the “early stopping” method96 to improve the fitting quality. The entire data set was divided into the training and validation sets, and the training procedure was inhibited when the overfitting occurs. We should have sufficient data points that cover all the important regions on a global PES, and thus a good sampling scheme of these data points is very important. An iterative fitting and sampling scheme proposed by Behler94 was employed in our work. Specifically, we first sample a small set of molecular configurations and do NN fittings to these points for several times. This is followed by extensive QCT calculations, by which a complete geometry space is sampled on one of the preliminary fittings. If a molecular configuration from a trajectory is located far from the existing data set, or the fitted energies for this configuration on several fittings are significantly different, we will perform ab initio calculations for this data configuration and add it to the data set for further NN fittings. We have to repeat the above iterative procedures for many times until the PES is converged with respect to the number of data points. A reliable PES should also be rigorously examined by dynamics calculations. The PESs developed in our group are always examined by the accurate quantum dynamics calculations, to investigate the final convergence of a PES. Usually, we perform the quantum dynamics calculations on a fitted PES when the RMSE is smaller than an acceptable value, only a small number of data points have high fitting errors, and the trajectories do not go to any artificial holes. After each fitting procedure, we save at least two PESs with least RMSEs and perform quantum dynamics calculations on them. If the dynamics results on these saved PESs and on the PESs fitted with fewer data points agree well, we believe the PES is finally converged. Otherwise, we will add more data points taken from the trajectories, again do the NN fitting for the new data set, and further carry out quantum dynamics calculations. 1. OH3. The H2 + OH ↔ H2O + H reactions and isotope analogues have become the prototypes for tetra-atomic reactions, which play an important role in both atmospheric chemistry and combustion.97,98 As three of the four atoms in the reactive system are hydrogens, the OH3 system is an ideal candidate for high-quality electronic structure calculation of the PES as well as for accurate quantum reactive scattering study. In the past decades, a number of PESs, including WDSE,99 OC,100,101 WSLFH,102 YZCL2,103 and XXZ23 PESs have been constructed for the reaction system using different methods. These PESs were either inaccurate or slow to be evaluated. Recently, we have reported an accurate PES using the NN

II. POTENTIAL ENERGY SURFACES II.A. Neural Networks. Neural networks are general fitting methods. The NN is extremely flexible, and in principle it can be used to fit any shape of function with very high accuracy. Over the past decades, many high-dimensional PESs have been fitted by the NN approach based mainly on density functional theory (DFT) or relatively low level ab initio calculations.37,88−94 However, it remains a problem to sample molecular configurations efficiently for ab initio calculations and the NN fitting. Recently, our group developed accurate global PESs for the polyatomic bimolecular reactions, such as OH3, HOCO, and CH5 reactive systems, using the NN fitting based on high-level ab initio data points. A systematic fitting procedure for selectively adding new ab initio points was proposed in our work. Next we will give a brief introduction for the NN fitting, with the OH3 system illustrative of the fitting procedure. We employed the feed forward NN to fit a PES from extensive ab initio calculations. The feed forward NN has two hidden layers, which connect the input layer and output layer. This is denoted as I−J−K−1 NN, with I nodes in the input layer, which corresponds to the number of input coordinates, i.e., the internuclear distances of the configuration. The output layer returns the potential energy of this configuration. As seen, there are J and K neurons in the two hidden layers, respectively. The output of jth neuron in the first hidden layer is ⎛ y1j = f 1 ⎜⎜b1j + ⎝

1 n

(3)

where xi (i = 1, ..., I) denote the internuclear distances of a molecular configuration. The ith and jth neurons of (l − 1)th and lth layers, respectively, were connected by the weights wlj,i. 2291

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation method based on high-level ab initio calculations.38 Roughly 17000 ab initio energies for the H2 + OH ↔ H2O + H reaction were calculated at UCCSD(T)-F12a/AVTZ level of theory and used in the NN fitting. For the OH3 system, the hydrogen atoms are permutated in advance to have the six input internuclear distances (ROH1, ROH2, ROH3,RH1H2, RH1H3, RH2H3) satisfying ROH1 ≤ ROH2 ≤ ROH3. To accelerate the NN fitting procedure, we divided the data points into three regions, i.e., OH + H2 product channel (I), H + H2O product channel (II) and interaction part (III), with some overlaps between them. We can use a much more efficient fitting scheme without any loss of accuracy by fitting the three regions separately, due to less neurons required to reach a desired RMSE. A global PES was obtained as a weighted sum of the three segmentally fitted parts

reactions, respectively, calculated on the three PESs for the total angular momentum J = 0. The three reaction probability curves are essentially identical, indicating the NN1 PES is well converged with respect to the number of data points as well as the fitting procedure. We denote the NN1 PES as the CXZ PES, with the overall RMSE of 1.61 meV. Compared with the existing YZCL2 and XXZ PESs which were constructed by the modified Shepard interpolation approach in this group, the CXZ PES is much faster on evaluating and substantially more smooth and accurate on fitting. In addition, the DCS of the D2 + OH → HOD + D reaction calculated by the TDWP approach on the CXZ PES achieves good agreement with the experiment.104 All of these indicated the CXZ PES represents the best available PES for the benchmark four-atom reaction. 2. HOCO. The OH + CO → H + CO2 reaction is considered as one of the most important reactions in combustion chemistry and atmospheric chemistry, due to its crucial role in the conversion of CO to CO2.105,106 This exothermic reaction proceeds via two pathways with the formation of a HOCO intermediate complex in trans and cis deep wells, which represents a prototype for complex-forming four-atom reactions. The further decomposition of the HOCO species over significant barriers leads to the production of CO2 and H. As the OH + CO reaction involves three heavy atoms and the long-lived intermediate complex, it presents a great challenge to both PES construction and quantum dynamics. There are some old PESs for this system,107−112 but they failed to provide a faithful representation of the complex potential topography, and are not sufficiently accurate for detailed dynamical studies. In 2012, Li et al. developed a global PES(CCSD- 1/d PES) and its modified version (CCSD- 2/d PES) using the PIP fit to roughly 35000 UCCSD(T)-F12a/ AVTZ data points.113 In 2013, we reported another global PES using the NN method,39 with data points spreading in a more complete configuration space and computed at essentially the same level of theory. It was found that the NN PES fits the ab initio points much better than the PIP PES. The NN fitting scheme used in constructing the PES of HOCO is similar to that of OH3, which has been introduced above. A total of 74400 data points were calculated at the UCCSD(T)-F12a/AVTZ level of theory and were included in the NN fitting, which were selected by iterative trajectory calculations and the scheme originally proposed by Behler.94 The number of ab initio data points we employed in the fitting is substantially larger than that for the CCSD-2/d PES by Li et al. Due to the large number of data points, we used the segmental fitting scheme by dividing the data points into HO + CO entrance part, HOCO interaction part, and H + CO2 exit part and performing the NN fitting separately for the three parts. Since the two deep wells supporting resonances exist on the HOCO PES, it is extremely difficult to obtain the converged quantum dynamics results to the level that we achieved for the OH3 reactive system. We used the neural network ensembles to make the overall convergence of the PES better. Explicitly, we performed NN fitting repeatedly for a given set of data points and roughly 40 PESs were generated, among which six PESs were picked up with the least fitting errors for average potential energy calculations. The final averaged PES results in the overall RMSE for data points below 1.5 eV is 4.2 meV. 3. CH5. The H + CH4 → H2 + CH3 reaction and its reverse are of crucial importance in CH4/O2 combustion chemistry.98

E = W1E1 + W2E2 + W3E3 W1 = logsig(50(R OH2 − 3.5)) W2 = (1 − W1) logsig(50(R OH3 − 6.4)) W3 = 1 − W1 − W2

where logsig is a logistic sigmoid function defined by the formula logsig(x) =

1 1 + e −x

Ei and Wi are the energy and weight of these three parts. The resulting global PES was denoted as NN1 PES. Another NN fitting was also performed to use the same data points as NN1, but with a larger number of neurons, to yield a global PES denoted as NN2 PES. Furthermore, to check the convergency, NN3 PES was fitted with 20% of data points removed randomly from the original data set. The quantum dynamics calculations were performed using the ISSWP method to check the convergence of fitted PES. Figure 1 shows the total reaction probabilities for the H2 + OH → H + H2O, and the exchange H + H′OH → H′ + H2O

Figure 1. Reaction probabilities for the OH + H2 →H + H2O reaction on NN1, NN2, and NN3 PESs as a function of collision energy. Reprinted with permission from ref 38. Copyright 2013 AIP Publishing. 2292

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation They represent benchmark six-atom reactions to advance our understanding of polyatomic chemical reactivity. Because five atoms involved are hydrogen atoms, it is amenable for highlevel ab initio calculations of the PES and quantum dynamics calculations. Substantial efforts have been devoted to constructing accurate PESs for this prototype six-atom reaction. The early semiempirical PESs are not considered quantitatively accurate.114,115 In 2004, Manthe and co-workers used the modified Shepard interpolation method to construct a full-dimensional PES based on CCSD(T)/cc-pVQZ level of theory.116 This PES is limited to the vicinity of the abstraction reaction barrier and is not reliable for global dynamics studies. Bowman and coworkers developed several versions of global ab initio PESs for this system,117,118 namely, ZBB1, ZBB2, and ZBB3, which are considerably more accurate than all previous PESs. These PESs were fitted by the PIP approach based on RCCSD(T)/AVTZ data points, which includes both the abstraction and exchange channels. Later we constructed another global PES (ZFWCZ) using the modified Shepard interpolation method based on UCCSD(T)/AVTZ data points,48 which has a comparable accuracy to ZBB3 PES. However, the evaluation of the PES is extremely slow due to the nature of the high-dimensional interpolation method. Recently, we developed a new global PES (XCZ) for this system using the NN fitting procedure based on UCCSD(T)F12a/AVTZ data points,54 resulting in a very small total RMSE. Starting from the ∼50000 molecular configurations incorporated in ZBB3 and ZFWCZ PESs, we first selected roughly 8000 configurations by the distance selecting method. Then iteratively, we carried out ab initio calculations for selected configurations, performed NN fitting, and ran trajectory calculations on one of the fittings to add more data points, as we employed in developing PESs of OH3 and HOCO systems. Finally a total of 47783 points were included in the NN fitting to get a converged PES. In order to accelerate the fitting procedure and improve the fitting accuracy, we divided the data points into the CH4 + H asymptotic part, CH5 interaction part, and CH3 + H2 asymptotic part, and fitted the three parts and checked the convergence properties separatively. The global PES was obtained by summing the three parts with a weighting factor. The equation showing the weighted sum with details on the weights given in ref 54, which was similar to that for the OH3 PES. To check the convergence of the PES with respect to the number of ab initio points, we randomly deleted 10% data points from the original 47783 data points and used the same fitting procedure to obtain another PES. Quantum scattering calculations were carried out on the two PESs. As shown in Figure 2, the total reaction probabilities calculated on the two PESs are essentially the same, indicating the PES is well converged with respect to the number of ab initio points. The XCZ PES is much faster on evaluating and more accurate on fitting compared with the earlier ZFWCZ PES, representing the best available PES for the CH5 system. II.B. PIP-NN. For a chemical reaction containing identical atoms, the potential energy should be permutation invariant with the permutations of identical atoms. Therefore, the permutation symmetry is an important factor which needs to be considered in the construction of PESs.35,94 Because the NN functions are not symmetric with respect to the permutation of identical atoms, the energy is not exactly continuous at a configuration with two or more equal internuclear distances,

Figure 2. Reaction probabilities for the H + CH4 → H2 + CH3 reaction on the XCZ PES and another PES fitted with 90% of data points included in fitting the XCZ PES. Reprinted with permission from ref 54. Copyright 2014 Chinese Physical Society.

which may introduce some errors in QCT calculations but can produce accurate results in quantum dynamics calculations. An exchange scheme to make the above NN PESs suitable for QCT calculations was proposed in this group, so that the conservation of the total energy can be substantially improved with the exchange scheme in the QCT calculations. Jiang and Guo proposed the PIP-NN approach to rigorously adapt the permutation symmetry in constructing the PESs of reactive systems.40,41 The PIP-NN approach employs a set of PIPs proposed by Bowman and co-workers instead of pairwise distances as input data in the NN. The PIP-NN approach first resolved the permutation problem in the three-atom reactions and was further successfully applied to four-atom and five-atom reactions. We also constructed PESs for the HOCO and CH5 systems using the PIP-NN method,119,120 with more ab initio data points included in the fitting. The new PIP-NN PESs are more suitable for QCT calculations, but the quantum dynamics results obtained on the PIP-NN PESs agree quite well with the NN PESs. II.C. FI-NN. The extra requirement of permutation symmetry in PES gives rise to the extra addition of functional basis used for fitting. This means the dimension of vector space with the transformed symmetrized basis should be equal to or larger than the original vector space before the transformation. Mathematically, the primary and secondary invariants can generate an invariant polynomial ring. The PIP method developed by Bowman and co-workers35 is generally based on this theorem. The invariant polynomials obtained from the monomial symmetrization was also developed by them.52 In the PIP-NN approach,40,41 the input pairwise internuclear distances of the NN is replaced by a set of PIPs, with all the polynomials truncated by a given degree. The degree for the truncation is the highest degree of the primary and secondary invariants. Nevertheless, the number of polynomials increases substantially with the degree bound. For example, the highest degree of secondary invariants for the A3B2 system is 11, resulting in 14984 invariant polynomials with degree 11. This number of invariant polynomials with the highest degree bound is huge, which is impossible to use in the input vector of a NN fitting. The number of polynomials can be reduced to 525 with a lower degree bound, such as six, but most of secondary polynomials 2293

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation will be discarded. Very recently, a new set of invariant polynomials was employed as the input vector of NN.42 Different from the invariant polynomials in PIP and also in PIPNN, the new set of polynomials contains the least number of invariants which can generate all the invariant polynomials. The new set of invariants is called fundamental invariants (FI),121 and the corresponding neural network is called a fundamental invariant neural network (FI-NN). For a molecular system AiBj···Xp, internuclear distances (r1, r2, ..., rn) are used as the variable set of the PESs in order to deal with the permutations of identical atoms. This choice of coordinates for permutational invariance was first accomplished by Bowman and co-workers.122 The permutation operator can be expressed as ĝ, which is an element of the direct product symmetric groups G = Si × Sj × ... × Sp, acting on the internuclear distances; here Sn is the symmetric group of degree n. The FIs can be calculated with King’s algorithm123 implemented in the computer algebra system called Singular.124 Here are the FI examples of molecular systems of A2B and A3B types. (1) For A2B molecules, the FIs calculated are f1 = x1 and

Table 1. Number of Fundamental Invariants (FIs), Invariant Polynomials Truncated by the Highest Degree of FIs, and the Speedup of FI-NN Relative to the Corresponding PIPNN system

FI

inv polynom

speedup

A2B A3 A2BC A2B2 A3B A4 A2BCD A2B2C A3BC A3B2 A4B A5

2,1 1,1,1 3,3 3,3,1 2,3,4 1,2,3,2,1 7,6 5,7,4 4,6,10 3,5,8,7,2,1 2,4,8,10,7 1,2,4,7,10,13,13,4,2

2,2 1,2,3 3,6 3,4,5 2,6,14 1,3,6,11,18 7,12 5,13,35 4,12,38 3,8,22,58,134,300 2,7,20,53,125 1,3,7,17,35,76,149,291,539

1.01 1.12 1.02 1.04 1.02 1.04 1.02 1.00 1.04 1.87 1.37 10.11

As discussed above, a total of nine FIs with a maximum degree of three are included for the A3B system. This number is considerably smaller than the number of PIP of 50 up to degree four in constructing a PIP-NN PES for the OH3 system. Therefore, a new PES for the OH3 system was constructed by the FI-NN method in this group. The same data set used in constructing the CXZ NN PES, i.e., a total of 16814 points computed at the UCCSD(T)-F12a/aVTZ level, was employed in fitting the FI-NN PES. The resulting overall RMSE and the maximal deviation in the FI-NN PES are 1.18 and 28.18 meV, respectively, while the RMSE and maximal deviation of the PIPNN PES are 1.4 and 36.3 meV. Since the fundamental invariants include three-body and four-body terms, the interaction between two fragments cannot exactly varnish in the reactant and product channels, which is also a numerical issue in the PIP-NN approach. We included many fragment data points in fitting these channels, resulting in very small fitting errors. The PESs in asymptote regions are reliable. A rigorous treatment is that the fit can be organized into a twobody fit; then one can discard the three-body terms unless at least three atoms are close and discard the four-body terms unless at least four atoms are close. This was done by Bowman and co-workers in several publications, where they ”purified” the fitting basis.126,127 We can further follow their approach to rigorously deal with the problem. In addition, we compared the fitting efficiency of the FI-NN and PIP-NN approaches. Table 2 lists the structures of model neural networks and the calculation time of 106 energy points on PIP-NN and FI-NN PESs. We can see from the speedup in the last column in Table 2 FI-NN PES is faster on evaluating than the corresponding PIP-NN PES, due to the smaller number of variables in the input layer. The weights in NN were randomly generated to guarantee the number of weights and neurons is nearly equal. No significant speedup was seen for small systems, but for the A3B2, A4B, and A5 systems, the speedup of FI-NN PESs is substantial, confirming the efficiency of the FI-NN approach. We anticipate the efficiency of the FINN approach can be further enhanced for larger systems. This FI-NN method is amenable to larger systems such as CH3CHO,49,128 and particularly important to floppy molecular systems with high permutation symmetry, such as the cases of H5+129 and H7+.130

f2 = x 2 + x3 f3 = x 2 2 + x32

(2) The FIs for the A3B system are given as follows. f1 = r1 + r2 + r4 f2 = r3 + r6 + r5 f3 = r12 + r22 + r42 f4 = r32 + r6 2 + r52 f5 = r1r3 + r2r3 + r2r6 + r4r6 + r4r5 + r1r5 f6 = r13 + r2 3 + r4 3 f7 = r33 + r6 3 + r53 f8 = r12r3 + r2 2r3 + r2 2r6 + r4 2r6 + r4 2r5 + r12r5 f9 = r32r4 + r1r6 2 + r2r52

The corresponding Fortran subroutines for the FIs of candidate systems up to five atoms are available from https://github.com/kjshao/FI. Table 1 gives the number of fundamental invariants, invariant polynomials truncated by the highest degree of FIs, and the speedup of FI-NN relative to the corresponding PIP-NN. Since the number of FI increases mildly, while the number of invariant polynomials increases very quickly with the degree of polynomial, the total number of FI is substantially smaller than that of invariant polynomials at the same degree, in particular for systems involving multiple identical atoms. The PIP-NN40,41 uses invariant polynomials as the input vector of NN, which are truncated by a highest degree bound to reduce the calculation.35,120,125 However, as the basis is incomplete, some extra error may be introduced into the PESs. The FI-NN employs FIs as the input vector of NN, which has the minimal number of polynomials. The network structure of FI-NN is much more flexible than PIP-NN, because it has less limitations in the input layer. 2294

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation Table 2. Comparison of Calculation Time (in Seconds) of 106 Energy Points on Model PIP-NN and FI-NN PESs system A2B A3 A2BC A2B2 A3B A4 A2BCD A2B2C A3BC A3B2 A4B A5

network structa

timeb

time of polynomc

time of NNd

3-30-30-1(1081) 4-30-29-1(1079) 3-30-30-1(1081) 6-45-15-1(1021) 7-30-30-1(1201) 10-25-35-1(1221) 7-30-30-1(1201) 12-25-35-1(1271) 9-30-30-1(1261) 22-20-40-1(1341) 9-30-30-1(1261) 22-20-40-1(1341) 13-23-37-1(1248) 19-20-40-1(1281) 16-30-30-1(1471) 53-15-40-1(1491) 20-40-30-1(2101) 54-20-50-1(2201) 26-40-30-1(2761) 525-4-75-1(2555) 31-40-40-1(2961) 207-10-70-1(2921) 56-27-33-1(2497) 1118-2-55-1(2459)

4.68 4.72 4.74 5.32 4.70 4.81 4.70 4.87 4.85 4.95 4.88 5.07 4.72 4.82 4.93 4.94 5.88 6.10 6.88 12.84 6.96 9.56 6.89 69.66

0.01 0.01 0.01 0.02 0.01 0.02 0.02 0.04 0.03 0.07 0.07 0.37 0.02 0.03 0.04 0.19 0.07 0.19 0.16 4.44 0.27 2.24 1.70 60.79

4.67 4.71 4.73 5.30 4.68 4.79 4.68 4.84 4.82 4.88 4.81 4.88 4.70 4.79 4.89 4.74 5.81 5.91 6.72 8.40 6.69 7.32 5.19 8.87

speedup 1.01 1.12

Figure 3. Jacobi coordinates (R,r1,r2,θ1,θ2,φ) for the diatom−diatom arrangement, and (R′,r′1,r′2,θ′1,θ′2,φ′) for the atom−triatom arrangement. Reprinted with permission from ref 77. Copyright 2012 AIP Publishing.

1.02 1.04 1.02 1.04 1.02 1.00 1.04 1.87

Figure 4. Comparison of (A) experimental and (B) theoretical differential cross-sections for the HD + OH → H2O + D reaction at the collision energy of 6.9 kcal/mol. From ref 23. Reprinted with permission from AAAS.

1.37 10.11

was the first time that quantum dynamics calculations reproduced the experimental DCSs for a four-atom reaction, indicating it is feasible to calculate complete dynamical information for some simple four-atom reactions without any dynamical approximation. This achievement represents a milestone in the quantum reactive scattering. In that study, the H2O product was found to be dominantly backward scattered with a large fraction of the available energy deposited into H2O internal excitation, which is consistent with a direct abstraction mechanism via a nearly collinear transition state. Subsequently, this RPD scheme based TDWP method was further applied to compute DCSs for some isotopically substituted reactions, HD + OH → H2O + D, D2 + OH → HOD + D, and H2 + OH → H2O + H.77,104,132,133 The reverse reaction H + H2O → H2 + OH has been widely studied as a prototype system for mode-specific chemistry. There are three vibrational modes of the water reagent, i.e., the symmetric stretching (ν1), bending (ν2), and asymmetric stretching (ν3) modes, respectively. It has proven to be an excellent candidate for demonstrating how the different vibrational modes of polyatomic reagents influence the reaction dynamics. A series of experimental studies on H + H2O and its isotopically substituted analogies reveal strong mode specificity and bond selectivity.134−140 Recently, the ISSWP method was employed to calculate the exact coupled-channel (CC) integral cross-sections (ICSs) for H + H 2O and H + HOD reactions141,142 on the YZCL2 and CXZ PESs, respectively. The reactivity enhancements from different initial vibrational states of H2O and HOD were obtained, including bending excited states, first and second stretching excited states, and simultaneous excitations of both bending and stretching modes, which suggest the strong mode specificity of H + H2O and bond selectivity of H + HOD as observed in experiment. The thermal rate constant and the contributions to this coefficient from individual vibrational states of H2O were also obtained

For each system, the first row is the network structure of FI-NN PES, the second row is the structure of PIP-NN PES. The values in parentheses are the number of weights for the corresponding network. b The total calculation time in the model PESs. cThe calculation time of polynomials in the model PESs. dThe calculation time of neural network in the model PESs. a

III. QUANTUM REACTIVE SCATTERING III.A. OH + H2↔ H + H2O. The H2 + OH ↔ H2O + H reactions and isotope analogues have become the prototypes for tetra-atomic reactions, which play an important role in both atmospheric chemistry and combustion.97,98 Since three of the four atoms in the reactive system are hydrogens, the OH3 system is an ideal candidate for high-quality electronic structure calculation of the PES as well as for an accurate quantum reactive scattering study. Much effort has been devoted to investigating the dynamics of this benchmark four-atom system at the state-to-state level.23,75−78,104,131−133 In 2011, Zhang and co-workers reported the full-dimensional DCSs for the HD + OH → H2O + D reaction on the XXZ PES, using the TDWP approach based on an improved version of the RPD scheme, i.e., multistep RPD.23,77 The basic strategy of the RPD scheme is to partition the full time-dependent (TD) wave function into a sum of reactant component (Ψr) and all product components (Ψp, p = 1, 2, 3, ...). Figure 3 shows the reagent Jacobi coordinates (R,r1,r2,θ1,θ2,φ) for the H2 + OH diatom−diatom arrangement and product Jacobi coordinates (R′,r′1,r′2,θ′1,θ′2,φ′) for the H + H2O atom−triatom arrangement, where r2 and r’2 share the same vector. Figure 4 shows the quantum dynamical DCSs and highresolution crossed molecular beam experimental results for the HD + OH → H2O + D reaction. Excellent agreement was achieved between the theoretical and experimental DCSs. This 2295

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation

nonreactive OH bond works well for the H + H2O → H2 + OH abstraction reaction and its reverse; however, both OH bonds should be treated as reactive bonds in order to accurately investigate the H′ + H2O → HOH′ + H exchange reaction, due to the saddle point close to a C3v geometry. In 2012, we reported a full-dimensional quantum dynamics study for the H + D2O → D + HOD reaction using the ISSWP approach on the YZCL2 PES, with both OD bonds in the D2O reactant treated as reactive bonds.145 This is the first time one can obtain CC results with two reactive bonds for a four-atom reaction. Due to the C3v minimum along the reaction path, a clear step-like feature is demonstrated in the CC cross-sections just above the threshold, presenting the existence of shape resonance in the reaction, as shown in Figure 6. The CC cross-

and compared with the available experimental data. The observed strong mode specificity is consistent with the extended Polanyi rules, indicating the vibrational excitation is supposed to be more efficient than the translational motion in promoting a reaction with a late barrier. Very recently, the H + H2O → H2 + OH reaction was investigated at the more detailed state-to-state quantum dynamical level. The full-dimensional state-to-state quantum calculation was carried out on the CXZ PES by using the multistep RPD scheme based TDWP method.143,144 The first full-dimensional quantum DCSs for this atom−triatomic reaction were reported with H2O initially in nonrotating (000), (100), and (001), where the influences of the initial vibrational excitation on the product state distribution and DCS were investigated.143 Figure 5 shows the total DCSs for the H +

Figure 5. (a) Differential cross-sections for the H + H2O → OH + H2 reaction for the three initial states at the collision energies of 0.8, 1.0, and 1.2 eV; (b) same as panel a, except at the total energies of 0.8, 1.1, and 1.4 eV. Reprinted from ref 143 with permission from the Royal Society of Chemistry.

H2O → H2 + OH reaction with H2O initially in the (000), (100), and (001) vibrational states at three collision energies and three total energies. As seen, the vibrational excitations dramatically enhance the DCSs, and one quantum in the symmetric or asymmetric stretching excitation of H2O has nearly identical effects on the DCS. The energy initially deposited in stretching vibrations is much more efficient than the translational energy in promoting the reaction but has effects on product angular distribution rather similar to those of the translational energy. In addition, the quantum dynamics results reveal that the reaction from (100) and (001) initial states of H2O produce vibrationally cold OH as the ground initial state, in agreement with the experimental observation of Zare and co-workers.139 This observation can be explained clearly by the local mode picture of H2O vibration and implies that the nonreacting OH does act as a spectator in the reaction. III.B. H′ + H2O → HOH′ + H. For the H + H2O → H2 + OH abstraction reaction and its reverse, the calculations were usually carried out on the basis that one of the OH bonds in the H2O reactant is basically a spectator bond that will not be cleaved or highly excited during the reaction. Therefore, one could use a very limited vibrational basis function for the bond, or sometimes just simply freeze it in some approximate studies in order to simplify the calculations. The treatment of

Figure 6. (a) CC and CS ICSs for the H + D2O → D + HOD exchange reaction; (b) same as panel a except in the low collision energy region. A distinct step-like feature is indicated in the CC curve, implying the signature of shape resonance. Reprinted from ref 145 with permission from the Royal Society of Chemistry.

sections show very good agreement with the experimental results and are a bit larger than the previous theoretical results, which reveals that the one reactive bond approximation and CS approximation together used in the previous study coincidentally did not lead to substantial errors in the ICSs. Similar behaviors for the CC cross-sections were found in later quantum dynamics studies for the H + H2O/HOD and D + H2O/HOD exchange reactions.146,147 III.C. OH + CO → H + CO2. Due to the involvement of heavy atoms of C and O, and HOCO intermediate complex, the OH + CO reaction presents a great challenge to quantum scattering calculations. Previous quantum dynamics was 2296

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation

experimental and theoretical interest for many years.152−155 Because five of the six atoms involved are hydrogens, this reaction has become a benchmark for developing and testing various theoretical methods for accurate studies of polyatomic chemical reactions. Due to the quantum nature of the reactive scattering problem, it is extremely difficult currently to treat such a sixatom reaction exactly in full dimension, although significant progress has been made by Manthe and co-workers on this direction.28,156−158 This group has been developing the TDWP methods to investigate the reaction of X + YCZ 3 type,24,26,29,48,159,160 employing an eight-dimensional (8D) model, which was originally proposed by Palma and Clary,161 by restricting the nonreacting CZ3 group in a C3v symmetry. Since the assumption holds very well in reality, the model has a quantitative level of accuracy. The eight degrees of freedom for the XYCZ3 system in the Jacobi coordinates (R, r, s, χ, θ1, φ1, θ2, φ2) are shown in Figure 8. Furthermore, the CH bond lengths in the nonreacting CH3

restricted to obtain the reaction probabilities for the total angular momentum J = 0.148 Full-dimensional ISSWP calculations were performed on the old Lakin−Troya− Schatz−Harding (LTSH) PES111 for many J partial waves using the CS approximation in this group.149 It was found the initial OH vibrational excitation substantially enhances the reactivity, while initial CO excitation has little effect on the reactivity. The CS approximation was shown to provide reasonably accurate total reaction probabilities for J > 0. The agreement with experimental rate constants is only qualitative, implying the possible inaccuracies of the LTSH PES. In 2011, the OH + CO → H + CO2 reaction was investigated at the state-to-state quantum mechanical level for the first time for total angular momentum J = 0. The fulldimensional dynamics calculation was performed on the LTSH PES for the ground and vibrationally excited initial states of OH and CO,76,150 using the multistep RPD method. It was found that the initial CO vibrational excitation essentially has the same effect on the product energy partition as the reagent translational motion, while the initial OH excitation leads to slightly more internal energy of CO2. Recently, the fulldimensional state-to-state quantum dynamics calculation for J = 0 on the most accurate PIP-NN PES was also carried out using the multistep RPD method.151 The quantum reaction probabilities for the ground initial rovibrational state on the PIP-NN PES are quite small in magnitude with many narrow and overlapping resonances. As indicated in Figure 7, the CO2

Figure 8. Eight-dimensional Jacobi coordinates for the X + YCZ3 system. Reprinted with permission from ref 48. Copyright 2011 AIP Publishing.

group can be fixed, and thus the number of degrees of freedom was reduced to seven. We carried out a seven-dimensional (7D) ISSWP study for the H + CD4 → HD + CD3 reaction on both the ZFWCZ and ZBB2 PESs, where good agreement was achieved between theory and experiment on the energy dependence of the ICS in a wide collision energy region.24 The influences of the CHD3 vibrational excitation on the H + CHD3 → H2 + CD3 reactivity were also investigated using the 7D ISSWP method on the ZFWCZ PES,159 which revealed that the CCH stretching excitation can promote the reaction dramatically. The thermal rate constants were obtained by taking into account the contributions from relevant initial vibration states. It was found that the ground initial state has a dominant contribution to the thermal rate constant at the lowtemperature region, while the relative contribution to the thermal rate constant from the vibrational excited states increases substantially with the increase of temperature. Furthermore, this group proposed a six-dimensional (6D) model to investigate the H + CH4 → H2 + CH3 reaction at the state-to-state quantum mechanical level,160 using the multistep RPD scheme on the XCZ PES. The 6D model was based on the quantum dynamics study of the H + CH4 → H2 + CH3 reaction that was carried out using the multistep RPD scheme on the XCZ PES,160 which was based on the 7D model by assuming that the CH3 group can rotate freely with respect to its C3v symmetry axis. The calculation indicated that the 6D treatment can produce essentially the same reaction probabilities as the 7D results. The product vibrational/rotational state distributions and product energy partitioning information were obtained for ground initial rovibrational state with the total angular momentum J = 0.

Figure 7. Product vibrational state distributions of CO2 for total angular momentum J = 0 at Ec = 0.1 (a) and 0.4 eV (b) obtained by the QM calculation for the OH + CO → H + CO2 reaction. Reprinted by permission from Springer-Verlag: Springer, Theor. Chem. Acc., ref 151, Copyright 2014.

products are dominantly excited in bending and symmetric modes, but essentially have no excitation in the antisymmetric stretching mode. Quite large population is also seen from the simultaneous excitations of combined bending and symmetric modes. In addition, the validity of the QCT method in describing the dynamics of this reaction was confirmed, except that the QCT calculations ignore important quantum effects, such as tunneling, zero-point energy, and resonances. III.D. H + CH4↔ H2 + CH3. The H + CH4 ↔ H2 + CH3 reaction and its reverse have been the subject of both 2297

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation

Very recently, Zhang and co-workers reported an accurate quantum dynamics study of the H + CH4 substitution reaction and its isotope analogues, employing the above 7D quantum model on the accurate ab initio LCZXZG PES.29 The calculations reveal that the reaction exhibits a strong normal secondary isotope effect on the ICSs measured above the reaction threshold, and a small but reverse secondary kinetic isotope effect at room temperature, as shown in Figure 10. The

The quantum dynamics study for this six-atom reaction without the CS approximation was first achieved in this group. The coupled channel (CC) results for the H + CHD3 → H2 + CD3 reaction using the 7D model was calculated on the XCZ PES. 162,163 It was found that the CS approximation considerably underestimates the CC results, and the initial rotational excitation of CHD3 up to J0 = 2 does not have any effect on the total cross-sections. Very recently, an 8D quantum dynamical approach based on Palma and Clary’s model was proposed for the H2 + CH3 → H + CH4 reaction,164,165 which is similar to that for its reverse reaction. The new 8D Hamiltonian was employed to treat the CH3 group, which is in a simple form and easy to implement. The reaction probabilities as well as the integral cross-sections were calculated for CH3 initially in the ground and various vibrationally excited initial states. The influences of vibrational excitations of both reagents were investigated for this reaction. It was found that the excitation of CH stretching has a negligible effect on the reactivity, and the 7D quantum model with the CH bond length fixed works very well for this reaction, as shown in Figure 9. The H2 vibrational excitation can

Figure 10. (a) ICSs for the T + CH4/CD4 substitution reactions as a function of collision energy; (b) same as panel a except for the H + CD4/CH3D → D + CD3H/CH4 reactions.

reaction proceeds along a path with a considerably higher barrier height than the static barrier, as the umbrella angle of the nonreacting CH3 group cannot change synchronously with the other reaction coordinates during the reaction owing to insufficient energy transfer from the translational motion to the umbrella mode. Those investigations provided unprecedented dynamical details for this simplest Walden inversion reaction and also shed valuable light on the dynamics of gas-phase SN2 reactions. III.F. Cl + CH4→ HCl + CH3. The Polanyi rules give the effects of different forms of energy in the atom−diatom reaction with different barrier locations. For a reaction with an early barrier, translational energy is more effective than vibrational energy in overcoming the barrier. On the other hand, vibrational excitation has a higher efficacy in enhancing the reactivity than the translational energy with a late barrier. Unlike the atom−diatom reaction, where only one vibration is invoked in the reactant, the reaction of polyatomic species involves multiple types of vibrational motion. It is unclear whether different vibrational modes are equivalent in their capacity to promote the reactivity. The extension of the Polanyi rules to reactions involving polyatomic molecules has drawn extensive interests in recent years, such as the Cl + CH4 → HCl + CH3 reaction. The Cl + CH4 → HCl + CH3 reaction has a product-like barrier (later barrier), and thus according to the extended

Figure 9. Reaction probabilities at the total angular momentum J = 0 for the H2 + CH3 → H + CH4 reaction, with reactants initially in different vibrational states (vH2, vu) and (vH2, vu,vs) for 7D and 8D, respectively. Reprinted with permission from ref 165. Copyright 2015 American Chemical Society.

promote the reaction but is less effective than the translational motion in low-energy regions. In contrast, the reaction threshold is largely reduced by exciting the first umbrella mode. The computed rate constants agree well with the available experimental results and previous theoretical calculations. III.E. H′ + CH4 → H + CH′H3. In addition to the H + CH4 → H2 + CH3 abstraction reaction, there exists a substitution reaction, H′ + CH4 → H + CH′H3, with a D3h transition state and a static barrier height of 1.6 eV. It is the simplest reaction proceeding through the back-side attack Walden inversion mechanism, very similar to the gas-phase SN2 reactions with central barriers, except that the SN2 reactions have pre- and postreaction wells arising from strong ion−dipole interaction between reagents and products. 2298

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation

more rotational states of product, they found the CH stretch excitation becomes more efficacious than the same amount of translational energy in promoting the reaction, which agrees reasonably well with the 7D quantum dynamics calculation.

Polanyi rules, the vibrational excitation is supposed to be more efficient than the translational motion in promoting the reaction. Surprisingly, crossed molecular beam experiments carried out by Liu and co-workers for the Cl + CHD3 reaction observed that the ICS for the CD3 (v = 0) product for the initial CH stretch excited state is smaller than that at low translational energies and becomes comparable at higher translational energies to the corresponding one for the ground initial state.166,167 Based on their estimation that most of CD3 product resulted in the ground vibrational state (v = 0), they revealed that the stretch excitation is no more effective or slightly more effective in promoting the late-barrier Cl + CHD3 reaction than an equivalent amount of translational energy, which is inconsistent with the Polanyi rules. In 2012, Zhang and co-workers performed a 7D quantum dynamics study based on the Palma−Clary model for the Cl + CHD3 reaction to investigate the influences of vibrational excitation on the reactivity.26 A high-quality, full-dimensional global PES developed by Czakó and Bowman (CB PES)168 for the reactive system based on accurate ab initio calculations was employed in the quantum dynamics calculations. The ICSs for different vibrational states were obtained based on the CS approximation on calculating the reaction probabilities for J > 0. The quantum dynamics results revealed that the vibrational enhancement by the C−H stretch excitation is actually very significant, except at very low collision energies, as indicated in Figure 11. After a careful reexamination of the experimental data, it was revealed that the discrepancy between the original experimental data and 7D quantum dynamical results arises from the incompleteness of the measurement.169 By including

IV. CONCLUSIONS We have reviewed methods developed in our group to obtain accurate potential energy surfaces and perform quantum reactive scattering calculations in high dimensions for polyatomic chemical reactions beyond three atoms in recent years. The TDWP methodologies as well as the construction schemes of highly accurate PESs were developed to achieve the quantum dynamical information and mechanism at the most fundamental level, which can be directly compared with the experiment. Illustrations of PESs and reaction dynamics were given to the H2 + OH, H + H2O, HO + CO, and H/Cl + CH4 bimolecular reactions. Particularly, we obtained the accurate state-to-state DCSs for the benchmark four-atom reactive system OH3, indicating quantum reactive scattering problems for many four-atom reactions can be finally resolved. In addition, accurate quantum dynamics calculations for some sixatom reactions of X + YCZ3 type were first accomplished in our group. Quantum reactive scattering calculations on reliable PESs have revealed and elucidated the dependence of collision energy, reactant rovibrational state, mode specificity, and isotope effects on detailed dynamical information, such as the reaction probability, product vibrational state distribution, and reactive resonances. Comparisons with experiments help us gain a deeper understanding of reaction dynamics.



AUTHOR INFORMATION

Corresponding Authors

*(B.F.) E-mail: [email protected]. *(D.H.Z.) E-mail: [email protected]. ORCID

Bina Fu: 0000-0003-1568-0259 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work was supported by the National Natural Science Foundation of China (Grant Nos. 21722307, 21673233, 21590804, 21433009, and 21688102), the Strategic Priority Research Program (Grant No. XDB17000000), and the Youth Innovation Promotion Association (Grant No. 2015143) of the Chinese Academy of Sciences.



REFERENCES

(1) Bowman, J. M.; Schatz, G. C. Theoretical studies of polyatomic bimolecular reaction dynamics. Annu. Rev. Phys. Chem. 1995, 46, 169− 195. (2) Zhang, J. Z. H.; Dai, J.; Zhu, W. Development of Accurate Quantum Dynamical Methods for Tetraatomic Reactions. J. Phys. Chem. A 1997, 101, 2746−2754. (3) Miller, M. H. Quantum and semiclassical Green’s functions in chemical reaction dynamics. J. Chem. Soc., Faraday Trans. 1997, 93, 685−690. (4) Clary, D. C. Quantum theory of chemical reaction dynamics. Science 1998, 279, 1879−1882. (5) Althorpe, S. C.; Clary, D. C. Quantum scattering calculations on chemical reactions. Annu. Rev. Phys. Chem. 2003, 54, 493−529.

Figure 11. (a) Total ICS and the cross-sections with product CD3(v = 0) for the ground initial state by QM calculations, in comparison with the standard QCT and ZPE constrained ICSs, for the Cl + CHD3 → HCl + CD3 reaction; (b) same as panel a, except for the initial CH excited state. The ICS for the ground initial state is also shown for direct comparison. Reprinted with permission from ref 26. Copyright 2012 American Chemical Society. 2299

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation

(28) Welsch, R.; Manthe, U. Loss of Memory in H + CH4 →H2 + CH3 State-to-State Reactive Scattering. J. Phys. Chem. Lett. 2015, 6, 338−342. (29) Zhao, Z.; Zhang, Z.; Liu, S.; Zhang, D. H. Dynamical barrier and isotope effects in the simplest substitution reaction via Walden inversion mechanism. Nat. Commun. 2017, 8, 14506. (30) Fu, B.; Shan, X.; Zhang, D. H.; Clary, D. C. Recent advances in quantum scattering calculations on polyatomic bimolecular reactions. Chem. Soc. Rev. 2017, 46, 7625−7649. (31) Pan, H.; Liu, K.; Caracciolo, A.; Casavecchia, P. Crossed beam polyatomic reaction dynamics: recent advances and new insights. Chem. Soc. Rev. 2017, 46, 7517−7547. (32) Zhao, B.; Guo, H. State-to-state quantum reactive scattering in four-atom systems. WIRES: Comput. Mol. Sci. 2017, 7, e1301−e1301. (33) Dawes, R.; Thompson, D. L.; Wagner, A. F.; Minkoff, M. Interpolating moving least-squares methods for fitting potential energy surfaces: a strategy for efficient automatic data point placement in high dimensions. J. Chem. Phys. 2008, 128, 084107. (34) Collins, M. A. Molecular potential-energy surfaces for chemical reaction dynamics. Theor. Chem. Acc. 2002, 108, 313−324. (35) Braams, B. J.; Bowman, J. M. Permutationally invariant potential energy surfaces in high dimensionality. Int. Rev. Phys. Chem. 2009, 28, 577−606. (36) Bowman, J. M.; Czakó, G.; Fu, B. High-dimensional ab initio potential energy surfaces for reaction dynamics calculations. Phys. Chem. Chem. Phys. 2011, 13, 8094−8111. (37) Raff, L. M.; Malshe, M.; Hagan, M.; Doughan, D. I.; Rockley, M. G.; Komanduri, R. Ab initio potential-energy surfaces for complex, multichannel systems using modified novelty sampling and feedforward neural networks. J. Chem. Phys. 2005, 122, 084104. (38) Chen, J.; Xu, X.; Xu, X.; Zhang, D. H. A global potential energy surface for the H2 + OH →H2O + H reaction using neural networks. J. Chem. Phys. 2013, 138, 154301. (39) Chen, J.; Xu, X.; Xu, X.; Zhang, D. H. Communication: An accurate global potential energy surface for the OH + CO →H + CO2 reaction using neural networks. J. Chem. Phys. 2013, 138, 221104. (40) Jiang, B.; Guo, H. Permutation invariant polynomial neural network approach to fitting potential energy surfaces. J. Chem. Phys. 2013, 139, 054112. (41) Li, J.; Jiang, B.; Guo, H. Permutation invariant polynomial neural network approach to fitting potential energy surfaces. II. Fouratom systems. J. Chem. Phys. 2013, 139, 204103. (42) Shao, K.; Chen, J.; Zhao, Z.; Zhang, D. H. Communication: Fitting potential energy surfaces with fundamental invariant neural network. J. Chem. Phys. 2016, 145, 071101. (43) Jiang, B.; Li, J.; Guo, H. Potential energy surfaces from high fidelity fitting of ab initio points: the permutation invariant polynomial - neural network approach. Int. Rev. Phys. Chem. 2016, 35, 479−506. (44) Bartlett, R. J.; Musiał, M. Coupled-cluster theory in quantum chemistry. Rev. Mod. Phys. 2007, 79, 291−352. (45) Werner, H.-J. Matrix-formulated direct multiconfiguration selfconsistent field and multiconfiguration reference configurationinteraction methods. Adv. Chem. Phys. 2007, 69, 1−62. (46) Schatz, G. C. The analytical representation of electronica potential-ennergy surfaces. Rev. Mod. Phys. 1989, 61, 669−688. (47) Hollebeek, T.; Ho, T. S.; Rabitz, H. Constructing multidimensional molecular potential energy surfaces from ab initio data. Annu. Rev. Phys. Chem. 1999, 50, 537−570. (48) Zhou, Y.; Fu, B.; Wang, C.; Collins, M. A.; Zhang, D. H. Ab initio potential energy surface and quantum dynamics for the H + CH4 →H2 + CH3 reaction. J. Chem. Phys. 2011, 134, 064323. (49) Fu, B.; Han, Y.-C.; Bowman, J. M.; Angelucci, L.; Balucani, N.; Leonori, F.; Casavecchia, P. Intersystem crossing and dynamics in O(3P)+C2H4 multichannel reaction: Experiment validates theory. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 9733−9738. (50) Yang, J.; Shao, K.; Zhang, D.; Shuai, Q.; Fu, B.; Zhang, D. H.; Yang, X. Trapped Abstraction in the O(1D) + CHD3 →OH + CD3 Reaction. J. Phys. Chem. Lett. 2014, 5, 3106−3111.

(6) Balint-Kurti, G. G. Time-dependent and time-independent wavepacket approaches to reactive scattering and photodissociation dynamics. Int. Rev. Phys. Chem. 2008, 27, 507−539. (7) Guo, H. Quantum dynamics of complex-forming bimolecular reactions. Int. Rev. Phys. Chem. 2012, 31, 1−68. (8) Yarkony, D. R. Nonadiabatic quantum chemistry: past, present and future. Chem. Rev. 2012, 112, 481−498. (9) Nyman, G.; Yu, H.-G. Quantum approaches to polyatomic reaction dynamics. Int. Rev. Phys. Chem. 2013, 32, 39−95. (10) Otto, R.; Ma, J.; Ray, A. W.; Daluz, J. S.; Li, J.; Guo, H.; Continetti, R. E. Imaging Dynamics on the F + H2O → HF + OH Potential Energy Surfaces from Wells to Barriers. Science 2014, 343, 396−399. (11) Liu, K. Crossed-beam studies of neutral reactions: State-Specific Differential Cross Sections. Annu. Rev. Phys. Chem. 2001, 52, 139−164. (12) Fernandez-Alonso, F.; Zare, R. N. Scattering resonances in the simplest chemical reaction. Annu. Rev. Phys. Chem. 2002, 53, 67−99. (13) Balucani, N.; Capozza, G.; Leonori, F.; Segoloni, E.; Casavecchia, P. Crossed molecular beam reactive scattering: from simple triatomic to multichannel polyatomic reactions. Int. Rev. Phys. Chem. 2006, 25, 109−163. (14) Yang, X. State-to-State Dynamics of Elementary Bimolecular Reactions. Annu. Rev. Phys. Chem. 2007, 58, 433−459. (15) Crim, F. F. Chemical dynamics of vibrationally excited molecules: Controlling reactions in gases and on surfaces. Proc. Natl. Acad. Sci. U. S. A. 2008, 105, 12654−12661. (16) Liu, K. Quantum dynamical resonances in chemical reactions: from A + BC to polyatomic systems. Advances in Chemical Physics 2012, 149, 1−46. (17) Guo, H.; Liu, K. Control of chemical reactivity by transitionstate and beyond. Chemical Science 2016, 7, 3992−4003. (18) Qiu, M.; Ren, Z.; Che, L.; Dai, D.; Harich, S. A.; Wang, X.; Yang, X.; Xu, C.; Xie, D.; Gustafsson, M.; Skodje, R. T.; Sun, Z.; Zhang, D. H. Observation of Feshbach Resonances in the F + H2 → HF + H Reaction. Science 2006, 311, 1440−1443. (19) Che, L.; Ren, Z.; Wang, X.; Dong, W.; Dai, D.; Wang, X.; Zhang, D. H.; Yang, X.; Sheng, L.; Li, G.; Werner, H.-J.; Lique, F.; Alexander, M. H. Breakdown of the Born-Oppenheimer Approximation in the F+ o-D2 → DF + D Reaction. Science 2007, 317, 1061−1064. (20) Wang, X.; et al. The Extent of Non−Born-Oppenheimer Coupling in the Reaction of Cl2(P) with para-H2. Science 2008, 322, 573−576. (21) Wang, T.; Chen, J.; Yang, T.; Xiao, C.; Sun, Z.; Huang, L.; Dai, D.; Yang, X.; Zhang, D. H. Dynamical Resonances Accessible Only by Reagent Vibrational Excitation in the F + HD→HF + D Reaction. Science 2013, 342, 1499−1502. (22) Yang, T.; Chen, J.; Huang, L.; Wang, T.; Xiao, C.; Sun, Z.; Dai, D.; Yang, X.; Zhang, D. H. Extremely short-lived reaction resonances in Cl + HD (v = 1) → DCl + H due to chemical bond softening. Science 2015, 347, 60−63. (23) Xiao, C.; Xu, X.; Liu, S.; Wang, T.; Dong, W.; Yang, T.; Sun, Z.; Dai, D.; Xu, X.; Zhang, D. H.; Yang, X. Experimental and Theoretical Differential Cross Sections for a Four-Atom Reaction: HD + OH → H2O + D. Science 2011, 333, 440−442. (24) Zhang, W.; et al. Depression of reactivity by the collision energy in the single barrier H + CD4 →HD + CD3 reaction. Proc. Natl. Acad. Sci. U. S. A. 2010, 107, 12782−12785. (25) Zhang, D. H.; Guo, H. Recent Advances in Quantum Dynamics of Bimolecular Reactions. Annu. Rev. Phys. Chem. 2016, 67, 135−158. (26) Zhang, Z.; Zhou, Y.; Zhang, D. H.; Czakó, G.; Bowman, J. M. Theoretical Study of the Validity of the Polanyi Rules for the LateBarrier Cl + CHD3 Reaction. J. Phys. Chem. Lett. 2012, 3, 3416−3419. (27) Liu, R.; Yang, M.; Czakó, G.; Bowman, J. M.; Li, J.; Guo, H. Mode Selectivity for a central Barrier Reaction: Eight-Dimensional Quantum Studies of the O3(P) + CH4 →OH + CH3 Reaction on an Ab Initio Potential Energy Surface. J. Phys. Chem. Lett. 2012, 3, 3776− 3780. 2300

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation (51) Shao, K.; Fu, B.; Zhang, D. H. A global full-dimensional potential energy surface and quasiclassical trajectory study of the O(1D) + CH4 multichannel reaction. Phys. Chem. Chem. Phys. 2015, 17, 24098−24107. (52) Xie, Z.; Bowman, J. M. Permutationally Invariant Polynomial Basis for Molecular Energy Surface Fitting via Monomial Symmetrization. J. Chem. Theory Comput. 2010, 6, 26−34. (53) Raff, L. M.; Komanduri, R.; Hagan, M.; Bukkapatnam, S. Neural networks in chemical reaction dynamics; Oxford University Press: New York, 2012. (54) Xu, X.; Chen, J.; Zhang, D. H. Global Potential Energy Surface for the H+CH4→H2+CH3 Reaction using Neural Networks. Chin. J. Chem. Phys. 2014, 27, 373−379. (55) Kuppermann, A.; Hipes, P. G. Three-dimensional quantum mechanical reactive scattering using symmetrized hyperspherical coordinates. J. Chem. Phys. 1986, 84, 5962−5964. (56) Pack, R. T.; Parker, G. A. Quantum reactive scattering in three dimensions using hyperspherical (APH) coordinates. J. Chem. Phys. 1987, 87, 3888−3921. (57) Skouteris, D.; Castillo, J.; Manolopoulos, D. ABC: a quantum reactive scattering program. Comput. Phys. Commun. 2000, 133, 128− 135. (58) Aquilanti, V.; Cavalli, S. The quantum-mechanical Hamiltonian for tetraatomic systems in symmetric hyperspherical coordinates. J. Chem. Soc., Faraday Trans. 1997, 93, 801−809. (59) Pogrebnya, S. K.; Echave, J.; Clary, D. C. Quantum theory of four-atom reactions using arrangement channel hyperspherical coordinates: Formulation and application to OH+H2→H2O+H. J. Chem. Phys. 1997, 107, 8975−8984. (60) Kosloff, R. Time-dependent quantum-mechanical methods for molecular dynamics. J. Phys. Chem. 1988, 92, 2087−2100. (61) Zhang, J. Z. H. Theory and Application of Quantum Molecular Dynamics; World Scientific, Singapore, 1998. (62) Neuhauser, D.; Baer, M.; Judson, R. S.; Kouri, D. J. The application of time-dependent wavepacket methods to reactive scattering. Comput. Phys. Commun. 1991, 63, 460−481. (63) Sun, Z.; Yang, W.; Zhang, D. H. Higher-order split operator schemes for solving the Schrodinger equation in the time-dependent wave packet method: applications to triatomic reactive scattering calculations. Phys. Chem. Chem. Phys. 2012, 14, 1827−1845. (64) Mandelshtam, V. A.; Taylor, H. S. A simple recursion polynomial expansion of the Green function with absorbing boundary conditions. Application to the reactive scattering. J. Chem. Phys. 1995, 103, 2903−2907. (65) Gray, S. K.; Balint-Kurti, G. G. Quantum dynamics with real wave packets, including application to three-dimensional (J = 0)D +H2D+H reactive scattering. J. Chem. Phys. 1998, 108, 950−962. (66) Chen, R.; Guo, H. Evolution of quantum system in order domain of Chebyshev operator. J. Chem. Phys. 1996, 105, 3569−3578. (67) Gray, S. K.; Balint-Kurti, G. G. Quantum dynamics with real wave packets, including application to three-dimensional (J = 0)D+H2 →DH+H reactive scattering. J. Chem. Phys. 1998, 108, 950−962. (68) Althorpe, S. C. Quantum wavepacket method for state-to-state reactive cross sections. J. Chem. Phys. 2001, 114, 1601−1616. (69) Yuan, K.; Cheng, Y.; Liu, X.; Harich, S.; Yang, X.; Zhang, D. H. Experimental and Quantum Dynamical Study on an Asymmetric Insertion Reaction: State-to-State Dynamics of O( 1 D)+HD(1Σg+,v’=0,j’=0)→OH(2Π,v’,N’)+D(2S). Phys. Rev. Lett. 2006, 96, 103202. (70) Lin, S. Y.; Guo, H. Quantum state-to-state cross sections for atom-diatom reactions: A Chebyshev real wave-packet approach. Phys. Rev. A: At., Mol., Opt. Phys. 2006, 74, 022703. (71) Hankel, M.; Smith, S. C.; Allan, R. J.; Gray, S. K.; Balint-Kurti, G. G. State-to-state reactive differential cross sections for the H+H2→ H2+H reaction on five different potential energy surfaces employing a new quantum wavepacket computer code: DIFFREALWAVE. J. Chem. Phys. 2006, 125, 164303.

(72) Gómez-Carrasco, S.; Roncero, O. Coordinate transformation methods to calculate state-to-state reaction probabilities with wave packet treatments. J. Chem. Phys. 2006, 125, 054102. (73) Sun, Z.; Guo, H.; Zhang, D. H. Extraction of state-to-state reactive scattering attributes from wave packet in reactant Jacobi coordinates. J. Chem. Phys. 2010, 132, 084112. (74) Peng, T.; Zhang, J. Z. H. A reactant-product decoupling method for state-to-state reactive scattering. J. Chem. Phys. 1996, 105, 6072− 6074. (75) Cvitaŝ, M. T.; Althorpe, S. C. State-to-state reactive scattering in six dimensions using reactant-product decoupling: OH + H2 →H2O + H (J = 0). J. Chem. Phys. 2011, 134, 024309. (76) Liu, S.; Xu, X.; Zhang, D. H. Communication: State-to-state quantum dynamics study of the OH + CO →H + CO2 reaction in full dimensions (J = 0). J. Chem. Phys. 2011, 135, 141108. (77) Liu, S.; Xu, X.; Zhang, D. H. Time-dependent wave packet theory for state-to-state differential cross sections of four-atom reactions in full dimensions: Application to the HD + OH →H2O + D reaction. J. Chem. Phys. 2012, 136, 144302. (78) Cvitas, M. T.; Althorpe, S. C. A Chebyshev method for state-tostate reactive scattering using reactant-product decoupling: OH + H2 →H2O + H. J. Chem. Phys. 2013, 139, 064307. (79) Welsch, R.; Huarte-Larranãga, F.; Manthe, U. State-to-state reaction probabilities within the quantum transition state framework. J. Chem. Phys. 2012, 136, 064117. (80) Manthe, U.; Welsch, R. Correlation functions for fully or partially state-resolved reactive scattering calculations. J. Chem. Phys. 2014, 140, 244113. (81) Miller, W. H. Quantum mechanical transition state theory and a new semiclassical model for reaction rate constants. J. Chem. Phys. 1974, 61, 1823−1834. (82) Zhao, B.; Sun, Z.; Guo, H. Calculation of state-to-state differential and integral cross sections for atom-diatom reactions with transition-state wave packets. J. Chem. Phys. 2014, 140, 234110. (83) Zhao, B.; Sun, Z.; Guo, H. State-to-State Mode Specificity: Energy Sequestration and Flow Gated by Transition State. J. Am. Chem. Soc. 2015, 137, 15964−15970. (84) Zhao, B.; Guo, H. Modulations of Transition-State Control of State-to-State Dynamics in the F + H2O →HF + OH Reaction. J. Phys. Chem. Lett. 2015, 6, 676−680. (85) Zhao, B.; Sun, Z.; Guo, H. A reactant-coordinate-based wave packet method for full-dimensional state-to-state quantum dynamics of tetra-atomic reactions: Application to both the abstraction and exchange channels in the H + H2O reaction. J. Chem. Phys. 2016, 144, 064104. (86) Zhao, B.; Sun, Z.; Guo, H. State-to-state mode selectivity in the HD + OH reaction: Perspectives from two product channels. J. Chem. Phys. 2016, 144, 214303. (87) Zhao, B.; Sun, Z.; Guo, H. A reactant-coordinate-based approach to state-to-state differential cross sections for tetratomic reactions. J. Chem. Phys. 2016, 145, 184106. (88) Blank, T. B.; Brown, S. D.; Calhoun, A. W.; Doren, D. J. Neural network models of potential energy surfaces. J. Chem. Phys. 1995, 103, 4129−4137. (89) Brown, D. F. R.; Gibbs, M. N.; Clary, D. C. Combining ab initio computations, neural networks, and diffusion Monte Carlo: an efficient method to treat weakly bound molecules. J. Chem. Phys. 1996, 105, 7597−7604. (90) Manzhos, S.; Wang, X.; Dawes, R.; Carrington, J. Tucker A Nested Molecule-Independent Neural Network Approach for HighQuality Potential Fits. J. Phys. Chem. A 2006, 110, 5295. (91) Behler, J.; Parrinello, M. Generalized Neural-Network Representation of High-Dimensional Potential-Energy Surfaces. Phys. Rev. Lett. 2007, 98, 146401. (92) Pukrittayakamee, A.; Malshe, M.; Hagan, M.; Raff, L. M.; Narulkar, R.; Bukkapatnum, S.; Komanduri, R. Simultaneous fitting of a potential-energy surface and its corresponding force fields using feedforward neural networks. J. Chem. Phys. 2009, 130, 134101. 2301

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation

Analysis of the Kinetic Isotope Effects Using Rectilinear and Curvilinear Coordinates. J. Phys. Chem. 1996, 100, 16561−16567. (116) Wu, T.; Werner, H.-J.; Manthe, U. First-Principles Theory for the H + CH4 → H2 + CH3 Reaction. Science 2004, 306, 2227−2229. (117) Zhang, X.; Braams, B. J.; Bowman, J. M. An ab initio potential surface describing abstraction and exchange for H+CH4. J. Chem. Phys. 2006, 124, 021104. (118) Xie, Z.; Bowman, J. M.; Zhang, X. Quasiclassical trajectory study of the reaction H+CH4(v3=0,1)→CH3+H2 using a new ab initio potential energy surface. J. Chem. Phys. 2006, 125, 133120. (119) Li, J.; Chen, J.; Zhang, D. H.; Guo, H. Quantum and quasiclassical dynamics of the OH + CO →H + CO2 reaction on a new permutationally invariant neural network potential energy surface. J. Chem. Phys. 2014, 140, 044327. (120) Li, J.; Chen, J.; Zhao, Z.; Xie, D.; Zhang, D. H.; Guo, H. A permutationally invariant full-dimensional ab initio potential energy surface for the abstraction and exchange channels of the H + CH4 system. J. Chem. Phys. 2015, 142, 204302. (121) Derksen, H.; Kemper, G. Computational Invariant Theory. Encyclopaedia of Mathematical Sciences, 2nd ed.; Springer: Berlin, 2015. (122) Brown, A.; Braams, B. J.; Christoffel, K.; Jin, Z.; Bowman, J. M. Classical and quasiclassical spectral analysis of CH5+ using an ab initio potential energy surface. J. Chem. Phys. 2003, 119, 8790−8793. (123) King, S. A. Minimal generating sets of non-modular invariant rings of finite groups. Journal of Symbolic Computation 2013, 48, 101− 109. (124) Decker, W.; Greuel, G.-M.; Pfister, G.; Schönemann, H. Singular 4-0-2: A computer algebra system for polynomial computations. http://www.singular.uni-kl.de, 2015. (125) Li, J.; Guo, H. Communication: An accurate full 15 dimensional permutationally invariant potential energy surface for the OH + CH4 → H2O + CH3 reaction. J. Chem. Phys. 2015, 143, 221103. (126) Conte, R.; Qu, C.; Bowman, J. M. Permutationally Invariant Fitting of Many-Body, Non-covalent Interactions with Application to Three-Body Methane-water-water. J. Chem. Theory Comput. 2015, 11, 1631−1638. (127) Yu, Q.; Bowman, J. M. Ab Initio Potential for H3O+ → H+ + H2O: A Step to a Many-Body Representation of the Hydrated Proton? J. Chem. Theory Comput. 2016, 12, 5284−5292. (128) Fu, B.; Han, Y.; Bowman, J. M. Three-state surface hopping calculations of acetaldehyde photodissociation to CH3 + HCO on ab initio potential surfaces. Faraday Discuss. 2012, 157, 27−39. (129) Xie, Z.; Braams, B. J.; Bowman, J. M. Ab initio global potentialenergy surface for H5+→H3++H2. J. Chem. Phys. 2005, 122, 224307. (130) Qu, C.; Prosmiti, R.; Bowman, J. M. MULTIMODE calculations of the infrared spectra of H7+ and D7+ using ab initio potential energy and dipole moment surfaces. Theor. Chem. Acc. 2013, 132, 1413. (131) Cvitaŝ, M. T.; Althorpe, S. C. Quantum Wave Packet Method for State-to-State Reactive Scattering Calculations on AB + CD → ABC + D Reactions. J. Phys. Chem. A 2009, 113, 4557−4569. (132) Zhao, Z.; Liu, S.; Zhang, D. H. State-to-state differential cross sections for a four-atom reaction: H2 + OH →H2O + H in full dimensions. J. Chem. Phys. 2016, 145, 134301. (133) Zhao, B.; Sun, Z.; Guo, H. State-to-state differential cross sections for D2 + OH →D + DOH reaction: Influence of vibrational excitation of OH reactant. J. Chem. Phys. 2016, 145, 134308. (134) Crim, F. F. Vibrational State Control of Bimolecular Reactions: Discovering and Directing the Chemistry. Acc. Chem. Res. 1999, 32, 877−884. (135) Sinha, A.; Hsiao, M. C.; Crim, F. F. Bond-selected bimolecular chemistry: H+HOD(4vOH)→OD+H2. J. Chem. Phys. 1990, 92, 6333. (136) Sinha, A.; Hsiao, M. C.; Crim, F. F. J. Chem. Phys. 1991, 94, 4928. (137) Metz, R. B.; Thoemke, J. D.; Pfeiffer, J. M.; Crim, F. F. Selectively breaking either bond in the bimolecular reaction of HOD with hydrogen atoms. J. Chem. Phys. 1993, 99, 1744.

(93) Handley, C. M.; Popelier, P. L. A. Potential Energy Surfaces Fitted by Artificial Neural Networks. J. Phys. Chem. A 2010, 114, 3371. (94) Behler, J. Neural network potential-energy surfaces in chemistry: a tool for large-scale simulations. Phys. Chem. Chem. Phys. 2011, 13, 17930. (95) Hagan, M. T.; Menhaj, M. B. Training feedforward networks with the Marquardt algorithm. IEEE Trans. Neural Netw. 1994, 5, 989. (96) Sarle, W. Stopped Training and Other Remedies for Overfitting. Proceedings of the 27th Symposium on the Interface of Computing Science and Statistics; Interface Foundation of North America: Fairfax Station, VA, USA, 1995; p 352. (97) Clary, D. C. Theory of Reactive Collisions at Low Temperatures. In Rate Coefficients in Astrochemistry; Millar, T. J., Williams, D. A., Eds.; Astrophysics and Space Science Library; Springer: Dordrecht, The Netherlands, 1988; Vol. 146, Chapter 1, pp 1−16. (98) Warnatz, J. In Combustion Chemistry; Gardiner, W. C., Jr., Ed.; Springer-Verlag: New York, 1985; Chapter 5. (99) Schatz, G. C.; Elgersma, H. A quasi-classical trajectory study of product vibrational distributions in the OH + H2 →H2O + H reaction. Chem. Phys. Lett. 1980, 73, 21−25. (100) Ochoa de Aspuru, G.; Clary, D. C. New Potential Energy Function for Four-Atom Reactions. Application to OH + H2. J. Phys. Chem. A 1998, 102, 9631−9637. (101) Pogrebnya, S. K.; Palma, J.; Clary, D. C.; Echave, J. Quantum scattering and quasi-classical trajectory calculations for the H2+OH→ H2O+H reaction on a new potential surface. Phys. Chem. Chem. Phys. 2000, 2, 693−700. (102) Wu, G.-S.; Schatz, G. C.; Lendvay, G.; Fang, D.-C.; Harding, L. B. A new potential surface and quasiclassical trajectory study of H +H2O→OH+H2. J. Chem. Phys. 2000, 113, 3150−3161. (103) Yang, M. H.; Zhang, D. H.; Collins, M. A.; Lee, S. Y. Ab initio potential-energy surfaces for the reactions OH+H2 ↔ H2O+H. J. Chem. Phys. 2001, 115, 174. (104) Liu, S.; Xiao, C.; Wang, T.; Chen, J.; Yang, T.; Xu, X.; Zhang, D. H.; Yang, X. The dynamics of the D2 + OH → HOD + D reaction: A combined theoretical and experimental study. Faraday Discuss. 2012, 157, 101. (105) Finlayson-Pitts, B. J.; Pitts, J. N., Jr. Chemistry of the upper and lower atmosphere; Academic Press: San Diego, 2000. (106) Miller, J. A.; Kee, R. J.; Westbrook, C. K. Chemical kinetics and combustion modeling. Annu. Rev. Phys. Chem. 1990, 41, 345−387. (107) Schatz, G. C.; Fitzcharles, M. S.; Harding, L. B. State-to-state Chemistry with Fast Hydrogen Atoms. Faraday Discuss. Chem. Soc. 1987, 84, 359−369. (108) Kudla, K.; Schatz, G. C.; Wagner, A. F. A quasiclassical trajectory study of the OH+ CO reaction. J. Chem. Phys. 1991, 95, 1635. (109) Bradley, K. S.; Schatz, G. C. A quasiclassical trajectory study of OH+ CO: Angular and translational distributions, and OH angular momentum alignment. J. Chem. Phys. 1997, 106, 8464. (110) Yu, H. G.; Muckerman, J. T.; Sears, T. J. A theoretical study of the potential energy surface for the reaction OH+ CO → H+ CO2. Chem. Phys. Lett. 2001, 349, 547−554. (111) Lakin, M. J.; Troya, D.; Schatz, G. C.; Harding, L. B. A quasiclassical trajectory study of the reaction OH+ CO→ H+ CO2. J. Chem. Phys. 2003, 119, 5848. (112) Valero, R.; van Hemert, M. C.; Kroes, G. J. Classical trajectory study of the HOCO system using a new interpolated ab initio potential energy surface. Chem. Phys. Lett. 2004, 393, 236−244. (113) Li, J.; Wang, Y.; Jiang, B.; Ma, J.; Dawes, R.; Xie, D.; Bowman, J. M.; Guo, H. Communication: A chemically accurate global potential energy surface for the HO+ CO→ H+ CO2 reaction. J. Chem. Phys. 2012, 136, 041103. (114) Jordan, M. J. T.; Gilbert, R. G. Classical trajectory studies of the reaction CH4+H→CH3+H2. J. Chem. Phys. 1995, 102, 5669−5682. (115) Espinosa-García, J.; Corchado, J. C. Recalibration of Two Earlier Potential Energy Surfaces for the CH4 + H →CH3 + H2 Reaction. Application of Variational Transition-State Theory and 2302

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303

Perspective

Journal of Chemical Theory and Computation

studying the H + CH4 →H2 + CH3 reaction on a neural network PES. J. Chem. Phys. 2015, 142, 064309. (159) Zhou, Y.; Wang, C.; Zhang, D. H. Effects of reagent vibrational excitation on the dynamics of the H + CHD3→H2 + CD3 reaction: A seven-dimensional time-dependent wave packet study. J. Chem. Phys. 2011, 135, 024313. (160) Liu, S.; Chen, J.; Zhang, Z.; Zhang, D. H. Communication: A six-dimensional state-to-state quantum dynamics study of the H + CH4 →H2 + CH3 reaction (J = 0). J. Chem. Phys. 2013, 138, 011101. (161) Palma, J.; Clary, D. C. A quantum model Hamiltonian to treat reactions of the type X+YCZ3→Y+CZ3: Application to O3(P)+CH4→ H+CH3. J. Chem. Phys. 2000, 112, 1859−1867. (162) Zhang, Z.; Chen, J.; Liu, S.; Zhang, D. H. Accuracy of the centrifugal sudden approximation in the H + CHD3 →H2 + CD3 reaction. J. Chem. Phys. 2014, 140, 224304. (163) Zhang, Z.; Zhang, D. H. Effects of reagent rotational excitation on the H + CHD3 →H2 + CD3 reaction: A seven dimensional timedependent wave packet study. J. Chem. Phys. 2014, 141, 144309. (164) Wang, Y.; Li, J.; Chen, L.; Lu, Y.; Yang, M.; Guo, H. Mode specific dynamics of the H2 + CH3 →H + CH4 reaction studied using quasi-classical trajectory and eight-dimensional quantum dynamics methods. J. Chem. Phys. 2015, 143, 154307. (165) Zhang, Z.; Chen, J.; Yang, M.; Zhang, D. H. Time-Dependent Wave Packet Study of the H2 + CH3 →H + CH4 Reaction. J. Phys. Chem. A 2015, 119, 12480−12484. (166) Yan, S.; Wu, Y.-T.; Zhang, B.; Yue, X.-F.; Liu, K. Do Vibrational Excitations of CHD3 Preferentially Promote Reactivity Toward the Chlorine Atom? Science 2007, 316, 1723−1726. (167) Yan, S.; Wu, Y.-T.; Liu, K. Tracking the energy flow along the reaction path. Proc. Natl. Acad. Sci. U. S. A. 2008, 105, 12667−12672. (168) Czakó, G.; Bowman, J. M. Dynamics of the Reaction of Methane with Chlorine Atom on an Accurate Potential Energy Surface. Science 2011, 334, 343−346. (169) Wang, F.; Lin, J.-S.; Cheng, Y.; Liu, K. Vibrational Enhancement Factor of the Cl + CHD3(v1 = 1) Reaction: Rotational-Probe Effects. J. Phys. Chem. Lett. 2013, 4, 323−327.

(138) Zare, R. N. Laser Control of Chemical Reactions. Science 1998, 279, 1875−1879. (139) Bronikowski, M. J.; Simpson, W. R.; Girard, B.; Zare, R. N. Bond-specific chemistry: OD:OH product ratios for the reactions H +HOD(100) and H+HOD(001). J. Chem. Phys. 1991, 95, 8647. (140) Bronikowski, M. J.; Simpson, W. R.; Zare, R. N. Effect of reagent vibration on the hydrogen atom + water-d reaction: an example of bond-specific chemistry. J. Phys. Chem. 1993, 97, 2194− 2203. (141) Fu, B.; Zhang, D. H. Mode specificity in the H + H2O →H2 + OH reaction: A full-dimensional quantum dynamics study. J. Chem. Phys. 2013, 138, 184308. (142) Fu, B.; Zhang, D. H. A full-dimensional quantum dynamics study of the mode specificity in the H + HOD abstraction reaction. J. Chem. Phys. 2015, 142, 064314. (143) Liu, S.; Zhang, D. H. A local mode picture for H atom reaction with vibrationally excited H2O: a full dimensional state-to-state quantum dynamics investigation. Chem. Sci. 2016, 7, 261−265. (144) Zhao, Z.-Q.; Liu, S.; Zhang, D. H. Differential Cross Sections for the H+D2O→D+OD Reaction: a Full Dimensional State-to-State Quantum Dynamics Study. Chin. J. Chem. Phys. 2017, 30, 16−24. (145) Fu, B.; Zhou, Y.; Zhang, D. H. Shape resonance in the H + D2O → D + HOD reaction: a full-dimensional quantum dynamics study. Chem. Sci. 2012, 3, 270−274. (146) Fu, B.; Zhang, D. H. Full-Dimensional Quantum Dynamics Study of the H + H2O and H + HOD Exchange Reactions. J. Phys. Chem. A 2012, 116, 820−825. (147) Fu, B.; Zhang, D. H. Full-dimensional quantum dynamics study of exchange processes for the D + H2O and D + HOD reactions. J. Chem. Phys. 2012, 136, 194301. (148) Ma, J.; Li, J.; Guo, H. Quantum Dynamics of the HO + CO → H + CO2 Reaction on an Accurate Potential Energy Surface. J. Phys. Chem. Lett. 2012, 3, 2482−2486. (149) Liu, S.; Xu, X.; Zhang, D. H. A full-dimensional timedependent wave packet study of the OH+CO → H+CO2 reaction. Theor. Chem. Acc. 2012, 131, 1068. (150) Wang, C.; Liu, S.; Zhang, D. H. Effects of reagent vibrational excitation on the state-to-state quantum dynamics of the OH+CO→H +CO2 reaction in six dimensions (J = 0). Chem. Phys. Lett. 2012, 537, 16−20. (151) Liu, S.; Chen, J.; Fu, B.; Zhang, D. H. State-to-state quantum versus classical dynamics study of the OH+CO→H+CO2 reaction in full dimensions (J = 0): checking the validity of the quasi-classical trajectory method. Theor. Chem. Acc. 2014, 133, 1558. (152) Camden, J. P.; Bechtel, H. A.; Zare, R. N. Dynamics of the Simplest Reaction of a Carbon Atom in a Tetrahedral Environment. Angew. Chem., Int. Ed. 2003, 42, 5227−5230. (153) Camden, J. P.; Bechtel, H. A.; Brown, D. J. A.; Zare, R. N. Effects of CH stretch excitation on the H+CH4 reaction. J. Chem. Phys. 2005, 123, 134301. (154) Camden, J. P.; Bechtel, H. A.; Ankeny Brown, D. J.; Martin, M. R.; Zare, R. N.; Hu, W.; Lendvay, G.; Troya, D.; Schatz, G. C. A Reinterpretation of the Mechanism of the Simplest Reaction at an sp3Hybridized Carbon Atom:H + CD4 →CD3 + HD. J. Am. Chem. Soc. 2005, 127, 11898−11899. (155) Hu, W.; Lendvay, G.; Troya, D.; Schatz, G. C.; Camden, J. P.; Bechtel, H. A.; Brown, D. J. A.; Martin, M. R.; Zare, R. N. H + CD4 Abstraction Reaction Dynamics: Product Energy Partitioning. J. Phys. Chem. A 2006, 110, 3017−3027. (156) Schiffel, G.; Manthe, U. Communications: A rigorous transition state based approach to state-specific reaction dynamics: Full-dimensional calculations for H+CH4→H2+CH3. J. Chem. Phys. 2010, 132, 191101. (157) Welsch, R.; Manthe, U. Communication: Ro-vibrational control of chemical reactivity in H+CH4→H2+CH3: Full-dimensional quantum dynamics calculations and a sudden model. J. Chem. Phys. 2014, 141, 051102. (158) Welsch, R.; Manthe, U. Full-dimensional and reduceddimensional calculations of initial state-selected reaction probabilities 2303

DOI: 10.1021/acs.jctc.8b00006 J. Chem. Theory Comput. 2018, 14, 2289−2303