Generalized Causal Set d'Alembertians

We introduce a family of generalized d'Alembertian operators in D-dimensional Minkowski spacetimes which are manifestly Lorentz-invariant, retarded, and non-local, the extent of the nonlocality being governed by a single parameter $\rho$. The prototypes of these operators arose in earlier work as averages of matrix operators meant to describe the propagation of a scalar field in a causal set. We generalize the original definitions to produce an infinite family of ''Generalized Causet Box (GCB) operators'' parametrized by certain coefficients $\{a,b_n\}$, and we derive the conditions on the latter needed for the usual d'Alembertian to be recovered in the infrared limit. The continuum average of a GCB operator is an integral operator, and it is these continuum operators that we mainly study. To that end, we compute their action on plane waves, or equivalently their Fourier transforms g(p) [p being the momentum-vector]. For timelike p, g(p) has an imaginary part whose sign depends on whether p is past or future-directed. For small p, g(p) is necessarily proportional to p.p, but for large p it becomes constant, raising the possibility of a genuinely Lorentzian perturbative regulator for quantum field theory. We also address the question of whether or not the evolution defined by the GCB operators is stable, finding evidence that the original 4D causal set d'Alembertian is unstable, while its 2D counterpart is stable.


Introduction
Causal set theory postulates that the fundamental structure of spacetime is that of a locally finite partially ordered set [1]. 1 Its marriage of discreteness with causal order implies that physics cannot remain local at all scales. To appreciate why this should be, let us consider how one might define a notion of "closeness" in a causal set, confining ourselves to causal sets C which are obtained by randomly selecting points from a Lorentzian manifold M 1 Characterized mathematically, this is a set C endowed with a binary relation ≺ such that for all x, y, z ∈ C the following axioms are satisfied: (1) transitivity: x ≺ y & y ≺ z ⇒ x ≺ z; (2) irreflexivity: x ⊀ x; (3): local finiteness: |{y ∈ C|x ≺ y ≺ z}| < ∞. Thus a causal set (causet) is in a certain sense both Lorentzian [in virtue of (1) and (2)] and discrete [in virtue of (3)]. and endowing the selected points with the causal relations inherited from the manifold. 2 Given such a causet, any intrinsically defined notion of closeness between two elements of C will reflect their Lorentzian distance in the embedding spacetime. But a small Lorentzian distance between two points of M does not mean that they are confined to a small neighbourhood within M . Rather, the second point can be "arbitrarily distant" from the first, as long as it is located near to the lightcone of the latter. Thus, an element of C will inevitably possess very many "nearest neighbours", no matter how that notion is formalized. In this manner, the concept of locality provided by the topology of a continuous spacetime manifold is lost.
This nonlocality manifests itself concretely when one seeks to describe the wave propagation of a scalar field on a causal set by defining a discrete counterpart of the d'Alembertian operator, . For the aforementioned reasons, it seems impossible to proceed in analogy with what one does when, for example, one discretizes the Laplacian operator in a Riemannian spacetime. Nevertheless, a non-local operator was suggested in [2] which on average reproduces in the appropriate continuum limit for 1 + 1 dimensional Minkowski space M 2 (i.e. for causets derived by sprinkling M 2 ). The expression introduced in [2] was generalized to D = 4 dimensions in [3] and recently to arbitrary D in [4].
We shall denote a discrete causal set d'Alembertian designed for M D by B (D) ρ , where ρ (dimensionally an inverse spacetime volume) is a volume-scale that controls the extent of the non-locality. In the case of causal sets which are well-approximated by D-dimensional Minkowski space M D , averaging B (D) ρ over all such causets (i.e. averaging over all sprinklings of M D in the sense of footnote 2) leads to a non-local and retarded continuum operator (D) ρ defined in M D . We shall refer to this operator as the continuum causal set d'Alembertian. Its crucial property is that it reproduces the usual d'Alembertian in the limit of zero non-locality scale: (D) ρ φ → φ as ρ → ∞ for test-functions φ of compact support.
Although the causet operator B (D) ρ is necessarily nonlocal, one might expect that the range of its nonlocality could be confined to the discreteness scale itself. In other words, one might expect that ρ ∼ −4 , being the -presumably Planckian -discreteness length. However, one can also cite reasons why one might need to have ρ −4 , leading to a more long-range nonlocality. 3 Although these reasons are not conclusive, let us accept them 2 This process is known as Poisson sprinkling: Given a spacetime M , let the discrete subset of points, C, be one particular realization of a Poisson process in M , and let the elements of C retain the causal relations they have when regarded as points of M . In order that the resulting precedence relation on C approximately encode the metric of M , one must exclude spacetimes with closed causal curves, for example by requiring M to be globally hyperbolic. 3 The issue here concerns the behavior of B (D) ρ for one particular sprinkling versus its behavior after averaging over all sprinklings. The latter converges to as ρ → ∞ but the former incurs fluctuations which grow larger as ρ → ∞ and which therefore will be sizable if ρ is the sprinkling density, −4 . Which behavior is relevant physically? In full quantum gravity some sort of sum over different causets will be involved, including in particular a sum over sprinklings. Such a sum differs from a simple average and might or might not damp out the fluctuations, or they might cancel in other ways. But if neither of these things happens, the only way out [2] would be to choose ρ small enough that the necessary averaging will occur within each individual causet. provisionally. A natural question then arises: might such a "mesoscopic" nonlocality show up at energy-scales accessible by current experiments?
Ideally, one would address this question in the fully discrete setting, but it seems much easier to begin with the continuum version of the same question by asking what changes when the local operator is replaced by the nonlocal operator ρ . In this paper, we make a start on answering this question by analysing the "spectral properties" (Fourier transform) of a family of continuum operators (D) ρ . In Section 2, we discuss the continuum operators corresponding to the original 2D [2] and 4D [3] causet d'Alembertians, and in Section 3 we generalize the discussion to an infinite family of operators parametrized by a set of coefficients, {a, b n }, for which we derive explicit equations that ensure the usual flat space d'Alembertian is recovered in the infrared limit. Based on the UV behaviour of these operators (which we determine for all dimensions and coefficients {a, b n }), we propose a genuinely Lorentzian perturbative regulator for quantum field theory (QFT). Finally, we address the question of whether or not the evolution defined by the (classical) equation (D) ρ φ = 0 is stable. We devise a numerical method to test for stability and present strong evidence that the original 4D causal set d'Alembertian is unstable in this sense, while its 2D counterpart is stable.

