Computer Simulation of Chloroform with a Polarizable Potential Model

Apr 24, 1997 - Computer Simulation Investigation of the Water−Benzene Interface in a Broad Range of Thermodynamic States from Ambient to Supercritic...
0 downloads 12 Views 374KB Size
Download PDF
J. Phys. Chem. B 1997, 101, 3413-3419

3413

Computer Simulation of Chloroform with a Polarizable Potential Model Tsun-Mei Chang and Liem X. Dang* EnVironmental Molecular Sciences Laboratory, Pacific Northwest National Laboratory, Richland, Washington 99352

Kirk A. Peterson Department of Chemistry, Washington State UniVersity, Richland, Washington 99352 ReceiVed: NoVember 18, 1996; In Final Form: February 26, 1997X

A series of molecular dynamics simulations were carried out to examine the thermodynamic and the structural properties of the liquid and liquid/vapor of chloroform. A polarizable potential model was used to describe the intermolecular chloroform-chloroform interactions. The computed liquid densities and the enthalpies of vaporization and their corresponding temperature dependence are in excellent agreement with experimental values. The radial distribution functions and the corresponding neutron scattering cross sections were critically evaluated against the experimental data. The agreement between both approaches was found to be quite reasonable. The equilibrium interfacial properties of the chloroform liquid/vapor were also evaluated. The computed density profile shows that the interface is not sharp at a microscopic level and has a thickness of 5.9 Å at 298 K. The calculated surface tensions as a function of temperature are in good agreement with the corresponding experimental data. Accurate ab initio calculations were also carried out on the chloroform dimer, and the result for the binding energy was in good agreement with the value obtained from the molecular dynamics model.

I. Introduction Past practices of the Department of Energy (DOE) have resulted in the discharge of chemical and radioactive material into the environment at many DOE sites, which has resulted in extensive contamination of the soils and ground water at these sites. Chlorinated hydrocarbons such as carbon tetrachloride (CCl4) and chloroform (CHCl3) have contributed significantly to the problem of nuclear waste remediation due to their use in the processing of radioactive materials. This is particularly true at the nearby DOE Hanford site, where 106 kg of CCl4 was injected into the soil over 40 years of plutonium production.1 Its interaction with oxides and ground water is considered a major problem in environmental restoration. Current data suggest that 21% of the total release of CCl4 has been lost to the atmosphere, 12% resides in the vadose zone, and 2% is dissolved in ground water. Some 65% is unaccounted for at this time.2 Many aspects of the environmental contamination from CCl4 are similar to that for other chlorocarbons in the environment (e.g., CHCl3) and are relevant to industrial contamination as well as to restoration of DOE sites. Computer simulation studies have contributed a great deal to progress in the understanding of chemical and physical processes. With the advent of powerful computers in the 1980s and 1990s, a wide range of applications connected to real problems became tractable with computer simulation methods such as Monte Carlo (MC) or molecular dynamics (MD).3 The field has been expanding since then, and applications have emerged in almost every field of biological, biophysical, chemistry, the engineering sciences, and many other disciplines. As part of our ongoing research on the development of a fundamental understanding of thermodynamic and structural properties of chlorinated hydrocarbon liquids, as well as the liquid/liquid interfaces of water-chlorinated hydrocarbons, we * Corresponding author. X Abstract published in AdVance ACS Abstracts, April 1, 1997.

S1089-5647(96)03855-2 CCC: $14.00

