Double copy for massive quantum particles with spin

The duality between color and kinematics was originally observed for purely adjoint massless gauge theories, and later found to hold even after introducing massive fermionic and scalar matter in arbitrary gauge-group representations. Such a generalization was critical for obtaining both loop amplitudes in pure Einstein gravity and realistic gravitational matter from the double copy. In this paper we elaborate on the double copy that yields amplitudes in gravitational theories coupled to flavored massive matter with spin, which is relevant to the problems of black-hole scattering and gravitational waves. Our construction benefits from making the little group explicit for the massive particles, as shown on lower-point examples. For concreteness, we focus on the double copy of QCD with massive quarks, for which we work out the gravitational Lagrangian up to quartic scalar and vector-scalar couplings. We find new gauge-invariant double-copy formulae for tree-level amplitudes with two distinct-flavor pairs of matter and any number of gravitons. These are similar to, but inherently different from, the well-known Kawai-Lewellen-Tye formulae, since the latter only hold for the double copy of purely adjoint gauge theories.


Contents
3 Double copy in three-point vertices

Introduction
With the discovery of gravitational waves by the LIGO and VIRGO experiments [1,2], we are now entering an era where observations and theory will join forces to explore gravitational physics at unprecedented scales. Future experiments are in the planning, and on the theory side calculation methods are being sharpened to handle the need for higher precision. Focusing on analytical methods that use scattering amplitudes, there are by now a number of impressive results and techniques relevant to black-hole scattering [3][4][5][6][7][8][9][10][11][12][13].
A common theme of this work is the observation that gravitational amplitudes are vastly simpler than any Lagrangian may suggest. Together with the general expectation that all physical information of a theory should be contained in its S-matrix, one may hope to circumvent traditional methods and extract the needed quantities directly from amplitudes. One of the most remarkable properties of gravitational scattering amplitudes is a double-copy structure that allows them to be factorized in terms of simpler gauge-theory amplitudes. It was first observed by Kawai, Lewellen and Tye (KLT) [14] in the context of closed-string amplitudes, and it became a systematic method applicable at loop level and for a multitude of theories with the introduction of the Bern-Carrasco-Johansson (BCJ) double copy [15,16]. The BCJ double copy identified the underlying gauge-theory property -a duality between color and kinematics [15,16]. This duality is present in many physically relevant gauge theories, such as pure Yang-Mills theory [15,[17][18][19], quantum chromodynamics (QCD) [20][21][22], Yang-Mills-Higgs models [23], and supersymmetric generalizations thereof [16,20,[23][24][25][26][27][28]. The color-kinematics duality provides a mechanism that extends gauge invariance of two gauge theories to a diffeomorphism symmetry of the resulting gravity [15,[29][30][31][32][33]. The mechanism lends support to the expectation that the double copy may be a general feature of gravitational theories [34,35]. Indeed, by now it has been observed to give correct amplitudes in a large set of gravitational theories, such as pure Einstein gravity [20], pure and matter-coupled supergravities [16,20,28,36,37], Einstein-Yang-Mills theories [23,34,35,38], and many others.
Apart from flat-space scattering amplitudes, the double-copy structure has been found in other gravitational settings, such as for non-trivial curved backgrounds [39,40]. Remarkably, a number of important classical solutions to Einstein's field equations -including the Schwarzschild and Kerr black holes -may be expressed in terms of solutions in the underlying gauge theory [41][42][43][44][45][46][47][48][49]. Such a classical double copy, introduced by Monteiro, O'Connell and White [41], relies on the metric parametrizations in which Einstein's equations effectively linearize. For genuinely non-linear problems, it has been perturbatively extended to problems involving classical sources and gravitational radiation [43,[50][51][52][53][54][55][56][57][58][59][60][61][62]. This gives hope that the double copy can become a method for state-of-the-art calculations related to gravitational waves and effective black-hole potentials. In fact, the amplitude-level double copy has already played a crucial role in the calculation of the effective potential between two non-spinning black holes to third order in the post-Minkowskian expansion [9].
The classical double copy is still less well understood than the original "quantum double copy" used for scattering amplitudes. The precise prescriptive details are still under development, and a number of puzzles remain. 1 For instance, there is the issue of eliminating the dilaton and axion from the double copy with sources, which has not been solved in an entirely satisfactory manner. In quantum theories without sources, these two massless particles can be systematically canceled from loop diagrams using matter that obey ghost statistics [20,28]. In the presence of massive sources, however, the mass parameter induces a linear coupling to the dilaton (for scalar sources) and axion (for sources with spin). Then it becomes non-trivial to consistently truncate these massless particles even at the classical level. Nevertheless, if one has access to the on-shell states, as in a unitarity cut, it should be possible to always consistently project out the unwanted states from the double copy [9]. Alternatively, one can resort to computing the much simpler Feynman diagrams that involve the undesired scalars, and subtract them from the double copy [54]. The latter approach resembles the ghost prescription at loop level, but is still unclear if it is part of a systematic and simple double-copy formalism.
These puzzles call for more systematic studies where classical perturbative results can be matched to the corresponding amplitude calculations, for which the double copy is better understood. An interesting study in this spirit is ref. [57] by Shen, where a five-point amplitude, directly related to the emission of gravitational radiation from binary system, was computed in the classical limit at the next-to-leading post-Minkowskian order. Notably, the calculation was simplified by mapping the one-loop two-source problem to a tree-level amplitude with three sources. It can thus be fruitful to improve our understanding of tree-level gravity with more than two sources as obtained through the double copy.
The classical double copy in the presence of multiple sources is indirectly related to earlier studies of a double copy for (super)gravity amplitudes coupled to distinctly flavored matter particles. In ref. [21] the double-copy formula M tree n,k = σ∈basis K(1, 2, σ)A(1, 2, σ) , (1.1) was introduced for general matter-coupled tree amplitudes in gravity. It takes as input gauge-theory partial amplitudes A(1, 2, σ) with distinctly flavored matter, such as QCD amplitudes, as well as kinematic factors K(1, 2, σ) that are dual to a basis of QCD color factors. This general construction bears strong resemblance to Shen's detailed calculation [57]. While it is applicable to tree-level gravitational processes for any number of sources, the formula has the drawback that the individual kinematic factors are gauge-dependent. The standard KLT formulae [14,64,65], which are manifestly gauge-invariant, cannot be used together with QCD amplitudes that contain multiple quark pairs, since the fundamental quarks alter the general amplitude properties, such as the BCJ relations [21]. Indeed, the KLT double copy implicitly assumes that the gauge theory has only adjoint particles, and if not, then the gauge-invariant double copy needs to take a different form. In ref. [66] a first attempt was made to find a gauge-invariant double copy, where a general formula applicable for QCD amplitudes was given in terms of an unknown KLT-like 1 Curiously, a recent classical double-copy calculation of the effective action of two massive spinless particles failed to match the known next-to-next-to-leading order result [63]. It remains to be seen if this is because of some inherent problem, or simply due to the ambiguities present in off-shell unphysical quantities.

Gauge theory
Gravitational theory = (Gauge theory) 2 YM Axiodilaton gravity = GR + dilaton + axion Table 1. Field content for pure Yang-Mills theory, self-conjugate QCD [74] and conventional QCD and their respective double copies matrix. Similarly to the adjoint case, it is related to the inverse [67] to the matrix of the doubly color-stripped amplitudes in a scalar theory with ϕ 3 interactions [60, [68][69][70]. The difference is that in the flavored case there are several types of scalars: one bi-adjoint and a family of bi-fundamental scalars [71]. Since the tree amplitudes of the cubic scalar theory are in principle straightforward to work out (see e.g. ref. [72] and section 5.3 here), the non-trivial problem is to invert the resulting matrix. Its rank is determined by the number of partial QCD amplitudes that are independent under the BCJ amplitude relations, which at multiplicity n, for k ≥ 2 sources, is given by (n − 3)!(2k − 2)/k! [21]. In section 5, we will determine this KLT-like matrix for certain sectors in the (n, k) space. In particular, we give all-multiplicity formulae for k ≤ 2, as well as lower-multiplicity results for k = 3.
Identifying the gravitational theory. In this paper, we set out to provide a general framework for massive double copies with spin. We use the double copy of QCD as our prime example. Picking a group representation of matter particles other than the adjoint (e.g. the fundamental representation of SU(N c ) for quarks) makes it possible to separate matter and gluons in the double copy [20,73]. While the construction is guaranteed to give diffeomorphism-invariant amplitudes [30], it is general a non-trivial problem to specify the gravitational theory that results from the double copy. The simplest step is, of course, to identify the spectrum of on-shell states. For massless quarks, the QCD double copy produces either massless scalars or massless vectors depending on how the spins of the quarks are aligned [20]. However, for more realistic matter we take the quarks to be massive, and then the double copy produces all the states in the tensor product of two spin-1/2 particles: a massive scalar and a massive abelian vector. If we consider a family of N f quarks, we obtain a family of flavored massive scalars and vectors in the gravitational theory. Note that each vector is by construction paired up with a scalar, and they have the same mass as the quark from which they originate. In addition to the massive matter, we also have a massless dilaton and an axion present in the gravitational theory, as reviewed in appendix A. See table 1 for a summary of the different constructions and state countings. A much harder problem is to understand how the gravitational matter interacts in the double-copy theory. For the massless and supersymmetric cases, one can identify the resulting gravitational theories with certain minimally-matter-coupled Maxwell-Einstein supergravities [31,34,38], such as the "generic Jordan family" (in 6D) [75], the "generic non-Jordan family" (in 5D) [76] and the Luciani model (in 4D) [77]. More precisely, the double copy (QCD) ⊗ (N = 2 SQCD) between two four-dimensional massless theories gives amplitudes in the N = 2 Luciani model [34,37]. Truncating away the N = 2 supersymmetry should give the double copy (QCD) ⊗ (QCD) albeit in the massless case.
In sections 2, 3 and 4, we work out much of the structure of the Lagrangian that results from double copy (QCD) ⊗ (QCD), by matching a Lagrangian ansatz with the amplitudes from the double copy. To construct the amplitudes, we exploit the little-group structure of the massive particles by using the massive spinor-helicity formalism [78]. The resulting Lagrangian inherits non-trivial matter interactions from the supersymmetric parent theory. However, the non-zero masses make the interactions considerably more complicated than in the Luciani model [77]. See refs. [35,79] for other recent massive double-copy construction for gauged supergravity theories which may offer some guidance to the massive Lagrangian.
The gravitational Lagrangian has to satisfy certain properties in order to be compatible with the double copy (QCD) ⊗ (QCD). First, for N f = 0 the Lagrangian should match Einstein gravity coupled to a dilaton and axion, called "axiodilaton gravity" henceforth. The 4D renormalizability of QCD imposes, by dimensional analysis, that the resulting gravitational Lagrangian have at most two derivatives in each term. Similarly, masses can appear at most quadratically, and each power of m lowers the number of derivatives by one. Moreover, the massive vector fields contain longitudinal modes, which in the massless limit behave as ε µ L ∼ p µ /m. Since they must produce well-behaved scalar interactions in the massless limit, each explicit vector field V µ must be multiplied by one power of m, or, alternatively, be part of a covariant field strength V µν = 2∂ [µ V ν] that is insensitive to V µ ∼ p µ . By dimensional analysis, this limits the vector fields to appear at most quadratically in the Lagrangian. Finally, one should expect to obtain interactions consistent with a truncation of the Luciani model [77] in the massless limit. However, this is quite a delicate matching, both because of the needed supersymmetric truncation and since the longitudinal vector mode should non-trivially combine with the scalar modes to form the complex scalars of the Luciani model. More work is required to work out the details of this last matching.
For convenience of the reader, below we quote all the determined terms in the gravitational Lagrangian obtained from the double copy of QCD, Here we follow the notation in which Lagrangians and actions are related as S = d 4 x √ −gL. We have explicitly checked through order κ 2 that it is consistent with the amplitudes given by the double copy.

Double copy in free theory
In this section we review the external wavefunctions for massive particles using the massive spinor-helicity formalism of ref. [78]. 2 As in the purely massless case, the idea of this formalism is to build up all of the scattering kinematics from basic SL(2, C) spinors such that all resulting formulae are covariant with respect to the little group of the each physical particle. For massive particles in four dimensions the little group is SU (2). As usual, massive momenta are constructed from two-spinors using the Pauli matrices: where the SU(2) little-group indices a, b, . . . = 1, 2 should be distinguished from the spinorial SL(2, C) indices α, β, . . . = 1, 2 andα,β, . . . = 1, 2, which represent the Lorentz group.
In that sense, the latter indices can be thought of as off-shell, whereas the little-group indices encode the on-shell spin degrees of freedom. In fact, any scattering amplitude can be represented [78] as having 2s symmetric indices for each s-spin massive state with momentum p, as carried by the helicity spinors λ a pα andλα a p .

Free spin 1/2
The most familiar massive particle with spin is the Dirac fermion. In the Weyl basis of the Clifford algebra, its external Dirac spinors can be built from the helicity two-spinors as as already explored in refs. [78,80]. We follow the textbook convention in that we normalize these spinors to 2m,ū This identity illustrates the general fact that the upper and lower little-group indices are related by complex conjugation, which will also be evident from numerous equations below. Let us use the Pauli-Lubanski (pseudo)vector to consider spin. More precisely, for m = 0 we use its mass-rescaled version as the spin vector where S µν is the intrinsic angular momentum operator that satisfies the Lorentz algebra In the Dirac case, it is explicitly where σ µν, β α andσ µν,αβ form separate left-and right-handed representations of the Lorentz algebra (2.5), that are selfdual and anti-selfdual, respectively. Rewriting the spin operator as S µ = − 1 4m [γ µ , p ]γ 5 , we can directly compute its one-particle matrix element and observe that it gives the exact same value as the familiar Dirac-spin operator γ µ γ 5 /2. We will rely on the associated textbook spin vector s µ (u s ) = 1 2mū s γ µ γ 5 u s -only we rewrite it using the little-group indices Here we have used that the two spin states of the Dirac fermion have opposite spin vectors, and picked one of them. Now recall that the little-group indices a = 1, 2 permit SO(3) rotations of the spin quantization axis. This axis can be set to the three-momentum of the particle, which corresponds to the definite-helicity spinor parametrization detailed in refs. [78,80]. For a momentum parametrized as p µ = (E, P cos ϕ sin θ, P sin ϕ sin θ, P cos θ), the definite-helicity spin vector is explicitly s µ p = 1 m (P, E cos ϕ sin θ, E sin ϕ sin θ, E cos θ). (2.9) However, even in a generic representation the one-particle expectation value of the spin operator is S µ a p =ū pa S µ u a p u pa u a p = s µ p /2, a = 1, −s µ p /2, a = 2. (2.10) Our conjugation conventions for the external spinors are consistent with the massless convention (λ pα ) * = sgn(p 0 )λ pα . Here we have used the charge-conjugation operator C, which plays a crucial role when one wishes to impose the Majorana reality condition on the fermionic fields As detailed in appendix B, the properties (2.11) imply that the same external spinors (2.2) may be used to construct amplitudes with external Majorana fermions. The translation from QCD to its self-conjugate version [74] is hence straightforward.

Free spin 1
Spins higher than one are represented by 2s symmetrized little-group indices [78]. The effect of this symmetrization first shows itself for massive vector particles. Their polarization vectors were first spelled out in refs. [8,81]. They are transverse, as required by the Proca equation: To see how the spin is encoded by the symmetric indices a, b = 1, we use the Lorentz generators in the vector representation They feed into the Pauli-Lubanski vector (2.4), and its one-particle expectation value is explicitly quantized in terms of the unit spin vector (2.8): Notice that the spacelike normalization in eq. (2.13) implies ε p11 · ε 11 p = ε p22 · ε 22 p = −1 but at the same time ε p12 · ε 12 p = −1/2. However, this turns out to be convenient: in the completeness relation, for instance, the factor of 1/2 is automatically compensated by summing over ε 12 p and its equal ε 21 p . One can also confirm that in the massless limit the longitudinal polarization vector ε µ12 p behaves as p µ /m, up to a numerical factor. There is a hint at the massive double-copy relation already in the above definition of the spin-1 wavefunction, which we can rewrite as a bi-spinor This is reminiscent of the well-known fact that a symmetric tensor product of two spinors gives a vector. The remaining antisymmetric combination of spinors can be regarded as a scalar via the spinor identities Of course, it is difficult to interpret the above equations as strict double copies at the level of external wavefunctions. For instance, the anti-aligned Dirac arrows can be shown to give rather perplexing identities: In section 3, however, we will show that the symmetric and antisymmetric combinations of QCD amplitudes with Dirac or Majorana spinors combine exactly into gravitational amplitudes with massive vector and scalar bosons, respectively.

Higher spins
Here let us remark here that the external wavefunctions for spins higher than one are built out of the lower-spin ones in a way that explicitly respects the double-copy construction. The integer-spin polarization tensors are simple symmetric tensor products of the vectors (2.13) [8,81] ε a 1 ...a 2s pµ 1 ...µs = ε (a 1 a 2 pµ 1 . . . ε a 2s−1 a 2s ) pµs . (2.20) The construction is similar for half-integer spins, and it may be instructive to consider massive Rarita-Schwinger particles here in more detail. The external wavefunctions for spin-3/2 particles considered e.g. in refs. [82,83] can be constructed as explicit double copies those for spin 1/2 and 1: and similarly for the v abc pµ . The angular momentum operator is constructed accordingly as a direct sum of the spin-1/2 and spin-1 ones: The expectation value of the Pauli-Lubanski operator gives the expected four states of the massive Rarita-Schwinger particle: (2.23) Similarly to Dirac fermions, we normalized the wavefunctions using the particle's mass: where the sign is due to them being spacelike in the mostly-minus metric convention. Again, the normalization of the lower-helicity states, such asū p112 =ū p121 =ū p211 , cancels their overcounting in state sums of the form u abc pµūpν abc . The formulae given here demonstrate that the massive-spinor helicity formalism is very useful if one wishes to consider a quantum field theory of higher-spin states. In this paper, however, we restrict ourselves to examples with lower-spin matter content, as present in QCD (s ≤ 1/2) and its double copy (s ≤ 1).

Double copy in three-point vertices
In this section we explore the double copy of the three-point matter amplitudes in QCD.
Note that at three points there is no need to distinguish between the KLT [14] and the BCJ double copy [15], as the full amplitude consists of one gauge-invariant diagram, and both constructions amount to where A 3 and A 3 are color-stripped coupling-stripped gauge-theory amplitudes. For instance, the gluonic self-interaction vertex amplitudes 3 square to the gravitational vertices The massive-matter results, summarized in section 3.4, are very simple and, up to minor modifications, satisfy the double-copy relation (3.1). For completeness, we go through the cases of s = 0, 1/2 and 1 both in gauge theory and in gravity and compute the threepoint amplitudes from explicit Feynman rules. For consistency purposes, we consider complex matter in a fundamental representation -without relying on specific details of the gauge group. This means that the representation can be easily switched, for example, to a real representation. Then the exact same Feynman rules as given below for complex matter may be used for self-conjugate matter, the least trivial case being the switch between Dirac and Majorana fermions, as reviewed in appendix B. As usual in that case, the necessary factors of 1/2 in the action are then compensated by the combinatoric factors due to self-conjugacy.

Three-vertex for spin 0
In this section we warm up with a scalar in gauge theory and gravity. The familiar Lagrangian implements both couplings minimally. Taking g µν = η µν and D µ = ∂ µ − igA µ gives the gauge-theory vertex where all the momenta are taken to be incoming. On the other hand, considering the massive scalar to be colorless and expanding the metric √ −gg µν = η µν − κh µν , we find the gravitational coupling vertex (3.6) where the second contribution in the bracket is a trace term due to √ −g = 1−κh/2+O(κ 2 ) from the spacetime measure. The corresponding three-point amplitudes are where we used momentum conservation p 1 + p 2 + p 3 = 0 for all states considered incoming, as well as the transverse and traceless polarization tensor ε µν 3 = ε µ 3 ε ν 3 of the graviton. The double-copy relation (3.1) is now evident.
This scalar example is convenient to explain the helicity variable x 3 , first introduced in ref. [78]. It follows from the on-shell conditions for the three-point kinematics that there exists a proportionality coefficient x 3 between the spinors |3 and |1|3]: where we have used an arbirary pair of spinors |q and |q] to write an explicit but nonlocal expression for x 3 . It is, however, independent of them as long as 3 The coefficient x 3 is designed to be dimensionless and carry a unit positive helicity with respect to momentum p 3 . The three-point amplitudes (3.7) are then proportional to the appropriate powers of this helicity factor:

Three-vertex for spin 1/2
In this section we consider a spin-1/2 particle minimally coupled to a gauge or gravitational field. For simplicity, we consider the complex case of the Dirac fermion with the projection to the Majorana case expained in appendix B.
In gauge theory, the three-point vertex is evidently We dress it with the external spinorsv a 1 and u b 2 from eq. (2.2) and use a Schouten identity to write the resulting amplitudes in the spinor-helicity form: Extracting the gravitational three-point vertex is somewhat more involved (its derivation from the Noether energy-momentum tensor is given in appendix C). It gives Plugging in the spinors and the polarization tensor is then no more challenging than in the case of QCD. We immediately get the amplitudes They obey a simple massive implementation of the double-copy relation (3.1) in the sense of multiplying the spin-1/2 amplitudes (3.13) with their scalar counterparts in eq. (3.10). Such gravitationally-interacting fermionic matter is, however, not part of the double copy of standard QCD, but it is relevant for supersymmetric versions of QCD, or for cases when some quarks are replaced by fundamental scalars.

Three-vertex for spin 1
In gauge theory, the interaction vertex for a massive vector V µ obtained from spontaneous symmetry breaking, is given by (see e.g. ref. [84]) Note that it does not coincide with what one would obtain from the free Proca Lagrangian (2.14) by a simple replacement ∂ µ → D µ [81]. The non-abelian-like interaction is required for obtaining a scattering amplitude that is well-behaved in the massless limit [78], which is compatible with spontaneous symmetry breaking. Using the massive polarization vectors introduced in eq. (2.13) and the three-point identities, we are able to arrive at the following expressions Here we have encoded the symmetrization of the little-group indices for particles 1 and 2 using a symmetrized tensor product operation, denoted by . It should be remembered that such a symmetrization is always applied to the spin indices of each massive particle independently.
While we do not include massive spin-1 bosons in the QCD Lagrangian, it is natural to consider double copies using spontaneously broken gauge theories, see refs. [23,35,79]. For our purposes, however, the double copy of QCD does contain gravitating massive abelian vectors. Let us then consider a covariantization of the Proca action (2.14), which amounts to replacing the Lorentz contractions with √ −gg µν = η µν − κh µν . Then the three-point vertex is where the last line is a trace term that vanishes in the resulting three-point amplitude Plugging in the spinor-helicity variables, we find 4

Three-point summary
The three cases for s = 0, 1/2 and 1 can be encapsulated into the following formulae. The gauge-theoretic amplitudes and their gravitational counterparts are all consistent with the arbitrary-spin amplitudes first proposed in ref. [78]. There, they were singled out of other possible spinor-helicity expressions by their tame behavior in the massless limit, which is said to define a notion of "minimal coupling" for massive matter to the gauge and gravitational fields. It now known to coincide with the textbook notion of minimal coupling only for lower-spin cases [81], on which we focus on in the present paper. The purpose of the subtle sign prefactors in eqs. (3.22) and (3.23) is to match the Feynman-rule calculations in the mostly-minus metric convention: they come from the standard requirement that the polarization tensors be spacelike. It is evident that the arbitrary-spin amplitudes above obey a simple extension of the double-copy relation (3.1) involving symmetrization of the little-group SU(2) indices, in which gravitational matter of spin s is constructed from gauge-theoretic matter with spins s 1 and s 2 : Below we consider in detail the three-point double copies that can be constructed from the QCD vertices (3.13). This includes special cases of the little-group-symmetrized double copy above, as well as less obvious little-group-antisymmetrized double copies.

Double copy of QCD three-vertex
As the simplest non-trivial example, let us consider the double copy of the quark-gluon amplitudes (3.13). First, we consider the little-group-symmetrized combinations 5 Here the first two amplitudes contain a graviton and match the minimal-coupling formulae (3.21). The last two amplitude is equal to −iκm 2 (ε 1 · ε 2 )/2 and can be seen to correspond to a non-trivial coupling of the axiodilaton scalar Z (orZ) to the vector mass term: (3.26) The scalar Z implements the long-known fact [85] that two massless scalar states are produced in the double copy of pure Yang-Mills theory in addition to the two graviton states. More precisely, the additional states correspond to a dilaton scalar φ and Kalb-Ramond two-form B µν . In appendix A we review the four-dimensional dualization of the latter field to a real pseudoscalar a, which can then be combined with the real scalar φ into a single complex axiodilaton scalar Z. The resulting Lagrangian for the massless sector of the double-copy theory is The fact that the vertices (3.26) stay proportional to m 2 even off shell (instead of developing additional momentum dependence) follows from inspecting the massless limit, in which the Z scalar should have a linearly realized U(1) symmetry (which is part of a larger U(N f +1) symmetry of the Luciani model [77]). Hence any U(1)-violating axiodilaton coupling must vanish in the massless limit. Now let us inspect the little-group-antisymmetrized double copies that produce the massive scalar ϕ. We see that it naturally couples both to the graviton and to the axiodilaton Z: 6 Here the last amplitude equals iκm 2 /2 and is due to the three-scalar coupling Again, that these vertices stay proportional to m 2 off shell follows from the fact that in the massless limit the scalar Z should have a conserved U(1) charge. 6 The following identities are helpful for considering antisymmetric double copies: Finally, we should not forget about the combination of litte-group symmetrization and antisymmetrization in the double copy: Here we see, as expected, that the gravitational interaction does not mix massive scalars and vectors, whereas the interaction with the axiodilaton does. The latter non-trivial interaction is consistent with the following vertices where the relative signs are dictated by the hermiticity of the resulting interactions with either complex or self-adjoint matter, hence the i prefactors in the double copies (3.31).

Double copy for four-point vertices
In this section we transition to four-point amplitudes in gauge theory and gravity. First we focus on Compton scattering, by which we understand a process with two matter particles 1 and 4 interacting with two massless bosons 2 and 3. The resulting four-point double copies are consistent with the precise off-shell forms of the three-point vertices given in the previous section. Then we consider the pure-matter amplitudes, for which we have to differentiate between the complex and self-conjugate matter, as well as discuss the introduction of multiple flavors of such matter.

Compton scattering in gauge theory
First of all, we use the BCFW approach [86,87] to derive the scattering amplitude for two massive particles and two gluons. Assuming for concreteness the distinct-helicity case, we shift the spinors of the gluons 2 + and 3 − to We also assume a vanishing behavior at z → ∞, which is guaranteed for s = 0 and 1/2 [88][89][90]. Instead of specializing to the latter case right away (also considered in ref. [80]), let us leave the spin s of the massive matter unspecified for the moment, which will make our discussion in section 6 more interesting.
There are two physical poles in z that come from the s 12 -and s 13 -channel propagators and localize at These poles imply the following kinematic identities: whereP is the cut momentum equal to either −(p 1 +p 2 ) or −(p 1 +p 3 ), respectively. Since on each of the poles the four-point amplitude factorizes into two three-point amplitudes of the type (3.22), we compute the residues as Here we have also taken care of putting a spin-dependent sign in the propagator. This is due to the fact that the completeness relations in the propagator numerators has a sign-flip pattern (−1) s under |−P e = −|P e , |−P e ] = |P e ]. Indeed, for spins 1/2 and 1 we have and so on for the higher-spin external wavefunctions given in eqs. (2.20) and (2.21). Therefore, the full color-dressed amplitude is given by the sum of eqs. (4.4a) and (4.4b): 7 The same method can be used to derive the identical-helicity amplitudes, that are even simpler. We find and A(1 {b}  ) is the same up to the bracket swap . . . ↔ [ . . . ]. Now let us apply the standard adjoint-representation color-ordering [91][92][93] to either of the above Compton amplitudes. For instance, the all-plus amplitude (4.7) decomposes into three ordered amplitudes We are now in a position to observe that, defined in this way for arbitrary spin, the colorordered amplitudes satisfy both the Kleiss-Kuijf (KK) relation [94] A(1, 2, 3, 4) + A(1, 3, 2, 4) + A(1, 3, 4, 2) = 0 (4.9) and the BCJ relation [21] (s 12 − m 2 )A(1, 2, 3, 4) = (s 13 − m 2 )A(1, 3, 2, 4). (4.10)

