Exact Density-Functionals with Initial-State Dependence and Memory

We analytically construct the wave function that, for a given initial state, produces a prescribed density for a quantum ring with two non-interacting particles in a singlet state. In this case the initial state is completely determined by the initial density, the initial time-derivative of the density and a single integer that characterizes the (angular) momentum of the system. We then give an exact analytic expression for the exchange-correlation potential that relates two non-interacting systems with different initial states. This is used to demonstrate how the Kohn-Sham procedure predicts the density of a reference system without the need of solving the reference system's Schr\"odinger equation. We further numerically construct the exchange-correlation potential for an analytically solvable system of two electrons on a quantum ring with a squared cosine two-body interaction. For the same case we derive an explicit analytic expression for the exchange-correlation kernel and analyze its frequency-dependence (memory) in detail. We compare the result to simple adiabatic approximations and investigate the single-pole approximation. These approximations fail to describe the doubly-excited states, but perform well in describing the singly-excited states.


I. INTRODUCTION
Time-dependent density-functional theory (TDDFT) [1,2] allows for an exact description of a many-body system in terms of an effective non-interacting system, known as the Kohn-Sham (KS) system. The external potential (known as the KS potential) in the non-interacting system is a functional of the density in such a way that the KS system has exactly the same density as the reference system.
The essential component in the KS construction is the exchange-correlation (xc) potential that contains all nontrivial many-body effects. It depends on the initial states of the interacting and the KS system (initial-state dependence) as well as the density at all previous times (memory). Both features of the xc potential are, however, not well understood and consequently virtually all commonly used approximations neglect them, which in important cases (doubly-excited states, molecular dissociation, charge transfer etc.) can lead to large errors in the calculated properties [1,2]. It is therefore highly desirable to have exact analytical functionals available for model systems that can serve as benchmarks and which can provide insight into how memory and initial state dependence can be incorporated into approximate functionals for real systems.
In this work we explicitly construct such exact analytic functionals that do incorporate initial-state dependence and memory for the case of a quantum ring (QR) with two particles in a singlet state. In Sec. II we will derive functionals with an explicit initial-state dependence for the case that the two particles are non-interacting. These functionals will then be used to construct an ex-plicit expression for the xc potential that connects two non-interacting systems. In Sec. III we will calculate the xc potential for two interacting particles at a specific density. For the same system we will then analytically construct the exact xc kernel of linear-response TDDFT and investigate its frequency-dependence. We conclude in Sec. IV.

II. FUNCTIONALS WITH INITIAL-STATE DEPENDENCE: NON-INTERACTING MODEL SYSTEM
The dynamical properties of many-electron systems, such as molecules or solids is well-described by the solution of the time-dependent Schrödinger equation (TDSE). If we restrict ourselves to external scalar potentials (such as laser fields in the dipole approximation) then the physical properties of an N -electron system evolving from a given initial state |Ψ 0 under the influence of an external scalar potential v(r, t) is determined by the Hamiltonian where ∇ i is the gradient with respect to the spatial coordinate r i and w(|r i − r j |) is the electron-electron interaction (usually chosen to be Coulombic). In molecules and solids the form of the kinetic energy operator and the two-body interactions is always the same, whereas the external potential v varies from system to system. For this reason we will treat v as a variable. Consequently the solutions of the corresponding time-dependent Schrödinger equation (TDSE) i∂ t |Ψ(t) =Ĥ(t)|Ψ(t) can be uniquely labeled by the initial state and the external potential, i.e. the quantum states |Ψ([Ψ 0 , v], t) depend functionally on the initial state and the external potential [25]. However, due to the large number of degrees of freedom of the many-body wave function, a numerical solution of the TDSE is only feasible for small systems. In the KS approach of density-functional theory the interacting many-body problem is mapped onto an effective non-interacting system which considerably reduces the computational effort. The effective potential in these equations is a functional of the density of the system which is defined as is the density operator. A number of observables of large interest, such as the optical absorption spectrum in linear response or the time-dependent dipole moment, are explicitly known as functionals of the density. The basic theorems of DFT actually guarantee that, at least in principle, all observables are a functional of the density. This is a consequence of the fact that the Runge-Gross (RG) theorem [3] and its generalizations [4,5] guarantee (under certain assumptions) that the full many-body wave function is uniquely determined by only knowing its initial state and the density, i.e. the wave function is a functional of the initial state and the density |Ψ([Ψ 0 , n], t) . As a consequence the knowledge of n[Ψ 0 , v] is enough to calculate all physical properties of a many-body system. In Sec. II A we will give an explicit example of this result by analytically constructing the wave-function functional |Ψ([Ψ 0 , n], t) for a specific system. How the density and the initial state determine the external potential of the Hamiltonian of Eq. (1) is then demonstrated in Sec. II B. This result is then employed to give an example of the KS scheme, which can be used to predict the density of a reference system by solving an auxiliary non-interacting problem, by explicitly constructing an initial-state dependent xc potential in Sec. II C.

