The Dirac–Frenkel Principle for Reduced Density Matrices, and the Bogoliubov–de Gennes Equations

. The derivation of eﬀective evolution equations is central to the study of non-stationary quantum many-body systems, and widely used in contexts such as superconductivity, nuclear physics, Bose–Einstein condensation and quantum chemistry. We reformulate the Dirac–Frenkel approximation principle in terms of reduced density matrices and apply it to fermionic and bosonic many-body systems. We obtain the Bogoliubov– de Gennes and Hartree–Fock–Bogoliubov equations, respectively. While we do not prove quantitative error estimates, our formulation does show that the approximation is optimal within the class of quasifree states. Furthermore, we prove well-posedness of the Bogoliubov–de Gennes equations in energy space and discuss conserved quantities.


Results on
The time evolution of the state ψ t of a system of N quantum particles is described by the time-dependent Schrödinger equation where ψ t ∈ L 2 a (R 3N ) for spinless fermions (where the wave function is totally antisymmetric under any exchange of particles) and ψ t ∈ L 2 s (R 3N ) for bosons (where the wave function is totally symmetric under any exchange of particles).In this generality, the Schrödinger equation models a vast range of physical systems, starting from nucleons in the atomic nucleus over electrons in semiconductors to stars for the fermionic theory, or Bose-Einstein condensates for the bosonic theory, depending on the choice of the one-particle Hamiltonian h (we think of h = −∆ + V ext (x) with some external potential V ext : R 3 → R; h i denotes this operator as acting in the variable x i ∈ R 3 ) and the pair interaction V : R 3 → R. Unfortunately, these systems also have an enormous number of degrees of freedom, making analytical and numerical solutions generally impossible.For this reason, there is a lot of interest in approximate theories (also called effective evolution equations), which contain fewer degrees of freedom and make analytical and numerical treatments possible.Of course, such theories do not achieve the broad validity of the Schrödinger equation and provide a good approximation only in specific physical regimes, which are mathematically modeled as scaling limits.In this paper, we discuss a geometric method for the derivation of effective evolution equations.This method, even though not the most convenient for proving quantitative error estimates for the obtained approximation, directly shows that the obtained equations are optimal as far as the available degrees of freedom permit.The method is also independent of any choice of scaling limit.
For this introduction we focus on (for simplicity of notation spinless) fermionic systems.The corresponding bosonic notions will be introduced in Sect. 4.Here we deal only with pure quasifree states.Only in Sect. 5 for the topic of well-posedness we also consider mixed states; we will always highlight explicitly when we talk about mixed states.
The most basic approximate theory of fermionic systems is obtained by restricting the Schrödinger equation to wave functions that are Slater determinants, ψ(x 1 , . . .x N ) = (N !) −1/2 det f i (x j ) N i,j=1 , also denoted as the antisymmetrized tensor product ψ = (N !) −1/2 f 1 ∧• • •∧f N , where the one-particle wave functions f j constitute an orthonormal set in L 2 (R 3 ).The corresponding effective evolution equation for the Slater determinant is given by the Hartree-Fock system of N nonlinear coupled PDEs for the one-particle wave functions: V * f i,t f j,t f j,t . (1.2) More conveniently, the Hartree-Fock equations can be formulated in terms of the oneparticle reduced density matrix.
The one-particle reduced density matrix of a state ψ ∈ L 2 (R 3N ) is defined as the non-negative trace-class operator γ on L 2 (R 3 ) obtained by taking the partial trace over N − 1 particles of the many-body density matrix |ψ ψ|, γ = N tr 2,...N |ψ ψ|, where we have chosen to normalize the one-particle reduced density matrix such that tr γ = N .If ψ = (N !) −1/2 f 1 ∧ . . .∧ f N , we find the rank-N projection γ = N j=1 |f j f j |.Conversely, every rank-N orthogonal projection specifies (uniquely up to a phase) a Slater determinant (just take the spectral decomposition of the projection to find the one-particle wave functions f j ).
If the one-particle wave functions have time dependence given by the Hartree-Fock equations (1.2), then γ t = N j=1 |f j,t f j,t | satisfies the equivalent equation where V * ρ γt (x) = dy V (x − y)γ t (y, y) is a multiplication operator called the direct term, and X V (γ t )(x, y) = V (x − y)γ t (x, y) is the integral kernel of an operator called the exchange term.
They satisfy the canonical anti-commutation relations (CAR), i. e.
{a(f ), a(g)} = 0, {a * (f ), a * (g)} = 0, and {a(f ), a * (g)} = f, g for all f, g ∈ L 2 (R 3 ).(The definition of the anti-commutator is {A, B} = AB + BA.) The vector Ω = (1, 0, 0, . ..) ∈ F a is called the vacuum state and is in the kernel of all annihilation operators, a(f )Ω = 0 for all f ∈ L 2 (R 3 ); it describes a system not containing any particles.It is convenient to introduce the operator-valued distributions a * x and a x with the defining property that (in a weak sense, within expectation values) They satisfy the formal canonical anti-commutation relations {a x , a y } = 0, {a * x , a * y } = 0, and {a x , a * y } = δ(x − y).The Hamiltonian is generalized to Fock space as where H n denotes the first quantized Hamiltonian as given in (1.1).The Hamiltonian H can also be represented in terms of creation and annihilation operators by 1 Restricted to N -particle states, this Hamiltonian agrees with H N .For the geometric considerations in this paper, the explicit form of H does not need to be specified, as long as it is a self-adjoint operator and conserves the number of particles.(A Hamiltonian conserves the number of particles if it commutes with the particle number operator N defined by In particular all operators of the form (1.4) conserve the particle number.)Only Sect. 5 refers to the particular Hamiltonian (1.5).
Using the creation and annihilation operators, the definition of the one-particle reduced density matrix γ is extended to states ψ ∈ F a by defining it to have the integral kernel γ(x, y) := ψ, a * y a x ψ . (1.6) We can now give a simple proof that γ ≤ 1: For all f ∈ L 2 (R 3 ), using the CAR, Slater determinants are a special case of a class of more general states in Fock space called quasifree states (see e. g. [57] for a very readable introduction).The defining property of quasifree states is that they are exactly those states ψ ∈ F a for which the Wick theorem holds, i. e. expectation values of creation and annihilation operators can be reduced to the sum of the expectation values of all possible pairings of just two operators (with the sign being the sign of the corresponding pairing); for example (Here we used the notation a ♮ j to denote an operator without specifying whether it is a creation or annihilation operator.)Notice that the Wick theorem allows us to express any expectation value of creation and annihilation operators in a quasifree state purely in terms of the one-particle reduced density matrix γ and the pairing density α(x, y) := ψ, a y a x ψ . ( Slater determinants are exactly those quasifree states for which the pairing density identically vanishes, α = 0. Notice that for a general quasifree state γ is not a rank-N projection; instead a quasifree state always satisfies γ 2 − γ = αα and αγ = γα. (1.8) Given γ and α satisfying (1.8), there is a (up to a phase unique) quasifree state ψ ∈ F a such that (1.6) and (1.7) hold.
There is a more compact way of writing the equations (1.8) by introducing the generalized one-particle reduced density matrix Γ.The generalized one-particle reduced density matrix Γ is an operator on (1.9) The characterization (1.8) of quasifree states is equivalent to Γ being an orthogonal projection, Γ 2 = Γ = Γ * .Any (not necessarily quasifree) generalized one-particle reduced density matrix has the property 0 ≤ Γ ≤ 1, and thus Γ 2 ≤ Γ. (1.10) In the theory of superconductivity, the pairing density is interpreted as the wave function of electrons that have formed Cooper pairs, which in many ways behave like bosons.These Cooper pairs are seen as the carriers of the superconducting current that has attracted so much attention for its technological applicability in the dissipationless transport of electricity.
In this paper, our focus lies on the effective evolution equation obtained by restriction of the many-body evolution to quasifree states with pairing.This system of effective evolution equations is known as the Bogoliubov-de Gennes equations with h HF (γ t ) as defined in (1.3) and the operator Π V (α t ) defined through its integral kernel Π as in [37,36] we use the generalized Hartree-Fock operator ).The Bogoliubov-de Gennes equations for fermionic systems are sometimes also called the generalized Hartree-Fock equations or fermionic Hartree-Fock-Bogoliubov equations; the usual Hartree-Fock equations correspond to the Bogoliubovde Gennes equations with α = 0.By restricting to physical regimes where direct and exchange term are negligible, and by including electron spin, one obtains the timedependent BCS equations (named after Bardeen, Cooper, and Schrieffer), which describe the dynamics of electrons and Cooper pairs in superconductors.
In the present paper, our goal is to formulate a systematic approximation principle by which we can obtain the Bogoliubov-de Gennes equations from many-body quantum theory.The approximation principle we establish is a reformulation of the Dirac-Frenkel principle in the space of reduced density matrices, and it yields the equations sometimes called the quasifree reduction principle.Applying the quasifree reduction principle to the Hamiltonian (1.5), one obtains the Bogoliubov-de Gennes equations.Afterward we study well-posedness and conserved quantities for the Bogoliubov-de Gennes equations.
While it is general knowledge that the quasifree reduction principle should be a consequence of the Dirac-Frenkel principle, we are not aware of a direct proof having appeared before; in particular the formulation of the Dirac-Frenkel principle in terms of reduced density matrices has not been given before.Among the advantages of our approach is that it shows that the obtained effective equations describe the optimal evolution possible within the approximation manifold.
Earlier results.The derivation of effective evolution equations for many-body systems has attracted a lot of attention in the community of mathematical physics and can be seen as a cornerstone of non-equilibrium statistical mechanics.Consequently, the literature is vast and we cannot claim to provide a complete overview.Let us say so much, that the geometric approximation principle on which we build in this paper goes back to the founding fathers of quantum mechanics [18,27].A rigorous mathematical discussion and highly valuable presentation has been given in [42].
The next step after the geometric derivation of the correct effective equation lies in the proof of convergence toward the effective equation and then the derivation of quantitative error bounds in given physical regimes modeled as scaling limits.This topic has attracted a lot of attention in recent years.For bosonic systems, many such results on the approximation of reduced density matrices have been proven for the meanfield model [6,26,53,16,48,40,56,1,39,29,2] and the Gross-Pitaevskii model [7,21,23,22,24,25,49,47].For fermionic systems, most results have only appeared in the last few years.The main regimes treated here are the mean-field regime on short time scales [5,28], the mean-field regime with slow variation of the effective interaction potential [4,46,45], and the combined mean-field/semi-classical limit for high-density systems [44,58,20,10,9,8,50].The papers cited in this paragraph do not use the Dirac-Frenkel principle but other methods that have been specifically developed for many-body systems, like the BBGKY hierarchy, coherent states, a Schwinger-Dyson expansion or counting the number of particles well-described by the effective evolution equation.Some applications of the Dirac-Frenkel principle with explicit error estimates can be found in [42]; another example is [30].
In the context of proving convergence toward effective equations for bosonic systems in mean-field and Gross-Pitaevskii scaling limits, equations of the Hartree-Fock-Bogoliubov-type appear when considering second-order corrections beyond the above results on approximation of reduced density matrices, that is approximations in the norm of the many-body Hilbert space [32,31,34,33,43,11].The difference to our discussion here is that we are just interested in the approximation of reduced density matrices, not in the norm of the many-body Hilbert space, but instead we focus on ensuring optimality of the derived effective equations.Moreover, while the scaling regimes are crucial for obtaining convergence and quantitative error estimates, in the qualitative geometric approach that we use, it is not necessary to specify a particular scaling limit.
For fermionic systems, the derivation of quantitative error bounds in appropriate scaling limits for the Bogoliubov-de Gennes equations or even the BCS equations remains an open problem.
Organization of the paper.In Sect.2.1 we shall recall the variational principle of Dirac and Frenkel, which is at the base of our paper.In Sect.2.2 we recall the principle of quasifree reduction, which is less fundamental and less general than the Dirac-Frenkel principle, but more convenient for explicit calculations.In Sect.3.1 we give a re-derivation of the Hartree-Fock equation as the simplest example to introduce our formulation of the Dirac-Frenkel principle in the space of one-particle reduced density matrices.In Sect.3.2 we use our formulation of the Dirac-Frenkel principle for systems with pairing, yielding the time-dependent Bogoliubov-de Gennes equations.In Sect.4, we present the analogous formulation of the Dirac-Frenkel principle for bosonic systems.Finally, in Sect. 5 we discuss well-posedness of the time-dependent fermionic Bogoliubov-de Gennes equations.