Compton scattering in gravity
Recall that the BCJ relations are a manifestation of the color-kinematics duality [15,16,20,21], and the duality implies the BCJ double copy [15]. The BCJ double copy is also known to be equivalent [95] to the tree-level KLT formulae [14] that are valid for external adjoint particles, or at most two non-adjoint ones [21]. Therefore, the fact that the colorordered amplitudes defined above obey the BCJ relation (4.10) means that the gravitational Compton amplitudes are given by the following gauge-invariant formula where again the symmetrized tensor products builds up gravitating matter of spin s = s 1 + s 2 from two copies of QCD matter of spin s 1 and s 2 . We have included the same sign prefactor that appears in the three-point double copy (3.24). This prescription implies The first of these expressions was earlier obtained in ref. [78] from matching to the residues of the three physical poles, while to our knowledge the second has not previously appeared in the literature. Refs. [78,81] go into some detail discussing the implications of the unphysical pole 3|1|2] that appears in eq. (4.12a) for s > 2. There is still a lot to be learned about massive higher-spin fields [96][97][98][99][100] from the on-shell perspective. We postpone the discussion of the unphysical-pole issue at higher spins to section 6. Instead, let us specialize to the case of s = 1 constructed from two spins 1/2. We find that for all helicity configurations the litte-group-symmetrized formulae (4.12) agree with the gravitational Compton amplitude obtained from the Feynman diagrams: Here we have used the four-point vertex implied by the minimal coupling of the massive vector to gravity, which is given explicitly in eq. (D.7) of appendix D. The corresponding scalar amplitude is obtained from SU (2) antisymmetrization which is a four-point analogue of eq. (3.29). The resulting amplitude matches the Feynmandiagram calculation involving the four-point minimal-coupling vertex (D.6).
Choosing the helicities of one of the gluons in the double copies (4.13) and (4.14) to be anti-aligned results in the amplitudes M(1 V * , 2 h , 3 Z , 4 V ) and M(1 ϕ * , 2 h , 3 Z , 4 ϕ ), which are consistent with the off-shell form of the three-point vertices (3.26) and (3.30), as well as their minimal four-point extensions which come exclusively from the nonlinearities in the h µν -graviton expansion.
One should not forget that the symmetric and antisymmetric matter double copies can be combined, which for the three-point amplitudes we did in eq. (3.31). Its generalization to four points is straightforward: If the gluon helicities are aligned to produce gravitons, these double copies continue to vanish, which reinforces our claim that the vector V and the scalar ϕ are only mixed by their interaction with the axiodilaton Z. Indeed, switching one of the gravitons to either Z orZ allows us to probe the following Feynman vertices, Going back to the double copies (4.13) and (4.14) and also switching both massless legs to scalars, we find that M(1 V * , 2 Z , 3 Z , 4 V ) and M(1 V * , 2Z, 3Z, 4 V ) each correspond to two trivalent Feynman diagrams, while M(1 V * , 2 Z , 3Z, 4 V ) calls for an additional quartic vertex All these vertices are consistent with the Lagrangian (1.2) stated in section 1.

