Volume 139, Issue 8, 28 August 2013
 COMMUNICATIONS


Communication: Analytic gradients in the randomphase approximation
View Description Hide DescriptionThe relationship between the randomphaseapproximation (RPA) correlation energy and the continuous algebraic Riccati equation is examined and the importance of a stabilizing solution is emphasized. The criterion to distinguish this from nonstabilizing solutions can be used to ensure that physical, smooth potential energy surfaces are obtained. An implementation of analytic RPA molecular gradients is presented using the Lagrangian technique. Illustrative calculations indicate that RPA with HartreeFock reference orbitals delivers an accuracy similar to that of secondorder Møller–Plesset perturbation theory.

Communication: Nonradiative recombination via conical intersection at a semiconductor defect
View Description Hide DescriptionLocalization of electronic excitations at moleculesized semiconductor defects often precedes nonradiative (NR) decay, and it is known that many molecules undergo NR decay via conical intersection. Herein, we report the direct simulation of fast and efficient NR decay via a conical intersection at a known semiconductor defect. It is suggested that this silicon epoxide defect may selectively quench photoluminescence (PL) in small silicon nanocrystals (band gap > ∼2.8 eV), and thus influence both the observed PL yield and PL energy of oxidized silicon nanocrystals.

Communication: Constructing an implicit quantum mechanical/molecular mechanics solvent model by coarsegraining explicit solvent
View Description Hide DescriptionTo avoid repeated, computationally expensive QM solute calculations while sampling MM solvent in QM/MM simulations, a new approach for constructing an implicit solvent model by coarsegraining the solvent properties over many explicit solvent configurations is proposed. The solvent is modeled using a polarizable force field that is parameterized in terms of distributed multipoles (electrostatics), polarizabilities (induction), and frequencydependent polarizabilities (dispersion). The coarsegraining procedure exploits the ability to translate these properties to the center of each coarsegraining cell and average them over many solvent configurations before interacting them with the solute. A single coarsegrained QM/MM calculation of the interaction between a formamide solute and aqueous solvent reproduces the much more expensive average over many explicit QM/MM calculations with kJ/mol accuracy.
 Top

 ARTICLES

 Theoretical Methods and Algorithms

ArisTaylor dispersion with drift and diffusion of particles on the tube wall
View Description Hide DescriptionA laminar stationary flow of viscous fluid in a cylindrical tube enhances the rate of diffusion of Brownian particles along the tube axis. This socalled ArisTaylor dispersion is due to the fact that cumulative times, spent by a diffusing particle in layers of the fluid moving with different velocities, are random variables which depend on the realization of the particle stochastic trajectory in the radial direction. Conceptually similar increase of the diffusivity occurs when the particle randomly jumps between two states with different drift velocities. Here we develop a theory that contains both phenomena as special limiting cases. It is assumed (i) that the particle in the flow can reversibly bind to the tube wall, where it moves with a given drift velocity and diffusivity, and (ii) that the radial and longitudinal diffusivities of the particle in the flow may be different. We derive analytical expressions for the effective drift velocity and diffusivity of the particle, which show how these quantities depend on the geometric and kinetic parameters of the model.

A boundary correction algorithm for metadynamics in multiple dimensions
View Description Hide DescriptionMetadynamics is an efficient method for simulation of the free energy of manyparticle systems. Over the last several years it has been applied to study a wide variety of systems, ranging from simple fluids to biological macromolecules. The method relies on uniform sampling along specified collective variables or order parameters. Such order parameters, however, are often bounded, and metadynamics algorithms as originally developed suffer from systematic errors at the corresponding boundaries. While several approaches have been proposed in the past to correct these errors for unidimensional systems, no method exists to fully correct these errors in multidimensional systems at points where multiple boundaries meet. Here we present a correction scheme that circumvents this limitation.

