Volume 136, Issue 14, 14 April 2012

Competition between π⋯πinteraction and halogen bond in solution has been investigated by using carbon nuclear magnetic resonance spectroscopy (^{13}C NMR) combined with density functional theory calculation. Both experimental and theoretical results clearly show that there are no C–Cl⋯π or C–Br⋯π halogen bonds and only the π⋯πinteractions exist in the binary liquid mixtures of C_{6}D_{6} with C_{6}F_{5}Cl and C_{6}F_{5}Br, respectively. The case is totally different for the binary liquid mixtures of C_{6}D_{6} with C_{6}F_{5}I in which the C–I⋯π halogen bonds not the π⋯πinteractions are present. The important role of entropy in the competition between π⋯πinteraction and halogen bond in solution was also discussed.
 ARTICLES

 Theoretical Methods and Algorithms

On the evaluation of the noninteracting kinetic energy in density functional theory
View Description Hide DescriptionThe utility of both an orbitalfree and a singleorbital expression for computing the noninteracting kinetic energy in density functional theory is investigated for simple atomic systems. The accuracy of both expressions is governed by the extent to which the Kohn–Sham equation is solved for the given exchange–correlation functional and so special attention is paid to the influence of finite Gaussian basis sets. The orbitalfree expression is a statement of the virial theorem and its accuracy is quantified. The accuracy of the singleorbital expression is sensitive to the choice of Kohn–Sham orbital. The use of particularly compact orbitals is problematic because the failure to solve the Kohn–Sham equation exactly in regions where the orbital has decayed to nearzero leads to unphysical behaviour in regions that contribute to the kinetic energy, rendering it inaccurate. This problem is particularly severe for core orbitals, which would otherwise appear attractive due to their formally nodeless nature. The most accurate results from the singleorbital expression are obtained using the relatively diffuse, highest occupied orbitals, although special care is required at orbital nodes.

Theory of binless multistate free energy estimation with applications to proteinligand binding
View Description Hide DescriptionThe weighted histogram analysis method (WHAM) is routinely used for computing free energies and expectations from multiple ensembles. Existing derivations of WHAM require observations to be discretized into a finite number of bins. Yet, WHAM formulas seem to hold even if the bin sizes are made arbitrarily small. The purpose of this article is to demonstrate both the validity and value of the multistate Bennet acceptance ratio (MBAR) method seen as a binless extension of WHAM. We discuss two statistical arguments to derive the MBAR equations, in parallel to the selfconsistency and maximum likelihood derivations already known for WHAM. We show that the binless method, like WHAM, can be used not only to estimate free energies and equilibrium expectations, but also to estimate equilibrium distributions. We also provide a number of useful results from the statistical literature, including the determination of MBAR estimators by minimization of a convex function. This leads to an approach to the computation of MBAR free energies by optimization algorithms, which can be more effective than existing algorithms. The advantages of MBAR are illustrated numerically for the calculation of absolute proteinligand binding free energies by alchemical transformations with and without softcore potentials. We show that binless statistical analysis can accurately treat sparsely distributed interactionenergy samples as obtained from unmodified interaction potentials that cannot be properly analyzed using standard binning methods. This suggests that binless multistate analysis of binding free energy simulations with unmodified potentials offers a straightforward alternative to the use of softcore potentials for these alchemical transformations.

Nonadiabatic Ehrenfest molecular dynamics within the projector augmentedwave method
View Description Hide DescriptionWe derive equations for nonadiabatic Ehrenfest molecular dynamics within the projector augmentedwave (PAW) formalism. The discretization of the electrons is timedependent as the augmentation functions depend on the positions of the nuclei. We describe the implementation of the Ehrenfest molecular dynamicsequations within the realspace PAW method. We demonstrate the applicability of our method by studying the vibration of NaCl, the torsional rotation of H_{2}C= in both the adiabatic and the nonadiabatic regimes, and the hydrogen bombardment of C_{40}H_{16}.

