We slightly extend prior results about the derivation of the Maxwell–Schrödinger equations from the bosonic Pauli–Fierz Hamiltonian. More concretely, we show that the findings from Leopold and Pickl [SIAM J. Math. Anal. **52**(5), 4900–4936 (2020)] about the coherence of the quantized electromagnetic field also hold for soft photons with small energies. This is achieved with the help of an estimate from Ammari *et al.* [arXiv:2202.05015 (2022)], which proves that the domain of the number of photon operator is invariant during the time evolution generated by the Pauli–Fierz Hamiltonian.

## I. INTRODUCTION

In this short paper, we derive the Maxwell–Schrödinger system of equations as an effective model describing a Bose–Einstein condensate of charged particles immersed in a coherent electromagnetic field. More precisely, we prove quantitatively that the Maxwell–Schrödinger system approximates well the many-body quantum evolution generated by the Pauli–Fierz Hamiltonian, provided that the total number of particles *N* is large, the particles are initially in a Bose–Einstein condensate, and the quantum nature of the field—quantified by the semiclassical parameter ℏ—is negligible. In particular, we focus on the combined regime $N\u223c1\u210f\u2192+\u221e$. Equivalently, the same effective dynamics approximates well *N* → ∞ bosons weakly interacting with a quantized electromagnetic field; see the discussion below.

This problem has already been studied by one of the authors, together with Pickl, in Ref. 1. The main focus here is to build on the results and techniques introduced there and to strengthen them by studying the convergence for the photons’ reduced density matrix. In the previous work, the quantum fluctuations around the coherent state of photons have been classified only by means of their energy. The extension to the reduced density matrix is physically relevant and mathematically nontrivial because the coherence of photons with small frequencies cannot be shown by the energy of the electromagnetic field due to its massless nature. We often refer to Ref. 1 throughout this paper, hopefully striking a good balance between being concise and being self-contained.

### A. The Maxwell–Schrödinger system of equations

The Maxwell–Schrödinger system of equations describes the wave function *φ* of a quantum particle (with a nontrivial charge distribution *κ*) interacting with the classical electromagnetic field, described by the vector potential **A** and the electric field $E=\u2212A\u0307$. We choose the Coulomb gauge

and this makes, indeed, **A** and **E** the only dynamical degrees of freedom of the field. Let us also preliminarily define the current

The Maxwell–Schrödinger system thus takes the form

where $V[\phi ]$ is an interaction term for the quantum particle and $1\u2212\u2207(\u2207\u22c5\Delta \u22121)$ is the Helmholtz projection onto the subspace of divergence free vector fields. The choice of gauge, Coulomb’s in this case, can be seen as a constraint, for it is preserved by the Maxwell–Schrödinger flow. A typical example for the particle interaction $V[\phi ]$ could be

where *W* is an external potential and ($v*$ |*φ*|^{2})*φ* is a nonlinear term, usually originating from a microscopic pair interaction. With this choice, (1.3) can used as a system of mean-field equations to model the dynamics of identical particles interacting via the potential *v* and the classical electromagnetic field which they generate. The Cauchy problem associated with (1.3) is obtained by fixing an initial datum (*φ*_{0}, **A**_{0}, **E**_{0}), subjected to the constraint ∇ · **A**_{0} = 0. In order to do so, it is convenient to introduce the complex scalar fields $\alpha 0(\u22c5,\lambda )\lambda =1,2$ by defining

where $\u03f5\lambda (k)\lambda =1,2$ are the polarization vectors satisfying

which implement the Coulomb gauge. In fact, there is such a unique decomposition for any time, i.e.,

which respects both the Coulomb gauge and $A\u0307=\u2212E$. This makes it possible to rewrite Maxwell’s equations (same as in Ref. 2, Chap. 1.C) in terms of *α*_{t} and to consider the equivalent system

with initial datum (*φ*_{0}, *α*_{0}(·, 1), *α*_{0}(·, 2)). Throughout this article, we use $F[f]$ to denote the Fourier transform of *f*. As it will be clarified shortly, *φ*_{t} and *α*_{t} appear naturally as the effective counterparts of the microscopic dynamical variables. Note that the energy functional of the Maxwell–Schrödinger system is given by