Seniority number in spinadapted spaces and compactness of configuration interaction wave functions
View Description Hide DescriptionThis work extends the concept of seniority number, which has been widely used for classifying Nelectron Slater determinants, to wave functions of N electrons and spin S, as well as to Nelectron spinadapted Hilbert spaces. We propose a spinfree formulation of the seniority number operator and perform a study on the behavior of the expectation values of this operator under transformations of the molecular basis sets. This study leads to propose a quantitative evaluation for the convergence of the expansions of the wave functions in terms of Slater determinants. The noninvariant character of the seniority number operator expectation value of a wave function with respect to a unitary transformation of the molecular orbital basis set, allows us to search for a change of basis which minimizes that expectation value. The results found in the description of wave functions of selected atoms and molecules show that the expansions expressed in these bases exhibit a more rapid convergence than those formulated in the canonical molecular orbital bases and even in the natural orbital ones.

Nonlinear optical response in solids from timedependent densityfunctional theory simulations
View Description Hide DescriptionA new method for computing nonlinear susceptibilities of periodic solids is presented. The electronic response and polarization current are obtained from timedependent Schrödinger equation dynamically coupled to the external electromagnetic field. Solid's polarization resulting from quasimonochromatic excitation is examined in frequency domain. The higher order susceptibilities are calculated nonperturbatively. The method is illustrated by examples of third harmonic generation in silicon and carbon diamond.

A comparative study of methods to compute the free energy of an ordered assembly by molecular simulation
View Description Hide DescriptionWe present a comparative study of methods to compute the absolute free energy of a crystalline assembly of hard particles by molecular simulation. We consider all combinations of three choices defining the methodology: (1) the reference system: Einstein crystal (EC), interacting harmonic (IH), or r −12 soft spheres (SS); (2) the integration path: FrenkelLadd (FL) or penetrable ramp (PR); and (3) the freeenergy method: overlapsampling freeenergy perturbation (OS) or thermodynamic integration (TI). We apply the methods to FCC hard spheres at the melting state. The study shows that, in the best cases, OS and TI are roughly equivalent in efficiency, with a slight advantage to TI. We also examine the multistate Bennett acceptance ratio method, and find that it offers no advantage for this particular application. The PR path shows advantage in general over FL, providing results of the same precision with 2–9 times less computation, depending on the choice of a common reference. The best combination for the FL path is TI+EC, which is how the FL method is usually implemented. For the PR path, the SS system (with either TI or OS) proves to be most effective; it gives equivalent precision to TI+FL+EC with about 6 times less computation (or 12 times less, if discounting the computational effort required to establish the SS reference free energy). Both the SS and IH references show great advantage in capturing finitesize effects, providing a variation in freeenergy difference with system size that is about 10 times less than EC. This result further confirms previous work for softparticle crystals, and suggests that freeenergy calculations for a structured assembly be performed using a hybrid method, in which the finitesystem freeenergy difference is added to the extrapolated (1/N→0) absolute free energy of the reference system, to obtain a result that is nearly independent of system size.

Quantum mechanical/molecular mechanical/continuum style solvation model: Timedependent density functional theory
View Description Hide DescriptionA combined quantum mechanical/molecular mechanical/continuum (QM/MMpol/C) style method is developed for timedependent density functional theory (TDDFT, including longrange corrected TDDFT) method, induced dipole polarizable force field, and induced surface charge continuum model. Induced dipoles and induced charges are included in the TDDFT equations to solve for the transition energies, relaxed density, and transition density. Analytic gradient is derived and implemented for geometry optimization and molecular dynamics simulation. QM/MMpol/C style DFT and TDDFT methods are used to study the hydrogen bonding of the photoactive yellow protein chromopore in ground state and excited state.

The time correlation function perspective of NMR relaxation in proteins
View Description Hide DescriptionWe applied over a decade ago the twobody coupledrotator slowly relaxing local structure (SRLS) approach to NMR relaxation in proteins. One rotator is the globally moving protein and the other rotator is the locally moving probe (spinbearing moiety, typically the 15N−1H bond). So far we applied SRLS to 15N−H relaxation from seven different proteins within the scope of the commonly used datafitting paradigm. Here, we solve the SRLS Smoluchowski equation using typical bestfit parameters as input, to obtain the corresponding generic time correlation functions (TCFs). The following new information is obtained. For actual rhombic local ordering and main ordering axis pointing along , the measurable TCF is dominated by the (K,K′) = (−2,2), (2,2), and (0,2) components (K is the order of the rank 2 local ordering tensor), determined largely by the local motion. Global diffusion axiality affects the analysis significantly when the ratio between the parallel and perpendicular components exceeds approximately 1.5. Local diffusion axiality has a large and intricate effect on the analysis. Modecoupling becomes important when the ratio between the global and local motional rates falls below 0.01. The traditional method of analysis − modelfree (MF) − represents a simple limit of SRLS. The conditions under which the MF and SRLS TCFs are the same are specified. The validity ranges of wobbleinacone and rotation on the surface of a cone as local motions are determined. The evolution of the intricate Smoluchowski operator from the simple diffusion operator for a sphere reorienting in isotropic medium is delineated. This highlights the fact that SRLS is an extension of the established stochastic theories for treating restricted motions. This study lays the groundwork for TCFbased comparison between mesoscopic SRLS and atomistic molecular dynamics.