Combining activespace coupledcluster methods with moment energy corrections via the CC(P;Q) methodology, with benchmark calculations for biradical transition states
View Description Hide DescriptionWe have recently suggested the CC(P;Q) methodology that can correct energies obtained in the activespace coupledcluster (CC) or equationofmotion (EOM) CC calculations, which recover much of the nondynamical and some dynamical electron correlation effects, for the higherorder, mostly dynamical, correlations missing in the activespace CC/EOMCC considerations. It is shown that one can greatly improve the description of biradical transition states, both in terms of the resulting energy barriers and total energies, by combining the CC approach with singles, doubles, and activespace triples, termed CCSDt, with the CC(P;Q)style correction due to missing triple excitations defining the CC(t;3) approximation.

The orbitalspecificvirtual local coupled cluster singles and doubles method
View Description Hide DescriptionWe extend the orbitalspecificvirtual tensor factorization, introduced for local MøllerPlesset perturbation theory in Ref. [J. Yang, Y. Kurashige, F. R. Manby and G. K. L. Chan, J. Chem. Phys.134, 044123 (2011)10.1063/1.3528935], to local coupled cluster singles and doubles theory (OSVLCCSD). The method is implemented by modifying an efficient projectedatomicorbital local coupled cluster program (PAOLCCSD) described recently, [H.J. Werner and M. Schütz, J. Chem. Phys.135, 144116 (2011)10.1063/1.3641642]. By comparison of both methods we find that the compact representation of the amplitudes in the OSV approach affords various advantages, including smaller computational time requirements (for comparable accuracy), as well as a more systematic control of the error through a single energy threshold. Overall, the OSVLCCSD approach together with an MP2 correction yields small domain errors in practical calculations. The applicability of the OSVLCCSD is demonstrated for molecules with up to 73 atoms and realistic basis sets (up to 2334 basis functions).

Analytical solution for the depolarization of hyperpolarized nuclei by chemical exchange saturation transfer between free and encapsulated xenon (HyperCEST)
View Description Hide DescriptionWe present an analytical solution of the Bloch–McConnell equations for the case of chemical exchange saturation transfer between hyperpolarized nuclei in cavities and in solvent (HyperCEST experiment). This allows quantitative investigation of host–guest interactions by means of nuclear magnetic resonance spectroscopy and, due to the strong HyperCEST signal enhancement, even NMR imaging. Hosts of interest can be hydrophobic cavities in macromolecules or artificial cages like cryptophaneA which was proposed as a targeted biosensor. Relevant system parameters as exchange rate and host concentration can be obtained from the monoexponential depolarization process which is shown to be governed by the smallest eigenvalue in modulus. For this dominant eigenvalue we present a useful approximation leading to the depolarization rate for the case of on and offresonant irradiation. It is shown that this rate is a generalization of the longitudinal relaxation rate in the rotating frame. We demonstrate for the free and cryptophaneAencapsulated xenon system, by comparison with numerical simulations, that HyperCEST experiments are precisely described in the valid range of this widely applicable analytical approximation. Altogether, the proposed analytical solution allows optimization and quantitative analysis of HyperCEST experiments but also characterization and optimal design of possible biosensors.

Distancedependent Schwarzbased integral estimates for twoelectron integrals: Reliable tightness vs. rigorous upper bounds
View Description Hide DescriptionA new integral estimate for fourcenter twoelectron integrals is introduced that accounts for distance information between the bra and ketcharge distributions describing the two electrons. The screening is denoted as QQR and combines the most important features of the conventional Schwarz screening by Häser and Ahlrichs published in 1989 [J. Comput. Chem.10, 104 (1989)10.1002/jcc.540100111] and our multipolebased integral estimates (MBIE) introduced in 2005 [D. S. Lambrecht and C. Ochsenfeld, J. Chem. Phys.123, 184101 (2005)10.1063/1.2079967]. At the same time the estimates are not only tighter but also much easier to implement, so that we recommend them instead of our MBIE bounds introduced first for accounting for chargedistance information. The inclusion of distance dependence between charge distributions is not only useful at the SCF level but is particularly important for describing electroncorrelation effects, e.g., within AOMP2 theory, where the decay behavior is at least 1/R ^{4} or even 1/R ^{6}. In our present work, we focus on studying the efficiency of our QQR estimates within SCF theory and demonstrate the performance for a benchmark set of 44 medium to large molecules, where savings of up to a factor of 2 for exchange integrals are observed for larger systems. Based on the results of the benchmark set we show that reliable tightness of integral estimates is more important for the screening performance than rigorous upper bound properties.

