A practical integral equation for the structure and thermodynamics of hard sphere Coulomb fluids

A closure for the Ornstein-Zernike equation is presented, applicable for fluids of charged, hard spheres. From an exact, but intractable closure, we derive the radial distribution function of nonlinearized Debye-Hückel theory by subsequent approximations, and use the information to formulate a new closure by an extension of the mean spherical approximation. The radial distribution functions of the new closure, coined Debye-Hückel-extended mean spherical approximation, are in excellent agreement with those resulting from the hyper-netted chain approximation and molecular dynamics simulations, in the regime where the latter are applicable, except for moderately dilute systems at low temperatures where the structure agrees at most qualitatively. The method is numerically more efficient, and more important, convergent in the entire temperature-density plane. We demonstrate that the method is accurate under many conditions for the determination of the structural and thermodynamic properties of homogeneous, symmetric hard-sphere Coulomb systems, and estimate it to be a valuable basis for the formulation of density functional theories for inhomogeneous or highly asymmetric systems.     

New evaluation of reconstructed spatial distribution function from radial distribution functions

Although three dimensional (3D) solvation structure is much more informative than one dimensional structure, its evaluation is difficult experimentally and theoretically. In our previous Communication [Yokogawa et al., J. Chem. Phys. 123, 211102 (2005)], we proposed a new method to present reconstructed spatial distribution function (RC-SDF) from a set of radial distribution functions (RDFs). In this article, we successfully extended the method more accurately with new basis sets. This new method was applied to two liquid solvation structures, methanol and dimethyl sulfoxide, as examples. Their RC-SDFs evaluated here clearly show that the former solvation structure is well defined while the latter one is broad, which agrees well with the SDFs calculated directly from molecular dynamics simulations. These results indicate that the method can reproduce well these 3D solvation structures in reasonable computational cost.     

Structure and phase behavior of Widom-Rowlinson model calculated from a nonuniform Ornstein-Zernike equation

The second-order integral-equation formalism of [Attard J. Chem. Phys. 91, 3072 (1989); 95, 4471 (1991)], applied previously to one-component hard spheres and Lennard-Jones fluids, as well as to their mixtures, is used to binary Widom-Rowlinson mixtures. Comparison with Monte Carlo simulations of the pair correlation functions and of the demixing phase diagram shows that this method is also quite accurate in the case of highly nonadditive mixtures. Moreover, the results of the second-order theory are compared with previous theoretical predictions. Our interest is also in the calculation of the bridge functions, i.e., parts of the radial distribution functions either not included or simply approximated in the usual theories.     

Pair correlation function of soft-sphere fluids

A closed-form analytic formula for the radial distribution function (RDF) or g(r) of inverse power fluids is proposed. The RDF is expressed as a sum of separate component functions, one monotonic and a series of exponentially damped oscillatory functions. Unlike previous treatments in the literature, this formula does not rely on patching different functional forms at arbitrary crossover distances. This expression, which we refer to as g(M)(r), yields the expected asymptotic behavior at large distance and reproduces the main features of the RDF generated by molecular dynamics (MD) simulations. The g(M) is applied to the soft n = 4 inverse power fluid, and it is shown that in this case seven or fewer terms are sufficient to represent accurately the MD-generated RDF over the entire fluid domain. The relative contributions of the separate terms of the g(M) as a function of density are analyzed and discussed. The key role played by the monotonic component function and two oscillatory terms is demonstrated. The origin of the crossover from the oscillatory to the monotonic behavior is shown to be the same as that recently proposed by Evans and Henderson [R. Evans and J. R. Henderson, J. Phys.: Condens. Matter 21, 474220 (2009)] for the dispersion interactions.     

A transferable coarse-grained model for hydrogen-bonding liquids

