Volume 131, Issue 15, 21 October 2009

A fully automated parallelized implementation of the incremental scheme for coupledcluster singlesanddoubles (CCSD) energies has been extended to treat molecular (unrelaxed) firstorder oneelectron properties such as the electric dipole and quadrupole moments. The convergence and accuracy of the incremental approach for the dipole and quadrupole moments have been studied for a variety of chemically interesting systems. It is found that the electric dipole moment can be obtained to within 5% and 0.5% accuracy with respect to the exact CCSD value at the third and fourth orders of the expansion, respectively. Furthermore, we find that the incremental expansion of the quadrupole moment converges to the exact result with increasing order of the expansion: the convergence of nonaromatic compounds is fast with errors less than 16 mau and less than 1 mau at third and fourth orders, respectively ; the aromatic compounds converge slowly with maximum absolute deviations of 174 and 72 mau at third and fourth orders, respectively.
 COMMUNICATIONS


Direct visualization of the H–Xe bond in xenon hydrides: Xenon isotopic shift in the IR spectra
View Description Hide DescriptionIR spectra of xenon hydrides (HXeCCH, HXeCC, and HXeH) obtained from different xenon isotopes ( and ) exhibit a small but detectable and reproducible isotopic shift in the absorptions assigned to H–Xe stretching (by ). To our knowledge, it is the first direct experimental evidence for the H–Xe bond in HXeY type compounds. The shift magnitude is in good agreement with quantumchemical calculations.

Tracing ultrafast hydrogen migration in allene in intense laser fields by tripleion coincidence momentum imaging
View Description Hide DescriptionUltrafast hydrogen migration in allene in intense laser fields was investigated by tripleion coincidence momentum imaging. The migrating proton covering the entire range of an allene molecule was visualized by the momentum correlation maps and by the geometrical structure of triply charged allene reconstructed from the observed momentum vectors of fragment ions. The extent of hydrogen migration was found to play a decisive role in breaking selectively one of the two initially equivalent C–C chemical bonds that become inequivalent in the course of the hydrogen migration.

Density scaling in viscous liquids: From relaxation times to fourpoint susceptibilities
View Description Hide DescriptionWe present numerical calculations of a fourpoint dynamic susceptibility, , for the Kob–Andersen LennardJones mixture as a function of temperature and density . Over a relevant range of and , the full dependence of and thus the maximum in , which is proportional to the dynamic correlation volume, are invariant for state points for which the scaling variable is constant. The value of the material constant is the same as that which superposes the relaxation time of the system versus . Thus, the dynamic correlation volume is a unique function of for any thermodynamic condition in the regime where density scaling holds. Finally, we examine the conditions under which the density scaling properties are related to the existence of strong correlations between pressure and energy fluctuations.
 Top

 ARTICLES

 Theoretical Methods and Algorithms

Automated calculation of anharmonic vibrational contributions to first hyperpolarizabilities: Quadratic response functions from vibrational configuration interaction wave functions
View Description Hide DescriptionQuadratic response functions are derived and implemented for a vibrational configuration interaction state. Combined electronic and vibrational quadratic response functions are derived using Born–Oppenheimer vibronic product wave functions. Computational tractable expressions are derived for determining the total quadratic response contribution as a sum of contributions involving both electronic and vibrational linear and quadratic response functions. In the general frequencydependent case this includes a new and more troublesome type of electronic linear response function. Pilot calculations for the FH, , , and pyrrole molecules demonstrate the importance of vibrational contributions for accurate comparison to experiment and that the vibrational contributions in some cases can be very large. The calculation of transition properties between vibrational states is combined with sumoverstates expressions for analysis purposes. On the basis of this some simple analysis methods are suggested. Also, a preliminary study of the effect of finite lifetimes on quadratic response functions is presented.