The Original 2D and 4D Causet d'Alembertians
In this Section we discuss the original continuum causet d'Alembertians for dimensions two [2] and four [3]. Let us start by establishing some terminology. Given any two elements x, y of a causal set C, we define the order interval Int(x, y) between them as the set of all elements common to the (exclusive) future of x and the (exclusive) past of y: Int(x, y) = {z ∈ C|x ≺ z ≺ y}. Notice that in our convention, Int(x, y) does not include x or y. An element y ≺ x is then considered a past nth neighbour of x if Int(y, x) contains n elements. For instance, y is a 0th neighbour of x if Int(y, x) is empty, a first neighbour if Int(y, x) contains one element, and so on (see Figure 1 for an example). We denote the set of all past nth neighbours of x by I n (x). Throughout the paper, we will only consider causal sets which are obtained by Poisson sprinklings of Minkowski space at density ρ.

2D
The original causet d'Alembertian for dimension 2, which we denote by B (2) ρ , acts on a scalar field Φ(x) on the causal set in the following way [2]: 1 Φ(y 1 ), respectively. The continuum limit, or rather average, of (B (2) ρ Φ)(x) can be understood as follows: fix the point x, keep sprinkling at density ρ and compute (B (2) ρ Φ)(x) for every sprinkling. The average of all these values is equal to (

