Ab Initio-Based Global Double Many-Body ... - ACS Publications

Article pubs.acs.org/JPCA

Ab Initio-Based Global Double Many-Body Expansion Potential Energy Surface for the First 2A″ Electronic State of NO2 V. C. Mota, P. J. S. B. Caridade, and A. J. C. Varandas* Departamento de Química, Universidade de Coimbra, 3004-535 Coimbra, Portugal S Supporting Information *

ABSTRACT: An ab initio-based global double many-body expansion potential energy surface is reported for the first electronic state of the NO2(2A″) manifold. Up to 1700 ab initio energies have been employed to map the full configuration space of the title molecule, including stationary points and asymptotic channels. The calculated grid energies have been scaled to account for the incompleteness of the basis set and truncation of the MRCI expansion and fitted analytically with a total root-mean-squared-deviation smaller than 1.0 kcal mol−1. The lowest point of the potential energy surface corresponds to the 2B1 linear minimum, which is separated from the Cs distorted minimum by a C2v barrier of ≈9.7 kcal mol−1. As usual, the proposed form includes a realistic representation of long-range interactions. Preliminary work indicates its reliability for reaction dynamics.

1. INTRODUCTION The open-shell NO2 species plays a central role in environmental chemistry. In the atmosphere it is involved in the diverse processes from the chemistry of ozone depletion1,2 to the mechanism that cools down the terrestrial thermosphere, which is believed to be essentially due to infrared emission from nitric oxide.3−5 Here, the major sources of vibrationally excited NO are the reactions6−10 2 3 N( 4S) + O2 (3Σ− g ) → NO( Π) + O( P)

(1)

2 3 N(2D) + O2 (3Σ− g ) → NO( Π) + O( P)

(2)