Implementation of the incremental scheme for oneelectron firstorder properties in coupledcluster theory
View Description Hide DescriptionA fully automated parallelized implementation of the incremental scheme for coupledcluster singlesanddoubles (CCSD) energies has been extended to treat molecular (unrelaxed) firstorder oneelectron properties such as the electric dipole and quadrupole moments. The convergence and accuracy of the incremental approach for the dipole and quadrupole moments have been studied for a variety of chemically interesting systems. It is found that the electric dipole moment can be obtained to within 5% and 0.5% accuracy with respect to the exact CCSD value at the third and fourth orders of the expansion, respectively. Furthermore, we find that the incremental expansion of the quadrupole moment converges to the exact result with increasing order of the expansion: the convergence of nonaromatic compounds is fast with errors less than 16 mau and less than 1 mau at third and fourth orders, respectively ; the aromatic compounds converge slowly with maximum absolute deviations of 174 and 72 mau at third and fourth orders, respectively.

An imagebased reaction field method for electrostatic interactions in molecular dynamics simulations of aqueous solutions
View Description Hide DescriptionIn this paper, a new solvation model is proposed for simulations of biomolecules in aqueous solutions that combines the strengths of explicit and implicit solvent representations. Solute molecules are placed in a spherical cavity filled with explicit water, thus providing microscopic detail where it is most needed. Solvent outside of the cavity is modeled as a dielectric continuum whose effect on the solute is treated through the reaction field corrections. With this explicit/implicit model, the electrostatic potential represents a solute molecule in an infinite bath of solvent, thus avoiding unphysical interactions between periodic images of the solute commonly used in the latticesum explicit solvent simulations. For improved computational efficiency, our model employs an accurate and efficient multipleimage charge method to compute reaction fields together with the fast multipole method for the direct Coulomb interactions. To minimize the surface effects, periodic boundary conditions are employed for nonelectrostatic interactions. The proposed model is applied to study liquid water. The effect of model parameters, which include the size of the cavity, the number of image charges used to compute reaction field, and the thickness of the buffer layer, is investigated in comparison with the particlemesh Ewald simulations as a reference. An optimal set of parameters is obtained that allows for a faithful representation of many structural, dielectric, and dynamic properties of the simulated water, while maintaining manageable computational cost. With controlled and adjustable accuracy of the multipleimage charge representation of the reaction field, it is concluded that the employed model achieves convergence with only one image charge in the case of pure water. Future applications to calculations, conformational sampling of solvated biomolecules and electrolyte solutions are briefly discussed.

Separating forward and backward pathways in nonequilibrium umbrella sampling
View Description Hide DescriptionUmbrella sampling enforces uniform sampling of steadystate distributions that are functions of arbitrary numbers of order parameters. The key to applying such methods to nonequilibrium processes is the accumulation of fluxes between regions. A significant difference between microscopically reversible and irreversible systems is that, in the latter case, the transition path ensemble for a reaction can be significantly different for “forward” and “backward” trajectories. Here, we show how to separately treat forward and backward pathways in nonequilibrium umbrella sampling simulations by working in an extended space. In this extended space, the exact rate (for equilibrium or nonequilibrium processes) can be calculated “for free” as a flux in phase space. We compare the efficiency of this rate calculation with forward flux sampling for a twodimensional potential and show that nonequilibrium umbrella sampling is more efficient when an intermediate is present. We show that this technique can also be used to describe steadystate limit cycles by examining a simulation of circadian oscillations. We obtain the path of the limit cycle in a space of 22 order parameters, as well as the oscillation period. The relation of our method to others is discussed.

Excitedstate reversible geminate recombination in two dimensions
View Description Hide DescriptionExcitedstate reversible geminate recombination with two different lifetimes and quenching is investigated in two dimensions. From the exact Green function in the Laplace domain, analytic expressions of twodimensional survival and binding probabilities are obtained at short and long times. We find that a new pattern of kinetic transition occurs in two dimensions. The longtime effective survival probabilities show a pattern of depending on the rate constants while the effective binding probabilities show .