examine here the thermodynamic and structural properties of chloroform in a widely varying environment, using MD techniques with a polarizable potential model. In view of the importance of CHCl3, there are several welltested chloroform potentials available in the literature.4-8 By applying the reference interaction site model theory, Hsu and Chandler proposed a five-site rigid chloroform potential that reasonably reproduced the static liquid structure of CHCl3 obtained from neutron scattering experiments.4 Using an effective pair potential, Evans5 analyzed the dynamical properties of chloroform and found that the rotation of chloroform is correlated with its center-of-mass translational motion. Several years later, Dietz and Heinzinger,6 as well as Bo¨hm and Ahlrichs,7 derived new all-atom potentials that yielded reasonable thermodynamics and static structures of liquid chloroform. More recently, Jorgensen and co-workers also developed a new set of CHCl3 potential parameters to study the relative partition coefficient for organic solutes.8 All of these studies have used pairwise, additive potentials, although in some cases the induced polarization effect has been incorporated implicitly in the potential parametrization. In consideration of the significance of CHCl3, along with the fact that CHCl3 has a high molecular polarizability of 8.53 Å3, it is essential to develop a polarizable potential model to describe chloroform interactions. In the present study, classical molecular dynamics are carried out to construct a polarizable CHCl3 model potential. The parametrization procedure involves fitting the computed properties to the available experimental measurements. The optimized parameters are obtained when the potential accurately reproduces experimental results such as the density and the heat of vaporization. The potential is then used in the investigation of the CHCl3 liquid/vapor interface. Further tests of this potential are made by comparing the calculated and experimental interfacial properties to assure the validity of this potential. The paper is organized as follows. In section II, we briefly describe the polarizable model potential used in this study and the © 1997 American Chemical Society

3414 J. Phys. Chem. B, Vol. 101, No. 17, 1997

Chang et al.

computational details of the MD simulations. Properties of the bulk liquid CHCl3 are presented and discussed in section III. Comparisons of these computed properties to those of experimental and previous theoretical work are also made. The simulation results of the CHCl3 liquid/vapor interface are summarized in section IV. The conclusions are presented in section V. II. Potential Model and Computational Details A. Potential Model. In the present study, a rigid five-site all-atom interaction potential is developed for chloroform. The CHCl3 molecule is kept at the experimental geometry with the C-Cl and C-H bond lengths of 1.76 and 1.07 Å and Cl-CCl and Cl-C-H bond angles of 111.2° and 107.6°, respectively.9 There are fixed atomic charges and Lennard-Jones potential parameters assigned to each atom to describe the intermolecular van der Waals and Coulombic interactions. In addition, there are point polarizabilities associated with each atom to account for the nonadditive induced polarization effect. The total interaction energy can be written as

Utot ) Upair + Upol

(1)

where the pairwise additive part of the potential is the sum of the Lennard-Jones and Coulomb interactions,

Upair )

( [( ) ( ) ] )

∑i ∑j 4

σij

σij

12

-

rij

qiqj

6

+

rij

rij

(2)

Here, rij is the distance between atom site i and j, q is the atom charge, and σ and  are the Lennard-Jones parameters. The nonadditive polarization energy is given by

∑i µi‚Ei0

Upol ) - 1/2

(3)

where Ei0 is the electric field at site i produced by the fixed charges in the system,

Ei0 )

∑ j*i

qjrij (4)

rij3

µi is the induced dipole moment at atom site i and is defined as

µi ) REi

(5)

where

Ei ) Ei0 +

Tij‚µj ∑ j*i

(6)

In the above, Ei is the total electric field at atom i, R is the point polarizability, and Tij is the dipole tensor

Tij )

(

)

1 3rijrij -1 rij3 rij2

(7)

During the molecular dynamics simulations, an iterative procedure is used to solve eqs 5 and 6. The iterative procedure has been discussed in details in a previous study.10 Convergence is achieved when the deviation of the dipole moment from two sequential iterations falls to within 0.000 01 D. B. Simulation Details. Classical molecular dynamics simulations and ab initio electronic structure calculations are performed to construct the nonadditive CHCl3 interaction potential. The partial charges associated with the carbon,

TABLE 1: Potential Parameters for CHCl3 Used in the MD Simulationsa atom type

σ (Å)

 (kcal/mol)

q (e)

R (Å3)

C Cl H

3.41 3.45 2.81

0.137 0.275 0.020

0.5609 -0.1686 -0.0551

0.878b 1.910b 0.135b

a σ and  are the Lennard-Jones parameters, q is the atomic charge, and R is the atomic polarizability. b Reference 12.

