Mathematical methods in quantum chemistry. Abstracts from the workshop held March 18–24, 2018. (English) Zbl 1409.00078

Summary: The field of quantum chemistry is concerned with the modelling and simulation of the behaviour of molecular systems on the basis of the fundamental equations of quantum mechanics. Since these equations exhibit an extreme case of the curse of dimensionality (the Schrödinger equation for \(N\) electrons being a partial differential equation on \(\mathbb R^{3N}\)), the quantum-chemical simulation of even moderate-size molecules already requires highly sophisticated model-reduction, approximation, and simulation techniques. The workshop brought together selected quantum chemists and physicists, and the growing community of mathematicians working in the area, to report and discuss recent advances on topics such as coupled-cluster theory, direct approximation schemes in full configuration-interaction (FCI) theory, interacting Green’s functions, foundations and computational aspects of density-functional theory (DFT), low-rank tensor methods, quantum chemistry in the presence of a strong magnetic field, and multiscale coupling of quantum simulations.


00B05 Collections of abstracts of lectures
00B25 Proceedings of conferences of miscellaneous specific interest
81Q05 Closed and approximate solutions to the Schrödinger, Dirac, Klein-Gordon and other equations of quantum mechanics
81V55 Molecular physics
65Y20 Complexity and performance of numerical algorithms
81-06 Proceedings, conferences, collections, etc. pertaining to quantum theory
Full Text: DOI


