J. Phys. Chem. 1994,98, 12506-12515
12506
Chemical Reaction Molecular Dynamics Simulation and the Energy-Transfer Mechanism in the Proton-Transfer Reaction of Formamidine in Aqueous Solution Masataka Nagaoka; Yoshishige Okuno,?,'and Tokio Yamabet Institute for Fundamental Chemistry, 34-4, Takano-Nishihiraki-cho,' Sakyo-ku, Kyoto 606, Japan Received: July 18, I994@
The energy-transfer mechanism in the proton-transfer reaction of formamidine in aqueous solution was investigated by using the chemical reaction molecular dynamics method. On the basis of our previous studies (J. Am. Chem. SOC.1991, 113, 769; Can. J . Chem. 1992, 70, 377), the reaction was described as the protontransfer reaction of a supermolecule consisting of formamidine and a reactive water molecule (i.e., the F W supermolecule) surrounded with medium water molecules. It was found that, in relaxation of the reaction energy, after fast damping of the potential energy of the FW supermolecule with a sudden increase in the vibrational kinetic energy of the supermolecule and a slight and gradual increase in the rotational energy, they decrease gradually, accompanied by dissipation of the solute energy into the medium. Some intemal vibrations of the FW supermolecule remained unrelaxed within the time duration of the present simulations, i.e., 0.5 ps. These qualitative features were in good agreement with the study of Gertner et al. (J. Am. Chem. Soc. 1991, 113, 74). Characteristically, the change of the interaction energy between the FW supermolecule and the medium water molecules was found to be small during the course of the reaction. This feature is contrasted with the case studied by Gertner et al. Since the present system is not involved in charge redistribution, the reorganization of medium water molecules is not essential for the reaction occurrence. In addition, the fact that there was no recrossing in the transition state region is also attributed to the same reason. Analysis of the work done on the FW supermolecule by medium water molecules has shown that the componential work vertical to the molecular plane contributes to the energy supply into the FW supermolecule although the energy relaxation process itself takes place. The medium water molecules that solvate the FW supermolecule by hydrogen bonding were also found to be the most influential in energy flow, whether they influenced the mechanism positively or negatively. In particular, it was remarkably found that most of the medium water molecules that provide the most and least work for the reaction progress are localized in the regions of those protons which locate far from the reaction site.
1. Introduction The synthesis of organic chemical substances has been usually performed in solution.' In the biosynthetic pathway, it is often important to consider the effect of the medium surrounding biochemical species.* Similarly, most of chemical reactions take place under the influence of medium or solvent. Although there are many chemists who want to understand these chemical reactions theoretically, it is not an exaggeration to say that the level of understanding of liquid-phase chemical reactions is, even now, still in its i n f a n ~ y ~in- contrast ~ to that of gas-phase chemical reactions. Remaining problems in the liquid phase are, more or less, associated with its many degrees of freedom. So far, much attention has been paid mainly to thermodynamic properties for understanding chemical reactions in s ~ l u t i o n , ~because -~ the properties of the liquids and the reactions can be characterized by only several thermodynamic variables. It was difficult for one to get a desire that the chemical reactions in solution are understood from the first principles of mechanics in the same way as the gas-phase reactions. However, recent rapid progress of computational facilities and remarkable developments of computer-simulation techniques have not only improved the circumstances but also changed the theoretical strategy for this subject. Recently,
* To whom correspondence
should be addressed. Also associated with the Department of Hydrocarbon Chemistry, Faculty of Engineering, Kyoto University, Sakyo-ku, Kyoto 606, Japan. Present address: Daicel Chemical Industries, Ltd., 1239, Shinzaike, Aboshi-ku, Himeji, Hyogo, 671-12, Japan. Abstract published in Advance ACS Absrracts, November 1, 1994.
*
@
0022-365419412098-12506$04.50/0
mechanical understanding has been obtained by some numerical simulations. In particular, the molecular dynamics (MD) method7**and the Monte Carlo methodg-'* have been theoretical tools used widely for the studies in this field. A liquid-phase chemical reaction can be regarded as such a process that reactants climb the potential energy barrier with the supply of energy from the surrounding medium. It is usually brought about by solvent fluctuation in equilibrium. In order for a reaction to occur, the solute potential energy must become large enough to surmount the reaction barrier involved. For this purpose, it is natural that various kinds of energy in solution should transform into the potential energy of the solute. It should be noted that in the gas phase the energy source is only the relative kinetic energy between colliding reactant species. Therefore, much attention in the gas-phase reaction has been paid often just to the potential energy surface, that is, the static feature in the reaction process, since the kinetic energy partitioning can be characterized easily by description of the potential energy surface.13-17 On the contrary, in solution, the energy to climb the barrier must be supplied from fluctuations of various kinds of energies: (i) the kinetic energies of the solute and solvent molecules; (ii) the intermolecular potential energy between the solute and solvent molecules; and (iii) that between solvent molecules. Therefore, it is natural and very important to include not only the potential energy but also the kinetic energy for understanding solution chemical reactions. Thus, complete understanding of energy flow through various kinds of fluctuations is a serious problem in solvent dynamics, and it is important to understand the role of solvent in liquid-phase 0 1994 American Chemical Society
Proton-Transfer Reaction of Formamidine in Aqueous Solution
J. Phys. Chem., Vol. 98, No. 48, 1994 12507
-7:
chemical reactions. Answering this problem is not so easy by experimental techniques as well. It is true that since the energy flow is really time-dependent, the theoretical treatment is also limited. Thus, the MD method, which solves directly a set of Y H classical equations of motion for molecules and is a timedependent theoretical tool, has no rival for this study. Ohminelg studied the energy dissipation mechanism of an optically excited ethylene molecule in argon or water liquid. He found that there is an important pathway for the fast energy dissipation in both argon and water solvent, which includes a Figure 1. Intrasupermolecular proton-transfer reaction. multistep collision process leading to energy flow from internal vibrations of the ethylene to solvent molecule motion through the driving force is induced and induces the proceeding reaction. ethylene rotation. Further, it was clarified that the energy decay Additionally, comparison with the previous supermolecule rate depends strongly on the nature of the solvent-solvent treatment^^^*^^-^^ will provide important and helpful information interaction and, especially in water solvent, ethylene motions in scrutinizing many kinds of theoretical approaches in liquidare coupled strongly with the large water kinetic energy phase chemical reactions. fluctuation which transfers efficiently the ethylene internal On the basis of our previous studies26,28where a water energy to the water librational kinetic energy. On the other molecule was found remarkably influential for potential energy hand, Benjamin et al.19 investigated the energy flow mechanism barrier reduction in the proton-transfer reaction of formamidine, for the sN2 reaction, c1 4-C12 c12 C1, in argon-like liquid the present system is made up of 1 supermolecule consisting as a weakly coupled system. They clarified the temporal order of a formamidine molecule and a reactive water molecule, and of excitations of various kinds of energy during the course of 61 medium water molecules in a cubic box, which is replicated the reaction. Similarly, Gertner et a1.20investigated the energythroughout space to form an infinite lattice in the CRMD transfer mechanism for a model sN2 reaction, C1- CH3C1simulation. Therefore, it is considered as a justifiable model ClCH3 C1-, in aqueous solution. In many respects, this for the proton-transfer reaction of a formamidine-water (FW) strongly coupled system is contrasted remarkably with the supermolecule (Figure 1) in aqueous solution. A method weakly coupled system. In particular, it is important that many starting trajectory calculations at the transition state (TS) of the water molecules were found active in the energy-transfer FW supermolecule was employed, following the techniques mechanism, whose situation is in contrast to that in the argon developed by several groups.32 Initial conditions of positions solvent. Li et aLZ1clarified that argon at moderate pressure of all the atoms were specified by the trajectory calculations, increases the initial instantaneous rate of dissociation of HOC1 assuming a canonical distribution under the condition that the to OH and C1 after the overtone-excitation of the OH vibration, FW supermolecule was fixed at the TS structure. Succeeding over that in the gas phase. The energy flow within the solute trajectory calculations were performed in the constant-energy was investigated also by analysis of the time-dependent velocity by using the velocities newly specified to obey power spectrum. Maxwell’s law of velocity distribution at 300 K. In these studies,18-21 however, the reacting species were We investigated not only the energy flow within the FW limited usually to the system consisting of a few atoms or the supermolecule but also the relaxation of the internal energy into simplified model system with a united atom under some the medium. The time-ordering of various kinds of energy flow constraints. However, for most realistic chemical reactions, in the relaxation process was clarified; after the vibrational solutes usually consist of many more atoms, and the intramokinetic energy within the FW supermolecule increased suddenly lecular energy redistribution can be expected to become more by the fast damping of its potential energy, the rotational energy complicated but important for the mechanism of energy flow. of the FW supermolecule increased slightly and then the total Furthermore, in such a system composed of a strongly coupled kinetic energy of the FW supermolecule gradually decreased solvent system and a solute that hardly involves any charge and was transformed into solvent energy, while some internal redistribution during the course of the reaction, there should vibrations of the FW supermolecule remained unrelaxed for a exist some simple principles which might be the subject we long time. The change of the interaction energy between the should clarify first, in contrast to the other studies.18-21 FW supermolecule and the medium water molecules was found Under the circumstances, we present here an investigation to be small during the course of the reaction because there was of the energy flow for a more realistic chemical reaction in no charge redistribution in the FW supermolecule. It was also solution, i.e., the proton-transfer reaction of formamidine in clarified that there was no recrossing in the TS region. This is aqueous solution, at the atomic level by using the chemical in good agreement with the previous result obtained by the reaction MD (CRMD) m e t h ~ d . ~Much ~ . ~ attention ~ has been supermolecule treatment.30q31Analysis of the work done on the paid to the reaction system especially in the context of F W supermolecule by the medium water molecules showed that, biochemical processes.24 Until now, we have investigated not while most of the reaction energy was dissipated into several only electronically but also dynamically the effect of the solvent directions within the molecular plane on which the reactive water on the reaction by supermolecule treatment.22%25-31 motion occurs, some amounts of energy were supplied, at the However, any investigation from the energy-transfer viewpoint same time, from the vertical direction to the FW supermolecular has not been done due to the complexity that many solvent plane. The reaction energy was found to flow mainly into the molecules should be involved in the mechanism in various medium water molecules solvating strongly with the FW complicated ways. Therefore, one important purpose of the supermolecule. It was characteristic that most of the influential present paper is to clarify the specific pictures of energy flow medium water molecules in each trajectory are coupled strongly in the reaction system and to establish the general pictures of to those protons which locate far from the reaction site. energy flow by comparison with other studies. The energy flow We clarify the energy flow process for the proton-transfer reaction of formamidine in aqueous solution, focusing on the must be classified into several types in accordance with those of chemical reactions. The other is to answer the question how difference of our result from other previous s t u d i e ~ . ~The ~.~~
-
+
+
+
12508 J. Phys. Chem., Vol. 98, No. 48, 1994
Nagaoka et al.
difference in energy flow mechanism is understood by comparing the qualitative feature of the present reacting system with others. Therefore, even though the results might be specific to such reacting systems as the present one, they are very instructive in understanding the general feature in liquid-phase chemical reactions. The following section provides the methodology for the CRMD simulations. Results and Discussion are given in section 3, followed by the Conclusion in section 4.
2. CRMD Method The energy-transfermechanism in the proton-transfer reaction of the FW supermolecule in aqueous solution was studied by executing CRMD simulations. One FW supermolecule and 61 medium water molecules were included in a cubic box with one side L = 12.43 A, and the periodic boundary conditions were implemented in order to overcome the problem of surface effects.33 A set of classical Hamiltonian equations of motion34 were integrated for all the atoms in the cubic box by the Verlet a l g ~ r i t h mwith ~ ~ ,a~time ~ step of 0.5 fs. A. Potential Functions and Periodic Boundary Conditions. In order to describe analytically the FW intrasupermolecular potential energy surface in the neighborhood of the intrinsic reaction ~ o o r d i n a t and e ~ to ~ facilitate ~ ~ ~ ~ ~the~ trajectory calculations, the analytic function obtained previously38by the heuristic potential function method39was used:
where rg is the distance between the atoms labeled i and j and uno are the fitted parameters. The numbering of the atoms is shown in Figure 1. Similarly, the interaction energy between the FW supermolecule and a medium water molecule was represented by using the intermolecular potential function fitted previously28by the least mean square method:
integration, we employed the potential truncation technique33 in order to make the force calculations manageable, and the cutoff distance was chosen to be 1.5 of the side L of the cubic box because the forces between the FW supermolecule and a medium water molecule were almost zero at this distance. The present functional form need not suffer from problems with respect to discontinuous forces at the potential cutoff distance owing to the fast falling-off of force intensities as the medium water molecule leaves from the FW supermolecule. For the interaction energy between medium water molecules, the Matsuoka-Clementi-Yoshimine (MCY) potential funct i ~ was n ~ used. ~ Since it adopts virtual sites, Le., nonatomic sites, as centers of interaction, the forces on the virtual sites were redistributed to the atoms in the medium water molecule in the simulations.44 Although the rigid-body appr~ximation~~ is assumed in the MCY potential function, the intramolecular motions of a water molecule have much higher frequencies than intermolecular motions, and thus we can adopt the rigid-body approximation for the present purpose. In fact, Gertner et aL20 showed that the intramolecular high-frequency vibrations of a water molecule do not directly contribute any of the energy required to drive the reactants up to the barrier. During the CRMD simulations, therefore, internal degrees of freedom of the medium water molecule were restricted to its stable state (SS) structure by using the SHAKE a l g ~ r i t h m . ~ ~ . ~ ~ Further, for numerical calculation of interaction energies between medium water molecules, the Ewald summation t e c h n i q ~ ewas ~ ~ .employed ~~ to take into account the long-range Coulomb interactions justifiably. The Coulomb interaction energy between medium water molecules under the periodic boundary conditions is, therefore, described as @W-W
2
- ep * w - w
+ @2w-w + (P3W-w)
erfc{Xln
(4)
+ rj* - rj*l} ( 5 )
N
with
[ x Z i s i n ( 2 n h ~ , * ) ] ~(6) } bl2ij
'0
(3)
where bmij are the fitted parameters. These functions were prepared so as to reproduce a set of energies calculated by the GAUSSIAN86 programa with the minimal STO-3G basis see1 at the Hartree-Fock level of the01-y.~~It should be noticed that the calculation by the STO-3G basis set was found previously to give correct barrier heights,28 compared with other studies using the MP4/6-31 lG(d,p) wave functions.24d Although the long-range Coulomb interaction energy should be in inverse proportion to the distance between a couple of charge centers, our intermolecular potential function does not include such a term proportional to rij-l. The reason why we selected the above functional form is that the potential function using a term proportional to rjj-l instead of r0-3 was found to poorly represent the correct potential energy surface. It is a direct result from the fact that the FW supermolecule is rather neutral even if a medium water molecule is polar. In fact, the present potential function reproduces satisfactorily even the long-range interaction energy between the FW supermolecule and a medium water molecule within a cubic box of 22 8,each side.28 For numerical
(7) where e is the elementary electric charge, Zi and Zj are the ionic valencies of atoms i and j , and N is the total number of atoms in medium water molecules. The sum over n is taken over all simple cubic lattice points (n,L, nyL, nJ), with n,, ny, and n, integers, ri* denotes a scaled length such as rdL where ri shows the position of the ith atom in the cubic lattice point n = (0, 0, 0). The reciprocal vector h is defined as (h,, hy,hJ with h,, h,, and h, integers. The prime indicates that we omit the case i = j for n = 0 or h = 0. x can be chosen arbitrarily and was taken presently as 6.0 with the condition that the summations were truncated at In/ = 0 and lhI2 = 108. These were determined by considering the balance between the necessity to reproduce the reliable interaction energy and the requirements of computational feasibility. B. Preparation of Initial Configurations and Velocities. The initial configurations of the system should be first deter-
Proton-Transfer Reaction of Formamidine in Aqueous Solution
J. Phys. Chem., Vol. 98, No. 48, 1994 12509
mined before starting trajectory calculations of the reaction. Since a chemical reaction event is very rare in the simulation time scale, the following technique was commonly used;32the trajectories are started from the top of the barrier with respect to the reaction coordinate, and all the other degrees of freedom are chosen from a thermal equilibrium distribution. This technique is based on the assumption that the distribution in the transition region can be described by the equilibrium distribution. The validity of this assumption has been discussed by Anderson.32c However, it is difficult to obtain numerically the equilibrium distribution on the dividing surface, because it is impossible for the reaction coordinate to be fixed correctly to the TS value in the FW supermolecule without any restriction of all the other degrees of freedom. The difficulty stems from the numerical error of the separation of overall rotation from internal vibrations of the N-body system47and anharmonicity of the potential energy surface in the neighborhood of the TS. Therefore, our equilibrium calculations were performed under the condition that all the internal degrees of freedom of the FW supermolecule were fixed to the TS structure. Strictly speaking, although the configurations obtained by this technique must be a little ill-balanced distribution, we believe that the trajectory calculations starting from these configurations should provide the essential and qualitative picture of the energy flow. Restriction of all the internal degrees of freedom of the FW supermolecule to their TS values was implemented also by using the SHAKE a l g ~ r i t h m . ~The ~ . ~constant-temperature ~ scheme2228 was employed for the initial calculations for equilibration in order to decide initial conditions from a canonical ensemble. The kinetic temperature of the FW supermolecule was set at 300 K by the velocity scaling m e t h ~ d . ~The ~ .initial ~ ~ configuration was first determined so that after the center of mass of the FW supermolecule was placed at the center of the box and the positions of the 64 medium water molecules were determined temporarily to form a uniform grid they were superimposed in the same cubic box. Then, the three medium water molecules nearest to any atoms of the FW supermolecule were eliminated from the cubic box. As a result, the system becomes composed of 1 FW supermolecule and 61 medium water molecules in the cubic box. It should be noted that the length of one side, L, in the cubic box was determined so that the density of the 64 water molecules in the cubic box corresponds to the density of real water at 300 K, Le., 1.0 g/cm3. The preparation of initial conditions consists of two steps. The fxst-step calculations were performed preliminarily for equilibration and were started with the velocities chosen from the Maxwell-Boltzmann distribution. Five configurations were recorded every 1 ps from 16 to 20 ps. These were used as seed-configurations for the second-step equilibration calculations. For each seed-configuration, independent 10 sets of velocities were generated randomly so that they should obey the Maxwell-Boltzmann distribution. The second-step equilibration calculations were carried out for 1.O ps. As a result, 50 initial conditions (=5 x 10) specifying the positions of all the atoms were obtained. C. CRMD Simulations. With each initial configuration obtained by the above equilibration calculations, 50 sets of trajectory calculations of the reaction were performed in the constant-energy scheme, Le., in the microcanonical way. All the internal degrees of freedom of the FW supermolecule, restricted initially to the TS values, were released to initiate the trajectory calculations. The initial translational and rotational velocities of the medium water molecules as well as those for all the atoms of the FW supermolecule were chosen from the Maxwell-Boltzmann distribution with the initial temperature To = 300 K. Each trajectory calculation was canied out for
0.5 ps. It should be, here, noted that any quantum effect is not considered in the present study. However, the tunneling effect for the present reaction system was discussed in our previous study,26where the tunneling probability for each of the 'H- and 2H-transfer reactions was compared. The following quantities were recorded every 0.5 ps: the total potential and kinetic energy; the internal vibrational, translational, and rotational energy of the FW supermolecule; the lengths of the 02-H4, 02-H5, H4-N6 and H5-N7 bonds; the kinetic energy of each atom of the FW supermolecule; the total accumulated work done on the FW supermolecule by the medium water molecules and its x, y, and z components (we will define and explain these in the following subsection); and the intermolecular potential energy between the FW supermolecule and the medium water molecules. D. Componential Work Analysis. According to the study by Ohmine,18 the accumulated work, AWj, done on the FW supermolecule by an atom j in medium water is defined as
A y ( t ) = L'Cf,.(t).v,(s) d t
(8)
V
where t is an instant time in the duration of the reaction, fvjis the force vector exerted on atom v of the FW supermolecule by atom j in the medium water, and vy is the velocity vector of atom v. In order to investigate the direction of energy flow in threedimensional space, we defined a body-fixed frame in the FW supermolecule and calculated the x, y, and z components of the accumulated work on the FW supermolecule done by the medium water molecules. This frame was set up initially for the TS structure of the FW supermolecule, its origin coinciding with the center-of-mass. The fundamental vector along the x axis is defined as
(9) where
r12= r2 - rl
(10)
rl and r2 are the position vectors of the C1 and O2 atoms, respectively (Figure 1). Those vectors along the y and z axis were defined similarly as e, =
r12 Ir12
r23
x
r23l
and
e, = e, x ey
(12)
respectively, where
r23= r3- r2 and r3 is the position vector of the H3 atom. We have selected these fundamental vectors of the body-fixed frame so that e, and e, should span the molecular plane, on which the reactive motion mainly occurs. The time-dependent rotation of the bodyfixed frame with respect to the laboratory-fixed frame was described conveniently by the procedure introduced previ0usl9~ in relation to the problem of separation of the energy of overall rotation. Additionally, the internal vibrational energy, the rotational energy, and the translational energy of the FW supermolecule were also calculated.
Nagaoka et al.
12510 J. Phys. Chem., Vol. 98, No. 48, I994
Figure 2. (a) Stable state geometry optimized by GAUSSIAN86 with the STO-3Gbasis set at the Hartree-Fock level of theory and (b) the geometry at 0.5 ps for a representative trajectory.
0.005.
.222.09 0.04
O.OO0 0.0
0.1
0.2
0.3
0.4
0.5
-222.09
A
8-
2"
2L a5l - mr a-
?
~
Time (DS) Figure 4. Internal vibrational energy curve (thick solid line), translational energy curve (thin solid line), and rotational energy curve (solid line) of the FW supermolecule. Energies were averaged over 50 trajectories.
.222.10
.2?2.10
a01
&:
-- -
-222.11
0
E0 ) 0 2
.222.11
n ; + .222 1 1
.222.12
'
I
0
01
02
0.3
0.4
0.01
0.5
l i m e (ps)
Figure 3. Potential energy curve (thin line), kinetic energy curve (thick solid line), and total energy curve (solid line) of the F W supermolecule. Energies were averaged over 50 trajectories.
Then, we define the component of work along an axis f (=x, y , or z ) as
wherefvic and vvc are the f components of the force vector fv, and the velocity vector vv, respectively.
3. Results and Discussion We observed the latter half of the whole reaction, Le., the relaxation process of the reaction energy. It should be pointed out that the former-half reaction could be obtained also by the time-reversal procedure. This is due to the time-reversal invariance of the MD simulation and the geometrical symmetry of the reactant and product of the present system. Therefore, we will discuss the energy-transfer mechanism in the latter half of the whole reaction. Most of the results were averaged over all 50 trajectories. The total time used to follow trajectories, 0.5 ps, is reasonable to discuss the present early time dynamics From TS because the present proton-transfer reaction of formamidine finishes completely at 0.5 ps of a representative trajectory (Figure 2). A. Energy Relaxation and Time-Scale Hierarchy. The relaxation of the potential energy of the FW supermolecule was found to be very rapid (Figure 3). The rapid increase of the potential energy at 5 fs might be an artifact due to the uncertainty of setting the initial configuration, discussed already in the previous section. Most of the energy required to reach the barrier top is relaxed in the first 0.05 ps, and 5 kcaVmol is lost in the potential energy of the FW supermolecule. This amount was transformed into the kinetic energy of the FW supermol-
ecule (Figure 3). Then, the longer-wavelength oscillations,in both the potential and the kinetic energy of the FW supermolecule continue for 0.15 ps mutually with the inverse phase. The strong correlation of these two oscillations is clear evidence of the mutual transformation of the two energy modes. As a whole, the total energy of the FW supermolecule decreases smoothly (Figure 3). It means that dissipation of the reaction energy into the medium is accompanied with the mutual transformation of these energy modes of the solute. The averaged slopes of these energy curves still remain negative at 0.5 ps. This indicates that the FW supermolecule still stores such energy to be dissipated. This excess energy should be stored almost in the intrasupermolecular vibrational modes because the time scale for vibrational kinetic energy relaxation is expected to be much longer than those for rotational or translational energy. This was verified by comparing the vibrational, rotational, and translational energies of the F W supermolecule (Figure 4). The negative slope was observed at 0.5 ps only in the vibrational kinetic energy curve. However, the excess vibrational kinetic energy is insufficient obviously to react again. It is important to point out that the kinetic energy of the FW supermolecule consists mainly of the internal vibrational kinetic energy. Some amount of the internal energy is transformed into the rotational kinetic energy which increases at 0.05 ps and shows larger fluctuation than the translational kinetic energy (Figure 4). The fluctuation of the rotational kinetic energy must be induced by the interaction between the FW supermolecule and the medium-water molecules. Some reaction energy is expected to be transformed into the medium energy through the solute-solvent interaction caused by the rotational motions of the FW supermolecule. On the other hand, the translational energy of the FW supermolecule hardly changes during the course of the reaction (Figure 4). However, this feature is not always the case for each trajectory because the averaging over 50 trajectories smooths over the large fluctuation seen in individual trajectories. B. TS Dynamics. Time evolution of the lengths of the 02H4, 02-H5, H4-N6, and H5-N7 bonds is shown in Figure 5. We have regarded the 02-H4 bond as the bond producing a water molecule in the product region. The fate of bond formation was decided immediately after passing through the
Proton-Transfer Reaction of Formamidine in Aqueous Solution
J. Phys. Chem., Vol. 98, No. 48, 1994 12511 281
-
24
-
- 02.H4 bond length
26
02-H4 bond length 02-H5 band length H5-N7 bond lenglh H4.N6 bond lenglh
24 h
s
22
5 0)
20
C
18
U
16
3 c 0 m
I
02-H5 bond length H5-N7 bond length H4-N6 bond length
14
12 10
t
_ _0 0
i 01
02
03
04
05
Time (PS) Figure 5. Time evolution of the lengths of the 02-H4 and 02-H5 bonds, and that of the lengths of the H5-W and H4-N6 bonds. Lengths were averaged over 50 trajectories.
TS region, and there was no recrossing at the barrier top. This should be explained in solute-solvent interaction terms. However, intuitively, since the charge distribution of the FW supermolecule hardly changes during the course of the reaction,28 it is easily understood that any configurational reorganization of the medium water molecules is not necessary for the reaction occurrence. Thus, any force against the forward motion along the reaction coordinate, which is induced by the solute-solvent interaction, should be substantially small. In other words, for the present reaction, the frictional effect is small. Then, the forward motion along the reaction coordinate should not be hindered by the opposite force, and any recrossing at the barrier top does not occur. This is in good agreement with previous s t u d i e ~ ,where ~ ~ . ~the ~ frictional influence was estimated for the reaction rate in the FW-water cluster system on the basis of the generalized Langevin equation for motion along the reaction coordinate. The present result is contrasted remarkably with other studies8s20where the solvent causes recrossing through the dividing surface by the strong solutesolvent interaction. Since the recrossing is caused by the change of the solute-solvent interaction during the course of reaction, no recrossing is a result intrinsic to such solution reactions that the solute-solvent interaction hardly changes during the course of reaction. To investigate the vibrational relaxation of the 02-H4, 02H5, H4-N6 and H5-N7 bonds, each trajectory was analyzed separately because fine structures of oscillatory behaviors in the time-evolution of bond lengths should disappear after averaging over all trajectories. The time evolution of the bond lengths for a representative trajectory is shown in Figure 6. The lengths of the bonds, H4-N6 and 02-H5, cleaved in the relaxation process, are retained around 1.7 and 2.0 A with a large oscillation after 0.05 ps, respectively. This separation is expected not to be sufficient for some medium water molecules to intervene between the reactive water molecule and the formamidine molecule. In fact, the reactive water molecule continued to be coupled with the formamidine molecule for a long time after reaching the product region from the barrier region. This represents clear cage formation by the surrounding medium water molecules which are hardly involved in proton exchange and reflects evidence that the characteristic time of diffusion in liquid water is ca. lo-" s48 and is long enough for
08' 00
I 01
02
0.3
05
04
Time (PSI Figure 6. Time evolution of the lengths of the 02-H4 and 02-H5 bonds, and that of the lengths of the H5-N7 and H4-N6 bonds for a representative trajectory. O.oO€
0.00:
0.W
0.003
0.002
0.001
OOOO
~
0.5
Time (ps)
Figure 7. Time evolution of the kinetic energy of the O2 atom. Energies were averaged over 50 trajectories.
the reaction time of the present proton-transfer reaction. This fact also justifies that supermolecule treatment can be employed to investigate the solvent effect as in our previous studies.22,26-31 C. Intrasupermolecular Energy Redistribution. The kinetic energies of the C1, 02,H4, H5, N6, and N7 atoms, which are involved in bond formation or cleavage (Figure 1), were found to increase at the incipient stage of energy relaxation. (The time evolution of the kinetic energy of the O2 atom is shown in Figure 7, because those of the other atoms are almost the same.) The increase of the kinetic energy of the O2 atom is prominent and contributes to the motion of water molecule formation accompanying H3 and H4 atoms. However, as explained in B, the distance between the reactive water and the formamidine molecules was found within 2.2 A during the course of the reaction (Figure 5 ) . The initial kinetic energy increase of these atoms at one-score fs is followed by the gradual decrease in the kinetic energy. On the other hand, the kinetic energy of the remaining atoms, H3, H8, H9,and HIo, increases steadily with oscillations. The combination tendencies lead to the equipartitioning of the kinetic energy among all the atoms
Nagaoka et al.
12512 J. Phys. Chem., Vol. 98, No. 48, 1994 1
00151
b
m
0 005
s3
0 000
5
.O 005
.0.010
t
-0015' 00
I 01
02
03
04
05
Time (ps)
Figure 8. Accumulated works done on the FW supermolecule along thex (thin solid line),y (solid line), and z (thick solid line) axes by all the medium water molecules. Works were averaged over 50 trajectories.
of the FW supermolecule. The time-lag in the energy redistribution may be caused by the vibrational frequency differences that prevent the intrasupermolecular energy from flowing rapidly. On the other hand, as shown already in Figure 3, the kinetic energy redistribution between the reacting system and the medium has not been accomplished by 0.5 ps. The energy flow within the FW supermolecule is a little faster than that between the reacting system and the medium. This indicates that the F W supermolecule itself plays a role as a temporary reservoir at the first stage of energy relaxation and it continues for a long time until the reaction energy relaxes into the medium completely. Finally, the equilibrium energy distribution is attained not only within the reacting system but also in the medium system. In general, the feature of energy relaxation should be dependent significantly on the number of atoms or the degrees of freedom of the reacting molecules. Although several simplified model reacting systems, where some intemal degrees of freedom are fixed, were employed often in previous works investigating energy one should treat all the intemal degrees of freedom of the reacting molecule for understanding more realistically and precisely the energy flow in chemical reactions in solution. D. Energy Flow. The accumulated work done on the F W supermolecule by all the medium water molecules is equivalent to the change in the total FW energy (Figure 3). The timedependent changes of the x, y, and z components of the accumulated work are shown in Figure 8. The x and y components become negative after the initial 0.1 ps. Since negative work corresponds to the dissipation of F W energy into the medium energy, the x and y components contribute to the relaxation of the reaction energy. The motion along the reaction coordinate occurs almost on the x-y plane, and it is natural that the x and y components of force exerted on the F W supermolecule by all the medium water molecules play such a role as a driving force for reaction occurrence. On the other hand, it is an interesting phenomenon that the z component acts inversely and provides positive energy to the FW supermolecule during the relaxation process. A part of the internal energy of the FW supermolecule is supplied from medium water molecules along the z axis, while the FW supermolecule dissipates its energy along the x and y axes. The balance between them leads
to the dissipation of the total F W energy into the medium, and the relaxation process is completed. To understand the energy flow mechanism in more detail, we have determined the spatial distribution of such medium water molecules that provide the most and least (Le., the most positive and most negative) work against the F W supermolecule. The medium water molecules providing the least (most negative) work in each trajectory are superimposed in Figure 9a. In other words, the medium water molecules providing the greatest contribution to energy dissipation in each trajectory are superimposed in the figure. Those medium water molecules are coupled with any of the 02,H3, and Hs atoms. Their interaction is found to contribute to energy stabilization because these couplings of the FW supermolecule and the medium water molecules correspond primarily to those in any of the previously optimized geometries of the FW(TS)-W(SS) complex (the geometries A, B, C, and D in Figure 2 in reference 28). The least work by the medium water molecules coupled strongly with the H8 atom can be explained by considering the electrostatic interaction between the oxygen atom of the medium water molecule and the H4 atom. As the H4 atom is dissociated from the N6 atom, the distance between the H4 atom and the oxygen atom of the medium water molecule becomes longer and its energy is dissipated into the medium water molecule through interaction between the H4 atom and the oxygen atom of the medium water molecule. On the other hand, the medium water molecules that provide the most work are localized mainly in positions coupled with the H9 atom (Figure 9b). It can be conjectured that the medium water molecules that have the most influence on the reaction, whether positively or negatively, should be localized in such regions as optimized for a medium water molecule by supermolecule treatment.28 Furthermore, it is remarkable that most of the medium water molecules providing the most and least work for the reaction progress are localized in the regions of those protons, i.e., the H8 and H9 atoms, which are far from the reaction site. It was observed that the interaction energy between the FW supermolecule and the medium water molecules decreases slightly in the relaxation process (Figure 10): Although this qualitative nature is similar to the study by Gertner et a1.,20the amount of decrease in the energy of the present system is much smaller than that in their study. The difference could be attributed to the fact that the charge distribution of the F W supermolecule hardly changes during the course of the reaction in the present reacting system.28 Therefore, solvent reorganization to solvate the solute is not necessary for the present reaction to occur. The present energy source of dissipation comes mostly from the relaxation energy. On the contrary, in the reacting system studied previously,2othe relaxation energy of the solutesolvent interaction contributes substantially to the energy source of dissipation. E. Summary of the Medium Effect. In our previous s t ~ d i e sa, medium ~ ~ ~ ~ water ~ molecule was found to have little influence on the potential energy barrier height because there is only a small change in the charge distribution between the TS and S S of the FW system. Then, one can anticipate that the surrounding medium waters hardly change the potential energy barrier itself. In fact, in the present simulation, it was found that only the small change has been observed in the interaction energy between the FW and the medium water during the course of the reaction (Figure 10). Thus, the static solvent effect, ie., the effect of the medium waters on the potential energy barrier, has been concluded to be small in the present reaction system.
Proton-Transfer Reaction of Formamidine in Aqueous Solution
J. Phys. Chem., Vol. 98, No. 48, 1994 12513
n
b
Figure 9. (a) Superimposed pictures over 50 trajectories for the position of the medium water molecule that did the least work on the FW supermolecule in each trajectory. The left and right pictures were drawn with the viewing direction along the z and y axes, respectively. (b) Similar pictures for the medium water molecule that did the most work on the FW supermolecule.
Many studies8 have pointed out that the frictional effect caused by the medium significantly reduces the rate constant. In the liquid-phase chemical reactions, the rate constant predicted by the transition state theory, ~TST,is modified by the transmission coefficient K : ’
where the transmission coefficient K, which is related to the frictional coefficient,8 is recognized as the reduction ratio of kTsT via the dividing surface recrossing. However, it has been shown in the preceding subsection B that there is no recrossing in the present reaction system. The result is in good agreement with our previous supermolecule study30*31 where the frictional effect was found to be quite small. Since the recrossing is caused by the change of the solute-solvent interaction during the course of reaction, no recrossing should occur in such a reaction that the charge distribution of the solute does not change during the course of reaction. One can certainly understand that the F W system has this characteristic.
In the present study, although there is no solvent effect which originates in the above two mechanisms, it is true that a role of the medium waters was significantly clarified on energy transfer in the relaxation. The electrostatic interaction of hydrogen bonds is strong and is involved in the present energy flow. Then, the energy transfers through the directions of the hydrogen bonds between the FW and the medium waters. 4. Conclusion
The mechanism of reaction dynamics, i.e., how the energy required to surmount the potential energy barrier is dissipated into the surrounding solvent molecules, has been clarified by analyzing the CRMD simulations for the proton-transfer reaction of formamidine in aqueous solution. Since the present reacting system was more complicated and realistic than those investigated previously by the MD method,I8-*I the study here should provide several essential pictures for understanding the liquidphase chemical reactions. In particular, we have focused on the time dependence of different energy modes within the FW supermolecule,the TS region dynamics, the spatial distribution
Nagaoka et al.
12514 J. Phys. Chem., Vol. 98, No. 48, 1994
Acknowledgment. We thank Professor Kenichi Fukui for helpful advice and discussions. The numerical calculations were performed by CONVEX C201 and C3210 at the computational facilities of LFC and by FACOM VP4OOE and VP2600 at the Data Processing Center of Kyoto University. This work was supported by a Grant-in-Aid for Scientific Research from the Ministry of Education, Science, and Culture in Japan.
.o 005
-0010'
-0015'
References and Notes %
P 0 W
-0.030
."
Yd.,
00
01
02
03
04
05
Time (ps)
Figure 10. Intermolecular potential energy curve between the F W supermolecule and the medium water molecule. averaged over 50 trajectories.
Energies were
of the influential medium water molecules, and the time evolution of the solute-solvent interaction. We have found that the potential energy of the FW supermolecule is transformed into the intemal vibrational kinetic energy at the first step of the relaxation. The next step is energy transfer of the vibrational and rotational kinetic energies into the medium. These energies were stored at the first step as excess energy of the supermolecule. On the other hand, the total FW energy decreases monotonously, and then the relaxation of the FW vibrational kinetic energy continues for a long time. This qualitative feature of energy flow is in good agreement with other studies. 19,20 Any recrossing at the TS region has never been found over all trajectories. This result supports the previous study30where the transmission coefficient which represents the rate of recrossing was determined by evaluating the frictional effect on the reaction in the cluster consisting of the FW supermolecule and a single medium water molecule. The validity of the supermolecule treatment has been confirmed. Characteristically, it has been found that, while the forces exerted on the FW supermolecule by all the medium water molecules contribute to the dissipation of the reaction energy, its component vertical to the FW molecular plane tends to hinder the relaxation. The balance between the energy flowed into from the vertical direction and that flowed out to any direction within the plane results in the total dissipation of the reaction energy into the medium. The most influential water molecule in energy flow has always been found to be in the neighborhood of any of the stable positions with respect to the FW supermolecule, and they are also localized in the regions coupled to those protons which are far from the reaction site. The change of the interaction energy between the FW supermolecule and the medium water molecules has been found to be small during the course of the reaction. This is because the F W charge distribution hardly changes during the reaction. Then, any reorganization of the medium water molecules is not induced. It is true that this result is specific to such reacting systems as the present one. Although our results may be specific to the present reaction system, the analyses and the discussion should provide some useful perspectives even in the other reactions. Further detailed investigation of the energy-transfer mechanism is in progress and will be discussed in the following article.49
(1) Mukaiyama, T. in Challenges in Synthetic Organic Chemistry; translated by Baldwin, J. E.; Clarendon Press: Oxford, 1990. (2) Zubay, G. Biochemistry; Macmillan Publishing Co.: New York, 1988. (3) Moreau, M.; Turq, P. in Chemical Reactivity in Liquids; Plenum Press: New York, 1988. (4) Popielawski, J., Ed. in The Dynamics qf Systems with Chemical Reactions; World Scientific: London, 1988. (5) Bagchi, B.; Krishnan, V., Eds. in Solvation Dynamics & Charge Transfer Reactions; World Scientific: Singapore, 1991. (6) Sathyamurthy, N. Reaction Dynamics; Narosa Publishing House: New Delhi, 1991. (7) van Gunsteren, W. F.; Weiner, P. K., Eds. in Computer Simulation of Biomolecular Systems; ESCOM Science Publishers B.V.: Leiden, 1989. (8) Hynes, J. T. in Theory of Chemical Reaction Dynamics; Baer, M., Ed.; CRC Press: Boca Raton, FL, 1985; Vol. IV and references therein. (9) Chandrasekhar, J.; Smith, S. F.; Jorgensen, W. L. J . Am. Chem. SOC.1984, 106, 3049. (10) Huang, J.-K.; King, G.; Creighton, S.; Warshel, A. J . Am. Chem. SOC.1988, 110, 5291. (11) Tapia, 0.;Lluch, J. M.; Cardenas, R.; Andres, J. J . Am. Chem. SOC.1989, 111, 829. (12) Huston, S. E.: Rossky, P. J.; Zichi, D. A. J . Am. Chem. SOC.1989, I l l , 5680. (13) Levine, R. D.; Bemstein, R. B. Molecular Reaction Dynamics and Chemical Reactivity; Oxford University Press: Oxford, 1987. (14) Fukui, K. J . Phys. Chem. 1970, 74, 4161. (15) Miller, W. H.; Handy, N. C.; Adams, J. E. J . Chem. Phys. 1980, 72, 99. (16) Marcus, R. A. J . Chem. Phys. 1966, 45, 4500. (17) Mezey, P. G. Potential Energy Hypersurfaces: Elsevier: New York, 1987. (18) Ohmine, I. J. Chem. Phys. 1986, 85, 3342. (19) Benjamine, I.; Gertner, B. J.; Tang, N. J.; Wilson, K. R. J . Am. Chem. SOC.1990, 112, 524. (20) Gertner, B. J.; Whitnell, R. M.; Wilson, K. R.; Hynes, J. T. J . Am. Chem. SOC.1991, 113,14. (21) Li, Y. S.; Whitnell, R. M.; Wilson, K. R.; Levine, R. D. J . Phys. Chem. 1993, 97, 3647. (22) Nagaoka, M.; Okuno, Y.; Yamabe, T. J . Am. Chem. SOC.1991, 113, 769. (23) Heidrich, D.; Kliesch, W.; Quapp, W. Properties of Chemically Interesting Potential Energy Surfaces (Lecture Notes in Chemistry, Vol. 56); Springer-Verlag: Berlin, 1991; p 176. (24) (a) Poirier, R. A.; Majlessi, D.; Zielinski, T. J. J . Compur. Chem. 1986, 7, 464. (b) Kinasiewicz, W.; Les', A,; Wawer, I. J. Mol. Struct. (Theochem.) 1988, 1, 168. (c) Kwiatkowski, J. S.; Bartlett, R. J.; Person, W. B. J . Am. Chem. SOC.1988,110,2353. (d) Nguyen, K. A,; Gordon, M. S.; Truhlar, D. G. J . Am. Chem. SOC.1991, 113, 1596. (e) Moog, R. S.; Maroncelli, M. J . Phys. Chem. 1991, 95, 10359. (f) Wang, X.; Nichols, J.; Feyereisen, M.; Gutowski, M.; Boatz, J.; Haymet, A. D. J.; Simons, J. J . Phys. Chem. 1991, 95, 10419. (g) Andres, J.; Beltran, A,; Carda, M.; Krenchl, J.; Monterde, J.; Silla, E. J . Mol. Struct. (Theochem.) 1992, 254, 465. (25) Yamashita, K.; Kaminoyama, M.; Yamabe, T.; Fukui, K. Theor. Chim. Acta. (Berlin) 1981, 60, 303. (26) Yamabe, T.; Yamashita, K.; Kaminoyama, M.; Koizumi, M.; Tachibana, A.; Fukui, K. J . Phys. Chem. 1984, 88, 1459. (27) Nagaoka, M.; Okuno, Y.; Yamabe, T. Chem. Phys. Lett. 1992,196, 197. (28) Nagaoka, M.; Okuno, Y.; Yamabe, T.; Fukui, K. Can. J . Chem. 1992, 70, 371. (29) Nagaoka, M.: Okuno, Y.; Yamabe, T. J . Chem. Phys. 1992, 97, 8143. (30) Nagaoka, M.; Okuno, Y.; Yamabe, T. (to be submitted). (31) Nagaoka, M.; Okuno, Y.; Yoshida, N.; Yamabe, T. Int. J . Quantum Chem. 1994, 51, 519. (32) (a) Keck, J. C. Discuss. Faraday SOC.1962, 33, 113; Adv. Chem. Phys. 1967, 13, 85. (b) Bennet, C. H in Diffusion in Solids; Burton, J. J., Ed.; Academic: New York, 1975; p 73. (c) Anderson, J. B. J . Chem. Phys. 1973, 58, 4684; Ibid. 1975, 62, 2246. (d) Smith, I. W. M. J . Chem. Soc.,
Proton-Transfer Reaction of Formamidine in Aqueous Solution Faraday Trans. 2 1981, 77,747. (e) Truhlar, D. G.; Garrett, B. C. Faraday Discuss. Chem. SOC. 1987, 84,464. (33) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: Oxford, 1987. (34) Goldstein, H. Classical Mechanics, 2nd ed.; Addision-Wesley: Reading, MA, 1980. (35) Verlet, L. Phys. Rev. 1967, 159, 98. (36) (a) Fukui, K. in The World of Quantum Chemistry; Daudel, R.; Pullman, B., Eds.; Reidel: Dordrecht, 1974; p 113; (b) Acc. Chem. Res. 1981, 14, 363. (37) Baker, J.; Gill, P. M. W. J . Comput. Chem. 1988, 9, 465. (38) Nagaoka, M.; Okuno, Y.; Yamabe, T.; Fukui, K. In?. J . Quantum Chem. 1992, 42, 889. (39) Swaminathan, S.; Whitehead, R. J.; Guth, E.; Beveridge, D. L. J . Am. Chem. SOC. 1977, 99, 7817. (40) Frisch, M. J.; Binkley, J. S.; Schlegel, H. B.; Raghavachari, K.; Melius, C. F.; Martin, R. L.; Stewart, J. J. P.; Bobrowicz, F. W.; Rohlfing, C. M.; Kahn, L. R.; Defrees, D. J.; Seeger, R.; Whiteside, R. A,; Fox, D.
J. Phys. Chem., Vol. 98, No. 48, 1994 12515 J.; Fleuder, E. M.; Pople, J. A. GAUSSIAN86 Camegie-Mellon Quantum Chemistry Publishing Unit: Pittsburgh, PA, 1984. (41) Hehre, W. J.; Stewart, R. F.; Pople, J. A. J . Chem. Phys. 1969,51, 2657. (42) Roothaan, C. C. J. Rev. Mod. Phys. 1951, 23, 69. (43) Matsuoka, 0.;Clementi, E.; Yoshimine, M. J . Chem. Phys. 1976, 64, 1351. (44) Berendsen, H. J. C.; van Gunsteren, W. F. in Molecular liquids-Dynamics and interactions; Barens, A. J., et al., Eds.; D. Reidel Publishing Company: Dordrecht, 1984; p 475. (45) Rychaert, J. P.;Ciccotti, G.; Berendsen, H. J. C. J . Comput. Phys. 1977, 23, 327. (46) Ewald, P. Ann. Phys. 1921,64, 253. Nijboer, B. R. A,; De Wette, F. W. Physica 1957, 23, 309. Ceperley, D. Phys. Rev. B 1978, 18, 3126. (47) Jellinek, J.; Li, D. H. Phys. Rev. Lett. 1989, 62, 241. (48) Eisenberg, D.; Kauzmann, W. in The Structure and Properties of Water; Clarendon: London, 1969. (49) Nagaoka, M.; Okuno, Y.; Yamabe, T. (to be published).
