Consistency Proof for Multi-time Schrödinger Equations with Particle Creation and Ultraviolet Cut-Off

For multi-time wave functions, which naturally arise as the relativistic particle-position representation of the quantum state vector, the analog of the Schrödinger equation consists of several equations, one for each time variable. This leads to the question of how to prove the consistency of such a system of PDEs. The question becomes more difficult for theories with particle creation, as then different sectors of the wave function have different numbers of time variables. Petrat and Tumulka (2014) gave an example of such a model and a non-rigorous argument for its consistency. We give here a rigorous version of the argument after introducing an ultraviolet cut-off into the creation and annihilation terms of the multi-time evolution equations. These equations form an infinite system of coupled PDEs; they are based on the Dirac equation but are not fully relativistic (in part because of the cut-off). We prove the existence and uniqueness of a smooth solution to this system for every initial wave function from a certain class that corresponds to a dense subspace in the appropriate Hilbert space.


Introduction
The quantum state of N particles is usually described by means of a wave function Ψ(t, x 1 , . . . , that is a function of time t ∈ R and the positions x j ∈ R 3 of the particles. The obvious relativistic generalization is a wave function of N space-time points (t j , x j ) ∈ M = R 4 (setting c = 1), called a multi-time wave function [4,8,29]. The relation between Ψ and Φ is straightforward from the fact that Ψ also refers to N space-time points, (t, x 1 ), . . . , (t, x N ), which are simultaneous relative to the chosen Lorentz frame; that is, the single-time wave function Ψ is recovered from Φ by setting all times equal, The domain of Φ is usually the set of spacelike configurations, S (N ) := (x 1 , . . . , x N ) ∈ M N : ∀j, k : (x j −x k ) μ (x j −x k ) μ < 0 or x j = x k , relative to the Minkowski metric diag(1, −1, −1, −1). The advantage of using Φ is that it is a covariant object, defined without reference to any hypersurfacenor in fact to coordinates, if we think of Φ as a function of N points in the space-time manifold M . The concept of a multi-time wave function Φ is closely related to that of associating with every spacelike hypersurface Σ a wave function Ψ Σ , as used by Tomonaga [50] and Schwinger [47]: Given Φ, we can define Ψ Σ on Σ N by simply setting Still, Φ is a more elementary concept than Ψ Σ , as Φ is simply a function of 4N variables. Moreover, as we will see, the multi-time wave function Φ can be adapted to the situation with an ultraviolet (UV) cut-off, while Ψ Σ cannot. In a quantum field theory (QFT), the particle-position representation of the quantum state vector |Ψ in Fock space, whenever that exists, naturally yields a multi-time wave function Φ(x 1 , . . . , x N ) := ∅|a(x 1 ) · · · a(x N )|Ψ , (6) where |∅ is the Fock vacuum and a(x) is the annihilation operator in the position representation at x = (x 0 , x) ∈ M . Since N is now variable, this Φ is a multi-time Fock function, i.e., a function on We take the approach of defining a QFT model in a manifestly covariant way (without ever choosing a Lorentz frame or spacelike hypersurface) by starting from Φ as a solution of a suitable multi-time variant of the Schrödinger equation. However, while this procedure is physically the ultimate goal, we have to introduce an UV cut-off for the sake of mathematical rigor, which breaks Lorentz invariance. That is why we here formulate a multi-time evolution law for Φ corresponding to a simplified QFT model with cut-off in a fixed Lorentz frame. In that setting, we can then prove existence and uniqueness of solutions for sufficiently regular initial conditions, which demonstrates that the approach makes sense.

Multi-time Evolution Laws
In order to discuss multi-time evolution laws, let us begin again with the case of a fixed number N of particles. As Φ depends on N time variables, its evolution can be governed by a system of N Schrödinger equations (setting = 1): with "partial Hamiltonians" H k . By (3) and the chain rule, the sum of the H k should be the Hamiltonian for Ψ at every simultaneous configuration.
As a consequence of ∂ tj ∂ t k Φ = ∂ t k ∂ tj Φ, we obtain a condition on the H k , known as the integrability condition or consistency condition for the system (8) [4,8,9,38]. If it is violated, then Eq. (8) cannot be expected to be simultaneously solvable, except perhaps for very special initial data. An initial datum means here the restriction of Φ to space-time configurations on the hyperplane {x 0 = 0} in space-time M or, in other words, the values of Φ when all time variables are set to 0. To prove consistency of (8) amounts to proving the existence and uniqueness of the solution for a sufficiently large set of initial data (such as a dense subset of Hilbert space). The commutator condition (9) is closely related to the integrability condition of the Tomonaga-Schwinger equation, which in turn is closely related to the causality axiom of the Wightman axioms [44, p. 65].
If the particle number is variable, then the wave function Φ on S consists (like a vector in Fock space) of sectors with different particle numbers, and thus with different numbers of time variables. Correspondingly, the terms H k Φ on the right-hand side of (8) can involve other sectors, and the consistency question becomes more involved. As a consequence, it is not obvious that the "consistency condition" (9) is actually sufficient for consistency; and even less obvious if the solution Φ is not required to exist on all space-time configurations but only on the set S of spacelike configurations. Thus, a consistency proof requires much more than just checking a commutator condition such as (9).

Model and Goals
Our model, adopted from [39] and inspired by models of Lee [19], Schweber [46,Sec. 12a], and Nelson [34], is a simple QFT in 1+3 space-time dimensions in which one species of particles, called x-particles in the following, can emit and absorb particles of another kind, called y-particles; in short, We take the x-particles to be fermions and the ys to be bosons, and we take both to be massive and to have spin 1 2 (although in nature bosons have integer spin). The free Hamiltonians are Dirac operators. A multi-time formulation of this model was given in [39], however in a UV divergent form. In order to enable a rigorous treatment, we introduce a UV cut-off; correspondingly, we need to slightly modify the definition of the domain S of Φ (i.e., to replace it by the set S δ of δ-spacelike configurations, see below). We formulate the multitime Schrödinger equations of our model with cut-off (Sect. 2) and prove the existence and uniqueness of a solution for every initial condition of sufficient regularity (Theorem 1 in Sect. 3). In particular, this result proves consistency of the multi-time equations; some steps toward a consistency proof of this kind were already taken in [37].
It would be desirable to treat more realistic models in the future. Here, the model is chosen for allowing a rigorous proof of the multi-time consistency for a variable number of time variables; in fact, our proof is the first such proof for a model in 1+3 dimensions.
The solutions we will construct are strong solutions, i.e., smooth functions that have derivatives in the classical sense. As a by-product of our analysis, we also prove smoothness of the single-time wave function Ψ as a function of t, the x k , and the y coordinates (Lemma 3 in Sect. 3.4); as far as we know, this fact was not in the literature before for any similar QFT model. Furthermore, we prove a statement (Lemma 2 in Sect. 3.4) on how fast the support of Ψ in R 3 can grow with time: It can only grow at the speed of light, except for an additional instantaneous growth by the cut-off length δ. A corresponding statement holds true for the multi-time wave function Φ: Its evolution is propagation local up to δ; it is also interaction local up to 2δ.

Motivation
Although in QFT one often focuses on the operators, there are of course also quantum state vectors. And quantum state vectors are something real, as various considerations such as the Pusey-Barrett-Rudolph theorem [42] indicate. By that it is meant that quantum states are things in nature and in that sense physical; they are not like probability distributions (representing the observer's knowledge). In view of that, and since the usual Fock space vectors refer to a chosen Lorentz frame, we are led to the question, how can the quantum state be represented in a covariant way? This is a conceptual, physically relevant question that naturally leads us to considering multi-time wave functions, as they provide a covariant particle-position representation of the quantum state.
Most approaches to QFT focus on the operators and put the link with space-time points into the field operators, while the quantum state is a functional on the algebra or an element of an abstract Hilbert space that is not directly related to space-time. In contrast, the approach we follow here puts a link with space-time into the quantum state. Indeed, just as in non-relativistic quantum mechanics, where the quantum state can be represented by a wave function depending on many space-time points, we describe the state of a quantum field theory as a function Φ of space-time points.
Multi-time wave functions thus offer a reformulation of QFTs and an alternative perspective on them. Richard Feynman [14] wrote in 1965: "I, therefore, think that a good theoretical physicist today might find it useful to have a wide range of physical viewpoints and mathematical expressions of the same theory (for example, of quantum electrodynamics) available to him." We think this is still true today.
Our perspective on QFT provides a kind of "pedestrian" approach that expresses a QFT as a system of covariant PDEs. And indeed, the perspective of the particle-position representation of the quantum state has proven useful, for example a few years ago through the realization that certain UV-divergent particle creation terms can be tamed by means of certain boundary conditions on the wave function, so-called interior-boundary conditions [49]. Another advantage of multi-time wave functions is that they are rather intuitive objects and governed by remarkably simple and natural equations, see (11). Moreover, this approach can be used as a derivation of the mathematical structure of QFT: Instead of starting from (say) the Wightman axioms for the field operators, we may start from non-relativistic quantum mechanics, ask for including particle creation and annihilation while observing Lorentz invariance, and are then led [29] to multi-time wave functions and evolution equations such as (11), which are presumably ultimately equivalent to operator-centered approaches to QFT.
Many physicists are worried that a particle-position representation may not be possible. Well, it is certainly of interest to explore to which extent it is possible, and the situation is often better than it may seem. For example, it is sometimes said that photons cannot be localized, but that boils down to difficulties with specifying the probability density in position space, not with specifying the wave function in position space [3]; a photon wave function is mathematically equivalent to a (complexified) Maxwell field and thus unproblematical.
Another advantage of multi-time wave functions is that they, like Tomonaga-Schwinger wave functions, provide a quantum state on every arbitrary (curved) Cauchy surface (at least in the absence of UV cut-offs), something not available through a representation of the Poincaré group as provided by the Wightman axioms. An advantage over Tomonaga-Schwinger wave functions is that multi-time wave functions can be considered with UV cut-offs. Finally, multi-time evolution equations may also be expected to generalize rather easily to curved space-time, although this question has not yet been fully explored in the literature. Now the possibility of multi-time wave functions is crucially linked to a mathematical problem: that of consistency of the multi-time evolution equations. And that is what this paper contributes to.
In this work, Lorentz invariance has a delicate status: On the one hand, Lorentz invariance is a main motivation for using multi-time wave functions; on the other, the wish for a rigorous discussion leads us to introducing an UV cut-off, which breaks the Lorentz invariance. Then again, its violation is limited to certain aspects, while other aspects of Lorentz invariance are retained. For example, although the propagation of wave functions is not confined to the light cone, it is confined to a neighborhood of the light cone with spatial radius given by the cut-off length δ > 0 (see Fig. 5 and Eq. (47) for more precise elucidation). It would be desirable to work without UV cut-off, but presently that is not an option in 1 + 3 dimensions.

Comparison to Prior Works
Our work is similar to recent work in [25,36]. In [25], a rigorous consistency proof is presented for a QFT model in 1+1 space-time dimensions involving one species ("x") of spin-1 2 particles with creation and annihilation according to x x + x instead of (10); instead of a UV cut-off, that model uses interiorboundary conditions. In [36], a rigorous consistency proof is presented for the QFT model of Dirac, Fock, and Podolsky [9] in 1+3 space-time dimensions with UV cut-off, in which x-particles interact through a quantized y-field; the main difference to our model is that here we give a separate time variable to every y-particle. Further consistency proofs for multi-time equations are contained in [38], where it was shown for any fixed number N of particles that the consistency condition (9) is necessary and sufficient for consistency of (8), provided that Φ is defined on all of M N (instead of S (N ) ) and that the partial Hamiltonians are bounded operators on L 2 (R 3N , C k ) that depend smoothly on t 1 , . . . , t N , or else are (possibly unbounded) self-adjoint operators and do not depend on t 1 , . . . , t N . Rigorous consistency proofs for multi-time equations governing N particles with zero-range interaction in 1+1 space-time dimensions were given in [17,20,22,24].
An overview of work on multi-time wave functions can be found in [27,29]. The idea of a multi-time wave function was conceived as early as 1929 [12,15,33]. Early examples of consistent multi-time evolutions with interaction were given in [5,10,11,51]. The approach is contrasted with the idea of multiple timelike dimensions in [28] and with multi-time equations in classical mechanics in [40]. The appropriate version of the Born rule for Φ is formulated and proved in [30] after pioneering work in [4,33].
A big difference between the mathematics of multi-time and single-time evolution concerns the role of Hilbert spaces. In the case of a fixed number N of particles, we cannot simply fix the values of all time variables t 1 , . . . , t N and consider the wave function as a function of the spatial variables x 1 , . . . , x N alone because for some values of the x j , ((t 1 , x 1 ), . . . , (t N , x N )) will not be a spacelike configuration and Φ will not be defined there. Nevertheless, Hilbert spaces play a role as a tool in our proof, as we will consider subsets of particles at equal time values. We will also make use of the fact that wave functions do not propagate faster than at the speed of light, so that Ψ t , when considered only in a given region of 3-space, is determined by initial data Ψ 0 in a suitably larger region. This paper is organized as follows: Section 2 is dedicated to the definition of the QFT model and the formulation of its multi-time Schrödinger equations. In Sect. 3, we describe our main results. In Sect. 4, we establish lemmas about the single-time evolution that we will need as tools for the multi-time analysis. In Sect. 5, we then construct the solution to the multi-time equations and prove uniqueness of the solution. In Sect. 6, we conclude.

The Physical Model
Our model is a toy QFT in which a fixed number M ∈ N of x-particles can emit and absorb y-particles. It is a UV-regularized version of the multi-time model of [39].

Original, UV Divergent Equations
The multi-time model of [39] is defined by the formal evolution equations: is a spacelike space-time configuration, H free x k is the free Dirac operator acting on particle x k , g ∈ C 4 is a fixed spinor playing the role of a coupling constant, s ∈ {1, 2, 3, 4} is the spin index of particle y , while most spin indices are not explicitly written, means omission, \ means to remove an entry, and G 0 is a Green function: It is the C 4 -valued distribution on R 4 that is the solution of with initial condition G 0,s (0, y) = g s δ 3 (y) .
Complete definitions tailored to our model are given in Sects. 2.2 and 2.3. (The superscripts (N ), (N + 1), etc., are actually unnecessary because the sector is determined by the argument of Φ; they are mentioned merely for easier readability.) In contrast to (8), we here encounter a q-dependent (and unbounded) number of PDEs, since the number N of y-particles is unbounded, and a separate PDE is assigned to each particle. In [39], non-rigorous arguments were given for the consistency of (11).
It follows from (11) that the single-time wave function Ψ evolves with the (UV-divergent) Hamiltonian at the spatial configuration q 3 := (x 3M , y 3N ) = (x 1 , . . . , x M , y 1 , . . . , y N ). Put in a different notation, where we have assumed a fermionic Fock space for the x-particles (although H will map its M -particle sector to itself). We have written dΓ(S) for the second quantization of the 1-particle operator S, H free y for the 1-particle Dirac operator, and a x and a y for the annihilation operators for x-and y-particles in the position representation.
In this paper, we formulate our proofs for the slightly more general possibility that the coupling coefficients g also act on the spin index r of the emitting or absorbing x-particle according to with the consequence that the "coupling constant" g is an element of C 4 ⊗ C 4 ⊗ C 4 . That is analogous to quantum electrodynamics, where the coupling coefficients are proportional to the 4-vector of Dirac gamma matrices γ μ with their 4 × 4 × 4 entries γ μ rr . The UV divergence of H manifests itself in (14) in the fact that δ 3 is not an L 2 function and thus HΨ / ∈ H ; the UV cut-off will consist in replacing δ 3 by an L 2 function ϕ. We now turn to defining mathematically the cutoff version, our model in this paper, first in the single-time formulation, then multi-time.

Single-time Formulation with Cut-Off
Although our goal is a multi-time formulation, we will use the single-time formulation as a tool, along with the Hilbert space and the Hamiltonian operator.
To begin with, the wave function of a single Dirac particle can be described in two different ways: Either we assign a C 4 -vector to each point of R 3 or we consider a C-valued function on R 3 × {1, 2, 3, 4} (i.e., 4 disjoint copies of R 3 , see Fig. 1): The corresponding spin-configuration spaces of M x-particles and N yparticles are In our model, M is fixed and N is variable, so An element of Q s3 , i.e., a spin-configuration in 3D, will be denoted by q s3 . In the analogous way, we define the "spin-free" configuration spaces Q 3 x , Q 3 y and Q 3 (where the superscript means an index, not a Cartesian power).
An equal-time quantum state is an element of the Hilbert space where M ∈ N is the (fixed) number of x-particles, ⊗ A denotes the antisymmetric and ⊗ S the symmetric tensor power; in particular, H y is the bosonic The free time evolution is given by the Dirac operators where α 1 , α 2 , α 3 and β are the Dirac alpha and beta matrices. Fermions are assumed to have a rest mass m x > 0, and bosons have a rest mass m y > 0. (Domains of operators will be specified in Sect. 3. Henceforth, we use bold face font for operators on Hilbert space, as well as still for vectors in dimensions 3, 3M, 3N, M, N .) Let e 1 , . . . , e 4 denote the standard basis in C 4 . The annihilation operator of a y-particle in the position representation at location x in the spin state e s will be denoted by a(e s δ 3 (· − x)). To implement the UV cut-off, we will replace δ 3 in the Hamiltonian by a smooth (i.e., infinitely often differentiable) cut-off function with compact support inside a ball B δ (0) of radius δ > 0 around the origin 0 ∈ R 3 . Physically, the radius δ should be small, but our mathematical results apply to any positive δ. The corresponding smeared-out annihilation operator a(e s ϕ(· − x)) will be abbreviated as a s (x). For the location x ∈ R 3 , we will need to insert the position operator of (say) the k-th x-particle, and the resulting annihilation operator will be denoted a s (x op k ). That is, The interaction Hamiltonian is of the form: where the k-th term represents emission and absorption by particle x k and is defined by The H free

Multi-time Formulation with Cut-Off
We now write down the system of multi-time equations corresponding to the Hamiltonian (27); the equations are a version of (11) with UV cut-off ϕ. As a preparation, we first need to • Define the set S δ of multi-time configurations q 4 = (x 4M , y 4N ) for which Φ(q 4 ) will be defined (Sect. 2.3.1); • Describe which functions on S δ are regarded as smooth, and which directional derivatives can be taken of them (Sect. 2.3.2).
Afterward, we will formulate the multi-time equations in Sect. 2.3.3.

Admissible Configurations.
We now define the set S δ of δ-spacelike configurations. A similar modification of S was used in Sections 1.5.3 and 4 of [38] for the purpose of using interaction potentials of range less than δ.
We write Q 4 for the set of all space-time configurations (4D configurations) q 4 = (x 4M , y 4N ), Sometimes, we want to talk about spin-configurations in 4D, (29) As a consequence of the cut-off ϕ over the distance δ > 0, an x-particle can instantaneously interact with (i.e., emit or absorb) a y-particle at a distance up to δ, and two x-particles can instantaneously interact with each other (if one emits a y and the other absorbs it) at a distance up to 2δ. Since this instantaneous interaction will generically conflict with the consistency of the multi-time evolution, we make the domain S of Φ slightly smaller and allow configurations involving two particles with such a small distance only if the time coordinates of the two particles are equal. So on each sector N , we define the set S (N ) δ of δ-spacelike configurations as follows: The set S δ is illustrated in Fig. 2. Whenever two particles satisfy the appropriate inequality in (30), we say that they keep their safety distance. The sets S , by the same conditions as in (30)] for spin-configurations in 4D.
A multi-time wave function is a mapping Note that the definition of S δ and S s δ depends on the frame of reference, as Lorentz invariance is broken by the cut-off ϕ.

Admissible Wave Functions.
Our considerations focus on smooth wave functions Φ; for them, derivatives can be understood in the classical sense. However, since S . Put briefly, we regard a function as smooth at q 4 ∈ S δ if it is smooth in the local number of dimensions of S δ at q 4 .
We now approach the detailed definition of smoothness, following [38]. To begin with, a function Ψ : Q s3 → C will be called smooth, Ψ ∈ C ∞ (Q s3 ), if its restriction to each sector is smooth. For any q 4 ∈ S δ , if two particles do not keep their safety distance (say, x j − x k < 2δ), then, by the definition (30), their times must be equal (x 0 j = x 0 k ) and must remain so in the vicinity of q 4 . Hence, we cannot vary the two time coordinates independently of each other, we can only increase both by the same amount. And hence, we cannot form the partial x 0 k derivative of Φ, we can only form the directional derivative in the direction e 0 xj + e 0 x k in (R 4 ) M +N (unless further particles do not keep their safety distance from either x j or x k ); here, e μ means the standard basis of R 4 ; see Fig. 3. That is why we proceed by grouping particles into families of equal time coordinate. To this end, let L (N ) be the M + N -element set of particle labels; that is, where x k and y are now not space-time points but names for the particles.
Let P (N ) denote the set of partitions of L (N ) , and let, as a unified notation, z μ x k = x μ k and z μ y = y μ , so z could represent any particle (x or y). For every P ∈ P (N ) , we define where we abbreviated the "safety distance" by In words, S P δ contains those configurations for which particles in the same family have equal time coordinate and particles in different families keep their ; that for different choices of P , the S P δ are not necessarily disjoint (because it is allowed that particles in different families have equal time coordinate); and that together, they cover S , we can group particles into the same family whenever they do not keep their safety distance; that is, P consists of the equivalence classes of the reflexive and transitive hull of the relation "do not keep their safety distance"; this is the finest partition with q 4 ∈ S P δ ). If q 4 ∈ S P δ , then we also write t j for the joint time variable of all particles in P j , and q j for the list of space coordinates of all particles belonging to P j ; using that notation, we also write so that S P δ can also be regarded as an open subset of R 3M +3N +J (while we continue to use the notation q 4 = (x 1 , . . . , x M , y 1 , . . . , y N )).
for every N ∈ N 0 . The set of all smooth complex-valued functions on S s δ is denoted by C ∞ (S s δ ), the set of smooth functions S P δ → C k by C ∞ (S P δ , C k ). Put differently, S P δ is a submanifold of M M +N of dimension 3M + 3N + J, and Φ counts as smooth on S P δ if it is smooth as a function on this submanifold. It then also follows that at q 4 ∈ S P δ , derivatives of Φ in any direction tangent to the submanifold can be taken, as well as higher-order derivatives. That is, Φ can be differentiated relative to any position coordinate and relative to any family time coordinate t j .