O + NO(v′≤18) → O + NO(v″ R c for all R2 values (10)

where Rc is the asymptotic limit of a 1 A″/2 A″ crossing seam,35 which at the current level of theory is Rc ≈ 2.8a0. In DMBE theory, the dissociation limits in eq 10 can be mimicked by the analytic expression 2

Vdc(R ) = −

∑ n = 6,8,10

2 V (R) = V (3)(R) + O2 (3Σ− g )(R1) + NO( Π)(R2)

+ NO(2 Π)(R3) + [N(2 D) − N(4S)]f (R)

(14)

where h(R1) and g(r1) are given by similar forms as in the previous treatment of NH2 case,58 namely

3. DMBE MODELING The dissociation scheme for the lowest state of the NO2(2A″) manifold is given by

2

(13)

(3) (2) where VEHF (R) and Vdc (R), account now for the extended Hartree−Fock and dynamic correlation two-body contributions (2) (R) is the asymptotic exchange to V(2)(R), whereas Vexc contribution. Finally, [N(2D)−N(4S)]f(R) is a pseudo onebody energy term controlled by the switching function58 f(R). As will be shown, such a scheme will lead to an accurate description within the single valued DMBE formalism of the N + O2 channel described by eq 10. 3.1. f(R) Switching Function. Following previous work,58 the switching function of the pseudo one-body energy term assumes the form

(9)

(2) The two-body scaling parameters FAB in eq 8 were chosen 55 according to DMBE-SEC method to reproduce the bond dissociation energy of the corresponding AB diatom. At the (2) 2 present level of theory it gives FO(2)2( 3Σ−g ) = 0.8502 and FNO ( Π) = (3) = 0.8337. For the three-body scaling factor, a value of FNO 2 0.8636 has been utilized so as to reproduce the experimental value of Te, the energy difference between the Cs distorted minimum of NO2(12A″) and the global minimum for the ground electronic state NO2(12A′), including the zero-point energy (ZPE) correction. Thus, Te = T0 + E0[12A′] − E0[12A″], where E0[12A′] and E0[12A″] are the ZPEs for NO2(12A′) and NO2(12A″), respectively. The value22 T0 = 16234 cm−1 has been used, with E0[12A′] obtained from the fundamentals56 ωb = 749.649 cm−1 (bending), ωs = 1319.794 cm−1 (symmetric stretching)and ωa = 1616.852 cm−1 (asymmetric stretching). In turn, E0[12A″] has been obtained from ref 19, but with the reassignment suggested in ref 31: ωb = 747 cm−1, ωs = 1690 cm−1, and ωa = 256 cm−1. Considering12 De[12A′] = −0.3618Eh, a value of De[12A″] = −0.2856Eh has been obtained. A slightly different scaling factor of FN(2D)+O2(3Σ−g )(3) = 0.8372 has been used so as to mimic the experimental excitation energy57 EN(2D)−N(4S) = 0.0875931Eh. Note that the experimental value of EN(2D)−N(4S) includes spin−orbit coupling, although this is minor (≈8.7 cm−1) and hence ignored. Tables 1 and 2 display the harmonic frequencies, geometries and relative energies of all stationary points of the PES.

⎧O (3Σ−) + N(2 D) ⎪ 2 g ⎪ NO2 (2A″) ⇒ ⎨O2 (3Δu) + N( 4S) ⎪ ⎪ NO(2 Π) + O(3P) ⎩

(12)

C χn(R ) nn R

asym Vexc(R ) = Vexc (R ) χ5(R )

(11) 3025

(17) (18)

dx.doi.org/10.1021/jp300031q | J. Phys. Chem. A 2012, 116, 3023−3034

The Journal of Physical Chemistry A

Article

Table 1. Stationary Points of the Title Systema energyb R2

R3

2.8040 2.850 2.786 2.829 2.7422 2.7534 2.7426

α1

ωs, ωb, ωa

2.2127 2.207 2.210 2.213 2.2225 2.2197 2.2227

109.75 109.73 110.1 109.7 109.77 109.67 109.65 ≈108

1641, 1735, 1014, 1483, 1733, 1668, 1663, 1491, 1690,

2.4103 2.412 2.409 2.397 2.41 2.4104 2.4054 2.3968

2.4103 2.412 2.409 2.397 2.41 2.4104 2.4054 2.3968

110.00 110.06 113.80 111.7 110.0 109.63 108.14 107.77

1312, 1300, 1082,

2.538 2.576 2.5363

2.538 2.576 2.5363

85.38 68.80 79.74

2.579 2.598 2.5836 2.5822 2.5821

2.579 2.598 2.5836 2.5822 2.5821

68.44 62.80 67.19 66.99 67.01

4.357 4.316 4.2672 4.3184 4.3183

4.357 4.316 4.2672 4.3184 4.3183

30.49 30.80 31.11 30.53 30.53

R2

R1

α3

4.078 4.002 3.9191 4.0219 4.0367 3.9822

2.289 2.289 2.2963 2.2998 2.2903 2.2808

112.97 109.30 111.71 111.12 109.44 108.84

≈2.52

level of theory

Cs Distorted Minimum (DM) 760, 357 768, 310 −110.29 740, 282 −22.76d 750, 275 820, 454 −113.95 780, 369 −113.95, −26.67d 786, 293 −113.95, −26.67d 747, 256 747, 256 C2v Transition State (TS1) 771, 556i 777, 432i −109.38 558, 826 −110.09 −21.52d

1320, 777, 558i −112.83 1454, 723, 45i −113.45, −26.17d 1686, 698, 663i −113.14, −25.86d C2v “Intersection” in the Insertion Path (TS2) −68.66 −81.08 −77.81 2 B1-Ring Form, C2v Minimum (BM) 1239, 722, 613 −82.79 1269, 506, 800 −81.40 1258, 777, 694 −83.81 1262, 790, 693 −86.09 1478, 944, 927 −85.59 Transition State for Insertion (TS3) 217i, 1629, 180 1.25 229i, 1411, 158i 1.43 405i, 1419, 174 2.66 477i, 1369, 210 2.95 944i, 6609, 214i 2.07 energyb

0.00 0.00 0.00 0.00 0.00 0.00 0.00

MRCI/AVQZ31 CASPT2/VTZ14 MRCI/VQZ40 MRCI(Q)/AVQZe DMBE-SECc DMBE PES DMBE+G PES exptf exptg

1.41 0.91

MRCI/AVQZ31 CASPT2/VTZ14 MBE PES14 MRCI/VQZ40 MRCI(Q)/AVQZh DMBE-SECc DMBE PES DMBE+G PES

1.25 1.48 1.12 0.50 0.81 42.59 36.14

CASPT2/VTZ14 MBE PES14 DMBE PES

28.87 27.86 28.36

CASPT2/VTZ14 MBE PES14 MRCI/AVQZc DMBE-SECc DMBE PES

115.34 116.91 116.03

CASPT2/VTZ14 MBE PES14 MRCI/AVQZc DMBE-SECc DMBE PES

27.50

ωs, ωb, ωa

level of theory

Transition State for Abstraction (TS4) 1534, 235, 164i 0.25 1543, 131, 174i 0.21 1437, 173, 248i 1.65 1448, 150, 239i 0.88 1446, 161, 253i 0.69 1554, 199, 200i 0.63

111.00 114.33 113.78 114.65 114.58

CASPT2/VTZ14 MBE PES14 MRCI/AVQZc MRCI(Q)/AVQZc DMBE-SECc DMBE PES

Distances in bohr, angles in degrees, frequencies in cm−1, and energies in kcal mol−1. bEnergies relative to the N + O2 limit and Cs minimum in the first and second subcolumns, respectively. cThis work, from a fit to a Taylor-series type expansion around the stationary point; see the text. dRelative to O + NO limit. eReference 48: geometries are ab initio optimized, whereas the vibrational energies stem from calculations in the LSN PES assuming a C 2v assignment. fGeometries are obtained22 by assuming a C2v assignment, whereas fundamentals utilize19 the Cs assignment. g Fundamentals from ref 19 after the reassignment suggested in ref 31. hValues from calculations in the LSN PES assuming C2v assignment.48 a

with the long-range damping functions χn(R) assuming the form ⎡ ⎤n ⎛ R R2 ⎞⎥ ⎢ ⎜ ⎟ χn(R ) = 1 − exp⎜ −A n − Bn 2 ⎟ ⎢⎣ ρ ρ ⎠⎥⎦ ⎝

dimensionless parameters for all isotropic interactions:49,50 α0 = 16.36606, α1 = 0.70172, β0 = 17.19338, and β1 = 0.09574. In turn, the scaling parameter ρ is defined as ρ = 5.5 + 1.25R 0

(19)

(20)

where R0 = 2(⟨rA ⟩ + ⟨rB ⟩ ) is the Le Roy parameter for onset the undamped R−n expansion, and ⟨rX2⟩ is the expectation value of the squared radius for the outermost 2 1/2

and An and Bn being auxiliary forms defined by An = α0n−α1 and Bn = β0 exp(−β1n), where α0, β0, α1, and β1 are universal 3026

2 1/2

61



Article

Table 2. Novel Stationary Points of the Title Systema energyb α1

R2

R3

2.262 2.3037 2.3050

2.262 2.3037 2.3050

180 180 180

2.39 2.3329 2.1475

2.39 2.3329 2.4761

131 138.08 135.13

5.3612 5.1160 5.1891

5.3612 6.0021 3.9822

24.55 21.85 24.58

R2

R1

α3

2.2342 2.5078

2.6771 2.6899

161.23 94.25

4.7509 5.8974

6.0441 6.4705

18.92 19.56

ωs, ωb, ωa

level of theory

Linear Global Minimum (LM) 1171, 957, 1827 −31.71c 1422, 849, 2351 −38.60c 1397, 857, 2250 −38.41c Pseudo Transition State (TS5) −3.46,c 28.25d −16.95,d 21.65d −6.16,d 32.25d vdW Structures in the N + O2 Channel 144, 1608, 48i −0.75 139, 1609, 63 −0.88 200i, 1554, 199 0.63 energya)

−8.95 −11.92 11.73

MRCI/VQZ40 DMBE PES DMBE+G PES

19.30 9.73 20.52

MRCI/VQZ40 DMBE PES DMBE+G PES

113.20 113.08 114.58

ωs, ωb, ωa

DMBE PES DMBE PES DMBE PES level of theory

Cs Transition States 703, 645i, 1551 −42.86 825, 769i, 1021 −46.56 vdW Structures in the O + NO Channel 127, 1908, 160 −1.57c 45, 1904, 76i −0.95c

71.10 67.39

DMBE PES DMBE PES

26.67 −25.73

DMBE PES DMBE PES

