Volume 136, Issue 19, 21 May 2012

We present a model for the local diffusionrelaxation dynamics of the C _{α}atoms in proteins describing both the diffusive shorttime dynamics and the asymptotic longtime relaxation of the position autocorrelation functions. The relaxation rate spectra of the latter are represented by shifted gamma distributions, where the standard gamma distribution describes anomalous slow relaxation in macromolecular systems of infinite size and the shift accounts for a smallest local relaxation rate in macromolecules of finite size. The resulting autocorrelation functions are analytic for any time t ⩾ 0. Using results from a molecular dynamics simulation of lysozyme, we demonstrate that the model fits the position autocorrelation functions of the C _{α}atoms exceptionally well and reveals moreover a strong correlation between the residue's solventaccessible surface and the fitted model parameters.
 COMMUNICATIONS


Communication: A minimal model for the diffusionrelaxation backbone dynamics of proteins
View Description Hide DescriptionWe present a model for the local diffusionrelaxation dynamics of the C _{α}atoms in proteins describing both the diffusive shorttime dynamics and the asymptotic longtime relaxation of the position autocorrelation functions. The relaxation rate spectra of the latter are represented by shifted gamma distributions, where the standard gamma distribution describes anomalous slow relaxation in macromolecular systems of infinite size and the shift accounts for a smallest local relaxation rate in macromolecules of finite size. The resulting autocorrelation functions are analytic for any time t ⩾ 0. Using results from a molecular dynamics simulation of lysozyme, we demonstrate that the model fits the position autocorrelation functions of the C _{α}atoms exceptionally well and reveals moreover a strong correlation between the residue's solventaccessible surface and the fitted model parameters.
 Top

 ARTICLES

 Theoretical Methods and Algorithms

Replica exchange with nonequilibrium switches: Enhancing equilibrium sampling by increasing replica overlap
View Description Hide DescriptionWe describe a replica exchange strategy where trial swap configurations are generated by nonequilibrium switching simulations. By devoting simulation time to the switching simulations, one can systematically increase an effective overlap between replicas, which leads to an increased exchange acceptance rate and less correlated equilibrium samples. In this paper, we derive our method for a general class of stochastic dynamics, and discuss various strategies for enhancing replica overlap through novel dynamical schemes and prudent choices of switching protocols. We then demonstrate our method on a model system of alanine dipeptide in implicit solvent, characterizing decreases in data correlations and gains in sampling efficiency.

A new multiscale modeling method for simulating the loss processes in polymer solar cell nanodevices
View Description Hide DescriptionThe photoelectric power conversion efficiency of polymersolar cells is till now, compared to conventional inorganic solar cells, still relatively low with maximum values ranging from 7% to 8%. This essentially relates to the existence of exciton and charge carrier loss phenomena, reducing the performance of polymersolar cells significantly. In this paper we introduce a new computer simulation technique, which permits to explore the causes of the occurrence of such phenomena at the nanoscale and to design new photovoltaic materials with optimized optoelectronic properties. Our approach consists in coupling a mesoscopic fieldtheoretic method with a suitable dynamic Monte Carlo algorithm, to model the elementary photovoltaic processes. Using this algorithm, we investigate the influence of structural characteristics and different device conditions on the excitongeneration and charge transport efficiencies in case of a novel nanostructured polymer blend. More specifically, we find that the disjunction of continuous percolation paths leads to the creation of dead ends, resulting in charge carrier losses through charge recombination. Moreover, we observe that defects are characterized by a low exciton dissociation efficiency due to a high charge accumulation, counteracting the chargegeneration process. From these observations, we conclude that both the charge carrier loss and the exciton loss phenomena lead to a dramatic decrease in the internal quantum efficiency. Finally, by analyzing the photovoltaic behavior of the nanostructures under different circuit conditions, we demonstrate that charge injection significantly determines the impact of the defects on the solar cell performance.

Equation of state of charged colloidal suspensions and its dependence on the thermodynamic route
View Description Hide DescriptionThe thermodynamic properties of highly charged colloidalsuspensions in contact with a salt reservoir are investigated in the framework of the renormalized Jellium model (RJM). It is found that the equation of state is very sensitive to the particular thermodynamic route used to obtain it. Specifically, the osmotic pressure calculated within the RJM using the contact value theorem can be very different from the pressure calculated using the KirkwoodBuff fluctuation relations. On the other hand, Monte Carlo simulations show that both the effective pair potentials and the correlation functions are accurately predicted by the RJM. It is suggested that the lack of selfconsistency in the thermodynamics of the RJM is a result of neglected electrostaticcorrelations between the counterions and coions.