Multi-time Equations.
We now formulate the analog of the multi-time equations (11) with UV cut-off.
Here, we have to face an issue that we touched upon already in Sect. 2.3.2 and Fig. 3: At configurations in which x j and x k have equal time coordinate but spatial distance < δ, we cannot vary x 0 j and x 0 k independently without leaving the set S δ on which Φ is defined. As a consequence, the partial time derivatives such as ∂ x 0 k that were specified in (11) (the original equations using Dirac delta functions) do not make immediate sense with UV cut-off. However, suitable linear combinations of these equations will specify derivatives such as ∂ x 0 j +∂ x 0 k , which can be understood as the derivative in the direction e 0 xj +e 0 x k , which does make sense. Correspondingly, we take these combinations as the precise evolution equations, as the way how (11) should be interpreted. Again, we group together particles that do not keep their safety distance (and thus have equal time coordinate) into families with common time t j and formulate equations for the ∂ tj -derivative; in particular, at q 4 ∈ S P δ with P = {P 1 , . . . , P J }, we specify J equations rather than M + N . Concretely, the evolution equations can be expressed as follows: For any N ∈ N 0 , P ∈ P (N ) , and q 4 ∈ S P δ , Here, x k ∈Pj means the sum over all x-particles in P j , etc.; the index r k is made explicit at some terms and not others, but occurs at each term. Note that the arguments of Φ (N +1) and Φ (N −1) lie in S δ again. Thus, the right-hand side of (37a) defines, for given N, P, j, an operator Note that the number of time variables t j involved in the system of Eq. (37a) is J = #P . Since q 4 can lie in S P δ for several P (say P and P ), one needs to check that Eq. (37a) from P and P are compatible with each other. Indeed, if P is a refinement of P , then the equations from P are linear combinations of the equations from P . Now for any given q 4 there is a coarsest partition (particles with equal time coordinate belong to the same family) and a finest partition (only particles that do not keep the safety distance belong to the same family), so the equations from either P or P are just linear combinations of the equations obtained from the finest partition. Thus, if Eq. (37a) hold for the finest partition, then they hold for every partition P such that q 4 ∈ S P δ . Equation (37a) is the system of equations that we require to hold for every N, P, j, q 4 ∈ S P δ . It is the system defining the multi-time evolution of Φ. The initial value problem amounts to solving (37a) for a given initial datum

