Perspective on Diabatic Models of Chemical Reactivity as Illustrated by the Gas-Phase S(N)2 Reaction of Acetate Ion with 1,2-Dichloroethane

Abstract

Diabatic models are widely employed for studying chemical reactivity in condensed phases and enzymes, but there has been little discussion of the pros and cons of various diabatic representations for this purpose. Here we discuss and contrast six different schemes for computing diabatic potentials for a charge rearrangement reaction. They include (i) the variational diabatic configurations (VDC) constructed by variationally optimizing individual valence bond structures and (ii) the consistent diabatic configurations (CDC) obtained by variationally optimizing the ground-state adiabatic energy, both in the nonorthogonal molecular orbital valence bond (MOVB) method, along with the orthogonalized (iii) VDC-MOVB and (iv) CDC-MOVB models. In addition, we consider (v) the fourfold way (based on diabatic molecular orbitals and configuration uniformity), and (vi) empirical valence bond (EVB) theory. To make the considerations concrete, we calculate diabatic electronic states and diabatic potential energies along the reaction path that connects the reactant and the product ion-molecule complexes of the gas-phase bimolecular nucleophilic substitution (S(N)2) reaction of 1,2-dichloethane (DCE) with acetate ion, which is a model reaction corresponding to the reaction catalyzed by haloalkane dehalogenase. We utilize ab initio block-localized molecular orbital theory to construct the MOVB diabatic states and ab initio multi-configuration quasidegenerate perturbation theory to construct the fourfold-way diabatic states; the latter are calculated at reaction path geometries obtained with the M06-2X density functional. The EVB diabatic states are computed with parameters taken from the literature. The MOVB and fourfold-way adiabatic and diabatic potential energy profiles along the reaction path are in qualitative but not quantitative agreement with each other. In order to validate that these wave-function-based diabatic states are qualitatively correct, we show that the reaction energy and barrier for the adiabatic ground state, obtained with these methods, agree reasonably well with the results of high-level calculations using the composite G3SX and G3SX(MP3) methods and the BMC-CCSD multi-coefficient correlation method. However, a comparison of the EVB gas-phase adiabatic ground-state reaction path with those obtained from MOVB and with the fourfold way reveals that the EVB reaction path geometries show a systematic shift towards the products region, and that the EVB lowest-energy path has a much lower barrier. The free energies of solvation and activation energy in water reported from dynamical calculations based on EVB also imply a low activation barrier in the gas phase. In addition, calculations of the free energy of solvation using the recently proposed SM8 continuum solvation model with CM4M partial atomic charges lead to an activation barrier in reasonable agreement with experiment only when the geometries and the gas-phase barrier are those obtained from electronic structure calculations, i.e., methods i-v. These comparisons show the danger of basing the diabatic states on molecular mechanics without the explicit calculation of electronic wave functions. Furthermore, comparison of schemes i-v with one another shows that significantly different quantitative results can be obtained by using different methods for extracting diabatic states from wave function calculations, and it is important for each user to justify the choice of diabatization method in the context of its intended use.

FIG. 1

Reactant and product VB structures and schematic correlation diagram for the acetate + DCE S_N2 reaction.

FIG 2

Geometric variables for the acetate + DCE reaction showing the formation of the O-C bond and the breaking of the C-Cl bond.

FIG. 3

Ground-state adiabatic potential energy profiles along the reaction coordinate computed with (a) MC-QDPT(4)//M06-2X (full line), M06-2X (dashed line), CDC-MOVB (dotted line), and B3LYP (dashed-dotted line), and (b) EVB. For each method, the zero of energy is taken as the ground-state adiabatic energy of the ion-molecule reactant complex. Note the change in the scale of the abscissa between parts (a) and (b).

FIG. 4

Definition of reaction energy and reaction barrier on the ground-state PES used in the present study.

FIG. 5

Contour plots of the ground-state adiabatic potential surface in the plane of the O-C and C-Cl distances. The levels of calculation used to construct the plots are (a) B3LYP, (b) M06-2X, and (c) EVB. “RS”, “TS”, and “PS” stand for the reactant ion-molecule complex, the transition state, and the product ion-molecule complex, respectively. Note the different ordinate between parts (a) and (b). The energy units used are kcal/mol.

FIG. 5

Contour plots of the ground-state adiabatic potential surface in the plane of the O-C and C-Cl distances. The levels of calculation used to construct the plots are (a) B3LYP, (b) M06-2X, and (c) EVB. “RS”, “TS”, and “PS” stand for the reactant ion-molecule complex, the transition state, and the product ion-molecule complex, respectively. Note the different ordinate between parts (a) and (b). The energy units used are kcal/mol.

FIG. 6

Comparison of adiabatic and diabatic potential profiles. (a) Adiabatic HF profile (black line) and adiabatic and diabatic MOVB potential energy profiles obtained by nonorthogonal VDC (red lines) and nonorthogonal CDC (blue lines) schemes. The panel in the center is a blowup of the central region of the larger panel. (b) Same as part (a) except MOVB results are in the orthogonal representation. In both (a) and (b) the zero of energy is chosen as the ground-state HF adiabatic energy of the ion-molecule reactant complex.

FIG. 7

Fourfold way (MC-QDTP(4)//M06-2X) adiabatic and diabatic energies along the reaction coordinate. The zero of energy is chosen as the ground-state adiabatic energy of the ion-molecule reactant complex.

FIG. 8

EVB adiabatic and diabatic energies along the reaction coordinate. The zero of energy is chosen as the ground-state adiabatic energy of the ion-molecule reactant complex.

FIG. 9

Diabatic coupling matrix elements as a function of the reaction coordinate. (a) Nonorthogonal MOVB resonance energy of eq 9a or the equivalent eq 10 (full line for the CDC method, and dashed line for the VDC method), (b) orthogonal MOVB $H_{12}^{\text{s}}$ (CDC is given as a 12 full line, and VDC is shown by a dashed line), and (c) fourfold way (MC-QDTP(4)//M06-2X) (curve with a peak) and EVB (straight line).

FIG. 9

Diabatic coupling matrix elements as a function of the reaction coordinate. (a) Nonorthogonal MOVB resonance energy of eq 9a or the equivalent eq 10 (full line for the CDC method, and dashed line for the VDC method), (b) orthogonal MOVB $H_{12}^{\text{s}}$ (CDC is given as a 12 full line, and VDC is shown by a dashed line), and (c) fourfold way (MC-QDTP(4)//M06-2X) (curve with a peak) and EVB (straight line).