Improved delayleaping simulation algorithm for biochemical reaction systems with delays
View Description Hide DescriptionIn biochemical reaction systems dominated by delays, the simulation speed of the stochastic simulation algorithm depends on the size of the wait queue. As a result, it is important to control the size of the wait queue to improve the efficiency of the simulation. An improved accelerated delay stochastic simulation algorithm for biochemical reaction systems with delays, termed the improved delayleaping algorithm, is proposed in this paper. The update method for the wait queue is effective in reducing the size of the queue as well as shortening the storage and access time, thereby accelerating the simulation speed. Numerical simulation on two examples indicates that this method not only obtains a more significant efficiency compared with the existing methods, but also can be widely applied in biochemical reaction systems with delays.

Incorporating a completely renormalized coupled cluster approach into a composite method for thermodynamic properties and reaction paths
View Description Hide DescriptionThe correlation consistent composite approach (ccCA), using the S4 complete basis set twopoint extrapolation scheme (ccCAS4), has been modified to incorporate the lefteigenstate completely renormalized coupled cluster method, including singles, doubles, and noniterative triples (CRCC(2,3)) as the highest level component. The new ccCACC(2,3) method predicts thermodynamic properties with an accuracy that is similar to that of the original ccCAS4 method. At the same time, the inclusion of the singlereference CRCC(2,3) approach provides a ccCA scheme that can correctly treat reaction pathways that contain certain classes of multireference species such as diradicals, which would normally need to be treated by more computationally demanding multireference methods. The new ccCACC(2,3) method produces a mean absolute deviation of 1.7 kcal/mol for predicted heats of formation at 298 K, based on calibration with the G2/97 set of 148 molecules, which is comparable to that of 1.0 kcal/mol obtained using the ccCAS4 method, while significantly improving the performance of the ccCAS4 approach in calculations involving more demanding radical and diradical species. Both the ccCACC(2,3) and ccCAS4 composite methods are used to characterize the conrotatory and disrotatory isomerization pathways of bicyclo[1.1.0]butane to trans1,3butadiene, for which conventional coupled cluster methods, such as the CCSD(T) approach used in the ccCAS4 model and, in consequence, the ccCAS4 method itself might fail by incorrectly placing the disrotatory pathway below the conrotatory one. The ccCACC(2,3) scheme provides correct pathway ordering while providing an accurate description of the activation and reactionenergies characterizing the lowestenergy conrotatory pathway. The ccCACC(2,3) method is thus a viable method for the analyses of reaction mechanisms that have significant multireference character, and presents a generally less computationally intensive alternative to true multireference methods, with computer costs and ease of use that are similar to those that characterize the more established, CCSD(T)based, ccCAS4 methodology.

Analytical evaluation of Fukui functions and realspace linear response function
View Description Hide DescriptionMany useful concepts developed within density functional theory provide much insight for the understanding and prediction of chemical reactivity, one of the main aims in the field of conceptual density functional theory. While approximate evaluations of such concepts exist, the analytical and efficient evaluation is, however, challenging, because such concepts are usually expressed in terms of functional derivatives with respect to the electron density, or partial derivatives with respect to the number of electrons, complicating the connection to the computational variables of the KohnSham oneelectron orbitals. Only recently, the analytical expressions for the chemical potential, one of the key concepts, have been derived by Cohen, MoriSánchez, and Yang, based on the potential functionaltheory formalism. In the present work, we obtain the analytical expressions for the realspace linear response function using the coupled perturbed KohnSham and generalized KohnSham equations, and the Fukui functions using the previous analytical expressions for chemical potentials of Cohen, MoriSánchez, and Yang. The analytical expressions are exact within the given exchangecorrelation functional. They are applicable to all commonly used approximate functionals, such as local density approximation (LDA), generalized gradient approximation (GGA), and hybrid functionals. The analytical expressions obtained here for Fukui function and linear response functions, along with that for the chemical potential by Cohen, MoriSánchez, and Yang, provide the rigorous and efficient evaluation of the key quantities in conceptual density functional theory within the computational framework of the KohnSham and generalized KohnSham approaches. Furthermore, the obtained analytical expressions for Fukui functions, in conjunction with the linearity condition of the ground state energy as a function of the fractional charges, also lead to new local conditions on the exact functionals, expressed in terms of the secondorder functional derivatives. We implemented the expressions and demonstrate the efficacy with some atomic and molecular calculations, highlighting the importance of relaxation effects.

