Auxiliary-field quantum Monte Carlo calculations of molecular systems with a Gaussian basis

We extend the recently introduced phaseless auxiliary-field quantum Monte Carlo (QMC) approach to any single-particle basis and apply it to molecular systems with Gaussian basis sets. QMC methods in general scale favorably with the system size as a low power. A QMC approach with auxiliary fields, in principle, allows an exact solution of the Schrodinger equation in the chosen basis. However, the well-known sign/phase problem causes the statistical noise to increase exponentially. The phaseless method controls this problem by constraining the paths in the auxiliary-field path integrals with an approximate phase condition that depends on a trial wave function. In the present calculations, the trial wave function is a single Slater determinant from a Hartree-Fock calculation. The calculated all-electron total energies show typical systematic errors of no more than a few millihartrees compared to exact results. At equilibrium geometries in the molecules we studied, this accuracy is roughly comparable to that of coupled cluster with single and double excitations and with noniterative triples [CCSD(T)]. For stretched bonds in H(2)O, our method exhibits a better overall accuracy and a more uniform behavior than CCSD(T).     

Generalizations for the potential of mean force between two isolated colloidal particles from Monte Carlo simulations

A substantial amount of experimental and numerical evidence has shown that the Derjaguin-Landau-Verwey-Overbeek theory is not suitable for describing those colloidal solutions that contain multivalent counterions. Toward improved understanding of such solutions, the authors report Monte Carlo calculations wherein, following Rouzina and Bloomfield, they postulate that, in the absence of van der Waals forces, the overall force between two isolated charged colloidal particles in electrolyte solutions is determined by a dimensionless parameter Gamma=z(2)l(B)/a, which measures the electrostatic repulsion between counterions adsorbed on the macroion surface, where z = counterion valence, l(B)=Bjerrum length, and a = average separation between counterions on the macroion surface calculated as if the macroion were fully neutralized. The authors find, first, that the maximum repulsion between like-charged macroions occurs at Gamma approximately 0.5 and, second, that onset of attraction occurs at Gamma approximately 1.8, essentially independent of the valence and concentration of the surrounding electrolyte. These observations might provide new understanding of interactions between electrostatic double layers and perhaps offer explanations for some electrostatic phenomena related to interactions between DNA molecules or proteins.     

Energy-consistent pseudopotentials for quantum Monte Carlo calculations

The authors present scalar-relativistic energy-consistent Hartree-Fock pseudopotentials for the main-group elements. The pseudopotentials do not exhibit a singularity at the nucleus and are therefore suitable for quantum Monte Carlo (QMC) calculations. They demonstrate their transferability through extensive benchmark calculations of atomic excitation spectra as well as molecular properties. In particular, they compute the vibrational frequencies and binding energies of 26 first- and second-row diatomic molecules using post-Hartree-Fock methods, finding excellent agreement with the corresponding all-electron values. They also show their pseudopotentials give superior accuracy than other existing pseudopotentials constructed specifically for QMC. Finally, valence basis sets of different sizes (VnZ with n=D,T,Q,5 for first and second rows, and n=D,T for third to fifth rows) optimized for our pseudopotentials are also presented.     

Parallel and distributed computing in problems of supercomputer simulation of molecular liquids by the Monte Carlo method

An effective strategy of molecular Monte Carlo simulation is proposed. The strategy is based on a combination of two key approaches to parallel computing. The advantage of spatial (domain) decomposition is the high scalability of computing algorithms by splitting “big tasks” into several simultaneously solvable subtasks. However, the domain size in this method can be reduced to a certain limit only. Particle decomposition (division of program loops into portions) is, by contrast, very efficient in the study of small and medium size objects, but is poorly scalable and quickly exhausts the computer system memory with increasing size of the model. The combination of the approaches helps neutralize their limitations and create efficient supercomputing programs for the study of molecular models consisting of hundreds of millions of atoms.     

Quantum Monte Carlo calculations of the dimerization energy of borane

Accurate thermodynamic data are required to improve the performance of chemical hydrides that are potential hydrogen storage materials. Boron compounds are among the most interesting candidates. However, different experimental measurements of the borane dimerization energy resulted in a rather wide range (-34.3 to -39.1) ± 2 kcal/mol. Diffusion Monte Carlo (DMC) simulations usually recover more than 95% of the correlation energy, so energy differences rely less on error cancellation than other methods. DMC energies of BH(3), B(2)H(6), BH(3)CO, CO, and BH(2)(+) allowed us to predict the borane dimerization energy, both via the direct process and indirect processes such as the dissociation of BH(3)CO. Our D(e) = -43.12(8) kcal/mol, corrected for the zero point energy evaluated by considering the anharmonic contributions, results in a borane dimerization energy of -36.59(8) kcal/mol. The process via the dissociation of BH(3)CO gives -34.5(2) kcal/mol. Overall, our values suggest a slightly less D(e) than the most recent W4 estimate D(e) = -44.47 kcal/mol [A. Karton and J. M. L. Martin, J. Phys. Chem. A 111, 5936 (2007)]. Our results show that reliable thermochemical data for boranes can be predicted by fixed node (FN)-DMC calculations.     