Approximation Principles
In this section we first introduce the Dirac-Frenkel principle, which can be seen as the fundamental principle for deriving optimal effective evolution equations.Afterwards we introduce the principle of quasifree reduction which does not exhibit the optimality but leads to the same results as the Dirac-Frenkel principle and is calculationally simpler.

The Dirac-Frenkel Variational Principle
In this section we discuss the abstract formulation of the Dirac-Frenkel variational principle for the approximation of the time-dependent Schrödinger equation by projection onto a submanifold.The Dirac-Frenkel variational principle is particularly interesting for its clear geometrical content, which shows that the obtained equation on the submanifold is the optimal choice.We follow [42,Chapter II].
We consider the Schrödinger equation as an evolution equation in a complex Hilbert space H.We use the convention that the scalar product is anti-linear in the first and linear in the second argument.The Hamiltonian H is a self-adjoint operator on H.The Schrödinger equation reads Now consider a smooth (typically infinite dimensional) submanifold M of H.The tangent space of M in the point u ∈ M is denoted by T u M; it consists of the derivatives u ′ 0 of all differentiable paths t → u t passing through u 0 = u.We are interested in approximating the solution ψ t of the Schrödinger equation with a path u t on the manifold M, assuming that initially ψ 0 = u 0 ∈ M. As pointed out already by Dirac, the path t → u t is to be chosen such that at every time t the derivative u ′ t ∈ T ut M is as close as possible to 1 i Hu t ; in other words, the path is determined by choosing its derivative as the orthogonal projection of 1  i Hu t onto the tangent space: P (u t ) being the orthogonal projection from H onto T ut M. From this formulation it is clear that the Dirac-Frenkel principle yields the effective evolution equation which at every infinitesimal time step is optimal.
In the case of fermionic many-body systems, one typically chooses H = L 2 a (R 3N ) (for the case of no pairing) and M as the set of N -particle Slater determinants.While this approach does yield the time-dependent Hartree-Fock equations (1.2), it is not expected to ever do so with controllable errors: it is a general fact that in many-body systems, the norm of many-body wave functions as a measure of distance has unfortunate behavior as the number of particles grows.In fact, the quantitative derivation of effective equations for many-body systems in appropriate scaling limits is typically proven in terms of the i Hu of the exact evolution is orthogonally projected into the tangent space T u M, yielding the optimal effective evolution in the approximation manifold M. Figure following [42].
trace norm or Hilbert-Schmidt norm of reduced density matrices (see the overview of results in Sect.1).For this reason, and to make the connection to the principle of quasifree reduction, in this paper we formulate the Dirac-Frenkel principle in the space of one-particle density matrices.This formulation will be calculationally very convenient when we consider quasifree states with pairing.