The role of the reference state in longrange random phase approximation correlation
View Description Hide DescriptionWe recently presented a combination of a shortrange density functional approximation with longrange random phase approximation (RPA) correlation [B. G. Janesko, T. M. Henderson, and G. E. Scuseria, J. Chem. Phys.130, 081105 (2009)]. Here we explore how this approximation’s performance is affected by the choice of reference state, i.e., the orbitals and orbital energy differences entering the RPA energy expression. Our previous results built the reference state using a nonlocal exchange potential. Rescaling the RPA correlation energy by an empirical factor gave very accurate results for a wide range of properties. We show here that reference states constructed from approximate local exchangecorrelation potentials give their best results with smaller rescaling factors . However, the tested potentials yield artifacts in some systems.

General formulation of pressure and stress tensor for arbitrary manybody interaction potentials under periodic boundary conditions
View Description Hide DescriptionThree distinct forms are derived for the force virial contribution to the pressure and stress tensor of a collection of atoms interacting under periodic boundary conditions. All three forms are written in terms of forces acting on atoms, and so are valid for arbitrary manybody interatomic potentials. All three forms are mathematically equivalent. In the special case of atoms interacting with pair potentials, they reduce to previously published forms. (i) The atomcell form is similar to the standard expression for the virial for a finite nonperiodic system, but with an explicit correction for interactions with periodic images. (ii) The atom form is particularly suited to implementation in modern molecular dynamics simulation codes using spatial decomposition parallel algorithms. (iii) The group form of the virial allows the contributions to the virial to be assigned to individual atoms.

Exact ground state Monte Carlo method for Bosons without importance sampling
View Description Hide DescriptionGenerally “exact” quantum Monte Carlo computations for the ground state of many bosons make use of importance sampling. The importance sampling is based either on a guiding function or on an initial variational wave function. Here we investigate the need of importance sampling in the case of path integral ground state (PIGS) Monte Carlo. PIGS is based on a discrete imaginary time evolution of an initial wave function with a nonzero overlap with the ground state, which gives rise to a discrete path which is sampled via a Metropolislike algorithm. In principle the exact ground state is reached in the limit of an infinite imaginary time evolution, but actual computations are based on finite time evolutions and the question is whether such computations give unbiased exact results. We have studied bulk liquid and solid with PIGS by considering as initial wave function a constant, i.e., the ground state of an ideal Bose gas. This implies that the evolution toward the ground state is driven only by the imaginary time propagator, i.e., there is no importance sampling. For both phases we obtain results converging to those obtained by considering the best available variational wave function (the shadow wave function) as initial wave function. Moreover we obtain the same results even by considering wave functions with the wrong correlations, for instance, a wave function of a strongly localized Einstein crystal for the liquid phase. This convergence is true not only for diagonal properties such as the energy, the radial distribution function, and the static structure factor, but also for offdiagonal ones, such as the onebody density matrix. This robustness of PIGS can be traced back to the fact that the chosen initial wave function acts only at the beginning of the path without affecting the imaginary time propagator. From this analysis we conclude that zero temperature PIGS calculations can be as unbiased as those of finite temperature path integral Monte Carlo. On the other hand, a judicious choice of the initial wave function greatly improves the rate of convergence to the exact results.

Ionspecific thermodynamics of multicomponent electrolytes: A hybrid HNC/MD approach
View Description Hide DescriptionUsing effective infinite dilution ionion interaction potentials derived from explicitwater molecular dynamics (MD) computer simulations in the hypernettedchain (HNC) integral equation theory we calculate the liquid structure and thermodynamic properties, namely, the activity and osmotic coefficients of various multicomponent aqueous electrolyte mixtures. The electrolyte structure expressed by the ionion radial distribution functions is for most ions in excellent agreement with MD and implicit solventMonte Carlo(MC) simulation results. Calculated thermodynamic properties are also represented consistently among these three methods. Our versatile HNC/MD hybrid method allows for a quick prediction of the thermodynamics of multicomponent electrolytesolutions for a wide range of concentrations and an efficient assessment of the validity of the employed MD forcefields with possible implications in the development of thermodynamically consistent parameter sets.

