# Research interests

My research work focuses on the numerical and theoretical analysis of numerical models and methods in quantum chemistry, whether they come from solid state physics or molecular chemistry. More precisely, they aim to estimate the errors that have several possible origins

1. modeling errors: the equations which govern quantum chemistry are known - the Schrödinger equation is a high-dimensional linear PDE - which can only be solved for systems reduced to a few electrons. The models (Hartree-Fock or Kohn-Sham) reduce the problem to solving a nonlinear system of PDEs in 3D.

2. discretization errors: for a given model, this error depends on the numerical scheme considered (plane waves, finite elements or finite differences).

3. algorithmic errors: at a fixed model and scheme, the chosen algorithm introduces an error which may be difficult to assess. The nonlinear PDEs in quantum chemistry are solved by fixed point methods. The underlying problem being in the vast majority of cases non-convex, the convergence with respect to a numerical tolerance is unclear.

## Tensor Networks

Tensor Networks (TN) have seen rapid and emerging interests. They aim to give a sparse representation of high-dimensional functions. The first instance of TN was Matrix Product States (MPS) or Tensor Trains that were introduced in the late 90s/early 2000s in physics as DMRG (Density Matrix Renormalization Group) to find the ground-state of statistical models (one-dimensional Hubbard or Heisenberg models). MPS has been used for the $N$-body electronic Schrödinger problem where the question of the optimal topology (or best ordering) is not clear.

In the contribution, we investigate the behaviour of the singular values of reshapes of a particular class of functions (Slater determinants). The decay of these singular values drives the approximability of the function by an MPS. For this class of functions, if we denote by $(\sigma_j)$ the singular values of the reshaped tensor $\Psi_{\mu_1 \dots \mu_{L/2}}^{\mu_{L/2+1}\dots \mu_L}$, then we show that there are $2^N$ positive singular values, where $N$ is the number of particles and they satisfy

$\sigma_1 \sigma_{2^N} = \sigma_2 \sigma_{2^N-1} = \sigma_3 \sigma_{2^N-2} = \cdots = \mathrm{Const}.$

The constant is moreover explicit and depends on the ordering of the modes of the tensor. This provides a natural criterion to optimze a TT approximation of the state.

## Linear response in quantum chemistry

Properties of a molecule in its ground-state perturbed by a time-dependent potential can be approximated in the linear regime (assuming that the perturbation is small) using the Kubo formula. To be more precise, if $\psi_0$ is the ground-state of $H = -\Delta + V$ an operator acting on $L^2(\mathbb{R}^d)$, for some observable $V_\mathcal{O}$ and $\psi(t)$ solution to

$i\partial_{t} \psi = H \psi + \varepsilon f(t) V_{\mathcal P} \psi , \quad \psi(0) = \psi_{0},$

the linear response function $K$ defined by

$\langle \psi(t), V_{\mathcal O} \psi(t) \rangle = \Big\langle \psi_{0}, V_{\mathcal O} \psi_{0}\Big\rangle + \varepsilon (K \ast f)(t) + \mathcal{O}(\varepsilon^2),$

is given by the Kubo formula

$K(\tau) = -i \theta(\tau) \Big\langle V_{\mathcal O} \psi_{0}, e^{-i (H-E_{0}) \tau} V_{\mathcal P} \psi_{0}\Big\rangle + {\rm c.c.}.$

The Fourier transform of $K$ is then

$\widehat K(\omega) = \lim_{\eta \to 0^{+}} \Big\langle \psi_{0}, V_{\mathcal O} \Big(\omega +i\eta - (H-E_0)\Big)^{-1} V_{\mathcal P} \psi_{0}\Big\rangle - \Big\langle \psi_{0},V_{\mathcal P} \Big(\omega +i\eta + (H-E_0)\Big)^{-1} V_{\mathcal O} \psi_{0}\Big\rangle.$

The above formula is in practice approximated by taking a positive $\eta$ and by truncating the space to a finite region $(-L,L)^d$. This necessarily yields a singular approximation of the linear response function $\widehat{K}$. In the paper, the regularity of the response function $\widehat{K}$ is determined, depending on the decay of the potential $V$. Moreover, the relationship between $\eta$ and $L$ guaranteeing the convergence of the approximation is shown, supported by numerical experiments.

## The particle-hole random phase approximation in density functional theory

The ground-state of a molecule is given by the lowest eigenvalue of the operator $H_N$ acting on $\bigwedge_{i=1}^N L^2(\mathbb{R}^3)$ with domain $H^2(\mathbb{R}^{3N})$

$H_N(v_\mathrm{ext},w) = \sum\limits_{i=1}^N \Big(-\frac{1}{2}\Delta_{r_i} + v_{\mathrm{ext}}(r_i) \Big) + \sum\limits_{1 \leq i < j \leq N} w(r_i-r_j),$

