J . Phys. Chem. 1993,97, 59-69
59
Cluster and Solution Simulation of Formaldehyde-Water Complexes and Solvent Effect on Formaldehyde l(n,r*) Transition Hiroo Fukunaga Ashigara Research Laboratories, Fuji Photo Film Co. Ltd.. Minamiashigara 250-01, Japan, and Institute for Molecular Science, Myodaiji, Okazaki 444, Japan
Keiji Morokuma’ Institute for Molecular Science, Myodaiji, Okaraki 444, Japan Received: July 13, 1992; In Final Form: October 12, 1992
On the basis of a b initio molecular orbital calculations, we have derived 12- 1-3 potential functions for interaction between formaldehyde, in both the ground and the n-A* singlet excited state, and one water molecule as well as the transition dipole moment functions in the n-A* transition of formaldehyde interacting with a water molecule. Using these functions, we carried out simulations of vibrational structures in the n-A* absorption spectrum of formaldehyde complexed with a few water molecules ( n = 1-3) and in water solution. Detailed analyses of structural and energetic features have been carried out. Formaldehyde in cluster interacts with the chain of hydrogen-bonding water molecules, while formaldehyde in water solution is surrounded by water molecules. Water molecules around the oxygen atom of formaldehyde interact attractively and those around the hydrogen atoms of formaldehyde interact both attractively and repulsively with formaldehyde. Though the transition energy shift in the H2CO-3(H20) cluster does not reach that in water solution, the bandwidth of the simulated absorption spectrum in the HzCO-3(H20) cluster is nearly the same as that in water solution. cluster size and temperature dependence of the absorption spectrum of a chromophore-rare gas cluster have also been To predict the electronic absorption spectrum of a dye, one studied. Bosiger et al. have studied the electronic spectra of would be required to calculate both the transition energy (the 0-0 carbazole-nAr clusters ( n = 1-14) by Monte Carlo simulation excitation energy) and the vibrational structure (or the band and shown that the spectral change due to size and temperature shape). These properties depend on the environment of the dye, is caused by the structural transition between wetting and e.g., whether the dye molecule is isolated, complexed with some nonwetting ~1usters.l~Mukamel et al. performed molecular solvent molecules (in cluster) or in solution. To understand these dynamics simulations of benzene-nAr clusters ( n = 1-21), properties at the electronic and molecular point of view, several showing that the observed spectral feature of inhomogeneous theoretical studies on this subject have been reported. broadening is attributed to the superposition of narrower As for the n-n* transition energy of an isolated formaldehyde subspectra with different shifts of various cluster configurations.14 molecule, several a b initio molecular orbital methods have been On the other hand, Herman and Berne have investigated the employed to reproduce the experimental The n-n* solvent effect on the vibronic absorption spectrum of Br2 in Ar transition energy shifts in the H2CO-n(H20) clusters ( n = 1-2) solvent, using a scheme based on the Born-Oppenheimer have been calculated using the supermolecule m e t h ~ d . ~ . ~ approximation, i.e., thevibrations and theotherdegreesof freedom The calculation of the n-n* transition energy shifts in the can be separated.I5 Using a hybrid method, Levy et al. have N02--n(H20)clusters ( n = 1-25) have been also carried out by carried out simulation of the 0-0excitation energy distribution Monte Carlo ~imulation.~ To evaluate the transition energy shift in the n-r* transition of formaldehyde in water solution and of a molecule in solution, two approaches have been applied. One reproduced the solvent energy shift and bandwidth of the is the solvent continuum model, in which the interaction between absorption ~ p e c t r u m . ~ a solute molecule and solvent molecules are considered as a Though as described above, the solvent effect on the transition continuous potential in the H a m i l t ~ n i a n . ~In, ~this model each energy shift and the band shape in an electronic absorption molecule is not treated explicitly, and therefore one cannot obtain spectrum has been extensively studied, we can find only few microscopic information concerning the behavior of solvent investigations which systematically discuss the difference of the molecules. The other is the statistical simulation, such as the absorption spectrum of a chromophore in cluster and in solution. statistical mechanical free energy perturbation method based on To clarify the differences between the absorption spectra in a molecular dynamics,s and the hybrid model of molecular dynamics cluster and in solution and understand their origins, we have thus and ab initio quantum mechanic^.^ In these methods, the performed the simulation of the absorption spectra in the n-n* transition energy shift is calculated as a difference of solvation transition of formaldehyde complexed with some water molecules energies between the ground and the excited states, and the ( n = 1-3) and in water solution. The reason we selected this structural features of the system can also be obtained. transition is that the absorption spectra of the n-n* transition of Investigations on the vibrational structure or the band shape formaldehyde both in the gas phase and in water solution have of the absorption spectrum have also been carried out. The been both experimentally and theoretically studied. Thespectrum vibrational structures for the n-n* transition of formaldehyde in in the n-n* transition of formaldehyde in the gas phase has been the gas phase were calculated by the ab initio molecular orbital well characterized experimentally, Le., the two progressions of method.10,” Molecular dynamics calculations has been used for 4;2;1 and 4:2;1 are characteristic in the spectrum.I6-*O Further the simulation of the line shape of the absorption spectrum for an I2 observed transition energies and the transition intensities have molecule in the gas phase.12 In both studies, the calculated been evaluated quantitatively by the ab initio molecular orbital absorption spectra reproduced the observed ones very well. The 1. Introduction
0022-3654/58/2091-0059%04.00/0 0 1993 American Chemical Society
Fukunaga and Morokuma
60 The Journal of Physical Chemistry, Vol. 97, No. I, 1993 method.10~11 Absorption spectra in several solvents have been also meas~red.2'-~3 In water solution, a structureless broad band with a large blueshift wasobserved. From theab initio molecular orbital calculations for the formaldehyde-water clusters, Morokuma et al.3 and Del Bene24had suggested that the origin of blue shift in the transition in water solution was attributed to breaking or weakening of the hydrogen bond in the n-A* excited state, while Taylor had deduced that the major component of the blue shift arose from the change induced in the component molecule (formaldehyde) geometries by hydrogen bonding.25 In this argument, the former explanation was supported from the recent studies based on the statistical simulation by Levy et aL9 and Kollman et ala8 In the simulation, we employed the combined method of quantum and classical mechanics which consists of two steps: one is the Monte Carlo step to generate the positions and orientations of water molecules around formaldehyde, and the other is to ensemble-average the band intensities evaluated for each configuration sampled by the M C simulation, using the vibrational wavefunctions and the energy levels of formaldehyde interacting with water. We at first derive the formaldehydewater potential functions and the transition dipole moment functions required in the simulation. We briefly describe ab initio calculations of interaction energies between formaldehyde in the ground and excited states and a water molecule and the transition dipole moments in the n-A* transition of the complexes in section 2A and explain potential functions and transition dipole moment functions in section 2B and 2C, respectively. In section 2D, we describe thesimulation method of absorption spectra. We present the fitted potential functions and transition dipole moment functions in section 3A and 3B, respectively. Further, on the basis of the MC simulations both in clusters and in solution, we discuss the structure of water molecules around formaldehyde, the interaction energies between formaldehyde and water and transition energy shifts, and the calculated absorption spectra in sections 3C, 3D, and 3E, respectively. In these sections, we emphasize the differences in the structure of water molecules around formaldehyde as well as in the formaldehyde-water interaction energies between cluster and solution.
2. Method (A) Ab InitioCalculationsof InteractionEnergiesandTransition Dipoles. Interaction energies and electric transition dipole moments of the 1:1 formaldehyde-water complex were calculated with the ab initio molecular orbital method. We employed the Hartree-Fock (HF) and electron-hole potential (EHP)26methods for the ground and the n-A* singlet excited stateof formaldehyde, respectively. The reason why we selected the E H P method, which can be considered a version of single-excitation configuration interaction (SECI), is as follows. The validity of E H P method has been ascertained in the calculation of both H2CO and the H2CO-H20 complex.3 Since the electron correlation is not taken into consideration in either the H F or the E H P method, they are expected to give energies in comparable accuracy to each other. In addition, transition dipole moments between the E H P excited states and the HFground state can be evaluated easily as transition dipole integrals over E H P MOs. We used the Huzinaga-Dunning contracted Gaussian double-{ (DZ) basis2'128 augmented with a set of d-type polarization functions (ad = 0.75 for C and 0.85 for 0)28 in the calculations. BSSE corrections were not considered. The calculations were mainly carried out with GAUSSIAN8229 except for EHP and transition dipole calculations which were done with IMSPACK.30 Since we are concerned with absorption spectra and in order to save computer time, we assume that formaldehyde and water molecules are fixed at their ground-state geometries in the gas phase. For this purpose, we employed the optimized geometry with the above-mentioned basis (RCO = 1.1886 A, RCH= 1.0934 A, and LHCH = 1 16.07') for formaldehyde and the experimental
\
I
i
'@
H
"To'
Y
4
-.
I
4
X
Figure 1. Definition of six intermolecular variables to represent the geometry of the 1:l H2CO-H20complex. The broken lineis the bisector of the water HOH angle.
one (ROH= 0.958 A and LHOH = 104.5')3' for water. Under this frozen geometry approximation, the configurations of the 1:l formaldehyde-water complex can be defined by six intermolecular variables as shown in Figure 1. We have carried out abinitiocalculations for both thegroundstate and the first excited singlet state for the 1:l complex at mesh points of interatomic configurations. For the C-0 distance, R , we selected four values for each configuration between 2.4 and 4.8 A, encompassing the region inside a minimum. For the other geometrical parameters, we chose the following values: ,'O 45', and 90' for 8; Oo, 45O, 90°,and 180' for I#J; ,'O 45', 90°,and 180' for a;,'O 90°,180°, and 270' for 8;,'O 45', 90', and 135' for y. More mesh points have been added in the vicinity of the global minimum (vide infra), for a better description of the potential functions in the statistically more important region. The total number of configurations is 530. (B) Two-BodyPotential Functions. Here we fit the calculated two-body H2CO-HZ0 interaction energies both for the ground and the excited singlet state to analytical potential functions. We express potential functions in the form of sums of atom (or pseudocenter)-atom (or pseudocenter) terms: on H2CO on H20
where the i's and j's denote atoms or pseudocenters on formaldehyde and water, respectively. The pseudocenters (labeled as L hereafter) are trigonally located around the oxygen atom in formaldehyde and tetrahedrally around the oxygen atom in water.32 We write each pairwise potential in a power series of the distance between sites i and j , ri,:
For bi,, we took C, H, and L centers as i and H and L a s j . For ai, and cij, we took C, H, and 0 as i and H and 0 as j . We considered three and six as a candidate for n, because 12-1-332 and 12-1-633functions had been used as potential functions for interaction of formaldehyde in the ground state and water. Since we have no foundation to select from, we have examined characteristics of both functions. The 12-1-6 function gave some negativevalues for aijresulting in phantom minima on the potential surface. On the contrary, the 12-1-3 function behaved rather smoothly without phantom minima. We decided therefore to employ n I= 3, Le., the 12-1-3 function as the potential function. We determined the distance between an oxygen atom and a pseudocenter, ROL,to reproduce the optimized geometries and the interaction energies calculated for the 1:l ground state formaldehyde-water complex by the ab initio method shown in Figure 2. Finally, we found 0.5 A as the optimum value of R ~ L , which was also employed for the potential function in the excited state. A choice of smaller value for RoL,'~giving an incorrect 1:l complex geometry, would influence the results of simulation for clusters. We performed linear least-squares fitting for the configurations with the interaction energy below + 10.0 kcal/mol for the ground
Formaldehyde-Water Complexes Y
The Journal of Physical Chemistry, Vol. 97, No. 1, 1993 61
Y
W e imagine that hfj is induced by each bond dipole of H 2 0 on each atom or pseudocenter of H2C0, and further assume the following form:
1.
0Id lO\H b
19
IC
II
II
11
19
lh
0
Figure 2. Stationary point geometries of the 1:l H2CGH20 complex. Values of geometrical parameters are in Table 11.
state (481 configurations) and below +15.0 kcal/mol for the excited state (476 configurations), using the SALS program.34 We used a weight function, W I = 1 + a exp(-(AEl - AEo)/kT),, in the fitting,32J3 where AEl is the integration energy in the Ith configuration, AEois -4.5 kcal/mol for the ground state and -2.0 kcal/mol for the excited state, kT is 0.6 kcal/mol, and a is 150 for both states. Furthermore, the neutral conditions35 for the parameters bij and ci, are used to guarantee a zero potential for large r:
+ bCH + 2bLL + 2bLH + 2bHL + 2bHH = 0 CCO + 2cCH + + 2cOH + 2cHO + 4cHH = 0
bCL
COO
(3)
(C) TransitionDipole Moment Functions. The transition dipole moment function for formaldehyde interacting with one water molecule is given as Mgf(Q,R) = (~*(q,Q,R)Ir(q)l~.r(Q,Q,R) )p (4) where agand ar are the electronic wavefunctions of the ground and excited state, respectively. r(q)is the electricdipole moment operator. q, Q, and R represent the electron and the nuclear coordinates of the solute H2C0, and the configuration (position and orientation) of the solvent molecule, respectively. The transition dipole moment M,,-(Q,R) can be written as the sum of three terms: Mn(Q),the transition moment of isolated formaldehyde, M*(Qo,R),the effect of the solvent molecule on formaldehyde in the equilibrium geometry Qo, and MnS(Q,R), the term coming from the interaction of vibration of formaldehyde with the solvent coordinate. However, we neglect the last term for practical reasons:
Since the transition moment for the n-r* singlet transition of H2CO in the equilibrium geometry is zero, one can consider Mn(Q) as the sum of Mn(Qk)'s, which were expressed as the polynomial of each normal coordinate, Qk, around the equilibrium geometry, Qo. For the solvent induced transition moment, MS(Qo,R), we justified the functional form upon the following c~nsideration.'~ We a t first assume that the induced transition moment can be written as a sum of atom and pseudocenter terms:
Here k represents the x, y , or z component. rj is the vector from the oxygen atom to t h e j t h center ( H I , H2, LI, and LZ f o r i = 1,2,3, and 4) within HzO, representing a bond dipole. Rij is the vector from the midpoint of t h e j t h dipole in H20 to the ith center (C, 0, HI, Hz, LI, and L2 for i = 1, 2, ..., 6) in H2C0, defined as
R, = Ri - 1/2rj
(8)
where Ri is the vector from the H20 oxygen to the ith H2CO center. W e have found that the first line only in eq 7 cannot fit the calculated changes in the transition moment and that the second and third lines are required. Tensors a,B, y, 6, t, and t are parameters to be determined by the fitting. Thecomponents generated by permutation among 1and m or 1, m and n, such as Ykijxy and Ykijyx, are equivalent each other, and such equivalent components were excluded in the summation. Some of the components are zero by symmetry. (D) Simulation of the Absorption Spectrum. To evaluate the transition energies and intensities between the vibronic states, we need thevibrational wave functions and energies in both the ground state and the singlet excited state of formaldehyde interacting with water molecules a t their given position and orientation, R. For this purpose, we assume the separation of time scale between the vibration of formaldehyde and the other degrees of freedom (the translation and the rotation of water).I5 Under this approximation, the vibrational wavefunction, +i(Q,R), and the energy, Ei(R), can be obtained by solving the following equation:
H(Q,R) $i(Q,R) = Ei(R) +i(QIR)
(9)
where V(Q,R) is the intramolecular potential energy perturbed by the interaction with water molecules, Le., the effective intramolecular potential. To reduce the problem of the six-dimensional vibrational equation into a one-dimensional problem, we consider the effective intramolecular potential as a sum of potentials along the normal coordinate, Qk:
In addition, we assume the normal coordinate of formaldehyde interacting with water molecules to be the same as that of formaldehyde in the gas phase, and approximate the components oftheeffective potential, V(Qk,R),asthesumof theintramolecular potential of formaldehyde in the gas phase, Ps(Qk), and the interaction energy between formaldehyde distorted along the normal coordinate, Qk, and water molecules, vntcr( &R):
on H K O on H,O
In the calculation of the effective intramolecular potential of
62
Fukunaga and Morokuma
The Journal of Physical Chemistry, Vol. 97, No. 1. 1993
Under the above approximations, the vibrational wavefunction,
TABLE I: Relative Intensities in the 2: Progression
$(Q,R), can be written as a product of the wave functions in the mode k , $(Qk,R), and the energy, E ( R ) , as a sum of the energies
calc n
exptl
harmonic
anharmonic
1 .oo 2.30 2.89 2.61 I .90 1.20 0.67 0.35 0.17
1 .oo
1 .oo
2.61 2.89 1.73 0.58
2.34 2.79 2.29 1.49 0.89 0.48 0.21 0.07
0.10 0.00 0.00 0.00
in the mode k, Ek(R):
k= I
formaldehyde in the excited state, we also used the same normal coordinate as that in the ground state. This approximation corresponds to neglecting the Dushinsky effect that the normal coordinates of the excited electronic states can be different from those of the ground state.37 The first term in eq 12 is taken from the observed harmonic frequencies3IJ8 for all modes except for mode 2 and 4 in the excited state. For mode 4 the double-well potential, represented as the sum of gaussian and quadratic terms,j9 is converted into a quartic polynomial for the later calculations, and for mode 2 the anharmonic term is included to reproduce the observed band intensities in the gas phase,20 as indicated in Table I:
The final equation to be solved for each mode is
T;’(Q,) = 0.167Q5’
We solved this equation by expressing the wavefunction, $(Qk,R), as a linear combination of the Hermit polynomials and by diagonalization of the Hamiltonian matrix for each sample of R . To sample the positions and orientations, R, of water molecules around formaldehyde, we used the Monte Carlo method for canonical ensemble based on the Metropolis algorithmS4O The geometry of formaldehyde were held at the equilibrium in the ground state during the simulation. The simulation has been done both for H2CO-n(H20) clusters (n = 1-3) at constant temperatures of 10, 50, 100, and 150 K, and for the system including one formaldehyde molecule and 100 water molecules at 298 and 323 K and at a density of 1 g/cm3 (as a model of solution). In the latter simulation, the periodic boundary conditions were employed for a cubic sample. After equilibration, trial moves of 50 000 cycles are executed in the simulation of both clusters and solution; in each cyclecontaining (n + 1) steps, every molecule was displaced in order by translation in all three directions and rotations about random axis. We used the 12-1-3 potential function determined in section 2B and the MCY potential function4’ to express formaldehyde-water and water-water interactions, respectively, as sums of two-body terms. Since we are concerned with the absorption spectrum at the Franck-Condon state, we also calculated the interaction energies between formaldehyde in the excited state and water for the sample of R determined above. All calculations were executed using the CARLOS program in the CCP5 package.42 The coordinates of formaldehyde and water molecules and the interaction energies were saved on a disk to be used both for analysis (e.g., energy and radial distribution functions) and for further simulation of the absorption spectrum. With the sample Ri ( i = 1, ...,N) coordinates generated by the Monte Carlo method based on the Metropolis algorithm, the band intensity, I( T,u), at a given wavenumber u and temperature T can be calculated using the following equation:
where all of the values in the equations are in atomic units. The normal coordinates were calculated by diagonalizing the massweighted force constant matrix obtained by the ab initio method with the above-mentioned basis set. For the second term in eq 12, we used the quartic polynomial obtained by the fitting of interaction energies between the water molecules at R and formaldehyde distorted along the normal coordinate Qk:
Here A E g P ~ J R i )and Mgp~y(Ri)are the transition energy and dipole moment, respectively, from thepth vibrational level in the ground state to the qth level in the excited state at the solvent coordinate R’, which was evaluated by the numerical integration of the following equation:
c:’(Q3) = 0.0426Q32