with ** A** being defined in analogy to (1.5). Global well-posedness for the Maxwell–Schrödinger system with $V[\phi ]=v*|\phi |2\phi $,

*κ*(

*x*) = e

*δ*(

*x*) (where e is the electric charge of the Schrödinger particle), and $v(x)=e2|x|$ has been proven in Refs. 3 and 4. We will also consider only the case $V[\phi ]=v*|\phi |2\phi $, but we may require the charge distribution

*κ*to be extended, in order to well-define the microscopic system, as discussed below. Typical examples of charge distributions that we will consider are of the form

representing a charged particle with total charge $e\u2208R$, distributed in a Gaussian fashion (“smoothed” spherical distribution of “diameter” *σ*), or

representing a sharp cutoff in momentum space with total charge $e\u2208R$. Let us remark that globally neutral particles can be considered as long as they have a nontrivial charge distribution: for example,

yields null total charge but nontrivial dipole and quadrupole interactions with the electromagnetic field.

If the charge distribution is not concentrated in a single point and the potential *v* represents an electrostatic mean-field self-interaction, then the form of the latter changes as well: a physically sensible choice would be $v=\kappa *1|\u22c5|*\kappa $. We will allow some liberty in the choices of *κ* and *v*; the specific requirements on the two will be made precise in Assumption 2.1. Concerning global well-posedness, let us remark that compared to the literature,^{3,4} our choices for *κ* and *v* will be at most “better” (i.e., more regular) and therefore do not affect the proof in any way since both *κ* and *v* act by convolution in the equation.

*(Ref.* *3 **). The Maxwell–Schrödinger system in Coulomb’s gauge* *(1.3)* *[with variables* (*φ*, **A**, **E**)*] is globally well-posed in* *H*^{1} × *H*^{1} × *L*^{2}*. More precisely, we have the following:*

*(rough solutions) For every*$(\phi 0,A0,E0)\u2208H1\xd7H1\xd7L2,$*there exists a unique global solution*$\phi ,A\u2208C0(R,H1)\xd7C0(R,H1)\u2229C1(R,L2)$*of the Cauchy problem associated with**(1.3)**, being the unique strong limit of a sequence of regular solutions in (1), whose initial data approximate the rough initial datum in**H*^{1}×*H*^{1}×*L*^{2}*.**(continuous dependence on initial data) The solutions*(*φ*,**A**)*in (2) continuously depend on the initial datum*(*φ*_{0},**A**_{0},**E**_{0}) ∈*H*^{1}×*H*^{1}×*L*^{2}*.*

For $m\u2208R$, let $hm$ denote the weighed $L2(R3)\u2297C2$-space with norm

Throughout this work, we will rely on the following statement, which results almost immediately from Proposition 1.1 (see the Appendix).

*Let* $\u22c5\u22121/2F[\kappa ]\u2208L2(R3,C)$*. For every initial datum* $(\phi 0,\alpha 0)\u2208H2(R3,C)\xd7h32$ *such that* $\u22c5\u22121/2\alpha 0\u2208L2(R3)\u2297C2$, *the Maxwell–Schrödinger system* *(1.10)* *has a unique global solution in* $H2(R3,C)\xd7h32$*.*

### B. The microscopic model: Pauli–Fierz Hamiltonian

The microscopic model corresponding to the Maxwell–Schrödinger system with mean-field self-interaction $v*$ |*φ*|^{2}*φ* consists of many identical nonrelativistic particles—obeying Bose–Einstein condensation—interacting among themselves by means of a weak pair potential and with a quantized electromagnetic field in Coulomb’s gauge. Contrarily to the “classical” case, the microscopic model is known to be well-defined only for extended charges. Let us start by defining a Hilbert space $H\u210f(N)$ depending on two parameters $N\u2208N,\u210f\u2208R+$ as follows:

where $Ls2(R3N)$ is the natural Hilbert space of *N* identical bosons (the subscript s indicates symmetry under the interchange of variables) and $\Gamma \u210f$ is the second quantization functor associating with any (pre-)Hilbert space $h=L2(R3)\u2297C2$ the corresponding Fock representation of the canonical commutation relations $[a\u210f(f),a\u210f*(g)]=\u210f\u27e8f,g\u27e9h,$ with ℏ being a semiclassical parameter measuring the degree of noncommutativity of the quantum field. The Fock representation is the natural one to describe noninteracting or regularized quantum field theories, the latter being the case here with $h\u2254L2(R3)\u2297C2$. With this interpretation, the limits *N* → ∞ and ℏ → 0 describe, respectively, the regimes in which the bosons are many and the quantum effects of the field are negligible. The time evolution is dictated by the Schrödinger equation,

where the Hamiltonian $HN,\u210f$, called Pauli–Fierz Hamiltonian, is given by

where $\mu N,\u210f$ describes the coupling strength between the particles and the field, *g*_{N} is the coupling strength between the particles, *κ* and *v* are the charge distribution and the pair potential introduced previously,

is the field’s kinetic energy, with $a\u210f\u266f(k,\lambda )$ being the polarized creation and annihilation operators satisfying the canonical commutation relations (CCR),

and^{37}

is the smeared quantized electromagnetic vector potential in Coulomb’s gauge. Let us remark that both $A\u0302\kappa (x)$ and *H*_{f} depend on ℏ through the creation and annihilation operators, which have ℏ-dependent CCRs. The Hamiltonian $HN,\u210f$ is self-adjoint on $D(HN,\u210f)=D(HN,\u210f(0))$, where $HN,\u210f(0)=HN,\u210f|\mu N,\u210f=gN=0$ whenever $\u22c5\u22121+\u22c51/2F[\kappa ]\u2208L2(R3)$ and *v* is Kato-infinitesimal with respect to −Δ.^{5–7}

### C. Scaling regime

Our aim is to prove that the Maxwell–Schrödinger system emerges in some limit *N* → ∞ and/or ℏ → 0 as an effective model of the microscopic Pauli–Fierz dynamics. This is true only if we couple the parameters *N*, ℏ suitably and choose the coupling constants $\mu N,\u210f$, *g*_{N} accordingly. A possible choice is given by *N* → ∞, $\u210f=1N$, $\mu N,\u210f=1$, and $gN=1N$. In this regime, the electromagnetic field becomes classical inverse proportionally to the increasing number of bosons. At the same time, the coupling between the particles and field is of order one, while the coupling between pairs of particles becomes weak (of order $1N$). Defining $N1/2aN\u22121(k,\lambda )$ and $N1/2aN\u22121*(k,\lambda )$ as new creation and annihilation operators leads to the mathematically equivalent but physically different choice *N* → ∞, ℏ = 1, $\mu N,\u210f=1N$, and $gN=1N$. Here, the physical interpretation is of many bosons that interact weakly both with the quantized electromagnetic field (coupling of order $1N$) and among themselves (pair coupling of order $1N$). Our result reads as follows.

*Provided that we choose an initial microscopic state that is “close enough” to a non-interacting state representing a complete condensate and a coherent field of minimal uncertainty, then at any time* *t* ≥ 0, *the evolution keeps the state “close” (in the same sense as above) to an analogous configuration in which the one-particle wave function and the argument of the coherent field have been evolved by the coupled Maxwell–Schrödinger equations.*

The described scaling regime has been considered in earlier works for the Nelson model with an ultraviolet cutoff,^{8–11} the renormalized Nelson model,^{12} and the Fröhlich model.^{13,14} In Ref. 15, the Nelson model with an ultraviolet cutoff has been studied in a limit of many weakly interacting fermions. The classical behavior of quantum fields has also been proven in different scaling regimes.^{16–31} We also would like to mention Ref. 32, which derives the Maxwell–Schrödinger equations in a nonrigorous manner by neglecting certain terms in the Pauli–Fierz Hamiltonian.

## II. MAIN RESULT