The Principle of Quasifree Reduction
The principle of quasifree reduction appears to be the computationally most accessible principle for deriving effective equations for many-body quantum systems.It applies to the particular case were the approximation manifold is given by a class of quasifree states.Typically it is formulated directly in the language of reduced density matrices.
The principle of quasifree reduction asserts that for fermionic systems the effective evolution equations are where ψ qf t ∈ F a is the quasifree state uniquely (up to a phase) assigned to γ t and α t .For bosonic systems the equations proposed as the quasifree reduction principle, including a condensate ϕ t ∈ L 2 (R 3 ), are the following [3] where ψ bog t in the bosonic Fock space is the Bogoliubov state associated with (ϕ t , γ t , α t ) (see (4.38) for the definition).
While this principle is easy to formulate, calculationally accessible and has been frequently used in many contexts (e. g. very recently to derive the Hartree-Fock-Bogoliubov equations for bosonic systems [3] and the Bogoliubov-de Gennes equations for fermionic systems [17]), it is not completely obvious that it is a consequence of the more fundamental Dirac-Frenkel principle.Maybe more severely, it is not at all clear whether the quasifree reduction principle yields the optimal approximation possible within the manifold of quasifree states.
In the present paper, we prove that the principle of quasifree reduction does follow directly from the Dirac-Frenkel principle, in particular showing that it yields the optimal effective evolution equations.

Derivation of the Quasifree Reduction Principle for Fermions
In this section we derive the principle of quasifree reduction from a re-formulation of the Dirac-Frenkel principle.We first sketch the instructive case of systems without pairing before generalizing to fermionic systems that also exhibit pairing.

Fermionic Systems without Pairing
Here we shall warm up with the case of no pairing (α = 0), i. e. giving a derivation of the standard Hartree-Fock equation.To this end we shall formulate the Dirac-Frenkel principle in terms of the one-particle reduced density matrix and from there, derive the principle of quasifree reduction.To our knowledge, this is the first formulation of the Dirac-Frenkel principle in terms of reduced density matrices.Notice that the derivation of the quasifree reduction principle does not make use of any particular form of the Hamiltonian; we assume only that it commutes with the particle number operator (i.e. the number of particles is conserved along the many-body evolution).
We consider the Schrödinger equation in the space L 2 a (R 3N ), i. e. describing a fermionic system of N -particles.The many-body evolution t → ψ t induces an evolution of the associated one-particle reduced density matrix γ ψt , which satisfies The one-particle reduced density matrix is a non-negative operator on L 2 (R 3 ).Since our system has a finite number of particles, tr γ ψt = N , the one-particle reduced density matrix is a Hilbert-Schmidt operator.We thus choose the ambient Hilbert space in which γ ψt lives to be Due to the condition of self-adjointness this is only a real-linear (instead of complexlinear) space; in the following all spaces are real-linear only.Corresponding to Slater determinants in the wave function picture, we choose our approximation manifold to be given by orthogonal projections, This is an infinite dimensional Hilbert submanifold of the self-adjoint Hilbert-Schmidt operators.The effective evolution equation is to be found in M.This is achieved in the optimal way by applying the Dirac-Frenkel principle reformulated in terms of the one-particle reduced density matrix.
Dirac-Frenkel Principle for Reduced Density Matrices.The effective evolution equation within the submanifold M is given by where ψ qf t is the state uniquely (up to the phase) associated with γ t , and proj(γ t ) : H → T γt M is the projection onto the tangent space in the point γ t .
We start by determining the tangent space and the projection onto the tangent space.
The orthogonal projection from H onto T γ M is given by Proof.Let A ∈ T γ M. By definition there exists a differentiable curve t → γ t in M such that γ 0 = γ and γ ′ 0 = A. (By definition of differentiability in the norm of the ambient Hilbert space H, A is a Hilbert-Schmidt operator.)Taking the derivative of the projection condition γ 2 t = γ t , we find Aγ + γA = A. Multiplying from the left and right by γ, we get 2γAγ = γAγ, so γAγ = 0. Furthermore, by multiplying it from the left and the right by and Hilbert-Schmidt, so e tB is a unitary.Now let γ t := e tB γe −tB .This is a curve of orthogonal projections with γ 0 = γ.Its derivative at zero is It is easy to see that proj(γ) is a projection and has the claimed image.
The manifold of orthogonal projections M has several connected components, corresponding to the value tr γ ∈ N. Clearly, any differentiable curve always stays within the same connected component, so we do not have to worry about this.Theorem 3.2 (Derivation of the Quasifree Reduction Principle, No Pairing).The Dirac-Frenkel principle (3.19), with H and M as chosen above, is equivalent to the quasifree reduction principle without pairing (α = 0): where ψ qf t is the Slater determinant uniquely (up to the phase) associated with γ t .(Equation (3.20) is the same as equation (2.16) but written without the use of operator-valued distributions.) We refer the reader to the more general proof of Theorem 3.4.
One may convince oneself that (3.20) indeed yields the Hartree-Fock equation (1.3) when evaluating the expectation value on the r.h.s. with the many-body Hamiltonian from (1.5), using the canonical anti-commutation relations and the Wick theorem.
Having (1.3) at hand, we can also obtain the quasifree reduction principle from the Dirac-Frenkel principle for reduced densities as follows.Clearly (1.3) implies that γ 2 t also satisfies the Hartree-Fock equation, So if we have a projection as initial data, γ 2 0 = γ 0 , assuming uniqueness, we conclude that γ 2 t = γ t for all times t.Alternatively, we could argue that the Hartree-Fock equation preserves the spectrum of γ t , which also implies γ 2 t = γ t for all times.Either way, we conclude that the derivative is in the tangent space of M, which makes the projection in the Dirac-Frenkel principle trivial and yields the quasifree reduction principle.However: • This argument uses (1.3) which is obtained by explicitly evaluating the quasifree reduction principle.Using only the equations of the quasifree reduction principle (2.16), there is no easy way to formulate (3.21); in fact, a direct verification that (2.16) stays within M seems complicated to us.
• The argument depends on the choice of the one-particle Hamiltonian h and regularity and decay of the interaction potential V and the initial data.For the initial value problem with pairing (α 0 = 0), uniqueness or conservation of the spectrum are by themselves non-trivial problems, see Sect. 5.
Our derivation does not require any specification of the Hamiltonian beyond its existence as a self-adjoint, particle number conserving, operator.Furthermore, our geometric approach makes it clear that the quasifree reduction principle is the optimal approximation within the set of quasifree states.(Also or ∂ t γ t = 0 would be an evolution in the manifold of quasifree states-but far from being the optimal approximation to the many-body problem.)