We present here a recent development of a generalized coarse-grained model for use in molecular simulations. In this model, interactions between coarse-grained particles consist of both van der Waals and explicit electrostatic components. As a result, the coarse-grained model offers the transferability that is lacked by most current effective-potential based approaches. The previous center-of-mass framework (P. A. Golubkov and P. Ren, J. Chem. Phys., 2006, 125, 64103) is generalized here to include arbitrary off-center interaction sites for both Gay-Berne and multipoles. The new model has been applied to molecular dynamic simulations of neat methanol liquid. By placing a single point multipole at the oxygen atom rather than at the center of mass of methanol, there is a significant improvement in the ability to capture hydrogen-bonding. The critical issue of transferability of the coarse-grained model is verified on methanol-water mixtures, using parameters derived from neat liquids without any modification. The mixture density and internal energy from coarse-grained molecular dynamics simulations show good agreement with experimental measurements, on a par with what has been obtained from more detailed atomic models. By mapping the dynamics trajectory from the coarse-grained simulation into the all-atom counterpart, we are able to investigate atomic-level structure and interaction. Atomic radial distribution functions of neat methanol, neat water and mixtures compare favorably to experimental measurements. Furthermore, hydrogen-bonded 6- and 7-molecule chains of water and methanol observed in the mixture are in agreement with previous atomic simulations.     

NMR relaxation parameters of a Lennard-Jones fluid from molecular-dynamics simulations

Ensembles of soft spheres or of Lennard-Jones atoms were studied by molecular dynamics at reduced temperatures from 0.8 to 3, and radial distribution functions, diffusion coefficients, and magnetic dipole-dipole correlation functions were measured as functions of system size. The expected relation between the values of the correlation functions at zero lag time and the integrals of the radial distribution was verified for each system. The measured correlation functions were compared with theoretical expressions derived by [Ayant et al., J. Phys. (Paris) 36, 991 (1975)] and by [Hwang and Freed, J. Chem. Phys. 63, 4017 (1975)]. It was shown that, in order to recover the long-time behavior characteristic of diffusion-controlled relaxation processes, the simulation must comprise at least 10 000 particles. By fitting the simulation results to the Hwang-Freed function, independent values of the diffusion coefficient were obtained, similar but not identical to those computed using the Green-Kubo formalism. The spectral densities of the dipole-dipole interaction were computed as Fourier transforms of the correlation functions. These quantities are less sensitive to model imperfections and reproduce quite well the values derived from theory. The dimensionless spin-lattice and spin-spin relaxation rates were derived from the spectral densities. It was shown that the spin-lattice (longitudinal) relaxation rate goes through a maximum as the temperature increases, while the spin-spin (transverse) rate decreases monotonously.     

A hybrid scheme for the resolution-of-the-identity approximation in second-order Møller-Plesset linear-r(12) perturbation theory

In the framework of second-order M?ller-Plesset linear-r(12) (MP2-R12) perturbation theory, a method is developed and implemented that uses an auxiliary basis set for the resolution-of-the-identity (RI) approximation for the three- and four-electron integrals. In contrast to previous work, the two-electron integrals that must be evaluated never involve more than one auxiliary basis function. The new method therefore scales linearly with the number of auxiliary basis functions and is much more efficient than the previous one, which scaled quadratically. A general formulation of MP2-R12 theory is presented for various ansatze, approximations, and orbitals (canonical or localized). The new method is assessed by computations of the valence-shell second-order M?ller-Plesset correlation energy of a few small closed-shell systems. The preliminary calculations indicate that the difference between the new and previous methods is about one order of magnitude smaller than the errors that occur due to basis-set truncations and RI approximations and under the assumptions of generalized and extended Brillouin conditions.     

Estimating statistical distributions using an integral identity

We present an identity for an unbiased estimate of a general statistical distribution. The identity computes the distribution density from dividing a histogram sum over a local window by a correction factor from a mean-force integral, and the mean force can be evaluated as a configuration average. We show that the optimal window size is roughly the inverse of the local mean-force fluctuation. The new identity offers a more robust and precise estimate than a previous one by Adib and Jarzynski [J. Chem. Phys. 122, 014114 (2005)]. It also allows a straightforward generalization to an arbitrary ensemble and a joint distribution of multiple variables. Particularly we derive a mean-force enhanced version of the weighted histogram analysis method. The method can be used to improve distributions computed from molecular simulations. We illustrate the use in computing a potential energy distribution, a volume distribution in a constant-pressure ensemble, a radial distribution function, and a joint distribution of amino acid backbone dihedral angles.     

