Monte Carlo simulations of the hydrophobic effect in aqueous electrolyte solutions

The hydrophobic interaction between two methane molecules in salt-free and high salt-containing aqueous solutions and the structure in such solutions have been investigated using an atomistic model solved by Monte Carlo simulations. Monovalent salt representing NaCl and divalent salt with the same nonelectrostatic properties as the monovalent salt have been used to examine the influence of the valence of the salt species. In salt-free solution the effective interaction between the two methane molecules displayed a global minimum at close contact of the two methane molecules and a solvent-separated secondary minimum. In 3 and 5 M monovalent salt solution the potential of mean force became slightly more attractive, and in a 3 M divalent salt solution the attraction became considerably stronger. The structure of the aqueous solutions was determined by radial distribution functions and angular probability functions. The distortion of the native water structure increased with ion valence. The increase of the hydrophobic attraction was associated with (i) a breakdown of the tetrahedral structure formed by neighboring water molecules and of the hydrogen bonds between them and (i) the concomitant increase of the solution density.     

Monte Carlo structure simulations for aqueous 1,4-dioxane solutions

Monte Carlo simulations in the NpT ensembles have been performed for the structure exploration of aqueous 1,4-dioxane solutions. Three different systems with all-atom dioxane:TIP4P water molar compositions of 2:500 (code:D2), 8:465 (D8), and 17:425 (D17) modeled solutions of 0.22, 0.88, and 1.86 mol/dm3 concentrations, respectively, at T = 298 K and p = 1 atm. The calculated solution densities increase from 0.992 to 1.002 g/cm3 with increasing dioxane concentration and approach the experimentally determined densities within 1%. This close agreement was achieved by utilizing RESP charges fitted to the in-solution IEF-PCM/B3LYP/6-31G* electrostatic potential of dioxane taken in its chair conformation and recently developed C, H steric parameters for ethers for calculations with a 12-6-1 all-atom potential. Solution structure analyses pointed out that the dioxane molecules arrange in the solutions with favorable distances of 4-8 angstroms for the ring symmetry centers. Within this range not only pairs of rings but triangular triads and tetrads have also been observed with center-center distances <8 angstroms. For the D8 system, about 25% of the sampled configurations included such a triad. In the case of the D17 model, two simulations starting from different solution configuration predicted different degrees for the dioxane aggregation in aqueous solution. In the more aggregated structure 3-21 triads are consistently maintained and 1-2 tetrads are formed in 58% of the configurations. Each dioxane oxygen forms about one hydrogen bond, on average, to a water molecule in the 0.22-1.86 molar range. The most likely O(dioxane)...H(water) hydrogen bond distance is 1.75-1.80 angstroms compared to the optimal distance of 1.72 angstroms in the isolated dimer. The optimal dioxane-water interaction energy of -5.65 kcal/mol indicates a remarkable hydrogen-bond acceptor character for dioxane.     

Monte Carlo study of the aqueous hydration of dimethylphosphate conformations

Monte Carlo computer simulations were performed on dilute aqueous solutions of the dimethylphosphate anion and the sodium dimethylphosphate ion pair, with the two phosphodiester torsional angles in the gauche–gauche, gauche–trans, and trans–trans conformations. The structural and energetic aspects of the aqueous hydration of each molecule were analyzed in terms of quasi component distribution functions based on the proximity criterion and partitioned into ionic, hydrophilic, and hydrophobic contributions to facilitate an understanding of the hydration pattern and conformational trends in these multifunctional solutes. Special attention was also paid to methodological issues affecting hydration, such as statistical uncertainty in the determined hydration indices, choice of partial atomic charges for the solute atoms, and solute–water interaction potentials adopted in the simulations. The results showed that gauche–trans and gauche–gauche forms are equally favorable for the dimethylphosphate anion with the trans extended form destabilized by hydration. The sodium dimethylphosphate ion pair hydration energetically favors the trans–trans conformation.     

Acceleration of Markov chain Monte Carlo simulations through sequential updating

