Holographic collisions in confining theories

We study the gravitational dual of a high-energy collision in a confining gauge theory. We consider a linearized approach in which two point particles traveling in an AdS-soliton background suddenly collide to form an object at rest (presumably a black hole for large enough center-of-mass energies). The resulting radiation exhibits the features expected in a theory with a mass gap: late-time power law tails of the form t^(-3/2), the failure of Huygens' principle and distortion of the wave pattern as it propagates. The energy spectrum is exponentially suppressed for frequencies smaller than the gauge theory mass gap. Consequently, we observe no memory effect in the gravitational waveforms. At larger frequencies the spectrum has an upward-stairway structure, which corresponds to the excitation of the tower of massive states in the confining gauge theory. We discuss the importance of phenomenological cutoffs to regularize the divergent spectrum, and the aspects of the full non-linear collision that are expected to be captured by our approach.


Introduction
The study of collisions and their outcomes is one of the most important ways of obtaining information about a theory and of testing it experimentally. This is true both in particle physics, where collision experiments have been dominant for a century now, and in gravitational physics, with the expected imminent detection of the gravitational radiation from collisions of black holes and neutron stars. The advent of gauge/gravity dualities brings about a merging of these two fields [1]: the collision of two high-energy particles in certain nonabelian gauge theories can be adequately described in terms of a dual gravitational collision in a higher-dimensional spacetime with negative curvature.
One phenomenon of current interest in this area is the collision at high energies of two objects (nuclei, nucleons, or partons) which, through the interactions of Quantum Chromodynamics (QCD), form a ball of quark-gluon plasma (QGP). Although the dual of QCD is not known, the analogous process in gauge theories with a gravity dual can be described via the collision of two objects of finite but small size that form a black hole in an asymptotically AdS spacetime. 1 The study of these collision processes is challenging because one must solve Einstein's equations in a dynamical setting, which generically must be done numerically. Several such studies have now been performed in cases in which the gauge theory is a Conformal Field Theory (CFT) [3,4,5,6]. 2 The goal of this paper is to give a first step towards extending this program to gravitational duals of confining gauge theories. For this purpose we will consider collisions in the so-called AdS-soliton [14,15].
One motivation for this extension is that a CFT has a continuous spectrum, so the result of the collision cannot be directly interpreted in terms of e.g. particle production. In contrast, we will see that, in a confining geometry, typical observables (e.g. the emitted radiation) admit an immediate particle interpretation. Another motivation is to explore the effects of confinement on the produced QGP. In a real heavy ion collision at RHIC or LHC the temperature of the produced QGP is roughly 2T c T 4T c , with T c ∼ Λ QCD the deconfinement temperature. This means that there is no hierarchical separation between the temperature of the fireball and the confinement scale, thus suggesting that the latter may play a role in the dynamics of the QGP. Similarly, in the range of temperatures above, the trace anomaly in QCD, which measures deviations from conformality, is still relatively sizable [16], again suggesting a possible role of Λ QCD .

The Zero-Frequency Limit framework
Our framework is very simple: we model the colliding objects as point particles moving along geodesics in a background spacetime, colliding instantaneously to form a single object at rest. The process amounts to specifying a conserved stress-energy tensor for point particles following these trajectories, and the gravitational field that they create is treated as a linearized perturbation of the background. Treating the collision in this approximation is well motivated, and this model is sometimes known in the literature as 'instantaneous collision framework' or 'Zero Frequency Limit' (ZFL) approximation [17,18,19,20,21,22]. The ZFL has been applied in a variety of contexts, including electromagnetism where it can be used to compute the electromagnetic radiation given away in β-decay (see for instance Chapter 15 in [23]). Wheeler used the ZFL to estimate the emission of gravitational and electromagnetic radiation from impulsive events [24]. In essence, we are reducing the full gravitational dynamics of the process to an effective theory of point particles interacting through a three-leg vertex, and then weakly coupling gravitational radiation to this system. This description is a reasonable one given what we know generically about black hole formation in similar collisions.
In asymptotically flat spacetimes, this simple approximation turns out to describe accurately all the main features of high-energy collisions of two equal-mass black holes [25,26,27,28]. These results are summarized in Figure 1, and refer to head-on collisions of two  [27]). The collision speed in the center-of-mass frame, β = v/c, is indicated in the legend. The energy spectrum is roughly flat (independent of frequency) up to the quasinormal mode (QNM) frequencies (marked by vertical lines), after which it decays exponentially. All quantities are normalized to the Arnowitt-Deser-Misner (ADM) mass of the system M ADM . The dashed horizontal lines are the ZFL prediction, obtained by a multipolar decomposition of (1.1).
equal-mass black holes, with center-of-mass energy parameterized by total energy M ADM and velocity v. The salient features of the ZFL analysis are: mass objects each with rest mass M/2, velocity v in the center-of-mass (CM) frame and Lorentz factor γ, the ZFL prediction [21,22] for the energy spectrum at an angle θ relative to the collision axis is (1.1) The independence of this spectrum on frequency ω follows from simple arguments: since we are working in the linearized approximation to gravity, the right-hand side of (1.1) must be proportional to M 2 . This fixes all the leading-order dependence on M . The absence of any other dimensionful parameter in the problem then forbids, on dimensional grounds, any possible dependence on ω. So the ZFL yields a flat spectrum, as shown in Fig. 1 for different CM velocities v. The nonlinear results are in good quantitative agreement and do show an approximately flat energy spectrum.
2. The need for an appropriate, physical cutoff. Because the spectrum is flat, estimates for the total radiated energy or time-domain signals formally diverge. The dependence on ω that would cutoff the spectrum is lost when we neglect non-linear effects and thus eliminate all the details of the interaction and the internal structure of the colliding and final objects. We can nevertheless reintroduce the cutoff in frequency (or momentum, via the dispersion relations) in a phenomenological way, which in asymptotically flat spacetimes is essentially uniquely determined (up to numerical, order-one factors). Although dimensional arguments do not fully fix the cutoff -besides the dimensionful scale M , there is a dimensionless parameter γ-we can expect that it is the size ∼ M γ of the final black hole that sets the cutoff: black holes absorb very efficiently frequencies that are larger than its lowest quasinormal mode (QNM) frequency, so we may expect that any higher frequencies will not be radiated away, and therefore ω cutoff ∼ ω QN M ∼ 1/(M γ) . (1.2) This implies that the frequency cutoff decreases as the velocity of the colliding particles increases, i.e. the quanta radiated are less energetic for larger CM energies -a characteristic property of collisions that involve black holes. However, since the ZFL spectrum (1.1) scales like γ 2 , the total radiated energy scales like γ and thus grows with v. Nonlinear simulations are in excellent agreement with this picture and show an exponential suppression of the spectrum for frequencies larger than the final black hole QNM frequency, as shown in Fig. 1.

3.
A "memory" effect in the signal, which is a consequence of the identity for the Fourier transformh(ω) of the time derivative of any metric perturbation h(t) (we omitted unimportant overall factors in the definition of the transform). Thus, the low-frequency spectrum depends exclusively on the asymptotic state of the colliding particles, which can be readily computed from their Coulomb gravitational fields. Because the energy spectrum is related toh(ω) via we immediately conclude that the energy spectrum at low-frequencies depends only on the asymptotic states [18,20,21,22,25].
4. Finally, the angular distribution of radiation is nearly isotropic at large collision energies (v → 1 in (1.1)), when higher multipoles become increasingly more relevant. This is also in agreement with nonlinear simulations.
The reason for the overall agreement of ZFL predictions with fully nonlinear simulations is not completely clear, one possibility being that nonlinearities are redshifted away. This simple model of particle collisions has provided useful benchmarking in the nonlinear simulations of asymptotically flat spacetimes [25,26,27,28]; we expect similar benefits in asymptotically anti-de Sitter (AdS) spacetimes where full-blown nonlinear evolutions are specially hard to perform.