Rejectionfree Monte Carlo scheme for anisotropic particles
View Description Hide DescriptionWe extend the geometric cluster algorithm [J. Liu and E. Luijten, Phys. Rev. Lett.92, 035504 (2004)], a highly efficient, rejectionfree Monte Carlo scheme for fluids and colloidal suspensions, to the case of anisotropic particles. This is made possible by adopting hyperspherical boundary conditions. A detailed derivation of the algorithm is presented, along with extensive implementation details as well as benchmark results. We describe how the quaternion notation is particularly suitable for the fourdimensional geometric operations employed in the algorithm. We present results for asymmetric LennardJones dimers and for the Yukawa onecomponent plasma in hyperspherical geometry. The efficiency gain that can be achieved compared to conventional, Metropolistype Monte Carlo simulations is investigated for rod–sphere mixtures as a function of rod aspect ratio, rod–sphere diameter ratio, and rod concentration. The effect of curved geometry on physical properties is addressed.

The splitting of atomic orbitals with a common principal quantum number revisited: np vs. ns
View Description Hide DescriptionAtomic orbitals with a common principal quantum number are degenerate, as in the hydrogen atom, in the absence of interelectronic repulsion. Due to the virial theorem, electrons in such orbitals experience equal nuclear attractions. Comparing states of severalelectron atoms that differ by the occupation of orbitals with a common principal quantum number, such as 1s ^{2} 2s vs. 1s ^{2} 2p, we find that although the difference in energies, ΔE, is due to the interelectronic repulsion term in the Hamiltonian, the difference between the interelectronic repulsions, ΔC, makes a smaller contribution to ΔE than the corresponding difference between the nuclear attractions, ΔL. Analysis of spectroscopic data for atomic isoelectronic sequences allows an extensive investigation of these issues. In the low nuclear charge range of pertinent isoelectronic sequences, i.e., for neutral atoms and mildly positively charged ions, it is found that ΔC actually reverses its sign. About 96% of the nuclear attraction difference between the 6p ^{2} P and the 6s ^{2} S states of the Cs atom is cancelled by the corresponding interelectronic repulsion difference. From the monotonic increase of ΔE with Z it follows (via the HellmannFeynman theorem) that ΔL > 0. Upon increasing the nuclear charge along an atomic isoelectronic sequence with a single electron outside a closed shell from Z _{ c }, the critical charge below which the outmost electron is not bound, to infinity, the ratio increases monotonically from to . These results should allow for a more nuanced discussion than is usually encountered of the crude electronic structure of manyelectron atoms and the structure of the periodic table.

Coupledmonomers in molecular assemblies: Theory and application to the water tetramer, pentamer, and ring hexamer
View Description Hide DescriptionWe present extensions to the localmonomer (LMon) Model, a general quantum method to describe coupled intramolecular vibrational modes of a molecular cluster consisting of a set of monomers[Y. Wang and J. M. Bowman, J. Chem. Phys.134, 154510 (2011)]10.1063/1.3579995, to incorporate monomermonomer coupling. A central aspect of the LMon model is a local normalmode analysis, done for each monomer, perturbed by all other mononers. Monomermonomer coupling is described by several approaches based on these normalmode analyses. Two are Hückeltype models, where coupling constants for each intramolecular mode are determined nonempirically from normalmode analyses. One model, the simple one, is limited to nearestneighbor interactions. The second and more general one determines monomermonomer couplings from the full and localmonomer Hessians, with no further assumptions. The simple approach is applied to the water tetramer, pentamer and ring hexamer. For the tetramer and ring hexamer cases, artificial degeneracies of the intramolecular energies in the LMon model, owing to the high symmetry of the cluster, are correctly lifted. The general approach to obtain coupling constants is illustrated for the ring hexamer, where new fundamental energies are reported. Other, more rigorous approaches are suggested but not implemented.