Fermionic Systems with Pairing
We now extend our formulation of the Dirac-Frenkel principle to derive the approximation of fermionic many-body systems by quasifree states with pairing, namely the Bogoliubov-de Gennes equations.As before, the derivation of the quasifree reduction principle does not require any particular form of the Hamiltonian; our only assumption is that it conserves the number of particles.
The geometry becomes very similar to the case of no pairing by using the generalized one-particle reduced density matrix Γ.Thus the manifold of quasifree states with pairing can be described by Γ 2 = Γ and the block structure (1.9), which however comes 'for free' since it is present in all generalized one-particle reduced density matrices (in particular also in the one derived from the many-body Schrödinger equation).
Let us be a bit more precise and define the involved spaces.First of all notice that, due to the form (1.9), we have the generalized one-particle density matrix is not a Hilbert-Schmidt operator.We remedy this problem by considering the generalized one-particle reduced density matrix as a point in an affine space, and the approximation manifold as a submanifold of this affine space.Let us denote by Γ vac := 0 0 0 1 the generalized one-particle reduced density matrix of the vacuum Ω ∈ F. Then any generalized one-particle reduced density matrix can be written as Every generalized one-particle reduced density matrix satisfies Γ 2 ≤ Γ, which implies γ 2 − αα ≤ γ (only for quasifree states we had equality here); thus which is twice the expected number of particles and as such assumed to be finite.The expected number of particles is trivially conserved along the many-body evolution since we assume the Hamiltonian to commute with the particle number operator; it is typically also conserved along the effective evolution, c. f.Lemma 5.4, so it is justified to take Γ as a Hilbert-Schmidt operator.Let us therefore introduce the affine space Similar to the no-pairing case, A is a real-linear space.
Now notice that the requirement of having the block structure of Γ in terms of γ and α as in (1.9) can be rewritten 2 as the condition ) being the anti-linear operator of complex conjugation.So we can think of the evolution of the many-body generalized one-particle reduced density matrix as living in the affine subspace of A given by (But not every Γ ∈ A − is the generalized one-particle reduced density matrix of a Fock space vector.) The approximation manifold is again given by the generalized one-particle density matrices corresponding to quasifree states: So compared to the no-pairing case not much has changed-the only additional complication is that we have to impose the block structure of Γ in terms of γ and α.Luckily, this block structure is present in any generalized one-particle reduced density matrix including the one of the many-body evolution.So the many-body evolution describes a curve in the affine subspace A − , of which M is a submanifold.
To provide a characterization of the tangent space, we also introduce as an auxiliary space the manifold of projections which do not necessarily have the block structure Notice that, since A is an affine space, T Γ A = A for any Γ ∈ A. (i) onto the tangent space of projection operators (ii) onto the tangent space of the affine subspace with the block structure

.26)
2 There is a subtlety here: The latter one however is not a Hilbert-Schmidt perturbation of Γvac and thus not a solution within A; in fact it corresponds to a state formally obtained from Γ by a particle-hole transformation replacing the vacuum by an infinite number of fermions filling up all the Hilbert space.
(iii) and onto the tangent space of quasifree states proj(Γ) : A → T Γ M.
Then, for Γ ∈ M, we have Proof.The projection onto the tangent space of projection operators (3.25) is known from Lemma 3.1.Since A − is an affine subspace, we can simply take the derivative of the defining equation to find It is easy to check that the formula (3.26) defines an orthogonal projection, maps into T Γ A − and is surjective onto ).Thus, using Lemma 3.1, Let P := proj − (Γ) proj aux (Γ).Obviously P is an orthogonal projection.It is easy to verify that it maps into T Γ M and is surjective onto T Γ M; therefore proj(Γ) = P .Now let A ∈ T Γ A − .Then proj aux (Γ)A = proj aux (Γ) proj − (Γ)A = proj(Γ)A, so (3.27) holds.
We can now derive the quasifree reduction principle from the Dirac-Frenkel principle.
Theorem 3.4 (Derivation of the Quasifree Reduction Principle, With Pairing).The effective equation for the generalized one-particle reduced density matrix Γ t obtained by applying the Dirac-Frenkel principle to the many-body evolution with M and A as chosen above yields the principle of quasifree reduction where ψ qf t is the quasifree state uniquely (up to its phase) assigned to Γ t .(Equation (3.28) is a compact way of writing (2.16), avoiding the use of operator-valued distributions by testing against F 1 and F 2 .) Proof.The proof uses some theory of Bogoliubov transformations, for which we recommend [57, Chapters 9 and 10] as a reference.(For the no-pairing case, the Bogoliubov transformation is a simple particle-hole transformations, see, e. g., [10,9].) Recall that the generalized one-particle reduced density matrix Γ ψ of a Fock space vector ψ is, avoiding the use of operator-valued distributions by testing against where the generalized creation and annihilation operators are and So for ψ t being the solution of the many-body Schrödinger equation, the associated generalized one-particle reduced density matrix satisfies Notice that the r.h.s., like the derivative of any differentiable curve of generalized oneparticle reduced density matrices, lies in T Γ A − .According to the Dirac-Frenkel principle, we have to project it onto the tangent space of quasifree states.We apply the projection as given by (3.27) and (3.25) to get where ψ qf t is the quasifree state uniquely assigned to Γ t .Comparing to the quasifree reduction principle (3.28), we see that we simply have to show that i H]ψ qf t = 0 and then also that Since ψ qf t is a quasifree state, it can be written in terms of an implementable Bogoliubov map Vt Ω (U Vt being the unitary implementation in Fock space).Take any ) of the Bogoliubov map and recalling (3.29), we calculate so we obtain Using this last identity we calculate that where we have introduced the notation V t −1 F i =: Fi =: fi gi for i ∈ {1, 2}.Now (Here we made use of the fact that g2 , g1 as a complex number commutes with everything, and of the fact that any annihilation operator applied to the vacuum gives zero.)Similarly, we find for the other diagonal block as well that it vanishes, Using the Wick theorem and the CAR, it is a simple calculation that the quasifree reduction principle (3.28), applied to the Hamiltonian (1.5), yields the time-dependent Bogoliubov-de Gennes equations (1.11).
Remark.The reader may wonder how it is possible that the many-body evolution gives rise to the equation and the effective evolution solves the seemingly identical equation yet the two evolutions in general differ even if they both start from quasifree initial data.
The answer is that (3.30) is not a well-posed initial value problem, simply because a general Fock space state has many more degrees of freedom than just the generalized oneparticle reduced density matrix; the r.h.s. is not a function of only Γ t .The equation (3.30)only makes sense if the r.h.s. is already prescribed by the Schrödinger equation (1.1).
On the other hand, (3.31) is a well-defined initial value problem because quasifree states in Fock space are (up to a phase) one-to-one with their generalized one-particle reduced density matrix.So the r.h.s. is a function only of Γ t here (alternatively think of the Wick rule which also shows that the r.h.s. can be expressed in terms of only Γ t ).
We provide the rigorous proof of well-posedness for a main class of physically relevant Hamiltonians and initial data in Sect. 5.