Figure 1 illustrates how B
(2) ρ is defined, given a Poisson sprinkling of 2D Minkowski space M 2 . The continuum operator (2) ρ is obtained by averaging B (2) ρ over all such Poisson sprinklings at density ρ : Here J − (x) denotes the causal past of x, and V (x − y) is the spacetime volume enclosed by the past lightcone of x and the future lightcone of y. Note that (2) ρ is a retarded operator, in the sense that (2.3) uses information only from the causal past of x.
The operator (2) ρ can be studied by analysing its action on plane waves. Due to translation symmetry of Minkowski space, 4 any plane wave e ip·x is an eigenfunction of (2) ρ (provided that the integrals in (2.3) converge, so that the left hand side is well defined): where p · x ≡ η µν p µ x ν and η µν = diag(−1, 1). Interestingly enough, g ρ (p) in this case can be expressed in closed form: 5 ρ −1 Im g (2) ρ (p) ρ (p), which (for real p) depends only on Z = ρ −1 p · p , and on sgn(p 0 ) when p is timelike. (b) The spectrum g (2) ρ (p) of the original 2D continuum causet d'Alembertian for real momenta p . For spacelike momenta (p·p > 0), g (2) (p) is real. For timelike momenta, it is complex with an imaginary part whose sign is opposite for past-directed and future-directed momenta.
where E 2 (z) is a generalized exponential integral function (see e.g. 8.19 of [5]) and (2.6) Here, as illustrated in Figure 2, E 2 (z) assumes its principal value, with a branch cut along the negative real axis. For real and spacelike momenta (Z > 0), g (2) is real. For real and timelike momenta (Z < 0), its value above/below the branch cut corresponds to past/future-directed momentum-vectors. There, g (2) ρ is complex and changes to its complex conjugate across the cut. That the spectrum is different for past and futuredirected momenta should come as no surprise, given that (2) ρ is retarded by definition. We will see in Section 3 that these features persist in all dimensions and for a much broader class of causet d'Alembertians.
The infrared (IR) and ultraviolet (UV) behaviours of g (2) ρ (p) are easily deduced from the asymptotic forms of E 2 (Z) (see e.g. 8.11.2, 8.19.1, and 8.19.8 of [5]): (2.8) The first of these two equations shows that the usual d'Alembertian is indeed reproduced in the limit of zero non-locality. The second equation, on the other hand, reveals a UV behaviour quite unlike that of the usual d'Alembertian; in Section 3.2 it will lead us to propose a new regularization scheme for quantum field theory.
An important question is whether the evolution defined by (2) ρ Φ = 0 is stable or not. To a large extent this is answered by the fact that the only zero of g (2) ρ (p) occurs at Z = ρ −1 p·p = 0. To demonstrate this, we note that g (2) ρ (p) has the following representation (see e.g. 8.19.1 and 8.6.4 of [5]): It therefore suffices to prove that f (Z) has no zeros when Z = 0. But the imaginary part of f (Z) is Because the integral that multiplies −Im(Z)/2 in (2.10) is strictly positive, Zf (Z) could vanish only for real Z. Obviously, it does vanish for Z = 0, but elsewhere on the real axis, it remains nonzero, as illustrated in Figure 2b. What we have just proven is that a plane wave solves the equation (2) ρ Φ = 0 iff it solves the equation Φ = 0. To the extent that the general solutions of these two wave equations can be composed of plane waves, they therefore share the same space of solutions. This, of course, is an important result in itself. But it also, a fortiori, answers the stability question in the affirmative, since we know that the evolution corresponding to is stable.
If there remains any doubt about stability or about the fact that both Φ = 0 and (2) ρ Φ = 0 yield the same evolution, it springs from a possible uncertainty about boundary conditions. In the usual situation (that of the ordinary d'Alembertian ), one understands how to relate a general solution to its initial data on an arbitrary Cauchy surface, and when Φ falls off suitably at infinity, its total energy is defined and conserved. From energy conservation, stability also follows -relative to the given choice of boundary conditions. On the other hand in the case of (2) ρ , a connection between solutions and Cauchy data remains to be found, as does a better understanding of appropriate falloff conditions. But absent some such boundary condition there is nothing to exclude complex momenta p that lead to exponential growth in time, e.g. an imaginary multiple of a real lightlike vector.
For these reasons, we would like to discuss stability from a slightly different angle, which also will be helpful when we come to deal with the 4D case. Quite generally, instabilities tend to be associated with exponentially growing "modes" (in this case plane waves). Let us then assume that we can take this as our criterion of (in)stability. And to exclude the kind of "fake instability" mentioned above, let us also require any putative unstable mode, Φ(x) = e ip·x , to be bounded at spatial infinity in at least one Lorentz frame. (Unfortunately we cannot say "in all Lorentz frames", since for a plane wave, exponential growth in time induces exponential growth in space via a Lorentz boost.) We might hope that the condition just formulated is equivalent to the following more natural one: consider only solutions of (2) ρ Φ(x) = 0 which have compact support on every Cauchy hypersurface (compact spatial support in every frame.) Be that as it may, if this criterion is accepted, then we can establish stability very simply in the present case, because an unstable mode, Φ(x) = e ip·x , is then precisely one such that p possesses a future-directed timelike imaginary part: p = p R + ip I with p I · p I < 0 and p 0 I > 0. This, however, is impossible for Z = 0, as one sees from the equation 0 = p · p = p R · p R − p I · p I + 2ip R · p I , whose right-hand side has a strictly positive real part when p I is timelike. For logical completeness, we should also observe that (2.5) is valid for all complex p whose imaginary parts are timelike and future-directed. (For more general complex momenta, the integral defining (2) ρ Φ might not converge, a circumstance that, depending once again on the choice of falloff conditions, might or might not impinge on the claimed identity between our solutions and those of the ordinary wave equation.)

4D
The causet d'Alembertian for dimension 4, has the same general form as that for M 2 , but with different coefficients [3] : The continuum average (4) ρ then also takes a similar form: We will show in Section 3.1 that the "spectrum" of (4) ρ , as defined by (4) ρ e ip·x = g (4) ρ (p)e ip·x , is given by where K 1 is a modified Bessel function of the second kind and All functions in (2.14) assume their principal values with branch cuts along the negative real axis. Many properties of the 2D function g (2) ρ (p) carry over to g (4) ρ (p) . For timelike p, the value of g (4) ρ (p) above/below the branch cut corresponds to past/future-directed momenta, and it changes to its complex conjugate across the cut. Also, g (4) ρ is real for spacelike momenta. Figure 3b shows the behaviour of g (4) ρ (p) for real momenta. The IR and UV behaviours of g (4) ρ (p), which are derived in Sections 3.2 and 3.3, are given by ρ | are plotted as a function of the real and imaginary parts of ρ (p) of the original 4D causet d'Alembertian for real momenta p. For spacelike momenta (p · p > 0), g (4) (p) is real. For timelike momenta, it contains also an imaginary part whose sign is opposite for past-directed and future-directed momentum-vectors.
Again, the IR behaviour confirms that the usual d'Alembertian is reproduced in the limit of zero non-locality. The UV limit has the form of a constant plus a term proportional to p −4 . The inverse of g (4) ρ (p), which defines the retarded Green's function in Fourier space, takes exactly the same form in the UV: (2.18) In any QFT based on (4) ρ , the propagator associated with internal lines in Feynman diagrams would presumably have the same UV behaviour. Subtracting the constant term from the propagator (which corresponds to subtracting a δ-function in real space) would then render all loops finite. This procedure could be the basis of a genuinely Lorentzian regularization and renormalization scheme for QFT. We will discuss these things more generally in Sections 3.3 and 3.4.
We have only been able to address the question of stability by numerical means in this case, and we refer the reader to Section 3.5. It turns out that g (4) ρ (p) does in fact have unstable modes in the sense that there exist complex momentum-vectors p which satisfy g (4) ρ (p) = 0, and whose imaginary parts are timelike and future-directed. Such a mode corresponds to a complex zero of g (4) ρ in the complex Z-plane, and Figure 3a shows one such zero (the other one being its complex conjugate).

The Generalized Causet Box (GCB) Operators
The key property of the causet d'Alembertians introduced in the previous Section is that they reproduce in the continuum-averaged (averaged over all sprinklings) and local (ρ → ∞) limit. In this Section, we explore a larger family of operators B Considering all these requirements, (B (D) Φ)(x) can be expressed in the following general form where {a, b n } are dimensionless coefficients and I n (x) is the set of all nth neighbours to the past of x (see beginning of Section 2). This is a straightforward generalization of (2.1) and (2.11), where we have now allowed ourselves up to L max neighbours. We will soon see that recovering requires keeping a minimum number of layers: e.g. L max ≥ 2 in 2D and L max ≥ 3 in 4D. The original 2D and 4D proposals are then the minimal cases in this sense.
The continuum-average acts on a scalar field Φ(x) in the following way: Here as before, J − (x) denotes the causal past of x, while V (x, y) is the spacetime volume enclosed by the past light cone of x and the future light cone of y.
The occurrence of the factor e −ρV in (3.2) shows that the parameter ρ (which dimensionally is an energy-density) functions as a kind of "nonlocality scale" controlling the distance over which the operator (D) ρ acts. As our definitions stand so far, this nonlocalityscale directly reflects the fundamental discreteness-scale, because (3.2) was derived under the assumption that ρ was the sprinkling-density in M D . However it turns out that one can decouple the two scales by tweaking the definition (3.1) in such a way as to produce a more general causet operator whose sprinkling-average reproduces the same continuum operator (3.2), even when ρ is smaller than the sprinkling density. With this operator, the nonlocality can extend over a much greater distance than that of the fundamental discreteness-scale. Although modifying B (D) ρ in this way has the disadvantge of introducing a second, independent length scale, it allows one to overcome a potential difficulty pointed out in [2], namely that (3.1) with fixed coefficients leads to fluctuations in (B (D) ρ Φ)(x) which grow with ρ, rather than diminishing. We have provided the definition of this "tweaked" operator and the derivation of its continuum average in Appendix D; but henceforth, we will concern ourselves exclusively with the continuum operator (D) ρ , without worrying about its relationship with any underlying discreteness. Correspondingly, ρ will henceforth denote a non-locality-scale with no necessary relation to any discreteness scale.

Spectrum
That any plane wave e ip·x is an eigenfunction of (D) ρ in M D follows from translational symmetry: V (x, y) = V (x − y). It can be shown in fact that 3) where V (y) = V (O, y) is the spacetime volume enclosed by the past light cone of y and the future light cone of the origin: .