hydrogen, and chlorine atoms were taken from ab initio calculations using the Gaussian 94 program.11 Specifically, the atomic charges are fit to the electrostatic potentials obtained from the MP2 CHCl3 wave functions optimized at a HartreeFock level with the 6-31+G* basis set; then the charges were scaled to reproduce the experimental monomer dipole moment of chloroform. The atomic site polarizabilities for carbon, chlorine, and hydrogen were taken from the work of Applequist et al.,12 which well reproduce the experimental CHCl3 molecular polarizability. By assuming that the cross-interactions between unlike atoms obey the Lorentz-Berthelot combining rules,3 this leads to only six Lennard-Jones parameters (σ and  for each atom) that remain to be determined. In order to derive these Lennard-Jones parameters, a series of molecular dynamics simulations were performed for liquid chloroform at 298 K. The structural and thermodynamic properties of CHCl3 including radial distribution functions, densities, and enthalpies of vaporization were evaluated from these simulations and compared to corresponding experimental measurements. The potential parameters were adjusted to reproduce the experimental properties of liquid CHCl3. In Table 1, the optimized potential parameters are summarized. The MD simulations were performed on a system consisting of 265 CHCl3 molecules in a cubic simulation cell with a linear dimension of roughly of 33 Å. The simulations were carried out in an isobaric-isothermal (NPT) ensemble at 1 atm with a time step for pressure relaxation of 0.2 ps. By coupling to a external thermal bath with a time constant of 0.1 ps, calculations of the liquid properties chloroform were carried out at several temperatures, T ) 283, 298, and 313 K. From a MaxwellBoltzmann distribution corresponding to the desired simulation temperature, the initial velocities of each atom were randomly assigned. During the dynamic simulations, the SHAKE13 algorithm was employed to constrain all the bond lengths to fix the CHCl3 internal geometry, and the periodic boundary conditions were applied in all three directions. A time step of 2 fs was used in integrating the equations of motion. In this study, the nonbonded (Lennard-Jones, Coulombic, polarization) interactions were truncated at a molecular separation of 15 Å. At all temperatures, the systems is equilibrated for 100 ps, followed by at least 200 ps of data collection for analysis. To further check the validity of the polarizable chloroform potential, additional simulations were carried out to study the CHCl3 liquid/vapor interface. The initial configuration of the liquid/vapor interface is obtained by replacing a box of 265 CHCl3 molecules that had been equilibrated previously in the bulk simulations into the middle of a simulation cell whose linear dimensions are 32.9 × 32.9 × 82.9 Å (see Figure 1). This creates two CHCl3 liquid/vapor interfaces that are chosen to be perpendicular to the z axis. The whole system is then equilibrated for 200 ps in a NVT ensemble at 1 atm and several temperaturess268, 283, and 298 Kswith periodic boundary conditions applied in all three directions. This is followed by a 200 ps simulation for data analysis. The interactions are truncated at a molecular center-of-mass separation of 15 Å.

Computer Simulation of Chloroform

J. Phys. Chem. B, Vol. 101, No. 17, 1997 3415

Figure 1. Schematic representation of the simulation box for the CHCl3 liquid/vapor interface. The liquid phase of CHCl3 is in the middle of the box. The z axis is chosen to be perpendicular to the interface.

Figure 3. Calculated (open circles) and experimental (filled circles) enthalpies of vaporization as a function of temperature.

Figure 2. Calculated (open circles) and the experimental (filled circles) liquid densities as a function of temperature.

III. Liquid CHCl3 Properties A. Temperature Dependence of Thermodynamic Properties. The thermodynamic results from the simulations are shown in Figures 2 and 3. The computed liquid densities as a function of temperature are depicted in Figure 2, along with the corresponding experimental liquid densities.14 It is clear that the agreement is quite good for the entire temperature range. The enthalpy of vaporization, ∆Hvap(T), can be evaluated from the MD simulations, according to the following equations

∆Hvap(T) ) -Epot(T) + RT

(8)