Four-matter vertices
In this section we turn to amplitudes with four matter particles. In QCD, the distinctflavor four-quark amplitude A(1 a i , 2 b  , 3 c k , 4 d l ) consists of the sole gauge-invariant Feynman diagram where the masses and gauge-group representations of the two quark lines may be different. The identical-flavor amplitude is then a combination of two such diagrams with the masses and representations taken equal (4.23) where the relative sign is due to Fermi-Dirac statistics. Now in the case of Majorana fermions, the distinct-flavor amplitude stays the same as eq. (4.22). That is, up to the switch to a real representation of the gauge group, which we take to be the adjoint representation for simplicity. Then the identical-flavor amplitude develops an additional pole that was forbidden by charge conservation in the Dirac case: 24) where the kinematic dependence is still given by permutations of eq. (4.22). Remarkably, a three-term identity for the kinematic numerators still holds with the relative signs of the Jacobi identity modified by the fermionic statistics. The BCJ relation however, is only implied by the color-kinematics duality (4.25) in the massless case, which is a well-known feature of supersymmetric Yang-Mills theories. Now in this paper we exclusively rely on the color-kinematics duality of amplitudes with distinctly flavored quarks [21,22]. Other QCD amplitudes may be then obtained as linear combinations of suitable relabellings thereof, as illustrated by eq. (4.23). In the same way, gravitational amplitudes with identically flavored matter can be defined through relabellings of those with distinctly flavored matter. So in the following we concentrate on how the latter are obtained as double copies of distinct-flavor QCD amplitudes.
The fully litte-group-symmetrized double copy of the QCD amplitude (4.22) gives the distinct-flavor four-vector amplitude that can be otherwise obtained from the three-point Feynman vertices (3.19) and (3.26) confirming that there is no four-vector vertex in the resulting gravitational theory. On the other hand, the fully little-group-antisymmetrized double copy results in such an amplitude for two massive scalar pairs 8 (4.28) which requires an additional four-scalar vertex Antisymmetrization of one or three pairs of little-group indices yields vanishing amplitudes with one or three scalars, in which only the axiodilaton may be exchanged: Note that in these amplitudes the two diagrams cancel due to the sign difference in the scalar-vector mixing vertices (3.32). Therefore, we see no need for vertices with more than two vector fields, in accord with our expectations explained in section 1. Finally, symmetrizing two pairs of SU(2) indices produces the following two-scalar two-vector gravitational amplitudes: where only the first one allows graviton exchange in addition to the axiodilaton exchange.
To coincide with the double copies above, these amplitudes require supplementary massinduced vertices: (4.32) 8 We use Mn as a shorthand for gravitational amplitudes without the trivial i(κ/2) n−2 prefactor.
This completes the list of vertices that could be determined from the three-and four-point double copy and enter the Lagrangian (1.2) advertized in section 1. We leave exploring the higher-point vertices implied by the double copy to future investigation.