Mimicking photoisomerisation of azomaterials by a force field switch derived from nonadiabatic ab initio simulations: Application to photoswitchable helical foldamers in solution
View Description Hide DescriptionA force field to induce isomerisation of photoswitchable azobenzene groups embedded in molecular materials has been developed in the framework of force field molecular dynamics simulations. A molecular mechanics switching potential has been tuned so as to reproduce both the correct photoisomerisation timescale and mechanism that has been generated by reference nonadiabatic ab initio molecular dynamics. As a first application, we present a force field molecular dynamics study of a prototype photoswitchable foldamer in acetonitrile as solvent. Our analyses reveal that the photoisomerisation of the azobenzene unit embedded in the foldamer occurs via the socalled NNtwist mechanism, and that there exist several distinct unfolding channels for the helix that could be exploited in novel applications of photoresponsive materials.

Assessment of mesoscopic particlebased methods in microfluidic geometries
View Description Hide DescriptionWe assess the accuracy and efficiency of two particlebased mesoscopic simulation methods, namely, Dissipative Particle Dynamics (DPD) and Stochastic Rotation Dynamics (SRD) for predicting a complex flow in a microfluidic geometry. Since both DPD and SRD use soft or weakly interacting particles to carry momentum, both methods contain unavoidable inertial effects and unphysically high fluid compressibility. To assess these effects, we compare the predictions of DPD and SRD for both an exact Stokesflow solution and nearly exact solutions at finite Reynolds numbers from the finite element method for flow in a straight channel with periodic slip boundary conditions. This flow represents a periodic electroosmotic flow, which is a complex flow with an analytical solution for zero Reynolds number. We find that SRD is roughly tenfold faster than DPD in predicting the flow field, with better accuracy at low Reynolds numbers. However, SRD has more severe problems with compressibility effects than does DPD, which limits the Reynolds numbers attainable in SRD to around 25–50, while DPD can achieve Re higher than this before compressibility effects become too large. However, since the SRD method runs much faster than DPD does, we can afford to enlarge the number of grid cells in SRD to reduce the fluid compressibility at high Reynolds number. Our simulations provide a method to estimate the range of conditions for which SRD or DPD is preferable for mesoscopic simulations.

An expanded calibration study of the explicitly correlated CCSD(T)F12b method using large basis set standard CCSD(T) atomization energies
View Description Hide DescriptionThe effectiveness of the recently developed, explicitly correlated coupled cluster method CCSD(T)F12b is examined in terms of its ability to reproduce atomization energies derived from complete basis set extrapolations of standard CCSD(T). Most of the standard method findings were obtained with augccpV7Z or augccpV8Z basis sets. For a few homonuclear diatomic molecules it was possible to push the basis set to the augccpV9Z level. F12b calculations were performed with the ccpVnZF12 (n = D, T, Q) basis set sequence and were also extrapolated to the basis set limit using a Schwenkestyle, parameterized formula. A systematic bias was observed in the F12b method with the (VTZF12/VQZF12) basis set combination. This bias resulted in the underestimation of reference values associated with small molecules (valence correlation energies <0.5 Eh) and an even larger overestimation of atomization energies for bigger systems. Consequently, caution should be exercised in the use of F12b for high accuracy studies. Root mean square and mean absolute deviation error metrics for this basis set combination were comparable to complete basis set values obtained with standard CCSD(T) and the augccpVDZ through augccpVQZ basis set sequence. However, the mean signed deviation was an order of magnitude larger. Problems partially due to basis set superposition error were identified with second row compounds which resulted in a weak performance for the smaller VDZF12/VTZF12 combination of basis sets.