The Zero-Frequency Limit in solitonic-AdS backgrounds
A collision in the AdS-soliton background is depicted in Fig. 2. These collisions differ from those in asymptotically flat spacetimes in several crucial respects. In particular, the process through which the initially-formed, highly-excited black hole radiates and relaxes to its equilibrium state is expected to be much more complex.
Scales in the collision and horizon evolution. In the Minkowski background, the entire collision and its evolution to a final state are characterized by the only scale in the problem, 3 namely the Schwarzschild radius ∼ M . In particular, this scale controls both the properties of the linearized field of point particles and the properties of the final equilibrium black hole. In contrast, in the AdS-soliton there are two additional scales, namely the AdS curvature radius L and the infrared length scale r 0 at which the geometry caps off smoothly. Actually, using coordinate reparametrizations they determine only one physical scale, which naturally may be taken as the gauge theory mass gap or confinement energy Λ QCD ∼ r 0 /L 2 . The presence of this scale besides the particle or black hole mass M implies that the linearized field of a point particle in a confining background need not give a good estimate of the size and shape of the horizon of an equilibrium black hole localized in the infrared. In fact, it is a very poor approximation to it when the mass is large, M Λ QCD [29]. The black hole is dual to a plasma ball [30] and is characterized by two very different scales: a thickness ∼ r 0 in the holographic direction, and a much larger proper extent ∼ L(M/Λ QCD ) 1/3 along the gauge-theory directions. Thus our approximation of the final state as a structureless point source is worse than in a Minkowski background.
Moreover, in a collision with M Λ QCD we can expect that the initially formed horizon will be largely insensitive to the scale r 0 . In this regime of energies the effects of the confinement scale can be neglected in this very early stage of the collision, but they will gradually appear in the relaxation to the final dual plasma ball. Thus we expect a richer evolution from collision to relaxation in a confining AdS background than in a Minkowski background, with different stages being characterized by different horizon scales.
ZFL spectrum cutoff. The linear approximation used in the ZFL neglects all the details of this process of horizon formation and relaxation, but some such information is nevertheless needed in order to specify the cutoff that renders finite the total radiated energy and timedomain signals. Here the differences between ZFL collisions in Minkowski and in AdS-soliton backgrounds can become significant in practice. In the former case, as we have seen, setting the frequency cutoff to be the final-state lowest QNM frequency is natural as this is the only scale in the problem. Moreover, in that background it does not matter whether we impose a frequency cutoff (the lowest QNM gives a characteristic time of horizon vibrations) or a momentum cutoff (implementing that a horizon only absorbs efficiently wavelengths shorter than its size).
The correct choice of a cutoff in a confining background is more convoluted. Attempting to find a cutoff on frequencies from the properties of the horizon formed in the collision is fraught with ambiguities: the horizon evolves through several complex stages and it is unclear to what extent these out-of-equilibrium horizons are well approximated by the properties of known stationary black holes. For instance, in a high-energy collision with M Λ QCD we might expect the initial horizon to be roughly similar to that of a large neutral black hole in global AdS, but even then it is not clear which kind of quasinormal modes would control the cutoff: 'fast' modes with frequencies ∝ M 1/4 , proportional to the black hole temperature, which correspond to the highest gauge-theory energies that the horizon presumably probes; or 'slow' modes, with much smaller frequencies ∼ 1/L, which approximate hydrodynamic modes for very large black holes. Moreover, the dual plasma ball that the system relaxes to also has 'slow' elastic modes associated to vibrations in the shape of the (dual) plasma ball. The complex time evolution of the horizon makes it unclear whether any of these QNMs can give a reliable frequency cutoff on the radiation.
The cutoff on the frequency of radiated AdS-soliton modes may be obtained more plausibly, via their dispersion relation, from a momentum cutoff. We would expect that the size of the horizon sets the largest wavelengths that it can efficiently absorb. For collision energies much larger than the confinement scale, at any stage the horizon will be larger (at least in some of its directions) than 1/Λ QCD , so it can absorb waves of momentum down to k cutoff Λ QCD . However, since the dispersion relation of the lowest AdS-soliton normal mode is ω 2 ≈ k 2 + Λ 2 QCD , the frequency is effectively cutoff by the mass gap. While these estimates seem reasonable, the correct choice requires truly non-linear information about the formation and evolution of the horizon, which at present is unknown to us. In order to deal with this uncertainty, we will not present our results in the time-domain, nor as an integrated total energy, but rather as frequency-domain signals, where the effect of the choice of different cutoffs is apparent and not obscured by integrating over the frequency spectrum. We will present some time-domain quantities only to illustrate their cutoff dependence.
It is important to note that in collisions at energies above the confinement scale, although it may be possible to neglect the confinement scale r 0 for certain aspects of the initial dynamics of the black hole formed, this scale is still crucial for the proper interpretation of the emitted radiation in terms of gauge theory particles, since it is responsible for the discreteness of the spectrum and the existence of a mass gap. These are universal, structure-independent features of the radiation produced in all such collisions that our model does capture.

Plan
The remainder of the paper is organized as follows. We begin by introducing the AdS-soliton geometry in Section 2, where we also collect its most relevant features for our study. In Section 3 we investigate, as a warm-up toy model, the scalar field radiation produced in the collision of two scalar-charged particles. This is a simpler problem that shares many features of the gravitational problem -except for the fact that in the former case the scalar charges in the initial particles simply add up to yield the final charge, whereas in the latter the kinetic energy of the colliding particles contributes to the final total charge (i.e. mass). Section 4 is concerned with the gravitational counterpart, the main objective of the paper. In this section we study the gravitational radiation resulting from colliding two point particles in the AdSsoliton background and we present results for the far-region behavior of the stress-energy tensor of the dual confining field theory, obtained via the AdS/CFT correspondence. Several technical details of the derivation are relegated to the Appendices. Finally, we summarize our results in Section 5.