Gauge-invariant double copy
In the previous sections we considered three-point and distinct-flavor four-point double copies, for which the transition between the BCJ and a KLT-like constructions was straightforward. Inspecting the amplitude properties was sufficient for identifying the necessary double-copy formulae, and so the main focus there was on the explicit gravitational expressions. In this section we switch focus to a more abstract approach for obtaining gaugeinvariant double-copy formulae at any multiplicity. This will reproduce the field-theory KLT relations [14], but will also include a more general setup that follows from the BCJ construction [15] using the color-kinematics duality and gauge-dependent numerators. The extension [20,21] of the latter to amplitudes involving massive and distinctly flavored matter in arbitrary gauge-group representations allows us to be very general in our discussion.

Generalities of BCJ double copy
Consider a tree-level QCD amplitude A n,k with (n − 2k) gluons and k distinctly flavored quark pairs. The group-theoretic content of QCD involves only structure constantsf abc and generators T a i , which build up all the color factors c j . All Feynman diagrams can be grouped into a distinct set of (2n − 5)!!/(2k − 1)!! purely trivalent graphs, dictated by the color factors and flavor assignments [21]. Now let us write the sum over these cubic graphs for three types of amplitudes: for the ϕ 3 theory of biadjoint and flavored bifundamental scalars (see the Lagrangian (5.22) below), for QCD with flavored quarks and for gravity with flavored matter: 9 where D j 's are the propagator denominators corresponding to physical poles, and n j are the kinematic numerators that captures the non-trivial kinematic dependence in the QCD amplitudes. We also introduce tildes to allow repeated color factors and numerators to originate from different theories. For the above connection between QCD and gravity to hold, the kinematic numerators n j must satisfy the same algebraic identities (Jacobi and commutation identities) as the color factors c j , i.e. they obey color-kinematics duality [15,21]. Hence the solution of the algebra of such identities is the same for c j and n j : where J j is an integer-valued rectangular matrix obtained by repeatedly applying Jacobi/commutation relations until a basis of irreducible color factors have been obtained. The rank of the solution matrix J j must equal (n − 2)!/k!, since this is the size of the Melia basis of color-ordered amplitudes in QCD [21,101,102]. Without going into details, we denote the independent color-ordered amplitudes using the abstract notation where the indices α can be understood as different cyclic orderings [101,102]. The Melia basis is independent under a subset of the Kleiss-Kuijf relations [94] that are compatible with the flavor assignments of the quarks. In ref. [21], a new color decomposition was given for a color-dressed QCD amplitude in terms of this basis of ordered amplitudes. 10 The corresponding color coefficients C α were given in terms of a explicit formula, and all the amplitudes in eq. (5.1) can now be written as The color coefficients C α are certain linear combinations of the cubic-graph color factors c j , which we can write schematically as C α = C α c. C α is an invertable square matrix that converts between C α 's and the basis of cubic-graph color factors c [21]. The gravitationalamplitude kinematic coefficientsK α are obtained from the color coefficients C α by swapping their constituent c j 's byñ j 's, namely Note that the above identity M n,k =K α A α is nothing but our initial double-copy formula (1.1) [21], rewritten in a more concise notation. Let us now indirectly define a propagator matrix P j α (see e.g. ref. [104]) for the colorordered QCD amplitudes, by factoring out the numerators, A α = P j α n j . Its rank must be (n − 2)!/k!, otherwise the Melia basis would not be a basis, so we can further reduce the numerators to the basis numerators n: Furthermore, inverting the color-factor matrix C α in eq. (5.6), we obtain In the purely gluonic case, the resulting square matrix P j α J j C −1 β turns out to be equal [67] to the sum over propagators of the graphs occuring in both color-ordered amplitudes with subtle signs: 11 10 Alternative amplitude bases and corresponding color decompositions can be found in ref. [103]. 11 The BCJ construction in this section defines A α|β with α and β being permutations in a given Melia basis. However, its amplitude interpretation (5.9) immediately extends the definition to a pair of arbitrary particle-label permutations α and β.
It also coincides with the doubly color-ordered amplitude for the theory of biadjoint and flavored bifundamental scalars. To prove eq. (5.9), it suffices to consider the zeroth-copy version of the double copy (5.5) where the last expression is a double color decomposition. An efficient way for computing this amplitude matrix is discussed below in section 5.3. The BCJ double copy in eq. (5.5) can therefore be rewritten as [21] M n,k = K βÃ β , where we also allowed the second copy to be in a different theory (with a matching flavor pattern), and β can in principle run over a different Melia-like basis.