Quantum Monte Carlo calculations on Be and LiH

Quantum Monte Carlo (QMC) methods, as recently developed for molecular systems, are applied to Be and LiH. The importance of the trial wavefunction, written as a product of a correlation factor and an orbital part, is emphasised. It is shown that significant improvements in the accuracy of the approach are achieved if multi-configuration wavefunctions are used in preference to self-consistent field wavefunctions. Various forms of the correlation factor are investigated.     

Maximum probability domains from Quantum Monte Carlo calculations

Although it would be tempting to associate the Lewis structures to the maxima of the squared wave function |Psi|2, we prefer in this paper the use of domains of the three-dimensional space, which maximize the probability of containing opposite-spin electron pairs. We find for simple systems (CH4, H2O, Ne, N2, C2H2) domains comparable to those obtained with the electron localization function (ELF) or by localizing molecular orbitals. The different domains we define can overlap, and this gives an interesting physical picture of the floppiness of CH5+ and of the symmetric hydrogen bond in FHF-. The presence of multiple solutions has an analogy with resonant structures, as shown in the trans-bent structure of Si2H2. Correlated wave functions were used (MCSCF or Slater-Jastrow) in the Variational Quantum Monte Carlo framework.     

Dissociation energy of the water dimer from quantum Monte Carlo calculations

We report a study of the electronic dissociation energy of the water dimer using quantum Monte Carlo techniques. We have performed variational quantum Monte Carlo and diffusion quantum Monte Carlo (DMC) calculations of the electronic ground state of the water monomer and dimer using all-electron and pseudopotential approaches. We have used Slater-Jastrow trial wave functions with B3LYP type single-particle orbitals, into which we have incorporated backflow correlations. When backflow correlations are introduced, the total energy of the water monomer decreases by about 4-5 mhartree, yielding a DMC energy of -76.428 30(5) hartree, which is only 10 mhartree above the experimental value. In our pseudopotential DMC calculations, we have compared the total energies of the water monomer and dimer obtained using the locality approximation with those from the variational scheme recently proposed by Casula [Phys. Rev. B 74, 161102(R) (2006)]. The time step errors in the Casula scheme are larger, and the extrapolation of the energy to zero time step always lies above the result obtained with the locality approximation. However, the errors cancel when energy differences are taken, yielding electronic dissociation energies within error bars of each other. The dissociation energies obtained in our various all-electron and pseudopotential calculations range between 5.03(7) and 5.47(9) kcalmol and are in good agreement with experiment. Our calculations give monomer dipole moments which range between 1.897(2) and 1.909(4) D and dimer dipole moments which range between 2.628(6) and 2.672(5) D.     

Quantum Monte Carlo calculations for ground and excited states

A brief overview of the diffusion quantum Monte Carlo method is given. We illustrate the application to ground‐state calculations by a study of the relative stability of carbon clusters near the crossover to fullerene stability, thereby determining the smallest stable fullerene. The application to excited states is illustrated via a study of excitonic states in small hydrogenated silicon clusters. © 2001 John Wiley & Sons, Inc. Int J Quantum Chem, 2001     

Reliable radical stabilization energies from diffusion Monte Carlo calculations

We assess the performance of variational (VMC) and diffusion (DMC) quantum Monte Carlo methods for calculating the radical stabilization energies of a set of 43 carbon-centered radical species. Even using simple single-determinant trial wavefunctions, both methods perform exceptionally well, with mean absolute deviations from reference values well under the chemical accuracy standard of 1 kcal/mol. In addition, the use of DMC results in a highly concentrated spread of errors, with all 43 results within chemical accuracy at the 95% confidence level. These results indicate that DMC is an extremely reliable method for calculating radical stabilization energies and could be used as a benchmark method for larger systems in future.     

Monte Carlo free energy calculations using electronic structure methods

The molecular mechanics-based importance sampling function (MMBIF) algorithm [R. Iftimie, D. Salahub, D. Wei, and J. Schofield, J. Chem. Phys. 113, 4852 (2000)] is extended to incorporate semiempirical electronic structure methods in the secondary Markov chain, creating a fully quantum mechanical Monte Carlo sampling method for simulations of reactive chemical systems which, unlike the MMBIF algorithm, does not require the generation of a system-specific force field. The algorithm is applied to calculating the potential of mean force for the isomerization reaction of HCN using thermodynamic integration. Constraints are implemented in the sampling using a modification of the SHAKE algorithm, including that of a fixed, arbitrary reaction coordinate. Simulation results show that sampling efficiency with the semiempirical secondary potential is often comparable in quality to force fields constructed using the methods suggested in the original MMBIF work. The semiempirical based importance sampling method presented here is a useful alternative to MMBIF sampling as it can be applied to systems for which no suitable MM force field can be constructed.     

Quantum Monte Carlo calculations of the dissociation energy of the water dimer