Derivation of the Quasifree Reduction Principle for Bosons
In this section we present the formulation of the Dirac-Frenkel principle for one-particle reduced density matrices of bosonic systems.This is slightly more complicated than for fermionic systems because the simple projection condition has to be replaced, and because we include a condensate, but can be treated by modifications of the previously developed geometric notions.
We start by reviewing some definitions for bosonic systems where they differ from the corresponding fermionic formulas.For a comprehensive introduction we refer to [57].Bosonic Fock space is defined in the same way as for fermionic systems, simply replacing antisymmetric by symmetric wave functions: Creation and annihilation operators a * (f ) and a(f ) (where f ∈ L 2 (R 3 ), a one-particle wave function) are defined as The bosonic creation and annihilation operators satisfy the canonical commutation relations (CCR), i. e.
[a(f ), a(g)] = 0, [a * (f ), a * (g)] = 0, and [a(f ), a * (g)] = f, g The corresponding operator-valued distributions satisfy the formal canonical commutation relations [a x , a y ] = 0, [a * x , a * y ] = 0, and [a x , a * y ] = δ(x − y).Quasifree states are defined as those states for which the Wick theorem holds, which only differs from the fermionic case by having all positive signs, e. g., ψ .The r. h.s. can be expressed in terms of γ(x, y) = ψ, a * y a x ψ and α(x, y) = ψ, a y a x ψ .
For any bosonic quasifree state, γ and α are related by conversely all γ and α satisfying these two equations define a (up to a phase) unique quasifree state in bosonic Fock space.The generalized one-particle reduced density matrix is defined as The relations (4.32) characterizing it as belonging to a quasifree state can be rewritten The generalized one-particle reduced density matrix Γ is a non-negative operator on As for fermionic systems, also bosonic quasifree pure states can be written in terms of a Bogoliubov transformation [57]: If ψ ∈ F s is quasifree, then there exists an implementable Bogoliubov map V : ), where the generalized creation/annihilation operators A(F ), A * (F ), F ∈ L 2 (R 3 ) ⊕ L 2 (R 3 ) are defined exactly the same way as for fermions.
With Γ vac = 0 0 0 1 the generalized one-particle reduced density matrix of the vacuum (identical to the fermionic case), we define the spaces (4.37) These take the role of: A the ambient affine space defining the scalar product, A + the affine subspace in which the many-body evolution can be found, and M the approximation manifold of generalized one-particle reduced density matrices of quasifree states.To see that M is indeed a submanifold of A + , notice that by (4.32) (or by (2Γ + S)S(2Γ + S) = S, which is equivalent to (4.34)), we can write every Γ ∈ M as and thus M as a graph.Alas!The computation of the tangent spaces of M from its graph representation involves the derivative of the operator square root around 1, which leads to Lyapunov equations of type {X, A} = B with A = √ 1 + αα and B given.There is no simple closed formula for the solution to this equation in operator form (one can only express X as a function of the eigenvalues and the eigenfunctions of A).We overcome this problem by noticing that it is sufficient to have a parametrization of the orthogonal complement of the tangent space.
Proof.Let us define the map ψ : Γ → P = −ΓS.We can explicitly write down its inverse: Γ = −P S since S 2 = 1.Let us specify domains and codomains.In parallel to the spaces A, A + and M we introduce (notice that P vac := ψ(Γ vac ) = Γ vac ) It is easy to check that ψ is an isometric isomorphism A → G (both sides with the Hilbert-Schmidt scalar product) and is also an isometric isomorphism Following the strategy of Lemma 3.1 and Lemma 3.3, with the condition SP * S = P taking the place of self-adjointness everywhere, we obtain The differential of ψ is given by D Γ ψB = −BS, which is also an isometric isomorphism T Γ A + → T P G + .In particular it conserves orthogonality, so it is also an isomorphism where the orthogonal complement is defined by the decomposition Noticing that T P G + = {B ∈ G : B + J BJ = 0}, the proof is complete.
Unlike fermionic states, bosonic states can exhibit condensation, so that for ψ ∈ F s it is possible that for some f ∈ L 2 (R3 ) we have the additional degree of freedom ψ, a(f )ψ = 0.
(For any quasifree state this is vanishing.)Let us define the Weyl operator The Weyl operator is unitary and W(ϕ) * = W(−ϕ); furthermore they also satisfy ).An ideal condensate is described by a coherent state Ψ = W(ϕ)Ω; we have and consequently Ψ, a(f )Ψ = f, ϕ .The expected number of particles in the coherent state is Ψ, N Ψ = ϕ 2 L 2 .Using the BCH formula 3 together with the CCR we find From the last formula we see that a coherent state is a linear combination of different particle numbers, where the probability to measure n particles is given by a Poisson distribution peaked at the value ϕ 2 L 2 .We now enlarge the class of quasifree states to the class of Bogoliubov states by including a condensate; more precisely, a Bogoliubov state4 is any state of the form where ϕ ∈ L 2 (R 3 ) (typically not normalized) and V is any implementable Bogoliubov map.Using the fact that any expectation value of an odd number of creation and annihilation operators in a quasifree state vanishes, Furthermore we find that the one-particle reduced density matrix is given by (4.40) Similarly, we find the pairing density to be α(x, y) = U V Ω, a y a x U V Ω + ϕ(x)ϕ(y) =: α(x, y) + ϕ(x)ϕ(y).So by first obtaining ϕ through (4.39) and then solving (4.40) and (4.41) for α, γ, we have a natural way of assigning a unique (ϕ, γ, α) to every quasifree state; conversely every triple (ϕ, γ, α) satisfying (4.42) defines a (up to a phase) unique Bogoliubov state in Fock space through (4.38).So as we just argued, Bogoliubov states are characterized by independently the condensate wave function ϕ ∈ L 2 (R 3 ) and the truncated expectations, i. e. Γ.We therefore introduce the manifold where M is the manifold of quasifree generalized one-particle reduced density matrices as determined before.Of course, the tangent space is given by then apply the projection onto the tangent space to ∂ t Γqf t .The projected derivatives describe the effective evolution.
Theorem 4.2 (The Quasifree Reduction Principle for Bogoliubov States).The Dirac-Frenkel principle, applied by projecting the curve of the many-body evolution from the space L 2 (R 3 )×A to the approximation manifold of Bogoliubov states M bog = L 2 (R 3 )×M yields the equations of the quasifree reduction principle (2.17).
As it was already the case for fermionic systems, the only difference between the Dirac-Frenkel principle and the principle of quasifree reduction is the projection onto the tangent space.So instead of really doing the projection onto the tangent space, we simply check that the right-hand sides of (2.17) are orthogonal to the orthogonal complement of the tangent space.
Proof.Recall (4.43): as far as ϕ is concerned, the tangent space is given by all of L 2 (R 3 ); i. e. the projection onto the tangent space is just the identity.Therefore we only have to take care of projecting the evolution of γ and α; more precisely we will check that the derivative of Γqf t already lives in the tangent space.For this, it is sufficient to show that A, ∂ t Γqf t S 2 = 0 for all operators A ∈ (T Γ M) ⊥ .Notice that a priori ∂ t Γqf t ∈ T Γ A + , so it is sufficient to consider the orthogonal complement as a subspace of T Γ A + instead of all of A. So by Lemma 4.1, we can write A = − (P * BP * + (1 − P * )B(1 − P * )) S for some operator B.
Since B is Hilbert-Schmidt, it has a singular value decomposition B = j λ j |ξ j ϕ j |, So it suffices that every such expectation value vanishes individually.Recall that we have ). Recall also that we can write ψ qf t = U V Ω, where the Bogoliubov map V satisfies Using these facts we find We now check that the second kind of expectation value vanishes.For all A similar calculation shows that also M, thus the projection on the tangent space is trivial, and the Dirac-Frenkel principle becomes the quasifree reduction principle.
The calculation to obtain the explicit evolution equations from the bosonic principle of quasifree reduction (2.17) was sketched in [3], based on the canonical commutation relations and the Wick theorem.The explicit equations are known as the Hartree-Fock-Bogoliubov equations, here written in terms of the truncated expectations αt and γt and the condensate wave function ϕ t , where h HFB differs only by the sign of the exchange term from the fermionic h HF : More compactly, (4.45) can be written symplectically [3, Eq. ( 41)] where Γt is the truncated generalized one-particle density matrix (4.42), and Γ t the (non-truncated) generalized one-particle density matrix.The generalized Hartree-Fock-Bogoliubov operator is 5 Well-Posedness of the Fermionic Bogoliubov-de Gennes Equations The well-posedness of the effective equation obtained from the Dirac-Frenkel principle is not automatic; however under reasonable assumptions on the interaction potential and the initial data it can be established by standard methods.Since to our knowledge there is no proof completely spelled out in the literature, in this section we give a detailed proof that the time-dependent fermionic Bogoliubov-de Gennes equations are well-posed.We consider only the case h = −∆ (in particular no external potential V ext is included), and we are interested in interaction potentials V including the Coulomb potential.
In this section we also consider mixed states as initial data for the Bogoliubov-de Gennes equation, i. e. generalized one-particle density matrices satisfying 0 ≤ Γ 0 ≤ 1 (whereas before in the derivation we considered only pure states, Γ 2 0 = Γ 0 ).Well-posedness for similar equations has been discussed before, e. g. in [36,41] for a relativistic system (which generally exhibits finite-time blow-up, and global well-posedness only for small initial data) and in [3] for the bosonic Hartree-Fock-Bogoliubov equations.They are generalizations of earlier work on the Hartree-Fock equations without pairing (α = 0) [14,15,12,38]; see also [2].All these works are applications of the abstract formalism developed by Segal [55].