Strict detailed balance is not necessary for Markov chain Monte Carlo simulations to converge to the correct equilibrium distribution. In this work, we propose a new algorithm which only satisfies the weaker balance condition, and it is shown analytically to have better mobility over the phase space than the Metropolis algorithm satisfying strict detailed balance. The new algorithm employs sequential updating and yields better sampling statistics than the Metropolis algorithm with random updating. We illustrate the efficiency of the new algorithm on the two-dimensional Ising model. The algorithm is shown to identify the correct equilibrium distribution and to converge faster than the Metropolis algorithm with strict detailed balance. The main advantages of the new algorithm are its simplicity and the feasibility of parallel implementation through domain decomposition.     

Acceleration of Monte Carlo simulations through spatial updating in the grand canonical ensemble

A new grand canonical Monte Carlo algorithm for continuum fluid models is proposed. The method is based on a generalization of sequential Monte Carlo algorithms for lattice gas systems. The elementary moves, particle insertions and removals, are constructed by analogy with those of a lattice gas. The updating is implemented by selecting points in space (spatial updating) either at random or in a definitive order (sequential). The type of move, insertion or removal, is deduced based on the local environment of the selected points. Results on two-dimensional square-well fluids indicate that the sequential version of the proposed algorithm converges faster than standard grand canonical algorithms for continuum fluids. Due to the nature of the updating, additional reduction of simulation time may be achieved by parallel implementation through domain decomposition.     

On the convergence of parallel tempering Monte Carlo simulations of LJ38

The convergence of parallel tempering Monte Carlo simulations of the 38-atom Lennard-Jones cluster starting from the Oh global minimum and from the C(5v) second-lowest-energy minimum is investigated. It is found that achieving convergence is appreciably more difficult, particularly at temperatures in the vicinity of the Oh --> C(5v) transformation when starting from the C(5v) structure. A strategy combining the Tsallis generalized ensemble and the parallel tempering algorithm is implemented and used to improve the convergence of the simulations in the vicinity of the Oh --> C(5v) transformation.     

Monte Carlo simulations of single crystals from polymer solutions

A novel "anisotropic aggregation" model is proposed to simulate nucleation and growth of polymer single crystals as functions of temperature and polymer concentration in dilute solutions. Prefolded chains in a dilute solution are assumed to aggregate at a seed nucleus with an anisotropic interaction by a reversible adsorption/desorption mechanism, with temperature, concentration, and seed size being the control variables. The Monte Carlo results of this model resolve the long-standing dilemma regarding the kinetic and thermal roughenings, by producing a rough-flat-rough transition in the crystal morphology with increasing temperature. It is found that the crystal growth rate varies nonlinearly with temperature and concentration without any marked transitions among any regimes of polymer crystallization kinetics. The induction time increases with decreasing the seed nucleus size, increasing temperature, or decreasing concentration. The apparent critical nucleus size is found to increase exponentially with increasing temperature or decreasing concentration, leading to a critical nucleus diagram composed in the temperature-concentration plane with three regions of different nucleation barriers: no growth, nucleation and growth, and spontaneous growth. Melting temperatures as functions of the crystal size, heating rate, and concentration are also reported. The present model, falling in the same category of small molecular crystallization with anisotropic interactions, captures most of the phenomenology of polymer crystallization in dilute solutions.     

Monte Carlo simulations of the hydration of substituted benzenes with OPLS potential functions

Intermolecular potential functions have been developed for use in computer simulations of substituted benzenes. Previously reported optimized potentials for liquid simulations (OPLS) for benzene and organic functional groups were merged and tested by computing free energies of hydration for toluene, p-xylene, phenol, anisole, benzonitrile, p-cresol, hydroquinone, and p-dicyanobenzene. The calculations featured Monte Carlo simulations at 25°C and 1 atm with statistical perturbation theory. The average difference between the computed results and experimental data for the absolute free energies of hydration is 0.5 kcal/mol. The AM1-SM2 method is also found to perform well in predicting the free energies of hydration for the substituted benzenes. In addition, the Monte Carlo simulations provided details on the hydration of the substituted benzenes, in particular for the solute–water hydrogen bonding. © 1993 John Wiley & Sons, Inc.     

Investigation of the salting out of methane from aqueous electrolyte solutions using computer simulations