Numerical integration of atomic electron density with double exponential formula for density functional calculation

In the preceding study, we reported an application of the double exponential formula to the radial quadrature grid for numerical integration of the radial electron distribution function. Three-type new radial grids with the double exponential transformation were introduced. The performance of radial grids was compared between the double exponential grids and the grids proposed in earlier studies by applying to the electron-counting integrals of noble gas atoms and diatomic molecules including alkali metals, halogens, and transition metals. It was confirmed that the change in accuracy of the quadrature approximation depending on atomic or molecular species is not significant for the double exponential integration schemes rather than the other integration schemes. In the present study, we further investigate the accuracy of the double exponential formula for the electron-counting integrals of all the atoms from H to Kr in the periodic table to elucidate the stable performance of the double exponential radial grids. The electron densities of the atoms are calculated with the Gauss-type orbital basis functions at the B3LYP level. The quadrature accuracy and convergence behavior of numerical integration are compared among the double exponential formula and the formulas proposed by Treutler et al. and by Mura et al. The results reveal that the double exponential radial grids remarkably improve the convergence rate toward high accuracy compared with the previous radial grids, particularly for heavy elements in the 4th period, without fine tuning of the radial grids for each atom.     

Molecular dynamics simulation of liquid trimethylphosphine

Structural and dynamical properties of liquid trimethylphosphine (TMP), (CH(3))(3)P, as a function of temperature is investigated by molecular dynamics (MD) simulations. The force field used in the MD simulations, which has been proposed from molecular mechanics and quantum chemistry calculations, is able to reproduce the experimental density of liquid TMP at room temperature. Equilibrium structure is investigated by the usual radial distribution function, g(r), and also in the reciprocal space by the static structure factor, S(k). On the basis of center of mass distances, liquid TMP behaves like a simple liquid of almost spherical particles, but orientational correlation due to dipole-dipole interactions is revealed at short-range distances. Single particle and collective dynamics are investigated by several time correlation functions. At high temperatures, diffusion and reorientation occur at the same time range as relaxation of the liquid structure. Decoupling of these dynamic properties starts below ca. 220 K, when rattling dynamics of a given TMP molecules due to the cage effect of neighbouring molecules becomes important.     

Liquid water simulations with the density fragment interaction approach

We reformulate the density fragment interaction (DFI) approach [Fujimoto and Yang, J. Chem. Phys., 2008, 129, 054102.] to achieve linear-scaling quantum mechanical calculations for large molecular systems. Two key approximations are developed to improve the efficiency of the DFI approach and thus enable the calculations for large molecules: the electrostatic interactions between fragments are computed efficiently by means of polarizable electrostatic-potential-fitted atomic charges; and frozen fragment pseudopotentials, similar to the effective fragment potentials that can be fitted from interactions between small molecules, are employed to take into account the Pauli repulsion effect among fragments. Our reformulated and parallelized DFI method demonstrates excellent parallel performance based on the benchmarks for the system of 256 water molecules. Molecular dynamics simulations for the structural properties of liquid water also show a qualitatively good agreement with experimental measurements including the heat capacity, binding energy per water molecule, and the radial distribution functions of atomic pairs of O-O, O-H, and H-H. With this approach, large-scale quantum mechanical simulations for water and other liquids become feasible.     

Auxiliary functions for molecular integrals with Slater‐type orbitals. I. Translation methods

Many types of molecular integrals involving Slater functions can be expressed, with the ζ‐function method in terms of sets of one‐dimensional auxiliary integrals whose integrands contain two‐range functions. After reviewing the properties of these functions (including recurrence relations, derivatives, integral representations, and series expansions), we carry out a detailed study of the auxiliary integrals aimed to facilitate both the formal and computational applications of the ζ‐function method. The usefulness of this study in formal applications is illustrated with an example. The high performance in numerical applications is proved by the development of a very efficient program for the calculation of two‐center integrals with Slater functions corresponding to electrostatic potential, electric field, and electric field gradient. © 2006 Wiley Periodicals, Inc. Int J Quantum Chem, 2006     