(3.6)
Evaluating χ(p, ρ) amounts to computing the Laplace transform of a retarded, Lorentzinvariant function, which has been done in [6]. It follows from their result that where K ν is the modified Bessel function of the second kind. All functions in (3.7) assume their principal values, with a branch cut along the negative real axis. This result is valid for all p whose imaginary part is timelike and future-directed, i.e. p I · p I < 0 and p 0 I > 0, where p = p R + ip I and the Lorentzian norm is given by p · p = p R · p R − p I · p I + 2ip R · p I . For momenta satisfying these conditions, the integral that defines χ(p, ρ), and consequently (D) ρ e ip·x , is absolutely convergent. Plugging (3.7) into (3.4) we find where Z is a dimensionless quantity defined by ρ (p) can be defined by first adding a small future-pointing and timelike imaginary part p I to p R , and then taking the limit as p I shrinks: When p R is timelike, this amounts to changing Z = ρ − 2 D p R · p R on the right hand side of (3.8) to Z + i for past-directed, and Z − i for future-directed p R . This is illustrated in Figure 2a. Because of the appearance of Z 1/2 in (3.8) and the fact that K ν (z) = K ν (z), it follows for timelike p that ρ (p) differs for past-and future-directed timelike p. This is to be expected, since requiring ρ (p) is real, as follows from the fact that K ν (z) is real when ν is real and ph(z) = 0 [5].