Spin densities from subsystem densityfunctional theory: Assessment and application to a photosynthetic reaction center complex model
View Description Hide DescriptionSubsystem densityfunctional theory(DFT) is a powerful and efficient alternative to Kohn–Sham DFT for large systems composed of several weakly interacting subunits. Here, we provide a systematic investigation of the spindensity distributions obtained in subsystem DFT calculations for radicals in explicit environments. This includes a small radical in a solvent shell, a πstacked guanine–thymine radical cation, and a benchmark application to a model for the special pair radical cation, which is a dimer of bacteriochlorophyll pigments, from the photosynthetic reaction center of purple bacteria. We investigate the differences in the spin densities resulting from subsystem DFT and Kohn–Sham DFT calculations. In these comparisons, we focus on the problem of overdelocalization of spin densities due to the selfinteraction error in DFT. It is demonstrated that subsystem DFT can reduce this problem, while it still allows to describe spinpolarization effects crossing the boundaries of the subsystems. In practical calculations of spin densities for radicals in a given environment, it may thus be a pragmatic alternative to Kohn–Sham DFT calculations. In our calculation on the special pair radical cation, we show that the coordinating histidine residues reduce the spindensity asymmetry between the two halves of this system, while inclusion of a larger binding pocket model increases this asymmetry. The unidirectional energy transfer in photosynthetic reaction centers is related to the asymmetry introduced by the protein environment.

Hybrid coupled cluster methods: Combining active space coupled cluster methods with coupled cluster singles, doubles, and perturbative triples
View Description Hide DescriptionBased on the coupledcluster singles, doubles, and a hybrid treatment of triples (CCSD(T)h) method developed by us [J. Shen, E. Xu, Z. Kou, and S. Li, J. Chem. Phys.132, 114115 (2010)10.1063/1.3359851;J. Shen, E. Xu, Z. Kou, and S. Li, J. Chem. Phys.133, 234106 (2010)10.1063/1.3518100;J. Shen, E. Xu, Z. Kou, and S. Li, J. Chem. Phys.134, 044134 (2011)10.1063/1.3541250], we developed and implemented a new hybrid coupled cluster (CC) method, named CCSD(T)qh, by combining CC singles and doubles, and active triples and quadruples (CCSDtq) with CCSD(T) to deal with the electronic structures of molecules with significant multireference character. These two hybrid CC methods can be solved with noncanonical and canonical MOs. With canonical MOs, the CCSD(T)like equations in these two methods can be solved directly without iteration so that the storage of all triple excitation amplitudes can be avoided. A practical procedure to divide canonical MOs into active and inactive subsets is proposed. Numerical calculations demonstrated that CCSD(T)h with canonical MOs can well reproduce the corresponding results obtained with noncanonical MOs. For three atom exchange reactions, we found that CCSD(T)h can offer a significant improvement over the popular CCSD(T) method in describing the reaction barriers. For the bondbreaking processes in F_{2} and H_{2}O, our calculations demonstrated that CCSD(T)qh is a good approximation to CCSDTQ over the entire bonddissociation processes.

Trapping photondressed Dirac electrons in a quantum dot studied by coherent two dimensional photon echo spectroscopy
View Description Hide DescriptionWe study the localization of dressed Diracelectrons in a cylindrical quantum dot(QD) formed on monolayer and bilayer graphene by spatially different potential profiles. Short lived excitonic states which are too broad to be resolved in linear spectroscopy are revealed by cross peaks in the photonecho nonlinear technique. Signatures of the dynamic gap in the twodimensional spectra are discussed. The effect of the Coulomb induced excitonexciton scattering and the formation of biexciton molecules are demonstrated.