We report diffusion quantum Monte Carlo (DMC) calculations of the equilibrium dissociation energy D(e) of the water dimer. The dissociation energy measured experimentally, D(0), can be estimated from D(e) by adding a correction for vibrational effects. Using the measured dissociation energy and the modern value of the vibrational energy Mas et al., [J. Chem. Phys. 113, 6687 (2000)] leads to D(e)=5.00+/-0.7 kcal mol(-1), although the result Curtiss et al., [J. Chem. Phys. 71, 2703 (1979)] D(e)=5.44+/-0.7 kcal mol(-1), which uses an earlier estimate of the vibrational energy, has been widely quoted. High-level coupled cluster calculations Klopper et al., [Phys. Chem. Chem. Phys. 2, 2227 (2000)] have yielded D(e)=5.02+/-0.05 kcal mol(-1). In an attempt to shed new light on this old problem, we have performed all-electron DMC calculations on the water monomer and dimer using Slater-Jastrow wave functions with both Hartree-Fock approximation (HF) and B3LYP density functional theory single-particle orbitals. We obtain equilibrium dissociation energies for the dimer of 5.02+/-0.18 kcal mol(-1) (HF orbitals) and 5.21+/-0.18 kcal mol(-1) (B3LYP orbitals), in good agreement with the coupled cluster results.     

Polymorphism in simple liquids: a Gibbs ensemble Monte Carlo study

We perform Gibbs ensemble Monte Carlo (GEMC) simulations of a one-component system of hard spheres with a repulsive shoulder and an attractive well. We show the existence of two distinct liquid-gas and liquid-liquid phase equilibria. The GEMC estimate of the critical parameters, as following from an interpolation of the binodal points, is only slightly influenced by finite size effects. The liquid-gas critical temperature and pressure are lower than those of the liquid-liquid phase separation. A discussion of our findings in comparison with those of previous numerical studies is also presented.     

Visualizing molecular wavefunctions using Monte Carlo methods

Using explicitly correlated wavefunctions and variational Monte Carlo we calculate the electron density, the electron density difference, the intracule density, the extracule density, two forms of the kinetic energy density, the Laplacian of the electron density, the Laplacian of the intracule density, and the Laplacian of the extracule density on a dense grid of points for the ground state of the hydrogen molecule at three internuclear distances (0.6, 1.4, 8.0). With these values we construct a contour plot of each function and describe how it can be used to visualize the distribution of electrons in this molecule. We also examine the influence of electron correlation on each expectation value by calculating each function with a Hartree–Fock wavefunction and then comparing these values with our explicitly correlated values. © 2008 Wiley Periodicals, Inc. Int J Quantum Chem, 2009     

Monte Carlo and modified Tanford-Kirkwood results for macromolecular electrostatics calculations

The understanding of electrostatic interactions is an essential aspect of the complex correlation between structure and function of biological macromolecules. It is also important in protein engineering and design. Theoretical studies of such interactions are predominantly done within the framework of Debye-Hückel theory. A classical example is the Tanford-Kirkwood (TK) model. Besides other limitations, this model assumes an infinitesimally small macromolecule concentration. By comparison to Monte Carlo (MC) simulations, it is shown that TK predictions for the shifts in ion binding constants upon addition of salt become less reliable even at moderately macromolecular concentrations. A simple modification based on colloidal literature is suggested to the TK scheme. The modified TK models suggested here satisfactorily predict MC and experimental shifts in the calcium binding constant as a function of protein concentration for the calbindin D(9k) mutant and calmodulin.     

Quantum Monte Carlo calculations of the potential energy curve of the helium dimer

We report results of two quantum Monte Carlo methods -- variational Monte Carlo and diffusion Monte Carlo -- on the potential energy curve of the helium dimer. In contrast to previous quantum Monte Carlo calculations on this system, we have employed trial wave functions of the Slater-Jastrow form and used the fixed node approximation for the fermion nodal surface. We find both methods to be in excellent agreement with the best theoretical results at short range. In addition, the diffusion Monte Carlo results give very good agreement across the whole potential energy curve, while the Slater-Jastrow wave function fails to bind the dimer at all.     

Ranking protein–protein docking results using steered molecular dynamics and potential of mean force calculations

Crystallization of protein–protein complexes can often be problematic and therefore computational structural models are often relied on. Such models are often generated using protein–protein docking algorithms, where one of the main challenges is selecting which of several thousand potential predictions represents the most near‐native complex. We have developed a novel technique that involves the use of steered molecular dynamics (sMD) and umbrella sampling to identify near‐native complexes among protein–protein docking predictions. Using this technique, we have found a strong correlation between our predictions and the interface RMSD (iRMSD) in ten diverse test systems. On two of the systems, we investigated if the prediction results could be further improved using potential of mean force calculations. We demonstrated that a near‐native (<2.0 Å iRMSD) structure could be identified in the top‐1 ranked position for both systems. © 2016 Wiley Periodicals, Inc.