From now on, we will keep *N* as the single parameter and choose ℏ = 1, $\mu N,\u210f=1N$, and $gN=1N$. We will use the notations $H(N)=H\u210f(N)$, $Fp=\Gamma 1L2(R3)\u2297C2$ with vacuum Ω, *H*_{N} = *H*_{N,1}, and Ψ_{N} = Ψ_{N,1}. Concerning the interaction potential and charge distribution, we will make the assumptions.

*The (repulsive) interaction potential*

*v*

*is a positive, real, and even function satisfying*

*The charge distribution*

*κ*

*with Fourier transform*$F[\kappa ]$

*satisfies*

In order to state our result, we define for $\Psi N\u2208H(N)$ the one-particle reduced density matrix of the charged particles $\gamma \Psi N(1,0):L2(R3)\u2192L2(R3)$ by

where Tr_{2,…,N} denotes the partial trace over the coordinates *x*_{2}, …, *x*_{N} and $TrFp$ is the trace over Fock space. In addition, we introduce the number of photon operators,

and the unitary Weyl operator $(f\u2208h)$,

satisfying

Our result is as follows.

*Let*

*v*

*and*

*κ*

*satisfy Assumption 2.1,*$(\phi 0,\alpha 0)\u2208H2(R3,C)\xd7h32$

*such that*$\u22c5\u22121/2\alpha 0\u2208h$

*and*$\phi 0L2(R3)=1$

*and*$\Psi N,0\u2208DHN\u2229DN1/2$

*such that*$\Psi N,0H(N)=1$

*. Define*

*Let*(

*φ*

_{t},

*α*

_{t})

*and*Ψ

_{N,t}

*be the unique solutions of*

*(1.10)*

*and*

*(1.17)*

*, respectively. Then, there exists a function*

*C*(

*s*)

*depending monotonically on the norms*$\phi sH2(R3)$

*,*$\u22c51/2\alpha sh$

*,*$vL2+L\u221e(R3)$,

*and*$\u22c5\u22121+\u22c51/2F[\kappa ]L2(R3)$

*such that*

*for any*

*t*≥ 0

*. In particular, for*$\Psi N,0=\phi 0\u2297N\u2297W(N\alpha 0)\Omega $,

*one obtains*

*Let*$\gamma \Psi N,t(0,1)$

*be the one-particle reduced density matrix of the photons with the integral kernel*

*By similar means as in Ref.*

*1*

*, Lemma 5.3, one obtains*

*from*

*(2.11)*

*and*

*for initial product states*$\Psi N,0=\phi 0\u2297N\u2297W(N\alpha 0)\Omega $

*from*

*(2.13)*

*.*

*In Ref.*

*1*,

*Theorem 2.2 was proven for the charge distribution*

*(1.13)*

*and with the number operator*$N$

*in*

*(2.8)*

*,*

*(2.11)*

*, and*

*(2.13)*

*being replaced by the field energy*

*H*

_{f}

*. Because of Markov’s inequality,*

*one can use the field energy to conclude that the quantum fluctuations around the coherent state are subleading for all photons with*$k\u2265I$

*. For sufficiently small*

*a*

_{N}

*,*

*b*

_{N},

*and*

*c*

_{N},

*one can choose*

*I*∼

*N*

^{−a}

*with*

*a*< 1

*. However, this choice does not provide information about the coherence of soft photons with frequencies below this threshold.*

## III. PROOF OF THE RESULT

The rest of the article outlines the Proof of Theorem 2.2. We will proceed as follows:

We define a functional

*β*[Ψ_{N},*φ*,*α*], which measures if the charges of the many-body state Ψ_{N}form a Bose–Einstein condensate with condensate wave function*φ*and if the photons are in a coherent state with mean photon number $N\alpha h2$.Next, we show that the domain of

*β*contains the solutions (*φ*_{t},*α*_{t}) of (1.10) and Ψ_{N,t}of (1.17) from Theorem 2.2 for all*t*≥ 0.Afterward, we compute the change of

*β*[Ψ_{N,t},*φ*_{t},*α*_{t}] in time.Finally, we control the growth of

*β*[Ψ_{N,t},*φ*_{t},*α*_{t}] with the help of Grönwall’s inequality. This concludes the proof.

### A. Definition of the functional

We define a functional that consists of three parts.