where Epot(T) is the total intermolecular potential energy of liquid CHCl3 at temperature T. The computed and experimental enthalpies of vaporization of chloroform are depicted in Figure 3. The predicted temperature dependence of the heats of vaporization are in good agreement with the experimental measurements with the largest deviation being 0.6% for the entire temperature range.14 The largest error in computed liquid densities and enthalpies of vaporization is 0.6% in comparison with experiments, which is an indication of the quality of the present chloroform model potential in describing the thermodynamic properties of liquid CHCl3 for a temperature range of 283-313 K. By examining the various components of the interaction energy, we find that the polarization energy constitutes roughly 1.5% of the total energy, which is on the same order as the Coulombic interaction, but small compared to the Lennard-Jones interactions. This leads us to expect that the pairwise additive

Figure 4. Gas-phase interaction potential for CHCl3 dimer as a function of C-C distance. Also shown is the minimum-energy structure of the CHCl3 dimer.

potential may be able to reproduce many properties of bulk liquid CHCl3. The average total dipole moment of a chloroform molecule in bulk liquid is 1.14 D, which is 13% larger than the gas phase dipole moment of 1.01 D. It has been established that, under certain circumstances such as liquid/liquid interfaces, it is essential to include the many-body polarization effect in the interaction potentials to properly describe the molecular properties. Because CHCl3 is a highly polarizable liquid, we believe that a polarizable potential model will give a more accurate and consistent description of chloroform molecules in complex heterogeneous environments. Shown in Figure 4 is the optimal geometry and gas phase interaction energy as a function of the C-C separation for the CHCl3 dimer. CHCl3 molecules are found to orient in such a way as to maximize the electrostatic interactions. The minimum interaction energy for the dimer is about -2.66 kcal/mol with a C-C separation of 4.18 Å. Accurate ab initio electronic structure calculations were also carried out on the CHCl3 dimer. The standard augmented correlation consistent polarized double-ζ (aug-cc-pVDZ) basis sets of Dunning and co-workers were used throughout.15-17 The aug-cc-pVDZ sets include diffuse s, p, and d functions, which have been shown to be important for accurately describing intermolecular interactions and mo-

3416 J. Phys. Chem. B, Vol. 101, No. 17, 1997

Chang et al.

Figure 6. C-C atomic radial distribution functions as a function of temperature from bulk CHCl3 simulations for T ) 283 (solid line), 298 (dashed line), and 313 K (dotted line).

Figure 5. Atomic radial distribution functions from bulk CHCl3 simulations at 298 K: (a, top) C-C (solid line), C-Cl (dashed line), and C-H (dotted line) radial distribution functions; (b, bottom) ClCl (solid line), Cl-H (dashed line), and H-H (dotted line) radial distribution functions.

lecular properties such as dipole moments and polarizabilities. The effects of electron correlation were computed at the MP2 level of theory. In addition, the standard function counterpoise (CP) method18 was used to correct for basis set superposition errors (BSSE). For the correction of the calculated equilibrium binding energies for zero-point vibrational effects, vibrational frequencies were obtained at the MP2/6-32+G** level. The calculations were carried out with both the Gaussian94 and Molpro programs.11,19 Two stable isomers of the chloroform dimer have been calculated. One of these (denoted isomer 1) consists of an endto-end type structure with C3V symmetry where the C-H end of one monomer is directed toward the C-Cl3 face to the other. The other isomer (denoted isomer 2), which is nearly isoenergetic, has an antiparallel C2h symmetry structure. At the MP2/ aug-cc-pVDZ level of theory, the counterpoise-corrected equilibrium binding energy of isomer 1 was computed to be -2.44 kcal/mol. Likewise, the analogous value for isomer 2 was calculated to be -2.41 kcal/mol. Addition of harmonic zeropoint vibrational corrections yield 0 K binding enthalpies of -2.28 and -2.06 kcal/mol for isomers 1 and 2, respectively. In the case of isomer 1 (end-to-end isomer), the MP2 C-C equilibrium distance was calculated to be 3.76 Å, while that of isomer 2 (antiparallel isomer) was 3.45 Å. It is not clear why the C-C distance obtained from the MD model is significantly larger than the corresponding C-C distance obtained from the electronic structure calculations. This is probably due to the