a Distances in bohr, angles in degrees, frequencies in cm−1, and energies in kcal mol−1. bEnergies relative to the N + O2 limit and Cs minimum in the first and second subcolumns, respectively. cRelative to O + NO limit. dRelative to the linear global minimum.

asym electrons of atom X (=A, B). In turn, the Vexc (R) term in eq 18 is written as

asym Vexc (R ) = (αexcR2 − βexcR )R1/2 exp( −σexcR )

3.3. Three-Body Dynamical Correlation Energy. For the three-body dynamical correlation, we have employed the semiempirical form63 (3) V dc (R) =

(21)

where αexc, βexc, and σexc are numerical values for the diatomics in the corresponding electronic state. As usual,62 the damping function χ5(R) in eq 19 has been utilized to mimic charge overlap effects in the repulsive asymptotic exchange energy contribution at inner configurations. The available parameters in eq 21 are then commonly chosen from a least-squares fits to ab initio data. For simplicity, all values in eqs 17−21 were taken from ref 11. The EHF-type energy term assumes the general form VEHF(R ) = −

D (1 + R

i

n

(22)

Cn(i)(R i ,γi) =

∑

CnL(R i) PL(cos γi)

L ∈ {0,2,4}

where γ = γ0[1 + γ1 tanh(γ2r )]

(24)

n

where gi = 1/2{1 − tanh[ξ(ηRi−Rj−Rk)]} is a convenient damping function, with similar expressions for gj and gk. Following the original proposition,64 the scaling parameters have been fixed at η = 6 and ξ = 1.0a0−1. Similarly,64 the damping functions χn(ri) were taken the same as in eq 19, but now using the center-of-mass separation for the relevant atom−diatom channel instead of R. To determine the respective parameter ρ, we have taken R0 in eq 20 as the average of the Le Roy parameters for NS and OP, which yields ρ = 12.83 a0. For n = 6, 8 and 10, Cn(i)(Ri,γi) are atom−diatom dispersion coefficients given by

∑ air i) exp(−γr) i=1

∑ ∑ gi(R) Cn(i)(R i ,γi) χn(ri)ri−n

(25)

where PL(cos γi) denotes the Lth term of the Legendre polynomial expansion. For every L component, the respective internuclear dependence CnL(Ri) has been estimated65 by using the generalized Slater−Kirkwood approximation66 jointly with the dipole polarizabilities, which have been calculated ab initio (see section 2). The atom−diatom dispersion coefficients so calculated were then fitted to the form

(23)

r = R − Re is the displacement from the equilibrium diatomic geometry, and D and ai (i = 1, ..., n) are adjustable parameters to be obtained as described elsewhere.50,59 For consistency, in this work we have determined the parameters in eqs 22 and 23 from calibration of these diatomic potentials exclusively to the DMBE-SEC energies (see section 2), instead of using the potentials from ref 11. A rmsd of ≲5 cm−1 has been achieved for both VNO(R) and VO2(R), with the final curves exhibited in Figure 1 and in the final set of linear parameters given in Table 1 of the Supporting Information.

CnA − BC(R ) = CnAB + CnAC 3

+ DM [1 +

3

∑ ai̅ r i] exp(− a1̅ r − ∑ bi̅ r i) i=1

i=2 (26)

3027



Article

where r = R − RM is the displacement relative to the position of the maximum and b̅1 ≡ a1̅ . Note that for the N−O2 interaction, only the polarizabilities of O2(3Σg−) have been considered, thus reducing the atom−diatom dispersion interaction on this channel exclusively to those of N(2D) + O2(3Σg−). Clearly, this is an approximation for R1 > Rc (see eq 10), which should, however, be of no great concern. The resulting parameters are gathered in Table 2 of the Supporting Information, and the internuclear dependencies of the dispersion coefficients are shown in Figure 2 of the Supporting Information.

where i, gi(R), Ri, ri, and γi have the same meaning as before. In turn, θa ,i is the angle that specifies the orientation of the atomic quadrupole, and ϕab ,i is the corresponding dihedral angle. For an A−BC interaction, C4(Ri,ri) and C5(Ri,ri) assume the form 3 Q DBC(R i) χ4(ri) and 2 A 3 C5(R i ,ri) = Q AQ BC(R i) χ5(ri) 4

C4(R i ,ri) =

(28)