Ionization energy of atoms obtained from GW selfenergy or from random phase approximation total energies
View Description Hide DescriptionA systematic evaluation of the ionization energy within the GW approximation is carried out for the first row atoms, from H to Ar. We describe a Gaussian basis implementation of the GW approximation, which does not resort to any further technical approximation, besides the choice of the basis set for the electronic wavefunctions. Different approaches to the GW approximation have been implemented and tested, for example, the standard perturbative approach based on a prior meanfield calculation (HartreeFock GW@HF or densityfunctional theoryGW@DFT) or the recently developed quasiparticle selfconsistent method (QSGW). The highest occupied molecular orbital energies of atoms obtained from both GW@HF and QSGW are in excellent agreement with the experimental ionization energy. The lowest unoccupied molecular orbital energies of the singly charged cation yield a noticeably worse estimate of the ionization energy. The best agreement with respect to experiment is obtained from the total energy differences within the random phase approximation functional, which is the total energy corresponding to the GW selfenergy. We conclude with a discussion about the slight concave behavior upon number electron change of the GW approximation and its consequences upon the quality of the orbital energies.

Potential flux landscapes determine the global stability of a Lorenz chaotic attractor under intrinsic fluctuations
View Description Hide DescriptionWe developed a potential flux landscape theory to investigate the dynamics and the global stability of a chemical Lorenz chaotic strange attractor under intrinsic fluctuations. Landscape was uncovered to have a butterfly shape. For chaotic systems, both landscape and probabilistic flux are crucial to the dynamics of chaotic oscillations. Landscape attracts the system down to the chaotic attractor, while flux drives the coherent motions along the chaotic attractors. Barrier heights from the landscape topography provide a quantitative measure for the robustness of chaotic attractor. We also found that the entropy production rate and phase coherence increase as the molecular numbers increase. Power spectrum analysis of autocorrelation function provides another way to quantify the global stability of chaotic attractor. We further found that limit cycle requires more flux and energy to sustain than the chaotic strange attractor. Finally, by detailed analysis we found that the curl probabilistic flux may provide the origin of the chaotic attractor.

Ab initio quantum dynamics using coupledcluster
View Description Hide DescriptionThe curse of dimensionality (COD) limits the current stateoftheart ab initio propagation methods for nonrelativistic quantum mechanics to relatively few particles. For stationary structure calculations, the coupledcluster (CC) method overcomes the COD in the sense that the method scales polynomially with the number of particles while still being sizeconsistent and extensive. We generalize the CC method to the time domain while allowing the singleparticle functions to vary in an adaptive fashion as well, thereby creating a highly flexible, polynomially scaling approximation to the timedependent Schrödinger equation. The method inherits sizeconsistency and extensivity from the CC method. The method is dubbed orbitaladaptive timedependent coupledcluster, and is a hierarchy of approximations to the now standard multiconfigurational timedependent Hartree method for fermions. A numerical experiment is also given.

Generalized Gibbs state with modified Redfield solution: Exact agreement up to second order
View Description Hide DescriptionA novel scheme for the steady state solution of the standard Redfield quantum master equation is developed which yields agreement with the exact result for the corresponding reduced density matrix up to second order in the systembath coupling strength. We achieve this objective by use of an analytic continuation of the offdiagonal matrix elements of the Redfield solution towards its diagonal limit. Notably, our scheme does not require the provision of yet higher order relaxation tensors. Testing this modified method for a heat bath consisting of a collection of harmonic oscillators we assess that the system relaxes towards its correct couplingdependent, generalized quantum Gibbs state in second order. We numerically compare our formulation for a damped quantum harmonic system with the nonequilibrium Green's function formalism: we find good agreement at low temperatures for coupling strengths that are even larger than expected from the very regime of validity of the secondorder Redfield quantum master equation. Yet another advantage of our method is that it markedly reduces the numerical complexity of the problem; thus, allowing to study efficiently largesized systemHilbert spaces.

Construction of highdimensional neural network potentials using environmentdependent atom pairs
View Description Hide DescriptionAn accurate determination of the potential energy is the crucial step in computer simulations of chemical processes, but using electronic structure methods onthefly in molecular dynamics (MD) is computationally too demanding for many systems. Constructing more efficient interatomic potentials becomes intricate with increasing dimensionality of the potentialenergysurface (PES), and for numerous systems the accuracy that can be achieved is still not satisfying and far from the reliability of firstprinciples calculations. Feedforward neural networks (NNs) have a very flexible functional form, and in recent years they have been shown to be an accurate tool to construct efficient PESs. Highdimensional NN potentials based on environmentdependent atomic energy contributions have been presented for a number of materials. Still, these potentials may be improved by a more detailed structural description, e.g., in form of atom pairs, which directly reflect the atomic interactions and take the chemical environment into account. We present an implementation of an NN method based on atom pairs, and its accuracy and performance are compared to the atombased NN approach using two very different systems, the methanol molecule and metallic copper. We find that both types of NN potentials provide an excellent description of both PESs, with the pairbased method yielding a slightly higher accuracy making it a competitive alternative for addressing complex systems in MD simulations.