*For*$\phi \u2208L2(R3)$,

*we define*$p1\phi :L2(R3N)\u2192L2(R3N)$

*by*

*and*$q1\phi :=1L2(R3N)\u2212p1\phi $

*. Now, let*$\Psi N\u2208D(HN)\u2229DN1/2$

*,*$\phi \u2208H1(R3)$

*, and*$\alpha \u2208h12$

*. Then,*

*and the functional*$\beta :D(HN)\u2229DN1/2\xd7H1(R3)\xd7h12\u2192R0+$

*is defined as*

*β*: =

*β*

^{a}+

*β*

^{b}+

*β*

^{c}

*.*

The functional *β*^{a} measures if the charges of the many-body state are in a Bose–Einstein condensate (we refer to Refs. 33 and 34 for a comprehensive introduction). Its relation to the trace norm distance of the one-particle reduced density matrix is given by (see, e.g., Ref. 1, Lemma 5.3)

The functional *β*^{b} quantifies the fluctuations of Ψ_{N} around the coherent state $W(N\alpha )\Omega $. Using property (2.6) of the Weyl operators, it can be written as

showing that it is the same quantity as usually considered in the coherent state approach.^{35} While *β*^{a} and *β*^{b} measure the deviation of Ψ_{N} from the product state $\phi \u2297N\u2297W(N\alpha )\Omega $, the functional *β*^{c} is introduced for technical reasons. It quantifies the fluctuations of the many-body energy per particle around the energy of the Maxwell–Schrödinger system.

In the original proof of Ref. 1, the functional *β* was considered with $\beta b\Psi N,\alpha $ being replaced by

This definition has the advantage that it can be defined for many-body states Ψ_{N} in the domain $DHN=H2(R3N,C)\u2297F\u2229DHf$, which is invariant under the time evolution $e\u2212iHNt$. The additional difficulties with respect to Ref. 1 actually originate from the fact that $DN1/2$ [in contrast to $DHf$] is not contained in the domain of the Pauli–Fierz Hamiltonian. On the contrary, $\beta \u0303b$ does not allow us to investigate the coherence of photons with small frequencies because the factor $k$ in the integral on the right-hand side of (3.5) suppresses contributions from photons with small energies.

### B. Invariance of the domain

Throughout the rest of this article, (*φ*_{t}, *α*_{t}) and Ψ_{N,t} denote the solutions of (1.10) and (1.17) from Theorem 2.2. In this section, we show that $(\Psi N,t,\phi t,\alpha t)\u2208D(HN)\u2229DN1/2\xd7H2(R3)\xd7h32$ for all *t* ≥ 0. The condition on the Maxwell–Schrödinger solutions is satisfied because of Corollary 1.2. While $DHN$ is invariant under the evolution of the Pauli–Fierz Hamiltonian, due to Stone’s theorem, the invariance of $DN1/2$ is less clear because the photon number is not conserved during the time evolution. The next statement, however, displays that the number of photons can be controlled by the energy of the system.^{38}

*Let*$\Psi N,0\u2208DHN1/2\u2229DN1/2$

*. Then, there exists a constant*

*C*

*(depending on*

*N*

*and the choice of*

*κ*

*) such that*

*This implies*$e\u2212iHNtDHN1/2\u2229DN1/2=DHN1/2\u2229DN1/2$

*.*

*δ*≥ 0, and consider the bounded operator $N\delta =Ne\u2212\delta N$. Using $\u2207j,A\u0302\kappa (xj)=0$, we obtain

*C*(

*N*,

*κ*) dependent on the number of particles and the choice of

*κ*such that $HN(0)+11/2\Psi \u2264C(N,\kappa )HN+11/2\Psi $ holds for all $\Psi \u2208DHN1/2$. This fact follows from $DHN(0)=DHN=H2(R3N,C)\u2297Fp\u2229DHf$ and the closed graph theorem (Ref. 5, Theorem 1.3 and Corollary 1.4). In that regard, note that $HN(0)1/2HN1/2+i\u22121$ is a closed operator. We consequently obtain

### C. Computing the change of *β* in time

