September 29 - October 3, 2008
The exact simulation of quantum mechanical systems on classical computers generally scales exponentially with the size of the system N. Using quantum computers, the computational resources required to carry out the simulation are polynomial. Our group has been working in the development and characterization of quantum computational algorithms for the simulation of chemical systems. We will give a tutorial on our algorithms for the simulation of molecular electronic structure, molecular properties and quantum dynamics, and will discuss the opportunities, open questions and challenges in the field of simulation of physical systems using quantum computers or dedicated quantum devices.
The application of conventional GGA, and meta-GGA, density functionals to van der Waals complexes is fraught with difficulties. Conventional functionals do not contain the physics of the dispersion interaction. To make matters worse, the exchange part alone can yield anything from severe over “binding” to severe over repulsion depending on the choice of functional. We rectify these problems by
a) adding a dispersion term with nonempirical C6, C8, and C10 dispersion coefficients (the Becke-Johnson dispersion model),
b) selecting a GGA exchange functional (PW86, also nonempirical) that gives excellent agreement with exact Hartree-Fock exchange repulsion curves.
The result is a simple GGA+dispersion theory giving excellent noble-gas pair interaction energies for He through Kr with only two adjustable parameters in the dispersion cutoff.
Studies of molecular dynamics and molecular spectroscopy generally
start from the Born-Oppenheimer approximation and require some form of
analytical potential energy surface fitted to ab initio electronic
structure calculations. We have used computational invariant theory
and the MAGMA computer algebra system as an aid to develop
representations for the potential energy and dipole moment surfaces
that are fully invariant under permutations of like nuclei. We
express the potential energy surface in terms of internuclear
distances using basis functions that are manifestly invariant. A
dipole moment is represented with use of effective charges at
positions of the nuclei, which must transform as a covariant, rather
than as an invariant, under permutations of like nuclei.
Malonaldehyde (CHOHCHCHO) provides an illustrative application. The
associated molecular permutational symmmetry group is of order 288
(4!3!2!) and the use of full permutational symmetry makes it possible
to obtain a compact representation for the surface.
We explored the relative stability, structure, and conductance of
crossed nanotube junctions with dispersion corrected density functional
theory. We found that the most stable junction geometry, not studied
before, displays the smallest conductance. While the conductance
increases as force is applied, it levels off very rapidly. This behavior
contrasts with a less stable junction geometry that show steady increase
of the conductance as force is applied. Electromechanical sensing
devices based on this effect should exploit the conductance changes
close to equilibrium.
Recent work in my group has focussed on the semiclassical
origins of density functional theory, and how much of modern DFT
can be understood in these terms, including the limitations of present
approximations. I will discuss this in detail for model systems,
describing a method that avoids DFT altogether.
This leads to a grand algorithmic challenge, whose solution could
revolutionize electronic structure calculations, by allowing much
larger numbers of electrons to be tackled.
Recent advances in treating electrons and nuclei (typically protons) quantum mechanically without the Born-Oppenheimer approximation using nuclear-electronic orbital (NEO) method and multicomponent density functional theory (MCDFT) is presented. Electron-proton dynamical correlation is highly significant because of the attractive electrostatic interaction between the electron and the proton. Inadequate treatment of electron-proton correlation produces nuclear densities that are too localized, resulting in abnormally high stretching frequencies, as well as inaccuracies in thermally averaged geometries and isotope effects. To address this problem, an explicitly correlated Hartree-Fock (NEO-XCHF) scheme has been formulated to include explicit electron-proton correlation directly into the nuclear-electronic orbital self-consistent-field framework. This approach is based on a general ansatz for the nuclear-electronic wavefunction in which explicit dependence on the electron-proton distance is incorporated into the total wavefunction using Gaussian-type geminal functions. A multicomponent density functional theory (MCDFT) has also been formulated by developing electron-proton functionals based on the explicitly correlated ansatz for the nuclear-electronic wavefunction. Benchmark calculations illustrate that these new methods significantly improve the description of the nuclear densities, thereby leading to more accurate hydrogen vibrational frequencies and vibrationally averaged geometries.
Density functional theory of electronic structure is widely and successfully applied in simulations throughout engineering and sciences. However, for many predicted properties there are spectacular failures that can be traced to the delocalization error and static correlation error
of commonly used approximations. These errors can be characterized and understood through the perspective of fractional charges and fractional
spins introduced recently. Reducing these errors will open new frontiers for applications of density functional theory.
Joint work with Reinhart Ahlrichs.
Density functional theory provides a systematic approach to the electronic structure of atoms, molecules and solids. It requires the repeated solution of single particle Schrodinger equations in a self consistent loop. Most techniques involve some sort of basis set, the most common ones being plane waves or Gaussians. In crystalline materials the most accurate solutions involve augmented basis sets. These combine numerical solutions of the Schrodinger equation in regions near the atomic nucleii with so called ‘tail functions’ in more distant regions. In the linear augmented plane wave (LAPW) method the tail functions are plane waves. This formulation has been incorporated into the WIEN2k code. With the current interest in nanoscale clusters, biomolecules, and other finite systems it is desirable to have a comparably accurate method for these. While it is always possible to build supercells, it is often convenient to have completely localized functions which eliminate interaction between periodic images. We recently proposed a finite cluster version of the linear augmented Slater-type orbital (LASTO) method [1]. STO’s have the correct behavior at large distances and possess an addition theorem – they can be re-expanded about other sites with analytic coefficients. We solve the Poisson equation by replacing the spherical part of the density near the nucleii with a smooth pseudo-density. The full potential, including the non-sphrical piece is then solved on a grid. Examples of small clusters and comparison with the Gaussian based program NWChem will be given.
[1] K. S. Kang, J. W. Davenport, J. Glimm, D. E. Keyes, and M. McGuigan, submitted to J. Computational Chemistry.
We have implemented relativisitic GVVPT2 using DKH integrals and ANO-RCC basis sets from Molcas package. It is done by developing an interface code accessing and transforming one- and two-electron integral array, making it available for any method implemented within the UNDMol package. The relativistic GVVPT2 is applied to calculations of ground and excited states potential energy curves of TiC and CrH.
Joint work with Alexander Gaenko and Jochen Autschbach.
We develop a time-dependence quasirelativistic density
functional ther based on the ZORA approximation for computing
frequency dependent linear response of molecules. Density fitting
was used for the calculation of complex components of the frequency dependent dipole-dipole polarizability. CPKS equations
based on 2-componenent ZORA response were derived. Using damping techniques excitation energy corresponding to the poles of the polarizability curves were calculated. We present the results of the calculations of complex dipole-dipole polarizability, of two and three dimensional gold clusters, and absorption spectra of heavy metal oxides.
Spectroscopic studies of the hydroxymethyl and 1-hydroxyethyl radicals have found an unusually large difference in their ionization energies (IE). The anticipated decrease in the IE of the latter radical due to its larger size does not fully account for the experimentally observed difference of 0.92 eV. Here we investigated the problem with the aid of electronic structure calculations. We found that the large drop in the IE of 1-hydroxyethyl radical is a result of the combined effects of the destabilization of its highest occupied molecular orbital and the stabilization of the corresponding cation due to hyperconjugation. This qualitative explanation agrees with a simple Huckel-like approach and is also consistent with Natural Bond Orbital calculation results.
We are interested in realizing the variational calculation with
the 2-order Reduced Density Matrix (2-RDM) for fermionic systems.
The known necessary N-representability conditions P, Q, G, T1, and
T2' are imposed, resulting in an optimization problem called semidefinite
programming (SDP) problem. The ground state energies for various small
atoms and molecules were calculated. Additionally, we show some results
of the one-dimensional Hubbard model for high correlation limit using
multiple precision arithmetic version of the solver.
We present a new linear scaling method for electronic structure computations in the context of Kohn-Sham density functional theory (DFT). The method is based on a subspace iteration, and takes advantage of the non-orthogonal formulation of the Kohn-Sham functional, and the improved localization properties of non-orthogonal wave functions. We demonstrate the efficiency of the algorithm for practical applications by performing fully three-dimensional computations of the electronic density of alkane chains.
This is joint work with Jianfeng Lu, Yulin Xuan, and Weinan E, at Princeton University.
The Resolution of the Identity (RI) is widely used in many-body algorithms. It expresses the completeness of a set of functions that possess the familiar orthonormality property.
In the first part of my lecture, I will discuss functions that possess an analogous property called Coulomb-orthonormality and which permit us to resolve the two-particle Coulomb operator into a sum of products of one-particle functions. Connections to Cholesky decomposition and Kronecker-product approximation will be made.
In the second part of the lecture, I will present and discuss numerical applications of Coulomb resolution in the context of electronic structure theory.
A new generation of density-functional methods is based on orbital-dependent
funtionals. With orbital-dependent functionals long-standing problems like
the occurence of unphysical Coulomb self-interactions or the qualitatively
wrong description of charge-transfer excitation in time-dependent
density-functional theory can be solved. Orbital-dependent functionals indeed
may represent the future of density-functional theory. In oder to determine
exchange-correlation potentials corresponding to orbital-dependent energy
functionals, however, it is necessary to solve a numerically very demanding
integral equation with the optimized effective potential (OEP) method.
Methods to handle this integral equation thus are required for the further
development of density-functional methods. The numerical stability of the OEP
integral equation is investigated and method to solve it are presented.
The green fluorescent protein (GFP) is widely used in biochemical and medical studies as a biomarker in living cells. Modeling properties of GFP is an essential step in the efforts to enhance efficiency of this in vivo marker and to expand the area of its applications. We apply the methods of electronic structure calculations, including the quantum chemistry methods and the combined quantum mechanical – molecular mechanical (QM/MM) approaches, to describe the structure, spectra and transformations of the GFP chromophore, 4-hydroxybenzylidene-imidazolinone, in the gas phase, solutions and in the protein matrix. We compare the results of calculations for the cis-trans chromophore isomerization in the ground electronic state in the gas phase by using the DFT approach PBE0/6-31+G** and the CASSCF(12,11)/cc-pVDZ approximation. The energy profiles computed with both methods are markedly different in the vicinity of the saddle point. The isomerization paths computed in the QM/MM approach for the chromophore buried inside the water cluster show that the CASSCF results are better consistent with the experimental observations than the DFT findings. We also report the results of QM/MM calculations with the DFT approximations in the QM subsystem for the geometry configurations of the chromophore binding pocket inside the protein by assuming various protonation states of the chromophore unit, anionic, neutral, zwitterionic, cationic. These structures are employed for calculations of the photoexcitation pathways in GFP.
First-principles molecular dynamics (FPMD) is a simulation method that combines molecular dynamics with the accuracy of a quantum mechanical description of electronic structure. It is increasingly used to address problems of structure determination, statistical mechanics, and electronic structure of solids, liquids and nanoparticles. The high computational cost of this approach makes it a good candidate for use on large-scale computers. In order to achieve high performance on terascale and petascale computers, current FPMD algorithms have to be reexamined and redesigned. We present new, large-scale parallel algorithms developed for FPMD simulations on computers including O(10^{3}) to O(10^{4}) CPUs. Examples include the problem of simultaneous diagonalization of symmetric matrices used in the calculation of Maximally Localized Wannier Functions (MLWFs), and the Orthogonal Procrustes problem that arises in the context of Born-Oppenheimer molecular dynamics simulations.
Supported by NSF-OCI PetaApps through grant 0749217.
We study convergence or divergence of the Moller–Plesset
perturbation series for systems with two electrons and a single nucleus of
charge Z > 0. This question is essentially to determine if the radius of
convergence of a power series in the complex perturbation parameter lambda
is greater than 1. We examine a simple one-dimensional model with delta
functions in place of Coulomb potentials and the realistic
three-dimensional model. For each model, we show rigorously that if the
nuclear charge Z is sufficiently large, there are no singularities for
real values of lambda between -1 and 1. Using a finite difference scheme,
we present numerical results for the delta function model.
Wave function-based quantum chemistry has two traditional lines of
development – one based on molecular orbitals (MO's), and the other on
valence bond (VB) theory. Both offer advantages and disadvantages for
the challenging problem of describing strong correlations, such as the
breaking of chemical bonds, or the low-spin (antiferromagnetic)
coupling of electrons on different centers.
Within MO methods, strong correlations can be viewed as those arising
within a valence orbital active space. One reasonable definition of
such a space is to supply one correlating orbital for each valence
occupied orbital. Exact solution of the Schrodinger equation in this
space is exponentially difficult with its size, and therefore
approximations are imperative. The most common workaround is to
truncate the number of orbitals defining the active space, and then
solve the truncated problem, as is done in CASSCF. An important
alternative is to systematically approximate the Schrödinger equation
in the full valence space, for example by using coupled cluster theory
ideas. I shall discuss progress in this direction.
Within spin-coupled VB theory, the target wave function consists of a
set of non-orthogonal orbitals, one for each valence electron, that
are spin-coupled together into a state of the desired overall
spin-multiplicity. The number of active orbitals is identical with
the valence space MO problem discussed above, though the problem is
not identical. Exact solution of the VB problem is exponentially
difficult with molecular size, and therefore approximations are
imperative. Again, the most common approach is to seek the exact
solution in a truncated valence orbital space, where other orbitals
are simply treated in mean-field. It is possible, however, to also
consider approximations that do not truncate the space, but rather
reduce the complexity. A new way of doing this will be introduced and
contrasted with the MO-based approaches.
We perform a rigorous mathematical analysis of the bending modes of
a linear triatomic molecule that exhibits the Renner-Teller effect.
Assuming the potentials are smooth, we prove that the wave functions
and energy levels have asymptotic expansions in powers of epsilon,
where the fourth power of epsilon is the ratio of an electron mass to the mass of a
nucleus. To prove the validity of the expansion, we must prove
various properties of the leading order equations and their
solutions. The leading order eigenvalue problem is analyzed in
terms of a parameter b, which is equivalent to the parameter
originally used by Renner. Perturbation theory and finite
difference calculations suggest that there is a crossing involving the ground bending vibrational
state near b=0.925. The crossing involves two states with
different degeneracy.
A recent algorithmic
revision of second order Generalized van Vleck perturbation theory (GVVPT2)
has proven to make the method efficacious for many challenging molecular
systems. 1 An extension to third order (GVVPT3) has been demonstrated to be
a close approximation to multireference configuration interaction including
single and double excitations (MRCISD). 2 To improve the computing
efficiency, new GVVPT codes have been developed to take advantage of
recently implemented configuration-driven configuration interaction (CI)
with unitary group approach (UGA).
Joint work with Jerzy Leszczynski,
Computational Center for Molecular Structure and Interactions,
Jackson State University, Jackson MS and Leonid Gorb,
US Army ERDC, Vicksburg, MS.
With the June announcement that RoadRunner supercomputer is the first
system to reach the petaflop level, the HPC community is entering a
realm of unprecedented computing power. More petascale computing
systems will soon be available to the scientific community. Recent
studies in the productivity of HPC platforms point to better software
as a key enabler to science on these systems.
The combination of computationally demanding electronic structure
methods with molecular dynamics is highly dependent on high-performance
computing resources. The availability of such applications constitutes
a big opportunity to evaluate both capabilities and limits of any HPC
system and software application within the framework of a real-life
feasibility study. The performance of benchmarks from the AIMD and
hybrid QM/MM simulations on two high performance computing platforms
will be discussed. Looking toward maximizing the computational
time/performance ratio, we analyzed performance data for the
Cray XT3/XT4 architectures available at ERDC.
The difficultly of approximate density functionals in describing the energetics of Diels-Alder reactions and dimerization of aluminum complexes is analyzed. Both of these reaction classes involve formation of cyclic or bicyclic products, which are found to be under-bound by the majority of functionals considered. We present a consistent view of these results from the perspective of delocalization error. This error causes approximate functionals give too low energy for delocalized densities or too high energy for localized densities, as in the cyclic and bicyclic reaction products. This interpretation allows us to understand better a wide range of errors in main-group thermochemistry obtained with popular density functionals. In general, functionals with minimal delocalization error should be used for theoretical studies of reactions where there is a loss of extended conjugation or formation of highly branched, cyclic, and cage-like molecules.
The bimolecular interaction potentials for various configurations of the ethylene dimer computed with coupled-cluster and spin-component scaled MP2 are reported. Of particular interest is any bias for particular orientations of the sigma- and pi-bonds introduced by the scaling of correlation components.
Joint work with Marat Valiev, Niri Govind, Peng-Dong Fan,
W.A. de Jong
(William R Wiley Environmental Molecular Sciences Laboratory and
Chemical Sciences Division,
Pacific Northwest National Laboratory
P.O. Box 999, MS K1-96, Richland, WA 99352) and Jeff R. Hammond
(The University of Chicago).
The coupled-cluster (CC) methodology has become a leading
formalism not only in gas-phase calculations but also in
modeling systems for which the inclusion of the surrounding
environment is critical for a comprehensive understanding of
complex photochemical reactions. At the same time it has been
proven that high-level CC formalisms are capable of providing
highly adequate characterization of excitation energies and
excited-state potential energy surfaces. With the ever
increasing power of computer platforms and highly scalable
codes, very accurate QM/MM calculations for large molecules
(defining the quantum region) can be routinely performed in the
foreseeable future even with iterative methods accounting for
the effect of triples ( CCSDT-n/EOMCCSDT-n). We will discuss
several components of recently developed and implemented CC
methodologies in NWChem. This includes: (1) Novel
iterative/non-iterative methods accounting for the effect of
triply excited configurations, (2) Massively parallel
implementations of the CC theories based on the manifold of
singly and doubly excited configurations. Several examples will
illustrate how these approaches can be used in multiscale
QM/MM framework.
Over the last decade the coupled-cluster (CC) methodology
played a dominant role in highly accurate predictions of
electronic structure. For this reason, the need for more
efficient parallel implementations is obvious. Currently, the existing
NWChem implementation can scale across thousand of CPUs and can be used
in correlating 250 electrons. Additionally, the CC codes can be used as
a quantum mechanical component of various multiscale approaches.
Range-separated (screened) hybrid functionals provide a
powerful strategy for incorporating nonlocal exact (Hartree-Fock-type)
exchange into density functional theory. Existing implementations
of range separation use a fixed, system-independent
screening parameter. Here, we propose a novel method that uses
a position-dependent screening. These locally range-separated
(LRS) hybrids add substantial flexibility for describing
diverse electronic structures and satisfy a high-density scaling
constraint better than the fixed screening approximation does.
Joint work with P.U. Manohar.
A non-iterative N7 triples correction for the equation-of-motion
coupled-cluster wave functions with single and double substitutions
(EOM-CCSD) is presented. The correction is derived by second order
perturbation treatment of the similarity-transformed CCSD Hamiltonian. The
spin-conserving variant of the correction is identical to the triples
correction of Piecuch and coworkers [Mol. Phys. 104, 2149 (2006)] derived
within method-of-moments framework and is not size-intensive. The
spin-flip variant of the correction is size-intensive. The performance of
the correction is demonstrated by calculations of electronic excitation
energies in methylene, nitrenium ion, cyclobutadiene, ortho-, meta-, and
para- benzynes, 1,2,3-tridehydrobenzene, as well as C-C bond-breaking in
ethane. In all cases except cyclobutadiene, the absolute values of the
correction for energy differences were 0.1 eV or less. In cyclobutadiene,
the absolute values of the correction were as large as 0.4 eV. In most
cases, the corrections reduced the errors against the benchmark values by
about a factor of 2 to 3, the absolute errors being less 0.04 eV.
The van der Waals density functional of Dion, Rydberg, Schroder, Langreth,
and Lundqvist [Phys. Rev. Lett. 92, 246401 (2004)] will be reviewed,
discussing implementations and applications by our group and others.
New results relevalent for hydrogen storage in metal-organic framework (MOF)
materials, as well for the intercalation of drug molecules in DNA will
be presented.
I will overview some open mathematical questions related to the
models and techniques of computational quantum chemistry. The talk is based
upon a recent article coauthored with E. Cances and PL. Lions, and
published in Nonlinearity, volume 21, T165-T176, 2008.
A number of exact relations are briefly discussed in terms of present developments and goals in DFT. In addition, conjectured relations are presented.
By means of rigorous thermodynamic limit arguments, we derive a new
variational model providing exact embedding of local defects in insulating or
semiconducting crystals. A natural way to obtain variational discretizations
of this model is to expand the perturbation of the periodic density matrix
generated by the defect in a basis of precomputed maximally localized Wannier
functions of the host crystal. This approach can be used within any
semi-empirical or Density Functional Theory framework. This is a joint work
with Eric Cancès and Amélie Deleurence (Ecole Nationale des Ponts et
Chaussées, France).
The quantal geometric phase [1-3] in a Born-Oppenheimer
(adiabatic) electronic wavefunction is a net phase change for
nuclear motion over a closed path. Effects of the quantal
geometric phase have been observed in theoretical studies of
the vibrational spectra of cyclic trinitrogen (N_{3}) molecule [4,
5]. Making a classical-quantum correspondence [6] relates the
quantal geometric phase to a classical one for cyclic nuclear
motion in N-body molecular dynamics. Each is described
differential geometrically as the holonomy of a connection [7],
physically in terms of the internal angular momentum, and with
examples. The classical geometric phase [8] is a net angle of
overall rotation in the center-of-mass frame. A net rotation
of 20 degrees has been observed experimentally in a triatomic
photodissociation and a net overall rotation of 42 degrees has
been observed computationally in protein dynamics. The
Hamiltonian operator in a generalized Born-Oppenheimer
Schrodinger equation for the electronic wavefunction is related to a classical
Hamiltonian for N-body molecular dynamics. Both the quantal
and classical geometric phases arise due to non-zero internal
angular momentum in N-body molecular dynamics.
References:
[1] C. A. Mead and D. G. Truhlar, J. Chem. Phys. 70, 2284
(1979).
[2] M. V. Berry, Proc. Royal Soc. London A 392, 45 (1984).
[3] B. Simon, Phys. Rev. Lett. 51, 2167 (1983).
[4] D. Babikov, B. K. Kendrick, P. Zhang, and K. Morokuma, J.
Chem. Phys. 122, 044315 (2005).
[5] D. Babikov, V. A. Mozhayskiy, and A. I. Krylov, J. Chem.
Phys. 125, 084306 (2006).
[6] F. J. Lin, Quantal and classical geometric phases,
2008.
[7] J. E. Marsden, R. Montgomery, and T. Ratiu, Memoirs of the
American Mathematical Society, Vol. 88, No. 436, American
Mathematical Society, Providence, RI, 1990.
[8] F. J. Lin, Discrete and Continuous Dynamical Systems,
Supplement 2007, 655 (2007).
The quantal geometric phase [1-3] in a Born-Oppenheimer
(adiabatic) electronic wavefunction is a net phase change for
nuclear motion over a closed path. Effects of the quantal
geometric phase have been observed in theoretical studies of
the vibrational spectra of cyclic trinitrogen (N_{3}) molecule [4,
5]. Making a classical-quantum correspondence [6] relates the
quantal geometric phase to a classical one for cyclic nuclear
motion in N-body molecular dynamics. Each is described
differential geometrically as the holonomy of a connection [7],
physically in terms of the internal angular momentum, and with
examples. The classical geometric phase [8] is a net angle of
overall rotation in the center-of-mass frame. A net rotation
of 20 degrees has been observed experimentally in a triatomic
photodissociation and a net overall rotation of 42 degrees has
been observed computationally in protein dynamics. The
Hamiltonian operator in a generalized Born-Oppenheimer
Schrodinger equation for the electronic wavefunction is related to a classical
Hamiltonian for N-body molecular dynamics. Both the quantal
and classical geometric phases arise due to non-zero internal
angular momentum in N-body molecular dynamics.
References:
[1] C. A. Mead and D. G. Truhlar, J. Chem. Phys. 70, 2284
(1979).
[2] M. V. Berry, Proc. Royal Soc. London A 392, 45 (1984).
[3] B. Simon, Phys. Rev. Lett. 51, 2167 (1983).
[4] D. Babikov, B. K. Kendrick, P. Zhang, and K. Morokuma, J.
Chem. Phys. 122, 044315 (2005).
[5] D. Babikov, V. A. Mozhayskiy, and A. I. Krylov, J. Chem.
Phys. 125, 084306 (2006).
[6] F. J. Lin, Quantal and classical geometric phases,
2008.
[7] J. E. Marsden, R. Montgomery, and T. Ratiu, Memoirs of the
American Mathematical Society, Vol. 88, No. 436, American
Mathematical Society, Providence, RI, 1990.
[8] F. J. Lin, Discrete and Continuous Dynamical Systems,
Supplement 2007, 655 (2007).
In this presentation I will give a review of the Cholesky Decomposition (CD) as it has
been implemented in the MOLCAS program package. These examples will include conventional CD, as
implemented for the HF, CASSCF, MP2, DFT, CASPT2 and CC methods, to the recent 1-center CD
approximation. In addition, the aCD abd acCD techniques for the on-the-fly generation of RI
auxiliary basis functions will be discussed. Analytic CD gradients will be introduced for
CD-HF, CD-DFT(pure and hybrid), and CD-CASSCF. If time allows I will briefly discuss the
use of CD technique for the fast evaluation of the exchange energy in CD-HF through the use
of CD localized orbitals.
We present a novel multiscale modeling approach that can simulate
multi-million atoms effectively via density functional theory. The
method is based on the framework of the quasicontinuum (QC) approach
with orbital-free density functional theory (OFDFT) as its sole
energetics formulation. The local QC part is formulated by the
Cauchy-Born hypothesis with OFDFT calculations for strain energy and
stress. The nonlocal QC part is treated by an OFDFT-based embedding
approach, which couples OFDFT nonlocal atoms to local region atoms.
The method - QCDFT- is applied to a nanoindentation study of an Al
thin film, and the results are compared to a conventional QC
approach. The results suggest that QCDFT represents a new direction
for the quantum simulation of materials at length scales that are
relevant to experiments.
We study the general problem of
mixing for ab-initio quantum-mechanical
problems. Guided by general mathematical principles and the underlying physics, we propose
a multisecant form of Broyden's second method for solving the self-consistent field
equations of Kohn-Sham density functional theory.
The algorithm is robust, requires relatively little fine-tuning and appears to
outperform the current state of the art, converging for cases that defeat
many other methods. We compare our technique to
the conventional methods for problems ranging from simple to nearly
pathological.
Optimization concepts will be reviewed with an eye on their proved or
potential application in Electronic Structure Calculations
and other Chemical Physics problems. We will discuss the role of
trust-region schemes, line searches, linearly and nonlinearly
constrained optimization, Inexact Restoration and SQP methods and the
type of convergence theories that may be useful in
order to explain the practical behavior of the methods. Emphasis will
be given on general principles instead of algorithmic details.
In the quantum mechanical treatment of molecules we use the Born-Oppenheimer (adiabatic) approximation, in which the motion of nuclei and electrons is separated. In this approximation the coupling between different electronic states is neglected and nuclei move on a single electronic potential energy surface. Nevertheless, non-adiabatic processes where the coupling between different electronic states becomes large and important. These processes are facilitated by the close proximity of potential energy surfaces, and especially by the extreme case where the potential energy surfaces become degenerate forming conical intersections. Modeling non-adiabatic processes requires accurate calculation of electronic structure states and their coupling. Methods for calculating excited states, the non-adiabatic couplings and conical intersections will be discussed.
Standard approximations for the exchange-correlation functional have been found to give big errors for the linearity condition of fractional charges, leading to delocalization error, and the constancy condition of fractional spins, leading to static correlation error. These two conditions are now unified for states with both fractional charge and fractional spin: the exact energy functional is a plane, linear along the fractional charge coordinate and constant along the fractional spin coordinate with a line of discontinuity at the integer. This sheds light on the nature of the derivative discontinuity and illustrates the need for a discontinuous functional of the orbitals or density. This is key for the application of DFT to strongly correlated systems.
The quantum mechanical - molecular mechanical (QM/MM) potential energy
functions are used in calculations of the
potential of mean force (PMF) following the conventional molecular
dynamics (MD) based procedure. Constant temperature
MD simulations, in particular, allowing for rigid-body MD algorithms, are
performed for the canonical (NVT) ensemble
in conjunction with the Nose-Poincare thermostat. The umbrella sampling
technique and the weighted histogram analysis
method are applied for PMF calculations. Two versions of the QM/MM method,
namely, the mechanically embedded cluster
technique and the flexible effective fragment approach are considered for
potential energy estimates. The QM subsystem
can be described by various electronic structure approximations including
the DFT and multiconfigurational CASSCF
methods. Conventional force field parameters can be employed in the MM
subsystem. The computer program utilizes the PC
GAMESS (A. Granovsky) and TINKER (J. Ponder) molecular modeling packages.
The first application considered the
mechanisms of proton conduction in the gramicidin A ion channel. The chain
of nine water molecules inside the channel
constituted the QM part described by the B3LYP/6-31G* approximation. The
peptide walls of the gramicidin channel and
two clusters of 20 water molecules placed at both ends of the channel
constituted the MM subsystem described by the
AMBER force field parameters. For the energy consuming stage of the water
file reorientation inside the channel we
calculated the activation free energy barrier of 7.7 kcal/mol at 300K as
compared to 6.5 kcal/mol in experimental
studies.
The exponential localization of Wannier functions in two or three dimensions is proven for all insulators that display time-reversal symmetry, settling a long-standing conjecture. The proof make use of geometric techniques, which also imply that Chern insulators cannot display exponentially localized Wannier functions.
Finally, a new algorithm to explicitly construct the exponentially localized Wannier functions is exhibited.
Principle Collaborator:
Natarajan Sukumar
(University of California, Davis)
Over the past few decades, the planewave (PW) pseudopotential method has established itself as the dominant method for large, accurate, density-functional calculations in condensed matter. However, due to its global Fourier basis, the PW method suffers from substantial inefficiencies in parallelization and applications involving highly localized states, such as those involving 1st-row or transition-metal atoms, or other atoms at extreme conditions. Modern real-space approaches, such as finite-difference (FD) and finite-element (FE) methods, can address these deficiencies without sacrificing rigorous, systematic improvability but have until now required much larger bases to attain the required accuracy. Here, we present a new real-space FE based method which employs modern partition-of-unity FE techniques to substantially reduce the number of basis functions required, by building known atomic physics into the Hilbert space basis, without sacrificing locality or systematic improvability. We discuss pseudopotential as well as all-electron applications. Initial results show order-of-magnitude improvements relative to current state-of-the-art PW and adaptive-mesh FE methods for systems involving localized states such as d- and f-electron metals and/or other atoms at extreme conditions.
This work performed under the auspices of the U.S. Department of Energy by Lawrence Livermore National Laboratory under Contract DE-AC52-07NA27344.
Simulation of ballistic transport in nanodevices are usually computationally intensive. To improve the efficiency of this simulation, the subband decomposition method approximates the solution $psi(x_1,x_2)$ by $sum_{i=1}^{n_e} varphi_i(x_1) xi_i(x_2;x_1)$, where $xi_i(x_2;x_1)$, $1 leq i leq n_e$ are solutions to an eigenvalue problem in the $x_2$-direction and parameterized by $x_1$; $varphi_i(x_1)$ are solutions to a Schr"{o}dinger equation with open boundary conditions; and $(x_1,x_2)$ denotes a point in the simulation domain $Omega$. However, determination of $xi_i(cdot;x_1)$ based on, say, finite element method, at all grid points in the $x_1$ direction can still be expensive. Here we propose approximating $xi_i(cdot;x_1)$ using the reduced basis method. By exploiting {em a posteriori} error estimators and greedy sampling algorithm, we can construct very efficient reduced basis approximation space for $xi_i(cdot;x_1)$. The computational cost only grows marginally when mesh spacing in the $x_2$-direction decreases, compared to exponential increase for a finite element approximation.
The exchange-correlation energy of Kohn-Sham density functional theory can be written as the integral over all space of an exchange-correlation energy density, which is a function of various density-dependent ingredients. Adding ingredients can produce approximations that satisfy more exact constraints, or fit data better. In our vision, the ladder has five rungs: (1) the local spin density approximation (LSDA), which uses only the local electron spin densities as ingredients, (2) the generalized gradient approximation (GGA), which adds the density gradients, (3) the meta-GGA, which adds the orbital kinetic energy densities, (4) the hyper-GGA, which adds the fully-nonlocal exact exchange energy densities, and (5) the generalized random phase approximation, which adds the unoccupied Kohn-Sham orbitals. All rungs but the fourth have been constructed without empiricism. Some recent developments will be sketched: (a) The bifurcation of the second rung into standard GGA's and GGA's for solids. (b) Possible refinement of the meta-GGA by recovering the gradient expansion for exchange over a wide range of density gradients, as in the PBEsol GGA for solids. Since meta-GGA is not much more expensive than GGA, and is potentially much more accurate for systems near equilibrium, an improved meta-GGA could replace LSDA or GGA in applications. (c) A hyper-GGA that interpolates between meta-GGA exchange
(in normal regions where the errors of meta-GGA exchange and correlation tend to cancel) and exact exchange (in abnormal regions where no such cancellation is possible). Extensive abnormal regions can occur in open subsystems connected by stretched bonds.
The first step in electronic-structure calculations is geometry
optimization: finding an atomic configuration that minimizes
the energy. A popular and successful model for the energy is the
total-energy functional from density-functional theory. However,
the evaluation of this functional at a given geometry is computationally
expensive. Therefore, standard minimization techniques may be costly
in practice. We propose a new approach for geometry optimization based
on surrogate modeling. We describe our approach and present preliminary
results.
Joint work with Filipp Furche.
We study the ability of TD-DFT to correctly predict the
splitting of
closely spaced singlet excited states in a system of two
spatially
separated chromophores. We find that functionals without at
least a
certain fraction of Hartree-Fock exchange kernel fare very
poorly at this
task because of a known problem these functionals have with the
underestimation of charge-transfer excitation energies.
This presentation will address our current efforts to develop
more accurate exchange-correlation forms
for density functional theory. There are two leading themes in
our current work: range separation and
local weights. On the first theme, I will present a three-range
hybrid functional and discuss the
rationale for the success of screened functionals like HSE and
LC-wPBE. On the second theme, the
emphasis will be on new metrics for local hybridization and
local range separation. Much of the focus
will be on the seemingly dissimilar needs between solids and
molecules, and on the computational
challenge of including nonlocal (Hartree-Fock type) exchange
efficiently in condensed systems.
A viable methodology for the exact analytical solution of the multiparticle Schrodinger and Dirac equations has
long been considered a holy grail of theoretical chemistry. Since a benchmark work by Torres-Vega and Frederick in the
1990's[1], the Quantum Phase Space Representation (QPSR) has been explored as an alternate method for solving
various physical systems, including the harmonic oscillator[2], Morse oscillator[3], one-dimensional hydrogen atom[4], and
classical Liouville dynamics under the Wigner function[5]. QPSR approaches are particularly challenging because of the
complexity of phase space wave functions and the fact that the number of coordinates doubles in the phase space
representation. These challenges have heretofore prevented the exact solution of the multiparticle equation in phase
space.
Recently, Simpao* has developed an exact analytical symbolic solution scheme for broad classes of differential
equations utilizing the Heaviside Operational Ansatz (HOA). It is proposed to apply this novel methodology to QPSR
problems to obtain exact solutions for real chemical systems and their dynamics. In his preliminary work, Simpao* has
already applied this method to a number of simple systems, including the harmonic oscillator, with solutions in agreement
to those obtained by Li [refs.2,3,4,5]. He has also demonstrated the exact solution to the radial Schrodinger Equation for
an N-particle system with pairwise Coulomb interaction**. In addition to the Schrodinger Equation, the HOA method is
capable of treating the Dirac equation*** as well as differential systems governing both relativistic and non-relativistic
particle dynamics.
Applying these methods would allow us to pursue further exploration of this methodology, starting with the exact
solution of multielectron atoms and moving toward complex molecules and reaction dynamics. It is believed that the
coupling of HOA with QPSR represents not only a fundamental breakthrough in theoretical physical chemistry, but it is
promising as a basis for exact solution algorithms that would have tremendous impact on the capabilities of computational
chemistry. As the theoretical foundation for spectroscopy is the Schrodinger equation, the significance of this discovery to
the enhanced analysis of spectroscopic data is obvious. For example, the analysis of the Compton line in momentum
spectroscopy necessitates the consideration of the momentum wavefunction for the molecular system under study. The
novel methods *,**,*** allow the exact determination of the momentum[and configuration] space wavefunction from the
QPSR wavefunction by way of a Fourier Transform. For example,the primariy focus of the PREPRINT ** is the pairwise
1/rij interaction in context of the radial equation in the nonrelativistic Schrodinger case. This application of the exact
solution ansatz developed above corresponds to the problem of n-particles with pairwise Coulomb interaction;scaling the
parameters and variables of the problem yields the exact solution of the QPSR Schrodinger equation for the first-principles
general polyatomic molecular Hamiltonian. Upon a straightforward slight adaptation of this non-relativistic Schrodinger
result, the QPSR Dirac equation addressed in *** immediately yields the relativistic counterpart for the first -principles
general polyatomic molecular Hamiltonian solution.
1. Torres-Vega, G. and J.H. Frederick, A quantum-mechanical representation in phase space. Journal of Chemical
Physics, 1993. 98(4): p. 3103-20.
2. Li, Q.S. and J. Lu, Rigorous solutions of diatomic molecule oscillator with empirical potential function in phase
space. Journal of Chemical Physics, 2000. 113(11): p. 4565-4571.
3. Hu, X.-G. and Q.S. Li, Morse oscillator in a quantum phase-space representation: rigorous solutions. Journal of
Physics A: Mathematical and General, 1999. 32(1): p. 139-146.
4. Li, Q.S. and J. Lu, One-dimensional hydrogen atom in quantum phase-space representation: rigorous solutions.
Chemical Physics Letters, 2001. 336(1,2): p. 118-122.
5. Li, Q.S., G.M. Wei, and L.Q. Lu, Relationship between the Wigner function and the probability density function in
quantum phase space representation. Physical Review A: Atomic, Molecular, and Optical Physics, 2004. 70(2):
p. 022105/1-022105/5.
* Electronic Journal of Theoretical Physics,1 (2004), 10-16
** PREPRINT Toward Chemical Applications of Heaviside Operational Ansatz:
Exact Solution of Radial Schrodinger Equation for Nonrelativistic N-particle System with Pairwise
1/rij Radial Potential in Quantum Phase Space[now published MAY 2008 Journal of Mathematical Chemistry]
*** Electronic Journal of Theoretical Physics, 3, No. 10 (2006) 239-247
http://www.springerlink.com/content/225x523327771420/?p=a04f3c2e1352400b87bde6ac7331e4b2π=9
Structures and binding in small water-benzene complexes (1-8 water molecules
and 1-2 benzene molecules) are studied using the general effective fragment
potential (EFP) method. The lowest energy conformers of the clusters were
found using a Monte-Carlo technique. The EFP method accurately predicts
structures and binding energies in the water-benzene complexes. Benzene is
polarizable and consequently participates in hydrogen bond networking of
water. Since the water-benzene interactions are only slightly weaker than
water-water interactions, structures with different numbers of water-water,
benzene-water, and benzene-benzene bonds often have very similar binding
energies. This is a challenge for computational methods.
Conventional density functional approximations for the
exchange-correlation energy fail to describe an important
class of systems formed by the van der Waals (vdW) interaction,
because they are unable to account for the long-range part of the
vdW interaction, while they may describe the short-range part well.
Here we first propose a density functional to simulate the coefficient
C_6 of the leading term of the long-range part. Then we
construct a nonempirical vdW-corrected meta-GGA functional by properly building
the long-range part into a sophiscated meta-generalized gradient
approximation (meta-GGA). Numerical tests on diverse atom pairs show that
the proposed C_6 model is remarkably accurate. Applications of the
vdW-corrected meta-GGA functional to rare-gas dimers show that the binding
energy curves and bond lengths obtained with the vdW-corrected meta-GGA
are well improved over those with the original meta-GGA, and agree fairly
well with experiments.
Full conﬁguration interaction (FCI) data are used to quantify the accuracy of approximate adiabatic connection (AC) forms in describing two challenging problems in density functional theory—the singlet ground state potential energy curve of H2 in a restricted formalism and the energies of the helium isoelectronic series, H− to Ne8+. For H2, an exponential-based form yields a potential energy curve that is virtually indistinguishable from the FCI curve, eliminating the unphysical barrier to dissociation observed previously with a [1,1]-Padé-based form and with the random phase approximation. For the helium isoelectronic series, the Padé-based form gives the best overall description, followed by the exponential form, with errors that are orders of magnitude smaller than those from a standard hybrid functional. Particular attention is paid to the limiting behavior of the AC forms with increasing bond distance in H2 and increasing atomic number in the isoelectronic series; several forms describe both limits correctly. The study illustrates the very high quality results that can be obtained using exchange-correlation functionals based on simple AC forms, when near-exact data are used to determine the parameters in the forms.
This lecture reports on work carried out in collaboration with
Yan Zhao.
We have developed a suite of density functionals. All four
functionals are accurate for noncovalent interactions and
medium-range correlation energy. The functional with broadest
capability, M06, is uniquely well suited for good performance
on both transition-metal and main group-chemistry; it also
gives good results for barrier heights. Another functional,
M06-L has no Hartree-Fock exchange, which allows for very fast
calculations on large systems, and it is especially good for
transition-metal chemistry and NMR chemical shieldings. M08-2X
and an earlier version, M06-2X, have the very best performance
for main-group thermochemistry, barrier heights, and
noncovalent interactions. M06-HF has no one-electron
self-interaction error and is the best functional for charge
transfer spectroscopy. A general characteristic of the whole
suite is the optimized inclusion of kinetic energy density and
higher separate accuracy of medium-range exchange and
correlation contributions with less cancellation of errors than
previous functionals [1-4]; for example, the functionals are
compatible with a range of Hartree-Fock exchange and, although
one or another of them may be more highly recommended for one
or another property or application, all four are better on
average than the very popular B3LYP functional. A few sample
applications, including catalytic systems [5,6] and
nanomaterials [7], will also be discussed. Recent work on
lattice constants, band gaps, and an improved version of M06-2X
will be summarized if time permits.
[1] "Design of Density Functionals by Combining the Method of
Constraint Satisfaction with Parametrization for
Thermochemistry, Thermochemical Kinetics, and Noncovalent
Interactions," Zhao, Y. ; Schultz, N. E.; Truhlar, D. G.; J.
Chem. Theory Comput. 2006, 2, 364-382.
[2] "A New Local Density Functional for Main Group
Thermochemistry, Transition Metal Bonding, Thermochemical
Kinetics, and Noncovalent Interactions," Zhao, Y.; Truhlar, D.
G. J. Chem. Phys. 2006, 125, 194101/1-18.
[3] “The M06 Suite of Density Functionals for Main Group
Thermochemistry, Thermochemical Kinetics, Noncovalent
Interactions, Excited States, and Transition Elements: Two New
Functionals and Systematic Testing of Four M06-Class
Functionals and 12 Other Functionals,” Zhao, Y.; Truhlar, D. G.
Theor. Chem. Acc. 2008, 120, 215-241.
[4] "Density Functionals with Broad Applicability in
Chemistry," Zhao, Y.; Truhlar, D. G. Acc. Chem. Res. 2008 41,
157-167.
[5] “Attractive Noncovalent Interactions in Grubbs
Second-Generation Ru Catalysts for Olefin Metathesis," Zhao,
Y.; Truhlar, D. G. Org. Lett. 2007, 9, 1967-1970.
[6] "Benchmark Data for Interactions in Zeolite Model Complexes
and Their Use for Assessment and Validation of Electronic
Structure Methods," Zhao, Y.; Truhlar, D. G. J. Phys. Chem. C
2008, 112, 6860-6868.
[7] "Size-Selective Supramolecular Chemistry in a Hydrocarbon
Nanoring," Zhao, Y.; Truhlar, D. G. J. Am. Chem. Soc.2007, 129,
8440-8442.
Constrained density functional theory methods (C-DFT) have been developed by several groups [1,2] to construct localized charge or spin states. At the same time, others have suggested ways of constructing variable-charge potential energy surfaces from certain valence-bond (VB) states [3,4]. In general, mainstream electronic structure codes do not generate the desired valence-bond states, and does not mix them in the manner needed to construct variable-charge potentials. C-DFT does provide estimates of the desired, individual valence-bond states and energies, but does not directly give the state mixing. While there are methods to overcome this limitation [5], here we explore an alternative.
That alternative, notionally suggested in Ref. 3, is to deduce state-mixing from C-DFT with fractional charges. One assumes that the normal integer-charge VB states remain the same when describing a state composed of atoms with fractional charges. For instance, we use C-DFT to calculate the energies for the three charge distributions for a geometrically symmetric water molecule. Two of the constrained VB states are integer-charge distributions, one with all atoms being neutral and a second with the oxygen integer being anionic and one of the hydrogens cationic. The third VB state places half of a negative and positive charge on the oxygen and one hydrogen, respectively. We present our initial efforts to calculate the variational energy surface for a water molecule at several bond-lengths along the symmetric stretch vibrational mode. To capture the fractional charge behavior of the energies, the B3LYP functional is employed [5-7].
[1] Q Wu and T Van Voorhis, Phys. Rev. A 72, 024502 (2005).
[2] J Behler, B Delley, K Reuter, and M Scheffler, Phys. Rev. B 75, 115409 (2007).
[3] SM Valone and SR Atlas, J. Chem. Phys. 120, 7262 (2004).
[4] SM Valone, J Li, and S Jindal, IJQC 108, 1452 (2008).
[5] Q Wu and T Van Voorhis, J. Chem. Phys. 125, 164105 (2006).
[6] AD Becke, J. Chem. Phys. 98, 5648 (1993).
[7] C Lee, W Yang, and RG Parr, Phys. Rev. B 37, 785 (1988).
Nearly all common density functional approximations fail to adequately
describe dispersion interactions responsible for binding in van der Waals
complexes. One of the most promising new methods is the nonlocal van der
Waals density functional (vdW-DF) of Ref. [1], which was derived from
first principles, describes dispersion interactions in a seamless fashion,
and yields the correct asymptotics. Recently we reported a self-consistent
implementation of vdW-DF with Gaussian basis functions [2]. Our code
includes analytic gradients of the energy with respect to nuclear
displacements, enabling efficient geometry optimizations. vdW-DF tends to
overbind molecular complexes, especially if used in combination with
Hartree-Fock exchange. We propose a slightly simplified construction of
the nonlocal vdW-DF correlation functional, for which we also derive a
semilocal gradient correction. This correction reduces the overbinding
tendency and improves the accuracy of vdW-DF. Adjusting an empirical
parameter in the semilocal part, vdW-DF can be made compatible with
different exchange approximations.
[1] M. Dion, H. Rydberg, E. Schroder, D.C. Langreth, and B.I. Lundqvist,
Phys.Rev.Lett. 92, 246401 (2004).
[2] O.A. Vydrov, Q. Wu, and T. Van Voorhis, J.Chem.Phys. 129, 014106 (2008).
A new method for calculating effective atomic radii within the generalized Born (GB) model of implicit solvation is proposed, for use in computer simulations of bio-molecules. First, a new formulation for the GB radii is developed, in which smooth
kernels are used to eliminate the divergence in volume integrals intrinsic in the model. Next, the Fast Fourier Transform (FFT)
algorithm is applied to integrate smoothed functions, taking advantage of the rapid spectral decay provided by the smoothing. The total cost of the proposed algorithm scales as O(N^3logN+M) where M is the number of atoms comprised in a molecule, and N is the number of FFT grid points in one dimension, which depends only on the geometry of the molecule and the spectral decay of the smooth kernel but not on M. To validate our algorithm, numerical tests are performed for three solute models: one spherical object for which exact solutions exist and two protein molecules of differing size. The tests show that our algorithm is able to reach the accuracy of other existing GB implementations, while offering much lower computational cost.
I will present a direct constrained minimization (DCM) algorithm
for solving the Kohn-Sham equations. The key ingredients of this
algorithm involve projecting the Kohn-Sham total energy functional
into a sequences of subspaces of small dimensions and seeking the
minimizer of total energy functional within each subspace. The
minimizer of a subspace energy functional not only provides a
search direction along which the KS total energy functional decreases
but also gives an optimal ``step-length" to move along this search
direction. I will provide some numerical examples to demonstrate
the efficiency and accuracy of this approach and compare it
with the widely used method of self-consistent field (SCF) iteration.
I will also discuss a few other numerical issues in algorithms
designed to solve the Kohn-Sham equations.