BCJ relations and KLT-like double copy
The rank of the matrix A α|β is in general smaller than the Melia basis size. Assuming that color-kinematics duality (5.2) holds, we need to ensure that the linear system A α|β K β = A α is not inconsistent. The Kronecker-Rouché-Capelli theorem in linear algebra then implies consistency conditions on the non-homogeneity of that linear equation which coincide with non-trivial BCJ relations. Indeed, we can write a BCJ relation [15,21] schematically as A α = Fᾱ α Aᾱ. It expresses an arbitrary color-ordered amplitude A α in the Melia basis in terms of the amplitudes Aᾱ in the BCJ basis with the kinematic coefficients Fᾱ α . For the BCJ basis itself we have simply Fγ α = δγ α . 12 The nontrivial part of the BCJ relations (5.13) is of course given by eq. (5.12): it resolves the non-BCJ-basis amplitudes in terms of the BCJ basis. Plugging the consistency condition (5.13) into our linear equation A α|β K β = A α we find A α|β K β = Fγ α Aγ = Fγ α Aγ |β K β . (5.14) As this linear system has not yet been solved for K β , they can be regarded as arbitrary, and we conclude that In other words, the BCJ-relation matrix (Fγ α − δ γ α ) spans the kernel subspace of A α|β . Now if we restrict the doubly color-ordered amplitude matrix A α|β to a square submatrix related to two sets of BCJ basis orderings, we can define the momentum kernel as its inverse, as was done already in ref. [66]. We have Sᾱβ ≡ A −1β|ᾱ = (A −1T )ᾱ |β ,ᾱ,β = 1, 2, . . . , (n−3)!(2k−2)/k! (BCJ basis). (5.16) 12 Here we assume that the BCJ basis is chosen as a subset of the Melia basis, but in principle we can use Kleiss-Kuijf relations [94] to extend the BCJ relation matrix Fβ α to more general basis choices.
Then multiplication of eq. (5.15) by the momentum kernel on the left and right gives the following formulae for the BCJ-relation matrix F : We can now go back to the double copy (5.11) and solve the linear equation Aᾱ |β K β = Aᾱ for a subset of Kβ corresponding to the second BCJ basis Inserting it into the double-copy relation M n,k = KβÃβ + K βÃ β , we see that the double copy takes a KLT-like form 13 [66] M n,k = AᾱSᾱβÃβ. (5.21) The double copy (5.21) has the nice feature that it is manifestly gauge-invariant, but as a consequence it does not make manifest the simple poles and crossing symmetry of the gravitational amplitude. This is in contradistinction to the original BCJ double copy (5.1), where the poles and crossing symmetry are both manifest but gauge (diffeomorphism) invariance is not manifest. Interestingly, the double copy in eq. (1.1) has both part of the gauge invariance manifest and the simple poles are also obvious from the gauge-theory amplitude, but crossing symmetry and gauge invariance of the kinematic factor is absent.