Duhamel Formula and Integral Form
The standard approach to show local well-posedness is through an application of the Banach fixed-point theorem to the integral equation obtained from the Duhamel formula (in the spirit of the Picard-Lindelöff Theorem).This is the strategy we also follow here.There is a small complication as it does not seem possible 5 to apply the Banach fixedpoint theorem in energy space (denoted Y below) when V has a Coulomb singularity.Instead we introduce an additional Banach space Z tailored to the Banach fixed-point theorem.We find local solutions in Z and we show afterward that these give rise to global energy space solutions. 5Let V (x) be the Coulomb potential 1  |x| and M the Fourier multiplier √ 1 − ∆.In energy space (5.50) the norm is given by γ . Consider the mild equation for γt: while we can control the term with [V(γs), γs] with sup 0≤s≤t γs 2 Y 1 , a homogeneity argument shows that we cannot control those with −ΠV (αs)αs and αsΠV (αs) with sup 0≤s≤t αs 2 Y 2 .Indeed we have tr ∂j t 0 e i∆(t−s) ΠV (αs)αse −i∆(t−s) ds ∂j = t 0 tr e i∆(t−s) ∂jΠV (αs)αs∂je −i∆(t−s) ds = t 0 tr ∂jΠV (αs)αs∂j ds, where we have used Fubini's theorem and the invariance of the trace under a unitary conjugation.This expression has contributions scaling like 1/[Length] 3 (the Coulomb potential, and one derivative contributing from each ∂j ) and is of degree 2 in αs.However, a bound in terms of αs 2 Y 2 can accommodate at most two derivatives, i. e. has contributions scaling at most like 1/[Length] 2 .The same holds for the term corresponding to αsΠV (αs).This provides us with reasonable doubts about the validity of the fixed-point scheme in energy space.To be controllable the two terms −ΠV (αs)αs and αsΠV (αs) cannot be separated, and the smoothing effect from the conjugation by e i∆(t−s) does not seem to help to estimate trace-class norms.
Recall the Bogoliubov-de Gennes equations from (1.11) (with h = −∆): (5.47) Using Duhamel's formula, we obtain the integral form6 of the Bogoliubov-de Gennes equations: where in the equation for α t , we are using the identification of the Hilbert-Schmidt operator α t with a two-particle wave function in L 2 (R 3 × R 3 ), denoting action of an operator on the i-th variable (i ∈ {1, 2}) by (•) i , and reading as a two-body Schrödinger operator.

Choice of Banach Spaces
For short we write S 1 = S 1 (L 2 (R 3 )) for the space of trace-class operators on L 2 (R 3 ), )) for that of Hilbert-Schmidt operators, and B = B(L 2 (R 3 )) for that of linear bounded operators (equipped with the operator norm • B ).For another Banach space X , the space of linear bounded operators on X will be written B(X ) with norm Let us introduce the Fourier multiplier M := (1 − ∆) 1/2 .We define the Banach space

usually called the energy space, by
(5.50) Here α T denotes the operator with integral kernel α T (x, y) = α(y, x).Here and below α(•, •) H 1 refers to the norm of α in H 1 (R 3 × R 3 ).We define the Banach space Z := Z 1 × Z 2 for the purpose of constructing local solutions by the Banach fixed-point theorem by For any operator A we have By (5.51), since the Hilbert-Schmidt norm is smaller than the trace norm, we also have To shorten notations, we sometimes identify an operator with its integral kernel and write α H 1 , hα and e ith α.Lemma 5.1 (Invariance under the Linear Evolution).Let V satisfy (5.53), and let the two-particle Schrödinger operator h : Then Y and Z are invariant under the group action where e ith acts on the integral kernel α(•, Proof.The action on γ is given by conjugation with the unitary Fourier multipliers e is∆ , conserving the self-adjointness and the S 1 -, Z 1 -and Y 1 -norms.For α ∈ H 2 (R 3 × R 3 ) = dom(h), the equality (hα) T = h(α T ) is straightforward, thus the action of e ith preserves the transpose symmetry (here α T = −α).The Hilbert-Schmidt norm is left invariant because h is self-adjoint.As V (x − y) is infinitesimally form-bounded w. r. t. −∆ x − ∆ y , the Sobolev space H 2 (R 3 × R 3 , dxdy) (being the domain of h) is invariant under e ith , and so is H 1 (R 3 × R 3 , dxdy) as its form domain.

Results on Well-Posedness
We now study the existence of solutions to the time-dependent Bogoliubov-de Gennes equations.As usual we expect them to conserve the number of particles tr(γ) and the energy of the system We need the following assumptions on the potential V : Observe that in fact the condition , where C 3 > 0 only depends on the dimension.This notation of C V will be used throughout this paper.Lemma 5.2 (Local Well-Posedness in Z).Assume that V satisfies (5.53).Consider a pair of initial data (γ 0 , α 0 ) ∈ Z. Then there exists a unique Z-continuous solution to the mild equations (5.48).
If [0, T ) is the maximal interval of existence of the solution, we have the usual blow-up alternative: either T = +∞ or lim t→T − ω t Z = ∞.Lemma 5.3 (Regularity of the Solution).If the pair of initial data (γ 0 , α 0 ) ∈ Z satisfies [−∆, γ 0 ] ∈ S 1 and α 0 (•, •) ∈ H 2 , then the mild solution of the previous lemma is a strong solution The reader might wonder why we do not simply refer to [55, Lemma 3.1].Transposed here, it states that if [−∆, γ 0 ] ∈ Z 1 and hα 0 ∈ Z 2 , that is α 0 ∈ dom S 2 (|h| 3/2 ), then the mild solution is a strong solution with time derivative in Z.When V is the Coulomb potential, even if a kernel α(•, •) is of Schwartz class, hα is generally not in H 1 (R 3 × R 3 ), unless the diagonal α(x, x) identically vanishes.Here, thanks to the transpose symmetry α T 0 = −α 0 , the result still has an important value, but Lemma 5.3 is in some sense optimal: it states regularity when the minimal requirements are satisfied.Lemma 5.4 (Conservation Laws).Let V satisfy (5.53) and let T > 0. Let (ω t ) 0≤t<T ∈ C([0, T ), Z) be a solution to (5.48).Then the expected particle number tr(γ t ) is conserved, there is a unitary propagator U (t, s) such that Γ t = U (t, 0)Γ 0 U (t, 0) * , and the spectrum of Γ t is conserved.
Remark.In the previous Lemma we have seen that the expected number of particles tr γ t is always a conserved quantity.Now let us introduce the quantity t − α t α) ≤ 2 tr γ t .It coincides with 2 tr γ t if and only if Γ t describes a pure quasifree state, and it is also a conserved quantity because Γ t is unitarily equivalent to Γ 0 and because there holds We can interpret (5.54) as a measure of the deviation from Γ t being pure quasifree.
In the bosonic case, the same role is played by tr (P t − P vac ) 2 + 2 tr γt = tr P 2 t − P t = 2 tr γt (1 + γt ) − αt αt , where P t = − Γt S, P vac = −Γ vac S, γt is the truncated one-particle density matrix, and Γt is the truncated generalized one-particle density matrix.To see that this is conserved, simply notice that tr P 2 t − P t = tr Γt S Γt + Γt S, and by [3, Lemma 3.10] Γ0 = U t Γt U * t for some symplectomorphism U t (meaning that U * t SU t = S = U t SU * t ).As a word of caution: for bosonic systems, tr γt is not conserved by itself, only the particle number tr γ t = tr γt + ϕ t 2 is conserved.
its spectral decomposition.By monotonicity of the square root, and writing down the spectral decomposition of a ± , we obtain For γ ∈ Z 1 self-adjoint, by splitting γ = γ + − γ − (positive and negative part), we find Now with (λ i ) i denoting the eigenvalues of M 1/2 γM 1/2 and (ϕ i ) i a corresponding orthonormal basis of eigenvectors, we thus get Using cyclicity of the trace in the last estimate, and afterward 4. Let us now prove the estimate on X V (γ)M −1 S 2 ( Π V (α)M −1 S 2 is estimated in the same way).We use the idea of [3,Lemma E.1].By cyclicity of the trace, for A ∈ S 2 , we have A S 2 = A * S 2 , and it suffices to show the boundedness of where M −2 (y − z) denotes the Yukawa potential at point y − z.For almost all x, we consider the L 2 -function g x (y) := γ(x, y)V (x − y).By operator monotonicity of the inverse, we have We obtain the upper bound: We introduce the "polarization" of K 1 and K 2 as a bilinear form (it is not necessary to give it a symmetrized form): (5.55) (5.56) Lemma 5.8.Let V satisfy (5.53).Then for the nonlinearities K 1 and K 2 seen as bilinear maps (as defined below in (5.55) and (5.56)) we have the estimates where the constant C depends only on C V .