Local CC2 response method based on the Laplace transform: Orbitalrelaxed firstorder properties for excited states
View Description Hide DescriptionA multistate local CC2 response method for the calculation of orbitalrelaxed first order properties is presented for ground and electronically excited states. It enables the treatment of excited state properties including orbital relaxation for extended molecular systems and is a major step on the way towards analytic gradients with respect to nuclear displacements. The Laplace transform method is employed to partition the eigenvalue problem and the lambda equations, i.e., the doubles parts of these equations are inverted onthefly, leaving only the corresponding effective singles equations to be solved iteratively. Furthermore, the state specific local approximations are adaptive. Densityfitting is utilized to decompose the electronrepulsion integrals. The accuracy of the local approximation is tested and the efficiency of the new code is demonstrated on the example of an organic sensitizer for solarcell applications, which consists of about 100 atoms.

Explicitly correlated plane waves: Accelerating convergence in periodic wavefunction expansions
View Description Hide DescriptionWe present an investigation into the use of an explicitly correlated plane wave basis for periodic wavefunction expansions at the level of secondorder MøllerPlesset (MP2) perturbation theory. The convergence of the electronic correlation energy with respect to the oneelectron basis set is investigated and compared to conventional MP2 theory in a finite homogeneous electron gas model. In addition to the widely used Slatertype geminal correlation factor, we also derive and investigate a novel correlation factor that we term YukawaCoulomb. The YukawaCoulomb correlation factor is motivated by analytic results for two electrons in a box and allows for a further improved convergence of the correlation energies with respect to the employed basis set. We find the combination of the infinitely delocalized plane waves and local shortranged geminals provides a complementary, and rapidly convergent basis for the description of periodic wavefunctions. We hope that this approach will expand the scope of discrete wavefunction expansions in periodic systems.

Efficient selfconsistent treatment of electron correlation within the random phase approximation
View Description Hide DescriptionA selfconsistent KohnSham (KS) method is presented that treats correlation on the basis of the adiabaticconnection dissipationfluctuation theorem employing the direct random phase approximation (dRPA), i.e., taking into account only the Coulomb kernel while neglecting the exchangecorrelation kernel in the calculation of the KohnSham correlation energy and potential. The method, denoted selfconsistent dRPA method, furthermore treats exactly the exchange energy and the local multiplicative KS exchange potential. It uses Gaussian basis sets, is reasonably efficient, exhibiting a scaling of the computational effort with the forth power of the system size, and thus is generally applicable to molecules. The resulting dRPA correlation potentials in contrast to common approximate correlation potentials are in good agreement with exact reference potentials. The negatives of the eigenvalues of the highest occupied molecular orbitals are found to be in good agreement with experimental ionization potentials. Total energies from selfconsistent dRPA calculations, as expected, are even poorer than nonselfconsistent dRPA total energies and dRPA reaction and noncovalent binding energies do not significantly benefit from selfconsistency. On the other hand, energies obtained with a recently introduced adiabaticconnection dissipationfluctuation approach (EXXRPA+, exactexchange random phase approximation) that takes into account, besides the Coulomb kernel, also the exact frequencydependent exchange kernel are significantly improved if evaluated with orbitals obtained from a selfconsistent dRPA calculation instead of an exact exchangeonly calculation. Total energies, reaction energies, and noncovalent binding energies obtained in this way are of the same quality as those of highlevel quantum chemistry methods, like the coupled cluster singles doubles method which is computationally more demanding.