IR Behaviour
We want to choose the coefficients a and b n so that the usual d'Alembertian operator is recovered in the limit of zero non-locality: This requirement is equivalent to demanding In Appendix A, we derive equations for a and b n which guarantee this behaviour for an arbitrary spacetime dimension D. We expand Z 2−D 4 K D 2 −1 (Z 1/2 s) on the right hand side of (3.8) about Z = 0, and arrange a, b n so that the terms which grow faster than Z vanish, while the coefficient of the term proportional to Z is −1. We state the main results here and refer the reader to Appendix A for the details.
In even dimensions, letting D = 2N + 2 with N = 0, 1, 2, . . . , the equations that need to be satisfied are where ψ(n) is the digamma function. Equations (3.14a) and (3.14c) determine b n , after which (3.14b) fixes a. The minimum number of terms required to solve these equations is determined by L max ≥ N + 2. In 2D and 4D in particular, keeping this minimum number of terms leads to the solutions (2.2) and (2.12), respectively. In odd dimensions, letting D = 2N + 1 with N = 0, 1, 2, . . . , the equation are (3.15c) Similarly to the even case, Equations (3.15a) and (3.15c) determine b n , after which (3.15b) fixes a. The minimum number of terms is determined by L max ≥ N + 1.

UV Behaviour and the Retarded Green's Function
The UV behaviour of g (D) ρ (p), as derived in Appendix B, is  ρ (p) has no zero in complex plane apart from at p · p = 0, this choice of contour ensures that G R is indeed retarded. As we will argue in the next section, the presence of such zeros implies that evolution defined by (D) ρ is unstable. Therefore, we shall ignore these cases for our current discussion.
The behaviour of G R (x − y) in the coincidence limit x → y is determined by the behaviour ofG R (p) at large momenta: ρ . The direction of integration is taken to be counter-clockwise.
Here we have assumed a = 0. When a = 0,G R (p) scales as p D for large momenta, a badly divergent UV behaviour. Therefore we will confine ourselves to cases where a = 0. The constant term 1 a represents a δ-function in real space. The other terms in the series have the form d D p p −nD , n = 1, 2, · · · , and it can be shown that they are all finite. It then looks like subtracting 1 a δ (D) (x − y) from ρ 2/D G R (x − y) must result in a completely smooth function in the coincidence limit, and we will now show this is indeed the case.
Although D = 4 is the dimension of greatest interest, the proof which we shall present is valid in all even dimensions. Let us define Our task is then to show G(x − y) is a smooth function at x = y. It follows from (3.18) that From this it follows for all x y that where the first equality is obtained by changing p → −p and then using (3.11), and the second equality is a direct consequence of (3.22) with x and y interchanged. Returning to (3.21), and subtracting zero in the form of (3.23), we obtain where the second equality is true because g ρ (p) is real for space-like momenta. (Note that the 1 a term contributes only when x = y.) In what follows, we let as given in the right hand side of (3.8).
The integral in (3.25) can be divided into two integrals over p 0 > 0 and p 0 < 0. For a fixed sign of p 0 , g (D) ρ (p) is only a function of p · p, making (3.25) the Laplace transform of a Lorentz-invariant function. Similarly to how we derived (3.7), we use the result of [6] to compute G(x − y): and is a small positive number which should be taken to zero at the end of calculations. When x − y is timelike and future-directed, we can let s 2 = −iτ xy where τ xy > 0. Using properties of Bessel functions (see e.g. 10.27.9 of [5]), (3.27) can be simplified into the following form for even D: (3.28) g. 10.2.2 of [5]) and the fact that Im g(−ξ 2 + i ) is exponentially damped for large ξ (see Appendix B.1), it can be verified that Thus G(x − y) approaches a constant in the coincidence limit. 6 Strictly speaking, the discussion above only analyzes the behavior of G(x − y) as τ xy approaches 0, and consequently it does not exclude the presence of terms which blow up discontinuously on the light cone, such as δ(τ 2 xy ). However, a similar treatment for the case where x − y = 0 is null rather than timelike removes this loophole.