Results
The main statement of this paper, Theorem 1 in Sect. 3.3, asserts that the system of multi-time equations (37a) has a unique solution for every sufficiently regular initial datum, and that the solution has the expected properties.

Supports
We also need to talk about the support of a wave function; in particular, it will be relevant to exploit the advantages of a compact support. However, we cannot expect Ψ to have compact support in Q s3 as that would imply a concentration of Ψ on finitely many sectors. We will thus additionally define the 3-support of Ψ, i.e., the region G ⊆ R 3 where particles can be encountered at all. We use again the notation q 3 := (x 3M , y 3N ).

Definition 2.
For Ψ ∈ H or Ψ : Q s3 → C let supp Ψ denote the essential support of Ψ in Q 3 , i.e., the smallest closed set G ⊆ Q 3 such that Ψ = 0 almost everywhere outside G. We set where the overbar means the closure in R 3 . Equivalently, supp 3x Ψ is the closed union of the projections of supp Ψ to each x k -variable in R 3 (and similarly for supp 3y Ψ).
The following definition will be convenient for expressing propagation locality.

Definition 3.
We define the grown set of G ⊆ R 3 as for any t ≥ 0, where B r (x) means the closed ball of radius r around x.

Admissible Initial Data
For every N, n ∈ N 0 , let M (N, n) be the set of multi-indices be the corresponding derivative.
Note that since the r-th Sobolev norm of Ψ (N ) is given by condition 3 is equivalent to the condition that for every m, r ∈ N 0 , the sum over N of N m times the square of the r-th Sobolev norm of Ψ (N ) is finite. Note also that H ∞ c is a dense subspace of H ; for example, it contains the dense subspace of all smooth functions Q s3 → C with compact support in Q s3 (in particular, which vanish outside finitely many sectors) that satisfy the fermionic and bosonic permutation symmetry. An alternative characterization of H ∞ c is given in Lemma 8 in Sect. 4.1.