The threeelectron harmonium atom: The lowestenergy doublet and quadruplet states
View Description Hide DescriptionCalculations of subμhartree accuracy employing explicitly correlated Gaussian lobe functions produce comprehensive data on the energy E(ω), its components, and the oneelectron properties of the two lowestenergy states of the threeelectron harmonium atom. The energy computations at 19 values of the confinement strength ω ranging from 0.001 to 1000.0, used in conjunction with a recently proposed robust interpolation scheme, yield explicit approximants capable of estimating E(ω) and the potential energy of the harmonic confinement within a few tenths of μhartree for any ω ⩾ 0.001, the respective errors for the kinetic energy and the potential energy of the electronelectron repulsion not exceeding 2 μhartrees. Thanks to the correct ω → 0 asymptotics incorporated into the approximants, comparable accuracy is expected for values of ω smaller than 0.001. Occupation numbers of the dominant natural spinorbitals and two different measures of electron correlation are also computed.

The multiscale coarsegraining method. VIII. Multiresolution hierarchical basis functions and basis function selection in the construction of coarsegrained force fields
View Description Hide DescriptionThe multiscale coarsegraining (MSCG) method is a method for determining the effective potential energy function for a coarsegrained (CG) model of a molecular system using data obtained from molecular dynamics simulation of the corresponding atomically detailed model. The coarsegrained potential obtained using the MSCG method is a variational approximation for the exact manybody potential of mean force for the coarsegrained sites. Here we propose a new numerical algorithm with noise suppression capabilities and enhanced numerical stability for the solution of the MSCG variational problem. The new method, which is a variant of the elastic net method [Friedman et al., Ann. Appl. Stat.1, 302 (2007)]10.1214/07AOAS131, allows us to construct a large basis set, and for each value of a socalled “penalty parameter” the method automatically chooses a subset of the basis that is most important for representing the MSCG potential. The size of the subset increases as the penalty parameter is decreased. The appropriate value to choose for the penalty parameter is the one that gives a basis set that is large enough to fit the data in the simulation data set without fitting the noise. This procedure provides regularization to mitigate potential numerical problems in the associated linear least squares calculation, and it provides a way to avoid fitting statistical error. We also develop new basis functions that are similar to multiresolution Haar functions and that have the differentiability properties that are appropriate for representing CG potentials. We demonstrate the feasibility of the combined use of the elastic net method and the multiresolution basis functions by performing a variational calculation of the CG potential for a relatively simple system. We develop a method to choose the appropriate value of the penalty parameter to give the optimal basis set. The combined effect of the new basis functions and the regularization provided by the elastic net method opens the possibility of using very large basis sets for complicated CG systems with many interaction potentials without encountering numerical problems in the variational calculation.

The multiscale coarsegraining method. IX. A general method for construction of three body coarsegrained force fields
View Description Hide DescriptionThe multiscale coarsegraining (MSCG) method is a method for constructing a coarsegrained (CG) model of a system using data obtained from molecular dynamics simulations of the corresponding atomically detailed model. The formal statistical mechanical derivation of the method shows that the potential energy function extracted from an MSCG calculation is a variational approximation for the true potential of mean force of the CG sites, one that becomes exact in the limit that a complete basis set is used in the variational calculation if enough data are obtained from the atomistic simulations. Most applications of the MSCG method have employed a representation for the nonbonded part of the CG potential that is a sum of all possible pair interactions. This approach, despite being quite successful for some CG models, is inadequate for some others. Here we propose a systematic method for including three body terms as well as two body terms in the nonbonded part of the CG potential energy. The current method is more general than a previous version presented in a recent paper of this series [L. Larini, L. Lu, and G. A. Voth, J. Chem. Phys.132, 164107 (2010)]10.1063/1.3394863, in the sense that it does not make any restrictive choices for the functional form of the three body potential. We use hierarchical multiresolution functions that are similar to wavelets to develop very flexible basis function expansions with both two and three body basis functions. The variational problem is solved by a numerical technique that is capable of automatically selecting an appropriate subset of basis functions from a large initial set. We apply the method to two very different coarsegrained models: a solvent free model of a two component solution made of identical LennardJones particles and a one site model of SPC/E water where a site is placed at the center of mass of each water molecule. These calculations show that the inclusion of three body terms in the nonbonded CG potential can lead to significant improvement in the accuracy of CG potentials and hence of CG simulations.