A Possible Regularization Scheme for Quantum Field Theoy
As was shown in the previous Section, changing the usual d'Alambertian to the nonlocal operator (D) ρ makes the coincidence limit more divergent, rather than smoothing it out as one might have initially expected. But it does so in an interesting way: all the divergences have now been absorbed into one δ-function at x = y. This feature has a natural application as a regularization tool for quantum field theory. In any QFT based on (D) ρ , one would expect the propagator associated with internal lines in Feynman diagrams to have the same UV behaviour as (3.19). Subtracting the constant term in (3.19) (which corresponds to subtracting a δ-function in real space) would then render all loops finite. This would be a genuinely Lorentzian regulator, with no need for Wick rotation. It would also be physically motivated, with the "UV completion" being understood as a theory on the causal set. It would be interesting to apply this technique to the renormalization of some well-understood scalar field theories.

Stability
Is the evolution defined by (D) ρ stable? As we discussed in Section 2.1, instabilities are in general associated with "unstable modes", and we agreed to use this as our criterion of instability for purposes of this paper. More specifically, we took such a mode to be a planewave Φ(x) = e ip·x satisfying the equation of motion (D) ρ Φ(x) = 0, with the wave-vector p possessing a future-directed timelike imaginary part (i.e. p = p R + ip I where p I · p I < 0 and p 0 I > 0). The necessary and sufficient condition for avoiding unstable modes is theñ whereg(Z) is defined in (3.26). Let us argue why this is the case. First observe that plane solutions of our wave-equation correspond exactly with zeros ofg(Z). If the above condition is verified, then the only such zero is at Z = 0, just as for the usual d'Alembertian. But we know (as is also easy to demonstrate ab initio) that there are no unstable modes in the usual case. Conversely, when the above condition is violated for some complex Z = 0, it is always possible to find a corresponding p with a timelike and future-directed 6 One can understand intuitively why GR(x − y) is the sum of a δ-function with a bounded remainder by noticing that (up to an overall numerical factor) our nonlocal d'Alembertian operator has the form 1 − S, where the '1' corresponds to the first term in (2.3) or (3.2) and the remainder S is given by an integralkernel which is both bounded and retarded. The inverse operator GR would then be GR = (1 + S) −1 = 1 + SGR = 1 + S + S 2 + S 3 · · · , a series that should converge sufficiently near to x = y. Since the operator 1 is represented by a term of δ(x − y) in GR(x − y), one sees that GR(x − y) is the sum of a δ-function with a term involving only smooth bounded functions.
imaginary part which satisfies p · p = ρ 2 D Z. To see this, we let p = p R + ip I and take p R =< π 0 R , π R > and p I =< π I , 0 > with π I > 0. This is always possible because p I is timelike and future-directed. The equations that need to be satisfied are Substituting for p I leads to This system of equations always has a solution. In fact, there is a whole family of such unstable modes parametrized by π I . Note however that the condition | π R | 2 > 0 puts an upper bound on the value of π I , and therefore on the growth rate of such an instability. We have thus reduced the question of whether or not has unstable modes to the question of whetherg(Z) has zeros other than Z = 0 in the complex plane. We can answer this question by counting the zeros ofg(Z) with the aid of the "argument principle" of complex analysis: 1 2πi where N and P are the number of zeros and poles, respectively, inside of the closed contour C, which we choose as shown in Figure 4b. The number of poles inside C is zero because all terms appearing ing(Z) are finite in that region (at least when L max is finite). As shown in Figure 4b, the path of integration C comprises four pieces: C 2 and C 4 run from −∞ to 0 a distance above and below the negative real axis respectively, C 3 is a semicircle of radius about the origin, and C 1 is (almost) a circle whose radius should be taken to infinity. For large Z we have from (3.16), and it follows that (We remind the reader of our standing assumption that a = 0. See the remarks following (3.19).) On the other hand the IR behaviour,g(Z) Also, becauseg(x + i ) =g(x − i ) for x < 0: Performing this last integral will allow us to determine whether (D) ρ has unstable modes or not.
Given a choice of the parameters a and b n , the last integral can be computed numerically. In the minimal 4D case discussed in Section 2.2, we find that (4) ρ has precisely two "unstable zeros". (Notice that becauseg(Z) =g(Z), if Z is a zero ofg(Z), so also isZ.) We have located these zeros numerically, as shown in Figure 3a. With different choices of the parameters {a, b n }, the number of zeros can change, but we have not been able to find any choice that would make (4) ρ stable. It would be interesting to find an analytical method to check for stability.