Theorem About Multi-time Evolution
The following summability property of the solutions Φ : S s δ → C will be relevant: For every J ∈ N and every choice of Here, is the set of all δ-spacelike configurations with N bosons in which only the times t 1 , . . . , t J occur; it has dimension 3M + 3N . We are now ready to formulate our main result. Theorem 1. For any initial datum Ψ 0 ∈ H ∞ c , the multi-time initial value problem (37a), (39) has a solution Φ ∈ C ∞ (S s δ ) satisfying the summability property (45). The solution is unique among the functions in C ∞ (S s δ ) satisfying (45) and has the following further properties: for any x k or y in q 4 .
Remark. 1. 4-supports The condition (47) for configurations outside the support of Φ can equivalently be re-expressed as: where Infl(A) means the domain of influence of the space-time set A ⊆ M , where future(A) is the union over x ∈ A of the closed future light cone of x, and correspondingly past(A). We can also express the situation in terms of the 4-support of Φ: Then, These sets are shown in Fig. 5. 2. Asymmetry Between x and y Patently, (47) is not symmetric in x and y. The reason is that the Hamiltonian H makes all x-particles create and annihilate y-particles, so even if at t = 0 there are no y-particles near an x-particle, this situation will change at t > 0 and supp 3y Ψ will immediately be extended by Gr(supp 3x Ψ 0 , δ). The same holds for t < 0. However, if there is a separated y-particle, H will not cause a creation of x-particles near it. So supp 3x Ψ does not grow by an additional amount depending on supp 3y Ψ. 3. Uniqueness The summability property (45) is relevant to the uniqueness of the solution; in C ∞ (S s δ ), the solution would not be unique, essentially because a further contribution to the wave function could come "from infinity on the N axis." This effect has nothing to do with the multiple time variables, it arises from considering smooth functions that solve the equations in the classical sense and is, in the usual one-time formulation, taken care of by choosing as the Hamiltonian a self-adjoint extension of the differential expression for H. This point is elucidated further in Sect. 3.5. 4. L 2 Solutions Usually, when considering solutions of the 1-time Schrödinger equation, one does not require classical differentiability. One allows weak (distributional) derivatives and thus wave functions from Sobolev spaces, and even more, one allows arbitrary L 2 functions as initial data, as the unitary time evolution operator exp(−iHt) is defined on the entire Hilbert space. In the same way, we can define a notion of L 2 solution of our multi-time equation (37a), and then allow for arbitrary Ψ 0 ∈ H . In fact, our construction of the solution Φ works in this way, and we have to invest further work for proving that for nice initial data Ψ 0 ∈ H ∞ c , the function Φ we construct is a classical solution. 5. Domain of Consistency As it happens, the simple model we are considering is consistent on an even larger set than S δ . This fact is presumably a curious artifact of our simple model that will not extend to more realistic models. This larger set S δ ⊃ S δ is defined by the same condition as for S δ between x-particles but no restriction on the y-particles, not even being spacelike separated, neither from x-particles nor among themselves. That is, and The physical reason is that the bosons do not interact with each other, nor with the fermions, except for being emitted and absorbed by the fermions. It is well known [29]  and S s δ ; in the same way as for S s δ , one defines C ∞ ( S s δ ). In Sect. 5.4, we show that for any initial datum Ψ 0 ∈ H ∞ c , the system of multi-time equations (56) below with initial values Ψ 0 as in (39) has a unique solution Φ ∈ C ∞ ( S s δ ) satisfying the summability condition (45). The multi-time equations are essentially given by (11) with the Green function G 0 replaced by a family (indexed by r, r ) of cut-off Green functions G rr , each in C ∞ (R 4 , C 4 ), viz., the solution of with initial condition G rr s (0, y) = g rr s ϕ(y) .
More precisely, since, as discussed in Sect. 2.3.3 for S δ , the derivative ∂ x 0 k is not admissible within S δ when two x-particles are closer than 2δ, we need to group together x-particles that do not keep their safety distance into a family with a common time. That is, for a given q 4 ∈ S δ let {P 1 , . . . , P J } be the finest partition of {x 1 , . . . , x M } such that particles that do not keep their safety distance belong to the same partition set. By the definition (52) of S δ , all particles in P j have the same time coordinate, now called t j . The multi-time equations (37a) can then be written as We note that every solution Φ ∈ C ∞ ( S s δ ) of (56), when restricted to S s δ , also solves (37a). For verifying this, the only real issue is that (56a) involves a sum over all y-particles ( runs from 1 to N ), whereas (37a) involves a sum over only the y-particles in the family P j of x k . However, this difference has no consequences in S δ , as there every y-particle outside the family of x k is δ-spacelike from x k , while G vanishes on all points δspacelike from x k .

Lemmas About Single-time Evolution
A few statements about the single-time evolution that we need as tools for proving Theorem 1 may be worth mentioning in their own right. (27) is essentially self-adjoint on H ∞ c ⊂ H .

Lemma 1. (Self-adjointness) The single-time Hamiltonian H defined in
By Lemma 1, H has a unique self-adjoint extension, which we will also denote by H and which provides the unitary single-time evolution operator for all t ∈ R. Lemma 2 yields information about the growth of 3-supports, and Lemma 3 the regularity of Ψ t .

Lemma 2.
(Support growth) Let Ψ 0 ∈ H and Ψ t = U (t)Ψ 0 for all t ∈ R. Supports will grow under time evolution at most as follows: in particular Remark. 6. Local Evolution Operator As a consequence, the wave function Ψ t on configurations with all particles in a region G ⊂ R 3 depends only on the initial data in Gr(G, |t|+δ). Thus, the single-time evolution defines an operator where H (G) means the Hilbert space associated with the region G, i.e., the subspace of L 2 (Q s3 (G)) with the appropriate permutation symmetry, where is the set of spin-configurations concentrated in G.
The following lemma ensures that the solution of the single-time evolution is smooth also in time; in fact, Ψ(t, q 3 ) will be smooth as a function of t and q 3 . q and t) For every Ψ 0 ∈ H ∞ c , one can choose a representative of each Ψ t in such a way that Ψ t (q s3 ) is a smooth function of t and q s3 .