The multiscale coarsegraining method. X. Improved algorithms for constructing coarsegrained potentials for molecular systems
View Description Hide DescriptionThe multiscale coarsegraining (MSCG) method uses simulation data for an atomistic model of a system to construct a coarsegrained (CG) potential for a coarsegrained model of the system. The CG potential is a variational approximation for the true potential of mean force of the degrees of freedom retained in the CG model. The variational calculation uses information about the atomistic positions and forces in the simulation data. In principle, the resulting MSCG potential will be an accurate representation of the true CG potential if the basis set for the variational calculation is complete enough and the canonical distribution of atomistic states is well sampled by the data set. In practice, atomistic configurations that have very high potential energy are not sampled. As a result there usually is a region of CG configuration space that is not sampled and about which the data set contains no information regarding the gradient of the true potential. The MSCG potential obtained from a variational calculation will not necessarily be accurate in this unsampled region. A priori considerations make it clear that the true CG potential of mean force must be very large and positive in that region. To obtain an MSCG potential whose behavior in the sampled region is determined by the atomistic data set, and whose behavior in the unsampled region is large and positive, it is necessary to intervene in the variational calculation in some way. In this paper, we discuss and compare two such methods of intervention, which have been used in previous MSCG calculations for dealing with nonbonded interactions. For the test systems studied, the two methods give similar results and yield MSCG potentials that are limited in accuracy only by the incompleteness of the basis set and the statistical error of associated with the set of atomistic configurations used. The use of such methods is important for obtaining accurate CG potentials.
 Gas Phase Dynamics and Structure: Spectroscopy, Molecular Interactions, Scattering, and Photochemistry

Fulldimensional quantum dynamics study of exchange processes for the D + H_{2}O and D + HOD reactions
View Description Hide DescriptionThe exchange processes of D + H_{2}O and D + HOD reactions are studied using initial stateselected timedependent wave packet approach in full dimension. The total reaction probabilities for different partial waves, together with the integral cross sections, are obtained both by the centrifugal sudden (CS) approximation and exact coupledchannel (CC) calculations, for the H_{2}O(HOD) reactant initially in the ground rovibrational state. In the CC calculations, small resonance peaks in the reaction probabilities and quick diminishing of the resonance peaks with the increase of total angular momenta J do not lead to clear steplike features just above the threshold in the cross sections for the title reactions, which are different in other isotopically substituted reactions where the hydrogen atom was included as the reactant instead of the deuterium atom [B. Fu, Y. Zhou, and D. H. Zhang, Chem. Sci.3, 270 (2012)10.1039/c1sc00684c; B. Fu and D. H. Zhang, J. Phys. Chem. A116, 820 (2012)10.1021/jp211096q]. It is interesting that the shape resonanceinduced features resulting from the reactiontunneling are significantly diminished accordingly in the reactions of the deuterium atom and H_{2}O or HOD, owing to the weaker tunneling capability of the reagent deuterium atom in the title reactions than the reagent hydrogen atom in other reactions. In the CS calculations, the resonance peaks persist in many partial waves but cannot survive the partialwave summations. The cross sections for the D^{′} + H_{2}O → D^{′}OH + H and D^{′} + HOD → D^{′}OD + H reactions are substantially larger than those for the D^{′} + HOD → HOD^{′} + D reaction, indicating that the D^{′}/H exchange reactions are much more favored than the D^{′}/D exchange.