deficiency in the MD model or the level of theory used in the electronic structure calculations was insufficient to describe the chloroform dimer interactions. Full details of the dimer electronic structure calculations will appear elsewhere.20 B. Structural Properties of Liquid CHCl3. The structure of liquid chloroform can be analyzed in terms of the atomic pair correlation functions and angular distribution functions. The six atomic radial distribution functions, gab(r), calculated from simulations at 298 K are shown in Figure 5. Well-defined peaks are observed for all radial distribution functions, suggesting that there is a structural order in liquid chloroform. In Figure 5a, we show the computed C-C, C-Cl, and C-H radial distribution functions. The C-C radial distribution function suggests that there is a long-range structural correlation between CHCl3 molecules beyond 15 Å, as indicated by the small oscillation at large C-C distance. The maximum of the first peak in gCC(r) is at 5.35 Å and the first minimum is at 7.6 Å, which are consistent with previous studies.5-7 By integrating over the first peak of gCC(r), it is found that roughly 12 neighboring carbon atoms reside in the first solvation shell of a central carbon atom. This result reflects that liquid CHCl3 resembles a closed-packed, monatomic liquid, similar to liquid CCl4. The first maxima in gCCl(r) and gCH(r) are located at 5.0 and 4.65 Å, which are slightly larger than those reported from earlier simulations.5-7 Again, weak oscillations are observed in these radial distribution functions at larger separations, implying that liquid CHCl3 is well structured. The calculated Cl-Cl, Cl-H, and H-H radial distribution functions are presented in Figure 5b. While gHH(r) exhibits weak oscillations at large H-H distance, the correlation of gClH(r) and gClCl(r) are much shorter ranged. This is contributed to more chlorine atoms being present in these simulations, and the averaging smoothes out the structural correlation. The first maxima of gHH(r), gClH(r), and gClCl(r) are positioned at 5.15, 3.25, and 3.70 Å, respectively. These locations are consistent with those reported in the literature with different CHCl3 interaction potentials.5-7 In Figure 6, the gCC(r) as a function of temperature is depicted. As expected, the peak positions of the radial distribution function shift to larger separations at higher temperature as a result of decreasing liquid density. The features of gCC(r) remain the same for the entire temperature range, indicating similar structural correlation for these liquids. The somewhat broadened features of the radial distribution functions are correlated with the larger thermal fluctuations at higher temperatures.

Computer Simulation of Chloroform

J. Phys. Chem. B, Vol. 101, No. 17, 1997 3417

Figure 8. Computed density profile (filled circles) from the CHCl3 liquid/vapor interface simulation at 298 K. The solid line is the best fit to eq 11. Figure 7. Cross section for liquid CHCl3. The solid line is from the MD simulation at 298 K, and the experimental results (circles) were obtained by Chieux et al. at 298 K.

In order to compare directly to the structural information from neutron scattering experiments,21,22 the atomic radial distribution functions are Fourier transformed into partial structure factors, Sab(k), as

∫0∞|gab(r) - 1|sink kr r dr

Sab(k) ) 4πF

(9)

where k is the scattering wave vector and r is the liquid density. The cross section measured in the neutron scattering experiment is simply a weighted sum of these structure factors

(dσ/dΩ)int(k) ) 0.33089SCC(k) + 2.86013SCCl(k) + 6.18045SClCl(k) + 0.10468SHH(k) - 1.60855SClH(k) 0.37129SCH(k) (10) We show in Figure 7 the scattering cross sections for liquid CHCl3 at 298 K calculated from the simulation and the experimental curves obtained by Chieux and co-workers at 298 K.22 The agreement between the simulations and experiments is satisfactory in general although small differences are observed. The discrepancy between experiment and our calculations at small k is most likely due to the finite size of the simulations, since it is well-known that long-range oscillations of g(r) dominate the behavior of S(k) at small k. IV. CHCl3 Liquid/Vapor Interface A. Density Profile. By computing the liquid densities in liquid slabs of 1 Å thickness parallel to the liquid/vapor interface, we obtain the CHCl3 density profile as a function of the z coordinate normal to the interface, which is presented in Figure 8. The sharp decrease in the liquid densities around z ) 25 and 58 Å indicates that there are two well-defined interfaces. Between these two interfaces exists a bulk CHCl3 region 33 Å wide. A hyperbolic tangent functional form is often used to describe the shape of the density profile23,24