Lemma 4. (Smoothness in
If, moreover, Ψ 0 depends smoothly on a parameter λ, then the solution Ψ(t, q 3 ) will also depend smoothly on λ. To formulate this statement precisely, let d ∈ N, and let H d denote the subspace of L 2 (R d × Q s3 ) of functions Ψ(λ, x 3M , r, y 3N , s) that are anti-symmetric against permutations of (x 3M , r) and symmetric against permutations of ( y 3N , s).
where M (N, n, d) is the set of multi-indices of degree n for the variables Then, Ψ is smooth in t, λ, q 3 , Ψ ∈ C ∞ (R × R d × Q s3 ).

Strong Single-time Solutions
Since we have chosen to consider strong solutions (smooth functions with classical derivatives) of our multi-time PDEs, and since we will apply singletime considerations to each time variable, we also need to consider strong solutions of our single-time evolution equation, which reads explicitly, with q s3 = (r, x 3M , s, y 3N ), t, x 3M , y 3N ,ỹ While Lemmas 1, 3, and 4 together guarantee that the function Ψ obtained by applying the exponential e −iH t of the self-adjoint operator H to Ψ 0 is a strong solution of (64), 1 we also need to verify that (64) has no further classical solutions with initial value Ψ 0 ; this is done in Lemma 7. We will also use that our bounds on support growth follow directly from (64), a fact established in Lemma 6.
The differences between Lemma 2 and Lemma 6 arise from their different methods of proof: Lemma 2 uses the Trotter product formula, Lemma 6 integration of the current; Lemma 2 concerns L 2 initial data and a time evolution defined by a self-adjoint operator H, Lemma 6 concerns a C ∞ solution of (64).

Lemma 7. (Uniqueness of strong single-time solutions) For every
Without demanding the summability properties Ψ t < ∞ and N 1/2 Ψ t < ∞, the solution would not be unique in C ∞ (R t × Q s3 ), with further contributions coming from N → ∞. This point can be illustrated by the following minimalistic variant of (64), in which spin is dropped and space R 3 is replaced by a single point, so that the configuration is fully described by the number of y-particles and no space variables appear any more, Ψ = (Ψ (0) , Ψ (1) , . . .), while Ψ (N ) (t) ∈ C is governed by the Schrödinger equation Then, an example of a smooth nonzero solution with initial condition Ψ 0 = (0, 0, . . .) is provided by with polynomials P N ,P N defined recursively by (the prime means the ν-derivative) with starting values P −1 = 0, P 0 = 1,P −1 = 0,P 0 = 1.
1 That is because Ψt ∈ H ∞ c by Lemma 3 and Ψt ∈ dom(H ) by Lemma 1, so t → Ψt is a differentiable curve in H , whose derivative agrees, on the one hand, with the strong t-derivative of the smooth function Ψ provided by Lemma 4 and, on the other hand, with the right-hand side of (64) by the definition of H on H ∞ c .

Proofs: Single-time Evolution
In the following, we will repeatedly use that for elements a i of a normed space.

Essential Self-adjointness
Proof of Lemma 1. We first verify that H is well defined as an operator is smooth and has compact support; k H free x k Ψ (N ) is anti-symmetric against x-permutations; and k H free because the α and β matrices are unitary. The first term is finite by (43) with m = 0 and n = 1, the second because Ψ ∈ H . In the same way, one can show that dΓ y (H free y )Ψ ∈ H , except that one has to apply (43) once with m = 1, n = 1 and once with m = 2, n = 0. Now we show that H int Ψ ∈ H . Since ϕ is smooth with compact support, H int Ψ is smooth and has compact 3-support, so each sector of it is square-integrable; it is clear that H int Ψ has the right permutation symmetry; it remains to show that N Ψ (N ) 2 N < ∞ forΨ = H int Ψ. We show this for Ψ = a(x op k )Ψ andΨ = a † (x op k )Ψ. In cases where the operators are annihilating or creating a fixed function f ∈ L 2 (R 3 ), we may bound them against the y-number operator as follows (e.g., [34]): In our case, the function f (y) = ϕ(y − x k ) depends on the position x k of the fermion, so f ∈ L 2 , but it is not fixed. However, we may write the Hilbert space H as a direct integral over all fermion configurations and decompose a(x op k ) and a † (x op k ) into fibers: (78) Applying the bound (76) on each fiber, we get the estimates The first term is finite by (43) with m = 1, n = 0, the second by Ψ ∈ H . This completes the proof that HΨ is well defined and ∈ H for Ψ ∈ H ∞ c . Now we turn to essential self-adjointness. Following [36] and [2], we use the following "commutator theorem" [

Then, H is essentially self-adjoint on D.
In our case, D = H ∞ c ; it is straightforward to verify that H is symmetric on H ∞ c . We choose as the comparison operator We first prove that R is essentially self-adjoint on H ∞ c , using standard arguments: To begin with, it is clearly defined on H ∞ c . The negative Laplacian −Δ is known to be positive and essentially self-adjoint on C ∞ c (R 3 ) (compactly supported smooth functions). Adding a bounded self-adjoint operator such as 1 has no effect on essential self-adjointness. It is known further [43, Thm. VIII.33 and Ex. 2 p. 302] that R (and thus also R) is essentially selfadjoint on C ∞ c (Q s3 ) (smooth functions that vanish on all but finitely many sectors and have compact support in each sector). Since the latter space is contained in H ∞ c and dense in H , R is also essentially self-adjoint on H ∞ c .
Positivity follows from positivity on each sector (using that R maps each sector to itself); in fact, R ≥ 1 because M ≥ 1. It remains to prove the two inequalities (81) and (82). For (81), we need to bound H free x k , H free y , and H int x k against R. The subsequent paragraphs deal with these terms separately.
H free x k : In the following, we focus on a single sector N and simply write Ψ (N ) = Ψ and · N = · . The free Dirac operator is rather easy to bound against the Laplacian: First, The left term can be bounded by −Δ k , using that the anti-commutator of the α-matrices is {α a , α b } = 2δ ab : By considering the Fourier transform of Ψ and using that |k| ≤ 1 + |k| 2 ≤ 1 + |k| 2 , we obtain that Hence, in every sector and thus also in H . H free y : The same reasoning shows that H int x k : is a finite sum of annihilation and creation operators a s (x op k ) and a † s (x op k ). They can be bounded using (79): [H int x k , −Δ k ]: vanishes unless k = k because H int x k does not involve x k . We drop the fermion index k, write −Δ = 3 a=1 i∂ a i∂ a and shift one derivative to the other side of the scalar product: What happens if a derivative i∂ a hits H int Ψ? The coordinate x k appears in each of a(x op k )Ψ (N ) and a † (x op k )Ψ (N ) twice, cf. (24): once in ϕ(y − x k ) and once in the factor Ψ (N +1) or Ψ (N −1) . Therefore, by the product rule, i∂ a H int Ψ will contain 2 terms: In the first term, ϕ is replaced by i∂ a ϕ; the annihilation and creation operators modified in this way will be denoted by a i∂aϕ (x op ) and a † i∂aϕ (x op ). The second term is just i∂ a directly hitting Ψ before creation or annihilation takes place. That is, This allows for further treatment of the commutator: with c 3 ∈ R + . By symmetry of (a+a † ), the second scalar product is real, so its imaginary part is 0. Using |Im z| ≤ |z|, the Cauchy-Schwarz (CS) inequality, and (69), we obtain that In the same way as in (86) but using |k| ≤ 1 + |k| 2 , we find that Hence, as desired.
[H int x k , −Δ ]: can be evaluated using the commutation relations which are perhaps most easily verified through a few lines of calculation by applying dΓ y (−Δ y ) to the defining equations (24) of a and a † ; for the first relation, integrate by parts twice inỹ; for the second, note that Ψ (N −1) does not depend on y . From (96), we obtain that with c 4 , d 2 ∈ R + . This concludes the proof of (82) and thus of Lemma 1.
In other words, the first and the third condition in the definition of H ∞ c can equivalently be replaced by Ψ ∈ dom(R p ) ∀p.
Proof. In terms of the (sectorwise) Fourier transform Ψ of Ψ, the third condition (43) can be expressed as: for all m, n ∈ N 0 , where P Nmn is the polynomial with d = 3M + 3N . Likewise, the condition Ψ ∈ dom(R p ) can be expressed as: where Q Np is the polynomial Every Q Np is bounded from above by a linear combination with N -independent coefficients of a finite number of P Nmn 's with the same N , m ≤ 2p, and n ≤ 2p; thus, H ∞ c ⊆ dom(R p ). Conversely, every P Nmn is bounded from above by a linear combination with N -independent coefficients of a finite number of Q Np 's with the same N and p ≤ (m + n + 1)/2. Suppose that Ψ ∈ dom(R p ) for all p; then, the weak derivatives of Ψ satisfy (43). Therefore, each Ψ (N ) lies in the Sobolev space +N ) ) for every n; by the Sobolev embedding theorem [1], Ψ (N ) possesses a smooth representative and satisfies (43) in the strong sense.

