Tanya Raschke

Hydration of Benzene and Cyclohexane

Raschke, T.M, Levitt, M. Detailed hydration maps of benzene and cyclohexane reveal distinct water structures. J. Phys. Chem. (submitted).

To investigate the origins of the hydrophobic effect, I have computed molecular dynamics simulations of two, small solutes, benzene and cyclohexane. While cyclohexane is considered a generic hydrophobic molecule, the aromatic nature of benzene increases its solubility in water. Both of these molecules are similar in size and shape, and the effects of shape can be distinguished from electrostatic effects by comparing the results from these two molecules.

Figure 1: Water Density Surrounding Benzene

Detailed 3-D Hydration Maps of Cyclohexane and Benzene

The density of the water O and H atoms surround the solutes was mapped in fine detail using a gridding technique. At each output time step, the system was rotated and translated to place the solute centroid at the center of the reference axes. The location of water O and H atom centers within a 20 Å cube surrounding the solute were placed in separate histograms, and the data normalized to the background density of a pure water simulation. This technique has the advantage that it can be applied to irregularly shaped solutes.

Several solute and solvent models were tested, including the ENCAD potential for benzene and cyclohexane along with the F3C flexible water model. In the GROMACS simulations the OPLS force field was used for the solutes, and the SPC, SPC/E, TIP3P and TIP4P water models were used. Surprisingly, all of the models give very similar results.

Figure2: Water Density Surrounding Cyclohexane

Water molecules form hydrogen bonds with the faces of the benzene ring.

Figure 1 shows the water O atom density surrounding benzene computed from the ENCAD simulations. Figure 2 shows the same data for cyclohexane. The distribution of water molecules around benzene and cyclohexane are clearly different. Both the O and H density maps show intense peaks centered above and below the faces of the benzene ring. The density along the edge of the benzene ring, however is very similar in intensity and volume to that surrounding cyclohexane.

The water O density maps reflect the average position of the water molecules surrounding the solutes, while the H density maps contain orientation information. The water molecules located in the peak regions above and below the benzene ring are in a hydrogen-bonding orientation, with the global H density peak located between the benzene center and the O density peak.

Different Solute/Solvent Models Give Remarkably Similar Distributions

Figure 3 shows more traditional, radial distribution functions computed from the solute centroid to either the O or H centers. The radial distribution functions have been averaged over the two volume regions: a cone lying 20° from the z-axis, termed the axial region, and the remaining volume, termed the equatorial region. All of the results from the different solute/solvent models are remarkably similar. The largest differences are seen in the ENCAD simulations, which use a fully-flexible water model and a short cutoff of 6 Å for the nonbonded interactions.

Figure3: Radial Distribution Functions For All Simulation Models

Panels E and F show a potential of mean force computed from gSol-O. PMF represents the energy of a water molecule at varying distance from the solute centroid, relative to that at infinite separation (i.e. in the bulk). We estimate the binding energy of water to the faces of the benzene ring at about -1 kcal/mol.

Aggregation of Simple Hydrophobic Solutes

Raschke, T.M., Tsai, J., Levitt M. Quantification of the hydrophobic interaction by simulations of the aggregation of small hydrophobic solutes in water. PNAS U.S.A., (2001) 98:5965-5969. [PDF]

The accuracy of our MD simulations for use in modeling the hydrophobic effect was assessed by studying the aggegation of small, hydrophobic solutes in water. The solutes studied were Methane, Butane, Isobutylene and Benzene. Many simulations were calculated with different numbers of solutes in the box (i.e. varying concentration) and at different box dimensions. All of the solutes formed fluctuating aggregates in the simulations. Figure 1A shows a snapshot from a simulation with 10 methanes grouped into 3 clusters. Figure 1B shows a snapshot from a simulation with 20 butanes, repeated 4 times to illustrate the peroidic boundary conditions. At very high concentrations, the solutes exhibited phase separation, forming cylindrical structures that spanned the boundaries of the simulation box.

Figure 1: Snapshots From MD Simulations Showing Clusters of Methane (A) and Butane (B)

The free energy of formation of the clusters was

Figure 2: Hydrophobic Energy Relates to Contact Surface