[ ]

z - z0 1 1 F(z) ) (FL + FV) - (FL - FV) tanh 2 2 d

(11)

where FL and FV are the liquid and vapor densities, respectively. Here, d gives an estimate of the interfacial thickness, and z0 denotes the position of the Gibbs dividing surface. The best

fit of the density profile gives a bulk liquid of 1.46 g/cm3, which is in good agreement with the experimental density of 1.47 g/cm3 at 298 K.14 The fit also yields an estimate of the so-called “1090” thickness of 5.9 Å. In Figure 9, snapshots of the CHCl3 liquid/vapor interface from the top and the side are shown. Obviously, the interface is very irregular and not sharp at a microscopic scale. B. Surface Tension. The present chloroform potential can be further inspected by evaluating the surface tension, γ, of the CHCl3 liquid/vapor interface. This will provide a consistent check on the validity of this model in describing chloroform interactions in heterogeneous environments. From molecular dynamics, the surface tension is simply the difference between the pressure components in the direction parallel and perpendicular to the interface25

γ)

(

)

1 pxx + pyy - pzz Lz 2 2

(12)

where pRR (R ) x, y, or z) is the RR element of the pressure tensor and Lz is the linear dimension of the simulation cell in the interface normal direction. The Rβ element of the pressure tensor can be calculated according to the virial equation as

pRβ )

1 V

(∑ N

i)1

miViRViβ +

1 N′

N′

)

∑ ∑Fi′j′Rrijβ

2i′)1j′)1

(13)

where N and N′ are the numbers of molecules and atoms, respectively. V is the total volume of the system, mi is the mass of molecule i, and ViR is the center-of-mass velocity of molecule i. In the above expression, Fi′j′R is the R component of the force exerted on atom i′ of molecule i due to atom j′ of molecule j, and rijβ is the β component of the vector connecting the center of mass of molecules i and j. In Figure 10, we show the calculated surface tension as a function of temperature, along with the experimental results.14 Clearly, as the temperature increases, the surface tension decreases. The good agreement between the computed and experimental surface tensions ensures us that the present polarizable chloroform potential indeed reasonably describes the intermolecular interactions of chloroform molecules at the liquid/ vapor interface. A summary of the calculated thermodynamic properties of CHCl3 is presented in Table 2. C. Structure at the Interface. It is of great importance to understand the structural changes of chloroform molecules at

3418 J. Phys. Chem. B, Vol. 101, No. 17, 1997

Chang et al.

Figure 10. Calculated (open circles) and experimental (filled circles) surface tension as a function of temperature.

Figure 11. C-C atomic radial distribution functions as a function of the z coordinate of the liquid slabs from the CHCl3 liquid/vapor simulations. The solid curves from the top to the bottom correspond to the liquid slab centered at z ) 34, 32, 30, 28, 26, and 24 Å.

TABLE 2: Thermodynamic Properties for Polarizable CHCl3 at 298 K and 1 atm g/cm3

Figure 9. Snapshots of the instantaneous structure of the CHCl3 liquid/ vapor interface. The light color corresponds to chloroform molecules at the surface, and the dark color corresponds to chloroform molecules in the bulk. The upper and lower panels show typical structures of the side and top of the interface, respectively.

the liquid/vapor interface. To characterize this, it is useful to define the atomic radial distribution functions, gab(r,z), as a function of the z coordinate of the CHCl3 molecule. We first divide the simulation box into liquid slabs of 2 Å thickness along the z direction normal to the interface. The number of neighboring atoms of type b that are in a spherical shell of thickness dr that is distance r away from a central atom a in a particular liquid slab is then recorded, regardless of the z position of the atoms b. The z-dependent radial distribution functions are simply defined as the number of neighbors divided by the total number of atoms of type a in the slab and the average number of atoms b inside the spherical cell in the bulk liquid. In other words, we use the bulk radial distribution functions as the references to the z-dependent radial distribution functions. In Figure 11, we show the computed C-C atomic radial