Simulation of fourwave mixing signals by a perturbative approach: Application to ultrafast twodimensional infrared spectroscopy
View Description Hide DescriptionWe propose an alternative method for the calculation of the phasematched contributions, which are responsible for the thirdorder optical signals measured in fourwave mixing experiments. In particular, we extend the strong field dissipation theory of Meier and Tannor [J. Chem. Phys.111, 3365 (1999)] to the case of a perturbative treatment with respect to the exciting laser fields. Our approach is based on an analytical expression of the thirdorder density matrix and hence it does not require to verify numerically the irrelevance of higher order terms or the calculation of a spatial Fourier transform. In order to illustrate this method, we simulate the experimental signal measured in femtosecond twodimensional infrared (2DIR) vibrational spectroscopy. We consider an intramolecular anharmonic vibrational mode modeled by a Morse potential and coupled to a dissipative bath of harmonic oscillators. We calculate the 2DIR correlation spectrum and we discuss the influence of the population decay on the line shapes. In particular, we compare two situations, one where only pure dephasing processes are considered, and another one where phase losses due to population relaxation are also taken into account. We show that the shape of the peaks observed in a 2DIR correlation spectrum differs in these two cases, and therefore this difference appears as a signature of population decay and gives information on the importance of pure dephasing processes in phase loss mechanisms.

The selfenergy beyond : Local and nonlocal vertex corrections
View Description Hide DescriptionIt is commonly accepted that the approximation for the electron selfenergy is successful for the description of the band structure of weakly to moderately correlated systems, whereas it will fail for strongly correlated materials. In the present work, we discuss two important aspects of this approximation: first, the “selfscreening error,” which is due to an incorrect treatment of induced exchange, and second, the atomic limit, in which, instead, correlation is directly responsible for the observed problem. Using the example of the removal of a particle from a box, we show that the selfscreening error stems from the use of test chargetestcharge screening and that it can be corrected by a twopoint vertex contribution to the selfenergy derived from timedependent density functional theory (TDDFT). We explain why the addition of a particle, instead, requires the use of a different approximate vertex. This illustrates why the general vertex function, valid both for valence and conduction states, must be a threepoint function. Moreover, we show that also the bad performance of in the atomic limit is due to the neglect of the vertex in the selfenergy; in that case, the TDDFTderived vertex correction is not sufficient in order to remove the error even qualitatively. We discuss the effects of the selfscreening error as well as the atomic limit using for the exactly solvable twosite Hubbard model.

Local hybrids as a perturbation to global hybrid functionals
View Description Hide DescriptionWe present new local hybrids of generalized gradient approximation exchange, designed to be small perturbations to the corresponding global hybrid. In general, local hybrids include a positiondependent admixture of nonlocal Hartree–Fock exchange. These new local hybrids incorporate a constant fraction of nonlocal exchange, plus additional nonlocal exchange contributions near nuclei. These functionals predict molecular thermochemistry and reaction barriers on average more accurately than their “parent” global hybrid.

Isochronal sampling in nonBoltzmann Monte Carlo methods
View Description Hide DescriptionNonBoltzmann sampling (NBS) methods are usually able to overcome ergodicity issues which conventional Monte Carlo methods often undergo. In short, NBS methods are meant to broaden the sampling range of some suitable order parameter (e.g., energy). For many years, a standard for their development has been the choice of sampling weights that yield uniform sampling of a predefined parameter range. However, Trebst et al. [Phys. Rev. E70, 046701 (2004)] demonstrated that better results are obtained by choosing weights that reduce as much as possible the average number of steps needed to complete a roundtrip in that range. In the present work, we prove that the method they developed to minimize roundtrip times also equalizes downtrip and uptrip times. Then, we propose a discreteparameter extension using such isochronal character as our main goal. To assess the features of the new method, we carry out simulations of a spin system and of lattice chains designed to exhibit folding transition, thus being suitable models for proteins. Our results show that the new method performs on a par with the original method when the latter is applicable. However, there are cases in which the method of Trebst et al. becomes inapplicable, depending on the chosen order parameter and on the employed Monte Carlo moves. With a practical example, we demonstrate that our method can naturally handle these cases, thus being more robust than the original one. Finally, we find an interesting correspondence between the kind of approach dealt with here and the committor analysis of reaction coordinates, which is another topic of rising interest in the field of molecular simulation.