A molecular H_{2} potential for heterogeneous simulations including polarization and manybody van der Waals interactions
View Description Hide DescriptionA highly accurate aniostropic intermolecular potential for diatomic hydrogen has been developed that is transferable for molecular modeling in heterogeneous systems. The potential surface is designed to be efficacious in modeling mixed sorbates in metalorganic materials that include sorption interactions with charged interfaces and open metal sites. The potential parameters are compatible for mixed simulations but still maintain high accuracy while deriving dispersion parameters from a proven polarizability model. The potential includes essential physical interactions including: shortrange repulsions, dispersion, and permanent and induced electrostatics. Manybody polarization is introduced via a pointatomic polarizability model that is also extended to account for manybody van der Waals interactions in a consistent fashion. Permanent electrostatics are incorporated using point partial charges on atomic sites. However, contrary to expectation, the best potentials are obtained by permitting the charges to take on values that do not reproduce the first nonvanishing moment of the electrostatic potential surface, i.e., the quadrupole moment. Potential parameters are fit to match ab initio energies for a representative range of dimer geometries. The resulting potential is shown to be highly effective by comparing to electronic structure calculations for a thermal distribution of trimer geometries, and by reproducing experimental bulk pressuredensity isotherms. The surface is shown to be superior to other similarly portable potential choices even in tests on homogeneous systems without strong polarizing fields. The present streamlined approach to developing such potentials allows for a simple adaptation to other molecules amenable to investigation by highlevel electronic structure methods.

Multielectron effects in high harmonic generation in N_{2} and benzene: Simulation using a nonadiabatic quantum molecular dynamics approach for lasermolecule interactions
View Description Hide DescriptionA mixed quantumclassical approach is introduced which allows the dynamical response of molecules driven far from equilibrium to be modeled. This method is applied to the interaction of molecules with intense, shortduration laser pulses. The electronic response of the molecule is described using timedependent density functional theory (TDDFT) and the resulting KohnSham equations are solved numerically using finite difference techniques in conjunction with local and global adaptations of an underlying grid in curvilinear coordinates. Using this approach, simulations can be carried out for a wide range of molecules and both allelectron and pseudopotential calculations are possible. The approach is applied to the study of high harmonic generation in N_{2} and benzene using linearly polarized laser pulses and, to the best of our knowledge, the results for benzene represent the first TDDFT calculations of high harmonic generation in benzene using linearly polarized laser pulses. For N_{2} an enhancement of the cutoff harmonics is observed whenever the laser polarization is aligned perpendicular to the molecular axis. This enhancement is attributed to the symmetry properties of the KohnSham orbital that responds predominantly to the pulse. In benzene we predict that a suppression in the cutoff harmonics occurs whenever the laser polarization is aligned parallel to the molecular plane. We attribute this suppression to the symmetryinduced response of the highestoccupied molecular orbital.

Photoelectron spectroscopy and theoretical studies of UF_{5} ^{−} and UF_{6} ^{−}
View Description Hide DescriptionThe UF_{5} ^{−} and UF_{6} ^{−} anions are produced using electrospray ionization and investigated by photoelectron spectroscopy and relativistic quantum chemistry. An extensive vibrational progression is observed in the spectra of UF_{5} ^{−}, indicating significant geometry changes between the anion and neutral ground state.FranckCondon factor simulations of the observed vibrational progression yield an adiabatic electron detachment energy of 3.82 ± 0.05 eV for UF_{5} ^{−}. Relativistic quantum calculations using density functional and ab initiotheories are performed on UF_{5} ^{−} and UF_{6} ^{−} and their neutrals. The ground states of UF_{5} ^{−} and UF_{5} are found to have C_{4v} symmetry, but with a large U−F bond length change. The ground state of UF_{5} ^{−} is a triplet state (^{3}B_{2}) with the two 5f electrons occupying a 5f_{z3 }based 8a_{1} highest occupied molecular orbital (HOMO) and the 5f_{xyz}based 2b_{2} HOMO1 orbital. The detachment cross section from the 5f_{xyz} orbital is observed to be extremely small and the detachment transition from the 2b_{2} orbital is more than ten times weaker than that from the 8a_{1} orbital at the photon energies available. The UF_{6} ^{−} anion is found to be octahedral, similar to neutral UF_{6} with the extra electron occupying the 5f_{xyz}based a_{2u} orbital. Surprisingly, no photoelectron spectrum could be observed for UF_{6} ^{−} due to the extremely low detachment cross section from the 5f_{xyz}based HOMO of UF_{6} ^{−}.