Summary and Remarks
We have defined an infinite family of scalar-field operators on causal sets which we dubbed Generalized Causet Box (GCB) operators. For causal sets made by sprinkling D-dimensional Minkowski space M D , these operators reproduce the usual d'Alembertian = ∇ µ ∇ µ when one averages over all sprinklings and takes the limit of infinite sprinkling-density ρ. If, on the other hand, one averages over all sprinklings while holding ρ fixed, one obtains an integral operator (D) ρ in M D which is manifestly Lorentz-invariant, retarded, and nonlocal, with the degree of nonlocality set by ρ. In the present paper, we have been concerned primarily with these continuum operators, whose nonlocality can be regarded as a "mesoscopic" residue of the underlying causal set discreteness.
The GCB operators B are parametrized by a set of coefficients, and we derived the equations in these coefficients which ensure that is recovered in the infrared limit. The minimal solutions of these equations turned out to reproduce the original operators proposed in [2]. We also computed the Fourier transform of (D) ρ , or equivalently its "spectrum of eigenvalues" obtained by applying it to an arbitrary plane wave. For spacelilke momenta the spectrum is real. For timelike momenta it contains also an imaginary part, which changes sign under interchange of past with future. The UV behaviour of the spectrum differs from that of in a way which led us to propose a genuinely Lorentzian, perturbative regulator for quantum field theory.
We also studied the question of whether the evolution defined by the continuumaveraged GCB operators is stable. This is of interest in relation to nonlocal field theories based on ρ . The continuum-average of the minimal 2D causal set d'Alembertian was shown to be stable by a direct proof. In 4D we did not settle the question analytically, but we devised a numerical diagnostic that applies to all the operators (D) ρ , and which disclosed a pair of unstable modes when applied to the minimal 4D causal set d'Alembertian. Are any of the continuum-averaged GCB operators stable in 3 + 1 dimensions? We were not able to find any, but there are an infinite number of such operators and a definitive search could only be conducted by analytical means. 7 Finally, it bears repeating that there might be more reliable indicators of instability than simply the existence of an exponentially growing plane-wave solution, which a priori tells us nothing about the behavior of solutions of limited spatial extent. For that reason, it would be worthwhile to analyze directly the late-time behavior of the Green function G R (x − y) which is inverse to (D) ρ . If it were bounded that would imply stability, and if it grew exponentially, that would imply instability.
Our results also suggest other problems for further work. It would be interesting, for example, to work out the continuum-averaged GCB operators in curved spacetimes. It was found in [3] that the minimal 4D operator has the following limit as ρ → ∞: where R is the Ricci scalar. (In fact one obtains the same limit in all dimensions D.) Would this feature persist for all of the GCB operators? This feature has also been used to define an action-functional for causal sets [3]. A final question then is whether the instability found above has any consequences for this causal set action?