[1] S. Kvaal, Manuscript in preparation (2018).
[2] J.S. Arponen, Variational Principles and Linked-Cluster exp S Expansions for Static and Dynamic Many-Body Problems, Ann. Phys., 151, pp. 311–382, (1983).
[3] P.-O. L¨owdin, On the stability problem of a pair of adjoint operators, J. Math. Phys., 24, pp. 70–87, (1982).
[4] K. Schm¨udgen, Unbounded Self-adjoint Operators on Hilbert Space, Graduate Texts in Mathematics 265, Springer, (2012).
[5] E. Zeidler, Nonlinear Functional Analysis and its Application II/B, Springer, (1990).
[6] A. Laestadius and S. Kvaal, Analysis of the extended coupled-cluster method in quantum chemistry, SIAM J. Numer. Anal., 56, pp. 660–683, (2018). Numerical Methods for Solving Coupled Cluster Equations Chao Yang (joint work with Jiri Brabec, Jinmei Zhang, Karol Kowalski and Edward Valeev) We examine two different approaches for accelerating numerical solutions of the couple cluster equations which are known to be nonlinear and highly complex. The first approach is based on the Jacobian-free Newton-Krylov method. We show that method tends to be more stable than the widely used DIIS method when the Jacobian of the coupled clustered equation is not dominated by a diagonal matrix consisting of the differences between the Hartree-Fock occupied and virtual orbital energies. We discuss the possibility of combining the Newton-Krylov method with DIIS method, although this combination does not seems to lead to a significant improvement in convergence rate, at least not for problems that are well behaved. We point out the importance of using level-shifting to regularize diagonal preconditioner consisting of HF orbital energy differences and stabilize the convergence of the Newton-Krylov method for problems that have nearly degenerate orbital energies. We also discuss the possibility of using alternative preconditioners that are based on the construction of a small active space approximation to the coupled cluster equation-of-motion Hamiltonian. In the second approach, we discuss methods that try to exploit sparsity and low rank structures of the coupled cluster amplitude. In our earlier work [1], we developed a practical numerical scheme to sparsify the correction to coupled cluster amplitudes in an inexact Newton iteration. We showed that 90
[7] J. Brabec, C. Yang, E. Epifanovsky, A.I. Krylov and E. Ng, Reduced-cost sparsityexploiting algorithm for solving coupled-cluster equations, Journal of Computational Chemistry, 37 (2016), pp. 1059–1067.
[8] M. Clement, J. Zhang, C. A. Lewis, C. Yang and E. F. Valeev, Assessment of Optimized Natural Orbitals for the Coupled Cluster Methods, submitted, (2018).
[9] C. Edmiston and M. Krauss, Pseudonatural Orbitals as a Basis for the Superposition of Configurations. I. He2+, J. Chem. Phys., 45, (1966), 1833.
[10] W. Meyer, PNO–CI Studies of electron correlation effects. I. Configuration expansion by means of nonorthogonal orbitals, and application to the ground state and ionized states of methane, J. Chem. Phys. 58, (1973), 1017.
[11] F. Neese, F. Wennmohs and A. Hansen, Efficient and accurate local approximations to coupled-electron pair approaches: An attempt to revive the pair natural orbital method., J. Chem. Phys., 130, (2009), 114108. Log-scaling tight-binding for material defects Christoph Ortner (joint work with Simon Etter) Intuitively, there is very limited “information content” in the electronic structure of a localised defect embedded in a homogeneous medium, independent of or slowly growing with the medium’s size. It should be possible to exploit this in order to develop highly efficient electronic structure methods for crystalline defects. To formalise this idea, we consider a model, developed in [1, 2, 3, 4], for point defects and straight dislocation lines embedded in a homogeneous host crystal, taking into account both mechanical relaxation and electronic relaxation within the tight-binding model. For the sake of simplicity, this note will only describe the case of an impurity in Zd, with a 2-centre s-orbital tight-binding model. More specifically, if (yℓ)ℓ∈Ω⊂ Rdis a finite collection of atoms then we assume that a tight-binding hamiltonian is given by H(y) = (h(rij))i,j∈Ωwith rij= |yi− yj| which is exponentially localised (h(r) . e−γrijfor some γ > 0). The potential energy landscape for the tight-binding model is then given by XN E(y) =f (εs),where H(y)ψs= εsψs, s=1 and f (ε) is an analytic function describing either the canonical or grand-canonical ensemble for the electrons, see [3] for more detail. Mathematical Methods in Quantum Chemistry643 It is then shown in [3, 2] that the problem of equilibrating E(y) under suitable boundary conditions has a limit as the domain Ω approaches Zdand that limit equilibrium displacements ¯u : Zd→ Rd, u = y− id satisfy (1) Dj¯u(ℓ) ≤ Cj|ℓ|1−d−j where Df (ℓ) = (f (ℓ + ei)− f(ℓ))di=1denotes a discrete gradient operator. Indeed, we may expect (though it is an open problem whether this is true) that the stronger estimate (2) Dju(ℓ)¯ ≤ Cj!|ℓ|1−d−j holds. Our aim is to exploit this a priori known regularity of the equilibrium displacements to construct highly efficient numerical schemes for solving the electronic structure problem. To that end, as a first step, we establish an analogous regularity result for the density matrix, given in the finite domain case by XN Γ :=fFD(εs) ψs⊗ ψs, s=1 where fFis the Fermi-Dirac function. In the thermodynamic limit it is more conveniently written as 1I 2πifFD(z)(z− H)−1dz, where the integration is taken along a contour that circles the spectrum of H(y) but does not circle any singularities of fFD. Establishing regularity of electronic structure can be reduced to regularity of Γ and hence of (z− H)−1. The following results use identities such as D(z− H)−1=−(z − H)−1(DH)S(z− H)−1, where (SekH)ij= Hi+ek,j+ekand DekH = SekH− H, as well as Coombe-Thomas estimates (see [4] and references therein) to bound off-diagonal decay of (z− H)−1. The regularity (1) then implies that (unpublished) (DkΓ)ij ≤ Ck′e−γ′krij|i|1−d−k+|j|1−d−k. Assuming the stronger regularity (2), we expect that (in preparation) (3) (DkΓ)ij ≤ Ck!e−γrij|i|1−d−k+|j|1−d−k. To exploit this regularity result, we propose to approximate each diagonal of the density matrix by a piece-wise polynomial (hp-FEM type approximation scheme). A little more precisely, if a uniform polynomial degree p is chosen and the finite element mesh size scales as h(x)∼ |x|, and is truncated at a suitable radius then under the regularity (1) it is relatively straightforward to establish a superalgebraic best-approximation rate kΓ − ΓhpkF≤ CtDOF−tfor all t > 0. 644Oberwolfach Report 13/2018 It turns out, however, that standard (continuum) polynomial spectral approximation error estimates do not immediately apply in the discrete setting. We aim to establish that (3) implies  (4)kΓ − ΓhpkF≤ C exp − cDOF1/r, for some c, r > 0. Proving this result is work in progress, however, there is already overwhelming theoretical and numerical evidence for its validity. The final ingredient, if (4) can be confirmed, is to develop an algorithm that efficiently computes Γhp. A Galerkin discretisation for Γhpcan achieve this in polynomial cost, which finally leads to a poly-logarithmic scaling result,  COST Γhp.logqkΓ − ΓhpkFfor some q > 0. References
[12] V. Ehrlacher, C. Ortner and A. V. Shapeev, Analysis of Boundary Conditions for Crystal Defect Atomistic Simulations, Arch. Ration. Mech. Anal., 222, pp. 1217–1268, (2016). doi:10.1007/s00205-016-1019-6
[13] H. Chen, F. Q. Nazar and C. Ortner, Geometry Equilibration of Crystalline Defects in Quantum and Atomistic Descriptions, http://arxiv.org/abs/1709.02770, (2017).
[14] H. Chen, J. Lu and C. Ortner, Thermodynamic Limit of Crystal Defects with Finite Temperature Tight Binding, to appear in Arch. Ration. Mech. Anal.
[15] H. Chen and C. Ortner, QM/MM Methods for Crystalline Defects. Part 1: Locality of the Tight Binding Model, Multiscale Model. Simul., 14, pp. 232–264, (2016). doi:10.1137/15M1022628 Current-density-functional theory for molecules in strong magnetic fields Andrew M. Teale (joint work with Tom J. P. Irons, James Furness, Erik I. Tellgren, Stella Stopkowicz and Trygve Helgaker) Recent progress in the theoretical development and implementation of non-perturbative current-density-functional theory (CDFT) using London atomic orbitals is reviewed. These calculations enable study of chemical systems in strong magnetic fields up to∼1 a.u. ≈ 235000 Tesla. Introducing a uniform magnetic field, B = ∇ × A described by a vector potential A leads to the electronic Hamiltonian (1) H(v, A) =p2i+|ri− rj|−1+v(ri)+A(ri)· pi+1XA2(r 22i). ii>jiii To setup the foundations of CDFT it is convenient to re-parameterize this Hamiltonian by collecting the scalar potential contributions as u = v +12A2, 1XXXX (2)H(u, A) =p2i+|ri− rj|−1+u(ri) +A(ri)· pi. 2 ii>jii Mathematical Methods in Quantum Chemistry645 This has the convenient consequence that the energy E(u, A) is concave in u and A (whereas E(v, A) is not concave) and so a formulation of CDFT analogous to Lieb’s formulation for DFT may be constructed [1]. A Kohn-Sham (KS) CDFT scheme can then be setup, with the KS equations,  (3)2p2+12{p, As} + us+ s· [∇ × As]ϕp= εpϕp where the scalar potential and vector potentials are 1 (4)2A2ext+ vJ+ vxcAs= Aext+ Axc. To provide an efficient computational implementation of CDFT we expand the KS orbitals in terms Gaussian type basis functions XKa (5)φa(r) = (x− Ax)ax(y− Ay)ay(z− Az)azdke−αk|r−A|2 k=1 and attach a field dependent phase factor to yield London atomic orbitals (LAOs), (6)ωa(r) = φa(r)e−2iB×(A−O)·r. The use of LAOs ensures gauge-origin independent energies are obtained. In this talk we will examine two key aspects to allow the practical application of KS-CDFT; (i) techniques for the efficient evaluation of molecular integrals over LAOs, (ii) the nature of approximate exchange-correlation functionals in CDFT, benchmarking them against ab initio methods such as coupled-cluster theory [2]. Implementation of these methods in London [3] and Quest [4] is outlined. Molecular Integrals The most time consuming step in practical calculations is the evaluation of the electron repulsion integrals (ERIs) over the atomic orbitals. For LAOs the importance of this step in determining the efficiency of the calculations is magnified by the fact that complex arithmetic is required and the permutational symmetry is reduced from eightfold (for real orbitals) to fourfold (for the complex LAOs). The ERIs over the LAOs take the form Z Zω∗ (ab|cd) =a(r1)ωb(r1)ωc∗(r2)ωd(r2)dr (7)|r1− r2|1dr2 and the most efficient implementations previously have used the McMurchieDavidson approach. This approach gives recursions for the evaluation of the ERIs involving two-centre intermediates of the form (see Ref. [5] for details of notation), [r + 1i](m)= PQi[r](m+1)+ ri[r− 1i](m+1) (8)[p|q] = (−1)q[p + q](0) 646Oberwolfach Report 13/2018 when applied to Gaussian type atomic orbitals of the type in Eq. (5). When applied to LAOs these recursions are significantly more complicated, [p + 1i|q](m)=− iχP,i[p|q](m)+ PQi[p|q](m+1) (9)+ pi[p− 1i|q](m+1)− qi[p|q − 1i](m+1) which leads to an intrinsic increase in computational cost of the McMurchieDavidson scheme using LAOs. In both cases the four centre ERI is then constructed via the recursion 1 (10)[a + 1ib p| = pi[a b p− 1i| + PAi[a b p| +[a b p + 1i| . 2ζ This intrinsic increase in computational effort for the algorithm coupled with the need to calculate more integrals when dealing with LAOs has led us to consider other approaches to molecular integral evaluation. Two of the most widely used schemes for evaluation of integrals over Gaussiantype atomic orbitals are the Head-Gordon-Pople (HGP) and Rys quadrature methods. For HGP the conventional recursion relations for the two-centre intermediates,  [e + 1i0|f0](m)= PAi[e0|f0](m)− PQi1[e0|f0](m+1) 2ζ 1  2ζ2ζi0|f0](m+1) 1 1 (11)+ fi[e0|f − 1i0](m+1) 2ζ2η carry through essentially unaltered for LAOs, meaning that the intrinsic complexity of the algorithm is not increased. Four centre ERIs are then obtained via, (12)(ab + 1i| = (a + 1ib| + ABi(ab| . For Rys quadrature the algorithm is also essentially unchanged (though the computation of the required roots and weights requires special care in the complex case—see Ref. [5]). In this scheme the integrands are determined recursively via  Ii(e + 1, f ; λ) =PAi−ηt2λPQ ζ + ηiIi(e, f ; λ) ηt2f t2 e (13)2ζ1−ζ + ηλIi(e− 1, f; λ) +2(ζ + η)λIi(e, f− 1; λ) and then summed over the Rys nodes XN (14)[e0|f0] =Ix(ex, fx; λ)Iy(ey, fy; λ)Iz(ez, fz; λ) wλ λ=1 and the final four-centre ERI is recovered using the horizontal recursion relation employed in the HGP scheme. Mathematical Methods in Quantum Chemistry647 As expected none of the above integral schemes is optimal for all types of ERI, in particular when integrals are classified according to their total angular momenta it becomes clear that the McMurchie-Davidson scheme remains most efficient for integrals over low angular momenta basis functions, whilst HGP is more efficient in the intermediate regime and Rys quadrature is most efficient for high angular momenta. I will discuss our recent implementation of a mixed scheme, illustrated in Fig. 1, in which the most appropriate algorithm is selected for each batch of integrals on the fly. The implementation of density-fitting (aka resolution of the identity) approaches will also be discussed. Figure 1.Schematic for selection of most appropriate integral evaluation scheme for LAO ERIs. The most efficient approach can be selected on-the-fly for batches of integrals classified according to their total angular momentum L. Exchange-Correlation functionals The outstanding challenge for realization of a practical scheme for CDFT calculations is the development of appropriate exchange-correlation functionals depending on the paramagnetic current density, iX (15)jp=−[φ∗i∇φi− φi∇φ∗i]. 2 i In previous work we have shown that simple vorticity based models fail catastrophically in strong magnetic fields [6]. However, we have recently demonstrated that meta-GGA class functionals, utilizing a modified kinetic energy density, (16)|jp|2, τ =⇒ ˜τ = τ − ρ show remarkably accurate performance in strong magnetic fields in comparison with high accuracy coupled-cluster results [7]. I will discuss how this performance 648Oberwolfach Report 13/2018 can be rationalized by considering the importance of 1-electron regions and visualized by an extended [8] definition of the electron localization function (ELF) 1[τ (r)− τvW(r)] (17)1 + α(r)2,α(r) =τUEG(r) with τvW(r) and τUEG(r) the von Weis¨acker and uniform electron gas definitions of the kinetic energy density, respectively. Finally, I will outline some outstanding challenges for CDFT calculations, including the accurate description of weak-field properties of general chemical interest such as nuclear magnetic resonance shielding constants [9] and the description of high spin ground states, which become prevalent as the strength of the magnetic field considered increases. References
[16] E. I. Tellgren, S. Kvaal, E. Sagvolden, U. Ekstr¨om, A. M. Teale and T. Helgaker, Choice of basic variables in current-density-functional theory, Phys. Rev. A, 86 (2012), 062506.
[17] S. Stopkowicz, J. Gauss, K. K. Lange, E. I. Tellgren and T. Helgaker, Coupledcluster theory for atoms and molecules in strong magnetic fields, J. Chem. Phys., 143 (2015), 074110.
[18] LONDON, a quantum-chemistry program for plane-wave/GTO hybrid basis sets and finite magnetic field calculations. By E. Tellgren (primary author), T. Helgaker, A. Soncini, K. K. Lange, A. M. Teale, U. Ekstr¨om, S. Stopkowicz, J. H. Austad and S. Sen. See londonprogram.orgfor more information.
[19] QUEST, a rapid development platform for Quantum Electronic Structure Techniques, 2017, http://quest.codes
[20] T. J. P. Irons, J. Zemen and A. M. Teale, Efficient Calculation of Molecular Integrals over London Atomic Orbitals, J. Chem. Theory Comput., 13 (2017), 2598.
[21] E. I. Tellgren, A. M. Teale, J. W. Furness, K. K. Lange, U. Ekstr¨om and T. Helgaker, Non-perturbative calculation of molecular magnetic properties within currentdensity functional theory, J. Chem. Phys., 140 (2014) 034101.
[22] J. W. Furness, J. Verbeke, E. I. Tellgren, S. Stopkowicz, U. Ekstr¨om, T. Helgaker and A. M. Teale, Current Density-Functional Theory using meta-Generalized Gradient Exchange-Correlation Functionals, J. Chem. Theory Comput., 11 (2015), 4169.
[23] J. W. Furness, U. Ekstr¨om, T. Hegaker and A. M. Teale, Electron localisation function in current-density-functional theory, Mol. Phys., 114 (2016), 1415.
[24] S. Reimann, U. Ekstr¨om, S. Stopkowicz, A. M. Teale, A. Borgoo, T. Helgaker, The importance of current contributions to shielding constants in density-functional theory, Phys. Chem. Chem. Phys., 17 (2015), 18834. Mathematical Methods in Quantum Chemistry649 Wannier functions for metals Antoine Levitt (joint work with Horia Cornean, Anil Damle, David Gontier, Lin Lin and Domenico Monaco) This talk summarized a recent effort towards the understanding of the localization properties of Wannier functions for metals [1, 2]. Wannier functions are a way to represent effectively the occupied subspace Span P where P = 1(H≤ εF) is the spectral projector of a crystal in the independent electron approximation H =−12∆ + Vper, where Vperis a d-dimensional periodic potential. By the Bloch theorem, this subspace is spanned by a set of Bloch functions, delocalized on the whole crystal. Wannier functions are a recombination of these Bloch functions to yield a set of localized functions that, together with their translates along the crystal lattice, form an orthogonal basis of the subspace. This is useful in applications to examine chemical bonding, interpolate band structures and speed up computations of exchange-correlation functionals [4]. By a Fourier-type duality, it is equivalent to find localized Wannier functions and to find a smooth and orthogonal basis of Span Pk, where Pk= 1(Hk≤ εF), Hkare the Bloch-Floquet fibers of H, and k runs on the Brillouin zoneB, a set having the topology of a d-dimensional torus. For an insulator, there is a gap at the Fermi level εF, and Pkis smooth as a function of k, and the existence of localized Wannier is well-known to be equivalent to the triviality of the fiber bundle Span Pk. For metals, Pkis not smooth as a function of k, and no localized Wannier functions in the classical sense can exist. However, a generalization of the concept is widely used in practice: instead of looking for an orthogonal basis of Span P , one looks instead for a larger set of orthogonal functions that span Span P . This is useful in reproducing band structure information, for instance in the context of Wannier interpolation [4]. In [1], we prove the existence of such a set. The main difficulty is the local topology (Chern number) carried by eigenvalue crossings. The construction is made possible by the fact that the sum on the Brillouin zone of the local topological numbers must vanish, by a Poincar´e-Hopf type argument. This implies the existence of almost-exponentially localized Wannier functions in the sense above for a large class of metals. In [2], we present a careful numerical study of Wannier functions for the simplest metal, the free electron gas. We show that while almost-exponentially localized functions do exist, the most widely used method to find them (the maximallylocalized Wannier functions scheme, see [4]) only yields poorly localized Wannier functions. 650Oberwolfach Report 13/2018 References
[25] H. Cornean, D. Gontier, A. Levitt and D. Monaco, Localised Wannier functions in metallic systems, submitted.
[26] A. Damle, A. Levitt and L. Lin, Variational formulation for Wannier functions with entangled band structure, submitted.
[27] G. Panati, Triviality of Bloch and Bloch-Dirac bundles, Annales Henri Poincar´e, (2007).
[28] N. Marzari, A. Mostofi, J. R. Yates, I. Souza and D. Vanderbilt, Maximally localized Wannier functions: Theory and applications, Rev. Mod. Phys., (2012). Numerical methods for Brillouin zone integration David Gontier (joint work with Eric Canc‘es, Virginie Ehrlacher, Antoine Levitt and Damiano Lombardi) In condensed matter, when studying numerically a periodic crystal, one needs to integrate some quantities over the Brillouin zone (a torus) T. This is the case forR the integrated density of states, of the formN (ε) :=T1(ε(k)≤ ε)dk, where ε(·) is a given smooth function over T, and for the energy per unit cell, of the formR E :=ε(k)1(ε(k)≤ εF)dk, where εFis the Fermi level, solution toN (εF) = N , T where N is a given number of electrons per unit cell. In this talk, we discussed the different numerical methods that have been proposed to discretise numerically these integrals. The main difficulty is that the error on the energy E depends not only on the discretisation, but also on the error made on the Fermi level εF, which itself depends on the error on the integrated density of statesN (·). In the case of insulators (existence of a spectral gap), we can rewrite the energy as the integral of a periodic smooth function. In this case, the approximation of the integral with a Riemann sum is enough to obtain an exponential rate of convergence with respect to the number of discretisation points (see e.g. [1]). In the case of metallic systems, the convergence is much slower, and two families of methods have been proposed. In interpolation methods, we interpolate the (smooth) function ε(·), and we compute the quantities of interest from the interpolation. In smearing methods, we replace the discontinuous step function 1 by a smooth one, so that the integrand becomes smooth. In this talk, we presented the results of [2], where we provided the rates of convergence of these methods. References
[29] D. Gontier and S. Lahbabi, Convergence rates of supercell calculations in the reduced Hartree-Fock model, M2AN (2016), pp. 1403–1424.
[30] ´E. Canc‘es, V. Ehrlacher, D. Gontier, A. Levitt and D. Lombardi, Numerical Quadrature in the Brillouin zone for periodic Schr¨odinger operators, in preparation. Mathematical Methods in Quantum Chemistry651 A unified approach to Wannier interpolation Anil Damle (joint work with Lin Lin) The construction of localized representations of electronic wavefunctions have a wide range of applications across quantum physics, material science, and chemistry. This is known as the Wannier localization problem, and the desired localized representations are known as Wannier functions [22, 13, 2]. For insulating materials with isolated eigenvalues, the problem is well studied mathematically [13, 2, 18, 12, 3, 10, 9, 19, 20, 17, 4, 7, 5] and good methods exist for computing Wannier functions in this setting [16, 6]. In contrast, when the eigenvalues are not isolated and become entangled far less is known. This situation arises naturally when considering metallic systems, but is also in insulating systems when considering a selected range of valance bands or conduction bands. In this setting the problem is significantly more difficult. Construction of Wannier functions now requires both identifying a subspace that contains a localized basis and constructing such a basis. A common methodology in this scenario is a disentanglement procedure [21]. Via two successive optimization procedures a subspace is identified and then a localized basis is computed. In this talk we propose a unified method to address this so-called Wannier localization problem for both the isolated and entangled cases. We avoid an explicit, initial disentanglement step and utilize a quasi-density matrix that “entangles” eigenfunctions of interest with the rest of them in a controlled manner. This allows us to simultaneously identify the desired subspace and construct the localized basis via a direct method, thereby avoiding any dependence on using a good initial guess. In the isolated setting, our approach reduces to the existing selected columns of the density matrix (SCDM) algorithm [6]. The key idea behind the SCDM methodology is the use of a column-pivoted QR factorization to select “good” columns of a density or quasi-density matrix from which the desired localized functions can be constructed directly. While the SCDM methodology extends cleanly from molecules to crystals, here we restrict our discussion to molecules for simplicity, though the crystal case will be addressed in the talk. Details of the generalization to crystals and corresponding numerical experiments may be found in the preprint [8] available online. We consider a single particle theory such as Kohn-Sham density functional theory (KSDFT) [11, 14], where the electronic wavefunctions are given by the eigenfunctions of a self-adjoint HamiltonianH. We denote the eigenfunctions as {ψi(r)} and are interested in those that satisfy (1)Hψi(r) = εiψi(r),εi∈ I for some energy windowI that indicates the eigenfunctions of interest. Mathematically, the Wannier localization problem corresponds to finding orthonormal and localized functions{wj} such that i∈I⊆ Vw:= span{wj}. 652Oberwolfach Report 13/2018 The functions{wj} are the desired Wannier functions. In this formulation, the aforementioned isolated case is characterized by an isolation condition (2)inf|εi− εi′| > 0. εi∈I,εi′∈I/ When (2) is violated, we consider the eigenvalues entangled. In the entangled setting, we construct the SCDM matrix through the use of a quasi-density matrix X (3)P =|ψiif(εi)hψi| = f(H), i for some smooth f (·) such that I is a subset of the support of f. The two most common scenarios for entangled eigenvalues corresponds to Wannier localization either below or around a certain energy level. In both cases we choose f to be large on the region of interest and smoothly decay to zero outsideI in a controlled manner. While formally the sum is over all the eigenfunctions, we can choose f to decay rapidly (e.g. as a complementary error function or Gaussian). By choosing f as an indicator function onI we can recover the density matrix for an isolated system. Given a smooth f , the kernel of the quasi-density matrix P (r, r′) decays rapidly [1, 15]. Therefore, to find NwWannier functions we select the Nw“most representative” and well conditioned columns of P and use them to construct Wannier functions. DefineE = diag [{εi}] ∈ RN ×Ncontaining all eigenvalues such that f (ε) is above some threshold and Ψ∈ CNg×Nbe the matrix containing the corresponding discretized eigenvectors. Computing the column-pivoted QR factorization (4)(f (E)Ψ∗) Π = QR we letC = {ri}Ni=1wdenote the Nwspatial points associated with the Nwcolumns permuted to the front by Π. Defining Ξ∈ CN ×Nwwith Ξi,i′= f (εi)ψ∗i(ri′), if the singular values of Ξ are bounded away from zero U = Ξ(Ξ∗Ξ)−12allows us to construct Wannier functions as XN (5)wi(r) =ψi′(r)Ui′,i i′=1 for i = i, . . . , Nw. In this setting U∈ CN ×Nwis a rectangular matrix with orthonormal columns refereed to as a gauge. The locality of the Wannier functions is a consequence of the locality of the quasi-density matrix (or the density matrix in the isolated case). Therefore, this SCDM methodology is able to construct Wannier functions for a wide range of systems directly, without relying on a good initial guess for the localized orbitals. One particularly important application of these localized functions is for band-structure interpolation of crystals. Examples of SCDM constructed Wannier functions in this setting may be found in [8]. Mathematical Methods in Quantum Chemistry653 References
[31] M. Benzi, P. Boito and N. Razouk, Decay properties of spectral projectors with applications to electronic structure, SIAM Rev., 55 (2013), pp. 3–64.
[32] E. I. Blount, Formalisms of band theory, Solid State Phys., 13 (1962), pp. 305–373.
[33] C. Brouder, G. Panati, M. Calandra, C. Mourougane and N. Marzari, Exponential localization of Wannier functions in insulators, Phys. Rev. Lett., 98 (2007), p. 046402.
[34] E. Canc‘es, A. Levitt, G. Panati and G. Stoltz, Robust determination of maximallylocalized Wannier functions, Phys. Rev. B, 95 (2017), p. 075114.
[35] A. Damle, L. Lin and L. Ying, Computing localized representations of the kohn-sham subspace via randomization and refinement, SIAM J. Sci. Comput., to appear.
[36] A. Damle, L. Lin and L. Ying, Compressed representation of Kohn–Sham orbitals via selected columns of the density matrix, J. Chem. Theory Comput., 11 (2015), pp. 1463– 1469.
[37] A. Damle, L. Lin and L. Ying, Scdm-k: Localized orbitals for solids via selected columns of the density matrix, J. Comput. Phys., 334 (2017), pp. 1–15.
[38] A. Damle and L. Lin, Disentanglement via entanglement: A unified method for wannier localization, arXiv:1703.06958, (2017).
[39] W. E, T. Li and J. Lu, Localized bases of eigensubspaces and operator compression, Proc. Natl. Acad. Sci., 107 (2010), pp. 1273–1278.
[40] F. Gygi, Compact representations of Kohn–Sham invariant subspaces, Phys. Rev. Lett., 102(2009), p. 166406.
[41] P. Hohenberg and W. Kohn, Inhomogeneous electron gas, Phys. Rev., 136 (1964), pp. B864–B871.
[42] E. Koch and S. Goedecker, Locality properties and Wannier functions for interacting systems, Solid State Commun., 119 (2001), p. 105.
[43] W. Kohn, Analytic properties of Bloch waves and Wannier functions, Phys. Rev., 115 (1959), p. 809.
[44] W. Kohn and L. Sham, Self-consistent equations including exchange and correlation effects, Phys. Rev., 140 (1965), pp. A1133–A1138.
[45] L. Lin, Localized spectrum slicing, Math. Comp., (in press).
[46] N. Marzari and D. Vanderbilt, Maximally localized generalized Wannier functions for composite energy bands, Phys. Rev. B, 56 (1997), p. 12847.
[47] J. I. Mustafa, S. Coh, M. L. Cohen and S. G. Louie, Automated construction of maximally localized Wannier functions: Optimized projection functions method, Phys. Rev. B, 92(2015), p. 165134.
[48] G. Nenciu, Dynamics of band electrons in electric and magnetic fields: rigorous justification of the effective hamiltonians, Rev. Mod. Phys., 63 (1991), pp. 91–127.
[49] V. Ozolin¸ˇs, R. Lai, R. Caflisch and S. Osher, Compressed modes for variational problems in mathematics and physics, Proc. Natl. Acad. Sci., 110 (2013), pp. 18368–18373.
[50] G. Panati and A. Pisante, Bloch bundles, Marzari-Vanderbilt functional and maximally localized Wannier functions, Commun. Math. Phys., 322 (2013), pp. 835–875.
[51] I. Souza, N. Marzari and D. Vanderbilt, Maximally localized wannier functions for entangled energy bands, Phys. Rev. B, 65 (2001), p. 035109.
[52] G. H. Wannier, The structure of electronic excitation levels in insulating crystals, Phys. Rev., 52 (1937), p. 191. 654Oberwolfach Report 13/2018 Magnetic fields, convexity, and gauge invariance in density-functional theory Erik I. Tellgren Density-functional theory is one of the most widely used electronic structure methods. Electrons interacting with a magnetic field are beyond the scope of the standard formulation of the theory, but two appropriate extensions are available [1, 2]. A novel extension can be obtained by using a magnetic Maxwell-Schr¨odinger model as the point of departure [3]. This extension in some respects generalizes Lieb’s seminal work on standard density-functional theory [4] and can be interpreted as a Moreau-Yosida regularization of a more conventional model. Denote the external electrostatic potential by v and the external magnetic vector potential by A. In atomic units, the standard Schr¨odinger Hamiltonian is given by XN (1)H(v, A) = T (A) +v(rl) + W l=1 where the kinetic energy and the electron-electron repulsion terms are given by 1XN2XN1 2− i∇l+ A(rl)and W =|r. l=1k<lk− rl| Following Lieb [4], densities ρ and electrostatic potentials v are taken to belongR to dual function spaces, with a pairing (ρ|v) =ρ(r) v(r) dr. Specifically, ρ∈ X := L1(R3)∩ L3(R3) and v∈ X∗:= L3/2(R3) + L∞(R3). A magnetic field B=∇ × A ∈ L2div(R3) is taken to be a divergence free L2(R3) function. The space of magnetic vector potentials is defined via the Fourier transform ˆAas those potentials that satisfy k2Aˆ(k)∈ L2. The Schr¨odinger ground state energy is obtained by minimization over all wave functions in a magnetic Sobolev space HA
[53] , (3)E(v, A) = infhΨ|H(v, A)|Ψi, Ψ∈HA Exploiting gauge invariance, E(v, A) = E(v, A +∇χ) for all suitable restricted χ, this energy can be reinterpreted as a function of B. Moreover, the energy can be expressed as a nested minimization over densities ρ and those wave functions Ψ7→ ρ that reproduce a given density,  (4)E(v, B) = inf(ρ|v) + FGH(ρ, B), ρ∈X with the Grayce-Harris functional given by (5)FGH(ρ, B) = infhΨ|H(0, A)|Ψi, Ψ7→ρ The energy functional E is known to be concave in v, but neither convex nor concave in B, and the formula above can be identified with a partial LegendreFenchel transformation of FGHwith respect to the first variable. Mathematical Methods in Quantum Chemistry655 Now let α be a vector potential that describes an internal magnetic field induced by the electrons. In the magnetic Maxwell-Schr¨odinger model, with the notation β=∇ × α,  kβk2+ inf α2µΨhΨ|H(v, A + α)|Ψi. The first term is the self-energy of the internal magnetic field, µ > 0 is a parameter with a physical interpretation as the magnetic permeability of vacuum, and the second term can be identified with the Schr¨odinger energy E(v, B + β). A change of variables b = B + β now yields  kB − bk2+ E(v, b), (7)Eµ(v, B) = inf b2µ which is a Moreau-Yosida regularization of the Schr¨odinger energy E(v, B) with respect to B [6, 7]. It can be remarked that this regularization is most natural to perform on a convex function, whereas E(v,·) is non-convex. Nonetheless, the regularization is beneficial in that it produces a paraconcave (“concave to within a square”) energy functional. More specifically, using Eq. (4) above, one finds (8)E¯µ(v, B) := Eµ(v, B)−kBk2= inf(ρ|v) −(b|B)+kbk2+ FGH(ρ, b), 2µρ,bµ2µ which is seen to be jointly concave in (v, B). Save for the unconventional factor −1/µ (which can be eliminated by a change of variables b′= b/µ), the above expression can be identified as a full Legendre-Fenchel transformation of the shifted Grayce-Harris functional F (ρ, b) = (2µ)−1kbk2+ FGH(ρ, b). I have not been able to establish that F (ρ, b) is jointly convex. However, a full Legendre-Fenchel transformation of the concave ¯Eµyields a functional (9)f¯µ(ρ, b) = supE¯µ(v, B)− (ρ|v) +(b|B) v,Bµ that is jointly convex by construction and also a lower bound ¯fµ(ρ, b)≤ F (ρ, b). Eq. (8) remains valid when F is replaced by ¯fµ. Hence, the shifted ground state energy ¯Eµ(v, B) and the universal density functional ¯fµ(ρ, b) constitute a pair of dual functionals that are each other’s Legendre-Fenchel transform (with suitable sign conventions—Eqs. (8) and (9) involve minimization and maximization, respectively). The interpretation in terms of regularization techniques can be taken further by introducing a second parameter λ≤ µ and defining the concave functional (10)E¯λµ(v, B) = Eµ(v, B)−kBk2= EkBk2 11 2λµ(v, B)−2µ+2µ−2λkBk2. The dual functional can be obtained as the full Legendre-Fenchel transformation f¯λµ(ρ, b) = supE¯λµ(v, B)− (ρ|v) +(b|B)= supe¯ (11)v,Bλvλµ(v, b)− (ρ|v). 656Oberwolfach Report 13/2018 Here, a type of partial Legendre-Fenchel transformation with respect to the magnetic field has been introduced, (12) ¯eλµ(v, b) = supE¯λµ(v, B)+(b|B)=kbk2+supEkB − bk2 Bλ2λBµ(v, B)−λ. The last term is a Lasry-Lions regularization of E, given by a double MoreauYosida regularization of E, first as a minimization problem with parameter µ and second as a maximization problem with parameter λ. In the case of non-convex functions on Hilbert spaces, and parameters λ < µ, this type of operation allows for stronger regularity results than simple Moreau-Yosida regularization [8, 9]. The formulation outlined above provides a gauge-invariant alternative to the Schr¨odinger-based convex formulation in Ref. [10], which featured the gauge dependent paramagnetic current density, Z (13)jp(r1) = ImΨ(r1, . . . , rN)∗∇Ψ(r1, . . . , rN)dr2· · · drN, as a variational parameter. The formulation also sheds light on attempts to formulate a density-functional theory that features the electron density ρ and the physical (gauge invariant) current density j as basic variables. (In the Schr¨odinger model, j = jp+ ρA, and in the Maxwell-Schr¨odinger model, j = jp+ ρ(A + α).) There are strong indications that such formulations are not compatible with the standard variational principle [11]. In the present setting, this problem is circumvented in an instructive way. The total magnetic field b = B + β can be equivalently represented as a current density k that is defined to satisfy µk =−∇ × b. The current density k is an independent variational parameter, with no general connection to the electronic wave function. However, the minimizing wave function Ψ and α in Eq. (6) are related—the stationarity condition reproduces a version of Ampere’s law, (14)µk = µj− ∇ × B, that relates a representation of the internal magnetic field (κ) to the physical current density. Hence, while j is not suitable as a variational parameter, k can play this role and the two quantities agree to within µ−1∇ × B at the minimum. References
[54] C. J. Grayce and R. A. Harris, Magnetic-field density-functional theory, Phys. Rev. A, 50(1994), pp. 3089–3095.
[55] G. Vignale and M. Rasolt, Density-functional theory in strong magnetic fields, Phys. Rev. Lett., 59 (1987), pp. 2360–2363.
[56] E. I. Tellgren, Density-functional theory for internal magnetic fields, Phys. Rev. A, 97 (2018), 012504.
[57] E. H. Lieb, Density functionals for Coulomb systems, Int. J. Quantum Chem., 24 (1983), pp. 243–277.
[58] E. H. Lieb and M. Loss, Analysis, American Mathematical Society (2001).
[59] R. T. Rockafellar and R. J-B Wets, Variational Analysis, Springer (2009).
[60] S. Kvaal, U. Ekstr¨om, A. M. Teale and T. Helgaker, Differentiable but exact formulation of density-functional theory, J. Chem. Phys., 140 (2014), 18A518. Mathematical Methods in Quantum Chemistry657
[61] J. M. Lasry and P. L. Lions, A remark on regularization in Hilbert space, Isr. J. Math. (1986), pp. 257–266.
[62] H. Attouch and D. Aze, Approximation and regularization of arbitrary functions in Hilbert spaces by the Lasry-Lions method, Ann. Inst. Henri Poincar´e (1993), pp. 289–312.
[63] E. I. Tellgren, S. Kvaal, E. Sagvolden, U. Ekstr¨om, A. M. Teale and T. Helgaker, Choice of basic variables in current-density-functional theory, Phys. Rev. A (2012), 062506.
[64] A. Laestadius and M. Benedicks, Nonexistence of a Hohenberg–Kohn variational principle in total current-density-functional theory, Phys. Rev. A (2015), 032508. The Future of the Electron-Correlation Problem: What about Full Configuration Interaction? J¨urgen Gauss (joint work with Janus J. Eriksen) The electron-correlation problem is at the heart of modern quantum chemistry, as its adequate treatment is essential for obtaining reliable and accurate results. While for so-called single-reference problems, i.e., cases where the wavefunction is well described by just a single Slater determinant, coupled-cluster theory has become the standard, a proper treatment of multireference cases still poses a challenge. To this day, a satisfactory formulation of a multireference extension of coupled-cluster theory still does not exist. The difficulties are here mainly due to the large set of requirements one tries to enforce, namely (1) size extensivity and consistency, (2) orbital invariance, (3) scaling of the cost independent of the size of the reference space, (4) proper number of parameters in the wavefunction ansatz, and (5) correct pole structure of the response function. Interest in multireference coupled-cluster theory is still high, but apparently there is currently a lack of new ideas. A possible conclusion might be instead to seek for radically different solutions to the electron-correlation problem, in particular, for solutions that avoid the complexity of coupled-cluster theory. In this context, we suggest to aim directly at the full configuration-interaction (FCI) ansatz, which provides the exact solution of the correlation problem within a given one-electron basis set. By doing this, one avoids all issues concerning the list of formal requirements given above. With respect to the enormous computational cost of FCI, we suggest to solve the FCI problem only in a numerically approximate manner, i.e., with an accuracy of about 0.01 to 0.1 mHartree. Such an accuracy is sufficient for most applications in chemistry, as µHartree accuracy is rarely needed. In addition, we recommend to use simple algorithms, as this might significantly facilitate (massive) parallelization. The latter is today required in order to exploit the existing computing power. There have actually been several suggestions for the solution of the electroncorrelation problem along the lines given above, e.g., the density-matrix renormalization group approach, stochastic FCI, as well as various variants of incremental or selected (full) CI. Our suggestion in this context is to use a many-body expansion (MBE) of the correlation energy [1]. This means that the FCI correlation energy EF CIis 658Oberwolfach Report 13/2018 obtained as nXoccnXoccnXocc EF CI=εi+∆εij+∆εijk+ . . . ii<ji<j<k where the first sum contains the correlation energies εiobtained with only one occupied spatial orbital correlated, the second sum accounts for corrections due to the simultaneous correlation of two occupied orbitals, the third for corrections due to the simultaneous correlation of three occupied orbitals, etc. The corrections are obtained in FCI calculations in which only the occupied orbitals referred to (indices i, j, k, . . .) are correlated together with the full set of virtual orbitals. This means that at the first level only up to double excitations appear in the FCI treatment, at the second level only up to quadruple excitations, etc. Our idea is not entirely new, as MBEs have often been used in quantum chemistry to treat large systems based on an expansion in terms of fragments. The MBE idea has also been used in the context of electron-correlation treatments; its first use goes back to work of Nesbet published in 1967. While, as seen in Figure 1a, the so-called MBE-FCI scheme indeed provides convergence to the FCI limit, an MBE with respect to occupied orbitals has only Figure 1.MBE-FCI results for water using a 6-31G basis set (H2O, 13 basis functions) and expansions with respect to (a) occupied and (b) virtual orbitals. >0(\(?@>&A``*&(``&()\)
[65] J. J. Eriksen, F. Lipparini and J. Gauss, Virtual orbital many-body expansions: a possible route towards the full configuration interaction limit, J. Phys. Chem. Lett., 8 (2017), pp. 4633–4639. Hartree-Fock Excited States Mathieu Lewin The Hartree-Fock model is an important nonlinear approximation of the N -particle Schr¨odinger linear equation which describes the N electrons in an atom or a molecule. Its successes and limitations for approximating the first eigenfunction of the N -body Hamiltonian (“ground state”) are well known, but its ability to describe excited states has probably been underestimated until recently. In a recent numerical work [1], a Hartree-Fock calculation has given very good predictions for 10 of the 11 first excited states of the H2molecule. These recent developments suggest that Hartree-Fock excited states need further mathematical and numerical investigation. In the talk I have explained a recent result from [4] where I constructed infinitely many excited states for neutral and positively charged atoms and molecules, with energy levels below the first energy threshold (the lowest energy of the same system with one electron removed), and above the true Schr¨odinger eigenvalues. To state the result, we introduce the N -particle operator describing the N electrons of a molecule XNXNXM H(N ) :=−∆xj−zm+X1 j=1j=1m=1|xj− Rm|1≤j<k≤N|xj− xk| with domain Ha2((R3× {↑, ↓})N)⊂ L2a((R3× {↑, ↓})N).The index a on Ha2 and L2adenotes the subspace of anti-symmetric functions, as is appropriate for fermions. A famous result of Zhislin and Sigalov [8, 9] asserts that when N < PM m=1zm+ 1, the operator H(N ) has infinitely many eigenvalues below its essential spectrum, denoted henceforth by λk(H(N )), k≥ 1. In the Hartree-Fock model, the set of wavefunctions is restricted to the special class of Slater determinants √−1 Ψ(x1, . . . , xN) =N !det(φj(xk)) where φ1, . . . , φNis any orthonormal system in L2(R3× {↑, ↓}). We call M the manifold of these special states and EHF(N ) = infhΨ, H(N)Ψi Ψ∈M Mathematical Methods in Quantum Chemistry661 the corresponding minimal energy. Minimizers were shown to exist in [5, 7] forP N <Mm=1zm+ 1 (neutral or positively charged molecules). Minimizers or other critical points of the quadratic form Ψ7→ hΨ, H(N)Ψi on M are found to solve the nonlinear coupled eigenvalue equations !N m=1|Rm− x|+k=1|φk|2∗| · |1φj−Xk=1(φkφj)∗| · |1φk= ǫjφj for j = 1, . . . , N . PM Theorem 1(Hartree-Fock excited states). Assume that N <m=1zm+1. Then there exists infinitely many critical points Ψ(k)of the energy Ψ7→ hΨ, H(N)Ψi on the Hartree-Fock manifoldM, such that for all k ≥ 1 (2)λk(H(H))≤ hΨ(k), H(N )Ψ(k)i < EHF(N− 1) < 0. In a nonlinear theory such as Hartree-Fock, there is no clear equivalent of the essential spectrum. But it has been proved in [2, 3] that the lowest Hartree-Fock energy EHF(N− 1) of N − 1 electrons satisfies similar mathematical properties as the bottom of an essential spectrum in the linear case, at least concerning the compactness of minimizing sequences. It is therefore natural to look for excited states with an energy below EHF(N− 1). In [6, 7], Lions was the first to construct infinitely many Hartree-Fock excited states, but those have an energy converging to 0, which does not correspond to the picture of the Zhislin-Sigalov theorem. He used a min-max method in the space H1(R3× {↑, ↓})Nof the N orbitals and Morse index bounds of min-maxing sequences to infer their compactness. His argument works the same for other similar nonlinear models, for instance when the exchange term (the last non-local term in (1)) is dropped. In our proof of [4] we construct different critical points by working in the N particle space. We therefore use that the Hartree-Fock model is obtained by restricting the true Schr¨odinger energy to the manifoldM and our argument does not work for other similar nonlinear models. We use min-max techniques which can also be applied to the linear problem, and obtain thereby the variational inequality on the left of (2), by a simple comparison principle. We think our critical points are better candidates for representing molecular excited states since their energies converge to EHF(N− 1) < 0 and not to 0. We did not use Morse index bounds to obtain the compactness of minimizing sequences, but actually proved that the Hartree-Fock energy satisfies the Palais-Smale condition below the first energy threshold EHF(N− 1), which allows us to use more classical techniques of nonlinear analysis. This is another property of Hartree-Fock theory which allows us to think of EHF(N− 1) as the bottom of a kind of nonlinear essential spectrum. This research has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement MDFT No 725528). 662Oberwolfach Report 13/2018 References
[66] G. M. J. Barca, A. T. B. Gilbert and P. M. W. Gill, Communication: Hartree-fock description of excited states of h2, The Journal of Chemical Physics, 141 (2014), p. 111104.
[67] G. Friesecke, The multiconfiguration equations for atoms and molecules: charge quantization and existence of solutions, Arch. Ration. Mech. Anal., 169 (2003), pp. 35–71.
[68] M. Lewin, Geometric methods for nonlinear many-body quantum systems, J. Funct. Anal., 260(2011), pp. 3535–3595.
[69] M. Lewin, Existence of Hartree-Fock excited states for atoms and molecules, Lett. Math. Phys., online first (2017).
[70] E. H. Lieb and B. Simon, The Hartree-Fock theory for Coulomb systems, Commun. Math. Phys., 53 (1977), pp. 185–194.
[71] P.-L. Lions, Sur l’existence d’´etats excit´es dans la th´eorie de Hartree-Fock, C. R. Acad. Sci. Paris S´er. I Math., 294 (1982), pp. 377–379.
[72] P.-L. Lions, Solutions of Hartree-Fock equations for Coulomb systems, Commun. Math. Phys., 109 (1987), pp. 33–97.
[73] G. M. Zhislin, Discussion of the spectrum of Schr¨odinger operators for systems of many particles. (in Russian), Trudy Moskovskogo matematiceskogo obscestva, 9 (1960), pp. 81– 120.
[74] G. M. Zhislin and A. G. Sigalov, The spectrum of the energy operator for atoms with fixed nuclei on subspaces corresponding to irreducible representations of the group of permutations, Izv. Akad. Nauk SSSR Ser. Mat., 29 (1965), pp. 835–860. On the approximation of wavefunctions by anisotropic Gauss functions Harry Yserentant (joint work with Stephan Scholz) The electronic Schr¨odinger equation describes the motion of N electrons under Coulomb interaction forces in a field of clamped nuclei. The solutions of this equation, the electronic wavefunctions, depend on 3N variables, three spatial dimensions for each electron. Approximating them is therefore inordinately challenging. I discussed the approximability of these wavefunctions by linear combinations of anisotropic Gauss functions, or more precisely Gauss-Hermite functions, products of polynomials and anisotropic Gauss functions in the narrow sense. Main result is that the original, singular wavefunctions can up to given accuracy and a negligibly small residual error be approximated with only insignificantly more such terms than their convolution with a Gaussian kernel of sufficiently small width and that basically arbitrary orders of convergence can be reached. This is a fairly surprising result, since it essentially means that by this type of approximation, the intricate hierarchies of non-smooth cusps in electronic wavefunctions have almost no impact on the convergence, once the global structure is resolved. References
[75] St. Scholz and H. Yserentant, On the approximation of electronic wavefunctions by anisotropic Gauss and Gauss-Hermite functions, Numer. Math. 136, pp. 841–874, 2017; also available as arXiv:1612.00360 [math.NA] Mathematical Methods in Quantum Chemistry663 X2C as an effective Hamiltonian Trond Saue A starting point for relativistic molecular quantum mechanical calculations is the Dirac equation   c (σ· ˆp)VeN− mec2ψS=ψψLSE, where VeNrepresents the interaction of electrons with the nuclear framework, c is the speed of light, methe mass of the electron, p the momentum operator and σ= (σx, σy, σz) is a vector collecting the Pauli spin matrices. The solutions are 4component complex vector functions, with the upper two components collected as the large components ψLand the lower two collected as the small components ψS. A troublesome feature of the Dirac equation, as well as effective one-electron equations such as the 4-component relativistic Hartree-Fock or Kohn-Sham equation, is the presence of solutions of negative energy. Contrary to classical mechanics such solutions can not be ignored, since there is always a finite probability of energy jumps, even across gaps of energy of the order of mec2. A major breakthrough in the field of relativistic molecular calculations was the development of the exact 2-component relativistic Hamiltonian (X2C),[4, 7, 10, 11, 13, 8, 5, 15, 16, 17, 14] which reproduces exactly the positive energy spectrum of the 4-component parent Hamiltonian. It is obtained by block diagonalization of the parent Hamiltonian  (2)U†hLLhLSU =h++0. hSLhSS0h−− In my contribution I point out that the X2C Hamiltonian is an example of an effective Hamiltonian, that is, that is an operator which reproduces exactly selected eigenvalues of the full Hamiltonian, but with approximate eigenvector restricted to a chosen model space. The machinery of effective Hamiltonians is nicely expressed in terms of principal angles and associated principal vectors relating the target space and the model space [6, 9]. I demonstrate that the X2C Hamiltonian is a canonical effective Hamiltonian[12] for which the effective eigenvectors are obtained by direct rotation of eigenvectors from the target space into the model space [18] under the constraint of least change [3]. An alternative 2-component relativistic Hamiltonian is obtained by the unnormalized elimination of the small component (UESC) [4]. This an example of a Bloch-Brandow effective Hamiltonian [1, 2]. In this case the effective eigenvectors are obtained by projection and the effective Hamiltonian is not self-adjoint. 664Oberwolfach Report 13/2018 References
[76] C. Bloch, Sur la th´eorie des perturbations des ´etats li´es, Nuc. Phys., 6, pp. 329–347, 1958.
[77] B. H. Brandow, Linked-Cluster Expansions for the Nuclear Many-Body Problem, Rev. Mod. Phys., 39, pp. 771–828, Oct 1967.
[78] L. S. Cederbaum, J. Schirmer and H. D. Meyer, Block diagonalisation of Hermitian matrices, J. Phys. A: Math. Gen., 22 (13), 2427, 1989.
[79] K. Dyall, Interfacing relativistic and nonrelativistic methods. I. Normalized eliminiation of the small component in the modified Dirac equation, J. Chem. Phys., 106, 9618, 1997.
[80] M. Filatov and K. Dyall, On convergence of the normalized elimination of the small component (NESC) method, Theor. Chem. Acc., 117, pp. 333–338, 2007.
[81] G. Golub and C. van Loan, Matrix Computations, Johns Hopkins, Baltimore, MD, 1989.
[82] M. Iliaˇs, H. J. Aa. Jensen, V. Kell¨o, B. O. Roos and M. Urban, Theoretical study of PbO and the PbO anion, Chem. Phys. Lett., 408, 210, 2005.
[83] M. Iliaˇs and T. Saue, An infinite-order two-component relativistic Hamiltonian by a simple one-step transformation, J. Chem. Phys., 126, 064102, 2007.
[84] A. V. Knyazev and M. E. Argentati, Principal Angles between Subspaces in an A-Based Scalar Product: Algorithms and Perturbation Estimates. SIAM Journal on Scientific Computing, 23 (6), pp. 2008–2040, 2002.
[85] W. Kutzelnigg and W. Liu, Quasirelativistic theory equivalent to fully relativistic theory, J. Chem. Phys., 123, 241102, 2005.
[86] W. Kutzelnigg and W. Liu, Quasirelativistic theory I. Theory in terms of a quasirelativistic operator, Mol. Phys., 104, pp. 2225–2240, 2006.
[87] S. Kvaal, Geometry of effective hamiltonians, Phys. Rev. C, 78, 044330, Oct 2008.
[88] W. Liu and D. Peng, Infinite-order quasirelativistic density functional method based on the exact matrix quasirelativistic theory, J. Chem. Phys., 125, 044102, 2006.
[89] W. Liu, Advances in relativistic molecular quantum mechanics, Phys. Rep., 537 (2), pp. 59– 89, 2014.
[90] W. Liu and W. Kutzelnigg, Quasirelativistic theory. II. Theory at matrix level, J. Chem. Phys., 126, 114107, 2007.
[91] W. Liu and D. Peng, Exact two-component hamiltonians revisited, J. Chem. Phys., 131, 031104, 2009.
[92] T. Saue, Relativistic Hamiltonians for Chemistry: A Primer, Chem. Phys. Chem., 12, 3077, 2011.
[93] J. H. Van Vleck, On σ-type doubling and electron spin in the spectra of diatomic molecules, Phys. Rev., 33, pp. 467–506, Apr 1929. An overview ofa posteriori estimation and post-processing methods for nonlinear eigenvalue problems Genevi‘eve Dusson (joint work with Eric Canc‘es, Yvon Maday, Benjamin Stamm and Martin Vohral´ık) Many mathematical models aiming at the determination of electronic structures give rise to nonlinear eigenvalue problems. This is for example the case of the Hartree-Fock model and the Density Functional Theory (DFT) models (e.g. KohnSham). Other nonlinear eigenvalue problems appear in different contexts, e.g. in nonlinear elasticity for the calculation of vibration modes of a mechanical structure, or the Gross-Pitaevskii equation which is used for the computation of the steady states of Bose-Einstein condensates. Mathematical Methods in Quantum Chemistry665 Controlling the errors arising in the resolution of such problems, and especially the errors due to (i)—the chosen discretization, and (ii)—the algorithm used to solve the eigenvalue problem (mostly iterative methods) is essential for at least two reasons. First, it allows to guarantee the accuracy of the approximation, and therefore to complement the computed output with an error bar. Second, it enables an optimization of the parameters used in the simulation (e.g. the number of degrees of freedom versus the number of iterations) to reach a given accuracy at a minimal computational cost. The control of this error is usually done by means of a posteriori error estimates. The main object of this work is to compare different error estimation methods for nonlinear eigenvalue problems in a unified framework. To do so, we consider the ground state of a generic relatively simple nonlinear eigenvalue problem arising from a minimization problem: Find (u, λ)∈ X × R, such that kuk = 1 and (1)Au + g(u2)u = λu, where A is a linear, bounded-below, self-adjoint operator on a separable Hilbert X with domain D(A) and compact resolvent. In our applications, g depends on the electronic density ρ, which is equal to u2in the case of only one eigenfunction. We focus on the lowest eigenvalue λ of the problem. A Galerkin discretization of this problem writes: Find (uNc, λNc)∈ XNc× R withkuNck = 1 such that (2)∀vNc∈ XNc,h A + g(u2Nc)uNc, vNci = λNchuNc, vNci. We compare and analyze here mainly four different works. These works provide error estimations for different linear and nonlinear eigenvalue problems. Note that these error estimations can be used in two different ways. They can be used to give a guaranteed bound on the error as explained above, but they can also be used to post-process the approximate solution, and therefore obtain a more accurate solution, hopefully at a low computational overcost. We first consider the error estimation and post-processing performed for the Hartree-Fock equations in 2003 [5]. In this article, the authors present an error estimation as well as a post-processing step for the discrete solution of the HartreeFock problem, which relies on a second-order Taylor expansion of the energy. Second, we present the a posteriori error estimation for the Gross-Pitaevskii equation presented in 2016 [4]. This estimation provides first a coarse guaranteed upper bound of the eigenvector error, based on a first-order Taylor expansion of the residual. The approach is then refined to propose a bound that is very close to the real error, valid only in the asymptotic regime where the error is small enough, but which can be checked in practice using the first guaranteed bound. Third, we introduce the post-processing method based on perturbation theory for the Kohn-Sham equations [2, 3]. This method provides an improvement for the energy as well as the orbitals in the case of a planewave discretization of these equations that is very cheap to compute. Finally, we mention the two-grid methods presented in 1999 for a linear eigenvalue problem [6], and in 2017 for a nonlinear eigenvalue problem [1]. In this 666Oberwolfach Report 13/2018 method, the full problem is solved on a coarse grid, and a boundary value problem (or a linear eigenvalue problem in the nonlinear case) is solved on a fine grid, to improve the solution on the coarse grid. In each case, it is possible to relate the error estimation or the post-processing to a first-order Taylor expansion of the error. More precisely, let us denote the residual of the equation by F , i.e. ∀v ∈ X, ∀µ ∈ R, F (v, µ) = Av + g(v2)v− µv. For the exact eigenpair (u, λ)∈ X × R, there holds F (u, λ) = Au + g(u2)u− λu = 0. Writing a first-order development of the residual around the coarse eigenpair (uNc, λNc)∈ XNc× R, where XNcis the discretization space gives an approximation of the error at first order by DF(u−1F (u Nc,λNc)Nc, λNc). These different methods providing error estimations or post-processing for the energy and the eigenfunction can all be derived from some approximation of this equation. This leads to different computational overcosts, and different precisions of the error estimates or the post-processed quantities. Presenting these methods in a same framework can help the practitioner to better choose the method corresponding to his/her needs and constraints in terms of computational cost and required accuracy. References
[94] E. Canc‘es, R. Chakir, L. He and Y. Maday, Two-grid methods for a class of nonlinear elliptic eigenvalue problems, IMA Journal of Numerical Analysis, 38 (2), pp. 605–645, 2017
[95] E. Canc‘es, G. Dusson, Y. Maday, B. Stamm M. Vohral´ık, A perturbation-method-based a posteriori estimator for the planewave discretization of nonlinear Schr¨odinger equations, Comptes Rendus Mathematique, 352 (11), pp. 941–946, 2014.
[96] E. Canc‘es, G. Dusson, Y. Maday, B. Stamm M. Vohral´ık, A perturbation-methodbased post-processing for the planewave discretization of Kohn–Sham models, Journal of Computational Physics, 307, pp. 446–459, 2016.
[97] G. Dusson and Y. Maday, A posteriori analysis of a nonlinear Gross–Pitaevskii-type eigenvalue problem, IMA Journal of Numerical Analysis, 37 (1), pp. 94–137, 2017.
[98] Y. Maday and G. Turinici, Error bars and quadratically convergent methods for the numerical simulation of the Hartree-Fock equations, Numerische Mathematik, 94 (4), pp. 739– 770, 2003.
[99] J. Xu and A. Zhou, A two-grid discretization scheme for eigenvalue problems, Mathematics of Computation, 70 (233), pp. 17–26, 1999. Mathematical Methods in Quantum Chemistry667 Density of states for optical spectra of molecules by low-rank approximation Boris Khoromskij (joint work with Peter Benner, Venera Khoromskaia, Sergey Dolgov and Chao Yang) The breaking through approach to low-parametric representation of multivariate functions and operators is based on the principle of separation of variables which can be realized by using nonlinear approximation in rank-structured tensor formats like canonical, Tucker and tensor train (TT) representations. The method of quantized TT (QTT) approximation is proven to provide the logarithmic datacompression on a wide class of discretized functions and operators [8]. The new range-separated tensor format [2] allows the efficient rank-structured representation of non-regular many-particle interaction potentials and scattered data. Thus bridging of the modern multilinear algebra and the nonlinear approximation theory allows to construct the novel tensor numerical methods with the linear complexity scaling in dimension, hence breaking the “curse of dimensionality” [9]. In this talk we provide the short introduction to tensor numerical methods and and briefly overview the most commonly used tensor formats. We then discuss how to construct the rank-structured tensor approximation of the basic two-electron integrals (TEI) tensor in the Hartree-Fock calculation [6] and in MP2 techniques [7], see the survey paper [12] on the related literature. The fourth order TEI tensor, B = [bµνκλ], is defined in the form ZZg (1)bµνκλ=µ(x)gν(x)gκ(y)gλ(y)dxdy, R3R3kx − yk where the size Nbof the basis set{gµ}, µ = 1, . . . , Nb(say, Gaussian type orbitals), chosen for approximating the Hartree-Fock equation, may vary from several tens to several hundreds and even larger. The large Nb2× Nb2TEI matrix B is approximated in the form of rank-RBdecomposition [6, 7] B := [bµν,κλ]≈ LLT,withL∈ RNb2×RB, RB∼ Nb, by using the truncated Cholesky decomposition of B with on-the-fly computation PR of the required column vectors B(:, j∗) =⊙3ℓ=1V(ℓ)Mk(ℓ)V(ℓ)(:, j∗)T. The latter k=1 is performed by using the precomputed factorization of B via the 1D “density fitting” scheme. We describe the reduced basis approach to the Bethe-Salpeter equation (BSE) (spectral problem describing excitation energies of molecules) introduced in [1]. The 2× 2-block representation of the BSE system matrix leads to the following eigenvalue problem,   yk−B∗−A∗yk= ωkxykk, 668Oberwolfach Report 13/2018 where the matrix blocks of size n× n with n = O(Nb2) are defined by (3)A = ∆ε + V− cW ,B = V− fW , and eigenvalues ωkcorrespond to the excitation energies. Here ∆ε is a diagonal matrix and V = [via,jb] is the rank-RBTEI matrix B≈ LLTprojected onto the Hartree-Fock molecular orbital basis [7, 1]. Matrices cW and fW can be constructed by using certain transforms of V . Our approach is based on the diagonal plus low-rank approximation of the matrix blocks A and B in BSE system matrix (by using some results in [11]) inherited from the low-rank Cholesky factorization of the TEI matrix in the molecular orbitals basis [7]. The structured eigenvalue solver, based on the use of the Sherman-Morrison-Wudbury formula for rank-structured representation of the matrix inverse, allows to accurately calculate several eigenvalues in the central part of the spectrum with the linear cost in the matrix size. The enhanced algorithm applies to the improved block-diagonal plus low-rank structure in the BSE matrix [1, 3]. The QTT tensor approximation of the matrix blocks allows to solve the partial structured eigenvalue problem in the logarithmic cost with respect to the large matrix size [3]. Finally, we present a new method to efficiently and accurately approximate the density of states (DOS) for large rank-structured matrices [4] with application to the diagonal block A in the BSE matrix (the so-called Tamm-Duncoff approximation). We use the simple definition of the DOS for symmetric matrices n 1X (4)φ(t) =δ(t− λj),t, λj∈ [0, a], n j=1 where δ is the Dirac delta and the λj’s are the eigenvalues of A = ATordered as λ1≤ λ2≤ · · · ≤ λn. The traditional methods for approximation of the DOS function [10], based on a polynomial or fractional-polynomial interpolation of the DOS regularized by Gaussians or Lorentzians, and computing traces of matrix resolvents, defined at a large set of interpolation points, are hardly applicable to the complicated spectral data of the large BSE system. Our approach is based on Lorentzian blurring, i.e. replacing each Dirac-δ by a Lorentzian, 1 1η1 πt2+ η2=πImt− iη, so that an approximate DOS can be written as n 1X nLη(t− λj). j=1 The regularized DOS takes the form (7) 1X11 φ(t)≈ φη(t) :=Im=Im Trace[(tI− A − iηI)−1], nπ(t− λj)− iηnπ j=1 Mathematical Methods in Quantum Chemistry669 or keeping real-valued arithmetics, 1Xnη1 nπ(t− λj)2+ η2=nπTrace[((tI− A)2+ η2I)−1]. j=1 In both cases the task of computing the approximate DOS φη(t) reduces to approximating the trace of the matrix resolvent (tI− A − iηI)−1or ((tI− A)2+ η2I)−1. Our method amounts to estimating the DOS by evaluating traces of the matrix resolvent represented in the block-diagonal plus low-rank form (the ShermanMorrison-Wudbury formula). The main contribution is to perform each function evaluation at low cost and to reduce the total number of function evaluations to O(log N ) in the case of fine representation grid of size N≫ n by using the QTT tensor interpolation via the adaptive cross approximation. The efficiency of presented tensor methods is illustrated by numerical examples. References
[100] P. Benner, V. Khoromskaia and B. N. Khoromskij, A reduced basis approach for calculation of the Bethe-Salpeter excitation energies using low-rank tensor factorizations, Molecular Physics, 114 (7–8), pp. 1148–1161, 2016.
[101] P. Benner, V. Khoromskaia and B. N. Khoromskij, Range-separated tensor formats for numerical modeling of many-particle interaction potentials, SIAM J. Sci. Comput., 40 (2), pp. A1034–A1062, 2018.
[102] P. Benner, S. Dolgov, V. Khoromskaia and B. N. Khoromskij, Fast iterative solution of the Bethe-Salpeter eigenvalue problem using low-rank and QTT tensor approximation, J. Comp. Phys., 334, pp. 221–239, 2017.
[103] P. Benner, V. Khoromskaia, B. N. Khoromskij and C. Yang, Computing the density of states for optical spectra by low-rank and QTT tensor approximation, available as epreprint: arXiv:1801.03852, 2018.
[104] V. Khoromskaia and B. N. Khoromskij, Tensor Numerical Methods in Quantum Chemistry: from Hartree-Fock to Excitation Energies, Phys. Chemistry Chem. Physics, 2015, 17, pp. 31491–31509.
[105] V. Khoromskaia, B. N. Khoromskij and R. Schneider, Tensor-structured calculation of two-electron integrals in a general basis, SIAM J. Sci. Comput., 35 (2), 2013, pp. A987– A1010.
[106] V. Khoromskaia and B. N. Khoromskij, Møller-Plesset (MP2) Energy Correction Using Tensor Factorizations of the Grid-based Two-electron Integrals, Comp. Phys. Comm., 185, pp. 2–10, 2014.
[107] B. N. Khoromskij, O(d log N )-Quantics Approximation of N -d Tensors in High-Dimensional Numerical Modeling, J. Constr. Approx. v. 34 (2), pp. 257–289 (2011).
[108] B. N. Khoromskij, Tensor Numerical Methods in Scientific Computing, Research monograph, De Gruyter Verlag, Berlin, 2018.
[109] L. Lin, Y. Saad and C. Yang, Approximating spectral densities of large matrices, SIAM Review, 58 (1), pp. 34–65, 2016. 670Oberwolfach Report 13/2018
[110] E. Rebolini, J. Toulouse and A. Savin, Electronic excitation energies of molecular systems from the Bethe-Salpeter equation: Example of H2molecule, in: Concepts and Methods in Modern Theoretical Chemistry (S. Ghosh and P. Chattaraj eds), vol 1: Electronic Structure and Reactivity, p. 367, 2013.
[111] S. Reine, T. Helgaker and R. Lindh, Multi-electron integrals, WIREs Comput. Mol. Sci. 2, pp. 290–303 (2012). A Quadratic Error Estimate for the Coupled-Cluster Method Tailored by Tensor Network States Fabian M. Faulstich (joint work with Andre Laestadius, Simen Kvaal, ¨Ors Legeza and Reinhold Schneider) The tailored Coupled-Cluster (TCC) method [4] is based on a very appealing ansatz, namely, the so-called split amplitude ansatz [6, 7, 1, 8]. Here, a CASFCI computation is performed on a subsystem yielding the solution ψCAS. This solution is then dynamically corrected, e.g., the CCSD scheme is such a correction capturing the dynamical correlation on the external system. Using this approach we are able to avoid the conventional MRCC methods and all their difficulties, but are still able to compute multi-configuration systems. A sever drawback of the TCC method is that it is based on the single-reference CC theory. Consequently, to cover the system’s static correlations we need a large CAS space, which scales exponentially in its dimensionality, i.e., it suffers from the curse of dimensionality. A correction to this was only recently achieved by approximating the FCI solution via the density matrix renormalization group (DMRG) method [13, 14], a high accuracy tool to compute static correlation, which fails, however, for dynamical correlation [14, 2]. Before presenting the TCC-energy error analysis performed in [3], we introduce the TCC method in a formal way: We start by defining the basis splitting for the considered N -electron problem. Let{χ1, . . . , χK} ⊆ H1(R3× {±12}) be a set of L2(R3× {±12})-orthonormal spin orbitals with K > N and φ0be the considered reference Slater determinant. We define BCAS={χ1, . . . , χN, χN +1, . . . , χk} , Bext={χk+1, . . . , χK} |{z}|{z}|{z} occupiedunoccupiedexternal and furthermoreBCAS={φ[µ1, . . . , µN] : µi∈ {1, . . . , k}, µ1<· · · < µN}. The corresponding FCI spaceHCASis then defined as the span ofBCAS. Analogously, we split the set of excitation-indicesJ describing the set of possible excitations, i.e.,JCAS={µ ∈ J : Xµφ0∈ HCAS} and Jext={µ ∈ J : Xµφ0∈ H/CAS}. An important observation is thatJextdoes not only contain excitations into states purely excited in Bextbut also into mixed states, i.e., for the virtual orbitals v(µ) = (A1, . . . , Ar) of the excitation index µ, there exists at least one l∈ {1, . . . , r} such that Al∈ {k + 1, . . . , K}. Note that the dynamic correction will note capture static correlation, hence, we needHCASto contain all Mathematical Methods in Quantum Chemistry671 of the system’s static correlation.Subsequently, we describe how this is ensured by using concepts from the field of quantum information theory. Let U be a low-rank DMRG solution on the FCI spaceHKrepresented in the TTformat [5]. From this solution we compute the one particle and two particle TT densities D[i] = U[i]U[i] and D[i, j] = U[i, j]U[i, j], where the matrisization 1,...,6mi,...,mK),(mi)∈ R2K−2×4is dei,...,6mj,...,mK)(mi,mj). From these densities we compute the mutual information I(i, j) = s(i) + s(j)− s(i, j), where s denotes the to D[i] resp. D[i, j] corresponding Von Neumann entropy, i.e., s(i) =−Tr(D[i] ln D[i]) resp. s(i, j) =−Tr(D[i, j] ln D[i, j]). For a more detailed description of the connection between these quantum information theory concepts and electronic structure theory, we refer the reader to [11, 12]. To describe the TCC method by means of a non-linear Galerkin scheme we consider K, N∈ N with K > N to be fixed, B={χ1, . . . , χK} to be a set of L2(R3× {±1})-orthonormal spin orbitals and φ 20 a reference state for an N -electron problem. Further, we assume BCASand Bextto be a partition of B and φCAS= eTCASφ0to be the FCI approximation on the complete active spaceHCAS. We define the TCC function f (· ; tCAS) :Vext→ Vext′ with t7→ f(t; tCAS), where for µ∈ Jext (f (t; tCAS))µ=hφµ, e−TCASe−THeTeTCASφ0i =hφµ, e−TCASexp(−XtνXν)H exp(XtνXν)eTCASφ0i . ν∈Jextν∈Jext In addition, let the TCC-energy functional be given by E(t; tCAS) =hφ0, e−TCASe−THeTeTCASφ0i . We now give a short descriptive summary of the energy-error analysis given in [3]. As a non-variational method, it is not immediate that the TCC method has a quadratic energy error bound. Besides deriving such a bound we subsequently present an error analysis solving the following three problems: First, note that the TCC method does not specify what method shall be used onHCASto compute ψCAS. For the sake of generality, we therefore want to derive an error bound that explicitly contains the error onHCAS, i.e., ∆ECAS=|hφ0, e−TCASP HP eTCAS− e−TFCICASP HP eTFCICASφ0i| . Second, we recall that for the single-reference CC method a quadratic error bound was derived by Schneider and Rohwedder [10, 9]. A major difference to the analysis presented in these articles is that we can not assume a HOMO-LUMO gap, as we are considering multi-configuration systems. Consequently, the in [10, 9] presented analysis needs to be generalized to the here considered case of a CAS-ext gap, i.e., ε0= λk+1− λk> 0 where{λk} denotes the spectrum of the Fock operator. Third, the TCC method will in general not converge to the FCI limit. This is due to its basis-splitting nature and shall be reflected in the energy-error bound. For the energy-error analysis let B ={χ1, . . . , χK} be a set of spin orbitals that are split into BCASand Bext. We denoteHKandHCASthe FCI space 672Oberwolfach Report 13/2018 corresponding to B resp. BCAS. Further, tCAS∗∈ VCAScorresponds to the FCI solution onHKprojected onHCAS, tCASFCI∈ VCASthe FCI solution onHCAS, and tCAS∈ VCASan approximation to tCASFCI. LetVext(d)⊂ Vextbe a subspace fulfilling d(t∗,Vext(d))≤γδ, γ + L where γ, L > 0 are the monotonicity and Lipschitz constants of f (· ; tCAS) on Bδ(t∗). This corresponds toVext(d)being a sufficiently good approximation space. Then there exists a unique solution td∈ Vext(d)of f (· ; tCAS)|V(d)= 0 that approxiext mates the solution t∗∈ Vextof f (· ; tCAS) = 0 onVext. Let (zd, z∗)∈ Vext(d)× Vext be the corresponding dual solutions of (td, t∗)∈ Vext(d)× Vext. Further, set ˜t∗∈ Vext the solution of f (· ; tCASFCI) = 0 onVextand text∗∈ Vextthe projection of the FCI solution onHKontoHCAS⊥. We now split the energy error ∆E =|E(td; tCAS)− E(text∗; tCAS∗)| |{z}|∗; tCAS)− E(t{z∗; tCASFCI)}| =∆ε=∆εCAS +|E(t∗; tCASFCI)− E(text∗; tCAS∗)| |{z} =∆ε∗CAS and derive upper bounds for each of the above ‘sub-errors’, i.e., ∆ε, ∆εCASand ∆ε∗CAS. By using a Taylor expansion and the fact that (text∗, tCAS∗) is a stationary point of the energy functionalE, we find (cf. Lemma 26 in [3]) ∆ε∗CAS=|E(t∗; tCASFCI)− E(text∗; tCAS∗)| . kt∗− text∗k2Vext+ktCASFCI− tCAS∗k22  ext∗− text∗k2Vext+ktCASFCI− tCAS∗k22. Furthermore, we can recast ∆ECAS, the energy error onHCASfrom the error expression ∆εCAS. We find (cf. Lemma 27 in [3]) ∆εCAS=|E(t∗; tCAS)− E(t∗; tCASFCI)| X .∆ECAS+kt∗− ˜t∗k2Vext+k(TCAS− TFCICAS)φ0k2H1+εµ(˜t∗)2µ. |µ|=1 Finally, it is straight forward to generalize the analysis presented in [10, 9] to hold for the above introduced CAS-ext gap assumption (cf. [3]). This yields ∆ε =|E(td; tCAS)− E(t∗; tCAS)| .ktd− t∗kVext(ktd− t∗kVext+kzd− z∗kVext) . Mathematical Methods in Quantum Chemistry673 Combining these estimates we obtain the following energy-error bound for the TCC method ∆E .ktd− t∗kVext(ktd− t∗kVext+kzd− z∗kVext) +kt∗− text∗k2Vext+kt∗− ˜t∗k2Vext X +ktCASFCI− tCAS∗k22+ktCAS− tCASFCIk22+εµ(˜t∗)2µ+ ∆ECAS. µ∈Jext |µ|=1 We see that if ∆ECASposses a quadratic error bound and the{˜t∗} amplitudes for single excitations are small, then we have a quadratic energy-error bound. For TNS-TCC methods we have such a quadratic bound for ∆ECASas the used tensor optimization algorithms are variational. Furthermore, we can assume without loss generality that the{˜t∗} amplitudes for single excitations are small, as a rotation into Brueckner orbitals ensures this property. Hence, TNS-TCC methods have the same asymptotic behavior as CI schemes. References
[112] L. Adamowicz, P. Piecuch and K. B. Ghose, The state-selective coupled cluster method for quasi-degenerate electronic states, Molecular Physics, (1998), pp. 225–234.
[113] G. K.-L. Chan and S. Sharma, The Density Matrix Renormalization Group in Quantum Chemistry, Annual Review of Physical Chemistry, (2011), pp. 465–481.
[114] F. M. Faulstich, A. Laestadius, S. Kvaal, ¨O. Legeza and R. Schneider, Analysis of The Coupled-Cluster Method Tailored by Tensor-Network States in Quantum Chemistry, arXiv:1802.05699 [math.NA].
[115] H. O. Kinoshita and R. J. Bartlett, Coupled-cluster method tailored by configuration interaction, The Journal of chemical physics 123 (7), (2005), 074106.
[116] I. V. Oseledets,Tensor-train decomposition, SIAM Journal on Scientific Computing (2011), pp. 2295–2317.
[117] P. Piecuch, N. Oliphant and L. Adamowicz, A state selective multireference coupled cluster theory employing the single reference formalism, The Journal of chemical physics, (1993), pp. 1875–1900.
[118] P. Piecuch and L. Adamowicz, State-Selective Multireference Coupled-Cluster Theory Employing the Single-Reference Formalism: Implementation and Application to the H8 Model System, The Journal of chemical physics, (1994), pp. 5792–5809.
[119] P. Piecuch, Active-space coupled-cluster methods, Molecular Physics, (2010), pp. 2987.3015.
[120] T. Rohwedder and R. Schneider, Error estimates for the Coupled Cluster method, ESAIM: Mathematical Modelling and Numerical Analysis, (2013).
[121] R. Schneider, Analysis of the projected coupled cluster method in electronic structure calculation, Numerische Mathematik, 113 (3), (2009), pp. 433–471.
[122] S. Szalay, G. Barcza, T. Szilv´asi, L. Veis and ¨O. Legeza, The correlation theory of the chemical bond, Scientific Reports, 7, (2016).
[123] S. Szalay et al., Tensor product methods and entanglement optimization for ab initio quantum chemistry, International Journal of Quantum Chemistry, (2015), pp. 1342–1391.
[124] L. Veis, A. Antal´ık, J. Brabec, F. Neese, ¨O. Legeza and J. Pittner, Coupled cluster method with single and double excitations tailored by matrix product state wave functions, The journal of physical chemistry letters, 7 (20), (2016), pp. 4072–4078.
[125] L. Veis, A. Antal´ık, ¨O. Legeza, A. Alavi and J. Pittner, Full Configuration Interaction Quantum Monte Carlo Benchmark and Multireference Coupled Cluster Studies of Tetramethyleneethane Diradical, arXiv:1801.01057 [physics.chem-ph]. 674Oberwolfach Report 13/2018 Singular Analysis of Coupled Cluster Equations and its Implication for Simulation Heinz-J¨urgen Flad (joint work with Gohar Flad-Harutyunyan) Modern first principles quantum chemistry relies to a large extend on models which make no explicit use of many-particle wavefunctions but instead are based on reduced quantities like density matrices, Green’s and response functions. For a system of particles interacting via singular Coulomb potentials these reduced quantities show a specific singular behaviour near diagonals. Whereas much effort has been spend to understand the singular behaviour of many-particle wavefunctions, cf. [8] much less is known about reduced quantities. We have started a research project to study the asymptotic properties of reduced quantities within the framework of coupled cluster (CC) theory, cf. [1]. This particular many-particle theory is presently considered as the ultimate benchmark in quantum chemistry. Despite its practical significance, a rigorous mathematical analysis of CC models is still in its infancy, see however recent progress on the subject by Schneider, Rohwedder, Kvaal and Laestadius [10, 12, 13, 14]. In order to achieve a detailed understanding of the reduced quantities mentioned before, it is first of all necessary to study the asymptotic behaviour of reduced pair and possibly higher order amplitudes which constitute the CC wavefunction via the exponential ansatz. This is the subject of the present work with focus on nonlinear CC models within the random phase approximation (RPA). Solutions of these models are commonly represented by series of a particular class of Goldstone diagrams so-called RPA diagrams. In Ref. [7], we presented a detailed asymptotic analysis of these RPA diagrams using techniques from singular analysis. In particular, we provided a connection between RPA diagrams and classical pseudo-differential operators which enables an efficient treatment of the linear and nonlinear interactions in these models. In order to give a brief outline of the problem and how we handled it, let us discuss a highly simplified nonlinear model problem for a pair amplitude τ (x1, x2), where x1, x2∈ R3denote the coordinates of an electron pair. This model contains various linear interaction terms of the pair amplitude given by fV(2)τ(x1, x2) := V(2)(x1, x2) τ (x1, x2),fV(1)τ(x1, x2) := V(1)(x1) τ (x1, x2), Z fV(2)◦τ(x1, x2) :=V(2)(x1, x3) τ (x3, x2) dx3 and a nonlinear term of the form Z fτ ◦V(2)◦τ(x1, x2) :=τ (x1, x3)V(2)(x3, x4)τ (x4, x2) dx3dx4 where V(1), V(2)correspond to one and two-particle interactions, respectively. The basic structure of the model can be represented by the nonlinear equation (1)Aτ =−V(2)− fV(2)τ+ fV(1)τ− fV(2)◦τ− fτ ◦V(2)◦τ. where A is an elliptic second order partial differential operator. Mathematical Methods in Quantum Chemistry675 It is common practice in quantum chemistry to solve (1) via a fixed point iteration scheme. First one solves the equation Aτ0=−V(2)with fixed right hand side. Calculate fV(2)τ0, fV(1)τ0, fV(2)◦τ0and fτ0◦V(2)◦τ0and in the next iteration step one solves (2)Aτ1=−V(2)− fV(2)τ0+ fV(1)τ0− fV(2)◦τ0− fτ0◦V(2)◦τ0. The last two steps can be repeated generating a iterative sequence of linear equations (3)Aτn+1=−V(2)− fV(2)τn+ fV(1)τn− fV(2)◦τn− fτn◦V(2)◦τn. which are solved in a consecutive manner until convergence of the sequence τ1, τ2, τ3, . . . has been achieved. In [7], we were mainly interested in the asymptotic behaviour of iterated pair amplitudes τi, i = 0, 1, . . ., near coalescence points of electrons. In order to extract these properties we applied methods from singular analysis [2, 9, 15] and solved (3) via an explicitly constructed asymptotic parametrix and corresponding Green operator, cf. [3, 4, 5, 6] for further details. In order to deal with the (non) linear interactions on the right hand side, we employed operator algebraic techniques in the framework of classical pseudo-differential calculus. In a nutshell, a parametrix P of a differential operator A is a pseudo-differential operator which can be considered as a generalised inverse, i.e., when applied from the left or right side it yields P A = I + Gl, AP = I + Gr, where the remainders Gl and Grare left and right Green operators, respectively. In contrast to the standard calculus of pseudo-differential operators on smooth manifolds, our calculus applies to singular spaces with conical, edge and corner type singularities as well, cf. the monographs [2, 9, 15] for a detailed exposition. In the smooth case, remainders correspond to compact operators with smooth kernel function, whereas in the singular calculus Green operators encode important asymptotic information which we want to extract. Acting on (3) with the parametrix from the left yields (4) τn+1=−P V(2)− P fV(2)τn+ P fV(1)τn− P fV(2)◦τn− P fτn◦V(2)◦τn− Glτn+1. The parametrix P maps in a controlled manner between functions with certain asymptotic behaviour which means that we can derive from the asymptotic properties of the terms on the right hand side of (3) its effect on the asymptotic behaviour of τn+1. Furthermore it is an essential property of Glthat the operator maps onto a space with specific asymptotic type. Therefore the asymptotic type of Glτn+1is fixed and does not depend on τn+1. Thus we have full control on the asymptotic properties of the right hand side of (4) and consequently on the asymptotic type of the iterated pair-amplitude τn+1. Let us briefly outline how one gets the singular asymptotic behaviour of individual Goldstone diagrams within this framework. The iteration scheme discussed in the previous paragraph can be actually further decomposed by taking into account the linearity of the differential operator A. Instead of solving the first iterated equation (2) as a whole, let us consider the decomposition Aτ1,1=−fV(2)τ0,Aτ1,2= fV(1)τ0,Aτ1,3=−fV(2)◦τ0,Aτ1,4=−fτ0◦V(2)◦τ0, 676Oberwolfach Report 13/2018 τττ a)b)c)d) τττττ e)f)g) Figure 1.Goldstone diagrams taken into account by our CC model. from which one recovers the first iterated solution via the sum (5)τ1= τ0+ τ1,1+ τ1,2+ τ1,3+ τ1,4, where each term actually corresponds to an individual Goldstone diagram. According to Goldstone’s theorem, cf. [11], each of these Goldstone diagrams is represented by a connected graph. In the next iteration step, one can further use the decomposition (5) to construct the interaction terms on the right hand side. For each new term on the right hand side obtained in such a manner one can again solve the corresponding equation which leads to the decomposition of the second iterated solution into Goldstone diagrams τ2= τ0+ τ1,1+ τ1,2+ τ1,3+ τ1,4+ τ2,1+ τ2,2+· · · This process can be continued through any number of successive iteration steps. Therefore from a diagrammatic point of view, iteration schemes correspond to the summation of an infinite series of Goldstone diagrams which represents a pairamplitude. Our main result concerning Goldstone diagrams is summarised in the following theorem where τRPAcan be any Goldstone diagram which appears in the RPA-CC model, depicted in Fig. 1. We want to emphasise that this result refers to the real RPA-CC model without requiring any of the simplifications discussed before. Theorem 1.Goldstone diagrams of RPA-CC pair-amplitudes1can be considered as kernel functions of classical pseudo-differential operators without logarithmic terms in their asymptotic expansion. Classical symbols corresponding to GoldstoneR diagrams, i.e., σRPA(x, η) :=e−izητRPA(x, z) dz, belong to symbol classes Sclp with p≤ −4. The asymptotic expansion of a Goldstone diagram τRPAwith symbol σRPA∈ Sclp, expressed in spherical coordinates (z, θ, φ), is given by X (6)τRPA(x, z)∼τp−j(x, z, θ, φ)modulo C∞(R3× R3), 0≤j 1 By abuse of notation, we refer to Goldstone diagrams τRPAeither with respect to (x1, x2) or (x, z) variables, with x := x1, z:= x1− x2, i.e., τRPA(x1, x2) ≡ τRPA(x, z). Mathematical Methods in Quantum Chemistry677 with j−p−3XXl τp−j(x, z, θ, φ) = zj−p−3gj,lm(x) Ylm(θ, φ), l=0m=−l j−p−l even where functions gj,lmbelong to C∞(R3). The symbol class of a diagram τRPAcan be determined in the following manner (i) Remove all ladder insertions in the diagram. (ii) Count the number of remaining interaction lines n. Then the corresponding symbol of the diagram τRPAbelongs to the symbol class Scl−4n. References
[126] R. J. Bartlett and M. Musial, Coupled-cluster theory in quantum chemistry, Rev. Mod. Phys., 79, (2007), pp. 291–352.
[127] Y. V. Egorov, B.-W. Schulze, Pseudo-Differential Operators, Singularities, Applications, Birkh¨auser, 1997.
[128] H.-J. Flad, G. Harutyunyan, R. Schneider and B.-W. Schulze, Explicit Green operators for quantum mechanical Hamiltonians. I. The hydrogen atom, manuscripta math., 135, (2011), pp. 497–519.
[129] H.-J. Flad, G. Harutyunyan and B.-W. Schulze, Singular analysis and coupled cluster theory, PCCP, 17, (2015), pp. 31530–31541.
[130] H.-J. Flad, G. Harutyunyan and B.-W. Schulze, Asymptotic parametrices of elliptic edge operators, J. Pseudo-Differ. Oper. Appl., 7, (2016), pp. 321–363.
[131] H.-J. Flad, G. Flad-Harutyunyan and B.-W. Schulze, Explicit Green operators for quantum mechanical Hamiltonians. II. Edge singularities of the helium atom, preprint: arXiv:1801.07552 [math-ph].
[132] H.-J. Flad and G. Flad-Harutyunyan, Singular analysis of RPA diagrams in coupled cluster theory, submitted to ESAIM: M2AN, preprint: arXiv:1801.07573 [math-ph].
[133] S. Fournais, M. Hoffmann-Ostenhof, T. Hoffmann-Ostenhof and T. Østergaard Sørensen, Analytic structure of many-body Coulombic wave functions, Comm. Math. Phys., 289, (2009), pp. 291–310.
[134] G. Harutyunyan and B.-W. Schulze, Elliptic Mixed, Transmission and Singular Crack Problems, EMS Tracts in Mathematics, Vol. 4, European Math. Soc., (2008).
[135] A. Laestadius and S. Kvaal, Analysis of the extended coupled-cluster method in quantum chemistry, preprint: arXiv:1702.04317v1.
[136] I. Lindgren and J. Morrison, Atomic Many-Body Theory, Springer, Berlin, (1986).
[137] T. Rohwedder, The continuous coupled cluster formulation for the electronic Schr¨odinger equation, ESAIM: M2AN, 47, (2013), pp. 421–447.
[138] T. Rohwedder and R. Schneider, Error estimates for the coupled cluster method, ESAIM: M2AN, 47, (2013), pp. 1553–1582.
[139] R. Schneider, Analysis of the projected coupled cluster method in electronic structure calculation, Num. Math., 113, (2009), pp. 433–471.
[140] B.-W. Schulze, Boundary Value Problems and Singular Pseudo-Differential Operators, Wiley, (1998). 678Oberwolfach Report 13/2018 Dynamical low-rank approximation Christian Lubich The talk reviewed differential equations on manifolds of matrices or tensors of low rank. They serve to approximate, in a low-rank format, large time-dependent matrices and tensors that are either given explicitly via their increments or are unknown solutions of high-dimensional differential equations, such as multi-particle time-dependent Schr¨odinger equations. Recently developed numerical time integrators are based on splitting the projector onto the tangent space of the low-rank manifold at the current approximation. In contrast to all standard integrators, these projector-splitting methods are robust with respect to the presence of small singular values in the low-rank approximation. This robustness relies on geometric properties of the low-rank manifolds: in each substep of the algorithm, the approximation moves along a flat subspace of the low-rank manifold. In this way, high curvature due to small singular values does no harm. The talk is based on work done intermittently over the last decade with Othmar Koch, Bart Vandereycken, Ivan Oseledets, Emil Kieri and Hanna Walach. References
[141] E. Kieri, C. Lubich and H. Walach, Discretized dynamical low-rank approximation in the presence of small singular values, SIAM J. Numer. Anal., 54, (2016), pp. 1020–1038.
[142] O. Koch and C. Lubich, Dynamical low-rank approximation, SIAM J. Matrix Anal. Appl., 29, (2007), pp. 434–454.
[143] O. Koch and C. Lubich, Dynamical tensor approximation, SIAM J. Matrix Anal. Appl., 31, (2010), pp. 2360–2375.
[144] C. Lubich, Time integration in the multiconfiguration time-dependent Hartree method of molecular quantum dynamics, Applied Mathematics Research eXpress, (2015), pp. 311–328.
[145] C. Lubich and I. V. Oseledets, A projector-splitting integrator for dynamical low-rank approximation, BIT Numer. Math., 54, (2014), pp. 171–188.
[146] C. Lubich, I. V. Oseledets and B. Vandereycken, Time integration of tensor trains, SIAM J. Numer. Anal., 53, (2015), pp. 917–941.
[147] C. Lubich, B. Vandereycken and H. Walach, Time integration of rank-constrained Tucker tensors, SIAM J. Numer. Anal., 56, (2018), to appear. Minimal Surface Hopping Caroline Lasser This talk aims at the mathematical underpinning of molecular quantum dynamics beyond the Born–Oppenheimer approximation. Born-Oppenheimer approximation describes the vibrational quantum dynamics of a molecule given the system’s electronic eigenvalue surfaces. Let us denote by x and y the nuclear and the electronic coordinates of a molecule, respectively, and by Em(x), m≥ 0, the electronic eigenvalue surfaces, that are the eigenvalues of Mathematical Methods in Quantum Chemistry679 the electronic Hamiltonian Hel(x). Assuming that all the nuclei have identical mass M > 0, the Born-Oppenheimer equations take the form  2 (1)i ∂tψm(x, t) =−∆xψm(x, t) + Em(x)ψm(x, t). 2M These equations yield an adequate description of vibrational quantum dynamics as long as the derivatives of the electronic eigenfunctions um(x) with respect to the nuclear coordinates x are not too large. In almost every molecular system, however, there are avoided or true crossings of electronic eigenvalue surfaces. That is, there are nuclear configurations x such that the gap gmn(x) =|Em(x)− En(x)| between different eigenvalue surfaces becomes small or even zero. The derivatives of the electronic eigenfunctions are intimately linked to the eigenvalue differences, since hum(x, y)| ∇xun(x, y)iy hum(x, y)| (∇xHel(x))un(x, y)iy. =− Em(x)− En(x) Hence, in the presence of small or even vanishing eigenvalue gaps the single surface equations (1) have to be replaced by coupled evolution equations. A particularly simple and widely popular way to deal with multiple surface dynamics is surface hopping. It combines the following two constitutive elements: a) classical equations of motion for nuclear positions along the electronic surfaces, M ¨x(t) =−∇xEm(x(t)), b) hopping between different electronic surface dynamics. The aim of the talk is to present a minimalist hopping mechanism and its mathematical justification according to [4, 3]. Minimal surface hopping monitors the time-evolution of a simple non-Born-Oppenheimer indicator as for the example the evolution of the surface gap along a classical trajectory, t7→ gmn(x(t)). Whenever a local minimum is attained, say at time t = t∗, the probability of switching to the other surface dynamics is determined according to the dynamic Landau-Zener rate s! πgmn(x(t))3, pmn= exp− 2 d2g dt2mn(x(t))t=t∗ see [1] for an application to the twelve-dimensional dynamics of the ammonia cation. The mathematical underpinning of minimal surface hopping draws from the theory of microlocal normal forms for Schr¨odinger systems with avoided and actual eigenvalue crossings, see for example [2]. In a nutshell, microlocal normal form 680Oberwolfach Report 13/2018 theorems construct invertible operatorsTlocand coupling parameters δlocsuch that  2E n(x)Tloc−1  ≈ −i ∂t+δtδloc loc−t The resulting ordinary differential equation tδ i  ˙ϕ(t) =locϕ(t) δloc−t is then used for a quantitative analysis of the hopping mechanism, while the almost explicit construction of the invertible operatorTlocallows to retransfer the key data back to the Schr¨odinger system. References
[148] A. Belyaev, W. Domcke, C. Lasser and G. Trigila, Nonadiabatic nuclear dynamics of the ammonia cation studied by surface hopping classical trajectory calculations, J. Chem. Phys., 142, (2015), 104307.
[149] C. Fermanian-Kammerer and P. G´erard, Mesures semi-classique et croisements des modes, Bull. Soc. Math. France, 130, (2002), pp. 123–168.
[150] C. Fermanian-Kammerer and C. Lasser, An Egorov Theorem for avoided crossings of eigenvalue surfaces, Comm. Math. Phys., 353, (2017), pp. 1011–1057.
[151] C. Lasser and S. Teufel, Propagation through Conical Crossings: an Asymptotic Semigroup, Comm. Pure Appl. Math., 58, (2005), pp. 1188–1230. The density to potential mapping in quantum dynamics: recent results Robert van Leeuwen (joint work with Michael Ruggenthaler, Markus Penz and Søren Nielsen) We address the question whether the particle density of a quantum many-particle system is uniquely determined by the external one-body potential and which properties this external potential must have to generate a certain density. The question is of particular importance to time-dependent density-functional theory [1]. The setting of the problem is the time-dependent Schr¨odinger equation (TDSE) for a system of N particles: i∂tΨ(x, t) = ˆH(t)Ψ(x, t)Ψ(x, t0) = Ψ0(x) where Ψ is the many-particle wave-function, ˆH(t) is a time-dependent Hamilton operator (dependent on a time variable t∈ R) specified below and x = (x1, . . . , xN) is a collection of N space-spin coordinates in which xi= (ri, σi) where ri∈ R3and σia discrete spin variable that assumes the values±1. The initial state at time t0is denoted by Ψ0. The wave-function is required to be antisymmetric under permutation of the N variables xi, i.e. ifPx = (xP(1), . . . , xP(N )) is a permutationP of the variables xithen Ψ(Px, t) = (−1)|P|Ψ(x, t) with|P| the Mathematical Methods in Quantum Chemistry681 number of transpositions out which the permutation is composed. The Hamilton operator is composed out of three terms H(t) = ˆˆT + ˆV (t) + ˆW in which ˆT represents the kinetic energy operator, ˆV (t) the time-dependent external potential and ˆW the inter-particle interactions. The operator ˆT is essentially given by the 3N -dimensional Laplacian operator ∆3N, explicitly we have XN T =ˆ−∇2i=−1∆ 223N i=1 where∇iis the gradient for variable ri. The operators ˆV (t) and ˆW are defined by XNXN V (t) =ˆv(ri, t),W =ˆw(ri− rj) i=1i>j and represent sums over one-body and two-body potentials respectively. The external potential v(r, t) is dependent on time and represents, for example, the attractive Coulomb potential of fixed external nuclear point charges plus an added laser field. The two-body interaction w(r) is typically taken to be a repulsive Coulombic interaction but we do not specify the function spaces for which the TDSE is well-defined here (see [2, 3] for mathematical details). We consider two relevant physical quantities, namely the density n(r, t) and the current density j(r, t) which are defined as XZ n(r, t) = NdN −1x|Ψ(x, t)|2 σ=±1 NZ 2idN −1x[Ψ∗(x, t)∇Ψ(x, t) − (∇Ψ∗(x, t))Ψ(x, t)] where dN −1xdenotes an integration over all variables ri(i = 2, . . . , N ) and a sum over the corresponding spin variables σiexcept for one variable, which we denote by x = (r, σ). The aim is to determine in which way the density is determined by the external potential, provided we fix the operators ˆT and ˆW as well as the initial state. We notice that changing v(r, t) by a purely time-dependent function v(r, t)→ v(r, t) + C(t) will only change the solution of the TDSE by a phaseR factor, Ψ(x, t)→ Ψ(x, t) exp (−ittd¯tC(¯t)) and leaves the density (as well as any 0 other physical observable) invariant. We therefore consider the equivalence class of potentials (gauge class), which we denote by [v(r, t)], of potentials that differ by a purely time-dependent function. The central theme of this presentation is the one-to-one correspondence between the density and initial state and the gauge class of external potentials: (n(r, t), Ψ0(x))↔ [v(r, t)] 682Oberwolfach Report 13/2018 In other words, the knowledge of the density and the initial state gives enough information to recover the external potential up to a purely time-dependent function. A physicist’s proof of this statement was first presented by Runge and Gross [4] although the one-to-one correspondence was already stated by Peuckert [5]. A rigorous mathematical justification of the Runge-Gross argument was given by Fournais et al. [6] but only for the class of potentials and initial states that are infinitely differentiable in space and time and which, in particular, excludes the Coulomb interaction between the particles. We have therefore in recent works [7, 8, 9, 2] pursued a different proof strategy to extend the domain of allowed potentials and initial states. This strategy is based on the quantum fluid formulation of the quantum many-body problem in the form: (1)∂tn(r, t) =−∇ · j(r, t) (2)∂tj(r, t) =−n(rt)∇v(r, t) − Q(r, t) in which the three spatial components of the vector Q(r, t) are given by X3 Qk(r, t) =∂iTik(r, t) + Wk(r, t) i=1 where ∂iis a spatial derivative and Tikand Wkare defined as NXZ Tik(r, t) =dN −1x[∂iΨ∗(x, t)∂kΨ(x, t) + (i↔ k)] −1∂ 22i∂k|Ψ(x, t)|2 σ=±1 Z Wk(r, t) = NdN −1x∂kw(r− r2)|Ψ(x, t)|2 Eqs. (1) and (2) represent a quantum fluid formulation of the many-body problem. Combining these two equations gives an equation that directly relates the density and the external potential: (3)−∇(n(r, t)∇v(r, t)) = q(r, t) − ∂t2n(r, t) where we defined q(r, t) =∇ · Q(r, t). This equation has been used in fixedpoint scheme to construct the potential from a given density and to extend the domain of potentials and initial states beyond the ones presented in [6]. The details of this procedure are given in [7, 8, 9, 2] and the procedure was numerically implemented [10, 11, 12]. An important recent theorem in this context was proven by M. Penz in [3] in the framework of non-autonomous evolution equations in Banach spaces which provides an important regularity theorem for the TDSE as a special case. The incorporation of this result in a rigorous treatment of the fixed-point proof based on Eq. (3) is subject of current research. References
[152] C. A. Ullrich, Time-Dependent Density-Functional Theory: Concepts and Applications, Oxford University Press, Oxford, (2012).
[153] M. Penz, The Density-Potential Mapping in Quantum Dynamics, PhD Thesis University of Innsbruck, arXiv:1610.05552, (2016). Mathematical Methods in Quantum Chemistry683
[154] M. Penz, Regularity for evolution equations with non-autonomous perturbations in Banach spaces, arXiv:1801.03361, (2018).
[155] E. Runge and E. K. U. Gross, Density-Functional Theory for Time-Dependent Systems, Phys. Rev. Lett., 52, 997, (1984).
[156] V. Peuckert, A new approximation method for electron systems, J. Phys. C.: Solid State Phys. 11, 4945, (1978).
[157] S. Fournais, J. Lampart, M. Lewin and T. Østergaard-Sørensen, Coulomb potentials and Taylor expansions in time-dependent density-functional theory, Phys. Rev. A, 93, 062510, (2016).
[158] M. Ruggenthaler and R. van Leeuwen, Global fixed-point proof of time-dependent density-functional theory, EPL, 95, 13001, (2011).
[159] M. Ruggenthaler, K. J. H. Giesbertz, M. Penz and R. van Leeuwen, Density-potential mappings in quantum dynamics, Phys. Rev. A, 85, 052504, (2012).
[160] M. Ruggenthaler, M. Penz and R. van Leeuwen, Existence, uniqueness, and construction of the density-potential mapping in time-dependent density-functional theory, Topical Review, J. Phys. Condens. Matter, 27, 203202, (2015).
[161] S. E. B. Nielsen, M. Ruggenthaler and R. van Leeuwen, Many-body quantum dynamics from the density, EPL, 101, 33001, (2013).
[162] S. E. B. Nielsen, M. Ruggenthaler and R. van Leeuwen, Quantum Control of ManyBody Systems by the Density, arXiv:1412.3794, (2014).
[163] S. E. B. Nielsen, M. Ruggenthaler and R. van Leeuwen, Numerical construction of the Runge-Gross mapping, Eur. Phys. J. B, to be published, (2018). Recent advances on the domain decomposition strategy for implicit solvation models Benjamin Stamm Continuum Solvation models (CSMs) [1, 2, 3, 4, 5] are nowadays part of the standard toolbox of computational chemists. Historically, in the quantum chemistry community, polarizable CSMs such as the Polarizable Continuum Model (PCM) [1] or the Conductor-like Screening Model (COSMO) [6] have been developed as a cheap, but physically sound way, to include solvation effects in the quantum mechanical (QM) description of a molecule and its properties. As the computational cost is usually dominated by the solution of the QM equations, the computational performance of CSM have not been historically taken into much consideration, as the setup and solution of the CSM equations, which are equivalent in some way to solving Poisson’s equation in a heterogeneous dielectric medium (see below), has always been assumed to be a negligible additional computational cost. Due to advances in hardware, more efficient implementations and the spread of linear scaling techniques within quantum chemistry, such an assumption started to be less and less true in the last decade. Further, the diffusion of multiscale methods such as quantum mechanics/molecular mechanics (QM/MM) has made large to very large systems accessible to computational chemists. For such systems, the computational cost associated with continuum solvation, which scales as the second or even third power of the size of the system, can easily become the real bottleneck of the calculation [9]. 684Oberwolfach Report 13/2018 From a modeling viewpoint, there are two main-ingredients for an implicit solvation model. First, the shape of the solute’s cavity Ω, or equivalently its surface ∂Ω, is introduced and determines the region where the implicit solvent is present, see Figure 1 (left). The second ingredient is the macroscopic description of the solvent that is used to model the (electrostatic) interaction between the solute and the continuum solvent region. More precisely, the electrostatic potential V generated by the solute’s charge ρ is then the solution of (1)−div(ε∇V ) + κ2V = 4πρ,in R3, and the electrostatic interaction energy is given by Z Es=1ρ V. 2R3 Here, the dielectric and Debye-H¨uckel constants of the solvent are given by εsand κsso that 0in Ω,0in Ω, (2)ε =κ = εsin Ωc,κsin Ωc. In this talk, we gave a review of the latest developments of the family of domain decomposition (dd) methods for polarizable continuum models. We started with the derivation of the domain decomposition method for COSMO [7, 8, 9, 10], i.e. with εs=∞ and κs= 0 in (1). Indeed, in this particular case, the problem is reduced to the bounded domain Ω (the solute’s cavity) which is then split into the union of overlapping and possibly scaled Van der Waals balls Ωiand equation (1) is solved inside each ball and coupled through appropriate boundary conditions. Numerical illustrations were underlying the efficiency of the present method compared to the state-of-the-art. This model can then be generalized to the PCM [11, 12] which corresponds to equation (1) with κs= 0 (no ionic screening) still using the Van der Waals cavity but for finite εs. The corresponding domain decomposition approach, but this time for an equivalent integral equation formulation as the problem is posed on the entire unbounded free space R3, was introduced. Recent results on larger molecules indicate the the model for the cavity consisting of the union of balls in form of the (possibly scaled) Van der Waals cavity does not provide phenomenologically accurate results. In turn, the Solvent Excluded Surface (SES) seems to be a favorable alternative. We highlighte the mathematical structure and recent analysis of the SES [14, 15], see also Figure 1 (right) for the SES of caffeine, and continued with a new domain decomposition approach for the PCM that uses the cavity the SES-cavity [13]. We discussed in detail the derivation of the numerical method and presented some numerical tests. Finally, the extension [16] to the full linearized Poisson-Boltzmann equation (1) was briefly discussed. Mathematical Methods in Quantum Chemistry685 Solvent Accessible Surface probe VdW-balls Solvent Excluded Surface Figure 1.An illustration of the different molecular surfaces (left) and the Solvent Excluded Surface (SES) for caffeine (right). References
[164] J. Tomasi, B. Mennucci and R. Cammi, Quantum Mechanical Continuum Solvation Models, Chem. Rev., 105 (8), (2005), pp. 2999–3093.
[165] B. Mennucci, Continuum Solvation Models: What Else Can We Learn from Them?, J. Phys. Chem. Lett., 1 (10), (2010), pp. 1666–1674.
[166] Ch. Cramer and D.G. Truhlar, Implicit Solvation Models: Equilibria, Structure, Spectra and Dynamics, Chem. Rev., 99 (8), (1999), pp. 2161–2200.
[167] B. Mennucci, Polarizable continuum model, WIREs Comput. Mol. Sci., 2 (3), (2012), pp. 386–404.
[168] F. Lipparini and B. Mennucci, Perspective: Polarizable continuum models for quantummechanical descriptions, J. Chem. Phys., 16 (144), (2016), 160901.
[169] A. Klamt and G. Sch¨u¨urmann, COSMO: a new approach to dielectric screening in solvents with explicit expressions for the screening energy and its gradient, J. Chem. Soc., Perkin Trans. 2, 5, (1993), pp. 799–805.
[170] E. Canc‘es, Y. Maday and B. Stamm, Domain decomposition for implicit solvation models, J. Chem. Phys., 139 (5), (2013), 054111.
[171] F. Lipparini, B. Stamm, E. Canc‘es, Y. Maday and B. Mennucci, Fast domain decomposition algorithm for continuum solvation models: Energy and first derivatives, J. Chem. Theory Comput., 9 (8), (2013), pp. 3637–3648.
[172] F. Lipparini, L. Lagard‘ere, G. Scalmani, B. Stamm, E. Canc‘es, Y. Maday, J.-P. Piquemal, M. J. Frisch and B. Mennucci, Quantum calculations in solution for large to very large molecules: A new linear scaling QM/continuum approach, J. Phys. Chem. Lett., 5 (6), (2014), pp. 953–958.
[173] F. Lipparini, G. Scalmani, L. Lagard‘ere, B. Stamm, E. Canc‘es, Y. Maday, J.-P. Piquemal, M. J. Frisch and B. Mennucci, Quantum, classical, and hybrid QM/MM calculations in solution: General implementation of the ddCOSMO linear scaling strategy, J. Chem. Phys., 18 (141), (2014), 184108.
[174] B. Stamm, E. Canc‘es, F. Lipparini and Y. Maday, A new discretization for the polarizable continuum model within the domain decomposition paradigm, J. Chem. Phys., 5 (144), (2016), 054101.
[175] P. Gatto, F. Lipparini and B. Stamm, Computation of Forces arising from the Polarizable Continuum Model within the Domain-Decomposition Paradigm, J. Chem. Phys., 22 (147), (2017), 224108. 686Oberwolfach Report 13/2018
[176] Ch. Quan, B. Stamm, and Y. Maday, Polarizable Continuum Model based on the Solvent Excluded Surface, Math. Models Meth. Appl. Sci, (2018), in press (available as e-print), doi.org/10.1142/S0218202518500331.
[177] Ch. Quan and B. Stamm Mathematical analysis and calculation of molecular surfaces, J. Comput. Phys., 322, (2017), pp. 760–782.
[178] Ch. Quan and B. Stamm, Meshing molecular surfaces based on analytical implicit representation, J. Mol. Graph. Model., 71, (2016), pp. 200–210.
[179] Ch. Quan, B. Stamm and Y. Maday A Domain Decomposition Method for the PoissonBoltzmann Solvation Models, Hal-preprint https://hal.archives-ouvertes.fr/hal-01489262, (2018). How to obtain the interacting Green’s functions of many-body quantum systems from a differential equation? Lucia Reining (joint work with Arjan Berger, Pierluigi Cudazzo, Stefano Di Sabatino, Matteo Guzzo, Giovanna Lani, Bernardo Mendoza, Pina Romaniello, Adrian Stan, Walter Tarantino, Marilena Tzavala and Jianqiang Sky Zhou) Many-body perturbation theory is a powerful approach to describe many properties of materials. Here we are interested in the one-body Green’s function of interacting electron systems, but much of the discussion can be generalized. To treat electrons in the presence of the Coulomb interaction, most often one solves a Dyson equation, namely an integral equation of the type (1)G(1, 2) = G0(1, 2) + G0(1, ¯3)Σ(¯3, ¯4)G(¯4, 2). Here G is the Green’s function of the interacting system, and G0is the Green’s function in absence of interaction. The self-energy kernel Σ contains all interaction effects, and it is in general unknown and has to be approximated. Arguments (i) stand for a space, spin and time argument: (i)R≡ (ri, σi, ti), and overlined arguments are integrated over: A(1, ¯3)B(¯3, 2)≡d3A(1, 3)B(3, 2). The state-of-the art approximation for Σ is Hedin’s GW approximation [1], where the self-energy is the product of the interacting Greens function and the Coulomb interaction, dressed by the linear response of all electrons to the propagation of a particle. Besides density functional theory, this is one of the most widely used approaches for electronic structure calculations in a broad range of materials. However, sometimes the GW approximation and related approaches are not sufficient, in particular when one is interested in structures that appear in excitation spectra exclusively because of the interaction (these structures are called satellites). It also breaks down in the case of strong coupling, where the quasi-particle picture is no longer adequate (i.e. when one cannot at all identify structures in an interacting spectrum with slightly modified peaks in the noninteracting spectrum). Mathematical Methods in Quantum Chemistry687 We explore an alternative route to the calculation of interacting electron Green’s functions. It is based on a set of functional differential equations (2) Gu(1, 2) = G0(1, 2) + G0(1, ¯3)[u(¯3) + vHu(¯3)]Gu(¯3, 2) + ivc(¯3, ¯4)δGu(¯3, 2). δu(¯4+) This set of equations [2] is created by applying a fictitious perturbing potential u to the system, which is then set to zero at the end of the derivation, such that limu→0Gu= G. Here vHuis the Hartree potential, i.e. the mean field potential due to the classical charge distribution of all electrons. All effects due to exchange and correlation are contained in the derivative. This set of equations can be used to generate the perturbation series. Here, instead, we are interested in solving the differential equations directly, in order to obtain non-perturbative expressions. Results along several lines have already been obtained, in particular, we have derived approximations close to the cumulant expansion for the one-body Green’s function [3, 4, 5], and, by projecting the equations onto a simple model that yields a scalar differential equation, we have discussed fundamental questions such as problems of multiple solutions in many-body perturbation theory [6, 7, 8, 9]. However, many important questions remain still open. The main issues brought to the audience of the workshop are the following: • What could be a general strategy (analytical and/or numerical) to solve the set of equations and to select the physical solution? • What is the meaning of the various solutions? Is there a link to phase transitions? • Could one derive a modified set of differential equations with only one, or at least less, solutions? • What is the impact of our findings on many-body perturbation theory (there seem to be problems in a regime of strong interaction)? • How could we restrict the domain of “physical” Green’s functions such that we can parametrize the solution, and derive a set of equations for these parameters? This would eliminate spurious solutions. The discussion concerning these questions has been prepared both using the simple model and the full equations. References
[180] L. Hedin, New Method for Calculating the One-Particle Green’s Function with Application to the Electron-Gas Problem, Phys. Rev., 139, A796, (1965).
[181] L. P. Kadanoff and G. Baym, Quantum Statistical Mechanics, W. A. Benjamin, New York, (1962).
[182] M. Guzzo et al., Valence Electron Photoemission Spectrum of Semiconductors: Ab Initio Description of Multiple Satellites, Phys. Rev. Lett., 107, p. 166401, (2011).
[183] M. Guzzo et al., Multiple satellites in materials with complex plasmon spectra: From graphite to graphene, Phys. Rev. B, 89, pp. 085425, (2014).
[184] J. Zhou et al. Dynamical effects in electron spectroscopy, J. Chem. Phys., 143, p. 184109, (2015). 688Oberwolfach Report 13/2018
[185] G. Lani, P. Romaniello, and L. Reining, Approximations for many-body Green’s functions: insights from the fundamental equations, New J. Phys., 14, p. 013056, (2012).
[186] J. A. Berger et al., Solution to the many-body problem in one point, New J. Phys., 16, p. 113025, (2014).
[187] A. Stan et al., Unphysical and physical solutions in many-body theories: from weak to strong correlation, New J. Phys., 17, p. 093045, (2015).
[188] W. Tarantino et al., Self-consistent Dyson equation and self-energy functionals: An analysis and illustration on the example of the Hubbard atom, Phys. Rev. B, 96, p. 045124, (2017). Ab Initio Dynamical Mean Field Theory Reinhold Schneider (joint work with Michael Lindsey and Lin Lin) Finite dimensional model: Full CI and discrete Fock spaces As an underlying model we consider the electronic Schr¨odinger equation in full CI spaces or discrete Fock spaces. Given a single particle space of dimension d Vd= span{(r, s) 7→ ϕi(r, s)∈ C : i = 1, . . . , d} which is spanned by orthonormal spin orbitals and is a finite dimensional subspace of the Sobolev spaceVd⊂ V = H1(R3,±12)⊂ L2(R3,±12), the corresponding N-particle full CI space is spanned by Slater determinants. It is the d-fold antisymmetric tensor product ofVd, denoted by ^N VF CIN:=Vd. i=1 In the sequel, we will consider the full CI spaces with N−1, N and N +1 particles. To unify the notation, we may embed these full CI spaces in the discrete Fock space Hd:=⊕dN =0VF CIN. The many body Schr¨odinger Hamiltonian in the discrete Fock-space, or in the corresponding full CI spaces is taken as the the underlying model, for further approximation or model reduction. In second quantization notation the electronic Schr¨odinger operator is a two particle Hamiltonian given by XdXd (1)H = h + HII=fpqa†pap+up,qr,sa†ra†sapaq p,q=1p,q,r,s=1 with the well known single and two electron integrals fqpand up,qr,s. Presently, we focus our interest on the following eigenvalue problems: HΨ0= E0Ψ0, E0= E0N∈ R , Ψ0∈ VNF CI, HΨN ±1p= EpN ±1ΨN ±1p, ΨN ±1n∈ VN ±1F CI. All operators H(z), with z∈ R, under consideration are supposed to be hermitian! Mathematical Methods in Quantum Chemistry689 Single Particle Greens Function The frequency depending Green’s functions z7→ G(z) = (gi,j(z)), as a function of a complex parameter z∈ C is defined by gi,j(z) :=hΨN0, a†i(zI + E0I + H)−1ajΨN0i + hΨN0, ai(zI + E0I− H)−1a†jΨN0i , for almost all z∈ C. We are interested in the poles z = E0− EpN −1of G(z). The self-energy for the Hamiltonian H = h + HII= a†h(z)a + HII, can be defined by (2)−Σ(z) := G−1(z)− G−10(z) = G−1(z)− ((z + E0)I− h) . Hence the poles of the Green’s function are the solutions of the following quasiparticle equations ((z + E0)I− (h + Σ(z)))xp= 0 These equations constitute non-linear eigenvalue problems,where the eigenvalue parameter z = EN− E0appears in a nonlinear form zI− Σ(z). Let us observe that the quasi-particle equations are low dimensional surrogates of the original linear eigenvalue problems, that were posed in extremely high-dimensional Full CI spaces. If we are interested in these solutions it seems to be advisable to approximate Σ(z) rather than G(z). Schur Complement The concept of Schur complement is of fundamental importance for the understanding of dynamical mean field theory. It shed already light into the nature of the exciton equations. −1 Cd:= Cn⊕ Cm=VA⊕ VB, G :=AB=S−1∗ CD∗∗ where the Schur complement S is given by S= A− ∆ = A − BD−1C, ∆ = BD−1C∈ Cn×n. It is obvious that, if D−1(z) exists: G−1(z)∈ Cd×dis invertible⇔ S(z) ∈ Cn×nis invertible. Which suggests that the Schur complement S(z)∈ Cn×non the subspace Cn provides a good surrogate for G−1(z). Fragmentation We decompose the system into K fragments Ak: k = 1, . . . K, by a decomposition of the basis set. We have in mind that this corresponds to a spatial decomposition if the orthonormal orbital basis functions ϕiare localized. We define the local single particle spaces and the decomposition of the single particle space: Vd=⊕Kk=1VAk=VA⊕ VAc, whereVAk:= span{ϕi: i∈ Ak} . 690Oberwolfach Report 13/2018 Given a matrix G(z) = (ga,b(z))a,b=1,...,dwith G(z) = (zI− (h + Σ(z)))−1for almost all z∈ C , we denote its restriction to a fragment subspace by GA(z) := PAG(z)PA:= (gi,j(z))i,j∈A, and hybridization ∆A(z) := G−1A(z)− PAG−1(z)PA For a given self energy z7→ Σ(z) and single particle operator hA, the corresponding z7→ ΣA(z), ∆A(z) are well defined. The Schur complement is given by G−1A(z) = ((z + E0)δi,j− (hi,j+ Σi,j(z) + ∆i,j(z)))i,j∈A Since PAG(z)PA= SA−1(z) on fragment A, the Schur complement SA(z) forms an ideal surrogate for the whole system zI− (h(z) + Σ(z)) on Vdand also for H on Fd Basic Assumption.The locality of the self-energy Σ is the basic assumption for the DMFT approximation. (Σi,j)i,j∈A=: ΣA, Σi,j(z) = 0 if i∈ A, j ∈ Ac. i.e. Σ is block-diagonal. Impurity problem For each fragment A∈ B, n = nA, the Schur complement covers the influence of the environment. We want to mimic this influence by a quantum mechanical model system. To this end, we consider an impurity system A∪ B with virtual bath B = BA, dimVB= m = mA, on the single particle levelVd=VA⊕ VBthis reads as  zIn+m− h − Σimp(z) =zIn− hA− ΣA(z)−T† −TzIm− d The corresponding Schur complement si,j(z) = zδi,j− (hi,j+ Σi,j(z) + ∆i,j(z)) , i, j∈ A , has a hybridization term ti,νtν,jX(Xµ)i,j cν− z=cµ− z ν=d+1µ where the matrices Xµ 0 are positively semi-definite and its rank µ is the multiplicity of cµ. This corresponds to a quantum mechanical system in the FockspaceFimp=FA⊗ FB≃ C2n+mwhere  H = HA⊗ I + a†hAt∗a:F tdimp→ Fimp and HA= HII|Ais the many particle operator at fragment A. The following theorem shows that the self-energy on the bath components vanishes. Mathematical Methods in Quantum Chemistry691 Lemma(M. Lindsey & L. Lin). Locality of self-energy  ⇒ Σimp:=ΣimpA0(z) 00 Given the impurity problem, the Green’s function z7→ GimpA(z) and self-energy term z7→ ΣimpA(z) can be computed for each fragment A = Ak, k = 1, . . . , K. Since we want that these quantities reflect the physics of the global system, we assume that GimpA(z) = GA(z), or equivalently ΣA+ ∆A= ΣimpA+ ∆impAto holds on each fragment. In a self-consistent iteration cycle, one assumes ΣA= ΣimpA, ∆A= ∆impA. It is not obvious that the self-consistency can be achieved, because of the particular structure of the hybridization. The following definitions and results of M. Lindsey provide an positive answer. Definition.Let M∈ N, z 7→ f(z) ∈ CM×Mis in PRPRP iff XK ∃K ∈ N, s.t. f(z) =Xkc1, c k=1k− zk∈ R and 0  Xk∈ CM×M, Theorem 1(M. Lindsey). There exist rational functions z7→ h(z) which cannot be represented or approximated by functions the class PRPRP. I.e. the class PRPRPis not universal. Theorem 2.The Greens functions z7→ G(z), z7→ G(z) ∈ PRPRP , z 7→ GA(z)∈ PRPRP Theorem 3.For the self-energy functions z7→ Σ(z), and z 7→ ΣA(z) there exists (a static) ˆh∈ CM×Ms.t. z7→ Σ(z) + ˆh ∈ PRPRP Theorem 4.The hybridization functions z7→ ∆A(z) z7→ ∆A(z)∈ PRPRP i.e. all are rational functions. Since z7→ g(z), gA(z), Σ(z), ΣA(z), ∆A(z), are rational functions, they are analytic outside the poles. Moreover, they are defined by finitely many values zi6∈ R, i = 1, . . . , ˜P . Choosing a sampling set Z := {zi6∈ R , i = 1, . . . , P } , (a popular choice is that zi= iyi, yi∈ R, are on the imaginary axis) one can compute ΣAk(zi), ∆Ak(zi) for all fragments A = Ak∈ B in a self consistent cycle Σ, ΣA⇒ ∆A= ∆impA⇒ HAimp⇒ ΣimpA= ΣA⇒ Σ To compute the parameters ti,kfor the impurity problem, on each fragment, we propose a new optimization formulation (Xi,j, cν) = argmin{kδℓ,j−Xνk2+ αXKkX cν− ziνk∗: 0 Xν} i=1ν=1ν=1 692Oberwolfach Report 13/2018 P wherekXk∗is the nuclear norm, and Xi,j=kti,ktj,k. For given poles cν∈ R, this is a convex optimization problem. Unique continuation principle and Hohenberg-Kohn theorem Aihui Zhou The Hohenberg-Kohn theorem addresses that the external potential of a manybody electronic system is a unique functional of the electronic density, apart from a trivial additive constant [3]. Since “all properties of the system are completely determined given only the ground state density” [11], the Hohenberg-Kohn theorem is the theoretical basis of density functional theory (DFT). We note that the mathematical version of the Hohenberg-Kohn theorem is presented in [10]. We understand that DFT has significantly impacted on modern science and engineering (see, e.g., [3, 11, 12, 13, 15, 21]). However, we observe from the literature (c.f., e.g., [5, 6, 7, 10, 14, 17, 18, 19] and references cited therein) that DFT is not entirely elaborated yet. For instance, the existing proofs of HohenbergKohn theorem usually assume directly or indirectly that the particle wavefunction does not vanish on a set of positive measure. We see that the particle wavefunction does not vanish on a set of positive measure is unclear in a real system (c.f. [14, 19]). We refer to [1, 2, 3, 5, 7, 8, 9, 10, 14, 17, 18] and references cited therein for discussions on the Hohenberg-Kohn theorem. In particular, we mention that a modified version of the Hohenberg-Kohn theorem that the external potential of electronic system is a unique functional of the electronic density for Coulomb potentials is presented in [18], in which a direct proof using the Fundamental Theorem of Algebra is provided. We see that no any constant is involved in the modified version [18]. In this talk, we present and prove the Hohenberg-Kohn theorem that if the two potentials differ by more than just a constant on an open set, then the associated densities must be different when the associated wavefunctions do not vanish on an open set [20], which is based a unique continuation principle (see, e.g., [16]). References
[189] P. W. Ayers, S. Golden, and M. Levy, Generalizations of the Hohenberg-Kohn theorem: I. Legendre transform constructions of variational principles for density matrices and electron distribution functions, J. Chem. Phys., 124 (2006), pp. 054101–054107.
[190] N. Hadjisavvas and A. Theophilou, Rigorous formulation of the Kohn-Sham theory, Phy. Rev. A., 30 (1984), pp. 2183–2186.
[191] P. Hohenberg and W. Kohn, Inhomogeneous gas, Phys. Rev., 136 (1964), pp. B864–B871.
[192] W. Kohn and L. J. Sham, Self-consistent equations including exchange and correlation effects, Phys. Rev. A, 140 (1965), pp. 4743–4754.
[193] E. S. Kryachko, On the original proof by reductio ad absurdum of the Hohenberg-Kohn theorem for many-electron Coulomb systems, Int. J. Quantum Chem., 103 (2005), pp. 818– 823.
[194] E. S. Kryachko, On the proof by reductio ad absurdum of the Hohenberg-Kohn theorem for ensembles of fractionally occupied states of Coulomb systems, Int. J. Quantum Chem., 106(2006), pp. 1795–1798. Mathematical Methods in Quantum Chemistry693
[195] E. S. Kryachko and E. V. Lude˜na, Density functional theory: Foundations reviewed, Phys. Rep., 544 (2014), pp. 123–239.
[196] M. Levy, University variational functionals of electron densities, first-order density matrices, and natural spin-orbitals and solutions of the v-representability problems, Proc. Nat. Acad. Sci. USA, 76 (1979), pp. 6062–6065.
[197] M. Levy, Electron densities in search of Hamiltonians, Phys. Rev. A, 26 (1982), pp. 1200– 1208.
[198] E. H. Lieb, Density functionals for Coulomb systems, Int. J. Quantum Chem., 24 (1983), pp. 243–277.
[199] R. Martin, Electronic Structure: Basic Theory and Practical Methods, Cambridge University Press, London, 2004.
[200] R. van Noorden, B. Maher and R. Nuzzo, The top 100 papers, Nature, 514 (2014), pp. 550–553.
[201] R. G. Parr and W. T. Yang, Density-Functional Theory of Atoms and Molecules, Oxford University Press, Oxford, 1989.
[202] R. Pino, O. Bokanowski, E. V. Lude˜na and R. L. Boada, A re-statement of the Hohenberg-Kohn theorem and its extension to finite subspaces, Theor. Chem. Account, 118 (2007), pp. 557–561.
[203] S. Redner, Citation statistics from 110 years of physical review, Physics Today, June 2005, pp. 49–54.
[204] M. Schechter and B.Simon, Unique continuation for Schrodinger operators with unbounded potentials, J. Math. Anal. Appl., 77 (1980), pp. 482–492.
[205] W. Szczepanik, M. Dulak and T. A. Wesolowski, Comment on “On the original proof by reductio ad absurdum of the Hohenberg-Kohn theorem for many-electron Coulomb systems”, Int. J. Quantum Chem., 107 (2007), pp. 762–763.
[206] A. Zhou, Hohenberg-Kohn theorem for Coulomb type systems and its generalization, J. Math. Chem., 50 (2012), pp. 2746–2754.
[207] A. Zhou, Some open mathematical problems in electronic structure models and calculations (in Chinese), Sci. Sin. Math., 45 (2015), pp. 929–938.
[208] A. Zhou, A mathematical aspect of Hohenberg-Kohn theorem, arXiv:1709.07118.
[209] http://www.nature.com/nmat/journal/v15/n4/pdf/nmat4619.pdf. The localization dichotomy for gapped periodic Schr¨odinger operators: ordinary vs Chern insulators Gianluca Panati (joint work with Domenico Monaco, Adriano Pisante and Stefan Teufel) The understanding of transport properties of quantum systems out of equilibrium is a crucial challenge in statistical mechanics. A long term goal is to explain the conductivity properties of solids starting from first principles, as e. g. from the Schr¨odinger equation governing the dynamics of electrons and ionic cores. While the general goal appears to be beyond the horizon, some results can be obtained for specific models, in particular for independent electrons in a periodic or random background. As a general paradigm, in this case the electronic transport properties are related to the spectral type of the Hamiltonian and to the (de-) localization of the corresponding (generalized) eigenstates. However, when periodic systems are 694Oberwolfach Report 13/2018 considered, the Hamiltonian operator has generically purely absolutely continuous spectrum2. Therefore, one needs a finer notion of localization, which allows for example to predict when a crystal, in the absence of any external magnetic field, exhibits a zero transverse conductivity, as it happens for ordinary insulators, and when a non-vanishing one, as in the case of the recently realized Chern insulators[1, 3] predicted by Haldane [5, 6]. Our main message is that such a finer notion of localization is provided by the rate of decay of composite Wannier functions (CWF) associated to the gapped periodic Hamiltonian operator. Equipped with this notion of localization, we are able to identify two different regimes: (i) whenever the system is time-reversal (TR) symmetric, there exist exponentially localized composite Wannier functions which are associated to the Bloch bands below the Fermi energy, assuming that the latter is in a spectral gap; correspondingly, the Hall conductivity vanishes; (ii) viceversa, as soon as the Hall conductivity is non-zero, as it happens for Chern insulators, the composite Wannier functions are delocalized. The relevant information to discriminate between the two regimes is of topological nature. It is provided by the the triviality of the Bloch bundle associated to the occupied states, that is, the vector bundle over the Brillouin torus whose fiber over k is spanned by the occupied Bloch states at fixed crystal momentum k. In a recent paper [9], we rigorously prove a Localization-Topology Correspondence. We consider a gapped periodic (magnetic) Schr¨odinger operator, and we assume that the Fermi projector corresponds to a non-trivial (magnetic) Bloch bundle, as it may happen when TR-symmetry is broken. For example, one might think of the operators modeling Chern insulators or Quantum Hall systems. The rate of decay of composite Wannier functions changes drastically in this case, from exponential to polynomial. We prove that the optimal decay for a system w = (w1, . . . , wm) of CWFs in a non-trivial topological phase is characterized by the divergence of the second moment of the position operator, defined as XmZ hX2iw≡|x|2|wa(x)|2dx. a=1Rd Heuristically, this corresponds to a power-law decay|wa(x)| ≍ |x|−α, with α = 2 for d = 2 and α = 5/2 for d = 3. The former exponent was foreseen by Thouless [14], who also argued that the exponential decay of the Wannier functions is intimately related to the vanishing of the Hall current. More precisely, we prove—under suitable technical hypothesis—the following statement: Localization-Topology Correspondence: Consider a gapped periodic (magnetic) Schr¨odinger operator. Then it is always possible to construct a system 2 A remarkable exception is the well-known Landau Hamiltonian. Notice, however, that if a periodic background potential is included in the model, one is generically back to the absolutelycontinuous setting. Mathematical Methods in Quantum Chemistry695 w = (w1, . . . , wm) of CWFs for the occupied states such that XmZ (1)|x|2s|wa(x)|2dx < +∞for every s < 1. a=1Rd Moreover, the following statements are equivalent: (a) Finite second moment: there exists a choice of Bloch gauge such that the corresponding CWFs w = (w1, . . . , wm) satisfy XmZ hX2iw=|x|2|wa(x)|2dx < +∞; a=1Rd (b) Exponential localization: there exists α > 0 and a choice of Bloch gauge such that the corresponding CWFsw = (ewe1, . . . ,wem) satisfy XZ e2β|x|| ewa(x)|2dx < +∞for every β∈ [0, α); a=1Rd (c) Trivial topology: the Bloch bundle associated to the occupied states is trivial. Our result can be reformulated in terms of the localization functional introduced by Marzari and Vanderbilt [7, 8], which with our notation reads XmZXmXdZ2 (2)FMV(w) =|x|2|wa(x)|2dx−xj|wa(x)|2dx a=1Rda=1j=1Rd =:hX2iw− hXi2w. In view of the first part of the statement, there always exists a system of CWFs satisfying (1) for fixed s = 1/2, so that the first momenthXiwis finite. Hence, the Marzari-Vanderbilt functional is finite if and only ifhX2iwis. By the second part of the Localization-Topology Correspondence, the latter condition is equivalent to the triviality of the Bloch bundle. The result is in agreement with previous numerical and analytic investigations on the Haldane model [13], and with more recent investigations on the subject [2]. As a consequence, the minimization of FMVis possible only in the topologically trivial case, and numerical simulations in the topologically non-trivial regime should be handled with care: we expect that the numerics become unstable when the mesh in k-space becomes finer and finer [2]. Future possible applications of the Localization-Topology Correspondence go beyond the realm of crystalline solids, including superfluids and superconductors [11, 16], and tensor network states [12]. In view of that, we hope that our results will trigger new developments in the theory of superconductors and of many-body systems, and possibly in other areas of physics. 696Oberwolfach Report 13/2018 References
[210] A. J. Bestwick, E. J. Fox, X. Kou, L. Pan, K. L. Wang and D. Goldhaber-Gordon, Precise quantization of the anomalous Hall effect near zero magnetic field, Phys. Rev. Lett., 114, 187201, (2015).
[211] E. Canc‘es, A. Levitt, G. Panati and G. Stoltz, Robust determination of maximallylocalized Wannier functions, Phys. Rev. B, 95, 075114, (2017).
[212] C. Z. Chang et al., High-precision realization of robust quantum anomalous Hall state in a hard ferromagnetic topological insulator, Nature Materials, 14, 473, (2015).
[213] D. Fiorenza, D. Monaco and G. Panati, Construction of real-valued localized composite Wannier functions for insulators, Ann. Henri Poincar´e, 17, pp. 63–97, (2016).
[214] F. D. M. Haldane, Model for a Quantum Hall effect without Landau levels: condensedmatter realization of the “parity anomaly”, Phys. Rev. Lett., 61, 2017, (1988).
[215] M. Z. Hasan and C. L. Kane, Colloquium: Topological Insulators, Rev. Mod. Phys., 82, pp. 3045–3067, (2010).
[216] N. Marzari and D. Vanderbilt, Maximally localized generalized Wannier functions for composite energy bands, Phys. Rev. B, 56, pp. 12847–12865, (1997).
[217] N. Marzari, A. A. Mostofi, J. R. Yates, I. Souza and D. Vanderbilt, Maximally localized Wannier functions: Theory and applications, Rev. Mod. Phys., 84, 1419, (2012).
[218] D. Monaco, G. Panati, A. Pisante and S. Teufel, Optimal decay of Wannier functions in Chern and Quantum Hall insulators, Commun. Math. Phys., 359, pp. 61–100, (2018).
[219] G. Panati and A. Pisante, Bloch bundles, Marzari-Vanderbilt functional and maximally localized Wannier functions, Commun. Math. Phys., 322, pp. 835–875, (2013).
[220] S. Peotta and P. T¨ormaSuperfluidity in topologically nontrivial flat bands, Nature Commun., 8, 8944, (2015).
[221] N. Read, Compactly-supported Wannier functions and algebraic K-theory, Phys. Rev. B, 95, 115309, (2017).
[222] T. Thonhauser and D. Vanderbilt, Insulator/Chern-insulator transition in the Haldane model, Phys. Rev. B, 74, 235111, (2006).
[223] D. J. Thouless, Wannier functions for magnetic sub-bands, J. Phys. C, 17, pp. L325–L327, (1984).
[224] D. J. Thouless, M. Kohmoto, M. P. Nightingale and M. den Nijs, Quantized Hall conductance in a two-dimensional periodic potential, Phys. Rev. Lett., 49, pp. 405–408, (1982).
[225] M. Tovmasyan, S. Peotta, P. T¨orma and S. D. Huber, Effective theory and emergent SU(2) symmetry in the flat bands of attractive Hubbard models, Phys. Rev. B, 94, 245149, (2016). A classical statistical mechanics approach to understanding Green’s function methods and the Luttinger-Ward formalism Michael Lindsey (joint work with Lin Lin) The Luttinger-Ward (LW) formalism [1] is an important component of Green’s function theories in quantum many-body physics. The LW functional Φ[G] provides the formal foundation for bold diagrammatic perturbation theory to all orders and is used to formally derive widely used numerical schemes such as the self-consistent Hartree-Fock approximation, the GW approximation [2], the dynamical mean-field theory (DMFT) [3, 4], and a number of its recent extensions such as the DMFT+GW method [5] and the dynamical cluster approximation [6]. Mathematical Methods in Quantum Chemistry697 Nonetheless, the very existence of the LW functional is currently under debate, with theoretical and numerical evidence favoring the contrary in the past few years for fermionic systems [7, 8, 9, 10]. Such failure has profound theoretical and practical implications. It suggests that many practically used Green’s function methods for computing static or dynamic properties might fail in unpredictable ways. In particular, even in the perturbative regime where bare diagrams converge, the bold diagrams may fail to converge or converge to the wrong quantity [7]. In this talk, we provide the first rigorously justified LW formalism, in the context of classical statistical mechanics, or Euclidean lattice field theory (such as the ϕ4theory [11, 12]). Due to an exact correspondence between the Feynman diagrammatic expansion of lattice field theory and that of quantum many-body physics [13, 11], Euclidean lattice field theory retains the valuable structural information of diagrammatic expansions. Meanwhile, it avoids a key theoretical challenge of the fermionic setup, in the sense that the Green’s function in the Euclidean lattice field theory, defined as a two-point correlator function, has a clearly defined domain, namely the set of positive definite matrices. Hence this work represents a key step towards understanding and potentially remedying the LW formalism and Green’s function methods for fermionic systems. Our theory also proves for the first time the widely used bold diagrammatic expansion, interpreted as an asymptotic series for approximating the LW functional. Independently, our adaptation of the LW formalism to a new setting may provide new insight into the study of Euclidean field theories. The material discussed in this talk is available in [14], and more detail from the mathematical perspective will be provided in future publications. References
[226] J. M. Luttinger and J. C. Ward, Ground-state energy of a many-fermion system. II, Phys. Rev., 118 (1960), p. 1417.
[227] L. Hedin, New method for calculating the one-particle Green’s function with application to the electron-gas problem, Phys. Rev., 139 (1965), p. A796.
[228] A. Georges, G. Kotliar, W. Krauth and M. J. Rozenberg, Dynamical mean-field theory of strongly correlated fermion systems and the limit of infinite dimensions, Rev. Mod. Phys., 68(1996), p. 13.
[229] M. Potthoff, Non-perturbative construction of the Luttinger-Ward functional, Condens. Matter Phys., 9 (2006), pp. 557–567.
[230] S. Biermann, F. Aryasetiawan and A. Georges, First-principles approach to the electronic structure of strongly correlated systems: Combining the GW approximation and dynamical mean-field theory, Phys. Rev. Lett., 90 (2003), p. 086402.
[231] P. Staar, T. Maier and T. C. Schulthess, Dynamical cluster approximation with continuous lattice self-energy, Phys. Rev. B, 88 (2013), p. 115101.
[232] E. Kozik, M. Ferrero and A. Georges, Nonexistence of the Luttinger-Ward Functional and Misleading Convergence of Skeleton Diagrammatic Series for Hubbard-Like Models, Phys. Rev. Lett., 114 (2015), p. 156402
[233] R. Elder, Comment on “Non-existence of the Luttinger-Ward functional and misleading convergence of skeleton diagrammatic series for Hubbard-like models”, arXiv:1407.6599, (2014). 698Oberwolfach Report 13/2018
[234] W. Tarantino, P. Romaniello, J. A. Berger and L. Reining, Self-consistent Dyson equation and self-energy functionals: An analysis and illustration on the example of the Hubbard atom, Phys. Rev. B, 96 (2017), p. 045124.
[235] O. Gunnarsson, G. Rohringer, T. Sch¨afer, G. Sangiovanni and A. Toschi, Breakdown of Traditional Many-Body Theories for Correlated Electrons, Phys. Rev. Lett., 119 (2017), p. 056402.
[236] D. J. Amit and V. Martin-Mayor, Field theory, the renormalization group, and critical phenomena: graphs to computers, World Scientific Publishing Co Inc., 2005.
[237] J. Zinn-Justin, Quantum Field Theory and Critical Phenomena, Clarendon Pr., 2002.
[238] A. L. Fetter and J. D. Walecka, Quantum theory of many-particle systems, Courier Corp., 2003.
[239] L. Lin and M. Lindsey, Variational structure of Luttinger-Ward formalism and bold diagrammatic expansion for Euclidean lattice field theory, Proc. Natl. Acad. Sci., 115 (2018), p. 2282. Generalized Kohn-Sham iteration on Banach spaces Andre Laestadius (joint work with Markus Penz, Erik I. Tellgren, Michael Ruggenthaler, Simen Kvaal and Trygve Helgaker) We here propose an infinite dimensional and general Kohn-Sham formulation on Banach spaces. Furthermore, a weak form of convergence by means of the Optimal Damping Algorithm is established. Let X be a real Banach space and X∗its dual. It is here assumed that X is strictly convex and X∗is uniformly convex. For x∗∈ X∗, define an energy functional E : X∗→ R and consider the variational problem of finding the system state x∈ X that corresponds to minimal energy of (1)E(x∗) = inf{ ˜F (x) +hx∗, xi | x ∈ ˜X}. The functional ˜F (x) stands for all energy contributions of internal effects, and hx∗, xi represents the potential energy that is seen as an external and controllable effect. The domain of ˜F is limited to a certain set of “physical” densities ˜X that by themselves usually do not form a linear space. It is here not given that ˜F is convex and weakly lower semicontinuous (l.s.c.). Instead (borrowing notation from density-functional theory) we introduce a universal Lieb-type functional [1], F : X→ R ∪ {+∞}, F (x) = sup{E(x∗)− hx∗, xi | x∗∈ X∗}. Now, such a functional is by construction convex and weakly l.s.c and has F (x) = +∞ whenever x ∈ X ˜X (corresponding to “unphysical” densities). It is key here that the principal problem stays the same when switching from ˜F to F . To find the state x of minimal energy in (1) by just relying on F we have to determine (2)arg min{F (x) + hx∗, xi | x ∈ X}. Mathematical Methods in Quantum Chemistry699 We do not assume that F is differentiable. The aim is to obtain the ground-state energy E(x∗) from (1) as well as a corresponding minimizer x from (2): −x∗∈ ∂F (x), x ∈ ∂E(x∗). To meet that end, we make use of the Moreau-Yosida regularization and follow [2] where this approach was applied to density-functional theory on Hilbert spaces. Set φ =12k · k2on both X and X∗and let J : X→ X∗be the duality map. For f convex and weakly l.s.c., define fε(x) = inf{f(y) + ε−1φ(x− y) | y ∈ X}. The infimum in the definition of fεis always uniquely attained. This unique minimizer for any given x gives rise to the definition of the so-called proximal mapping, proxεf(x) := arg min{f(y) + ε−1φ(x− y) | y ∈ X}. Applying the Moreau-Yosida regularization to F gives a Fεthat is Fr´echet differentiable and∇Fε(x) = ε−1J(x− proxεf(x)). We moreover set Eε(x∗) = inf{Fε(x) +hx∗, xi | x ∈ X}. It then holds E(x∗) = Eε(x∗) + εφ(x∗). The Kohn-Sham idea [3] is based on the study of a simpler problem. Assume the existence of a reference functional ˜F0: ˜X→ R that belongs to a so-called Kohn-Sham system that captures parts of the real system’s internal physics. This leads to a variational problem that is supposedly easier to solve. Just like with ˜F we derive the regularized functionals Fε0and Eε0and set up the reference problem in analogous fashion −xKS∗∈ ∇Fε0(x)⇐⇒ x ∈ ∂E0ε(xKS∗), −x∗∈ ∇Fε(x)⇐⇒ x ∈ ∂Eε(x∗). Note that the Kohn-Sham system and the real system yield the same state, i.e., the minimizer state x is the same in both cases. It is necessary to choose a different and at this stage completely undetermined potential for the reference system, the so-called Kohn-Sham potential xKS∗. We remark that the problem of “v-representability” plays no role in the regularized setting: By differentiability of Fεand Fε0, potentials always exist xKS∗=−∇Fε0(x),x∗=−∇Fε(x). To determine the Kohn-Sham potential we take the difference of xKS∗=−∇Fε0(x) and x∗=−∇Fε(x), which yields the first step in a self-consistent iteration scheme. The next iteration step towards xKS∗is then given by x∗i+1= x∗+∇Fε(xi)− ∇Fε0(xi). 700Oberwolfach Report 13/2018 Using the Optimal Damping Algorithm [4] we find the steepest descent direction of Fε(x) +hx∗, xi and minimize Fε(x) +hx∗, xi along this direction. Our main result is the following: Theorem 1(Regularized Kohn-Sham iteration). Set x∗1= x∗and solve x1∈ ∂E0ε(x∗1). Iterate i = 1, 2, . . . according to: (a) x∗i+1= x∗+∇Fε(xi)− ∇Fε0(xi) stop if x∗i+1=−∇Fε0(xi) = xKS∗, (b) select ˜xi+1∈ ∂Eε0(x∗i+1), (c) choose ti∈ (0, 1] maximally such that for xi+1= xi+ ti(˜xi+1− xi) one has h∇Fε(xi+1) + x∗, ˜xi+1− xii ≤ 0. The descending sequence{Fε(xi) +hx∗, xii}iconverges as a sequence of real numbers to eε(x∗) = inf{Fε(xi) +hx∗, xii} ≥ Eε(x∗). i Thus, the shifted eε(x∗) + εφ(x∗) is an upper bound for the ground-state energy E(x∗). References
[240] E. H. Lieb, Density functionals for coulomb systems, Int. J. Quant. Chem., 24 (3), 1983.
[241] S. Kvaal, U. Ekstr¨om, A. Teale and T. Helgaker, Differentiable but exact formulation of density-functional theory, J. Chem. Phys., 140, 18A518 (2014).
[242] W. Kohn and L. J. Sham, Self-Consistent Equations Including Exchange and Correlation Effects, Phys. Rev., 140, A1133 (1965).
[243] E. Canc‘es, SCF algorithms for Kohn-Sham models with fractional occupation numbers, J. Chem. Phys., 114, 10616 (2001). Tensor-based modelling of long range electrostatic potentials in many-particle systems Venera Khoromskaia The novel tensor-structured numerical methods appeared as a result of bridging of the algebraic tensor decompositions originating from chemometrics and of the nonlinear approximation theory on separable representation of multivariate functions and operators [9]. The prior results on theory of the low-rank tensor-product approximation of multivariate functions and operators based, in particular, on sinc-quadrature techniques [4] were a significant background for the advanced approaches in scientific computing. The tensor numerical techniques are based on the representation of d-variate functions and operators on large n⊗dgrids in the rank-structured tensor formats which provide O(dn) complexity of numerical calculations instead of O(nd) by conventional methods. A starting point was the Hartree-Fock solver based on the tensor-structured calculation of the two-electron integrals, and of the Laplace and nuclear potential operators using a Gaussian-type basis, discretized on n× n × n 3D Cartesian Mathematical Methods in Quantum Chemistry701 grids [10, 7]. The tensor approach enables calculation of the 3D convolution integral operators in O(n log n) complexity, with the representation of the convolution kernelkxk1, x∈ R3by a low-rank canonical tensor [8], XR (1)PR=p(1)q⊗ p(2)q⊗ p(3)q∈ Rn×n×n. q=1 Thus for all operators in the non-linear 3D integro-differential Hartree-Fock equation the 3D analytical integration is completely avoided, since it is substituted by the grid-based tensor algorithms in 1D complexity. Recently the rank-structured tensor approach suggested a progress in the numerical treatment of the long-range electrostatic potentials in many-particle systems. Let us consider a sum of single Coulomb potentials on a finite L× L × L lattice in a volume box Ω0= [−b/2, b/2]3, XZ k1,k2,k3∈Kkx − a1(k1, k2, k3)k,x∈ ΩL∈ R3. The assembled tensor summation method applies to the lattice potentials defined on n×n×n 3D Cartesian grid. It reduces the calculation of the collective potential sum over a rectangular 3D lattice, XXXR PcL=Wν(k)Pe=W(k)(pe(1)q⊗ ep(2)q⊗ ep(3)q)∈ Rn×n×n, k∈K3k1,k2,k3∈Kq=1 to the summation (assembling) of shifted directional vectors of the canonical tensor representation for a single Newton kernel [5], XRXXX (3)PcL=(W(k1)ep(1)q)⊗ (W(k2)pe(2)q)⊗ (W(k3)pe(3)q). q=1k1∈Kk2∈Kk3∈K Here eP(similar to (1)) represents the single Newton kernel on a twice larger grid and Wν(k)= Wk1⊗ Wk2⊗ Wk3is the shift-and-windowing (onto ΩL) separable transform along the k-grid. For rectangular finite 3D lattices, the rank of the resulting sum is proven to be the same as for the R-term canonical reference tensor Perepresenting a single Coulomb potential. For lattices with multiple vacancies the tensor rank is increased by a small factor [6]. The numerical cost and storage size for 3D lattices are bounded by O(RLn) and O(Rn), respectively, where n is the univariate grid size. Figure 1 presents the assembled canonical vectors for the cluster of 16× 8 × 2 Hydrogen nuclei, with a distance 1.4 bohr between nuclei and the cross-section of its collective electrostatic potential. The method of grid-based assembled tensor summation of the electrostatics potentials on L×L×L finite lattices provides the computational complexity of the order of O(L), and calculation of the interaction energy in O(L2) [7]. Note that the traditional approaches, like Ewald-type summation [3] exhibit O(L3) complexity. 702Oberwolfach Report 13/2018 0.15 0.25 0.2 0.1 0.15 0.1 0.05 0.05 00 -20020-20-1001020 0.1 0.08 0.06 0.04 0.02 0 -10-50510 Figure 1.Assembled x-, y- and z-axis canonical vectors for a cluster of 16×8×2 Hydrogen atoms and the collective electrostatic potential. Recent range-separated (RS) tensor format introduced in [1] applies to many particle systems with rather generally located potentials. These can be the electrostatic potentials of large atomic systems like bio-molecules or the multidimensional scattered data modeled by radial basis functions. The main advantage of the RS tensor format is that the partition into the long- and short-range parts is performed just by sorting skeleton vectors in the tensor representation of the generating Newton kernel, PR= PRs+ PRl. It was proven [1] that the sum of long range contributions from all particles in the collective potential is represented by a low-rank canonical/Tucker tensor with a rank which only logarithmically depends on the number of particles N . The representation complexity of the short range part is O(N ) with a small constant independent on the number of particles. The RS-canonical tensor format specifies the class of d-tensors A∈ RPn1×···×ndwhich can be represented as a sum of a rankR canonical tensor U =Rk=1ξku(1)k⊗ · · · ⊗ u(d)k∈ Rn1×···×ndand a (uniform) PN cumulated canonical tensor bU=ν=1cνUνgenerated by U0with rank(U0)≤ R0, XRXN (4)A=ξku(1)k⊗ · · · ⊗ u(d)k+cνUν. k=1ν=1 Mathematical Methods in Quantum Chemistry703 Figure 2.Left: the free space electrostatic potential of a small biomolecule computed by RS tensor format. Right: the shortrange part of the RS tensor. The basic tools here are the canonical-to-Tucker algorithm and the reduced higher order SVD (RHOSVD) introduced in [10]. The error of the RS tensor representation is defined by the ε-truncation threshold for the Tucker-to-canonical algorithm in the rank reduction scheme. Figures 2 and 3 illustrate the cross-sections of the collective free space electrostatic potential for the 379 atomic biomolecule [2] computed by the RS tensor representation. Figure 3.The long-range part of the collective free space electrostatic potential of a small biomolecule (left) and the error of RS representation (right). References
[244] P. Benner, V. Khoromskaia and B. N. Khoromskij, Range-Separated Tensor Format for Many-Particle Modeling, SIAM J. Sci. Comp., 40 (2), pp. A1034–A1062, 2018.
[245] P. Benner, V. Khoromskaia, B. N. Khoromskij, C. Kweyu and M. Stein, Regularization of the Poisson-Boltzmann equation by using the range-separated tensor format, Manuscript, 2018. 704Oberwolfach Report 13/2018
[246] P.P. Ewald, Die Berechnung optische und elektrostatischer Gitterpotentiale, Ann. Phys 64, 253 (1921).
[247] I. P. Gavrilyuk, W. Hackbusch and B. N. Khoromskij, Data-Sparse Approximation to Operator-Valued Functions of Elliptic Operator, Math. Comp. 73, (2003), pp. 1297–1324.
[248] V. Khoromskaia and B. N. Khoromskij, Grid-based Lattice Summation of Electrostatic Potentials by Assembled Rank-structured Tensor Approximation, Comp. Phys. Comm., 185, 2014, pp. 3162–3174.
[249] V. Khoromskaia and B. N. Khoromskij, Fast tensor method for summation of long-range potentials on 3D lattices with defects, Numer. Lin. Algebra Appl., 23, pp. 249–271, 2016.
[250] V. Khoromskaia and B. N. Khoromskij, Tensor Numerical Methods in Quantum Chemistry: from Hartree-Fock to Excitation Energies, Phys. Chemistry Chem. Physics, 2015, 17, pp. 31491–31509.
[251] B. N. Khoromskij, Fast and Accurate Tensor Approximation of a Multivariate Convolution with Linear Scaling in Dimension, J. Comp. Appl. Math., 234 (2010), pp. 3122–3139.
[252] B. N. Khoromskij, Tensor Numerical Methods in Scientific Computing, Research monograph, De Gruyter Verlag, Berlin, 2018.
[253] B. N. Khoromskij and V. Khoromskaia, Multigrid Tensor Approximation of Function Related Arrays, SIAM J. Sci. Comput., 31 (4), pp. 3002–3026, 2009. New directions in coupled-cluster theory Thomas Bondo Pedersen (joint work with Gustav Baardsen, Audun Skau Hansen, Karl R. Leikanger, Elisa Rebolini, Lorenzo Maschio and Simen Kvaal) While not the most widely used electronic-structure theory—a position held firmly by Kohn-Sham density-functional theory, coupled-cluster (CC) theory is the highaccuracy method of choice for systems dominated by a single electronic configuration and is colloquially referred to as the “gold standard” in the quantum chemistry community. The application of CC theory, however, has been restricted to stationary electronic states of finite molecular systems. It is perfectly reasonable to expect CC theory to excel at describing stationary states of extended systems and for electron dynamics. Although CC theory exhibits a benign polynomial-scaling computational cost relative to the factorial scaling of full configuration interaction, it is still too expensive for straightforward application to extended systems. We have recently generalized the linear-scaling divide-expand-consolidate (DEC) approach to ground-state correlation energies developed for molecules by Jørgensen and coworkers [1, 2, 3] to extended systems with periodic boundary conditions [4]. The extended DEC (X-DEC) framework attempts to exploit the short-ranged nature of electron correlation in gapped systems using a basis of well-localized (generalized) Wannier functions spanning both the occupied and the unoccupied single-particle spaces. Assigning each Wannier function, unoccupied as well as occupied, to an atomic site A, the CC correlation energy per unit cell may be written as a sum of atomic site energies EAin the central unit cell and interaction energies ∆EABbetween Mathematical Methods in Quantum Chemistry705 the central unit cell sites and all other sites: XX (1)Ecorr=EA+∆EAB AB<A XX (2)EA=Eijab ij∈Aab  XX X (3)∆EAB=+Eijab i∈Ai∈Bab j∈Bj∈A (4)Eabij= tabij− taitbj(ia|jb) ZZ (5)(ia|jb) =ψi(r)ψa(r)|r − r′|−1ψj(r′)ψb(r′) dr′dr where taiand tabijare single- and double-excitation amplitudes. While this expression is exact, the X-DEC approximation attempts to exploit the fast asymptotic decay of the electron repulsion integrals ((ia|jb) ∼ O(R−6)) to truncate summations and thus make calculations feasible. Error control is obtained by solving amplitude equations in systematically increased orbital subspaces for each site A until EAis stable within a predefined tolerance. The converged orbital subspaces are then used for each significant pair A, B to solve for pair amplitudes required to compute ∆EAB. Since all fragment calculations are independent, this not only leads to linear-scaling complexity but also allows for an embarrassingly parallel implementation. Intuitively, one would expect the performance of the X-DEC scheme to improve with increased locality of the Wannier functions. In [4], however, we observed that the orbital subspaces used for solving the atomic site amplitude equations are surprisingly large, limiting the computational performance of the method. To investigate this problem further, we have used the local-correlation algorithm implemented in the Cryscor program [5] to test the convergence of the secondorder Møller-Plesset correlation energy with respect to excitation domain size, comparing unoccupied Wannier functions with projected atomic orbitals (PAOs). The PAOs are obtained by projecting out the occupied states from the atomicorbital basis functions, thus creating a linearly dependent nonorthonormal set that spans the unoccupied orbital space but having significantly greater spread (less locality) than the linearly independent orthonormal Wannier functions spanning the same space. Preliminary results for bulk LiH show that the correlation energy is converged to within 1 mHa using the PAOs on the 6 nearest-neighbor sites only. With unoccupied Wannier functions, however, a domain size comprising 32 atoms (up to fourth-nearest neighbor sites) is required to achieve milli-Hartree accuracy of the correlation energy. From a local correlation perspective, therefore, the requirement that Wannier functions decay exponentially at infinite distance is necessary but not sufficient. Defining a localization functional that leads to a set of unoccupied Wannier functions with exponential asymptotic decay and, simultaneously, captures as much of the correlation energy with as few functions as possible is a subject for further research. 706Oberwolfach Report 13/2018 Going beyond the calculation of stationary electronic states, solving the timedependent CC equations paves the way for highly accurate studies of atto-second processes in strong laser pulses as well as (nonperturbative) treatments of conventional molecular response properties in the frequency domain through Fourier transformations. Relatively little work has been done in this direction [6, 7, 8, 9] and the first step is the integration of the time-dependent CC equations, which are complex analogues of Hamilton’s equations of classical mechanics [10], ∂H (6)˙t = −i(t) ∂λ (7)˙λ = i∂H(t, λ) ∂t where t and λ are amplitude vectors. With the linear parameterization chosen here, λ becomes the (time-dependent) Lagrange multipliers of Helgaker and Jørgensen [11]. The HamiltonianH is nonseparable and due to the high computational cost of evaluating the derivatives on the right-hand side, we have tested the explicit symplectic second-order integrator recently proposed by Tao [12], aiming at long-time stability with (relatively) large time steps. For the He atom exposed to a uniform electric-field with maximum amplitude at 72 as, however, we find that the Tao integrator diverges with time steps above approximately 2.4 as. For comparison, the fourth-order Runge-Kutta method shows a less dramatic drift in energy. −0.0025−3.243e1Tao−0.0025−3.243e1RK4 −0.0030−0.0030 )−0.0035)−0.0035 HaHa E−0.0040 (mE−0.0040 −0.0045−0.0045 −0.0050−0.00500500010000150002000025000 010002000300040005000 time (as)time (as) References
[254] M. Zi´o lkowski, B. Jans´ık, T. Kjærgaard and P. Jørgensen, Linear scaling coupled cluster method with correlation energy based error control, J. Chem. Phys., 133 (2010), 014107.
[255] K. Kristensen, M. Zi´o lkowski, B. Jans´ık, T. Kjærgaard and P. Jørgensen, A Locality Analysis of the Divide–Expand-Consolidate Coupled Cluster Amplitude Equations, J. Chem. Theory Comput., 7 (2011), pp. 1677–1694.
[256] I.-M. Høyvik, K. Kristensen, B. Jans´ık and P. Jørgensen, The divide-expandconsolidate family of coupled cluster methods, J. Chem. Phys., 136 (2012), 014105.
[257] E. Rebolini, G. Baardsen, A. S. Hansen, K. R. Leikanger and T. B. Pedersen, DivideExpand-Consolidate Second-Order Møller-Plesset Theory with Periodic Boundary Conditions, J. Chem. Theory Comput. In press (2018). DOI: 10.1021/acs.jctc.8b00021 Mathematical Methods in Quantum Chemistry707
[258] C. Pisani, M. Sch¨utz, S. Casassa, D. Usvyat, L. Maschio, M. Lorenz and A. Erba, Cryscor: a program for the post-Hartree-Fock treatment of periodic systems, Phys. Chem. Chem. Phys., 14 (2012), pp. 7615–7628.
[259] S. Kvaal, Ab initio quantum dynamics using coupled-cluster, J. Chem. Phys., 136 (2012), 194109.
[260] D. R. Nascimento and A. E. DePrince III, Linear Absorption Spectra from Explicitly Time-Dependent Equation-of-Motion Coupled-Cluster Theory, J. Chem. Theory Comput., 12(2016), pp. 5834–5840.
[261] D. R. Nascimento and A. E. DePrince III, Simulation of Near-Edge X-ray Absorption Fine Structure with Time-Dependent Equation-of-Motion Coupled-Cluster Theory, J. Phys. Chem. Lett., 8 (2017), pp. 2951–2957.
[262] T. Sato, H. Pathak, Y. Orimo and K. L. Ishikawa, Time-dependent optimized coupledcluster method for multielectron dynamics, J. Chem. Phys., 148 (2018), 051101.
[263] T. B. Pedersen and H. Koch, On the time-dependent Lagrangian approach in quantum chemistry, J. Chem. Phys., 108 (1998), pp. 5194–5204.
[264] T. Helgaker and P. Jørgensen, Calculation of geometrical derivatives in molecular electronic structure theory, Methods in Computational Molecular Physics, S. Wilson and G. H. F. Diercksen (Eds.), Plenum, New York (1992), pp. 353–421.
[265] M. Tao, Explicit symplectic approximation of nonseparable Hamiltonians: Algorithm and long time performance, Phys. Rev. E, 94 (2016), 043303. Optimization of the first bands of a Hill’s operator Virginie Ehrlacher (joint work with Athmane Bakhta, David Gontier) The aim of this talk was to present some new considerations on the problem of optimizing the first bands of a periodic Hill’s operator, motivated by materials design applications. The question can be stated as follows: given M functions defined on (−1/2, 1/2), denoted by b1, . . . , bM, is it possible to find an optimal periodic potential V which minimizes the functional XZ1/2 bm(q)− ǫVm,q2dq m=1−1/2 where (−1/2, 1/2) ∋ q 7→ ǫVm,qis the mthband of the 2π-periodic Hill’s operator −12∂xx+ V . It is possible to rigorusly prove an existence result for the optimization of the first band (M = 1) provided that the optimal potential V is searched for in a particular class of singular potentials (Borle measures). Several numerical methods can be proposed for the resolution of the above mentioned optimization problem in practice. 708Oberwolfach Report 13/2018 A Numerical Linear Algebra Perspective of Projection-based Embedding Theory Leonardo Zepeda-N´u˜nez (joint work with Lin Lin) Solving the many body Schr¨odinger equation, even in its approximate form, can be very challenging when the system size becomes large. Under certain circumstances, large quantum systems can be partitioned into two parts: the “system” (or “defect‘”) part containing the degrees of freedom of interest, which need to be treated accurately, and a “bath” (or “environment”) part containing the remaining degrees of freedom that can be treated approximately. In these cases, it naturally becomes desirable to have a numerical method that can restrict the calculation to the system part only, resulting in a much lower computational cost. This intuitive idea has led to the development of various “quantum embedding theories” [3, 2] in which different models have been coupled together resulting in cheaper methods. We focus on the quantum embedding theory within the context of Kohn-Sham density functional theory (KSDFT). In particular, we focus on the recently developed Projection-based Embedding Theory (PET) [1], which we reformulate within the formalism of linear algebra, and extend using a perturbation argument. After proper discretization, KSDFT can be represented as the following nonlinear eigenvalue problem (1)H[P ]Ψ = ΨΛ,P = ΨΨ∗, where H[P ]∈ CN ×Nis a Hermitian matrix. To solve this problem, we need to compute the algebraically lowest Neeigenvalues encoded in the diagonal matrix Λ∈ RNe×Neand the associated eigenvectors Ψ∈ CN ×Ne, which satisfy the orthonormality condition Ψ∗Ψ = INe. P is a spectral projector, usually called the density matrix. The Hamiltonian H[P ] depends nonlinearly on the density matrix P , and (1) needs to be solved self-consistently. We are mainly interested in the accurate computation of the block of P corresponding to the system part, thus sharply reducing the computation cost. To reformulate the PET within the context of linear algebra, we start with the linear case. Within the linear regime, H[P ] = H∈ CN ×Nis a Hermitian matrix, whose eigenvalues are ordered non-decreasingly as λ1≤ λ2≤ · · · ≤ λNe< λNe+1≤ · · · ≤ λN, and we assume that there is a positive energy gap between the Ne-th and (Ne+ 1)-th eigenvalues. Then (1) is the Euler-Lagrange equation of the following energy minimization problem (2)E =infE[P ], where E[P ] := Tr[HP ], P2=P,P∗=P TrP =Ne where the condition P = P2requires P to be a projector. Mathematical Methods in Quantum Chemistry709 In addition, let H0∈ CN ×Nbe a reference Hamiltonian, and P0to be the known minimizer of (3)E0=infTr[H0P ]. P2=P,P∗=P TrP =Ne0 Within this context, PET efficiently computes an approximation to P by leveraging knowledge of P0, under the assumption that the density matrices P and P0 do not differ much in the bath. This can be encoded by splitting P0as (4)P0= P0,b+ P0,s. where P0,band P0,sare both projectors and orthogonal to each other, and by assuming that P can be approximately split in a similar fashion. These assumptions lead to the PET ansatz to the density matrix: (5)P≈ PPET= P0,b+ Ps, where Ps2= Psis also a projector and it is orthogonal to P0,b. Since the rank of P0,bis already Nb, the rank of Psis thus equal to Ns:= Ne− Nb. Thus for the linear problem, PET solves the following constrained minimization problem with respect to Ps: (6)EPET=infTr[H(Ps+ P0,b)]. Ps2=Ps,Ps∗=Ps P0,bPs=0,TrPs=Ns Compared to (2), we find that PET restrains the domain of the density matrices only to those satisfying the ansatz (5). Hence by the variational principle EPET≥ E provides an upper bound for the energy. Computationally, PET reduces the number of eigenpairs that need to be computed in the eigenvalue problem from Neto Ns. Note that the dimension of H0and H must be the same, but Ne0and Ne can be different, and thus the ranks of P0,sand Pscan also be different. This is necessary in the context of KSDFT, where the system part in H can involve different numbers and/or types of atoms from that in H0. In order to compute the minimizer to (6), we define the deflated Hamiltonian, H|B⊥= (I−P0,b)H(I−P0,b), to be the restriction of H to the subspace orthogonal 0 to the range of P0,b, which we denote byB0, and assume that there is a positive gap between the Ns-th and (Ns+ 1)-th eigenvalue of H|B⊥. Then the variational 0 problem (6) has a unique minimizer, which is given by the solution to the following linear eigenvalue problem (7)H|B⊥Ψs= ΨsΛs,Ps= ΨsΨ∗s, 0 where (Ψs, Λs) are the lowest Nseigenpairs of H|B⊥. 0 We show that PET can be understood as an eigenvalue problem of a matrix whose off-diagonal blocks are neglected, resulting in solving two decoupled eigenvalue problems: one involving only the bath, and the other involving the system, as shown in (7). In order to uncover the structure of the problem, we conceptually 710Oberwolfach Report 13/2018 define a unitary N× N matrix W := [Ψ0,b, Ψs, Ψu], where Ψ0,bis an orthogonal matrix in the range of P0,bsuch that Ψ0,bH0Ψ0,b= Λ0,band P0,b= Ψ0,bΨ∗0,b, Ψs is the solution to (7), and Ψuare the unoccupied eigenvectors in (7). The matrix representation of H with respect to the basis W , denoted by HW, can be written as  Ψ∗0,bHΨ0,bΨ∗0,bHΨsΨ∗0,bHΨu HW= W∗HW =Ψ∗sHΨ0,bΨ∗sHΨsΨ∗sHΨu , Ψ∗uHΨ0,bΨ∗uHΨsΨ∗uHΨu  Ψ∗0,bHΨ0,bΨ∗0,bHΨsΨ∗0,bHΨu =Ψ∗sHΨ0,bΛs0 . Ψ∗uHΨ0,b0Λu where we used the orthogonality of Ψsand Ψu, and Λuis a diagonal matrix with the unoccupied eigenvalues of (7) in its diagonal. Computing the eigendecomposition of HWwould yield the exact answer to (2). Instead, PET computes an eigen-decomposition of an auxiliary matrix HWPET. In particular, PET performs two approximations: first, it discards all the off-diagonal blocks, and second, in the first diagonal term it approximates H = H0+ ∆H, by H0, which leads to  Λ0,b00Ψ∗0,bH0Ψ0,b00Ψ∗0,bHΨ0,b00 HWPET=0Λs0 =0Λs0 ≈0Λs0 . 00Λu00Λu00Λu PET can be recast as solving an optimization problem using HWPET, whose minimizer, under certain gap conditions, is given by  INb00 PWPET=0INs0 . 000 When rotated back to the standard basis, the density matrix becomes PPET= W PWPETW∗= Ψ0,bΨ∗0,b+ ΨsΨ∗s= P0,b+ Ps. Using this formalism we can define a perturbation of P = PPET+ δP , where δP is computed from the perturbation Ψ∗0,b∆HΨ0,bΨ∗0,bHΨsΨ∗0,bHΨu δHW= HW− ∆HW=Ψ∗sHΨ0,b00 . Ψ∗uHΨ0,b00 Using the structure of δHW, the Dyson equation, and the Cauchy integral formula, we can deduce an expression of δP , which is given by: (8)δP = δΨ0,bΨ∗0,b+ h.c. where δΨ0,b∈ CN ×Nbsatisfies the equation (9)Q (λi;0,bI− H) Qδψi;0,b= Q(Hψi;0,b),Qδψi;0,b= δψi;0,b. Mathematical Methods in Quantum Chemistry711 Here the projector Q = I− (Ps+ P0,b) = ΨuΨ∗u. λi;0,bis the i-th diagonal element of Λ0,b. ψi;0,b, and δψi;0,bare the i-th columns of Ψ0,b, and δΨ0,b, respectively. h.c. represents for the Hermite conjugation of the first term. Even though the theory for the linear problem seems to be rather complete, its extension to the non-linear case remains mostly unexplored. From a algorithmic point of view, formulating the PET as a variational problem makes the extension of the algorithm to the nonlinear case seamless. Consider a more general energy functional such as (10)E[P ] = Tr[P HL] + EHxc[P ], where HLis a given matrix derived from the discretized Laplacian operator and the electron-nuclei interaction potential. EHxc[P ] is a nonlinear functional of the density matrix P , and it consists of the Hartree, and exchange correlation energy. In such case, the minimization problem can be stated as (2), whose EulerLagrange equation is precisely (1), which in this case takes the form (11)H[P ]Ψ = (HL+ VHxc[P ])Ψ = ΨΛ,P = ΨΨ∗, where VHxc[P ] =δEHxcδP[P ]is called the exchange-correlation potential. The same procedure can be utilized as in the linear case. PET evaluates the modified variational problem by restricting the domain as (12)EPET=infE[Ps+ P0,b], Ps2=Ps,Ps∗=Ps P0,bPs=0,TrPs=Ns which is solved by the following nonlinear eigenvalue problem (13)H[P ]|B⊥Ψs= ΨsΛs,Ps= ΨsΨ∗s,P = Ps+ P0,b. 0 Numerical results for real chemical systems, using KSSolve [4], indicate that the PET approach can provide reasonably accurate approximation to the density matrix and the energy, even when the system is relatively small. References
[266] F. R. Manby, M. Stella, J. D. Goodpaster and T. F. Miller III, A simple, exact density-functional-theory embedding scheme, J. Chem. Theory Comput., 8 (2012), pp. 2564– 2568.
[267] T. Nguyen, A. A. Kananenka and D. Zgid, Rigorous Ab Initio Quantum Embedding for Quantum Chemistry Using Green’s Function Theory: Screened Interaction, Nonlocal SelfEnergy Relaxation, Orbital Basis, and Chemical Accuracy, J. Chem. Theory Comput., 12 (2016), pp. 4856–4870.
[268] Q. Sun and G. K.-L. Chan, Quantum embedding theories, Acc. Chem. Res., 49 (2016), pp. 2705–2712.
[269] C. Yang, J. C. Meza, B. Lee and L. W. Wang, KSSOLV—a MATLAB toolbox for solving the Kohn-Sham equations, ACM Trans. Math. Software, 10 (2009), 110. 712Oberwolfach Report 13/2018 Cubic scaling algorithm for RPA correlation energy using interpolative separable density fitting Kyle Thicke (joint work with Jianfeng Lu) In Kohn-Sham Density Functional Theory (KS DFT), the ground state energy of a system can be written as (1)E = Ts+ Uext+ UH+ Exc, where Tsis the kinetic energy of the non-interacting system, Uextis the energy from the external potential, UHis the Hartree energy, and Excis the so-called exchange-correlation energy. The first three terms are well known, and given the Kohn-Sham orbitals, it is straightforward to compute them. A simple explicit expression for the exchange-correlation energy is not known, so it is approximated in practice [2]. In this work, we consider the Random Phase Approximation (RPA) for the correlation energy [3]. Z∞ (2)EcRPA=1trln(1− χ0(iω)vc) + χ0(iω)vcdω, 4π−∞ where XoccXvirψ∗ (3)χ0(x, y; iω) =j(x)ψk(x)ψk∗(y)ψj(y)+ c.c., ǫj− ǫk− iω jk vcis the Coulomb kernel, and{ψj} are the Kohn-Sham orbitals. We are considering the zero temperature case where the first Noccorbitals are occupied and rest are unoccupied (virtual). Under this approximation, the bottleneck in the computation is the construction of χ0, which takes O(N4) time, where N is a parameter describing the system size. All other steps take only cubic time in N . In order to compute (3) in cubic time, we split up the dependence of j and k in the denominator by using Cauchy’s integral formula, Z! (4)χ0(x, y; iω) =1Xoccψj∗(x)ψj(y)Xvirψk(x)ψ∗k(y)dλ + c.c., 2πiCλ− ǫj+ iωλ− ǫk jk whereC is a negatively oriented curve that encircles all {ǫk}k∈vir, but not any {ǫj± iω}j∈occ. The strategy to obtain (4) has been utilized in [8]. (An alternative strategy, the so-called Laplace transform method, is utilized in [6, 7].) The major contribution of this work is to lower the prefactor in front of the N3in the computational cost. Note that if we use (4), then the computation of (2) will still take time proportional to Ngrid3, where Ngridis the number of grid points. In actual calculations, Ngridis typically very large. Therefore, it would be advantageous for us to write the problem into a smaller basis set. We do this using the newly developed interpolative separable density fitting [4, 5]. Mathematical Methods in Quantum Chemistry713 The traditional density fitting method factors the orbital pair functions as NXaux (5)ψ∗j(x)ψk(x)≈CjkµPµ(x), µ=1 where{Pµ} are user-chosen basis functions and the coefficients Cjkµare found by solving least squares problems (which takes O(N4) time). The interpolative separable density fitting method factors the orbital pair functions as NXaux (6)ψj∗(x)ψk(x)≈ψj∗(xµ)ψk(xµ)Pµ(x), µ=1|{z} coefficient where{Pµ} are the new basis functions (which are computed by the algorithm) and ψ∗j(xµ)ψk(xµ) are the coefficients in front of the basis functions. The splitting of the j and k dependence in the coefficient is crucial for the application to our problem. Additionally, computing the factorization only takes O(N3) time. By writing the problem into this new basis, we reduce the computational complexity of the cubic scaling algorithm from O(Ngrid3) to O(Naux2Ngrid). This is a significant savings since Nauxhas been shown empirically to be much smaller than Ngrid[4, 1]. References
[270] J. Lu and K. Thicke, Cubic scaling algorithms for RPA correlation using interpolative separable density fitting, Journal of Computational Physics, 351, (2017), pp. 187–202.
[271] J. Perdew and K. Schmidt, Jacob’s ladder of density functional approximations for the exchange-correlation energy, AIP Conference Proceedings, Vol. 577, No. 1, AIP, 2001.
[272] X. Ren et al., Random-phase approximation and its applications in computational chemistry and materials science, Journal of Materials Science, 47 (21), (2012), pp. 7447–7471.
[273] J. Lu and L. Ying, Compression of the electron repulsion integral tensor in tensor hypercontraction format with cubic scaling cost, Journal of Computational Physics, 302, (2015), pp. 329–335.
[274] J. Lu and L. Ying, Fast algorithm for periodic density fitting for Bloch waves, Annals of Mathematical Sciences and Applications, 1 (2), (2016), pp. 321–339.
[275] M. Kaltak, J. Klimeˇs and G. Kresse, Low scaling algorithms for the random phase approximation: Imaginary time and Laplace transformations, Journal of Chemical Theory and Computation, 10 (6), (2014), pp. 2498–2507.
[276] M. Kaltak, J. Klimeˇs and G. Kresse, Cubic scaling algorithm for the random phase approximation: Self-interstitials and vacancies in Si, Physical Review B, 90 (5), (2014), 054115.
[277] J. Moussa, Cubic-scaling algorithm and self-consistent field for the random-phase approximation with second-order screened exchange, The Journal of Chemical Physics, 140 (1), (2014), 014107. 714Oberwolfach Report 13/2018 molsturm: Modular electronic structure theory framework Michael F. Herbst (joint work with James E. Avery, Guido Kanschat and Andreas Dreuw) The talk gives an overview of our recent efforts to simplify the investigation of novel types of basis functions for quantum-chemical calculations [1]. Motivated by Coulomb-Sturmians as basis functions for quantum-chemistry simulations the design of the flexible and light-weight quantum-chemical method development framework molsturm [2, 3] is outlined, where new discretisation methods for electronic structure theory calculations can be easily implemented and tested. First a list of desirable properties for a basis in the context of quantum chemistry is reviewed. Amongst other aspects an ideal basis would be able to (1) accurately represent the physics of a chemical system, and (2) allow for systematic improvements aided by error estimates, whilst it (3) gives rise to numerically feasible discretised problems. The predominant basis function type in electronic structure theory, namely contracted Gaussian-type orbitals (cGTO) [4, 5], leads to comparatively simple discretised problems as well as an acceptable accuracy for most applications. These functions do, however, not satisfy the first two aforementioned criteria perfectly. For example, they are not able to properly describe the region of a state where electron and nucleus are close, see figure 1. Our recent research has looked into so-called Coulomb-Sturmians (CS) as an alternative. These exponentially decaying functions are isoenergetic, analytical solutions to a partial differential equation, which is related to the Schr¨odinger equation for a single-electron system by scaling the nuclear attraction potential [6, 7]. They form a complete basis for H1(R3) and are able to correctly reproduce the physical features of the wave function, i.e. both the nuclear cusp as well as proper exponential decay [8]. Figure 1 displays the local energy for the hydrogen atom at small electron-proton distances if selected cGTO and CS basis sets are used for the discretisation. For cGTO discretisations the large fluctuations of the local energy at small distances are indicative of deviations of the obtained approximate eigensolution from being a true eigensolution. This is not improved much for distances less than 0.1 Bohr if a larger cc-pV6Z basis is used. In contrast for CS-based discretisations the behaviour is more uniform and the larger (4, 1, 1) basis has local energies closer to−0.5 than (3, 1, 1) over the full depicted range. Ideally one would not need to restart development from scratch for each new type of basis function, but would be able to utilise already existing quantumchemistry programs as much as possible. A direct implementation of a new basis function type into an existing quantum chemistry programs is, however, typically challenging, since the numerical properties differ substantially between discretisation methods. To give an illustration, figure 2 shows the sparsity as well as the size of the Fock matrices resulting from different discretisation approaches applied to a Hartree-Fock (HF) self-consistent field (SCF) calculation of beryllium. In some cases, like CS-based discretisations, unusual angular-momentum selection rules Mathematical Methods in Quantum Chemistry715 Figure 1.Plot of local energy EL(r) versus radial distance of electron and proton for the hydrogen atom. Shown is the region where electron and proton are close. The exact local energy is−0.5 and deviations can be understood as the relative error of the discretised solution to being a true eigensolution of the Schr¨odinger equation [1, 7]. Q2finite elements,contracted Gaussians,Coulomb-Sturmians adaptively refined 3D gridpc-2 basis set(3, 2, 2) basis with k= 2.0 Figure 2.Examples of Fock matrices in a Hartree-Fock selfconsistent field calculation for beryllium, discretised using finite elements, contracted Gaussians or Coulomb-Sturmians [1, 3]. can reduce computational effort for obtaining the Fock matrix elements further, provided that appropriate evaluation schemes can be used, see [1]. Such deviating numerical demands of different discretisations can be supported by so-called contraction-based methods [1, 3]. In this well-established approach one avoids to store the problem matrices in memory and instead uses matrix-vector 716Oberwolfach Report 13/2018 product expressions as the algorithmic basis. Employing a concept called lazy matrices, where all matrix operations are subject to lazy evaluation, molsturm [2, 3] features an SCF algorithm, which is completely separated from the code dealing with the discretisation. In molsturm it is thus possible to add a new type of discretisation or basis function in a plug and play fashion, i.e. just by implementing the link to the code performing the integral evaluations. Conversely a new algorithm — like an SCF scheme — only needs to be programmed once and immediately can be used with all basis function types available. Thus molsturm facilitates employing new mathematical approaches in practice and comparing them on a common setting to the methods already implemented. Lazy matrices are not necessarily limited to the SCF step or the context of quantum chemistry. The library lazyten [9] allows them to be used in the context of other physical problems as well. It is explicitly not the goal of molsturm to become yet another complete quantum-chemistry package. Much rather the package focuses on solving SCF problems in the most general way possible and handing the results over to thirdparty packages or other user-provided algorithms building on top. For this reason molsturmoffers a readily applicable python interface, where all relevant quantities can be easily obtained for post-processing. Due to the basis-function-independent nature of the SCF algorithm any novel basis function type implemented in molsturmimmediately becomes available to algorithms built on top of this framework. Therefore molsturm can be thought of as a mediator between low-level developments with respect to novel algorithms and discretisation schemes for solving the HF problem on the one hand and high-level quantum-chemical method development: Both ends are abstracted from each other, but can still benefit from mutual progress. References
[278] M. F. Herbst, Development of a modular quantum-chemistry framework for the investigation of novel basis functions, Ph.D. thesis, Heidelberg University (2018). Dissertation submitted.
[279] M. F. Herbst and J. E. Avery, https://molsturm.org, A modular electronic structure theory code. Accessed on 10th April 2018.
[280] M. F. Herbst, A. Dreuw and J. E. Avery, molsturm: A light-weight quantum-chemistry framework for rapid method development not restricted to a particular type of basis function. In preparation.
[281] F. Jensen, Atomic orbital basis sets, Wiley Interdisciplinary Reviews: Computational Molecular Science, 3, 273 (2013).
[282] J. G. Hill, Gaussian basis sets for molecular applications, International Journal of Quantum Chemistry, 113, 21 (2013).
[283] J. Avery and J. Avery, Generalized Sturmians and Atomic Spectra, World Scientific (2006). Mathematical Methods in Quantum Chemistry717
[284] P. E. Hoggan, How Exponential Type Orbitals Recently Became a Viable Basis Set Choice in Molecular Electronic Structure Work and When to Use Them, pp. 199–219, SpringerVerlag, Dordrecht (2009).
[285] J. S. Avery, J. E. Avery, V. Aquilanti and A. Caligiana, Atomic Densities, Polarizabilities, and Natural Orbitals Derived from Generalized Sturmian Calculations, Advances in Quantum Chemistry, 47, 157 (2004).
[286] M. F
This reference list is based on information provided by the publisher or from digital mathematics libraries. Its items are heuristically matched to zbMATH identifiers and may contain data conversion errors. It attempts to reflect the references listed in the original paper as accurately as possible without claiming the completeness or perfect precision of the matching.