Support Growth
Proof of Lemma 2. For any number M + N of free Dirac particles, it is well known [7,48], [6,Thm. 2.20], [38,Lemma 14] 2 that the wave function propagates no faster than light, in particular that In fact, if we took the Hamiltonian of the M + N particles to be only without the H free y , then the 3y-support would be invariant, supp 3y Ψ t = supp 3y Ψ 0 , while supp 3x Ψ would grow according to (102). Conversely, if we took only H free y as the Hamiltonian, then the 3x-support would be invariant, while the 3y-support grows according to (103).
Let us turn to the case with particle creation. Following a strategy of [48,Theorem 3.4], we now prove that supp 3x Ψ t is not altered if we change the free evolution by including H int .
First, we consider the evolution with only H int x k as the Hamiltonian, without the free Dirac operators, and claim that then the 3x-support is invariant, Indeed, since H int x k can be decomposed into fibers as in (78) (the continuous analog of being block diagonal), so can exp(−iH int x k t) [45,Thm. XIII.85(c)]. Considering Ψ 0 in this fiber decomposition (i.e., as a function of x 3M with values in H y ), it follows that Ψ t = exp(−iH int x k t)Ψ 0 vanishes in those fibers where Ψ 0 does. So the 3M -x-support is invariant, and therefore also the 3xsupport.
Second, under H int Indeed, due to the fiber decomposition over x 3M just mentioned, it suffices to consider just one fiber x 3M and show that on it, for any set G ⊆ R 3 with B δ (x k ) ⊆ G and any vector Ψ 0 ∈ H y with 3y-support in G, also exp(−iH int x k t)Ψ 0 has 3y-support in G. To see this, note that for disjoint sets . Now the claim follows, and with it (106). 3 Third, back at the full time evolution, we can now decompose H into a part H free Mx which makes supp 3x Ψ grow at most at the speed of light and a sum of terms leaving supp 3x Ψ invariant, All contributions to H are now put together using Trotter's product formula [43,Thm. VIII.31]: (108) We claim that for all k = 0, . . . , n, 3 An alternative argument might go as follows but would require further work to be made rigorous. On any fiber x 3M and for any set G ⊆ R 3 with B δ (x k ) ⊆ G, H int x k commutes with 1 G . Therefore, also exp(−iH int x k t) commutes with 1 G . As a consequence, a vector Ψ 0 ∈ Hy with 3-support in G evolves to one with 3-support in G. supp 3y U k t,n Ψ 0 ⊆ Gr supp 3y Ψ 0 , k n |t| ∪ Gr supp 3x Ψ 0 , k n |t| + δ .
We proceed by induction along k. For k = 0, the claim is trivially true. For the induction step k → k + 1, we use that Gr(Gr(G, s), t) = Gr(G, s + t) and A ⊆ B ⇒ Gr(A, t) ⊆ Gr(B, t) and conclude that each factor e −i t n H free Mx makes supp 3x Ψ grow by |t|/n, while all other factors in U t,n leave it invariant. Thus, Likewise, using Gr(A ∪ B, t) = Gr(A, t) ∪ Gr(B, t), What we need is the 3y-support after U k+1 t,n , and this operator comprises in addition the interaction factors exp(−i t n H int x k ), see (108). The effect of these factors on supp 3y is that it will at most be joined with the δgrown set of supp 3x of the wave function obtained at that point, i.e., of e −i t n dΓy(H free completing the induction and thus the proof of (109) and (110). Now we set k = n and let n → ∞. If a vector Ψ t ∈ H is a limit of a sequence Ψ t,n → Ψ t (such as Ψ t,n := U n t,n Ψ 0 ) and supp 3x Ψ t,n ⊆ G for all n, then supp 3x Ψ t ⊆ G; likewise for supp 3y . Thus, we have proved (58) and therefore also (59).

Smoothness Conditions
We now prove Lemma 3. The proof is very similar to [36, Thm. 7, Lemma 8] and uses a commutator theorem due to Min-Jei Huang [16, Theorem 2.3] that we will use in the following specialized form: Let R be a strictly positive, self-adjoint operator and H a self-adjoint operator. Suppose that, for every n ∈ N 0 , is densely defined and bounded. Then, for every p ∈ N, Proof of Lemma 3. Due to Lemmas 2 and 8, it suffices to show that for the strictly positive comparison operator R given by (83) and every p ∈ N. By Huang's theorem, it suffices to show that every Z n is densely defined and bounded.
To begin with, Z n is defined on the dense subspace dom(R 2 ). Indeed, for Ψ ∈ dom(R 2 ), R −n Ψ ∈ dom(R n+2 ). Whenever Ψ ∈ dom(R p ) for p ≥ 1, then RΨ exists and lies in dom(R p−1 ), and HΨ exists by virtue of (81) and lies in dom(R p−1 ) as well. Thus, [H, R]R −n Ψ exists and lies in dom(R n ), so Z n Ψ exists. Now we show that Z n is bounded.
First, we consider n = 1 and show that Z 1  H int x k , −Δ k is evaluated using the product rule as in (91). We obtain that A similar equality holds true for a † . So, with suitable c 3 , c 4 , c 5 ∈ R + , which is what we wanted to show.
H int x k , dΓ y (−Δ y ) is bounded as in (97): So we have obtained the desired inequality (123). We now turn to arbitrary n ∈ N and show how the bounds generalize. Since R ≥ 1, clearly R −n ≤ 1, so it remains to show for a suitable d n ∈ R + . To this end, we write Bounding these terms works similarly to the n = 1 case. The N -sector of the first term applied to Ψ is a linear combination of terms of the form: where a # means either a or a † and β, γ are multi-indices with |β + γ| ≤ 2n (with β acting on 3 variables and γ on 3M + 3N ). By (79), each such term is bounded by c(N + 1) 1/2 R |γ|/2 Ψ (N ) , and they are jointly bounded by c R n Ψ . The second term in (128) is bounded by d 1 R n Ψ , as follows from (123) by inserting Ψ → R n−1 Ψ, which still lies in H ∞ c .
Proof of Lemma 4. Let T > 0 be arbitrary but fixed, and let H T be the subspace of L 2 ((−T, T ) ×Q s3 ) of functions with the appropriate fermionic and bosonic permutation symmetry. We can identify H T with L 2 ((−T, T ), H ). The function t → Ψ t (Ψ for short) belongs to H T because U (t) is unitary, in fact Ψ H T = √ 2T Ψ 0 H . We now show that Ψ ∈ H p ((−T, T ), H ) for all p ∈ N 0 , i.e., that Ψ possesses weak time derivatives of any order that are square integrable over t.
It follows further that also weak time derivatives of Ψ exist and are given by (−iH) p Ψ , and that t → ∂ p t Ψ t 2 is differentiable for every p, therefore continuous for every p, and thus integrable over any (−T, T ). The upshot is that (−iH) p Ψ ∈ H T , which shows that Ψ ∈ H p ((−T, T ), H ). Now we show that Ψ ∈ H p ((−T, T ) × Q s3 ). The weak t derivative in L 2 ((−T, T ), L 2 (Q s3 )) gets translated to the weak t derivative in L 2 ((−T, T ) × Q s3 ). Since R q Ψ 0 ∈ H ∞ c , also ∂ p t R q Ψ ∈ L 2 ((−T, T )×Q s3 ) for every q ∈ N 0 . It follows that Ψ ∈ H p ((−T, T ) × Q s3 ). By the Sobolev embedding theorem, Ψ ∈ L 2 ((−T, T )×Q s3 ) has a smooth representative, and by continuity of t → Ψ t , this smooth function is also a representative of Ψ t for every t ∈ (−T, T ). Since T was arbitrary, Ψ is smooth on R × Q s3 .
We define the probability current j, a vector field on R × Q 3 , in analogy to the Dirac current; on the n-y-particle sector, it is defined by , q 3 , r, s)  *  Ψ(t, q 3 , r, s) j x k ,a (t, q 3 ) := r ,s,r k Ψ(t, q 3 , r k ) * (γ 0 γ a ) r k ,r k Ψ(t, q 3 , r k )    where a = 1, 2, 3, and not all spin indices are always made explicit. The current vector field has the property that for each k and each , (j 0 , j x k ,1 , j x k ,2 , j x k ,3 ) and (j 0 , j y ,1 , j y ,2 , j y ,3 ) are future-causal (i.e., future-timelike or futurelightlike). It follows that j · n ≥ 0 (with · the Euclidean inner product in 3M + 3n + 1 dimensions) on ∂C lat and ∂C top .
Let C (n) := C ∩ (R × Q 3,(n) ) be the n-y-particle sector of C. By the Ostrogradsky-Gauss integral theorem (divergence theorem), The integral over ∂C (n) bot is 0 by hypothesis, those over ∂C (n) lat and ∂C (n) top are non-negative, so the left-hand side must be non-negative. We will show that It then follows that each summand must vanish, so the right-hand side of (137) vanishes, and in particular the integral over ∂C (N ) top vanishes, which is what we wanted to show. So it remains to prove (138).
In Lemmas 13 and 14 of [38], the (3M + 3n + 1-dimensional) divergence of j vanished everywhere in C. In our situation, this is not the case, due to creation and annihilation terms in the Hamiltonian. However, as we will show, the integral of the divergence of j over C still vanishes because the creation and annihilation terms transfer probability to other places in C but not outside of C. Indeed, Consider any x 3M ⊂ G := Gr X 3M , T − t + ε . From the definition of the annihilation operator a s (x), cf. (24), and the fact that B δ (x k ) ⊂ Gr X 3M , T − t + δ + ε , we obtain that Thus, by the definition (26) of H int x k in terms of a and a † , the left-hand side of (140) equals a telescopic sum with partial sum up to N 0 given by G M dx 3M Im c * N0 , whose modulus is using the Cauchy-Schwarz inequality, (69), (79), and the hypotheses Ψ t < ∞ and N 1/2 Ψ t < ∞. This shows that (140) indeed converges to zero (in fact absolutely). This concludes the proof.

Proofs: Multi-time Evolution
In this section, we prove Theorem 1; in particular, we prove existence and uniqueness of solutions of the multi-time equations (37a) with given initial data (39). (45) to the multi-time equations (37a) with initial data Ψ 0 as in (39).