A. Wave-Function Functional
In this Section we give a non-trivial analytical realization of the wave-function functional |Ψ([Ψ 0 , n], t) with explicit initial-state dependence, for the case of two noninteracting particles on a QR (a one-dimensional system with periodic boundary conditions) of length L.
We assume the non-interacting wave function |Φ(t) to be in a spin-singlet configuration. In a position-spin basis we then make an orbital product Ansatz for the spatial part of the resulting wave function where x and y are the spatial coordinates of the particles along the ring. The full position-spin dependence is obtained by multiplication with the usual anti-symmetric singlet spin-function. Here the orbital ϕ(x, t) satisfies the one-dimensional Schrödinger equation with periodic boundary conditions on the interval [0, L] and starting from the initial state ϕ 0 (x) = ϕ(x, t 0 ) (We adopt the convention that an external potential belongs to a non-interacting system if we use the subindex s).
The periodic boundary conditions on the orbital ϕ(x, t) then correspond to periodic boundary conditions on the norm |ϕ| and quasi-periodic boundary conditions on the phase S, i.e.
for some integer m. Note that the initial orbital ϕ 0 (x) = |ϕ 0 (x)| exp(iS 0 (x)) determines the choice of m since S 0 (x) = S(x, t 0 ) must obey condition (4). To proceed, we use that the density and current of the non-interacting system, are connected by the continuity equation which expresses the local conservation of particles. This is a Sturm-Liouville equation [5] depending parametrically on the time t and thus the density determines the phase function S(x, t) for a given set of boundary conditions (4) and (5). More precisely S is determined uniquely up to a purely time-dependent constant C(t), since the constant function is eigenfunction of the Sturm-Liouville operator in Eq. (8) with eigenvalue zero and also satisfies the boundary conditions (4) and (5). Physically this freedom amounts to the gauge freedom in the potential. Following similar derivations as in reference [5] we find that where we defined η(x, t)η(y, t) Note, these functions are defined only within the interval [0, L] but can be extended periodically outside of it. At t = t 0 this equation determines S 0 (x) in terms of n(x, t 0 ), ∂ t n(x, t 0 ) and m, up to an overall constant and therefore for a given choice of m the density completely determines the initial state ϕ 0 (x) up to a global phase factor e iα . Thus, if we restrict ourselves to the product Ansatz of Eq. (2), there is only a countably infinite number of physically different initial states possible for any given time-dependent density. Obviously these initial states all share the same initial density and time-derivative of the density, but their phases differ. If we compare the resulting currents given by Eq. (7) as functionals of the density and the initial state m we find with the help of Eq. (9) that i.e. the currents differ only by a time-dependent constant. Accordingly the integral of the local velocity fields v(x, t) = j(x, t)/n(x, t) differ exactly by 2π(m − m ′ ). So the density rotates differently around the QR for the different values of m, but in such a way as to yield the same density.
The resulting density-functional for the orbital (and with this the full wave function) is then given by ϕ([m, n], x, t) = n(x, t) 2 exp (iS([m, n], x, t)) .
This is an explicit realization of the RG result. As pointed out before, a direct consequence is that we can calculate all observables of the particles in terms of the density and the choice of initial state only. For instance, the kinetic-energy functional in this case becomes n(x, t) where the first term on the right hand side is the famous Weizsäcker kinetic-energy functional. The second term is an initial-state dependent correction that together with the Weizsäcker term constitutes the exact kinetic energyfunctional.