Optimal coherent control of coherent antiStokes Raman scattering: Signal enhancement and background elimination
View Description Hide DescriptionThe ability to enhance resonant signals and eliminate the nonresonant background is analyzed for coherent antiStokes Raman scattering(CARS). The analysis is done at a specific frequency as well as for broadband excitation using femtosecond pulseshaping techniques. An appropriate objective functional is employed to balance resonant signal enhancement against nonresonant background suppression. Optimal enhancement of the signal and minimization of the background can be achieved by shaping the probe pulse alone while keeping the pump and Stokes pulses unshaped. In some cases analytical forms for the probe pulse can be found, and numerical simulations are carried out for other circumstances. It is found that a good approximate optimal solution for resonant signal enhancement in twopulse CARS is a superposition of linear and arctangenttype phases for the pump. The wellknown probe delay method is shown to be a quasioptimal scheme for broadband background suppression. The results should provide a basis to improve the performance of CARSspectroscopy and microscopy.

Improved constraint satisfaction in a simple generalized gradient approximation exchange functional
View Description Hide DescriptionThough there is fevered effort on orbitaldependent approximate exchangecorrelation functionals, generalized gradient approximations, especially the PerdewBurkeErnzerhof (PBE) form, remain the overwhelming choice in calculations. A simple generalized gradient approximation (GGA) exchange functional [A. Vela, V. Medel, and S. B. Trickey, J. Chem. Phys.130, 244103 (2009)] was developed that improves substantially over PBE in energetics (on a typical test set) while being almost as simple in form. The improvement came from constraining the exchange enhancement factor to be below the LiebOxford bound for all but one value of the exchange dimensionless gradient, s, and to go to the uniform electron gas limit at both s = 0 and s → ∞. Here we discuss the issue of asymptotic constraints for GGAs and show that imposition of the large s constraint, , where F xc (n, s) is the enhancement factor and n is the electron density, upon the VelaMedelTrickey (VMT) exchange functional yields modest further improvement. The resulting exchange functional, denoted VT{8,4}, is only slightly more complicated than VMT and easy to program. Additional improvement is obtained by combining VT{8,4} or VMT exchange with the LeeYangParr correlation functional. Extensive computational results on several datasets are provided as verification of the overall performance gains of both versions.

A comparison of methods for melting point calculation using molecular dynamics simulations
View Description Hide DescriptionAccurate and efficient prediction of melting points for complex molecules is still a challenging task for molecular simulation, although many methods have been developed. Four melting point computational methods, including one free energybased method (the pseudosupercritical path (PSCP) method) and three direct methods (two interfacebased methods and the voids method) were applied to argon and a widely studied ionic liquid 1nbutyl3methylimidazolium chloride ([BMIM][Cl]). The performance of each method was compared systematically. All the methods under study reproduce the argon experimental melting point with reasonable accuracy. For [BMIM][Cl], the melting point was computed to be 320 K using a revised PSCP procedure, which agrees with the experimental value 337–339 K very well. However, large errors were observed in the computed results using the direct methods, suggesting that these methods are inappropriate for large molecules with sluggish dynamics. The strengths and weaknesses of each method are discussed.

Relativistic explicit correlation: Coalescence conditions and practical suggestions
View Description Hide DescriptionTo set up the general framework for relativistic explicitly correlated wave function methods, the electronelectron coalescence conditions are derived for the wave functions of the DiracCoulomb (DC), DiracCoulombGaunt (DCG), DiracCoulombBreit (DCB), modified DiracCoulomb (MDC), and zerothorder regularly approximated (ZORA) Hamiltonians. The manipulations make full use of the internal symmetries of the reduced twoelectron Hamiltonians such that the asymptotic behaviors of the wave functions emerge naturally. The results show that, at the coalescence point of two electrons, the wave functions of the DCG Hamiltonian are regular, while those of the DC and DCB Hamiltonians have weak singularities of the type with ν being negative and of . The behaviors of the MDC wave functions are related to the original ones in a simple manner, while the spinfree counterparts are somewhat different due to the complicated electronelectron interaction. The behaviors of the ZORA wave functions depend on the chosen potential in the kinetic energy operator. In the case of the nuclear attraction, the behaviors of the ZORA wave functions are very similar to those of the nonrelativistic ones, just with an additional correction of to the nonrelativistic cusp condition. However, if the Coulomb interaction is also included, the ZORA wave functions become close to the largelarge components of the DC wave functions. Note that such asymptotic expansions of the relativistic wave functions are only valid within an extremely small convergence radius R _{ c } of . Beyond this radius, the behaviors of the relativistic wave functions are still dominated by the nonrelativistic limit, as can be seen in terms of direct perturbation theory (DPT) of relativity. However, as the two limits α → 0 and r _{12} → 0 do not commute, DPT is doomed to fail due to incorrect descriptions of the smallsmall component Ψ^{ SS } of the DC wave function for r _{12} < R _{ c }. Another deduction from the possible divergence of Ψ^{ SS } at r _{12} = R _{ c } is that the DC Hamiltonian has no bound electronic states, although the last word cannot be said. These findings enrich our understandings of relativistic wave functions. On the practical side, it is shown that, under the nopair approximation, relativistic explicitly correlated wave function methods can be made completely parallel to the nonrelativistic counterparts, as demonstrated explicitly for MP2F12. Yet, this can only be achieved by using an extended nopair projector.
 Gas Phase Dynamics and Structure: Spectroscopy, Molecular Interactions, Scattering, and Photochemistry