where $N$ is the number of electrons, $v_\mathrm{ext}$ accounts for the interaction with the nuclei and $w$ is the electron-electron interaction. Solving the high-dimensional eigenvalue problem is quickly intractable for systems with more than a few electrons, but thankfully, the ground state energy can be computed using density functional theory (DFT). DFT states the existence of a universal functional $F$ of the electronic density $\rho(r) = N\int_{\mathbb{R}^{3(N-1)}} |\Psi_0(r,\overline{r})|^2 \, \mathrm{d}\overline{r}$ such that

$E_0 = \min_{\|\Psi\|=1} \langle \Psi, H_N \Psi \rangle = \min_{\rho \in \mathcal{I}_N} F(\rho) + \langle v_{\mathrm{ext}}, \rho\rangle,$

where $\mathcal{I}_N$ is the admissible set of electronic densities.

It is a tremendous reduction of dimensionality, at the expense that the functional $F$ is unknown... Over the past 60 years, great efforts have been put in the design of more and more accurate density functionals, that have been used for systems with tens of thousands electrons, but it appeared that the simple dissociation of the H$_2$ molecule, i.e. a two-electron molecule remains a challenging system for most density functionals. At the dissociation limit, when the hydrogen atoms are pulled away to infinity, we expect that the energy of the system at the limit is twice the energy of a single hydrogen atom. This fails for example for restricted Hartree-Fock as the spin restriction introduces a spurious self-interaction energy that can only be lifted by breaking the spin symmetry. In other words, the correct depiction of the dissociation requires a superposition of two determinantal states, hence it cannot be captured at the level of restricted Hartree-Fock.

phRPA, also simply called RPA, is one of the few ansatz for the correlation energy that is able to dissociate H$_2$ correctly. In the paper with Kyle Thicke, we show that by taking the RPA correlation energy, which can be derived from an adiabatic connection argument that is explained, on top of a restricted Hartree-Fock functional, we have the exact dissociation of the H$_2$ molecule.

## Anderson acceleration

Anderson acceleration, also known as DIIS (Direct Inversion of the Iterative Space) in quantum chemistry, aims to accelerate the convergence of fixed-point iteration $x_{k+1} = g(x_k)$. Such fixed-point problems are ubiquitous, and in quantum chemistry come from the Hartree-Fock or Kohn-Sham equations. DIIS is an extrapolation method where the next iterate is built using the knowledge of previous iterates. The DIIS sequence is given by

$(c^{(k)}_i)_{0 \leq i \leq m} = \argmin_{\sum c_i = 1} \left\|\sum\limits_{i=0}^m c_i f(x_{k-i})\right\|_2 \quad \text{and} \quad x_{k+1} = g\left(\sum_{i=0}^m c^{(k)}_i x_{k-i}\right),$

where $f$ is an error function that locally vanishes only at the solution of the fixed-point problem. A usual choice is $f(x) = x-g(x)$.

An important parameter in the efficiency of the method is the width $m$ of the history/number of iterates kept. Two variants are considered in our publication:

• one where the history increases until the least-square problem is ill-conditioned.

• another where teh size of the history is controlled by the ratio of the residuals $\frac{\|f(x_i)\|}{\|f(x_j)\|}$ for $k-m\leq i,j \leq k$. The heuristic behind this variant is to keep the iterates that are the most relevant to build the next one.

For both variants, we prove local linear convergence and superlinear convergence by tuning the DIIS numerical parameters.

## Solid state physics

In solid-state physics, the periodic boundary conditions motivate a plane-wave discretisation. However because of the Coulomb singularities at the nuclei, the eigenfunctions of the Hamiltonian $-\Delta + V$ have cusps which impedes the convergence. Pseudopotential methods are widely used, where the Coulomb singularities are smoothened. In the PAW method (Projector Augmented-Wave), introduced by Blöchl, the pseudopotential is introduced using an invertible operator acting locally around the nuclei. The PAW equations involve infinitely many terms that are truncated in the numerical treatment which causes a modelling error. This modelling error has been analysed for the one-dimensional model $-\frac{\mathrm{d}^2}{\mathrm{d}x^2} - Z_0 \sum_{k \in \mathbb{Z}} \delta_k - Z_a\sum_{k \in \mathbb{Z}} \delta_{a+k}$ for which explicit formulas of the eigenpairs are available.

Another approach consists in introducing an invertible operator $T$ acting locally around the nuclei. This operator is defined such that it maps smooth functions to functions with cusps satisfying the Kato cusp condition. This operator acts as a preconditonner for a plane-wave discretisation. This method is called the Variational Projector Augmented-Wave and has been analysed for one-dimensional and three-dimensional Hamiltonians.