B. Potential Functional
The basic theorems of TDDFT further establish the uniqueness and existence of a density-potential mapping, i.e. for a given initial state there is a one-to-one correspondence between the external potentials and the densities. This allows for the determination of the external potential that produces a given density by propagation of an initial state, i.e. the external potential is a functional of the initial state and the density v([Ψ 0 , n], r, t). This fact forms the basis of the KS construction, which allows us to determine the density of an interacting system by solving an auxiliary non-interacting problem.
Here we will give an explicit example for the functional v[Ψ 0 , n]. We will rely upon our previous results of two non-interacting particles on a QR. The external potential v s can readily be expressed in terms of the orbital by inverting the Schrödinger Eq. (3) and we find [1,6,7 The last term on the right hand side vanishes as a consequence of the continuity Eq. (8) and we thus find using |ϕ| = n/2 that  [5]. Thus we have analytically defined a density-potential mapping which is also explicitly initial-state dependent. We stress that the periodic boundary conditions on the wave function were essential in deriving Eq. (10). This excludes, for instance, the example of a homogeneous electric field on a ring of constant density given in reference [8] [26].

C. Exchange-Correlation Functional
The functional v[Ψ 0 , n] plays a central role in TDDFT. In practice, however, we are usually not directly interested in this mapping. We are rather interested in the density of a particular system that has a specific external potential v ext . For example, in the case that we want to describe a single molecule in a laser field, the potential v ext is simply given by the Coulombic attraction of the atomic nuclei in the molecule with the addition of the laser field. For a given choice of v ext every observable we want to know is then determined by solution of the TDSE for the given initial state |Ψ 0 . In particular we can calculate the density of the system, which, for future reference, we denote by n 0 (r, t). However, the full solution of the TDSE is usually not feasible in practice, due to the large degrees of freedom that we need to consider. The main idea of the KS construction in TDDFT is to reduce the complexity by mapping the interacting manybody problem to a non-interacting many-body problem with the same density. This leads to one-particle equations that are computationally much easier to deal with. The price we pay for this simplification is that the functional v[Ψ 0 , n] now appears implicitly as part of the xc potential v xc in the KS equations. Below we will present an analytic example of an xc potential for our QR system. However, we will start with a brief description of the KS method and define the KS and xc potentials.
The existence of a density-potential mapping v[Ψ 0 , n] does not depend on the chosen two-body interaction. Specifically this means that we have a density-potential mapping for interacting as well as non-interacting systems. For the case of a non-interacting system this mapping is called v s [Φ 0 , n]. Since in this case we have no twobody interactions the Hamiltonian is then simply given byĤ The initial state |Φ 0 of the non-interacting system is usually chosen to be a single Slater determinant of orbitals ϕ i (r). This allows us to reduce the TDSE for the non-interacting system of 2N electrons to single-orbital equations of the form can be reconstructed from solving Eqs. (11) and (12). In particular, if n = n 0 is the density of an interacting system with external potential v ext and initial state |Ψ 0 then, provided that we chose |Φ 0 such that Φ 0 |n(r)|Φ 0 = Ψ 0 |n(r)|Ψ 0 , the potential v s [Φ 0 , n 0 ] reproduces the density n 0 of the interacting system in a non-interacting system. However, it is clear that the Eqs. (11) and (12) can not predict the density n 0 (r, t) of interest since they contain no information on the interacting system that we are trying to solve. To set up a predictive scheme we need to connect the interacting and the non-interacting system. To do this we introduce the KS potential [28] v If we assume full knowledge of the functionals v[Ψ 0 , n] and v s [Φ 0 , n] then the set of equations does have a unique solution [4,5,9] for a self-consistent density n sc . By definition of v s [Φ 0 , n] the self-consistent density n sc is exactly attained whenever which according to Eq. (13) is precisely satisfied when In turn, this is exactly true when n sc = n 0 as there is a unique potential producing a given density. We therefore see that the set of Eqs. (14) and (15) has exactly a self-consistent solution at the density n 0 of the interacting system with initial state |Ψ 0 and external potential v ext . To make the scheme practical we need to know the (13) or at least have a reasonable approximation for it. The first non-trivial approximation to this expression is given by the classical electrostatic potential of the electrons, i.e. the Hartree potential Usually this approximation is made explicit and the rest is then called the xc potential The KS potential may thus also be written as Therefore, the fundamental approximation in TDDFT is that of the xc potential v xc [Ψ 0 , Φ 0 , n] and the results thus only depend on the quality of this approximation.
However, the xc potential v xc [Ψ 0 , Φ 0 , n] is still a complicated functional that depends on the initial states of both the interacting and non-interacting system (initialstate dependence) and the density at all previous times (memory).
Let us now give an example for the KS and xc potentials for our model system. The construction of these functionals requires the knowledge of the functional v[Ψ 0 , n], which is not explicitly known. However, if the reference system is also non-interacting then v[Ψ 0 , n] = v s [Ψ 0 , n] and v H [n] = 0, and we find that For our case of a QR with two particles in a single-orbital singlet state the functional v s [Φ 0 , n] is given by Eq. (10), and we find where ∂ x S([0, n], x, t) is defined only in terms of n and ∂ t n and corresponds to the spatial derivative of the first term on the right hand side of Eq. (9). Note, the integers m and m ′ play the role of the initial state |Ψ 0 respectively |Φ 0 . The corresponding KS equations are thus . This equation determines the density n(x, t) of the reference system when we prescribe v ext . We note that the xc potential is given only in terms of n and ∂ t n. In contrast, the functional v s ([m, n], x, t) of Eq. (10) that reproduces a prescribed density via propagation of the KS equation also contains a second-order time-derivative of the density (in the term ∂ t S as can be seen with the help of Eq. (9)). We therefore can explicitly see that the second order time-derivative of the density vanishes if we connect the two systems. This is an important fact which sometimes is overlooked in the literature and can lead to misunderstandings about the KS approach [1, 10, 11].

III. FUNCTIONALS WITH MEMORY: INTERACTING MODEL SYSTEM
In the previous section we have constructed functionals that depend only on the density at one time. Although also time-derivatives of the densities appear in the expressions we call these functionals time-local and accordingly they do not exhibit memory. At this point it is useful to give a more precise definition of memory. We first define the xc kernel as the functional derivative of the xc potential, i.e.
Any approximation to the xc potential that depends only locally on the density and its time-derivatives gives rise to an xc kernel that is proportional to time-derivatives of the delta function δ(t − t ′ ). These functions vanish for t = t ′ and therefore have zero memory depth. If the xc kernel is non-zero for t = t ′ we will say that the xc potential has memory. We can find another useful characterization of memory in the case that the functional derivative of Eq. (16) is evaluated at a ground state density. Due to the time translation invariance of the ground state Hamiltonian the kernel f xc will then only depend on the time-arguments through the combination t − t ′ , i.e. f xc (r, t, r ′ , t ′ ) = f xc (r, r ′ , t − t ′ ). We can therefore by means of a Fourier transform define a frequency-dependent xc kernel by f xc (r, r ′ , ω) = dτ e iωτ f xc (r, r ′ , τ ).
In this case memory is characterized by a non-polynomial frequency dependence of f xc (since the Fourier transform of the n-th time derivative of a delta function gives a frequency dependence proportional to ω n ).
We now address the question whether for our QR system we can construct an xc potential with memory. We have seen that the xc potential that arises in the modeling of a non-interacting system by another non-interacting system with a different initial state has no memory (at least not for the product Ansatz used). One way to induce memory is to introduce many-body interactions. However, in the case that the reference system is interacting we do not know v xc [Ψ 0 , Φ 0 , n] as we do not know v[Ψ 0 , n]. However, if we can determine n[Ψ 0 , v ext ] for some |Ψ 0 and a specific external potential v ext , we can still calculate v xc [Ψ 0 , Φ 0 , n, v ext ] as a function of space and time for this density, since v[Ψ 0 , n] is then known, i.e. v[Ψ 0 , n] = v ext . Here we will do this for the case of two particles on a QR of length L with external potential v ext = 0 and which interact via a squared cosine potential, i.e. for the Hamiltonian where λ is the strength of the interaction. In Sec. III B we will then construct the resulting xc potential which takes the simple form since in our example v ext = 0. To give an explicit expression of functionals with memory we will further construct the xc kernel of TDDFT in Sec. III C for this system. The frequency-dependence (memory) of the xc kernel will then be investigated in detail in Sec. III D. Finally, in Sec. III E we will test the validity of the single-pole approximation for this model system. However, to do all these things, it will prove helpful to first explicitly construct all the eigenstates of the Hamiltonian of Eq. (17).

A. Spectrum of the Model System
The eigenfunctions of the Hamiltonian of Eq. (17) can be written as the product of a spatial wave function Ψ(x, y) and a spin-function. We have a spin-singlet (spintriplet) configuration if Ψ(x, y) is (anti)-symmetric with respect to an interchange of x and y, i.e.
Ψ(x, y) = ±Ψ(y, x) (19) where + refers to the singlet state and − to the triplet state. We further have the periodic boundary conditions with similar conditions on the spatial derivatives. It is convenient to introduce the center-of-mass coordinate R = (x + y)/2 and the relative coordinate r = x − y. In terms of these coordinates the Hamiltonian of Eq. (17) attains the form The eigenstates Φ(R, r) = Ψ(x, y) in the new coordinates then satisfy the equivalent property of Eq. (19) Φ(R, r) = ±Φ(R, −r), and the periodic boundary conditions and similarly for the spatial derivatives. With the Ansatz Φ(R, r) = f (R)g(r) the Schrödinger equation can be separated. The periodic boundary conditions on f and g become g(r + L) = ±g(r), and similarly for the spatial derivatives, where the signs on the right hand side of these equations must be the same for f and g in order to fulfill Eq. (21). The equation for the center-of-mass coordinate R becomes a free particle Schrödinger equation which has the eigenstates (up to normalization) where the boundary conditions with ± in Eq. (23) correspond to k being even and odd respectively. The energy eigenvalue is ǫ = (kπ/L) 2 . After changing coordinates to z = rπ/L the Schrödinger equation in the relative coordinate becomes where we defined M (z) = g(Lz/π). We further defined with E the eigenenergy of the full Hamiltonian of Eq. (17). The boundary condition of Eq. (22) then becomes M (z + π) = ±M (z). Eq. (24) is the wellknown Mathieu equation [12]. The solutions are given by the Mathieu-sine and Mathieu-cosine functions denoted by SE(l, q, z) and CE(l, q, z) where l is a nonnegative integer labelling certain discrete values a l for the constant a in Eq. (24). In the limit λ → 0 (noninteracting case) we simply have CE(l, 0, z) = cos(lz) and SE(l, 0, z) = sin(lz) and a l = l 2 . We thus see that the ± signs in the boundary conditions Eq. (22) correspond to the case that l is even and odd respectively. From Eq. (20) we see that the singlet and triplet case corresponds to the symmetry g(r) = ±g(−r) or equivalently M (z) = ±M (−z) for the Mathieu functions. This means that the singlet solution corresponds to the Mathieucosine function and the triplet to the Mathieu-sine function. The full solution of the problem is therefore given by where + and − refer to the singlet and triplet cases respectively and N ± l is a normalization factor. In both cases k and l need to be both even or both odd. The associated energy eigenvalues are where a ± l (q) are the characteristic values for the Mathieucosine and Mathieu-sine function respectively [12]. For q = 0 the characteristic values obey a + 0 (q) < a − 1 (q) < a + 1 (q) < a − 2 (q) < ... , while in the non-interacting case a + l (0) = a − l (0) = l 2 . We thus nicely see how the twoparticle interaction splits the degeneracy of the spinsinglet and spin-triplet states. In this noninteracting limit the wave functions attain the simple orbital product form where φ n (x) = e inπx/L and where k±l is always even. For any interaction strength the ground state of the QR is the spin-singlet state Ψ + 00 (x, y). We see from Eq. (25) that all states with |k| = l correspond to doubly excited states relative to the ground state which are notoriously difficult to describe by adiabatic functionals. We will return to this issue in Sec. III D. For large values of q the Mathieu functions become localized around z = π/2 (and hence r = L/2) corresponding to the strongly correlated limit of well-localized electrons on opposite parts of the ring. The limit L → ∞ corresponds to q → ∞ and to a limit where the density goes to zero. This limit corresponds to the famous Wigner crystal [13].

B. Exchange-Correlation Potential
We now start to construct the exact xc potential for a specific density n that corresponds to a solution of the time-dependent Schrödinger equation with the Hamiltonian of Eq. (17). For such a density the xc potential is given by Eq. (18). The xc potential can be further split into an exchange (x) and a correlation (c) part v xc = v x + v c where, for our two-electron system, the x potential is simply given by [1] v We choose the density n to come from a freely propagating superposition of two normalized eigenstates of our QR which is a solution to the time-dependent Schrödinger equation. This wave function is properly normalized whenever C 2 0 + C 2 1 = 1. Note that both eigenstates have a constant density. If the constant C 0 is almost 1 (or 0), the density of the system only deviates slightly from being homogeneous. If we look at a small QR, e.g. L = 1 and different interaction strengths λ, we find that even for small deviations from homogeneity the c potential is at least of the same order of magnitude as the x potential. In this case, increasing the density variations by changing C 0 makes the correlation potential v c the dominant contribution to v xc . A notable exception is an initial KS state that has approximately the right initial angular momentum (in the case of λ = 100 and L = 1 this is the state m ′ = 1 as can be seen in Fig. 1). For this case the c potential plus the x potential mainly needs to cancel the Hartree potential. The KS orbital would travel around the ring in approximately the right manner if there were no external perturbations. Besides the initial-state dependence one also clearly sees the non-locality of the c potential in time (memory) and space, as it has in general no obvious simple relation to the local density (see the c potential for m ′ = 0 in Fig. 1). If we go to larger QRs, e.g. L = 2π, the x potential becomes the dominant contribution to the xc potential. This seems counterintuitive since for this case the value of q is larger, corresponding to a more correlated state. It should, however, be remembered that the relation between the density profile (and hence the shape of v s ) and the electronic correlations is rather indirect. For example, the ground state density and KS potential of the QR are spatially constant, independent of the interaction strength. To get more insight into the influence of interactions, it is therefore more useful to study a two-point function. We will therefore now construct the (equilibrium) xc kernel for this problem, which is defined to be the first functional derivative of v xc with respect to the density n, evaluated at the ground state density. We will be able to do so because the ground-state density of the system is homogeneous irrespective of the interaction. Therefore the λ = 0 case is the KS system for any interaction strength λ.

C. Exchange-Correlation Kernel
The xc kernel is the central object of interest in linearresponse TDDFT from which one can determine the perturbative dynamics of the quantum system and its excitation energies. We start by calculating how the ground state spin-density reacts to small external perturbations, i.e.
(see e.g. in Refs. [1,2]), where with ǫ > 0 an infinitesimal,n(xσ) the usual spin-density operator and −∞ ≤ k ≤ ∞ and 0 ≤ l ≤ ∞ (k and l are always either both even or both odd). Here with p = − we refer to the triplet state with spin function (δ σ,↑ δ σ ′ ,↓ + δ σ ′ ,↑ δ σ,↓ )/ √ 2 only, since the spin-triplet functions orthogonal to this one give a zero contribution in the sum. In a first step we can deduce using the periodicity of the solutions that We note that the Mathieu-cosine and Mathieu-sine are real and thus we have D ± (k, l) * = D ± (−k, l). Further we note that D ± (0, l) = 0 for l = 0. After some manipulations of the general expression for the linear-response kernel we end up with where where the sum runs over all even values of l if k is even and over all odd values if k is odd. In the noninteracting case we find due to |D ± (k, l)| 2 → δ |k|,l /(2L 2 ) and ν 0,+ k (ω) = ν 0,− k (ω) the simple expressions and µ 0,− k (ω) = 0. Thus the non-interacting linear response kernel χ 0 has non-zero contributions only from excited states with |k| = l. As discussed below Eq. (25) the states with |k| = l are exactly the singly-excited states of the non-interacting system. We therefore recovered the well-known fact that the non-interacting response function χ 0 has only poles at singly-excited states.
In linear-response (spin) TDDFT the interacting response function χ is expressed in terms of the response function of a non-interacting system with the same density. In our case, since the ground-state density is homogeneous irrespective of the interaction strength λ, the KS system is the one with λ = 0 and the corresponding KS response function is χ 0 . Therefore we can express where the Hartree-exchange-correlation (Hxc) kernel is defined as and integration as well as summation over reoccurring position-spin variables is implied. With the inverse kernels of Eq. (28) we find that The xc kernel is then trivially found by subtracting the interaction potential, i.e. f xc = f Hxc − w. The restriction to k = 0 in the sum is a consequence of the fact that the response functions are only invertible in the space of functions orthogonal to constant function, since a constant potential variation gives no density change. For the Hxc kernel this amounts to the freedom of adding any function of the form g(xσ, x ′ σ ′ , ω) = g 1 (xσ, ω) + g 2 (x ′ σ ′ , ω), since it is always constant either in x ′ σ ′ or xσ when integrating over the internal degrees of freedom in Eq. (29). Therefore adding a function g to f Hxc does not change the linear response kernel χ [14]. We have now fully characterized the behavior of the interacting particles on a QR in terms of the KS system for weak external perturbations. The xc kernel exhibits a strong frequency dependence as it needs to shift the poles of χ 0 and generate new poles in order to have the correct density response of the correlated system. If we Fourier-transformed the kernel from frequency to time, the frequency-dependence would translate to a dependence on previous times, i.e. the frequencydependence corresponds to memory. Therefore we have constructed the first exact density-functional with memory.

D. Frequency-Dependence of the Exchange-Correlation Kernel
In a next step we investigate the frequency dependence of the Hxc kernel in more detail. Such considerations are of importance for developing frequency-dependent approximations to the Hxc kernel [15][16][17][18][19][20], since even advanced approximation schemes can result in unphysical behaviour [21]. To simplify the forthcoming discussion a little we will restrict ourselves to spin-independent linearresponse theory, i.e. we only allow for spin-independent pertubations δv(x, ω) in Eq. (27) and are interested in δn(x, ω) = σ δn(xσ, ω). Therefore we can straightaway sum over all spin-degress of freedom in Eq. (28), leading to Accordingly we no longer couple to the spin-triplet states, and of the whole physical spectrum only the spin-singlet transitions ∆E + kl show up in our linear-response calculations. If we further note that ζ k (x) = x|k is a spatial basis (for square-integrable functions) we can expresŝ It is now interesting to compare the exact expression to some standard approximations for the Hxc kernel. We first note that Therefore the Hartree-exchange approximation (Hx) reads with1 = k |k k| since according to Eq. (26) it is simply obtained by functional differentiation 1/2 of the Hartree term for the case of a two-particle spin-singlet state. The local-density approximation (LDA) together with the Hartree (H) term amounts tô where ǫ ′′ QR is determined by the second functional derivative of the xc energy functional of the (homogeneous) ground-state density [1], i.e., We approximate it from Fig. 2, where we employ a value of ǫ ′′ QR ≃ 0.5 such that we on average reproduce the exact f k Hxc = k|f Hxc |k for |k| > 1. We compare the exact expression for the Hxc kernel to the (frequencyindependent) approximations for k = 1 and 2 in Fig. 2.
The Hx approximation has only a contribution for k = 1 while the HLDA approximation has a contribution for every value of k. We further see in Fig. 2 that for |k| > 1 and ω → ∞ the HLDA approximation and the exact kernel become identical. In order to understand how the frequency-dependence that is missing in the above approximations works, it is useful to express the interacting kernelχ in a different form. From Eq. (29) we find that The task of the denominator1 −χ 0fHxc is two-fold: it shifts the existing poles ofχ 0 and it generates poles that are missing in the bare KS kernel. In order to do so, the denominator has to become zero at the values of the physical resonance frequencies ∆E + kl . This condition reads as Therefore k ′ |χ 0fHxc (∆E ± kl )|k ′′ = δ k ′ k ′′ . If we interpret χ 0fHxc as an infinite-dimensional matrix (in the above basis set), then at the resonance frequencies it has only entries in the diagonal and is zero otherwise. We see that in our case the matrix expression off Hxc is already diagonal for any frequency. This does also not change if we multiply by the matrix expression forχ 0 and find k|.
We immediately see that whenχ has a pole (ν + k (ω) → ∞) then ν 0,+ k (ω)/ν + k (ω) → 0.   This behaviour is nicely visible in Fig. 3, where for k = 1 and k = 2 we have zeros at the first four eigenfrequencies of the interacting system (∆E + 11 = 22.5, ∆E + 20 = 39.5, ∆E + 22 = 79.5 and ∆E + 13 = 99.0). As explained below Eq. (25) the excitations of the form ∆E + (±l)l are singly-excited states, whereas the excitations of the form ∆E + kl with |k| = l correspond to doublyexcited states which do not generate poles in the noninteracting response function. Indeed, we see that the Hxc kernel generates new poles (∆E + 20 and ∆E + 13 ) in the exact response function corresponding to doubly-excited states which are missing entirely in the KS kernelχ 0 . In contrast, the frequency-independent approximations can only shift the already existing poles. The exact kernel does this by canceling the bare KS poles since ν 0,+ k (ω) → ∞, and generates the corresponding physical poles (∆E + 11 and ∆E + 22 ). We note that the Hx approximation does only change the position of the first KS resonance but leaves all others unmodified.

E. Single-Pole Approximation
In practice, even if one had the exact Hxc kernel, calculations are often performed employing certain approximations [1,2]. One of the most important approximations used to determine excitation energies from a linearresponse TDDFT calculation is the so-called single-pole approximation (SPA) [1,2]. This approximation can be derived from our previous considerations on the Hxc kernel, from where we know that if1 −χ 0fHxc = 0 for some frequency ω then it corresponds to a physical resonance. This equation can be rewritten in terms of a generalized eigenvalue problem called the Casida equation [22]. By expanding in terms of KS frequencies ω k [23] one finds that the spin-singlet and spin-triplet excitation energies are perturbatively given by where Φ 0k (x) = ϕ 0 (x)ϕ k (x) and ϕ k (x) the k-th KS orbital as well as In the case at hand we can perform these integrals analytically and find by using that 1/µ 0,+ k (ω k ) = 0 the simple expression In Fig. 3 we have indicated the first two shifted (singletsinglet) eigenvalues calculated with the SPA (Ω + 1 = 25.2 and Ω + 2 = 79.5). While the first shift of the eigenvalues is overestimated (from the bare KS resonance ω 1 = 19.7) the second resonance frequency is extremely well reproduced (with the bare KS resonance being ω 2 = 78.9). Still the calculation of the SPA results includes a sum over (infinitely many) l values. In order to more easily investigate the behaviour of the SPA for a large set of resonances and cases we make a further approximation to find a closed expression for Ω ± k . To do so we note, that the term |D ± (k, l)| 2 gives its main contribution for k = l and falls off rapidly. Therefore it seems a reasonable approximation to employ the value of the non-interacting case, i.e. |D ± (k, l)| 2 → δ |k|,l /(2L 2 ). This leads to the simple explicit expression With this explicit expression [29] we can now easily investigate properties of the SPA. We do so by comparing the relative error of the SPA ∆ ± k = (∆E ± kk − Ω ± k )/∆E ± kk with the relative error of the bare KS eigenvalues The results are displayed in Figs. 4 and 5.  We first note, that the bare KS response does not have a spin-triplet component. Nevertheless, we take the bare KS resonances as a zeroth order guess for both, the spinsinglet and spin-triplet transitions. Therefore in general ω k has a different relative error with respect to ∆E ± kk . This can be nicely seen in Fig. 5 in the difference between δ + k and δ − k for small values of k. For large values of k these differences as well as the relative error of the different approximations becomes small. For the bare KS resonance frequencies this can be explained by the fact, that the higher lying states are dominated by the kinetic energy and are less influenced by the interaction. On the other hand, since the SPA corresponds to the first term in a Laurent expansion with respect to ω k [23], it becomes more accurate when the bare KS excitation energy is closer to the true resonance frequency. Hence the SPA inherits the high k behaviour of the bare KS values. For the lower lying states we find that the SPA usually strongly improves upon the bare KS resonances (see Figs. 4 and 5). In our example the singly-excited states are well-separated from the doubly-excited states. In this case the SPA describes well the singly-excited states.
The adiabatic approximations based on LDA and the exchange-only kernel do not generate new poles and, not surprisingly, fail to describe the doubly-excited states of the system. In order to describe doubly-excited states a frequency-dependent xc kernel is required. We are not aware of any existing simple approximation for a memory kernel that would be able to reproduce the doublyexcited states of our model system. Some approximate kernels exist based on electron gas models, but such kernels would lead to artificial complex excitation energies. In any case, our exact expression for the xc kernel will serve as a useful benchmark to study future density functionals.

IV. OUTLOOK
In this work we have presented analytical expressions of exact density-functionals with initial-state dependence and memory. The functionals were used to give explicit examples of the otherwise very abstract concepts of the xc potential and the KS construction. We demonstrated how one can calculate the exact xc potential for an interacting model system and how one can construct the corresponding exact, frequency-dependent xc kernel. We have then shown how these analytical examples can be used to investigate the basic properties of time-dependent density-functionals and to test approximations.
These results will help to understand the properties of time-dependent density-functionals in more detail. For instance, they have already been used to investigate the Floquet-approach to TDDFT [24]. Further, these exact functionals show how initial-state dependence and memory have to be incorporated into more accurate functional approximations. As such one can employ these exact expressions as benchmarks for the development of new and more reliable functional approximations in TDDFT.