F, ∆Hvap, kcal/mol surface tension, dyn/cm diffusion, 10-5 cm2/s dimer energy, kcal/mol C-C dimer distance, Å a

this work

expta

1.473 ( 0.005 7.47 ( 0.01 26.9 ( 0.3 3.5 ( 0.2 -2.66 4.18

1.473 7.47 26.6 2.8

Reference 14.

distribution functions as a function of the z coordinates of the chloroform molecules. Here, small z indicates the interfacial region and large z indicates the bulk region. No significant changes are observed in the peak positions of these radial distribution functions, suggesting that the local structures of CHCl3 remain unchanged from the bulk to the interface. The decrease in the peak height of the radial distribution functions as the CHCl3 approaches the surface simply reflects the departure of CHCl3 molecules into the vapor side of the interface. In principle, the magnitude of the molecular radial distribution functions at the interface should be half that of the bulk radial distribution functions for a planar interface. This is demonstrated by the gab(r,z) centered around the Gibbs dividing surface at z ) 26 Å. It has been established that although the local structural arrangement of molecules remain unchanged at the liquid/vapor

Computer Simulation of Chloroform

J. Phys. Chem. B, Vol. 101, No. 17, 1997 3419 a weak preferred orientational order is observed for CHCl3 molecules at the interface which diminished in the bulk liquid. It is also observed that the calculated dipole moments of chloroform are close to their gas-phase value near the interface, while CHCl3 molecules have dipole moments corresponding to their bulk value far from the interface.

Figure 12. Probability distributions for the angle between the z axis and the intramolecular C-H bond. The data points correspond to the liquid slab centered at z ) 24 (squares), 28 (circles), and 32 Å (crosses). A straight line describes a uniform angular distribution.

interface, the interface may induce a certain orientational order. This situation is especially true for molecules with directional interactions such as the case of the H2O liquid/vapor interface.26,27 For spherical molecules such as CCl4, there is little orientational order observed.21 For chloroform molecules, we expect a weak orientational order because of the competition of the high molecular symmetry and the weak directional interaction. To verify this, we examine the normalized probability distribution for angles between the intramolecular C-H bond and the interfacial normal direction as a function of the z coordinates of the CHCl3 molecules. These results are shown in Figure 12. Here, a straight horizontal line indicates a uniform angular distribution. Clearly, while the angular distribution function deviates more from the straight line in the interfacial region (small z coordinate), the angular distribution functions are approaching the uniform distribution in the bulk region (larger z value). These results suggest that CHCl3 has a weak preferred orientational order at the interface as expected. V. Conclusion With the use of statistical mechanical molecular dynamics techniques, we have developed a new interaction model potential for chloroform that incorporates many-body polarization effects. From a series of molecular dynamics simulations, the computed liquid densities and enthalpies of vaporization are found to be in excellent agreement with experimental results for temperatures ranging from 283 to 303 K. The static structure of liquid chloroform is also found to be consistent with the experimental measurements. The equilibrium interfacial properties of the chloroform liquid/ vapor interface is also examined using this potential model. The computed density profile shows that the interface is not smooth at a microscopic level and has a thickness of 5.9 Å at 298 K. The calculated surface tensions of the chloroform interface agree well with the corresponding experimental data as a function of temperature. This suggests that the present chloroform potential gives a good description of the CHCl3 interactions at the interface. We observed that while the local structure of chloroform remains unchanged from the interface into the bulk,

Acknowledgment. This work was performed under the auspices of the Division of Chemical Sciences, Office of Basic Energy Sciences, US Department of Energy, under Contract DEAC06-76RLO 1830 with Battelle Memorial Institute, which operates the Pacific Northwest Laboratory, a multiprogram national laboratory. Computer resources were provided by the Division of Chemical Sciences and by the Scientific Computing Staff, Office of Energy Research, at the National Energy Research Supercomputer Center (Berkeley, CA). References and Notes (1) Haag, W. R.; Yao, C. C. D. EnViron. Sci. Technol. 1992, 26, 1005. (2) U.S. Department of Energy, Office of Environmental Management, FY 1995 Technology Development Needs Summary, 1994; p 2-17. (3) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford University Press: Oxford, 1987. (4) Hsu, C. S.; Chandler, D. Mol. Phys. 1979, 37, 299. (5) Evans, M. W. J. Mol. Liq. 1983, 25, 211. (6) Dietz, W.; Heinzinger, K. Ber. Bunsen-Ges. Phys. Chem. 1985, 88, 543. (7) Bo¨hm, H. J.; Ahlrichs, R. Mol. Phys. 1985, 54, 1261. (8) Jorgensen, W. L.; Briggs, J. M.; Contreras, M. L. J. Phys. Chem. 1990, 94, 1683. (9) Bartell, L. S.; Brockway, L. O.; Schwendeman, R. H. J. Chem. Phys. 1955, 23, 1854. (10) Chang, T.-M.; Peterson, K. A.; Dang, L. X. J. Chem. Phys. 1995, 103, 7502. (11) Gaussian 94 (Revision A.1): Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Gill, P. M. W.; Johnson, B. G.; Robb, M. A.; Cheeseman, J. R.; Keith, T. A.; Petersson, G. A.; Montgomery, J. A.; Raghavachari, K.; AlLaham, M. A.; Zakrzewski, V. G.; Ortiz, J. V.; Foresman, J. B.; Cioslowski, J.; Stefanov, B. B.; Nanayakkara, A.; Challacombe, M.; Peng, C. Y.; Ayala, P. Y.; Chen, W.; Wong, M. W.; Andres, J. L.; Replogle, E. S.; Gomperts, R.; Martin, R. L.; Fox, D. J.; Binkley, J. S.; Defrees, D. J.; Baker, J.; Stewart, J. P.; Head-Gordon, M.; Gonzalez, C.; Pople, J. A. Gaussian Inc., Pittsburgh, PA, 1995. (12) Applequist, J.; Carl, J. R.; Fung, K.-K. J. Am. Chem. Soc. 1972, 94, 2952. (13) Ryckaert, J.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327. (14) Yaws, C. L. Thermodynamic and Physical Property Data; Golf Publishing: Houston, TX, 1992. (15) Dunning, T. H., Jr. J. Chem. Phys. 1989, 90, 1007. (16) Kendall, R. A.; Dunning, T. H., Jr.; Harrison, R. J. J. Chem. Phys. 1992, 96, 6796. (17) Woon, D. E.; Dunning, T. H., Jr. J. Chem. Phys. 1993, 98, 1358. (18) Boys, S. F.; Bernardi, F. Mol. Phys. 1970, 19, 553. (19) MOLPRO is a package of ab initio programs written by H.-J. Werner and P. J. Knowles with contributions from J. Almlo¨f, R. D. Amos, M. J. O. Deegan, S. T. Elbert, C. Hampel, W. Meyer, K. A. Peterson, R. M. Pitzer, A. J. Stone, P. R. Taylor, and R. Lindh. (20) Peterson, K. A. Manuscript in preparation. (21) Bertagnolli, H.; Leight, D. O.; Zeidler, M. D.; Chieux, P. Mol. Phys. 1978, 36, 1769. (22) Bertagnolli, H.; Chieux, P. Mol. Phys. 1984, 51, 617. (23) Alejandre, J. D.; Tildesley, J.; Chapela, G. A. J. Chem. Phys. 1995, 102, 4574. (24) Lie, G. C.; Grigoras, S.; Dang, L. X.; Yang, D.-Y.; Mclean, A. D. J. Chem. Phys. 1995, 102, 4574. (25) Kirkwood, J. G.; Buff, F. P. J. Chem. Phys. 1949, 17, 338. Salomons, E.; Mareschal, M. J. Phys.: Condens. Matter 1991, 3, 3645. (26) Toxvaerd, S. Faraday Symp. Chem. Soc. 1981, 16, 159. (27) Croxton, C. A. Physica 1981, 106A, 239.