A pair natural orbital implementation of the coupled cluster model CC2 for excitation energies
View Description Hide DescriptionWe demonstrate how to extend the pair natural orbital (PNO) methodology for excited states, presented in a previous work for the perturbative doubles correction to configuration interaction singles (CIS(D)), to iterative coupled cluster methods such as the approximate singles and doubles model CC2. The original scaling of the PNO construction is reduced by using orbitalspecific virtuals (OSVs) as an intermediate step without spoiling the initial accuracy of the PNO method. Furthermore, a slower error convergence for chargetransfer states is analyzed and resolved by a numerical Laplace transformation during the PNO construction, so that an equally accurate treatment of local and chargetransfer excitations is achieved. With statespecific truncated PNO expansions, the eigenvalue problem is solved by combining the Davidson algorithm with deflation to project out roots that have already been determined and an automated refresh with a generation of new PNOs to achieve selfconsistency of the PNO space. For a large test set, we found that truncation errors for PNOCC2 excitation energies are only slightly larger than for PNOCIS(D). The computational efficiency of PNOCC2 is demonstrated for a large organic dye, where a reduction of the doubles space by a factor of more than 1000 is obtained compared to the canonical calculation. A compression of the doubles space by a factor 30 is achieved by a unified OSV space only. Moreover, calculations with the still preliminary PNOCC2 implementation on a series of glycine oligomers revealed an early break even point with a canonical RICC2 implementation between 100 and 300 basis functions.

Derivation of a true (t → 0_{+}) quantum transitionstate theory. II. Recovery of the exact quantum rate in the absence of recrossing
View Description Hide DescriptionIn Paper I [T. J. H. Hele and S. C. Althorpe, J. Chem. Phys.138, 084108 (Year: 2013)]10.1063/1.4792697 we derived a quantum transitionstate theory (TST) by taking the t → 0+ limit of a new form of quantum fluxside timecorrelation function containing a ringpolymer dividing surface. This t → 0+ limit appears to be unique in giving positivedefinite Boltzmann statistics, and is identical to ringpolymer molecular dynamics (RPMD) TST. Here, we show that quantum TST (i.e., RPMDTST) is exact if there is no recrossing (by the realtime quantum dynamics) of the ringpolymer dividing surface, nor of any surface orthogonal to it in the space describing fluctuations in the polymerbead positions along the reaction coordinate. In practice, this means that RPMDTST gives a good approximation to the exact quantum rate for direct reactions, provided the temperature is not too far below the crossover to deep tunnelling. We derive these results by comparing the t → ∞ limit of the ringpolymer fluxside timecorrelation function with that of a hybrid fluxside timecorrelation function (containing a ringpolymer flux operator and a MillerSchwarzTromp side function), and by representing the resulting ringpolymer momentum integrals as hypercubes. Together with Paper I, the results of this article validate a large number of RPMD calculations of reaction rates.

On the uniqueness of t → 0_{+} quantum transitionstate theory
View Description Hide DescriptionIt was shown recently that there exists a true quantum transitionstate theory (QTST) corresponding to the t → 0+ limit of a (new form of) quantum fluxside timecorrelation function. Remarkably, this QTST is identical to ringpolymer molecular dynamics (RPMD) TST. Here, we provide evidence which suggests very strongly that this QTST (≡ RPMDTST) is unique, in the sense that the t → 0+ limit of any other fluxside timecorrelation function gives either nonpositivedefinite quantum statistics or zero. We introduce a generalized fluxside timecorrelation function which includes all other (known) fluxside timecorrelation functions as special limiting cases. We find that the only nonzero t → 0+ limit of this function that contains positivedefinite quantum statistics is RPMDTST.

Simulations of nanocrystals under pressure: Combining electronic enthalpy andlinearscaling densityfunctional theory
View Description Hide DescriptionWe present an implementation in a linearscaling densityfunctional theory code of an electronic enthalpy method, which has been found to benatural and efficient for the ab initio calculation of finite systems underhydrostatic pressure. Based on a definition of the system volume as that enclosed within anelectronic density isosurface [M. Cococcioni, F. Mauri,G. Ceder, and N. Marzari, Phys. Rev. Lett.94, 145501 (2005)], it supports bothgeometry optimizations and molecular dynamics simulations. We introduce an approach for calibratingthe parameters defining the volume in the context of geometry optimizations and discuss theirsignificance. Results in good agreement with simulations using explicit solvents are obtained, validating ourapproach. Sizedependent pressureinduced structural transformations andvariations in the energy gap of hydrogenated silicon nanocrystals areinvestigated, including one comparable in size to recent experiments. A detailed analysis of thepolyamorphic transformationsreveals three types of amorphous structures and theirpersistence on depressurization is assessed.