Uniqueness Lemma 9. (Uniqueness of solutions) For every
Proof. For any q 4 ∈ Q 4 , let J(q 4 ) denote the number of different values of time variables that occur in q 4 ; for q 4 ∈ S δ , this is the minimal number of families. Let S δJ = {q 4 ∈ S δ : J(q 4 ) ≤ J} .
(149) We proceed by induction along J.
Everywhere on S δ1 , Φ is a function of just 1 distinct time variable, and the restriction of Φ to S δ1 obeys the single-time equation (64). By (45) for J = 1, Ψ t < ∞ and N 1/2 Ψ t < ∞, and by Lemma 7, Φ is unique on S δ1 .
The induction assumption asserts that Φ is unique on S δJ−1 , and we need to prove uniqueness on S δJ . Fix Q 4 ∈ S δJ and sort the families in increasing order of the time variables, T 1 < . . . < T J−1 < T J . As in (36), we write Q 4 = (T 1 , Q 1 , . . . , T J , Q J ) with Q j = (X j , Y j ) the configuration of all particles with time variable T j . Let φ be the function obtained from Φ on S δJ by inserting t 1 = T 1 , . . . , t J−1 = T J−1 while keeping t J variable, and inserting the configurations of the families 1 through J − 1, This set includes (T J , Q J ) because Q 4 ∈ S δ . Since Φ is a solution of (37a) for j = J, φ is a solution of (64) with time variable t = t J , M J rather than M fermions, and initial data given by Φ where t J = t J−1 . By (45), φ t < ∞ and N 1/2 φ t < ∞ with the norm · taken over the configurations in the bracket in (151). By Lemma 7, φ is uniquely fixed on its domain (151), in particular at (T J , Q J ), as claimed.
Proof. It now plays a role that we evolve negative t j only toward the past and positive ones only toward the future. So let J + (q 4 ) denote the number of positive time values, J − (q 4 ) that of negative ones, and let J 0 (q 4 ) be 1 or 0 depending on whether 0 occurs as a time value. Let S δJ+J− = {q 4 ∈ S δ : We proceed by induction, first along J + , then along J − .
On S δ10 and S δ01 , the statement is provided by Lemma 6. If true on S δ,J+−1,0 , it follows on S δ,J+,0 in the same way as in the previous proof. For J − > 0, order the negative time variables so that the least comes last. Then, the induction step from S δ,J+,J−−1 to S δ,J+,J− works in the same way as the previous proof but in the opposite time direction.

Commutator Conditions
The commutator condition (9) arises heuristically as the consistency condition. It will also play a role in our proof of the existence of solutions, specifically for proving that the function obtained by solving one of the multi-time equations also solves the others. To this end, we verify the commutator condition in this section. The appropriate condition concerns the Hamiltonians H P j corresponding to a family with a common time variable t j as in (37a).
To check it explicitly, we begin with other partial Hamiltonians, associated with individual particles and defined by H x k = H free x k + H int x k and H y = H free y . In the following, we will understand these Hamiltonians as acting on functions of q 4 instead of q 3 and then write them as H z instead of H z , z ∈ {x 1 . . . x M , y 1 . . . y N }. In fact, we will assume for H z the slightly more general form corresponding to (56) instead of (37a) (using cut-off Green functions and applicable also outside of S δ ), In order to convince oneself that the H z can be applied to smooth functions on Q s4 or on S s δ or on S s δ , one needs to verify the following three properties for each of these sets: (i) It contains with every q s4 also a spacelike neighborhood of q s4 . (ii) It contains with every q s4 also all other spin components of configurations with the same space-time positions. It now follows that the derivatives in H free are well defined. We still need certain configurations to describe creation and annihilation: (iii) It contains with every q s4 also all configurations with a y-particle added or removed in the 3d δ-neighborhood of any x-particle. Thus, all configurations needed for defining H z are included.
Likewise, corresponding to a partition of q 4 ∈ S δ into families as in Sect. 2.3.3, we regard the Hamiltonian H P j of family j defined in (37b) as acting on functions of q 4 . We set (adopting notation from [39,Sec. 5.3]). This operator can also be defined as acting on C ∞ (Q s4 ) by Since K P j Φ(q 4 ) depends only on values of Φ in a neighborhood of q 4 and with a y-particle added or removed, the action of K P j at q 4 ∈ S δ does not depend on whether we regard it as an operator on C ∞ (Q s4 ) or C ∞ (S s δ ).
This follows from Lemma 12. On C ∞ (Q s4 ), the commutator K P j , K P j vanishes at every q 4 ∈ S P δ . Even more, at every q 4 ∈ Q 4 , and whenever x k and x k keep their safety distance.
Proof. The yy commutator (157a) can easily be seen to vanish, as the free Hamiltonians are time-independent and commute everywhere. The xy commutator (157b) yields (with r indices acting on Φ). Since supp 3 ϕ ⊆ B δ (0) and supp 3 G rr s (t, ·) ⊆ B δ+|t| (0), the last expression will vanish if the first safety distance condition in (30) holds true, i.e., if (In the special cases that g rr s = δ rr g s or g rr s = h rr g s with self-adjoint matrix h, (159) also vanishes when x 0 k = x 0 k , but we do not use this fact.)

Existence
We now construct a solution Φ ∈ C ∞ (S s δ ) to the multi-time equations (37a) from initial data Ψ 0 ∈ H ∞ c . The construction proceeds in a way similar to the proof of Lemma 9 (and to the construction in Sections 5.3 and 5.5 of [39]). We define Φ on S δJ by induction over J. Put briefly, to obtain Φ at a configuration with J time values t 1 < t 2 < . . . < t J , we consider the families defined by a common time value: first evolve all particles to time t 1 , then those belonging to families 2 and up to time t 2 , and so on; see Fig. 10. Since the number M J of x-particles in family J may be less than M , we will have to use the appropriate versions of the Hilbert space H and Lemmas 1 through 5.
For J = 1, we obtain Φ ≡ Ψ on S s δ1 from the single-time evolution; by Lemmas 3 and 4, it is smooth on S δ1 ∼ = R × Q 3 and belongs to H ∞ c for every fixed t ∈ R.
The cases J = 2 and J = 3 allow a particularly simple construction that we want to describe first as it will play a role also for J > 3; we describe it for J = 2. (The same strategy can be applied if the number M of x-particles is ≤ 3.) Suppose that two regions G 1 , G 2 ⊂ R 3 and times t 1 ≤ t 2 ∈ R are such that {t j } × G j are 2δ-spacelike separated, i.e., for every x 1 ∈ G 1 and x 2 ∈ G 2 , x 1 − x 2 > 2δ + |t 1 − t 2 |. Then, for any t ∈ [t 1 , t 2 ], the sets G j := Gr(G j , |t j − t| + δ) are disjoint, and Φ on configurations concentrated in and in fact given by Φ(t 1 , ·, t 2 , ·) = W 1,t1−t ⊗ W 2,t2−t Ψ t (161) with W t as in (60) and W j,t acting on H Mj (G j ). Since, as explained before (132), for every configuration in S δ2 there is ε > 0 so that the union of the εballs around each particle are still 2δ-spacelike separated, Φ can be determined from Ψ via (161) everywhere on S δ2 . For J > 3 time values, this strategy cannot be directly applied because for some q 4 ∈ S δJ , there is no t at which the G j = Gr(q j , |t j − t| + δ) would all be mutually disjoint; for example, q 4 = (x 1 . . . x 4 ) with x 1 = (0, 0, 0, 0), x 2 = (1, 2, 0, 0), x 3 = (3, 5, 0, 0), x 4 = (4, 7, 0, 0).
We now turn to the strategy for general J by induction with anchor J = 1. The induction hypothesis asserts that Φ ∈ C ∞ (S s δJ ) is well defined and satisfies (37a) and (45) as well as a further condition that we will formulate in (162) below. We assume it for J and prove it for J + 1.
On S δJ+1 \S δJ , let us label the time variables so that t 1 < t 2 < . . . < t J+1 . The strategy is to solve the multi-time equation (37a) for j = J + 1 in Figure 10. Time evolution for a partition P = {P 1 , P 2 , P 3 } into three sets. At each t j , the time evolution for particles in P j generated by H j is switched off. The shaded regions indicate where the demand of being spacelike forbids particle coordinates of the multi-time configuration the variable t J+1 (while keeping t 1 , . . . , t J unchanged) from initial data given by Φ on S δJ , i.e., for t J+1 = t J . In fact, we solve it for 4 M1+N1+···+MJ +NJ functions, the components of Φ for different values of the spin indices for all particles in the families 1, . . . , J, where (M j , N j ) are the particle numbers in family j. Since the indices of other families are not acted upon in (37a), (37a) can be solved separately for each choice of values for those indices. By propagation locality, the solution on S (t 1 . (46) is determined by initial data on S (t 1 . . . t J ). We want to use Lemma 5 to conclude that Φ exists and is smooth where t 1 < t 2 < . . . < t J+1 . To this end, we regard t J+1 as the variable t of Lemma 5; we consider the evolution separately for every fixed choice of M 1 , N 1 , . . . , M J , N J ; we regard the t 1 , . . . , t J and the configurations q 1 , . . . , q J of the families 1 through J as the parameters λ in Lemma 5, so d = J + 3M 1 + 3N 1 + · · · + 3M J + 3N J . We know that the time evolution preserves the compactness of the 3-support in all space variables; however, in order to be able to apply Lemma 5, we need compact support in all λ variables, including t 1 , . . . , t J . Since wave functions do not have compact support on the time axis, we need to cut off the time dependence; this does not cause any harm because the solution provided by Lemma 5 is obtained by solving the 1-time evolution for every value of λ separately. That is, we apply Lemma 5 to the function f (t 1 ) · · · f (t J ) Φ(t 1 , q 1 , . . . , t J , q J , q J+1 ) instead of Φ(t 1 , q 1 , . . . , t J , q J , q J+1 ), where f : R → R is a smooth function such that f = 1 on [−T, T ] and f = 0 outside [−2T, 2T ]; for the desired result at q 4 ∈ S δ,J+1 , the value of T must be chosen larger than all absolute time values in q 4 . To fulfill the hypotheses of Lemma 5, we need that the smoothly cut off wave function lies in H ∞ cd . We will assume a little more as part of the induction hypothesis: (162) (We talk about extension because Φ(t 1 , q s3 1 , . . .) is not defined for every q s3 1 but only those δ-spacelike from the other families.) Then, Lemma 5 applies, and we obtain the desired function Φ on S δ,J+1 where all |t j | < T , and on all of S δ,J+1 by letting T → ∞. It remains to verify four things: (i) that Φ satisfies the summability condition (45); (ii) that Φ satisfies the multi-time equations (37a) for all j, not just j = J + 1; (iii) that Φ is smooth (the issue here is whether the transition from t J < t J+1 to t J > t J+1 is smooth); and (iv) that (162) holds on S δ,J+1 .
(i) Follows from the unitarity of the single-time evolution and Lemma 3. (ii) Since Φ was constructed on S δ,J+1 using the unitary time evolution in t J+1 , it is a solution of (37a) for j = J + 1 where t 1 < . . . < t J < t J+1 (in the strong sense for the same reasons as discussed after (64), i.e., K J+1 Φ = 0). Now we show that it is also a solution of (37a) for every j ≤ J, i.e., that K j Φ = 0. By Lemma 11, K J+1 K j = K j K J+1 on S δ,J+1 , so By induction hypothesis, K j Φ = 0 on S δJ . By (162), (62), and Lemma 3, K j Φ also satisfies (45) with J + 1 instead of J on S δ,J+1 . By Lemma 7, K j Φ = 0 on S δ,J+1 . (iii) By Lemma 5, Φ is smooth where t 1 ≤ . . . ≤ t J < t J+1 , but it is not obvious that Φ is smooth on the set where t J = t J+1 . To be sure, Lemma 5 provides a smooth function of t for t > 0, t = 0, and t < 0, but when we reduce the time variable of the family P J+1 to values less than t J we reach configurations for which P J and P J+1 would be labeled P J+1 and P J according to our prescription that the times are labeled increasingly. Still, for a fixed choice of the partition (in particular P J and P J+1 ) and in a neighborhood of a 4-configuration in which t J = t J+1 , while q J and q J+1 keep the safety distance, the constructed function Φ agrees with (161) with t J , t J+1 (the time variables of families P J , P J+1 ) playing the roles of t 1 , t 2 , in fact for t J > t J+1 as well as for t J < t J+1 . But by Lemma 5, Φ given by (161) is a smooth function. Hence, the constructed function Φ is smooth everywhere in S δ,J+1 . (iv) Statement (162) for J + 1 involves splitting q J+1 into two families, q J+1 = (q J+1 ,q J+2 ), considering them at time t J+1 , regarding t J+1 and q J+1 as parameters, and introducing a further factor f (t J+1 ). To check the definition of H ∞ cd , we know that f · · · f Φ is smooth, has compact support in t J+1 because of f , compact support in t 1 , . . . , q J by assumption, and compact 3-support in q J+1 because, by Lemma 2, the 3-support can grow (from that of Φ| S δJ ) by at most 2T + δ; so it remains to check (62): Indeed, for any given degree n of differentiation, at most the first n derivatives of f can occur (independently of N ), all of which have finite L 2 norm; moreover, t J+1 -derivatives in (62) (without loss of generality the rightmost derivatives in ∂ α ) can be replaced by −iH J+1 as in (56); now each term in H J+1 Φ (N ) involves either a spatial derivative of Φ (N ) (which obeys (62) by Lemma 3) or ϕ Φ (N +1) or ϕ Φ (N −1) , whose contributions to (62) involve derivatives of ϕ of up to n-th order and a factor of √ N ≤ N (so m → m + 1), so they remain finite since Φ itself obeys (62) at every t = t J+1 by Lemma 3. Thus, (162) holds on S δ,J+1 .
This completes the construction of Φ on S δ and the proof that Φ is smooth and solves the multi-time equations (37a).
Remark. 7. One might be tempted to think that the following procedure yields an alternative construction of the solution. The idea is to construct a function Φ on Q 4 which agrees with the desired solution on S δ (or S δ ) by writing the multi-time equations in the form: r k ,sN+1 x 4M , y 4N , (0,ỹ) and to integrate them in a particular order. Specifically, solve first the equation for x 1 from 0 to any desired x 0 1 , keeping all other x 0 k and y 0 at 0; then, solve the equation for x 2 , and then for x 3 , . . . , x M , keeping all y 0 at 0; finally, solve (165) for all ys. However, the function Φ thus obtained will not be a solution of our multi-time equations (37a) on S δ ; indeed, Φ does not even agree with Ψ on the set S δ1 of simultaneous configurations. That is because when all x 0 k = t = y 0 , then Φ = e −iH y t T e −i t 0 H x M (s) ds · · · T e −i t 0 H x 1 (s) ds Ψ 0 whereas (166) with T e the time ordered exponential, H y = H y1 + · · · + H yN on the Nsector and H x k (x 0 k ), H x k , and H y the right-hand sides of (164), (152), and (165). That the expressions are not the same is strongly suggested by the facts that the H x k (s) do not commute with each other (and e A e B = e A+B when AB = BA) and H x k = H x k .
Proof of Theorem 1. Theorem 1 follows by putting together the statements obtained about the existence of solutions with Lemmas 9 and 10. It remains to prove the permutation (anti-)symmetry of Φ. Since all Hamiltonians and time evolution procedures are invariant under relabeling, Φ with permuted labels is the solution with the correspondingly permuted initial condition. Thus, if permuting two x's changes the sign of the initial data, then the corresponding permutation in Φ changes the sign of Φ; likewise, if permuting two y's leaves Ψ 0 invariant, then the corresponding permutation in Φ leaves Φ invariant, as claimed.

Proof of Remark 5
Proof. Given Φ on S δ , it can be extended to Φ on S δ by freely varying all y , i.e., The resulting Φ is smooth, for example, by Lemma 5. Since the H y (and the H y ) commute pairwise, Φ satisfies the multi-time equations (56). Since the strong solution of the 1-particle Dirac equation is unique and Φ is unique by Theorem 1 and (45), Φ is unique.

Conclusions
In this paper, we have provided a rigorous study of a system of multi-time equations (37a) for a model quantum field theory with UV cut-off. We have proved the existence and uniqueness of solutions on the set S δ of δ-spacelike configurations (and even on a larger set S δ allowing arbitrary points for bosons), and thus the consistency of the multi-time equations. This result supports the viability of multi-time wave functions as a covariant expression of the quantum state in the particle-position representation, although the model considered here is not fully covariant, partly because of the UV cut-off.
Our proof is one of the first rigorous consistency proofs for multi-time formulations of quantum field theories. Only two other results of this kind are known to date: First, for a similar set of equations, first proposed by Dirac, Fock, and Podolsky [9] and involving a fixed number of time variables, consistency was recently proved in [36]. And second, for a model with a variable number of time variables in 1+1 dimensions and with a cut-off in the particle number, consistency was recently proved in [25].
Since we considered solutions in the classical sense (i.e., differentiable functions, rather than weak derivatives), we also had to prove smoothness of the solutions. Furthermore, our proof establishes in particular that the nonrigorous arguments for consistency in [39] also apply rigorously. In fact, it turns out that all considerations of [39] are rigorously valid if formulated appropriately in view of the UV cut-off. In particular, the model satisfies, up to a tolerance of the size δ of the UV cut-off, the conditions "propagation locality" and "interaction locality" that played important roles for the derivation of Born's rule on arbitrary Cauchy surfaces in [30].
For the future, it would be of interest to move toward more realistic models of quantum field theory and to obtain a multi-time formulation of quantum electrodynamics (QED). For example, we believe that a consistent multi-time formulation of the Landau-Peierls model of QED [18] is possible and will have the advantage of allowing to switch easily between electromagnetic field tensors and vector potentials, which is not possible in the single-time formulation used by Landau and Peierls.