In order to estimate the emerging correlations between the particles and the photons during the evolution of the system, we compute the change of $\beta \Psi N,t,\phi t,\alpha t$ in time.

*Let*(

*φ*

_{t},

*α*

_{t})

*and*Ψ

_{N,t}

*be the solutions of*

*(1.10)*

*and*

*(1.17)*

*from Theorem 2.2 and*

*β*

^{b}

*be defined as in Definition 3.1. Then,*

*t*from $DN1/2$ to $H(N)$ with (see Ref. 36, Lemma 3.1)

From Ref. 1, Secs. 6.2 and 6.4 and Duhamel’s formula, we immediately obtain the following.

*Let*(

*φ*

_{t},

*α*

_{t})

*and*Ψ

_{N,t}

*be the solutions of*

*(1.10)*

*and*

*(1.17)*

*from Theorem 2.2 and*

*β*

^{a}

*,*

*β*

^{c}

*be defined as in Definition 3.1. Then,*

*Moreover,*

*due to energy conservation.*

### D. Controlling the growth of *β* in time

In this section, we classify the growth of $\beta \Psi N,t,\phi t,\alpha t$ in time. The first two inequalities of Theorem 2.2 then follow from (3.3) and (3.4) and the statement below. By similar estimates as in Ref. 1, Chap. 7, one obtains (2.12) and (2.13).

*Let*(

*φ*

_{t},

*α*

_{t})

*and*Ψ

_{N,t}

*be the solutions of*

*(1.10)*

*and*

*(1.17)*

*from Theorem 2.2. Then, there exists a monotone increasing function*

*C*(

*s*)

*of the norms*$\phi sH2(R3)$

*,*$\u22c51/2\alpha sh$

*,*$vL2+L\u221e(R3)$,

*and*$\u22c5\u22121/2+\u22c5\u22121F[\kappa ]L2(R3)$

*such that*

*holds for any*

*t*≥ 0

*.*

*η*with Fourier transform,

*F*^{♯}with

*♯*∈ {−, +} and use

## ACKNOWLEDGMENTS

N.L. would like to thank Peter Pickl for many fruitful discussions within the project.^{1} M.F. acknowledges the support from “Istituto Nazionale di Alta Matematica (INdAM)” through the “Progetto Giovani GNFM 2020: Emergent Features in Quantum Bosonic Theories and Semiclassical Analysis.” N.L. acknowledges the support from the Swiss National Science Foundation through the NCCR SwissMap and funding from the European Union’s Horizon 2020 Research and Innovation Programme under Marie Skłodowska-Curie Grant Agreement No. 101024712.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**Marco Falconi**: Writing – original draft (equal); Writing – review & editing (equal). **Nikolai Leopold**: Writing – original draft (equal); Writing – review & editing (equal).

## DATA AVAILABILITY

Data sharing is not applicable to this article as no new data were created or analyzed in this study.

### APPENDIX: PROPERTIES OF (1.10)

*φ*

_{0},

**A**

_{0},

**E**

_{0}) ∈

*H*

^{2}×

*H*

^{2}×

*H*

^{1}. The existence of a unique global solution $\phi ,\alpha $ with $\phi t\u2208H2(R3)$ and $\u22c51/2+\u22c53/2\alpha t\u2208h$ then follows from Proposition 1.1. In order to see that $\alpha t\u2208h$, we bound the integral version of (1.10) by

## REFERENCES

To simplify the notation we assume $F[\kappa ](k)\u2208R$ for all $k\u2208R3$. Theorem 2.2 equally applies if $F[\kappa ]$ is complex valued. In this case, $A\u0302\kappa (x)=\u2211\lambda =1,2\u222bd3k12k\u03f5\lambda (k)F[\kappa ](k)\u0304eikxa\u210f(k,\lambda )+F[\kappa ](k)e\u2212ikxa\u210f*(k,\lambda )$.

Inequality (3.6) was originally proven by Hiroshima and appeared in a slightly different form (for the second instead of the first moment of the number operator) in Ref. 16, Proposition 3.11. We would like to thank Hiroshima for sharing his notes with us. The proof is presented again for the convenience of the reader.