We calculate the excess chemical potential of methane in aqueous electrolyte solutions of NaCl using Monte Carlo computer simulations. In a recent work [Docherty et al. J. Chem. Phys. 2006, 125, 074510], we presented a new potential model for methane in water which is capable of describing accurately the excess chemical potential of methane in pure water over a range of temperatures, a quantity that can be related to the solubility and which is commonly used to study the hydrophobic effect. Here, we use the same potential model for the water-methane interactions and investigate the effect of added salt on the chemical potential of methane in the solution. The methane molecules are modeled as single Lennard-Jones (LJ) interaction sites, and the water molecules are modeled with the TIP4P/2005 model. A correcting factor of chi = 1.07 for the energetic Berthelot (geometric) combining rule of the methane-water interaction is also used, which mimics the polarization of methane in water. We consider NaCl as the salt and treat the ions with the Smith and Dang model (i.e., as charged LJ interaction sites). Ion-water, ion-ion, and ion-methane interactions are treated using Lorentz-Berthelot combining rules. In addition, the Coulombic potential is used to model charge-charge interactions which are calculated using the Ewald sum. We have carried out isobaric-isothermal (NpT) simulations to determine the equilibrium densities of the solutions. The simulation data is in excellent agreement with experimental densities of aqueous NaCl solutions of different concentration. Hydration numbers are also obtained and found to be in agreement with reported data. Canonical (NVT) simulations at the averaged densities are then performed using the Widom test-particle insertion method to obtain the excess chemical potential of methane in the saline solutions. An increase in the chemical potential of methane, corresponding to a salting out effect, is observed when salt is added to the solution. We investigate different concentrations and ion sizes. An overprediction of the salting out effect as compared with experimental data is observed, which we believe is due to the polarizing effect of the ions in the solution, which is not taken into account by the model. We also find a direct correlation between the increase in the chemical potential and the packing fraction of the solution and argue that the main cause of the observed salting out effect (as represented by an increase in the excess chemical potential) is the increase in the packing fraction of the solutions due to the added salt. Together, with this, we put forward an argument toward explaining the anomalous Hofmeister effect of Li(+).     

Titration of hydrophobic polyelectrolytes using Monte Carlo simulations

The conformation and titration curves of weak (or annealed) hydrophobic polyelectrolytes have been examined using Monte Carlo simulations with screened Coulomb potentials in the grand canonical ensemble. The influence of the ionic concentration pH and presence of hydrophobic interactions has been systematically investigated. A large number of conformations such as extended, pearl-necklace, cigar-shape, and collapsed structures resulting from the subtle balance of short-range hydrophobic attractive interactions and long-range electrostatic repulsive interactions between the monomers have been observed. Titration curves were calculated by adjusting the pH-pK(0) values (pK(0) represents the intrinsic dissociation constant of an isolated monomer) and then calculating the ionization degree alpha of the polyelectrolyte. Important transitions related to cascades of conformational changes were observed in the titration curves, mainly at low ionic concentration and with the presence of strong hydrophobic interactions. We demonstrated that the presence of hydrophobic interactions plays an important role in the acid-base properties of a polyelectrolyte in promoting the formation of compact conformations and hence decreasing the polyelectrolyte degree of ionization for a given pH-pK(0) value.     

A new technique for metropolis Monte Carlo simulation to capture aggregate structures of fine particles: Cluster-moving Monte Carlo algorithm

We present a new Monte Carlo simulation procedure which is capable of capturing aggregate structures in a suspension where fine particles are dispersed. The algorithm we call the “cluster-moving” Monte Carlo algorithm involves moving aggregates (clusters) as unitary particles at every certain Monte Carlo step. We discuss here the theoretical background of the cluster-moving Monte Carlo algorithm and the availability of the algorithm for simulations of systems where fine particles aggregate. The results of simulations for two model systems, magnetic fluids and colloidal dispersions, have shown that the new algorithm produces much more rapid convergence than the conventional one for unstable dispersion systems and reproduces physically reasonable aggregate structures of fine particles.     

Optimization of 89Zr production using Monte Carlo simulations