Nonuniqueness of magnetic fields and energy derivatives in spinpolarized density functional theory
View Description Hide DescriptionThe effect of the recently uncovered nonuniqueness of the external magnetic field corresponding to a given pair of density and spin density on the derivative of the energy functional of spinpolarized density functional theory, and its implications for the definition of chemical reactivity descriptors, is examined. For ground states, the nonuniqueness of implies the nondifferentiability of the energy functional with respect to . It is shown, on the other hand, that this nonuniqueness allows the existence of the onesided derivatives of with respect to . Although the electron ground state can always be obtained from the minimization of without any constraint on the spin number , the Lagrange multiplier associated with the fixation of does not vanish even for ground states. is identified as the left or rightside derivative of the total energy with respect to , which justifies the interpretation of as a (spin) chemical potential. This is relevant not only for the spinpolarized generalization of conceptual density functional theory, the spin chemical potential being one of the elementary reactivity descriptors, but also for the extension of the thermodynamical analogy of density functional theory for the spinpolarized case. For higherorder reactivity indices, ’s nonuniqueness has similar implications as for , leading to a split of the indices with respect to into onesided reactivity descriptors.

Making the random phase approximation to electronic correlation accurate
View Description Hide DescriptionWe show that the inclusion of secondorder screened exchange to the random phase approximation allows for an accurate description of electronic correlation in atoms and solids clearly surpassing the random phase approximation, but not yet approaching chemical accuracy. From a fundamental point of view, the method is selfcorrelation free for oneelectron systems. From a practical point of view, the approach yields correlation energies for atoms, as well as for the jellium electron gas within a few kcal/mol of exact values, atomization energies within typically 2–3 kcal/mol of experiment, and excellent lattice constants for ionic and covalently bonded solids (0.2% error). The computational complexity is only , comparable to canonical secondorder Møller–Plesset perturbation theory, which should allow for routine calculations on many systems.

Resonating valence bond wave function with molecular orbitals: Application to firstrow molecules
View Description Hide DescriptionWe introduce a method for accurate quantum chemical calculations based on a simple variational wave function, defined by a single geminal that couples all the electrons into singlet pairs, combined with a real space correlation factor. The method uses a constrained variational optimization, based on an expansion of the geminal in terms of molecular orbitals. It is shown that the most relevant nondynamical correlations are correctly reproduced once an appropriate number of molecular orbitals is considered. The value of is determined by requiring that, in the atomization limit, the atoms are described by Hartree–Fock Slater determinants with Jastrow correlations. The energetics, as well as other physical and chemical properties, are then given by an efficient variational approach based on standard quantum Monte Carlo techniques. We test this method on a set of homonuclear (, , , , , and ) and heteronuclear (LiF and CN) dimers for which strong nondynamical correlations and/or weak van der Waals interactions are present.

A gradientdirected Monte Carlo method for global optimization in a discrete space: Application to protein sequence design and folding
View Description Hide DescriptionWe apply the gradientdirected Monte Carlo (GDMC) method to select optimal members of a discrete space, the space of chemically viable proteins described by a model Hamiltonian. In contrast to conventional Monte Carlo approaches, our GDMC method uses local property gradients with respect to chemical variables that have discrete values in the actual systems, e.g., residue types in a proteinsequence. The local property gradients are obtained from the interpolation of discrete property values, following the linear combination of atomic potentials scheme developed recently [M. Wang et al., J. Am. Chem. Soc.128, 3228 (2006)]. The local property derivative information directs the search toward the global minima while the Metropolis criterion incorporated in the method overcomes barriers between local minima. Using the simple HP lattice model, we apply the GDMC method to proteinsequence design and folding. The GDMC algorithm proves to be particularly efficient, suggesting that this strategy can be extended to other discrete optimization problems in addition to inverse molecular design.