The AdS-soliton background
As discussed in Section 1, we will be interested in colliding point particles in a gravitational dual to a four-dimensional gauge theory in a confined phase. The prototype for such a geometry is the six-dimensional AdS-soliton [14,15]. We will take this spacetime as a background on top of which we then consider a linear analysis of perturbations induced by the point particles.
The AdS-soliton is a vacuum solution of Einstein's equations with a negative cosmological constant, which is given in terms of the AdS curvature radius by −5/L 2 . It asymptotes to AdS with one of the spatial coordinates periodically identified, thus forming an S 1 . This circle smoothly shrinks to zero size at a finite radial coordinate and consequently the geometry is regular everywhere, without possessing an event horizon. The soliton metric is given by and y has periodicity We will group the coordinates as The S 1 shrinks to zero size as r → r 0 and the boundary lies at r → ∞. The gauge theory lives on Mink 1,3 × S 1 . As we will see, the Kaluza-Klein scale associated to the compact direction, also sets the confining scale (hence our choice of notation), reflecting the well-known fact that these two scales cannot be decoupled within the gravity approximation. 4 Of some interest for us is the eikonal limit of massless field propagation, which is generically associated to spacetime geodesics. Let us focus on radial geodesics along the holographic direction r and along the flat directions x i . These satisfy where = −1, 0 for timelike and null particles respectively, and E t is a conserved (dimensionless) energy parameter defined by The quantity defines the other conserved parameters and we set E y = 0. Note that dx i /dt = E x i /E t and therefore the projection of null geodesics onto flat (constant-r) slices follows straight lines. Equation (2.8) shows that c x represents the speed of light projected along the x-space. In particular, light propagates along constant-r slices with constant speed c x = 1, but a non vanishing component along the holographic direction implies c x < 1.
There are turning points at r = r 0 and at r = r i ≡ LE t 1 − c 2 x . In general the motion is bounded in the r−direction, and periodic. For null or high-energy timelike particles, the (coordinate) time it takes for a roundtrip from r = r 0 to r = ∞ and back is with Γ(x) a factorial-Gamma function. This implies a characteristic frequency of This frequency is also the frequency of high-energy timelike particles. On the opposite end we have low-energy particles, which oscillate between r = r i and r = r 0 with a period (for and a frequency Note that the oscillation period does not depend on the amplitude of the oscillation, in the small amplitude limit.

A toy model: scalar interactions
Although our final goal is to study gravity, let us first take a look at a simplified problem, that of scalar fields in the background (2.1). We will see later that this problem shares many common features with the more involved gravitational case. We will study the stability of the spacetime against scalar field perturbations, compute the field of static scalar charges, and finally let these collide. The action of our generic setup is whereŻ 2 = g abŻ aŻb ,Ż a = dZ a /ds, and we have adopted the somewhat arbitrary 1/8π normalization for the scalar. This theory describes a point particle of mass M and worldline Z a (s) minimally coupled with strength (charge) Q to a massless scalar field Φ. Note that the factor of −Ż 2 in the last term is necessary to make this term invariant under worldline reparametrizations. One important consequence of this is that the coupling to the scalar field vanishes in the ultra-relativistic limitŻ 2 → 0. We will see manifestations of this fact in our results below.

Stability and normal modes
We start by understanding the stability of the vacuum spacetime when the particle source is absent. The evolution of the scalar is then described by the massless Klein-Gordon (KG) equation This equation separates under the ansatz with n y = 0, 1, 2... due to periodicity y ∼ y + 2π/Λ QCD , and yields with k 2 = k 2 1 + k 2 2 + k 2 3 . Using a new variable ρ = r/r 0 the above can be written in manifestly dimensionless format, or equivalently where the dimensionless quantitỹ is the invariant four-dimensional mass ω 2 − k 2 measured in units of the confinement scale Λ 2 QCD .
Generically, the point ρ = 1 is a regular singular point of the ODE, and close to this point the dominant asymptotic behavior is of the form Ψ ∼ (ρ − 1) ±ny/2 . Specializing to n y = 0 modes, the equation simplifies to ρ(ρ 5 − 1) Ψ + (6ρ 5 − 1) Ψ + ρ 2ω2 Ψ = 0 , (3.9) At ρ = 1 the solutions behave as At infinity (for generic n y ) they behave as We require as boundary conditions that the solution be regular in the infrared, i.e. that C r 0 1 = 0, and that it be normalizable near the boundary, i.e. that C ∞ 1 = 0. Equation (3.6) together with these boundary conditions defines a Sturm-Liouville eigenvalue problem. The eigenvalues satisfyω 2 n > 0 and physically they characterize the mass spectrum of scalar excitations in the gauge theory in the limit in which their possible mixing with higher-spin excitations is neglected, since we have ignored their possible mixing with e.g. gravitational perturbations. We will denote the dimensionful eigenvalues in (3.4) as In the case n y = 0, which will be our focus later, the corresponding eigenfunctions Ψ n satisfy 13) and are orthonormal with respect to the scalar product ∞ r 0 dr r 2 L 2 Ψ n (r) Ψ m (r) = δ mn . (3.14) By direct integration, with two independent codes, we find the modes in the first column of Table 1. The first eigenfunctions are shown in Fig. 3. Note that the spectrum is discrete and gapped despite the fact that the radial direction is infinite, due to the fact that the AdS-soliton geometry acts like a 'box'. As we will see, by expanding any function of the radial direction in the complete basis provided by these eigenfunctions, one may 'Kaluza-Klein reduce' along the radial direction and reduce the problem to one dimension lower. The mass of the scalar mode with n = 0, n y = 0 determines the mass gap of the gauge theory in the scalar channel as   Table 1: Resonances for Klein-Gordon modes (first two columns) and for vector and scalar gravitational modes (last three columns). Hereω 2 ≡ L 4 (ω 2 − k 2 )/r 2 0 . Notice that for large overtone the spacing is roughly constant and equal toω n −ω n−1 ∼ 2.506 for all fields. As we show in Sec. 4, the second family of vector gravitational modes (vector II) is described by the same equation as in the Klein-Gordon case and therefore these two sets of modes exactly coincide. In the particle collisions we consider in Sec. 4.3, only the gravitational scalar and vector-II modes are excited. The spectrum of scalar gravitational modes agrees within numerical precision with the spectrum of scalar, 0 ++ glueballs found in [32] for the AdS 6 soliton. mode of Φ in the geometry, as it would imply an imaginary wavenumber. For very low frequencies the field is exponentially suppressed with distance; we will explicitly show this in Section 3.3. These results are similar to those of waveguides in classical electromagnetism, as is the physical setting [23].
In Table 1 we also present vector and scalar gravitational modes which are discussed in the next sections. Our results for scalar perturbations agree perfectly with those of Ref. [33] (see also [34]). For all families of modes and for large overtone number n, we find that the resonant frequencies are of the formω n ∼ω 0 + 2.506 n, withω 0 a field-dependent constant. Notice that at large overtone, where the frequency is large and where the eikonal limit is valid, the geometrical optics regime of the wave equation should go over to the geodesic equation. It is therefore pleasing to notice that the spacing of modes in this regime is exactly (to within numerical precision) the same as the one predicted by a geodesic analysis. This agreement can be made more formal by using the variable χ = 1/ρ, and the wavefunction   Table 1.
The wave equation (3.9) is now brought to the form A standard WKB analysis gives the asymptotic modes satisfying Dirichlet boundary conditions as [35]ω The quantity q is to be evaluated at largeω. We then get a spacing of A comparison with (2.10) gives finally a spacing equal to 2π/P roundtrip , and the asymptotic relation is established.

A static scalar charge in the AdS-soliton background
Let us consider now the solution of the KG equation with source at x i = X i , r = b described by the action (3.1). Notice that, when b > r 0 , this is in fact not a point particle, but rather a ring of matter that extends along the y-direction. Since this preserves the rotational symmetry along this direction, only the n y = 0 mode will get sourced. The equation of motion is By Fourier expanding the field, we get the following equation for the field of a static scalar charge, 24) or equivalently There are several distinct but equivalent ways of solving this equation. The standard procedure uses variation of parameters. Let us introduce two linearly independent solutions, Ψ I and Ψ II , of the homogeneous equations and their Wronskian W ≡ Ψ I Ψ II − Ψ II Ψ I . For our case with r 1 an integration constant. Now, all we have to do is define the homogeneous solutions Ψ I , Ψ II such that The quantity aω is generically a function ofω 2 , and we provide a detailed characterization of it in Appendix A. For the present, we are focusing on the case ω 2 = 0 and thus in the static case aω = aω(k). Note that Ψ I is regular at r = r 0 , that Ψ II is normalizable, and that Ψ II is singular at r = r 0 unless aω = 0. In this particular case Ψ II is a regular and normalizable solution of the homogeneous equation. Recall that such solutions exist only for real and strictly positive values of ω 2 − k 2 . Since for ω = 0 this translates into strictly negative values of k 2 , it follows that aω(k) has zeros at strictly imaginary values of k.
The solution can be written as , and we can write Let us rewrite the KG equation as Performing an integration on both sides from b − to b + we get Finally, This result is finite and continuous everywhere, except in the limit r → b → r 0 , since in this case we have (see Appendix A) and therefore For k = 0 the homogeneous equation can be solved exactly, with the result Thus aω(ω = 0) = −1. Furthermore, we find that at large k the function aω increases exponentially (in absolute value), aω →−∞ ∼ −6.4 × 10 −3 e −1.156ω . A solution of the inhomogeneous equation (3.24) can also be obtained by expanding Ψ in the normal modes Ψ n . Setting substituting into (3.24) and using (3.13) we get Multiplying both sides by Ψ m (r), integrating over r and using the orthonormality conditions (3.14) we obtain (3.40) We thus see that the coefficients in the expansion (3.38) are proportional to the support of the corresponding wave function at the location of the particle. We also see that these coefficients have poles at (imaginary) values of k determined by the mass spectrum of the normal modes.

An equivalent, matrix-valued Green function approach
Here we discuss an equivalent approach, based on Green function techniques for coupled systems of ordinary differential equations (see e.g. Ref. [36]). This approach is advantageous because, as we shall discuss later, it can be directly extended to the gravitational case. Any second-order system of coupled ODEs can be written in a first-order form, where Y and S are generically n dimensional vectors and V is a n × n matrix. We define the n × n matrix X whose mth column contains the mth solution of the homogeneous system i , where the j index denotes a solution of the homogeneous system and i is the vector index. The matrix X constructed in such a way is also a solution of the associated homogeneous system, in the sense that In order to solve (3.41), we impose the ansatz Y = XΞ, where Ξ is a vector to be determined. Substituting this ansatz into the inhomogeneous system and using Eq. (3.42) we find 43) and the solution to (3.41) formally reads Let us now apply this method to Eq. (3.25). In this case, Y ≡ (Ψ, Ψ ) and Furthermore, using the same notation as in the previous section, the matrix of the homogeneous system reads Evaluating the first component of Eq. (3.44) we obtain (3.47) where the limits of integration were chosen in order to have the correct boundary conditions. Finally, evaluating the expression above when r < b and when r > b, we recover the same result as in Eq. (3.34). In this and in some subsequent plots we work with dimensionless units by setting r 0 = L = 1 and we have chosen, without loss of generality, Q = 1.