where QA represents the quadrupole of atom A, DBC(Ri) and QBC(Ri) are the dipole and quadrupole moments of the diatom BC as a function of the interatomic distance Ri, and χ4(ri) and χ5(ri) are the appropriate dumping functions defined as in eq 19. Following ref 67, the terms (DQ and (QQ in eq 27 read67 (DQ (θa,γ,ϕab) = cos γ(3 cos2 θa − 1) + 2 sin θa sin γ cos θa cos ϕab

(29)

and (QQ (θa,γ,ϕab) = 1 − 5 cos2 θa − 5 cos2 γ + 17 cos2 θa cos2 γ + 2 sin 2 θa sin 2 γ cos2 ϕab + 16 sin θa sin γ cos θa cos γ cos ϕab Figure 2. T-shaped profile. Key: ●, points used in DMBE fit; ○, points used in DMBE+G fit; , DMBE PES; ---, MBE PES;13 ···, atom−diatom dissociation limit; −·−, DMBE+G PES; −··−, DM energy.

(30)

As usual in the DMBE theory, the classical optimized quadrupole (COQ) model69−73 is utilized to eliminate the angle θa ,i in the electrostatic potential of eq 27. In physical terms, it is assumed that at any molecular arrangement of the atom the orientation of the quadrupole will instantaneously adapt the orientation such as to maximize the involved electrostatic interaction. One then obtains for the present case

As pointed out elsewhere,63 eq 24 causes an overestimation of the dynamical correlation energy at the atom−diatom dissociation channels. To correct such a behavior, we have multiplied the two-body dynamical correlation energy for the ith pair by ∏l≠i[1 − gl(R)]; correspondingly for channels j and k. This ensures63 that the only two-body contribution at the ith channel is that of jk. 3.4. Three-Body Electrostatic Interaction Energy. The dominant forces involved in long-range atom−diatom interactions are known to be of electrostatic nature. Within the multipole expansion treatment,67 they can be approximated by the interactions involving the atomic quadrupoles and diatomic dipoles and quadrupoles. For the O−NO channel, the O(3P) atom is known to exhibit a nonzero quadrupole momentum, thus allowing for the presence of all such interactions. For the N−O2 channel of the 12A″ state, the nitrogen electronic state will be N(2D) for R1 ≤ Rc as shown in eq 10, implying in a nonzero quadrupole moment and hence corresponding interactions to those for the O−NO channel. Of course, this differs from what is observed for the same channel of the NO2 ground state PES, because the electronic state of atomic nitrogen is then 4S. In DMBE theory, the total electrostatic interaction in eq 12 is written as68

∑ gi(R)⎢ i=1

+

/{[256C52 sin γ2 cos γ2 +2C4C5r cos γ(36 − 25 sin γ2) + 361C52 sin γ4 + 9C4 2r 2 + 144C52]1/2 − (5C4 2r 2 + 456C52) sin γ2 ∓ 3C4r cos γ ± C5(19 sin γ2 − 12)}}

⎢⎣

ri 4

C5(R i ,ri) (QQ (θa ,i,γi,ϕab ,i) ⎤ ⎥ ⎥⎦ ri 5

(31)

with C4 and C5 accounting respectively for C4(Ri,ri) and C5(Ri,ri) in eqs 28. Note that eq 31 is valid when both quadrupole−dipole and quadrupole−quadrupole interactions are present. Because O2(3Σg−) has a vanishing dipole moment, eq 31 will be valid only for the O−NO channel. For the N−O2 channel, one must use instead the expression from refs 70−72. Note also that for R1 ≥ Rc in eq 10 the nitrogen atom is in the N(4S) state, and hence the inclusion of the electrostatic terms for such regions of configuration space is a mere approximation. For the atomic quadrupoles in eqs 28 we have used the values QN(2D) = −0.0998ea02 and QO(3P) = −0.9449ea02, which were obtained via ab initio calculations. Because QN(2D) is significantly smaller than QO(3P), one expects the electrostatic interaction at the N−O2 channel to be of somewhat smaller importance. However, for completeness, they have all been included in our modeling. To model QNO(R), we have performed ab initio calculations relative to the center-of-mass

⎡ C4(R i ,ri) (DQ (θa ,i,γ ,ϕ ) i ab ,i

3

(3) (R) = V ele

θa = ∓arctan{2 sin γ(8C5 cos γ + C4r )

(27) 3028



Article

where S2a2 = Q 22 + Q 32, S2b2 = Q 22 + Q 32, and S33 = Q 33 + 3Q 22Q 3, and Q i (i = 1−3) are symmetry-adapted coordinates defined as

of the diatomic molecule for a given set of interatomic distances. The resulting values have then been fitted to the form73 Q (R ) = A(1 +

3

3

i=1

i=1

∑ aĩ r i) exp(− ∑ bĩ r i) +

M6 R6

χ8(R ) + Q ∞

⎛ 0⎞ ⎛Q1 ⎞ ⎛ 1/3 1/3 1/3 ⎞⎜ R1 − R1 ⎟ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎜Q 2 ⎟ = ⎜ 0 1/2 − 1/2 ⎟⎜ R2 − R20 ⎟ ⎟ ⎜⎜ ⎟⎟ ⎜ ⎟ ⎜ ⎟⎜ ⎝Q 3 ⎠ ⎝ 2/6 − 1/6 − 1/6 ⎠⎝ R3 − R30 ⎠

(32)

where r = R − Rref, with Rref denoting the distance at which the maximum of QNO arises, which is determined through the fitting. In turn, Q∞ corresponds to the value of the diatomic quadrupole at the separate-atoms limit. Finally, the function DNO(R) has been modeled in a similar fashion by fitting the ab initio data to the form D(R ) = A(1 +

3

3

i=1

i=1

(3) To model VEHF (R) in eq 34, we have started as usual by choosing convenient geometries where to center a given set of polynomials. Those were then calibrated via clj and γij, with the quality of the resulting least-squares fit determining the addition of further polynomials. The final set of values for all parameters in eqs 34 and 35 are displayed in Table 4 of the Supporting Information. Table 3 shows the corresponding

∑ aĩ r i) exp(− ∑ bĩ r i) (33)

with r and Rref defined as for Q(R) in eq 32. Figure 3 of the Supporting Information shows the final fits for QNO(R) and DNO(R), together with the respective ab initio calculations, whereas Table 3 of the Supporting Information gathers the numerical values of the least-squares parameters. For the N + O2 channel, we have used for QO2(R) the functional form from ref 74, which has been obtained by fitting eq 32 to MRCI/AVQZ calculations by Lawson and Harrison.75 3.5. Three-Body Extended Hartree−Fock Energy. For a given geometry, by removing from eq 11 the sum of the two-body energy terms plus the three-body [N(2D)−N(4S)]f(R) term, modeled in sections 3.1 and 3.2, one gets the total three-body energy. By subtracting now the three-body dynamical correlation plus electrostatic interaction obtained in sections 3.3 and 3.4, one then gets the three-body extended Hartree−Fock energy contribution shown in eq 12. To model the latter, we have utilized the three-body distributed-polynomial scheme76 by writing (3) V EHF (R) =

n

3

j=1

i=1

Table 3. Stratified Root-Mean-Squared-Deviation, rmsda

∑ P(j)(Q1,Q 2 ,Q 3) × ∏ {1 − tanh[γij(R i − R ij ,0)]}

ΔEb

no. of pts

max dev

rmsd

10 20 30 40 50 60 70 80 90 100 150 250 500 1000

287 354 379 524 550 587 659 710 740 778 1508 1646 1670 1681

2.982 2.982 4.317 4.317 4.317 4.558 4.558 4.558 4.558 4.558 5.963 5.963 5.963 5.963

0.517 0.543 0.669 0.600 0.652 0.766 0.794 0.872 0.896 0.934 0.878 0.983 0.986 0.988

a Energies and rmsd in kcal mol−1. bEnergies with respect to the linear minimum.

(34)

stratified rmsd. Note that although all ab initio data has been used in the fit, some points belonging to regions of strong nonadiabatic coupling have been given a vanishing weight in the calculation of the rmsd, as indicated in Figure 2.

with the polynomials P(j)(Q1,Q2,Q3) in the many-body expansion46,47 including up to sixth-order terms: P(Q1,Q 2 ,Q 3) = c1 + c 2Q1 + c3Q 3 + c4Q12 + c5S2a2 + c6Q1Q 3

4. RESULTS AND DISCUSSION Figures 2−8 depict the final fit, and Tables 1 and 2 gather the attributes of the resulting PES, hereinafter denoted as DMBE PES whereas the one from ref 14 will be addressed as MBE and the seminumerical one from ref 40 as SN. The global minimum of the DMBE PES is predicted to be a 2 B1 linear structure (LM), whose properties are given in Table 2. It is separated from the Cs distorted minimum (DM) by a pseudo-first-order transition state (TS5), which appears solely as a consequence of the single-valued formalism here utilized to model the 2A2/2B1 crossing seam that separates LM and DM (Figures 3−5). Note that the MBE PES does not describe either LM or TS5, as also seen from panels d and b of Figures 3 and 4, respectively. Although based on ab initio calculations10,26 to characterize LM, it is not clear whether such structures are described by the PESs from refs 41 and 10. Thus, the DMBE PES can only by compared in this region with the SN PES. Indeed, a comparison of Figure 2 from ref 40 with Figure 5 of this work shows that they have similar topographies, although Table 2 reveals that their attributes can be quantitatively

+ c 7S2b2 + c8Q13 + c 9Q1S2a2 + c10S33 + c11Q12Q 3 + c12Q1S2b2 + c13Q 3S2a2 + c14Q14 + c15Q12S2a2 + c16S2a 4 + c17Q1S33 + c18Q13Q 3 + c19Q12S2b2 + c 20Q1Q 3S2a2 + c 21Q 3S33 + c 22S2a2S2b2 + c 23Q15 + c 24Q13S2a2 + c 25Q1S2a 4 + c 26Q12S33 + c 27S2a2S33 + c 28Q14Q 3 + c 29Q13S2b2 + c30Q12Q 3S2a2 + c31Q1Q 3S33 + c32Q1S2a2S2b2 + c33Q 3S2a 4 + c34S2b2S33 + c35Q16 + c36Q14S2a2 + c37Q12S2a 4 + c38Q13S33 + c39Q1S2a2S33 + c40S2a6 + c41S36 + c42Q15Q 3 + c43Q14S2b2 + c44Q13Q 3S2a2 + c45Q12Q 3S33 + c46Q12S2a2 S2b2 + c47Q1Q 3S2a 4 + c48Q1S2b2S33 + c49Q 3S2a2S33 + c50S2a 4S2b2

(36)

(35) 3029



Article

Figure 3. Relaxed contour plots. Contours, in hartree, start at (a) −0.104001, (b) −0.243091, (c) −0.103990, and (d) −0.243074. Lines: , ---, and ··· are shown at intervals of +0.000835, −0.00835, and −0.000183, respectively, except for panel b, where ··· corresponds to −0.0003. Limits of optimization: 2.0 ≤ R1/a0 ≤ 2.5, ΔR1/a0 = 0.001 in panels a and c; 2.1 ≤ R2/a0 ≤ 2.6, ΔR2/a0 = 0.001 in panels b and d. The ab initio points in panels a and b refer to R1 = 2.28a0 and R2 = 2.17a0, respectively.

Figure 5. Contour plots which will allow a direct comparison with Figure 2 of ref 40. Contours start at −8.30 eV separated by 0.25 eV, with solid and dashed lines for positive and negative contours, respectively. Geometries: (a) x = R2 and y = R3, forming an angle of γxy/deg = 109.67; (b) relaxed plot for x = R2, y = γxy, and 2.1 ≤ R3/a0 ≤ 2.5 with ΔR3/a0 = 0.001; (c) same as panel a but with γxy/deg = 180.

Figure 4. Contour plots for T-shaped geometries (in hartree). Contours start at (a) −0.104001 and (b) −0.103990 [with intervals of 0.018375 (≈0.5eV)]. Solid and dashed lines indicate positive and negative energy spacings, respectively.

different. The DMBE PES is seen to predict LM to be 6.89 kcal mol−1 lower than the corresponding SN prediction, with the ωa and ωs values also showing significant differences (up to 251 and 524 cm−1, respectively). Such a difference on the energetics is likely due to the use of scaled energy points to built the DMBE PES. In turn, the differences in the normal modes frequencies show that the DMBE PES exhibits narrower wells for both the symmetric and asymmetric stretching modes

than the SN PES. Moreover, the height of TS5 as predicted by the DMBE PES with respect to LM and DM is, respectively, 6.6 and 9.57 kcal mol−1 lower than the corresponding predictions from the SN PES. Note, however, that TS5 is an artifact of the single-sheeted approaches devised for both the DMBE and SN PESs, which makes any rigorous comparison unworthy. 3030



Article

the Cs minimum. Regarding the DMBE PES description of the DM, we recall that the ab initio data employed on its calibration have been scaled so as to reproduce the experimental Te value (see section 2) relative to the ground state. Because both 12A′ and 12A″ states dissociate into atoms in the same electronic state, the DMBE PES also reproduces the experimental well depth. Table 1 includes the properties of the local semi numerical (LSN PES) form of ref 48, a suggested improved description of the DM region as given by the SN PES. Note that, with Te = 16502 cm−1, LSN PES predicts DM to be ≈0.77 kcal mol−1 above the experimental Te value, and hence also above the DMBE PES by the same value. The DM minima show as two symmetrical isomers separated by a C2v transition state (TS1). As indicated in Table 1, the barrier height predicted for TS1 by the DMBE PES (relative to DM) is lower than the value predicted by the scaled data by only 0.62 kcal mol−1. Similarly, the asymmetric stretching frequency is lower than the ab initio prediction by 513i cm−1. On the other hand, the predicted geometries show good agreement with the data from the scaled ab initio energies. Note that, as sketched in Figure 4 of ref 10, TS1 lies near to the 2 A2/2B1 crossing seam, which makes hard any description of this region of configuration space via a single-sheeted formalism. A similar situation occurs for the so-called 2B1-ring form C2v minimum10,14,28 (BM), also lying close to the 2A2/2B1 crossing seam. In single-sheeted models, such a crossing seam gives rise to a so-called14 peaked insertion (TS2), which arises along the insertion path. Its properties are shown in Table 1. However, the BM description has been less affected than TS1′s description, with geometries, normal modes frequencies, and relative energies in good agreement with the direct predictions from the scaled ab initio data. Also shown in Table 1 are the transition states for insertion (TS3) and abstraction (TS4) mechanisms. TS3 can be visualized in panel a of Figures 3 and 4. Although the scaled grid data for TS4 is satisfactorily well described by the DMBE PES, qualitative and quantitative differences are observed for TS3. The DMBE PES describes TS3 as a second-order transition state, whereas for both ab initio and scaled data it is a first-order transition state. Also, the ωs and ωb values differ substantially from the scaled data predictions. Note that similar difficulties have been reported by González et al.14 This may be due to a combination of two factors: first, in this energy range (≈116 kcal mol−1 above the Cs minimum) the rmsd is ≈0.84 kcal mol−1 (Table 3), which is ≈41% of the height of TS3 with respect to the N + O2 limit; second, as shown in ref 35, at least two asymptotic evolving crossing seams approach the N + O2 atom−diatom limit. Figure 7 compares the first term of the Legendre expansion of the DMBE PES for the O + NO channel with other results. Unlike the MBE PES, the DMBE curve agrees well with realistic representations of this channel, as would be expected given the careful treatment devoted to the long-range parts of the PES. We have also compared the DMBE PES with the MRCI/VQZ calculations employing spin−orbit (SO) coupling for O + NO interactions in Figures 9 and 10 of ref 40. Although lower in energy by typically ≈200 cm−1, the DMBE PES is practically parallel to the SO calculations. Regarding the N + O2 channel, we highlight the more realistic representation provided by the DMBE PES in comparison with the MBE PES, as shown in Figure 4. This is due primarily to the use of the switching function formalism from ref 58 for this channel. The atom− diatom interactions can be further visualized from panels a and b in Figure 3, which show the van der Walls (vdW) structures

To improve further the description of TS5, we have calculated an additional grid of ≈200 DMBE-SEC points covering this region of configuration space, and fitted them by adding to the DMBE PES the following localized term: G(R) =

1 2

n

3

j=1

i=1

∑ P(j)(Q1,Q 2 ,Q 3) ∏ exp[− γij(R i − R ij ,0)2 ]

× {1 + tanh[d0(r12 − d1)]}

(37)

with d0 = 5a0−2 and d1 = 0.4a02, where all coordinates are defined as in sections 2 and 3.5, and P(j)(Q1,Q2,Q3) is given by eq 35. Note that the Gaussian forms have been chosen such as to essentially restrict G(R) to the region of TS5. In turn, the extra tanh form has been used to preserve as much as possible the properties of LM. Table 5 of the Supporting Information gathers the G(R) final parameters, and Figure 2 illustrates the corresponding fit. Note that the open circles (which have not been used to fit the original DMBE PES) are the ones utilized to calibrate G(R). In the following we refer to the DMBE PES plus the extra G(R) term as DMBE+G PES and stress that both the DMBE PES and DMBE +G PES are independent forms that may be employed separately according to convenience. Tables 1 and 2 illustrate the attributes of the DMBE+G PES whenever they differ from the ones of the DMBE PES. As expected, the major changes observed occur for TS5. Indeed, its geometry shows a slight bent arrangement, with unequal RNO′ and RNO″ distances, whereas its energetics is much closer to the SN PES predictions. The DM attributes did not suffer any significant change, although LM and TS1 have been slightly shifted in energy (by 0.19 and 0.31 kcal mol−1, respectively). Although the shift in LM may not be relevant, for TS1 it represents an increase of 62% of the barrier high toward the DMBE-SEC prediction, thus being of potential implications for the isomerization dynamics (see text below). Figure 6 illustrates the DMBE+G PES.

Figure 6. Contour plots of DMBE+G PES for T-shaped geometries. Key for contours as in panel a of Figure 4.

The DMBE and SN PESs are also the only reported singlesheeted global forms that predict the DM. In fact, although the ab initio calculations utilized by González et al.14 describe the Cs minimum, its MBE PES shows the minimum to be of C2v symmetry. In turn, Duff and Sharma41 have used calculations26 that identify the geometry of the PES minimum as having C2v symmetry, whereas Braunstein and Duff10 make no reference to 3031



Article

Figure 7. Spherically averaged component of the O(3P)−NO(2Π) interaction potential.

superimposed on a relaxed triangular plot79 of the PES utilizing hyperspherical coordinates. As shown, the reaction is of the oxygenatom exchange type, illustrating the quenching of NO prepared in the first vibrational state. As also seen, the trajectory samples the deep Cs minima as well as the collinear one. A full report on the O + NO(v) quenching process will be published elsewhere.

exhibited in Table 2. The radial dependence is shown with further detail as inserts in panels a and d of Figure 4 of the Supporting Information. Also visible in panels a and b of Figure 3 are linear barriers for the N−OO and O−NO attacks. The N−OO barrier lies ≈7.7 kcal mol−1 above the N−O2 dissociation energy, and its coordinates are R1 ≈ 2.28a0, R2 ≈ 5.94a0, and R3 ≈ 3.66a0. In turn, the O−NO barrier lies ≈30.7 kcal mol−1 above the ONO dissociation energy, with a geometry defined by R1 ≈ 5.62a0, R2 ≈ 3.38a0, and R3 ≈ 2.24a0. Finally, the DMBE PES predicts two additional transition states; see Table 2. Although ab initio data covering this region have been included in the fit, the presence of strong couplings with further excited states of 2A″ symmetry suggests that the realism of such stationary states and their relevance for reaction dynamics requires further analysis. To test of the reliability of the present DMBE+G PES in dynamics studies, we have run many thousands of quasiclassical trajectories78 (QCT). Figure 8 shows one such a trajectory

5. CONCLUSIONS Extensive high-level ab initio calculations have been performed to map the entire configuration space of NO2(12A″), including all known stationary points as well as the atom−diatom channels. To account for the incompleteness of the basis set and truncation of the MRCI expansion, the data has been subsequently scaled using the DMBE-SEC method. For the title molecule, the scaling factor has been chosen so as to reproduce the experimental energy difference between the NO2(12A″) Cs distorted minimum and the NO2(12A′) global minimum. In turn, the diatomic fragments have been scaled such as to mimic the corresponding experimental dissociation energies. Using the scaled data, a single-sheeted DMBE PES for NO2(12A″), which reproduces all the stationary states and atom−diatom channels of the title system has been modeled. The final form provides also a realistic description of long-range forces, an asset of DMBE theory. In building the single-sheeted DMBE PES, a major source of concern has been a description as accurate as possible of regions nearby the various nonadiabatic structures. These are determinant for the DMBE PES along the insertion path, particularly near the linear (global) minimum, Cs, and B2 ringform minima. Also, the existence of asymptotic crossing seams for the N + O2 channel has demanded a careful treatment of the latter. We believe that such difficulties have been properly circumvented. Further improvements on the description of the insertion path and the N + O2 channel, within an analytic model, will therefore demand the construction of a multisheeted DMBE PES for the 2A″ manifold. This work as well as the results of trajectories run on the above DMBE+G PES are currently in progress and will be reported elsewhere.

Figure 8. Relaxed contour plot79 of DMBE+G PES. Contours start at 0.2650 Eh, spaced by 0.01 Eh. The curly solid line indicates a QCT trajectory; see text. Key for symbols: shaded circle, DM; shaded square, TS1; solid square, TS2; solid diamond, BM; solid pentagon; TS3; solid circle, TS4; shaded pentagon, LM; shaded diamond, TS5. 3032



■

Article

(27) Jackels, C. F.; Davidson, E. R. J. Chem. Phys. 1976, 64, 2908. (28) Jackels, C. F.; Davidson, E. R. J. Chem. Phys. 1976, 65, 2941. (29) Handy, N. C.; Goddard, J. D.; Schaefer, H. F. III. J. Chem. Phys. 1979, 71, 426. (30) Burton, N. A.; Yamaguchi, Y.; Alberts, I. L.; Schaefer, H. F. III. J. Chem. Phys. 1991, 95, 7466. (31) Bera, P. P.; Yamaguchi, Y.; Schaefer, H. F. III; Daniel Crawford, T. J. Phys. Chem. A 2008, 112, 2669. (32) Helgaker, T.; Jørgensen, P.; Olsen, J. Molecular ElectronicStructure Theory; John Wiley & Sons: Chichester, U.K., 2000. (33) Dunning, T. H. Jr. J. Chem. Phys. 1989, 90, 1007. (34) Kendall, R. A.; Dunning, T. H. Jr.; Harrison, R. J. J. Chem. Phys. 1992, 96, 6796. (35) Mota, V. C.; Caridade, P. J. S. B.; Varandas, A. J. C. Int. J. Quantum Chem. 2011, 111, 3776. (36) Köppel, H.; Domcke, W.; Cederbaum, L. S. Adv. Chem. Phys. 1984, 57, 59. (37) Sharma, R. D.; Duff, J. W. Geophys. Res. Lett. 1997, 24, 2407. (38) Caridade, P. J. S. B.; Mota, V. C.; Mohallem, J. R.; Varandas, A. J. C. J. Phys. Chem. A 2008, 112, 960. (39) Galvão, B. R. L.; Corzo-Espinoza, J. A.; Caridade, P. J. S. B.; Varandas, A. J. C. Int. J. Chem. Kinet. 2011, 43, 345. (40) Ivanov, M. V.; Zhu, H.; Schinke, R. J. Chem. Phys. 2007, 126, 054304. (41) Duff, J. W.; Sharma, R. D. J. Chem. Soc., Faraday Trans. 1997, 93, 2645. (42) Parr, C. A.; Truhlar, D. G. J. Phys. Chem. 1971, 75, 1844. (43) Katagiri, H.; Kato, S. J. Chem. Phys. 1993, 99, 8805. (44) Andersson, K.; Malmqvist, P. A.; Roos, B. O.; Sadlej, A. J.; Wolinski, K. J. Phys. Chem. 1990, 94, 5483. (45) Krishnan, R.; Binkley, J. S.; Seeger, R.; Pople, J. A. J. Chem. Phys. 1980, 72, 650. (46) Varandas, A. J. C.; Murrell, J. N. Faraday Discuss. Chem. Soc. 1977, 62, 92. (47) Murrell, J. N.; Carter, S.; Farantos, S. C.; Huxley, P.; Varandas, A. J. C. Molecular Potential Energy Functions; Wiley: Chichester, U.K., 1984. (48) Schinke, R. J. Chem. Phys. 2008, 129, 124303. (49) Varandas, A. J. C. J. Mol. Struct.: THEOCHEM. 1985, 120, 401. (50) Varandas, A. J. C. Adv. Chem. Phys. 1988, 74, 255. (51) Varandas, A. J. C. Lecture Notes in Chemistry; Laganá, A., Riganelli, A., Eds.; Springer: Berlin, 2000; Vol. 75, p 33. (52) Varandas, A. J. C. Advanced Series in Physical Chemistry; World Scientific Publishing: Singapore, 2004; Chapter 5, p 91. (53) Varandas, A. J. C. Chem. Phys. Lett. 1992, 194, 333. (54) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M.; et al. MOLPRO, a package of ab initio programs, version 2010.1; see http://www.molpro.net, Cardiff, U.K., 2010. (55) Varandas, A. J. C. J. Chem. Phys. 1989, 90, 4379. (56) Delon, A.; Jost, R. J. Chem. Phys. 1991, 95, 5686. (57) Ralchenko, Yu.; Kramida, A. E.; Reader, J.; NIST ASD Team. NIST Atomic Spectra Database (ver. 4.1.0), [Online]. Available: http:// physics.nist.gov/asd [2012, February 10]. National Institute of Standards and Technology: Gaithersburg, MD, 2011. (58) Varandas, A. J. C.; Poveda, L. A. Theor. Chem. Acc. 2006, 116, 404. (59) Varandas, A. J. C.; Silva, J. D. J. Chem. Soc., Faraday Trans. 2 1992, 88, 941. (60) Varandas, A. J. C. Mol. Phys. 1987, 60, 527. (61) Le Roy, R. J. Spec. Period. Rep. Chem. Soc. Mol. Spectrosc. 1973, 1, 113. (62) Varandas, A. J. C.; Voronin, A. Asian J. Spectrosc. 1997, 1, 135. (63) Varandas, A. J. C. J. Chem. Phys. 1996, 105, 3524. (64) Varandas, A. J. C.; Rodrigues, S. P. J J. Phys. Chem. A 2006, 110, 485. (65) Varandas, A. J. C.; Rodrigues, S. P. J Chem. Phys. Lett. 1995, 245, 66. (66) Matias, M. A.; Varandas, A. J. C. Mol. Phys. 1990, 70, 623. (67) Buckingham, A. D. Adv. Chem. Phys. 1967, 12, 107.

ASSOCIATED CONTENT

S Supporting Information *

Tables 1−5 of the Supporting Information gather the fitted parameters obtained in sections 2−4 of this work. In turn, Figure 1 of the SI contains a description of the utilized coordinate sets, and Figures 2 and 3 the fits of the dispersion coefficients and permanent electric moments obtained in sections 3.3 and 3.4. Additionally, Figure 4 shows selected profiles of the DMBE and MBE+LR PESs together with the corresponding scaled ab initio points. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work has the support of European Space Agency under ESTEC Contract No. 21790-/08/NL/HE, and Fundaçaõ para a Ciência e Tecnologia, Portugal, under contracts PTDC/QUIQUI/099744/2008 and PTDC/AAC-AMB/099737/2008.

■

REFERENCES

(1) Rowland, F. S. Philos. Trans. R. Soc. B 2006, 361, 769. (2) Crutzen, P. J.; Oppenheimer, M. Climatic Change 2008, 89, 143. (3) Kockarts, G. Geophys. Res. Lett. 1980, 7, 137. (4) Zachor, A. S.; Sharma, R. D.; Nadile, R. M.; Stair, A. T. J. Geophys. Res. 1985, 90, 9776. (5) Sharma, R. D.; Roble, R. G. Chem. Phys. Chem 2002, 3, 841. (6) Rusch, D. W.; Gerard, J. C.; Fesen, C. G. J. Geophys. Res. 1991, 96, 11331. (7) Barth, C. A. Planet. Space Sci. 1992, 40, 315. (8) Sharma, R. D.; Kharchenko, V. A.; Sun, Y.; Dalgarno, A. J. Geophys. Res. 1996, 101, 275. (9) Gerard, J. C.; Bisikalo, D. V.; Shematovich, V. I.; Duff, J. W. J. Geophys. Res. 1997, 102, 285. (10) Braunstein, M.; Duff, J. W. J. Chem. Phys. 2000, 113, 7406. (11) Varandas, A. J. C.; Voronin, A. I. Mol. Phys. 1995, 85, 497. (12) Varandas, A. J. C. J. Chem. Phys. 2003, 119, 2596. (13) González, M.; Miquel, I.; Sayós, R. Chem. Phys. Lett. 2001, 335, 339. (14) González, M.; Miquel, I.; Sayós, R. J. Chem. Phys. 2001, 115, 8838. (15) Mahapatra, S.; Köppel, H.; Cederbaum, L. S.; Stampfuβ, P.; Wenzel, W. Chem. Phys. 2000, 259, 211. (16) Kurkal, V.; Fleurat-Lessard, P.; Schinke, R. J. Chem. Phys. 2003, 119, 1489. (17) Schinke, R.; Grebenshchikov, S. Y.; Zhu, H. Chem. Phys. 2008, 346, 99. (18) Reignier, D.; Stoecklin, T.; Halvick, P.; Voronin, A.; Rayez, J. C. Phys. Chem. Chem. Phys. 2001, 3, 2726. (19) Shibuya, K.; Terauchi, C.; Sugawara, M.; Aoki, K.; Tsuji, K.; Tsuchiya, S. J. Mol. Struct.: THEOCHEM 1997, 413, 501. (20) Weaver, A.; Metz, R. B.; Bradforth, S. E.; Neumark, D. M. J. Chem. Phys. 1989, 90, 2070. (21) Shibuya, K.; Kusumoto, T.; Nagai, H.; Obi, K. J. Chem. Phys. 1991, 95, 720. (22) Aoki, K.; Hoshina, K.; Shibuya, K. J. Chem. Phys. 1996, 105, 2228. (23) Kaldor, U. Chem. Phys. Lett. 1991, 185, 131. (24) Huzinaga, S. J. Chem. Phys. 1965, 42, 1293. (25) Dunning, T. H. Jr. J. Chem. Phys. 1970, 53, 2823. (26) Gillispie, G. D.; Khan, A. U.; Wahl, A. C.; Hosteny, R. P.; Krauss, M. J. Chem. Phys. 1975, 63, 3425. 3033



Article

(68) Varandas, A. J. C.; Brandão, J. Mol. Phys. 1986, 57, 387. (69) Varandas, A. J. C.; Brandão, J.; Quintales, A. M. J. Phys. Chem. 1988, 92, 3732. (70) Varandas, A. J. C. J. Mol. Struct.: THEOCHEM 1988, 166, 59. (71) Varandas, A. J. C.; Pais, A. A. C. C. Mol. Phys. 1988, 65, 843. (72) Varandas, A. J. C.; Rodrigues, S. P. J J. Chem. Phys. 1997, 106, 9647. (73) Rodrigues, S. P. J; Varandas, A. J. C. Phys. Chem. Chem. Phys. 2000, 2, 435. (74) Rodrigues, S. P. J; Sabín, J. A.; Varandas, A. J. C. J. Phys. Chem. A 2002, 106, 556. (75) Lawson, D. B.; Harrison, J. F. J. Phys. Chem. A 1997, 101, 4781. (76) Martı ́nez-Núñes, E.; Varandas, A. J. C. J. Phys. Chem. A 2001, 105, 5923. (77) Harding, L. B.; Stark, H.; Troe, J.; Ushakov, V. G. Phys. Chem. Chem. Phys. 1999, 1, 63. (78) Hase, W. L.; Duchovic, R. J.; Hu, X.; Komornicki, A.; Lim, K. F.; Lu, D.; Peslherbe, G. H.; Swamy, K. N.; Linde, S. R. V.; Varandas, A. J. C. QCPE Bull. 1996, 16, 43. (79) Varandas, A. J. C. Chem. Phys. Lett. 1987, 138, 455.

3034


Ab Initio-Based Global Double Many-Body ... - ACS Publications

Recommend Documents