Recursion for bi-colored scalar amplitudes
The gauge-invariant form of the double copy (5.21) involves the kinematic matrix Sᾱβ, defined in eq. (5.9) as the inverse of the doubly color-ordered amplitude in the theory of biadjoint and flavored bi-fundamental scalars. The Lagrangian for this theory is explicitly There is a good reason why the gauge-invariant double copy (5.21) looks like a symmetric bilinear form.
A basis change from one pair Aᾱ,Ãβ of BCJ bases to another such pair Aγ,Ãδ can be described by a generalized BCJ-relation matrix, that allows for Kleiss-Kuijf relations as well. The gravity amplitude is invariant with respect to such basis changes: It only produces trivalent Feynman diagrams that encodes the same color and propagator structure as QCD. An efficient way to compute amplitudes in this theory is via the Berends-Giele recursion for the doubly color-ordered current: where we use ⊕ to denote concatenation of two shorter particle-label permutations. As indicated, the current vanishes unless α is a permutation of β. The base of recursion is simply J i|j = δ ij . This recursion is different from the formula written in ref. [72] (for the purely adjoint case) by the θ factor, the purpose of which is to eliminate illegal flavor assignments. 14 The doubly color-ordered amplitudes are obtained from currents by removing the final denominator and putting the corresponding external momentum on shell: Therefore, the recursion relation given above provides an efficient way for numerical evaluation of gauge-invariant double copies. In eqs. (5.23) and (5.26) above and in the rest of this paper, we use the massless Mandelstam invariants s α = i∈α p i 2 for simplicity. The correct way to extend them to the massive case is by switching to a mass-deformed notation defined for a permutation α as a sum of flavor symbols of each particle, f (α) = i∈α f (i). Assuming that, as in ref. [21], the first k pairs of particle labels are assigned to the distinctly flavored matter particles and the remaining (n − 2k) labels correspond to the massless interaction particles, we specify f (i) = 0, i = (2k + 1), (2k + 2), . . . , n. For instance, f (1, 2, 5, 4) = −f1 + f1 + 0 + f2 = f2 is a single term, whereas f (1, 5, 4) = −f1 + f2 is not.
the Melia-type bases have a simple form of A(1, 2, α), with the relative ordering of 3 and 4 inside α fixed to either 3 ← 4 or 4 → 3 [101,102]. A BCJ basis is then imposed by requiring that one of the previously unfixed quarks be adjacent to the first pair, as in e.g. A(1, 2, 3, α), where α is an arbitrary permutation of the remaining particle labels [21]. In section 5 we assumed the same BCJ basis for both sides of the gauge-invariant double copy (5.21). Such a choice, however, leads to non-localities in the momentum kernel. This is actually true even in the flavorless case, i.e. for A(1, 2, 3, α) ×Ã (1, 2, 3, β). In fact, even the purely adjoint KLT double copy is realized non-locally for the BCJ basis combinations of the form A(1, 2, α, 3, β) ×Ã(1, γ, 2, 3, δ) with fixed sizes of α, β, γ and δ, unless at least one of them is set to zero. This can be explained by the fact all the involved primitives miss some physical poles, such as s 13 , which have to occur in the gravity amplitude and therefore must be introduced by the momentum kernel. A similar analysis implies that local combinations of flavored BCJ bases require that the gluon-free slots be adjacent to the quark pair fixed by the Melia basis, but it cannot be simply A(1, 2, 3, α) × A(1, 2, 3, β), as then neither side contains the s 2n -channel pole (where n ≥ 5).
The simplest gauge-invariant double copy that we managed to find is where b i is defined as the position of α i within β, i.e. β b i = α i . In fact, this basis combination corresponds to a momentum kernel that mimics precisely the unflavored one. For instance, it satisfies the same recursion relation that was previously noticed in ref. [105]: where the gluon set γ is defined as a complement of the ordered subset γ * , which picks some of the gluons trapped inside bracket {3 . . . 4} both in α and β: where auxiliary gluon subsets γ * (ordered) and γ (unordered) are defined as and the prefactor θ 31 equals These three k = 2 formulae have been checked through eight points. Needless to say, the latter two reveal a much richer structure than seen in the flavorless case. The seemingly most complicated one, eq. (5.31), actually contains much fewer nonzero terms than the other two. Indeed, the prefactor (5.31c) vanishes for a growing fraction of permutations, starting from 1/4 for n = 5 and reaching 63% for n = 10. This means that using this double copy to evaluate gravitational amplitudes must be numerically more efficient at higher points.