Local
To establish local existence we show that the nonlinearities (when simply taking the norm inside the integral and neglecting the unitaries-for the second equation this is possible because h is infinitesimally operator bounded w. r. t. the Laplacian) are locally Lipschitz; then ♣ on a sufficiently short time interval I is a contraction, and a solution is found by the Banach fixed-point theorem.We refer to [55, Theorem 1] for details.
To prove local Lipschitz continuity of the quadratic terms K 1 and K 2 , it suffices to show continuity of the corresponding bilinear maps K 1 (•, •) and K 2 (•, •) introduced in (5.55)- (5.56).Indeed, we can write down the difference of the quadratic terms in terms of the polarization according to the following formula: Similarly, using the pull-through formula for Π V (α), we obtain Combining everything we get M K 1 (ω 1 , ω 2 ) S 1 ≤ C ω 1 Z ω 2 Z .We can estimate the term K 1 (ω 1 , ω 2 )M S 1 in the same way.We conclude that To prove continuity of the nonlinearity K 2 , we use (5.51) as an upper bound for K 2 2 H 1 , and then proceed by the same method as for K 1 .We obtain 5.6 Regularity of the Solution (Proof of Lemma 5.3) We now assume that the initial data satisfy also the additional regularity conditions [γ 0 , −∆] ∈ S 1 and α 0 (•, •) ∈ H 2 (or equivalently α 0 (•, •) ∈ dom S 2 (h)).It suffices to adapt Segal's result [55,Lemma 3.1] to ensure that the solution is a strong solution: instead of the Z-norm, we apply the same argument to the S 1 × S 2 -norm.
Identifying γ t and α t with their integral kernel, the derivatives γt , αt are well-defined as bounded linear operators from H 2 (R 3 ) to its dual H −2 (R 3 ), and they are equal to −i times the r.h.s. of (5.47).This establishes the equations, and it just remains to prove that ω t is Fréchet-differentiable in S 1 × S 2 .
We construct the putative derivatives by the Banach fixed-point theorem.We start by formally differentiating the Bogoliubov-de Gennes equations to see that the derivatives should satisfy the equations where ∂ ∂ω K 1 (ω t ) and ∂ ∂ω K 2 (ω t ) denotes the Fréchet derivative of the nonlinearities.As initial data for the fixed-point problem of the derivatives we have (as given by the Bogoliubov-de Gennes equations) We now write the equations for γt and αt in mild form: (5.59) As K 1 and K 2 are quadratic functions of ω, we have for all δω ∈ Z where the K j (•, •)'s are the (un-symmetrized) polarizations as defined in (5.55) and (5.56).By Lemma 5.8, their extensions ∂ ∂ω K j (ω t ) : S 1 × S 2 → S j to S 1 × S 2 ⊃ Z are continuous with norm controlled by C V ω t Z .So we can apply the Banach fixed-point theorem to (5.59) in the Banach space S 1 × S 2 .We obtain a unique local solution v t = (g t , a t ) t∈[0,T ) .
It is easy to see through a simple Grönwall argument that v t S 1 ×S 2 growths at most like exp Ct sup s∈[0,t] ω s Z , hence the maximal interval of existence of v t is the same as that of ω t .
Following the proof of [55, Lemma 3.1], we show that w ε (t) := ε −1 (ω t+ε − ω t ) − v t converges to 0 in S 1 × S 2 as ε → 0 by another Grönwall argument.We fix 0 < T 1 < T and work on [0, T 1 ].It is convenient to write the mild equations in terms of ω s and v t : By the Hille-Yosida theorem the first line tends to 0 uniformly in t ∈ [0, T ) (the operators e tA are unitary in S 1 × S 2 ).By the change of variable εu = s the second line is By dominated convergence it converges to 0 uniformly in t ∈ [0, T 1 ].Let us now deal with the third line.As K : Z → Z is • Z -differentiable and that (ω t ) 0≤t<T is Z-valued, we have where T s,ε : Z → Z is a convex combination of the Fréchet derivatives in the integrand.Here K is quadratic, hence T s,ε = 1 2 ( ∂ ∂ω K(ω s+ε ) + ∂ ∂ω K(ω s )).By Lemma 5.8, T s,ε B(S 1 ×S 2 ) is uniformly bounded on [0, T 1 ], and since the convergence lim ε→0 ω s+ε − ω s Z = 0 holds point-wise, so does lim ε→0 T s+ε − ∂ ∂ω K(ω s ) B(S 1 ×S 2 ) = 0. We decompose the third line: By dominated convergence, the second integral converges to 0 as ε → 0, uniformly in t ∈ [0, T 1 ].By Lemma 5.8, the norm of the first integral is bounded by Putting everything together, we obtain the integral inequality for t ∈ [0, T 1 ]: , and the solution ω t is differentiable in S 1 × S 2 on [0, T 1 ].As 0 < T 1 < T was arbitrary, this shows that the same holds on the whole interval [0, T ).