Single particle force distributions in simple fluids

The distribution function, W(F), of the magnitude of the net force, F, on particles in simple fluids is considered, which follows on from our previous publication [A. C. Bran?ka, D. M. Heyes, and G. Rickayzen, J. Chem. Phys. 135, 164507 (2011)] concerning the pair force, f, distribution function, P(f), which is expressible in terms of the radial distribution function. We begin by discussing the force on an impurity particle in an otherwise pure fluid but later specialize to the pure fluid, which is studied in more detail. An approximate formula, expected to be valid asymptotically, for W(F) referred to as, W(1)(F) is derived by taking into account only binary spatial correlations in the fluid. It is found that W(1)(F) = P(f). Molecular dynamics simulations of W for the inverse power (IP) and Lennard-Jones potential fluids show that, as expected, W(F) and P(f) agree well in the large force limit for a wide range of densities and potential forms. The force at which the maximum in W(F) occurs for the IP fluids follows a different algebraic dependence with density in low and high density domains of the equilibrium fluid. Other characteristic features in the force distribution functions also exhibit the same trends. An exact formula is derived relating W(F) to P(x)(F(x)), the distribution function of the x-cartesian components of the net force, F(x), on a particle. W(F) and P(x)(F(x)) have the same analytical forms (apart from constants) in the low and high force limits.     

New translation method for STOs and its application to calculation of two-center two-electron integrals

The new translation method for Slater-type orbitals (STOs) previously tested in the case of the overlap integral is extended to the calculation of two-center two-electron molecular integrals. The method is based on the exact translation of the regular solid harmonic part of the orbital followed by the series expansion of the residual spherical part in powers of the radial variable. Fair uniform convergence and stability under wide changes in molecular parameters are obtained for all studied two-center hybrid, Coulomb, and exchange repulsion integrals. Ten-digit accuracy in the final numerical results is achieved through multiple precision arithmetic calculation of common angular coefficients and Gaussian numerical integration of some of the analytical formulas resulting for the radial integrals. © 2000 John Wiley & Sons, Inc. Int J Quant Chem 79: 91–100, 2000     

Structural features and molecular assembly of amorphous phosphazenic materials in the bulk--combined theoretical and experimental techniques: Tris-(2,2'-dioxy-1,1'-binaphthyl)cyclotriphosphazen