Quantum Monte Carlo for the xray absorption spectrum of pyrrole at the nitrogen Kedge
View Description Hide DescriptionFixednode diffusion Monte Carlo (FNDMC) is used to simulate the xray absorptionspectrum of a gasphase pyrrole molecule at the nitrogen Kedge. Trial wave functions for coreexcited states are constructed from groundstate KohnSham determinants substituted with singly occupied natural orbitals from configuration interaction with single excitations calculations of the five lowest valenceexcited triplet states. The FNDMC ionization potential (IP) is found to lie within 0.3 eV of the experimental value of 406.1 ± 0.1 eV. The transition energies to antibonding virtual orbitals match the experimental spectrum after alignment of IP values and agree with the existing assignments.

Timedependent wave packet theory for statetostate differential cross sections of fouratom reactions in full dimensions: Application to the HD + OH → H_{2}O + D reaction
View Description Hide DescriptionTimedependent wave packet method has been developed to calculate differential cross section for fouratom reactions in full dimension, utilizing an improved version of reactantproductdecoupling scheme. Differential cross sections for the title reaction were calculated for collision energy up to 0.4 eV. It is found that the differential cross sections for the reaction are all peaked in the backward direction. The majority of H_{2}O is produced in the first stretch excited state, with a large fraction of available energy for the reaction going into H_{2}O internal motion. As compared in a previous report by Xiao et al. [Science333, 440 (2011)]10.1126/science.1205770, the differential cross section at E_{ c } = 0.3 eV and the differential cross section at the backward direction as a function of collision energy agree with experiment very well, indicating it is possible now to calculate complete dynamical information for some simple fouratom reactions, as have been done for threeatom reactions in the past decades.

Quasiclassical trajectories study of Ne_{2}Br_{2}(B) vibrational predissociation: Kinetics and product distributions
View Description Hide DescriptionThe vibrational predissociation of the Ne_{2}Br_{2}(B) van der Waals complex has been investigated using the quasiclassical trajectory method (QCT), in the range of vibrational levels v ^{′} = 16–23. Extensive comparison is made with the most recent experimental observations [Pio et al., J. Chem. Phys.133, 014305 (2010)]10.1063/1.3456550, molecular dynamics with quantum transitions simulations [Miguel et al., Faraday Discuss.118, 257 (2001)]10.1039/b009222n, and preliminary results from 24dimensional Cartesian coupled coherent state (CCCS) calculations. A sequential mechanism is found to accurately describe the theoretical dynamical evolution of intermediate and final product populations, and both QCT and CCCS provide very good estimates for the dissociation lifetimes. The capabilities of QCT in the description of the fragmentation kinetics are analyzed in detail by using reduceddimensionality models of the complexes and concepts from phasespace transport theory. The problem of fast decoupling of the different coherent states in CCCS simulations, resulting from the high dimensionality of phase space, is tackled using a reexpansion scheme. QCT rovibrational product state distributions are reported. Due to the weakness of the van der Waals couplings and the low density of vibrational states,QCT predicts a larger than observed propensity for Δv′ = −1 and −2 channels for the respective dissociation of the first and second Ne atoms.