Zirconium-89 (⁸⁹Zr) can be produced in a cyclotron by focusing the proton beam on an yttrium-89 (⁸⁹Y) foil target. Optimal combination of beam energy and target assembly configuration enables maximum production of ⁸⁹Zr while minimizing the formation of contaminant nuclides such as ⁸⁸Zr and ⁸⁸Y to allow efficient and effective radiopharmaceutical labeling. Accurate modeling of the proton beam and the target is therefore an essential step to assure the best beam and target specification. We used the radiation transport code MCNPX to simulate the transport of protons through the irradiation assembly and the nuclear reaction code TALYS to obtain the production cross sections of various nuclides from proton-⁸⁹Y reactions. Results from simulating the irradiation of 14 mm diameter targets with aluminum (Al) degrader at 19.8 mA for 1 h suggested that the 0.15 mm thick one would produce 227 MBq while the 0.3 mm thick one would produce 413 MBq of ⁸⁹Zr with less than 1 % uncertainty. These results show excellent agreement with experimental work in literature. This work provides the basis for further experimental and theoretical assessments of the use of ⁸⁹Zr as radiopharmaceutical labels.     

Hybrid Monte Carlo simulations of vertical electronic transitions in acetone in aqueous solution

A simulation of the absorption and the fluorescence of acetone in aqueous solution is reported. The model has an explicit solvent representation with an effective ab initio treatment of the solute. The model attempts to balance quantum chemistry, intermolecular interactions and statistical thermodynamics. It includes a non-electrostatic perturbation on the solute which models the solute–solvent exchange repulsion and the restriction put on the electronic structure of the solute by the antisymmetry to the solvent. The solvent shift to the absorption transition is found to be between 0.16 and 0.21 eV; the shift to the fluorescence transition is found to be between 0.02 and 0.05 eV. The simulation supports the conclusion that the first peak in the fluorescence spectrum of acetone is from a single molecule in equilibrium with the solvent, not from an excimer.     

Stochastic potential switching algorithm for Monte Carlo simulations of complex systems

This paper describes a new Monte Carlo method based on a novel stochastic potential switching algorithm. This algorithm enables the equilibrium properties of a system with potential V to be computed using a Monte Carlo simulation for a system with a possibly less complex stochastically altered potential V. By proper choices of the stochastic switching and transition probabilities, it is shown that detailed balance can be strictly maintained with respect to the original potential V. The validity of the method is illustrated with a simple one-dimensional example. The method is then generalized to multidimensional systems with any additive potential, providing a framework for the design of more efficient algorithms to simulate complex systems. A near-critical Lennard-Jones fluid with more than 20,000 particles is used to illustrate the method. The new algorithm produced a much smaller dynamic scaling exponent compared to the Metropolis method and improved sampling efficiency by over an order of magnitude.     

Flexible docking simulations: Scaled collective variable Monte Carlo minimization approach using Bezier splines,and comparison with a standard Monte Carlo algorithm

An algorithm for docking a flexible ligand onto a flexible or rigid receptor, using the scaled‐collective‐variables Monte Carlo with energy minimization approach, is presented. Energy minimization is shown to be one of the best techniques for distinguishing between native‐ and nonnative‐generated conformations. Incorporation of this technique into a Monte Carlo procedure enables one to distinguish the native conformation directly during the conformational search. It avoids the generation of a large number of ligand conformers for which more sophisticated energy evaluation tools would have had to be applied to identify the nativelike conformations. The efficiency of the Monte Carlo minimization was greatly improved by incorporating a new grid‐based energy evaluation technique using Bezier splines for which the energy function, as well as all of its derivatives, can be deduced from the values at the grid points. Comparison between our ECEPP/3‐based algorithm and the Monte Carlo algorithm presented elsewhere (Hart, T. N.; Read, R. J. Prot Struct Funct Genet 1992, 13, 206–222) has been made for docking NH₂ D Phe Pro Arg COOH, the noncovalent analog of NH₂ D Phe Pro Arg chloromethylketone (PPACK), onto the active site of human α‐thrombin. ©1999 John Wiley & Sons, Inc. J Comput Chem 20: 244–252, 1999     

Three-dimensional Monte Carlo simulations of the dynamics of macromolecular particles in solutions flowing in mesopores