Six-matter double copies at low multiplicity
All-multiplicity gauge-invariant double-copy formulae for more than two matter particle pairs have proven highly non-trivial to derive. In

Summary and Discussion
In this paper, we have elaborated on a general framework for double copies applied to gravitational theories with flavored matter and spin. The framework is based on the results and lessons learned from the extension of the color-kinematics duality to fundamental matter [20], and in particular to QCD with massive quarks [21]. However, we expect it to be more generally applicable to gravitational amplitudes obtained from a pair of gauge theories that obey color-kinematics duality and that have massive matter transforming in non-adjoint representations of the gauge group. Due to this generality, it should also be useful for more systematic studies related to classical scattering of multiple gravitational sources and associated gravitational wave emission, as demonstrated by Shen [57]. See also refs. [6, 7, 9, 12, 13, 51, 54-56, 58, 59, 62, 63, 107, 108]. Specifying the properties that identify the gravitational theory given by the double copy is in general a non-trivial problem. The exception is when one can rely on uniqueness properties such as (super)symmetries and other easily identifiable traits of the underlying gauge theories. A major focus of the current work is to consider the double copy obtained from QCD with with N f massive quarks, which is a less straightforward task than previous identifications and classifications (see e.g. refs. [34,35,38,38,109]). The resulting gravitational theory contains N f massive vector bosons minimally coupled to gravity. Such non-self-interacting massive vectors are useful models for spinning black holes, which in the classical limit of scattering amplitudes can be used to obtain the lower-multipole corrections to gravitational observables [7,12]. It is an automatic feature of our QCD construction that the resulting double-copy gravitational theory can at most be quadratic in the massive vectors.
For every massive vector there is a companion scalar of the same mass, which interacts highly non-linearly with itself and with other matter field such as the (massless) axiodilaton scalar. This feature is an echo of the supersymmetric cousins of the same construction [34], which has non-linear scalar interactions parametrized by supergravity moduli spaces. This makes it harder to use a brute force method to construct the full scalar dependence in the Lagrangian. Nevertheless, we have through the double copy constrained the gravitational Lagrangian up to quartic couplings; the result is given in eq. (1.2). It would be interesting to pinpoint the remaining scalar interactions in theory, such that one obtain closed form for functions that parametrize the scalar space.
From the perspective of using the current construction for classical calculations of relevant observables, the massive scalars would most likely obfuscate any physical process. One way to avoid dealing with such scalars is to consider a double copy where the spinning matter belongs to only one gauge theory, and the other gauge theory has only scalar matter, as done in ref. [12]. The general framework of this paper is well suited for such a construction, and in particular problems with multiple sources may strongly benefit from the considerations in this work. Having the massive spin degrees of freedom come exclusively from one of the two gauge-theories should also make the problem of identifying the gravitational theory more straightforward.
A second major part of this paper was to elaborate on general features of the double copies with flavored matter -agnostic of the type of theories being considered -and to give manifestly gauge-invariant expressions for the double-copy formulae. Doing so, we found new all-multiplicity formulae for a gauge-invariant double copy for amplitudes with two flavored pairs of matter particles, as well as partial results for three matter pairs. In particular, we provide a double-copy formula [106] for a seven-point amplitude with three matter pairs. In the context of classical gravitational radiation, such a seven-point amplitude contains the information of the same post-Minkowskian order as considered in ref. [57]. A remaining task is to find the general formulae for the gauge-invariant double copy, with three or more flavors. While these general formulae would be useful for classical calculations, they might contain more information than what is necessary for practical calculations, as the classical limit usually involves placing various particles on shell, where the double copy and the corresponding momentum kernel factorizes into smaller pieces. Understanding this more systematically should be a promising research direction.
Double-copy tension at higher spins. Coming back to the double copy of non-zero massive spins, we recall the higher-spin issue [78] that we encountered in section 4. Namely, in the general-spin formulae (4.6) and (4.12a) the power of the factor 3|1|2] becomes negative for s > 1 and s > 2, respectively. In that case the amplitudes contain unphysical poles. If so, the BCFW-shift variable must not fall off fast enough at infinity for s ≥ 3/2 in gauge theory. The 3|1|2] pole can be removed if the neglected boundary contribution is recovered by other means. Refs. [78,81] treated this as an invitation to find an additional amplitude contribution that vanishes on all of the physical poles and subtracts out the unphysical pole from the BCFW result (4.6). Such a contribution, however, spoils the highenergy/massless limit as the mass parameter m must appear explicitly in the denominator to compensate for the engineering dimension of the interaction. Indeed, as argued in ref. [78], the singular high-energy limit of such amplitudes is consistent with the expected composite nature of higher-spin particles.
Let us now add a new perspective to this puzzle. At higher spins there seems to also be a tension between the double copy and locality. Indeed, consider the double copy (4.11) of the gauge-theory amplitudes with spins s 1 = 1 and s 2 = 1/2. They both contain only physical poles and satisfy color-kinematics duality, so we have a consistent double copy that gives a local gravitational amplitude for a massive Rarita-Schwinger field (which has the interpretation of a massive gravitino in a gravitational theory with broken supersymmetry [35,79]). At the same time, let us apply the double copy (4.11) to gauge theories with spins s 1 = 3/2 and s 2 = 0, respectively. In order to reproduce the above amplitude via 2) the gauge-theory amplitude with the spin-3/2 particle must contain a spurious pole: If we follow refs. [78,81] and attempt to cancel out the pole by additional contributions free of physical poles, this will likely non-trivially interfere with color-kinematics duality and therefore the double copy. Even if one manages to find a color-dual gauge-theory amplitude for massive spin-3/2 particles, it must introduce a correction to the double copy (6.2) such that it becomes seemingly impossible to manifest the trivial identity s = 1 + 1 = 3/2 + 1/2 = 2 + 0 by various partitioning of the spins on the two sides of the double copy.
Having pointed out the above potential problems for higher-spin double copies, let us end on a positive note. We know that the double copy is an intrinsic feature of string theory. Not only does it follow from the work of KLT [14], but also more recent results show that the double-copy structure permeates all types of strings: the bosonic, heterotic, closed and open superstrings [109][110][111][112][113][114][115]. As string theory contains a multitude of higherspin massive particles we can conclude that some version of a double copy must exist for these particles. While the infinite tower of massive modes may non-trivially conspire and cause technical challenges, there are good reasons to expect that one can isolate individual modes by investigating the residues of massive poles in the tree-level amplitudes. Indeed, it would be interesting to initiate a systematic study using string theory and see if the modes can be isolated without disturbing the double-copy structure. We leave this question and those raised above for the future.