The structure and the assembly of tris-(2,2'-dioxy-binaphthyl)cyclotriphosphazene [(+)-[NP3(O2C20H12)3], DBNP, in the solid amorphous state was studied using molecular dynamics (MD) including ad hoc quantum mechanically derived force field (FF) parameters, in combination with the energy dispersive X-ray diffraction (EDXD) technique. The atom-atom radial distribution function (RDF) curve obtained through the EDXD experiment revealed low intensity peaks not attributable to the intramolecular distances of the single molecule, but clearly featuring a low energy state of long-distance three-dimensional assembly. The radial distribution functions (RDF) were calculated for various models of DBNP submitted to theoretical MD simulations. Based on the comparison of theoretically calculated RDFs and those obtained from the EDXD experiment, the predominant structural motif of the material in the bulk was found to have DBNP molecules laid one upon the other to form tubular nanostructures. These contain eight DBNP units each (length ca. 46 A) with two and three of these units aligned in parallel and held together. The material can be represented as a bulk of tubular snake-like chains undergoing distortions with a step of eight DBNP units. The bending angles, that vary randomly, attain limited values sufficient to induce disorder and thus nonperiodic structure. The present application of MD simulations combined with EDXD data appear to be a general approach to solve for the first time otherwise intractable issues concerning structural features and assembly of molecular materials in the bulk.     

Molecular mass dependence of point-to-set correlation length scale in polymers

We use a recently proposed metric, termed the point-to-set correlation functions, to probe the molecular weight dependence of the relevant static length scales in glass-forming oligomeric chain liquids of 4, 5, 8, and 10 repeat units. In agreement with the results for simple, monatomic fluids, we find that static length scales of the oligomers increase monotonically when the temperature is lowered towards the glass transition temperature of the fluid. More interestingly, the static length scale increases with increasing chain length. Within the bounds of error in our simulations, the static length scale appears to scale as the radius of gyration of the oligomer, but with a prefactor, which is much larger than unity and which grows with the temperature. The preceding behavior contrasts with the length scales extracted from the radial distribution function of the oligomer system, which is practically independent of the chain length.     

Analytical coarse-grained description for polymer melts

Starting from the Ornstein-Zernike equation the authors derive an analytical theory, at the level of pair correlation functions, which coarse grains polymer melts into liquids of interacting soft colloidal particles. Since it is analytical, the presented coarse-graining approach will be useful in developing multiscale modeling procedures to simulate complex fluids of macromolecules. The accuracy of the theory is tested by its capacity to reproduce the liquid structure, as given by the center-of-mass intermolecular total pair correlation function. The theory is found to agree well with the structure predicted by molecular dynamics simulations of the liquid described at the united atom level as well as by molecular dynamics simulations of the liquid of interacting colloidal particles. The authors perform simulations of the liquid of interacting colloidal particles having as input the potential obtained from their analytical total pair correlation function by enforcing the hypernetted-chain closure approximation. Tests systems are polyethylene melts of chains with increasing degrees of polymerization and polymer melts of chains with different chemical architectures. They also discuss the effect of adopting different conventional approximations for intra- and intermolecular monomer structure factors on the accuracy of the coarse-graining procedure, as well as the relevance of higher-order corrections to their expression.     

A novel algorithm for creating coarse-grained, density dependent implicit solvent models

Implicit solvent simulations are those in which solvent molecules are not explicitly simulated, and the solute-solute interaction potential is modified to compensate for the implicit solvent effect. Implicit solvation is well known in Brownian dynamics of dilute solutions but offers promise to speed up many other types of molecular simulations as well, including studies of proteins and colloids where the local density can vary considerably. This work examines implicit solvent potentials within a more general coarse-graining framework. While a pairwise potential between solute sites is relatively simple and ubiquitous, an additional parametrization based on the local solute concentration has the possibility to increase the accuracy of the simulations with only a marginal increase in computational cost. We describe here a method in which the radial distribution function and excess chemical potential of solute insertion for a system of Lennard-Jones particles are first measured in a fully explicit, all-particle simulation, and then reproduced across a range of solute particle densities in an implicit solvent simulation.     

一维函数方法求解原子和分子薛定谔方程

对薛定谔方程的严格数值求解, 尤其是发展标准方法之外的、 包含新功能的解法, 一直是物理学研究的基本关注点. 本文介绍一种近些年发展的一维函数近似解方法, 该方法通过对波函数的不同坐标分量进行处理来求解原子和分子体系的薛定谔方程. 电子的试探波函数被离散化到实空间均匀格点上, 因此可以通过残差矢量校正的方法对其进行改进. 一维函数方法本身的特征决定其非常利于数值积分, 避免了很多由常规的多电子、 多中心势分子积分所带来的问题. 计算中, 最终能量是从严格的能量上限逐渐收敛所获得, 计算出的两电子薛定谔波函数呈现出常规单电子近似方法所含有的电子关联效应. 不同于密度泛函理论及Hartree-Fock的单电子解法, 本方法对电子-电子排斥能的多体效应的处理更加精确.     

Integral equation theory for correcting truncation errors in molecular simulations

Various strategies for correcting structural and energetic artefacts of molecular simulations with truncated potentials based on integral equation theory are described and applied to liquid water. The performance of the methods is examined for a range of cutoff distances and different shifted-force potentials. With the recently enhanced damped Coulomb potential (D. Zahn, B. Schilling, S.M. Kast, J. Phys. Chem. B, 106 (2002) 10725), parameterised and corrected by integral equation theory, radial distribution functions and excess internal energy very close to the Ewald simulation limit are obtained from a simulation with a cutoff distance of only 6 Å.