A IR Behaviour of the GCB Operators: Details
Here we will derive the equations that the constants a and {b n } should satisfy in order for (D) ρ to have the desired IR behaviour (3.12), or equivalently (3.13), which in turn is equivalent tog as given in the right hand side of (3.8).

A.1 Even Dimensions
Let D = 2N + 2 where N = 0, 1, 2, . . . . Theñ In order to examine the behaviour ofg(Z) as Z → 0, we need to expand (Z 1/2 s) −N K N (Z 1/2 s) in this regime. From the power series expansion of K N (see e.g. 10.31.1 and 10.25.2 of [5]), it follows that where ψ(n) is the digamma function. Because we need the leading behaviour of ρ − 2 Dg (Z) to be −Z, we have only considered terms up to this order. All the terms in (A.4a) and (A.4b) diverge as Z → 0, forcing us to pick the b n such that none of them contribute tõ g(Z) in the Z → 0 limit. The contribution of the term (A.4d) is also unwanted and should be made to vanish by choosing b n appropriately. This leads us to the following series of equations: The integration over s can be performed (see e.g. 5.9.1 of [5]) to give us the condition reproduced above as equation (3.14a): Lmax n=0 b n n! Γ(n + k + 1 N + 1 ) = 0, k = 0, 1, . . . , N + 1.
Requiring the contribution of the constant term (A.4c) to vanish yields We can perform the integral over s by using the formula (see e.g. 5.9.19 and 5.9.1 of [5]) Finally, requiring the contribution of (A.4e) to reproduce the desired −Z behaviour leads to As before, we have only kept track of terms up to Z. The contributions of all the terms in (A.13a) should be made to vanish; this leads to the equation

B UV Behaviour of the GCB Operators: Details
Here we derive the UV behaviour of (D) ρ . We will make use of the following identity [6], which holds for arbitrary natural number m: Notice that both these terms are real for both positive and negative Z, because D/2 is an integer when D is even. In order to produce the sub-leading terms, one can continue integrating by parts in (B.9). The sub-leading terms are thus also real, whence the imaginary part ofg(Z) must, for even D, decay faster than any power of Z for Z → ∞. This behavior can be seen in Figures 2b and 3b.