A Axiodilaton gravity in four dimensions
Here we review the massless gravitational theory obtained as a double copy of pure Yang-Mills theory. The string-theoretic origin of the KLT double copy [14] implies that this theory (sometimes referred to as N = 0 or NS-NS gravity [116]) is Einstein's gravity coupled to a dilaton scalar φ and a Kalb-Ramond two-form B µν . Its Lagrangian is given e.g. in ref. [117]: where H λµν = ∂ λ B µν +∂ µ B νλ +∂ ν B λµ is the fully antisymmetric field strength for B µν . The double-copy relation of this theory to pure Yang-Mills in a classical setting was investigated in refs. [43,51].
In four dimensions we it is convenient to dualize the two-form B µν to an axion pseudoscalar a through a type of the Legendre transform where the Levi-Civita tensors are defined through the Levi-Civita symbol as E λµνρ = √ −g λµνρ , E λµνρ = g λα g µβ g νγ g ρδ E αβγδ = λµνρ √ −g , 0123 = − 0123 = 1. (A.3) Now the field equation for a sets λµνρ ∂ λ H µνρ to zero, which means that H is a closed and therefore exact three-form. Having thus guaranteed H = dB on the equations of motion, we can now consider H λµν as a separate field. Its field equations are Then for both terms in the axion Lagrangian (A.2) we have e −4φ H λµν H λµν = 4e 4φ λµνρ λµνσ ∂ ρ a ∂ σ a = −24e 4φ ∂ µ a ∂ µ a, (A.5a) −4H λµν E λµνρ ∂ ρ a = −8e 4φ λµνρ λµνσ ∂ ρ a ∂ σ a = 48e 4φ ∂ µ a ∂ µ a, (A.5b) so that the full Lagrangian (A.1) becomes equivalent to At this point, we could rescale the dilaton and axion by κ/(2 √ 2) to fix the normalization of their kinetic terms to 1/2 and deal with a perturbation theory of two real scalars. The four-dimensional double copy, however, makes it preferable to introduce a complex scalar field [118,119]  which we refer to as the axiodilaton. This change of variables produces the following Lagrangian (see e.g. [120]) where all energies p 0 = p 2 + m 2 are taken positive. For instance, in the complex case the expectation value 0|Ψ(x)|ψ( p, a) with an electron/quark state vanishes, but in the real case of a neutrino state it givesv a p e −ip·x , as shown above. Under the Majorana projection, the Dirac Lagrangian (3.11) becomes 3) The resulting three-point vertex is composed of two possible contractions where for simplicity we assumed the adjoint representation of the fermion and used the antisymmetry of the structure constants.

C Gravitational coupling of Dirac fermion
In section 3.2 we used the three-point vertex (3.14) with which a Dirac particle interacts gravitationally at lowest order in κ. Let us here derive its on-shell part, assuming that the end result must look like Here the linearized stress-energy tensor T µν should be symmetric and conserved on the equations of motion, ∂ µ T µν ∂ ν T µν 0. We wish to obtain it from the Noether energymomentum tensor in flat space which is not symmetric. A known remedy for that is to add a Belinfante contribution [122,123] T µν = T µν N + ∂ λ B λµ ν , B λµ ν = 1 2 S λ µν that is constructed using the Noether intrinsic spin tensor The needed properties of the resulting tensor T µν follow from the conservation of the total angular momentum tensor J λ µν N = x µ T λν N − x ν T λµ N + S λ µν N . Explicitly, the stress-energy tensor becomes Coupling it to the linearized graviton via eq. (C.1), we arrive at the interaction vertex (3.14).

D Cheung-Remmen gravitational vertices
In this appendix we summarize the perturbation theory for general relativity formulated by Cheung and Remmen in ref. [124]. They discovered that after a suitable gauge fixing the Einstein-Hilbert action can be rewritten in a form, in which the infinite tower of gravitational vertices is generated by a single geometric series: Similarly to the Landau-Lifshitz formulation of general relativity [125,126], the main variable here is the perturbation of the "gothic inverse metric" where the Lorentz indices are raised and lowered with the flat background metric η µν . At the linearized level, h µν corresponds to the trace-reversed version h µν − hη µν /2 of the more conventional graviton κh µν = g µν − η µν . The graviton propagator that follows from the action (D.1) is h αβ K µσ 1 νσ 1 µσ 2 νσ 2 µσ 3 νσ k (p 1 , p 2 )η νσ 3 µσ 4 η νσ 4 µσ 5 . . . η νσ k−1 µσ k , K αβ γδ µν (p, q) = − p α q γ − q α p γ − 1 2 q 2 η αγ η βµ η δν − 1 2 (p − q) α q β η γµ η δν . (D. 4) where every pair of indices (µ i , ν i ) should in principle be symmetrized, but this symmetrization can be left implicit as long as the symmetric propagators (D.3) and wavefunctions ε µν = ε µ ε ν are used. As a consistency check, we verified that the above vertices reproduce the three-, four-and five-graviton amplitudes obtained by the standard KLT relations (mimicked by eq. (5.28)). Moreover, this choice of variables is convenient for considering gravitationally covariantized matter Lagrangians. Namely, the scalar kinetic term only contrubutes to the three-point vertex (3.6), so vertices with two scalar fields and more than one graviton come exclusively from the mass term, which involves For vector bosons it is the mass term that is only contributes to the three-point vertex (3.19), and the kinetic term generates higher-point vertices, such as the quartic vertex