Numerical simulations are developed to calculate the dynamic equilibrium probability distribution functions (PDF) for macromolecular rod-like particles suspended in a fluid under hydrodynamic flow inside mesopores. The simulations take into account the effects of Brownian and hydrodynamic forces acting on the particles, as well as diffusive collisions of the particles with the solid surface boundaries. An algorithm is developed for this purpose based on Jeffery’s equations for the dynamics of ellipsoidal objects in bulk fluids, and on a mechanism of restitution for the diffusive collisions. The results are presented with a focus on the depletion layer next to two types of solid boundaries, ideally flat and rough. They demonstrate the significance of numerical simulations in 3D compared to previous results based on a 2D approach. In particular, we are able to obtain a complete topography for the PDFs segmented as a hierarchy in the depletion layer.     

Computational methods in coupled electron-ion Monte Carlo simulations.

In the last few years, we have been developing a Monte Carlo simulation method to cope with systems of many electrons and ions in the Born-Oppenheimer approximation: the coupled electron-ion Monte Carlo method (CEIMC). Electronic properties in CEIMC are computed by quantum Monte Carlo rather than by density functional theory (DFT) based techniques. CEIMC can, in principle, overcome some of the limitations of the present DFT-based ab initio dynamical methods. The new method has recently been applied to high-pressure metallic hydrogen. Herein, we present a new sampling algorithm that we have developed in the framework of the reptation quantum Monte Carlo method chosen to sample the electronic degrees of freedom, thereby improving its efficiency. Moreover, we show herein that, at least for the case of metallic hydrogen, variational estimates of the electronic energies lead to an accurate sampling of the proton degrees of freedom.     

Macroion solutions in the cell model studied by field theory and Monte Carlo simulations

Aqueous solutions of charged spherical macroions with variable dielectric permittivity and their associated counterions are examined within the cell model using a field theory and Monte Carlo simulations. The field theory is based on separation of fields into short- and long-wavelength terms, which are subjected to different statistical-mechanical treatments. The simulations were performed by using a new, accurate, and fast algorithm for numerical evaluation of the electrostatic polarization interaction. The field theory provides counterion distributions outside a macroion in good agreement with the simulation results over the full range from weak to strong electrostatic coupling. A low-dielectric macroion leads to a displacement of the counterions away from the macroion.     

A Monte Carlo simulation study of the aqueous hydration of d(CGCGCG) in Z form

Summary Monte Carlo computer simulation is described for the hexamer d(CGCGCG) in the Z form together with 910 water molecules at an environmental density of 1 gm/cc in a cubic cell under periodic boundary conditions. Water-water interactions were treated using the TIP4P potential and the solute water interactions by TIP4P spliced with the non-bonded interactions from the AMBER 3.0 force field. The simulation was subjected to proximity analysis to obtain solute coordination numbers and pair interaction energies for each solute atom. Hydration density distributions partitioned into contributions from the major groove, the minor groove, and the sugar-phosphate backbone were examined, and the probabilities of occurrence for one- and two-water bridges in the simulation were enumerated. The results are compared with observations of crystallographic ordered water sites from the X-ray diffraction studies.     

Reactive Monte Carlo and grand-canonical Monte Carlo simulations of the propene metathesis reaction system

The influence of silicalite-1 pores on the reaction equilibria and the selectivity of the propene metathesis reaction system in the temperature range between 300 and 600 K and the pressure range from 0.5 to 7 bars has been investigated with molecular simulations. The reactive Monte Carlo (RxMC) technique was applied for bulk-phase simulations in the isobaric-isothermal ensemble and for two phase systems in the Gibbs ensemble. Additionally, Monte Carlo simulations in the grand-canonical ensemble (GCMC) have been carried out with and without using the RxMC technique. The various simulation procedures were combined with the configurational-bias Monte Carlo approach. It was found that the GCMC simulations are superior to the Gibbs ensemble simulations for reactions where the bulk-phase equilibrium can be calculated in advance and does not have to be simulated simultaneously with the molecules inside the pore. The confined environment can increase the conversion significantly. A large change in selectivity between the bulk phase and the pore phase is observed. Pressure and temperature have strong influences on both conversion and selectivity. At low pressure and temperature both conversion and selectivity have the highest values. The effect of confinement decreases as the temperature increases.