Existence of Unitary Propagator and Conservation of Spectrum
In this section we prove that the solution Γ t at time t of the Bogoliubov-de Gennes equation is related to Γ 0 by conjugation with a unitary propagator.This implies that the spectrum of Γ t is time-independent.
For regular initial data We start with the case where [−∆, γ t ] ∈ S 1 and α(•, •) ∈ H 2 (R 6 ).We split F Γt into unbounded time-independent and bounded time-dependent part as Using the estimates of Lemma 5.7, we find that t → B(∂ t Γ t )M −1 is continuous (w.r. t. the operator norm).
Regularization Since the kinetic part of the energy functional is not Z-continuous, we cannot use the same strategy as for particle number conservation.Instead, we introduce a regularization for which the conservation of energy holds.As before P Λ denotes 1 (−∆<Λ) , and to shorten notations, we also denote P Λ ⊗ 1 C 2 (acting on L 2 (R 3 ) 2 ) by P Λ .We regularize both the equation and the functions: for any Λ > 0, we consider the solution (Γ (Λ) t ) t to: i∂ t S t = P Λ F St P Λ , S t , S 0 = P Λ Γ 0 P Λ . (5.63) Above, P Λ F St P Λ denotes the bounded operator: The Duhamel form of the equation is similar to (5.48).The infinitesimal generators of the free evolutions of the one-body and the pairing densities are are obtained from the original ones K 1 , K 2 by replacing V(γ) and Π V (α) by P Λ V(γ)P Λ and P Λ Π V (α)P Λ respectively.
We can apply the Banach fixed-point theorem to the regularized equations with the • Z -norm: the estimates are the same, and the interval of existence [0, T ) only depends on the initial data Γ 0 , it does not depend on the cutoff Λ > 0. For any Λ > 0, we thus obtain a solution (Γ for the corresponding objects of the regularized solution.
As the operator Z ), mild solutions to (5.63) are also strong solutions as we can differentiate in the Duhamel formulas for γ and α.Indeed, by the Hille-Yosida theorem the integrand is point-wise differentiable and by dominated convergence [19,Theorem III.6.16]we can differentiate inside the integral.Furthermore, by [19,Theorem III.6.20]we can pull the bounded operators [−∆P Λ , •] and h Λ outside the integral, where h Λ denotes the operator h Λ α := {−P Λ ∆, α} + P Λ Π V (α)P Λ , α(•, •) ∈ H 2 . (5.64) Since the regularized equation has the same structure as the original one, conservation of the spectrum still holds (by the argument we gave before).In particular we have: 0 ≤ Γ (Λ) t ≤ 1.We show consecutively the following four points.(5.65) • For the second technical issue, we introduce a second level of cutoff Λ ′ > 0. For 0 ≤ s ≤ T 1 , the operator K 2 (ω s ) is compact (and its integral kernel is in H 1 ), hence we have point-wise in s: On [0, T 1 ], the norm K 2 (ω s ) S 2 is uniformly bounded by sup t∈[0,T 1 ] K 2 (ω t ) < ∞.
Hence by dominated convergence, we obtain Consider now the second line (5.69).Taking Λ ′ = Λ, we get that it converges to 0 as Λ → ∞.
Consider now the fourth line (5.71).The first term converges to 0 since we have lim Λ→0 P Λ α 0 P Λ − α 0 S 2 = 0.For its second term, splitting α 0 into two: The second summand vanishes as Λ ′ → ∞.At fixed Λ ′ , we then have This follows from the convergence in the strong operator topology of e −ith Λ and the fact that P Λ ′ α 0 P Λ ′ ∈ dom S 2 (h) = H 2 which gives by functional calculus the crude estimate where t 1 , t 2 ∈ R and where h denotes h or h Λ .Hence by an ε/2-argument, the fourth line converges to 0 as Λ → ∞.
Therefore the third line (5.70)tends to 0 as Λ → ∞ (by an additional ε/2-argument used to choose the auxiliary cutoff level Λ ′ > 0 at the very beginning).
As in Lemma 5.8, by using Lemma 5.7 we obtain the following estimate: where Γ i denotes the generalized density matrix corresponding to ω i .Point (ii) and Point (iii), namely the energy conservation of Γ  Notice also that by assumption γ ≥ 0, so M 1/2 γM 1/2 S 1 = tr M 1/2 γM 1/2 .The first term in the integral (the direct term) can be estimated using Lemma 5.
To estimate M 1/2 γM 1/2 S 1 we employ the spectral decomposition γ = i λ i |ψ i ψ i | (where all λ i ≥ 0) to get The third term term is of the same type as the previous one, so we easily find By assumption αα * + γ 2 ≤ γ, and thus α 2 From the estimates above (and adjusting the choice of δ) we conclude that Finally, by using the symmetry α T = −α and going to Fourier space, we get which is estimated as above.This concludes the proof of the first bound.
For the second bound, by the above estimates, choosing δ = 1, we obtain 5.9 Global Well-Posedness (Proof of Theorem 5.6) Observe that the conservation laws (Lemma 5.4) together with Lemma 5.5 imply that the maximal interval of existence for Γ t is [0, ∞).By energy conservation, the solution lies in Y. Since we have conservation of the spectrum of Γ t we have in particular γ t ≥ 0 for all times t.Thus γ t Y 1 = tr((1 − ∆)γ t ).
The conservation of the energy E(Γ t ) and of the particle number tr γ t together with the Z-continuity of E pot give tr((1 − ∆)γ t 0 ) = lim n→∞ tr((1 − ∆)γ tn ).This concludes existence and continuity in Y for positive times.By time-reversal symmetry (defined as in Sect.5.7.3), we also obtain the solution for negative times.

Figure 1 :
Figure1: The Dirac-Frenkel principle: Consider u ∈ M. At every "time step", the tangent 1 i Hu of the exact evolution is orthogonally projected into the tangent space T u M, yielding the optimal effective evolution in the approximation manifold M. Figure following[42].

Lemma 3 . 1 (
Tangent Space, No Pairing).The tangent space in a point γ ∈ M is
. More compactly, γ t and α t satisfy (1.11) if and only if the generalized one-particle density matrix Γ t satisfies i∂ t Γ t = [F Γt , Γ t ]; (1.12) 43))So we can now formulate the Dirac-Frenkel principle for the condensate wave function and the generalized reduced density matrix of bosonic Bogoliubov states: Calculate the derivative of the condensate wave function evolving by the many-body Hamiltonian H in the Bogoliubov state associated with ϕ t and Γqf then apply the projection onto the tangent space to ∂ t ϕ t .Calculate the derivative of the generalized one-particle reduced density matrix evolving by the many-body Hamiltonian H in the quasifree state associated with Γqf t , f, ∂ t ϕ t = ψ bog t , [a(f ), t , and in particular unitary (see the proof of [51, Theorem X.71]).Let us now verify the C 1 -condition.Since we have assumed regular initial data, the mild solution Γ t is continuously differentiable.Since ϕ ∈ H 2 (R 3 ), we can insert 1 = M −1 M and obtain d dt t := B(Γ t ) prescribed by the solution Γ t of the nonlinear Bogoliubov-de Gennes equations (5.47))i∂ t U (t, s) = A + B t U (t, s), U (s, s) = 1,provided that the domain D(A + B t ) is independent of t and the function t → (A + B t )ϕ is continuously differentiable for every ϕ ∈ D(A + B t ).The solution U (t, s) then is a propagator (i.e. U (t, s)U (s, r) = U (t, r) for all r, s, t ∈ R), Λ (K 2 (ω s ) − P Λ ′ K 2 (ω s )P Λ ′ ) S 2 ds −→ T 1 ] is compact, and the bilinear map K 2 (•, •) : Z 2 → S 2 is continuous, the map t ∈ [0, T 1 ] → K 2 (ω t ) is • S 2 -equicontinuous.At fixed s ∈ [0, T 1 ] and Λ ′ > 0, the integral kernel of P Λ ′ K 2 (ω s )P Λ ′ is in H 2 ,hence a similar estimate to (5.72) holds for this operator.By an ε/3-argument, we get the uniform estimate sup s,t∈[0,T 1 ] × S 2 imply the convergence of E pot (Γ t ) to E pot (Γ t ).