Yukawa-like potential at large distances
Let us now consider a point particle located at r = b = r 0 and X i = 0 and compute explicitly the large-R behavior. From equation (3.29), we get that the Wronskian between Ψ I , Ψ II is W = aω(k)/(r − r 0 ). Using (3.34) we then find that for r > r 0 , where the function aω depends only on k = k 2 1 + k 2 2 + k 2 3 . At r ∼ r 0 , recall we have from Eq. (3.29) that Ψ II ∼ cω + aω log(r − r 0 ) and so we get Thus, at leading order, the solution is localized in the radial direction. It is not our priority here, but it would be interesting to extract the coefficient cω defined in (3.29). This coefficient would presumably dictate the small R, small r 0 dependence of the field away from the source location.
Let us focus instead on the r → ∞ regime, where Ψ II → 1/r 5 . In this case, we get The integral above can be simplified, and we have used the fact that aω(k) is an even function of k. Note that the values of k corresponding to the normal modes of the system correspond to poles of the integral in the complex plane. A contour plot of the scalar field strength Φ(r, R) in the R − r plane is shown in Fig. 4, whereas the the integral (3.51) is shown in Fig. 5 as a function of R. Note that the potential is finite as R → 0 and it decays exponentially for large values of R.
Our numerical results are consistent with a Yukawa-like decay, I ∼ e −µR /R (cf. right panel of Fig. 5). Indeed, the integral above can also be computed by deforming the integration contour in the complex plane and using the residue theorem. We shall discuss this technique in detail in the following sections; here we just give the final result. Using the fit (A.11) to approximate the function aω(k) close to the poles, we obtain a sum of Yukawa-like potentials, 5 where µ n are the modes listed in the first column of Table 1 and the coefficientsc n read (cf. where P = 2 × 2.506. In Sec. 3.4 we shall confirm these results by obtaining the leading terms of Eq. (3.52) via a completely independent approach. In Fig. 5 we compare the numerical results with the superposition of normal mode solutions (3.52) (cf. also Section 3.4 below) truncated at N = 5. An asymptotic analysis of the integral I in (3.51) at large R confirms this Yukawa-like behavior: using a stationary phase approach, the relevant function to study is f (k) = ka −1 ω at k = 0 [37]. The asymptotic behavior at large R is strongly dependent on the behavior of the even derivatives of f (k) and in particular the existence of power-law behavior seems to be connected to non-zero even derivatives of this function. As we show in Appendix A, the function aω is an even function of its argument for small enough argument. Together with our numerical data, which is consistent with f (2n) = 0 at least for n = 0, 2, 4, this kind of asymptotic analysis also predicts what we find numerically, i.e. a Yukawa-like suppression.
In fact, this behavior can be proven exactly by using the solution in the form (3.38) with the coefficients (3.40). Through the Fourier transform (3.23), the poles in the coefficients produce precisely the Yukawa-like terms: where d n is proportional to the k-independent part of c n . Note that this expression is valid for all values of r, R.
The Yukawa-like decay is not a peculiarity of having a point source. Our results generalize to any compact distribution along the flat directions x i , as long as the source is localized in the holographic direction. In fact, the large distance behavior is controlled by the small-ω asymptotics, which is independent of how the source is distributed, as long as it is compact. In subsection 3.4 below, we exhibit vacuum solutions not localized along the holographic direction and that also display a Yukawa-like decay.

High-energy collisions of point particles
We now consider the collision of two scalar particles of equal mass moving towards each other at velocity v along the flat directions, depicted in Fig. 2. For concreteness, we focus on two particles flying at each other along the x 1 axis and colliding at t = 0. Thus, we model the problem as where γ = 1/ √ 1 − v 2 is the relativistic boost factor, Θ(t) is the Heaviside function, and BH will be defined momentarily. This is the Instantaneous Collision Framework or Zero Frequency Limit (ZFL) approximation described in the Introduction. It is well-motivated for the high-energy collision of two objects and is known to work well to describe classical processes in electromagnetism [23,22]. Two objects flying at close to the speed of light barely feel each other's field, and therefore the interaction takes place right at the moment of collision. Since we are eventually trying to describe the formation of a single black hole from the collision of two objects, we let the final particle be at rest, as described by the second term on the right-hand side of (3.55).
Notice that already at this level a choice of the final state is crucial: if the final state is charge-conserving then BH = 1. If instead the final state is a black hole and black holes in this theory continue to have no hair, then the scalar charge of the final black hole is presumably zero, and BH = 0 in this case. This technical detail yields a difference between a radiation output that scales as γ 0 if BH = 1 and a radiation that scales γ −1 if BH = 0. The physical intuition behind this is that the radiation is the dislocation in the field created by the change in the source at t = 0. For large γ this change is of order unity if BH = 1 and of order 1/γ if Let us proceed by Fourier analyzing the fields. We expand any function Z as In Fourier space, (3.55) yields and we have exhibited the appropriate ε-prescription, which we may not show explicitly in subsequent expressions.
With the same procedure as before, we define two homogeneous solutions Ψ I , Ψ II by Eqs. (3.27-3.29). The solution to the inhomogeneous problem can then be written as (3.61)

Reduction to four dimensions
We expand the five-dimensional field in the basis of normal modes as (3.62) Substituting into (3.55) and using (3.13), eq. (3.55) becomes Multiplying by (r 2 /L 2 )Ψ m (r), integrating over r and using (3.14) we arrive at an infinite set of independent equations, one for each mode: Thus the five-dimensional problem reduces to an infinite set of identical four-dimensional problems, each of them consisting of the determination of the massive scalar field generated by a source proportional to J(t, x i ). Each of these problems is a classical bremsstrahlung problem in which two particles moving in opposite directions collide and come to a complete stop, thus emitting radiation into the massive scalar field that they couple to.

Cutoffs
As we discussed in the introduction, in our approximation the spectrum needs to be cut off at high frequencies in order to yield sensible results for the total radiated energy and timedomain signals. While in the case of gravitational interactions one expects that non-linear effects, through black hole formation, dynamically introduce such a cutoff, in the case of scalar interactions the origin of the cutoff is less well defined; it would naturally depend on properties of the final object such as, possibly, its size scale. If this is larger than 1/Λ QCD , the argument in sec. 1.2 extends to this case and the frequency cutoff is set by Λ QCD . At any rate, these details are of not of much interest to us as we are taking this calculation as a toy model for gravitational interactions and, if needed, we may simply borrow the cutoff from a gravitational estimate.

High-energy collisions, conserved scalar charge
In the large-velocity limit and assuming that the final state has conserved charge ( BH = 1), the source function S does not depend on k 1 , and the process is spherically symmetric in the flat directions at leading order. This is a manifestation of the fact that the source for the scalar field vanishes in the ultra relativistic limit, as mentioned above. Thus in this limit one is left with the spherically symmetric field sourced by the 'sudden appearance' at t = 0 of a point-like source at rest, namely by the second line of (3.55). For large holographic coordinate r the field takes the form (see Eq. (3.50)) Our final goal is to extract the stress-energy tensor of the dual gauge theory from the asymptotic behavior of the metric perturbations; we shall tackle that problem in Section 4. For now we are considering only a massless scalar field on a non-dynamical AdS-soliton geometry. To complete this warm-up exercise we can determine the expectation value of the dual scalar operator, which is proportional to the coefficient of the r −5 term [38]: (3.67) We will use contour integration to perform the integral, and for that we employ an extension of the standard approach to compute massive field propagators [39,40]. We start by estimating the residues of 1/aω. Using the fit (A.11), we find close to the poles of this function that Let us focus on the k-integral, which has poles at where has the same sign as . Because we intend to close the contour on the upper half plane, only poles in the upper half plane contribute. These are located at It is gratifying to recover the expected result that frequencies smaller than the effective mass are exponentially suppressed at large distances R. On the other hand, for frequencies above the effective mass the field displays oscillatory behavior for large R. The wavelength of these modes decreases as ω increases but each time the frequency crosses above a modeω n a new oscillatory term appears, with a wavelength that is shorter than the previous mode.
In Fig. 6 we show the quantity | O(ω, R) |ωR (modulo a constant factor proportional to Q) as a function of the frequency. The different effective masses show up as plateaux in the spectrum, and the zero-frequency limit of the spectrum vanishes at large radii. This is a universal solid prediction coming from this model.
In Fig. 7 we show the full dependence of | O(ω, R) | on ω and R. As expected, the details of O depend on the number N of modes included in Eq. (3.73). As N increases, the waveform displays a complicated behavior due to the superposition of several massive modes. The various cutoffs for ω > ω n are also evident in the right panel of Fig. 7.
With the frequency-domain quantities under control, we can Fourier-transform back to the time-domain and discuss the time evolution of the scalar operator. Numerically, this is achieved by evaluating the inverse-Fourier amplitude of a generic functionψ(ω): whereψ is evaluated at a fixed spatial position, ω k = k∆ω with 2∆ω = ω max /(N − 1) and we assume a frequency domain [−ω max , ω max ] discretized in N equidistant points. The resolution in time is given by 2π/ω max , so that the larger the frequency domain the more refined is the resolution of the time evolution.
In Fig. 8 we show several snapshots of O in the time domain as a function of R. Several pulses (corresponding to modes with different masses) propagate with different velocities. As expected from causality, the waveform must vanish when R > t. This is consistent with our results to within our numerical accuracy. Note that the details of the waveform depend on the frequency cutoff, but the qualitative behavior is generic. The larger the frequency cutoff, the larger the number of massive modes that can be excited and the waveform displays some beating effects. The full time-domain dependence of | O(t, R) | is shown in Fig. 9 for different cutoff frequencies.
Finally, there is a nontrivial dispersion relation, and different frequencies propagate at different speeds. In other words, Huygens' principle is not satisfied and there is propagation inside the entire light cone [41,42]. This peculiarity gives rise to a wake behind the main pulse, which dies off at late times as [41] O where O static is the static, Yukawa-like potential to which O asymptotes at late times. The t −3/2 fall-off can be proven analytically from the properties of the retarded Green's function of a massive scalar field in four dimensions. This takes the form The delta function only contributes on the light-cone, whereas the Bessel function contributes inside the light-cone. Because of this, the field generated by a particle at a point p is the integral of the Green's function along the world line of the particle from the remote past to the latest time t ret (t) from which the particle could causally affect p. For a particle that has been sitting at x = 0 forever, we can write the resulting field at time t schematically as This of course yields the static Yukawa potential. In contrast, the first term on the left-hand side of (3.75) is only sourced from t = 0, so the difference in that equation is We thus see that for t → ∞ this difference is generated at times that are in the far past of the point of observation. Consequently, they are controlled by the fall-off at large σ of the Bessel function, which is  Figure 10 shows that even when the signal includes several modes, the tail of the waveform in time is precisely t −3/2 as predicted by the formula above.

High-energy collisions, final state with no scalar charge
In this case, BH = 0, and subleading terms have to be taken into account. At large distance r and for b ∼ r 0 we have where we defined k 2 ρ = k 2 2 + k 2 3 , and where E t = γr 0 /L. Finally, introducing cylindrical coordinates 6 we find, for large r, with J 0 (x) a Bessel function of the first kind [43]. Using Eq. (3.68), close to the poles we have aω ∼ − f n (k 1 − k 1,n )k 1,n L 4 ω n r 2 0 , (3.84) 6 Here, ρ = x 2 2 + x 2 3 should not be confused with the dimensionless holographic coordinate r/r0 used in Sec. 3.1.
By first integrating over k 1 we get where we have defined andãω is computed for k 2 1 = ω 2 /v 2 . Therefore Ψ(ω, x 1 , k ρ ) has poles at coming from the first term of the equation above, and at  where and H (1) n is the Hankel function of the first kind. Let us solve the integral above in the ultrarelativistic limit, v → 1. We split it into two contributions, I 2 = I n the branch points are on the real axis. In this case we can still integrate numerically, 7 but the branching points must be suitably excluded from the integration domain. The quantity ω ρ 2 + x 2 1 O(ω, x 1 , ρ) obtained by integrating numerically and summing the two contributions is shown in Fig. 11 (modulo a coefficient proportional to Q/E t ).
We find the same qualitative features observed in the spherically symmetric case. When the frequency is smaller than the fundamental mode, ω < 4.062r 0 /L 2 , the spectrum is exponentially suppressed at large ρ. As ω increases, several mass barriers can be overcome and

Sources extended in the holographic direction
In the previous sections, we constructed the field of point-like scalar charges, showed it decays exponentially fast in the gauge theory directions and collided them along one of the flat directions. What are the effects of finite-size on the previous results? To investigate this, we now construct a general class of solutions which are not localized in the holographic direction, and which still display Yukawa-type asymptotics. The formalism below is quite general and can handle any finite sized object in the holographic direction, reducing the problem to a Minkowski space evolution of massive fields.
Up to now we have kept factors of r 0 and L in most of the equations. In this section and in the remainder of the paper we will set r 0 = L = 1 to reduce cluttering of the equations. This can always be accomplished by appropriately rescaling the holographic coordinate, together with the fields. However, we shall explicitly reinstate such factors in the main results.
Consider then the Klein-Gordon equation with the general source where for the moment ρ(r) is an arbitrary regular function. If we look for spherically symmetric solutions by introducing a radial coordinate R = x 2 1 + x 2 2 + x 2 3 we get the following Separable solutions to this problem exist, and we can study them by using the following decomposition: ρ(r) = n a n Ψ n (r) , where Ψ n are the vacuum eigenfunctions studied in Section 3.1, with eigenvalueω 2 n = −k 2 n . We then get the following ODE for Z n (R), This is nothing but the equation for the field of a point-like particle coupled to a (massive) Klein-Gordon field in flat spacetime, whose solutions have the classical Yukawa-like form: Thus, a generic distribution in the holographic direction can be understood as the sum of the field generated by point-like particles in a Minkowski background and carrying a massive interaction. In the generic case, the field is a superposition of such solutions and the coefficients a n are evaluated as an overlap of different eigenfunctions with weight ∼ r 2 , a n = ∞ 1 dr r 2 ρ(r)Ψ n (r) .
(3.100) Suppose as a first example that a n = δ n n 0 , so that the stress-tensor profile in the holographic direction coincides with one of the normal modes of the field. We then get and at large distances we find where η n 0 is a normalization constant, which for the fundamental mode is η 0 ≈ 3.4. Let us recover as a final example the point particle results within this approach. In this case, ρ = δ(r − 1) and we find for the fundamental mode contribution (an infinite tower of modes contribute to the point particle field) This can be evaluated to be, at large holographic distances, which agrees to within 0.01% with a fit of our numerical results for the point particle calculation (which yields η 0 ≈ 12.213, to be compared with the analytical approximate results, η 0 = c 0 /(10π 2 ) ≈ 11.516, cf. Eqs. (3.50)-(3.53)). We also get c 1 /c 0 ∼ −5.4, in rough agreement with Eq. (3.53). 8 It is interesting to note that at large distances the field generated by the same total charge Q does depend on the charge distribution. In fact, the field generated by a point particle is roughly 3 times stronger than that created by a smooth ρ(r) distribution identical to the fundamental mode.
Collisions of these non-pointlike configurations can also be studied with well-known methods. Using the same notation as above, our equation now becomes a flat-space massive field equation: For BH = 1 and in the E t → ∞ limit, the solution in the space domain reads Generically, the collision of extended particles is quantitatively different but qualitatively identical to the collision of point particles. The output can be quantitatively the same by correcting only the static profile. For instance, the collision of a point particle with a spectrum cut at the fundamental mode results in a spectrum 3 times larger than the collision of a smooth extended distribution along r, with ρ(r) = ψ 0 and with the same total charge.

Gravitational interactions
We now turn to the main interest of this paper: the study of gravitational perturbations resulting from the head-on collision of two point particles in the AdS-soliton background. We tackle this problem step by step, as we did for the scalar toy model, first addressing the normal modes of the spacetime, then investigating the gravitational field created by a static particle and finally considering the collision process. From our results we will be able to infer the behavior of the stress-energy tensor of the dual gauge theory. Recall we are now setting r 0 = L = 1 for simplicity. We will explicitly display such factors only in the main results.

Stability and normal modes of the gravitational waveguide
Let us start by studying the gravitational normal modes of the AdS-soliton. Gravitational perturbations have more degrees of freedom, but we will be mainly interested in perturbations that keep some of the symmetries of the background intact. In particular, our goal here is not to perform a full perturbative decomposition of the gravitational field, so we now focus on the type of perturbations which are more directly relevant for the physics we wish to understand. In Appendix B we show the existence of a special type of vector-like gravitational perturbations, which are not excited by colliding objects head-on. These vector-type perturbations, are in principle excited in other, more generic situations; we show in the appendix that their spectrum shares the same main features as the ones we discuss below.
We focus on a subset of gravitational perturbations, appropriate for the symmetries we want to consider. In (t, r, x 1 , x 2 , x 3 , y) coordinates the metric reads where the perturbation quantities are defined as Here, an integral over ω and k i is implicit, as well as a summation over i = 1, 2, 3. We have singled out the coordinate x 1 to be aligned with the collision axis and for notational convenience we are setting x 1 ≡ x. The transverse directions x 2 and x 3 are on an equal footing and will be denoted indistinctly by x ⊥ . In the configurations that we consider the stress-energy tensor of the particles does not have any components T tx ⊥ . Thus, when working in transverse gauge, the components h tx ⊥ obey homogeneous equations with trivial boundary conditions and must vanish. From the viewpoint of the gauge theory this may seem counterintuitive, since we expect that the collision generates radiation with momentum in the x ⊥ directions. We will see that this apparent puzzle is actually resolved by the fact that the change of coordinates from (4.1) (which are a natural gauge choice from the bulk viewpoint) to Fefferman-Graham coordinates (appropriate for boundary observables) generates the expected components T tx 2 , T tx 3 of the boundary stress tensor. Inserting the ansatz (4.2) into the Einstein equations we get four algebraic equations for h tx , h yy , h ⊥ and h tr : In addition, the function h rx satisfies a homogeneous second order equation. We completely fix the gauge by requiring h rx = 0. Although the denominator of h tr above vanishes as r → 1, it is easy to show that in the same limit the numerator vanishes with the same power of (r − 1), cf. Table 2 in Appendix E for detail. 9 Therefore, h tr is regular at r ∼ 1 if the other functions are also regular. Finally, we get three dynamical equations for h xx , h rr and h tt . The latter can be simplified by introducing two new variables z − (r) and z + (r) such that (4.7) The perturbation equations in these new variables read (4.10) Here, we used the identity r 2 F = 2(3F − rF ) to avoid the explicit appearance of the second derivative of the metric function F (r). It is worth noting that our ansatz (4.2) includes, but it is not limited to, gravitational scalar modes. The latter are defined as those that, in their rest frame, transform as scalars under the little group SO(3). This is equivalent to the condition h xx = h ⊥ at k 2 = 0 in the ansatz (4.2). Using this condition and the perturbation equations above, it is easy to show that Eq. (4.8) is identically satisfied. Therefore, the gravitational scalar sector is entirely described by the system (4.9)-(4.10). On the other hand, Eq. (4.8) is decoupled from the other two and it describes a subsector of the gravitational vector perturbations. Notice that Eq. (4.8) is equivalent to the scalar equation (3.9), so that this subset of gravitatonal vector modes coincides with the scalar spectrum, as shown in Table 1. In Table 1 we refer to these modes as "vector II" family, to distinguish them from the "vector I" family presented in Appendix B which, however, is not excited in the collision discussed in Sec. 4.3.
Close to the bottom r = 1, we impose regularity of the perturbation functions: where the coefficients A can all be written in terms of just two independent parameters, namely A  can all be expressed in terms of C where the expansion coefficients A ∞ can be written in terms of four independent parameters and the coefficients C (j) ∞ can be written in terms of two independent parameters. We guarantee that all metric perturbations decay at infinity by fixing By imposing the conditions above, the asymptotic behavior of the metric functions reads , (a, b) = (t, r) (4.14) and the explicit form of the constants A ij and B ij is given in Appendix F.1. There we show that the large-distance behavior only depends on the three parameters We have integrated Eqs. (4.9)-(4.10) imposing the expansion (4.11) at r = 1 and requiring Eqs. (4.13) at infinity. With these boundary conditions, all metric components are guaranteed to be regular at r = 1 and to vanish as r → ∞. We have searched for the eigenvalues using two different methods, one of which we now describe. An alternative method based on Frobenius expansions is outlined in Appendix C.  Figure 13: Eigenfunctions corresponding to the fundamental mode and the first few overtones of scalar-type gravitational perturbations (left panels) and of vector-type II gravitational perturbations (right panels) as listed in Table 1.
The most efficient method is a standard technique to deal with matrix-valued eigenvalue problems (cf. Ref. [44] for a review). First, we perform two integrations starting from r = 1 with (A b ), respectively. The latter also correspond to two sets of solutions, (h I rr , z I + ) and (h II rr , z II + ). Finally, the eigenfrequencỹ ω is obtained by searching for the roots of det S. We find a discrete set of modes, which are listed in Table 1. Curiously, the modes of the system (4.9)-(4.10) include, but they are not limited to, the modes of Eq. (4.8). We stress that the latter coincide with the scalar modes previously discussed.
In Fig. 13 we show the eigenfunctions corresponding to the first three gravitational scalartype and first three vector-II type modes in Table 1. Within this direct integration method, the eigenfunctions are defined as are constants obtained after the determinant of the matrix S in Eq. (4.17) has been minimized. By construction, since det S = 0, the functions above are automatically eigenfunctions. In Fig. 13, the classification of modes into two different families is manifest. Finally, we show in Appendix D that ω = k = 0 is not a regular solution of the problem, and therefore that no zero modes exist in this background. In the gauge theory this means that the glueball spectrum is gapped, as expected.

A static point particle in the AdS-soliton background
In this section, we investigate the gravitational field generated by a static point-like source located at the tip r = r 0 = 1. In (t, r, x 1 , x 2 , x 3 , y) coordinates the metric reads as in Eqs. (4.1)-(4.2) but with all off-diagonal terms set to zero, in addition to ω = 0 and h ⊥ = h xx . Inserting this ansatz into the Einstein equations we get two coupled, second-order differential equations for h xx and h rr : whereT tt (r) is the Fourier transform (defined as in Eq. (3.58)) of the only nonvanishing component of T µν ,T tt (r) = µ δ(r − 1). Note that the Einstein equations allow a source term of this form, which is consistent with the fact that in this spacetime a point-particle can be static only at the tip r = 1. The remaining perturbation functions h yy and h tt are algebraically related to h xx , h rr , and their derivatives: The series expansions in Table 2 then determine the large (holographic) distance behavior: 10 The dimensionless functions a k (k 2 ), b k (k 2 ), c k (k 2 ) and d k (k 2 ) are related to the behavior of the gravitational perturbations at r ∼ r 0 . They can be constructed by a numerical integration of the homogeneous system, cf. Appendix E for details. Finally, in the space domain at leading order in r we obtain  Table 1. We consider the fundamental mode and the first three overtones including both families listed in Table 1. Right: same for the quantity R I(R).
so that all the information about the metric perturbations is encoded in the following integral The functions b k and D k = a k d k − b k c k have behaviors qualitatively similar to that of the function aω in the negative x−axis shown in Fig. 18. In order to evaluate the integral (4.28), we proceed as in the scalar case: The integral above is shown in Fig. 14, where we again compare the numerical results with a superposition of normal-mode solutions analogous to Eq. (3.52), but where now µ n are only the first four gravitational modes listed in Table 1 (we considered the fundamental mode plus three overtones including both scalar-type and vector-II type modes).

High-energy collision of particles
In this section, we compute the linear gravitational emission during the head-on collision of two point-particles with mass m boosted with speed v, each following a straight geodesic along the x 1 direction in the hyperplane defined by r = r 0 . We consider Einstein's equations with the stress-energy tensor where γ = 1/ √ 1 − v 2 and, at first order, the mass of the final (static) object is equal to the total energy in the initial system, M = 2m L r 0 γ. Accordingly, the spacetime velocity vectors of the three particles are 0, 0, 0, 0, 0) . (4.31) Recall we are using coordinates x µ = (t, r, x 1 , x 2 , x 3 , y) but, to avoid cluttering subsequent formulas, we shall adopt the notation x 1 ≡ x as in Section 4.1. In the final results we revert to the original notation.
The Fourier transform, defined as in Eq. (3.58), of the stress-energy tensor (with covariant indices) may be expressed as follows: Note two qualitative differences with respect to the scalar case: (i) because kinetic energy gravitates, the stress-energy tensor depends on v also in the ultrarelativistic limit v → 1 through the γ term; (ii)T µν explicitly depends on k 1 even in the ultrarelativistic limit, so that the source term is not spherically symmetric in the k-space. This is analogous to the case of black hole formation with no scalar charge, BH = 0, discussed in Section 3.3.4. Note also that the nonvanishing components ofT µν are related to each other byT xx = ω 2T tt /k 2 1 and T tx = −k 1Txx /ω. By using these relations, in the following we shall write the perturbation equations in terms ofT tt only.
A consistent ansatz for the metric perturbation of two point particles boosted along the x 1 direction is given by Eq. (4.2). The metric perturbations h tx , h yy , h ⊥ and h tr read as in Eqs. (4.3)-(4.6), whereas the dynamical variables h xx , h tt and h rr satisfy an inhomogeneous system of equations. The latter takes a simpler form after introducing the functions z ± as defined in Eq. (4.7). We obtain one decoupled inhomogeneous equation for z − , and a system of two coupled inhomogeneous equations for z + and h rr , (4.36) Note that no source term appears in Eq. (4.36), which reads as in the vacuum case. Our final goal is to compute the expectation value of the holographic stress-energy tensor in the dual theory. This is computed explicitly in Appendix F.3 and the entire computation is detailed in Appendix F. Here we only give the final result. In Fourier space and in the v → 1 limit, we get and αω, βω, γω and δω are related to the behavior near r 0 of the metric functions (see Table  3). The other five nonvanishing components, T x 3 x 3 , T tx 3 , T x 1 x 2 , T x 1 x 3 and T x 2 x 3 , can be obtained by using symmetry arguments and the tracelessness and divergence-free conditions where η mn stands for the Minkowski metric in the five-dimensional space covered by coordinates x m = (t, x 1 , x 2 , x 3 , y). As a check on our calculations, we have computed all components of T mn and checked that the conditions above are satisfied by virtue of the Einstein equations.
A relevant quantity is the energy flux across a sphere of radius R, where we have introduced spherical coordinates (R = x 2 1 + x 2 2 + x 2 3 , θ, ϕ) and T tR is the R−t component of the holographic stress-energy tensor. By performing a change of coordinates, the latter reads where, by symmetry, T tx 3 can be obtained from Eq. (4.41) by replacing k 2 ↔ k 3 . Even though T tx 2 and T tx 3 have no particular symmetry, one can show that T tR in Fourier space reads where the functions have cylindrical symmetry in k-space, and recall k 2 ρ ≡ k 2 2 + k 2 3 .

Numerical results for the stress-energy tensor of the dual theory
The results (4.37)-(4.42) were obtained in Fourier space. By inverse-Fourier transforming from k-space, we can obtain the spatial dependence of the operators. Note that Eqs. (4.37)-(4.39) and Eq. (4.47) explicitly depend on k 2 and k 2 1 only and the inverse Fourier transform can be performed in cylindrical coordinates as shown in the scalar case (cf. Eq. (3.83)). The term (4.40) only depends on k 2 and its transform can be evaluated as in the spherically symmetric scalar case (cf. Eq. (3.66)). Finally, the last two terms (4.41) and (4.42) are less symmetric because they explicitly depend on (k 1 , k 2 , k 3 ).
As in the scalar case, some of the integrals above have to be performed numerically. In Fig. 15 we show the operator |T yy (ω, R)|, which is spherically symmetric and qualitatively similar to the scalar case shown in Fig. 6. In Fig. 16 we show the operators that are cylindrically symmetric. In this case the spectrum is qualitatively similar to the cylindrically symmetric emission in the scalar case, cf. Fig. 11.
Finally, we can compute the energy flux. Let us first compute the inverse-Fourier transform of T tR . From Eq. (4.47), we obtain where cylindrical coordinates are related to spherical ones via ρ = R sin θ, x 1 = R cos θ and there is no explicit dependence on ϕ. Therefore, the energy flux reads where, using the symmetries of the problem, we integrate over half of the cos θ-space. The energy flux as a function of the frequency is shown in Fig. 17. As expected, the flux is vanishing for frequencies smaller than the first normal mode, ω <ω 1 ∼ 2.52.  Figure 16: The operators |T tt (ω, x 1 , ρ)| x 2 1 + ρ 2 , |T xx (ω, x 1 , ρ)| x 2 1 + ρ 2 and |T tx (ω, x 1 , ρ)|(x 2 1 + ρ 2 ) (modulo a coefficient proportional to ξ) for N = 11. Left panels: x 1 = 1, Right panels: x 1 = 10. The spectrum is qualitatively similar to the cylindrically symmetric emission in the scalar case, cf. Fig. 11.

Discussion and Conclusions
Gravity-dominated high-energy collisions are a fascinating topic: the efficiency for gravitational wave emission is huge, and these processes typically give rise to the largest known luminosities. In addition, fine-tuned collisions of this kind provide useful tests of Cosmic Censorship. In four-dimensional, asymptotically flat spacetimes the simulations of such events took several decades to perform and turned out to be 'mostly linear' in that waveforms are smooth and simple perturbative models capture most of the physics [25,26,27,45].
By contrast, full non-linear simulations of collisions in AdS spacetime are in their infancy [3,4,5,6]. Validation of any such simulations requires benchmarking with perturbative results, either for the final ringdown stage or for intermediate stages of the process. We have explored a simple and compelling model for such collisions. Some of the important physical observables, such as total radiated energy and time-dependence of the stress-energy tensor depend quantitatively on the magnitude of the cutoff, which in turn depends sensitively on the (unknown) final state. Nevertheless, we also obtained what we expect are universal features to be seen in any simulation and experiments: Yukawa-type potentials for static particles and power-law decay of perturbations at late times, characteristic of massive fields. The energy distribution of the particles produced during such events is also cutoff-independent.
As shown in Figs. 6, 11, 15 and 16, the spectrum is exponentially suppressed for frequencies smaller than the fundamental mode of the AdS-soliton, i.e. smaller than the mass gap in the gauge theory. Using Eq. (1.3), this result implies that the waveform shows no memory effect in this spacetime. This property is most likely related to the fact that the AdS boundaries are timelike and can be reached by the emitted radiation in a finite time.
Furthermore, for larger frequencies the spectrum shows a peculiar upward-stairway structure, which is formed by various plateaux corresponding to the excitation of various normal modes with increasing overtone number. From the dual theory perspective, as more massive states becomes available in the confining gauge theory, the energy density dE/dω for a given frequency grows monotonically with the energy. For example, for the scalar charge-conserving collisions studied in Section 3.3.3, Eq. (3.73) predicts a distribution dE/dω ∼ ω 2.53 at large energies. For the gravitational emission, our numerical results are less accurate, but the behavior is also consistent with a w Υ dependence, with 2.5 Υ 3 for T yy , T XX and T tX .
These qualitative features are robust because they follow from the fact that the fields in a confining geometry can be decomposed into a discrete set of massive four-dimensional fields propagating in Minkowski spacetime. From the gauge theory viewpoint these are just the different glueball states in the theory. This reduced description in Minkowski spacetime allows us to give a very simple and intuitive picture of the radiation field coming out of the collision, along the same lines as in electromagnetism. Imagine that we regularize the problem by assuming that the two particles come to rest in a small but finite amount of time δt, so that they start slowing down at t = 0 and come to a full stop at t = δt. After this time, the solution has three regions. The field sufficiently far from the particles, at R t, is just the sum of the two boosted fields created by the incoming point particles. In the opposite limit, R t, the field is the spherically symmetric solution created by the resulting particle at rest. In between these regions there is a thin shell of thickness ∼ δt in which the field smoothly connects these two solutions. This 'dislocation' is the gravitational wave, and it is not spherically symmetric because it must connect a spherically symmetric solution (at small R) to a non-spherically symmetric one (at large R).
Note that the only difference between the present situation and that in electromagnetism is that in our case the effective four-dimensional fields are massive. This is not an essential difference though, since the front wave of a massive field still propagates at the speed of light even if it is followed by slower modes. As in electromagnetism, the regularized picture illustrates the fact that the radiated energy diverges in the limit δt → 0 and hence explains the need for a cutoff. Indeed, the derivative of the field across the shell of thickness δt is of order v/δt, and hence the total radiated energy scales as δt × (v/δt) 2 .
Although our model has no internal information about the colliding objects, presumably this is not a serious limitation if the collision is sufficiently energetic: horizon formation will cloak any multipolar strucutre of the colliding particles and presumably a point-particle approximation is just as good as any other [28,46,47].
A more important limitation is the fact that our linear approximation cannot describe strong-gravity effects, in particular the formation of a black hole and its subsequent relaxation. This means that the part of the gravitational radiation that is accurately captured by our approximation is that near the future lightcone of the collision point. This pulse will be followed by radiation emitted in the relaxation process to the final, equilibrium state. The crudest features of the final state can be accounted for by introducing appropriate cutoffs, as we have explained. However, a more precise determination of the relaxation dynamics to this final state will require a non-linear analysis, which we leave for future work.

A The function aω
Because it plays such an important role in our analysis, in this appendix we present the main properties of the function aω introduced in (3.29). The overall behavior is shown in Figure 18. h,2 = log The solution which is regular everywhere is h,1 . With this procedure, what we have accomplished is to maintain the singularity behavior of Ψ (0) II at r 0 , at the expense of changing the normalization at infinity. This way, the constant aω needs to be redefined. In particular, we get that close to the AdS boundary we have As an example, let's work to second order inω. We find that close to the tip of the soliton, r = r 0 , the singular solution behaves like Ψ II = (1 + ω 2 ) log r 5 Thus, to summarize, we have The numerical results are in perfect agreement with this prediction: we find aω = −1 + 0.13229ω 2 from fitting the data to a parabola dependence. The procedure can be easily extended to higher orders.
Large momentum behavior. At large negative values ofω 2 , we find an exponential increase The WKB analysis can be used to understand this scaling, giving a large k expansion of the wavefunction [35] Ψ(ρ) ∼ e α(ρ)k , (A.10) Here 2 F 1 is an hypergeometric function, while F was defined in (2.1) as a function of r = r 0 ρ. Close to the tip, ρ = 1, we find α(ρ = 1) 1.25, in good agreement with a numerical fit of our data.
Large frequency behavior. Finally, at large positive values ofω 2 , we find an oscillating power-law decay of the form Notice that the period of this ringing pattern is roughly twice as large as the spacing of the resonant modes in Table 1, which is a good consistency check.

B Vector-I gravitational perturbations
In the main text we studied gravitational modes which are excited by the axisymmetric collision described in this work. These modes were referred to as gravitational scalar and vector-II modes in Table 1. However, other types of perturbations are excited more generically. Here we briefly show the existence of at least one other perturbation mode, which is vector-type (with respect to the (t, r)−subspace) and takes the form (4.1) with the perturbation quantities defined as where an integral over ω and k i (i = 1, 2, 3) is implicit. Inserting the ansatz above into the Einstein equations we get E xy : (4r 5 + 1)h ry + ir 2 ωh ty + r(r 5 − 1)h ry = 0 , where primes stand for derivatives in r and note we are setting r 0 = L = 1. One can solve E ry for h ry and its derivatives,  Table 1 as "vector I" family. Note that for ω = k there is an exact solution regular at infinity, h ty = r −3 . However, this solution is not compatible with the exact solution for (B.5) regular at r = 1, h ty = (r 5 − 1)/r 3 . Thus, in particular there are no zero modes.
At large overtone numbers the spacing is again consistent with the geodesic approximation. In fact, the substitution r ≡ 1/η leads the equation above to the following WKB-amenable form, Following the scalar-field analysis, with the same largeω behavior, one finds the same asymptotic expression.

C Determination of eigenmodes by Frobenius expansions
In this appendix we describe a method to obtain the eigenmodes of a boundary value problem defined by a system of coupled ODEs. This procedure is an alternative to the method discussed in Section 4.1.
Here we adopt a series solution to the problem, by first defining the wavefunctions H rr and Z + according to which are then expanded in a power series around η = 1/r = 1: The radius of convergence of this series is at least as large as the distance to the closest singular point, thus it should converge on the entire interval 0 < η < 1, corresponding to 1 < r < ∞. The coefficients a q , b q can be obtained by direct substitution into the two coupled ODEs. Because it is linear, we can choose a 0 = 1, and all coefficients will be functions of ω and b 0 . This expansion satisfies the boundary conditions at the tip η = 1; to satisfy the boundary conditions at infinity (η = 0) we require that This results in two conditions for two quantities,ω and b 0 . The series has to be truncated at some value; we typically need around 30 terms in the expansion to get an accuracy of 1%. This second method yields values in very good agreement with the direct integration procedure, but computationally it seems more costly. To get higher overtones, one needs to keep more terms in the series to get a good convergence rate. This method is an extension of a Frobenius expansion used originally in Refs. [48,49,50] (see also Ref. [44] for a review).

D Absence of a zero mode
Our numerical investigations in Section 4.1 did not reveal the existence of any zero mode (ω = 0). Here we explicitly rule out this mode in the static limit, which then dictates generic Yukawa-type decay at large distances R. For ω = k = 0, equations (4.20)-(4.23) simplify to The solution does have B ∞ = 0, and thus it corresponds to a spacetime with a deformed boundary, because both h tt and h yy asymptote to a constant. We should exclude these solutions. Therefore, no zero-mode solution exists in this spacetime. E Green's function analysis for a static point particle in the AdS-soliton background In this appendix we discuss a Green's function approach to solve the system (4.20)-(4.21 (E.1) For clarity of notation, we restore all factors of r 0 and L in this appendix.
In order to construct the matrix X (see Section 3.2.1), we need four independent solutions of the homogeneous system which satisfy the correct boundary conditions. Close to r 0 , the general solution reads We require the fields to be regular at r = r 0 , so a k = c k = 0 and the correct asymptotic behavior reads where c xx , by solving the equation in the r → r 0 limit perturbatively. Then, two independent solutions can be found by imposing these two parameters to be (1, 0) and (0, 1). We shall denote these solutions by Y (r 0 ,1) and Y (r 0 ,2) , respectively.
Likewise, at infinity we have the following general behavior: Table 2: Schematic asymptotic behavior of the homogeneous solutions of the system (4.20)- (4.21). Recall that the superscripts (1,2) in h rr and h xx denote the choice (1, 0) or (0, 1), respectively, for the couple of independent parameters of the expansion at infinity. Making similar choices for the independent parameters of the expansion near r = r 0 and integrating out to infinity gives the asymptotic behavior of h where with W ≡ det(X). In writing this we are localizing the particle at r = b but in the end we want to take the limit b → r 0 , as in Section 3.2. The functions C (i) ± depend on the solutions of the homogeneous system. For completeness, their expressions read If r > b → r 0 , then the solution reads , , (E. 20) and, as in the scalar case with n y = 0, this limit is finite. Finally, we obtain Eqs. (4.24) and (4.25) in the main text.
F Computation of the holographic stress-energy tensor for high-energy particle collisions in the bulk

F.1 Asymptotic behavior of the metric perturbation
The asymptotic behavior at large holographic distance r of the metric functions defined in Eq. (4.2) reads (recall we are choosing a gauge such that h rx = 0) which straightforwardly give Therefore, the coefficients A ∞ and C ∞ , which shall be crucial for our analysis, are related to the dominant and subdominant terms of z + and to the dominant term of z − , respectively.
Asymptotic behavior of the solutions of the homogeneous system. In order to apply the Green's function technique discussed in Sect. F.2 below, we need the asymptotic behaviors (at r ∼ r 0 and at infinity) of the solution of the homogeneous field equations, i.e. Eqs. (4.8)-(4.10) without the source terms. The Green function method requires two independent solutions of the homogeneous system, the first one being regular at r ∼ r 0 and (for generic values of the frequency ω) irregular at infinity; the second solution is regular at infinity and generically irregular at r ∼ r 0 . We shall denote these solutions as X (r 0 ) and X (∞) , respectively, where X collectively denotes any perturbation variable. By analyzing Eqs. (4.8)-(4.10) at infinity and at r ∼ r 0 , we obtain the behavior of the relevant metric functions given in Table 3.
with (i = 1, 2) and W ≡ det(X). The functions C (i) ± , which have lengthy expressions and so we avoided presenting explicitly, are similar to the functions C (i) ± defined in Appendix E for the static case. If r > b → r 0 , then the solution reads Finally, using the expressions above and the asymptotic behaviors shown in Table 3, we obtain which, together with Eq. (F.13), constitute the main results of this section. Note that, while A ∞ and B ∞ do not depend on k 1 in the ultrarelativistic limit, C ∞ still has explicit dependence on the longitudinal component of the wavevector.

F.3 Holographic stress-energy tensor
In this part we shall obtain the stress-energy tensor T ab of the holographic boundary theory, as determined by the asymptotic expansion of the bulk gravitational field via holographic renomalization.
The metric we consider is ds 2 = r 2 L 2 (η ab + h ab ) dx a dx b + (1 + h rr ) dr 2 F (r) + 2h ar drdx a + F (r)(1 + h yy )dy 2 , (F.27) Here x a = (t, x 1 , x 2 , x 3 ). The h µν are linear perturbations, which we require to vanish at infinity, and we only consider cases where Asymptotics and stress-energy tensor. We will obtain the stress-energy tensor by writing the metric near infinity in the Fefferman-Graham gauge ds 2 = L 2 r 2 dr 2 +r 2 L 2 η mn + 1 r 5 Here m, n = t, x 1 , x 2 , x 3 , y.

(F.57)
This leaves five coefficients, e.g. a rr , b rr , b tt , a 1r , b 11 as independent coefficients. Two of them can be eliminated. If we choose the gauge h 1r = 0, as in Section 4.1, then a 1r = 0. Other non-dynamical Einstein equations correspond to constraints due to gauge choices in directions other than the radial one. This is the case with the algebraic equations (4.22) which imply that b rr = − A Since in our case we have h rr (r 0 ) = h yy (r 0 ), the y circles in the solution and in the background will match if we take the same parameter r 0 for both of them. Then, in order to subtract the background contribution from the stress-energy tensor we simply have to remove the terms ∝ r 5 0 from T mn (this subtraction could equally well be done by adding an appropriate counterterm to the action).
Final result. After performing the background subtraction we obtain The other five non-vanishing components, T x 3 x 3 , T tx 3 , T x 1 x 2 , T x 1 x 3 and T x 2 x 3 , can be obtained by using the tracelessness and divergence-free conditions. As a check on our calculations, we have computed these components explicitly and checked that the tracelessness and divergencefree conditions hold as a consequence of the asymptotic behavior shown in Eqs.