Scattering of Massless Particles: Scalars, Gluons and Gravitons

In a recent note we presented a compact formula for the complete tree-level S-matrix of pure Yang-Mills and gravity theories in arbitrary spacetime dimension. In this paper we show that a natural formulation also exists for a massless colored cubic scalar theory. In Yang-Mills, the formula is an integral over the space of n marked points on a sphere and has as integrand two factors. The first factor is a combination of Parke-Taylor-like terms dressed with U(N) color structures while the second is a Pfaffian. The S-matrix of a U(N)xU(N') cubic scalar theory is obtained by simply replacing the Pfaffian with a U(N') version of the previous U(N) factor. Given that gravity amplitudes are obtained by replacing the U(N) factor in Yang-Mills by a second Pfaffian, we are led to a natural color-kinematics correspondence. An expansion of the integrand of the scalar theory leads to sums over trivalent graphs and are directly related to the KLT matrix. We find a connection to the BCJ color-kinematics duality as well as a new proof of the BCJ doubling property that gives rise to gravity amplitudes. We end by considering a special kinematic point where the partial amplitude simply counts the number of color-ordered planar trivalent trees, which equals a Catalan number. The scattering equations simplify dramatically and are equivalent to a special Y-system with solutions related to roots of Chebyshev polynomials.


Introduction and Summary of Results
In 2003, Witten proposed a formula for the tree-level S-matrix of Yang-Mills in four dimensions, as an integral over the moduli space of certain rational maps from a n-punctured sphere to twistor space [1]. Shortly after, Roiban, Spradlin and Volovich (RSV) studied the formula in momentum space and gave non-trivial evidence for its validity [2]. In 2012, an analogous construction for gravity in four dimensions was found [3,4]. A natural question is whether similar constructions exist in arbitrary dimension.
In recent work [5] we presented compact formulas for the complete tree-level S-matrix of Yang-Mills and gravity theories in any dimension. In their simplest form, both formulas can be written in a unified manner as Tr(T a 1 T a 2 · · · T an ) (σ 1 − σ 2 ) · · · (σ n − σ 1 ) + . . .

2−s
Pf Ψ s (1.1) with s = 1 for Yang-Mills and s = 2 for gravity. In this formula T a denotes the U (N ) color group generators, and the ellipsis means a sum over all permutations of labels modulo cyclic ones. Here s ab = (k a + k b ) 2 and by defining e ab = ( a + b ) 2 and d ab = ( a + k b ) 2 with the understanding that the polarization vectors 's are null, one can write Ψ, which is a 2n × 2n antisymmetric matrix, as (1.2) for 1 ≤ a, b ≤ n, and the block Ψ a,b+n follows from the antisymmetry of the matrix. Pf Ψ ≡ (−1) i+j σ i −σ j PfΨ ij ij for 1 ≤ i < j ≤ n is called the reduced Pfaffian of Ψ, where Ψ ij ij denotes the non-singular matrix obtained by removing columns {i, j} and rows {i, j} from Ψ 1 . The meaning of the symbol and illustrations on how to use the formula explicitly are reviewed in section 2.
Here we would like to consider (1.1) not only as a convenient way to write Yang-Mills and gravity scattering matrices in a unified way but also as a definition of the S-matrix for particles of spin s. This means that the case s = 0 should correspond to a scalar theory. In this paper we show that this is indeed the case.
In order to make the claim more precise, recall that the formula for pure gravity can be slightly generalized [5] by replacing the integrand by the product of two independent Pfaffians, each with its own choice of gauge for polarization vectors (Pf Ψ( , k, σ)) 2 → Pf Ψ( , k, σ) × Pf Ψ(˜ , k, σ). (1.3) As proven in [6], the corresponding formula is obtained by applying the Kawai-Lewellen-Tye (KLT) relations to two copies of the Yang-Mills formula with polarizations and˜ . It is well known that the result gives amplitudes with gravitons coupled to dilatons and B-fields [7]. This suggests that when we set s = 0 in (1.1) to get a scalar theory, the integrand can also be generalized as Tr(T a 1 T a 2 · · · T an ) σ 12 σ 23 · · · σ n1 + . . .
(1.4) where σ ab denotes σ a − σ b (this notation is used here and in the rest of the paper in order to keep formulas more compact). Note that while the original squared factor depends on a single color group U (N ), the new factor has, in general, a different color group U (Ñ ) with T b as its generators.
This naturally motivates the study of a theory of scalars in the adjoint of the product of two different color groups U (N ) × U (Ñ ). The simplest possibility is the theory with only cubic interactions of the form, − f abcfa'b'c' φ aa' φ bb' φ cc' (1.5) where f abc andf a'b'c' are structure constants of U (N ) and U (Ñ ) respectively. In this note we show that M (0) n defined as Tr(T a 1 T a 2 · · · T an ) σ 12 σ 23 . . . σ n1 + . . . Tr(T b 1T b 2 · · ·T bn ) σ 12 σ 23 · · · σ n1 + . . .
(1.6) gives the full tree-level S-matrix of such scalar theory in any dimension.
All of the above also leads to the natural conclusion that the factors C U (N ) ≡ Tr(T a 1 T a 2 · · · T an ) σ 12 σ 23 . . . σ n1 + . . . and E ≡ Pf Ψ( ) (1.7) are interchangeable as they give rise to physical theories in the process. This is a colorkinematics correspondence which is valid for individual solutions to the scattering equations. More precisely, the formulas with s = 0, 1, 2 are closely related. Starting from the formula for the scalar theory, with integrand C U (N ) × C U (Ñ ) , if we replace C U (N ) (or C U (Ñ ) ), by the Pfaffian E (or E˜ ), we get the Yang-Mills formula with color group U (Ñ ) (or U (N )); if we further replace C U (Ñ ) (or C U (N ) ) in the Yang-Mills formula by another copy of the Pfaffian E˜ (or E ), we arrive at the gravity formula. We summarize the relations by the following diagram Tr(T b α(1)T b α(2) · · ·T b α(n) )M (0) n (α(1), α(2), . . . , α(n)), (1.8) where M Tr(T a β(1) T a β(2) · · · T a β(n) )m (0) n (α|β), (1.9) into what we call 'double-partial' amplitudes, denoted as m (0) n (α|β). There are two natural specializations of our formula which give rise to the two cases of partial amplitudes. The first case is (for simplicity we write the explicit formula with the canonical ordering α = I), (1.10) By definition M (0) n must be given by the sum over trivalent U (Ñ )-color-ordered Feynman diagrams where each vertex is dressed with a structure constant f abc of U (N ).
The second one gives double-partial amplitudes . (1.11) In section 3, we show that m (0) n (α|β) is given by the sum over all trivalent graphs which have two planar embeddings, one consistent with the α ordering and the other with the β ordering. In the special case α = β, the formula gives rise to the sum over all planar (ordered) trivalent graphs with weights given by the product of all scalar propagators. In section 5 we focus on the kinematic regime where each trivalent graph evaluates to unity. This means the double-partial amplitude with n particles evaluates to the total number of planar trivalent graphs which is C n−2 , the (n − 2) th Catalan number. We find that on this special kinematics the scattering equations simplify dramatically and are reduced to a special Y-system whose solutions are related to the roots of Chebyshev polynomials.We make a conjecture for what the integrand of our formula (1.11) gives when evaluated on each solution.
Having a direct connection to scalar trivalent graphs and the fact that the color and kinematic factors in (1.7) are interchangeable indicates a link to the color-kinematics duality discovered by Bern, Carrasco and Johansson (BCJ) in 2008 [8]. In fact, we are able to show that an expansion of the kinematic factor E analogous to that of the color factor C U (N ) exists. This leads to a formula for the Yang-Mills amplitude as a linear combination of double-partial amplitudes, m (0) n (α|β). Using their expansion in terms of trivalent graphs, our formula leads directly to the BCJ color-kinematics duality. Furthermore, by using our transformation from Yang-Mills amplitudes to gravity amplitudes one finds a new proof for the BCJ double-copy relations [8]. These facts will be discussed in Section 4. Section 6 is devoted to consistency checks of our formula by showing that (1.6) has correct soft limits and factorization limits. There we also give a proof for the structure of the double-partial amplitudes (1.11). We end in section 7 with conclusions and discussions.

Details of the Formula and Examples
Let us discuss the precise definition of all the elements entering (1.11) and then show how explicit computations are carried out. Let us denote by {σ 1 , σ 2 , . . . , σ n } the position of n punctures in the complex plane. Below we will see that our formula has an SL(2, C) invariance which means that we are dealing with n punctures on CP 1 . One of the main ingredients of the formula are the scattering equations [6,9] b =a which connect the space of kinematic invariants defined by scalar products of the momenta of external particles {k µ a } and the puncture locations {σ a }. Using the fact that in any physical process momentum is conserved and all particles are on-shell (i.e. k 2 a = 0), one can show that the equations are SL(2, C) invariant and only n − 3 are linearly independent. Here SL(2, C) acts as usual A simple way of imposing the support of the scattering equations is by noticing that is independent of the choice {i, j, k} and hence permutation invariant [6].
In Section 1, we discussed color decomposition by using traces of products of color group generators. However, in practice it is more convenient to use an alternative color basis proposed in [10] which contains only (n − 2)! elements. It fixes the position of two particles, e.g. 1 and n while permuting the remaining n − 2 labels. The traces are then replaced by where α ∈ S n−2 , and similarly forc α . C U (N ) and C U (Ñ ) can be decomposed in terms of c α andc α respectively. In this color basis, partial amplitudes M (0) n (α) and double-partial amplitudes m (0) n (α|β) are identical to those in the trace basis, but with the position of 1 and n fixed.
Combining these with the Faddeev-Popov Jacobian obtained by fixing the SL(2, C) redundancy acting on the σ s one finds that the scalar amplitude, (1.6), becomes where in the new color basis α(1) = β(1) = 1, α(n) = β(n) = n, and the sum is over all the solutions to the scattering equations and (2.6) is a corank 3 matrix. To get the reduced determinant det Φ, one removes any three rows {i, j, k} and any three columns {p, q, r} to get a reduced matrix, whose determinant, which we denote as |Φ| ijk pqr , is non-vanishing. Then This matrix was first encountered in [3] as a natural analog of Hodge's MHV gravity formula [11]. Finally, it is important to mention that an inductive algorithm for solving the scattering equations is known [6]. The total number of solutions for generic kinematics is (n − 3)!. In section 5 we show that in some special situations the equations simplify dramatically and all solutions can be obtained analytically.

Examples
The simplest example is the three particle amplitude. In this case one finds that no equations have to be solved as all σ variables can be fixed using the SL(2, C) invariance. Moreover, det Φ = 1/(σ 12 σ 23 σ 31 ) 2 which cancels the Parke-Taylor like factor squared in the integrand of (2.5) to give M which is the correct answer for a cubic scalar theory. Next we compute the four particle amplitude. The scattering equations become a single equation for one variable. Solving the scattering equations with σ 1 = 0, σ 2 = 1, σ 3 = ∞ gives σ 4 = s 24 /s 34 . Let us define s 12 = s, s 23 = t, s 13 = u, and for four particles, we adopt the standard convention for labelling the color factors using the s, t, u channels: and similarly forc s ,c t ,c u . We also denote the ordering (1324) as P (the canonical ordering is denoted as I). Computing det Φ = |Φ| 123 234 /(σ 12 σ 2 23 σ 31 σ 34 σ 42 ) and plugging the solution into the amplitude one gets as expected for a color-dressed cubic theory amplitude. In the last equality we have used the Jacobi identities c s − c u − c t =c s −c u −c t = 0. We have also checked explicitly that formula (1.6) gives the correct five point amplitude. The full amplitude can be written as a term proportional to the partial amplitude M (0) 5 (1,2,3,4,5), plus five terms related by permutations, and each partial amplitude can be decomposed into six double-partial amplitudes, (2,3,4), (2.11) where we have denoted the orderings as I = P 0 , (13245) = P 1 , (12435) = P 2 , (14325) = P 3 , (13425) = P 4 , (14235) = P 5 . Let us first consider how to compute these double-partial amplitudes. In arbitrary dimensions the scattering equations give rise to an irreducible quadratic polynomial. We have checked that by summing over solutions as dictated by (2.5) one reproduces the desired answer for the pair of two canonical orderings, m By plugging them into the combination inside the bracket of (2.11) and after using Jacobi identities repeatedly, one gets the expected partial amplitude M (0) 5 (1, 2, 3, 4, 5). By summing over permutations, (2.11) gives the full five-point amplitude.
In four dimensions, something special happens and when using spinor variables the quadratic polynomial from the scattering equations factorizes. In the terminology of Yang-Mills theory, one solution gives rise to the MHV sector while the other gives the MHV sector. In our computation we have to add up both solutions. It is easy to find the explicit formula for an arbitrary choice of permutations, m  Finally, we have checked that our formula also reproduces the sum of planar trivalent diagrams up to eight particles. Just as in previous cases, it is crucial to sum over all (n − 3)! solutions to find the amplitude.
In the next section we give an interpretation of each m  Plugging this into (2.14) immediately gives m

Double-Partial Amplitudes
Double-partial amplitudes are the building blocks of scalar amplitudes. In section 4, we show that they are also building blocks of Yang-Mills and gravity amplitudes. In this section, we show that double-partial amplitudes can be expanded as sums of trivalent graphs. Also, using a property of scattering equations called KLT orthogonality, we prove that the matrix with double-partial amplitudes as entries equals the inverse of the KLT matrix.

Trivalent Graph Expansion
Here we study double-partial amplitudes as generating functions of sums of scalar diagrams. We first recall the formula defining double-partial amplitudes The sum is over all (n − 3)! solutions to the scattering equations.
The simplest examples are those for n = 3, 4 and n = 5 given in the previous section. Let us recall some of the results in order to motivate the proposal for the meaning of m From these examples it is easy to see that when both permutations in m (0) n (α|β) are the same then the answer is a sum over all color-ordered trivalent graphs, each contributing the product of its propagators. When the two permutations are different it gives a subset of terms appearing in the formula for m(α|α). In extreme cases like m (0) 5 (I|13524) the subset is the empty set. A closer look at the relation between the diagrams that contribute to a given case straightforwardly motivates the following proposition.
Proposition: The function m (0) n (α|β) computes the sum of the collection of all trivalent scalar diagrams that can be regarded both as α-color-ordered and β-color-ordered, where each diagram's contribution is given by the product of its propagators. In other words, only diagrams that belong to the intersection of both sets contribute to the double-partial amplitude.
More explicitly, let T (α) denote the set of α-color-ordered diagrams and T (β) the set of β-color-ordered ones. Then where the integer n flip (α|β) is defined below and s e = P 2 e where P e is the momentum flowing along the edge e in the set of edges, E(g), for the Feynman diagram g. In particular, whenever We give a proof of this proposition in section 6 by showing that m (0) n (α|β) has the same soft limits and factorizations properties as the sum over the corresponding Feynman diagrams. In the rest of this section we present an efficient description of the diagrams that appear in a particular amplitude as well as a procedure to determine n flip (α|β).
Consider a particular double-partial amplitude m (0) n (α|β), without loss of generality take α to be the identity permutation, i.e., α = I. Start by drawing a disk with n nodes sitting on the boundary in the ordering α. Then link the n nodes together with a loop of line segments according to the ordering β. Generally the segments intersect each other in the middle of the graph. It is convenient to introduce the following terminology: A subset of all n points which are consecutive with respect to the α (or β) ordering will be said to be α-consecutive (or β-consecutive).
The way to compute m n (α|α) is given as above). If the set can be extended, e.g., by adding i − 1 while still remaining β-consecutive then take its maximal extension.
Assuming {i, i + 1, . . . , i + r} is maximal already, redraw the graph by moving all points in the set along the boundary of the disk, until they are close to each other 2 . The other external points must be kept fixed. If the lines coming out of i and i + r intersect then give a name to identify the intersection point, e.g. R (If the lines do not intersect then go back to the original graph and move on to the next set of both α-and β-consecutive external points).
Assuming the lines intersect, note that {i, i + 1, . . . , i + r, R} form a convex polygon. Now remove the polygon from the graph, bring R to the boundary of the disk and treat the new graph, which has R as an external point, as a new problem and repeat the procedure. If at any given point one fails to find a set which forms a polygon with at most one internal point after trying all of them, then m This procedure comes to an end when one finds a graph where both orderings agree completely and therefore gives a single polygon with only 'external' edges.
Finally, after completing all iterations and obtaining a list of all polygons found in the process one computes for each polygon its corresponding sub-amplitude. These sub-amplitudes can be computed just as regular amplitudes, where both orderings coincide, because for scalar particles it is straightforward to go off-shell.
The double partial amplitude, m n (α|β), is then given by the product of all sub-amplitudes times a propagator for each internal point that was obtained in the process of removing polygons.
In figure 1 we present m and propagators Putting all of them together gives, up to an overall sign, Now we give the rule to determine the overall sign of the results. First define the orientation of the disk by the ordering α (Figure 2 (a)), and define the orientation of the loop of segments by the ordering β (Figure 2 (b)), which induces an orientation in every convex polygon. The rule is as follows: (1) each polygon with odd number of vertices contributes a plus sign if its orientation is the same as that of the disk and a minus sign if opposite, (2) each polygon with even number of vertices always contributes a minus sign, and (3) each intersection point contributes a minus sign (Figure 2 (c)). Then the product of all these signs determines the overall sign of the double-partial amplitude relative to its corresponding scalar diagrams.
In our example one finds that the three external polygons give a minus sign while the internal one gives a plus sign. There are three propagators, and each gives a minus sign (see figure 2). Altogether one finds six minus signs and therefore the overall sign is plus. This means that (3.7) is indeed the right answer. , and think of each cubic vertex as a triangle whose vertices are glued either to a node or to a vertex of another triangle. Then the sign can be obtained by the same rules as above. This is because in replacing each convex polygon by a specific diagram, the number of triangles and propagators are fixed by the polygon, and the triangles therein inherit the same orientation. Whenever a triangle picks up a minus sign, the ordering of its three vertices are flipped in β as compared to α 3 . In this way, the sign of the double-partial amplitude can be interpreted as the number of propagators n − 3 and the number of ordering flips n flip (α|β) in the cubic vertices of any specific diagram contained therein. This explains the factor
In order to state what KLT orthogonality is let us define a function defined for any two given solutions to the scattering equations {σ I } and {σ J } evaluated on the I th solution while U (J) β stands for the following expression evaluated on the J th solution.
As proven in [6], the following holds for any I and J, This is known as KLT orthogonality. For our purposes, only two facts, also explained in detail in [6], are necessary. The first is that S[α|β] is only a function of the kinematic invariants s ab . The second is that (3.14) It is interesting to note that the dimension of space of permutations of labels 2, . . . , n−2 (for which both U and V vectors form a basis, evaluated on any solution), and the dimension of the solution space, are both (n−3)!. This allows us to define the following (n−3)! × (n−3)! matrices based on U, V vectors normalized with respect to the inner product (3.10) where the rows (superscripts) are labeled by solutions and columns (subscripts) by permutations (or orderings), and KLT orthogonality is the simple statement that the product of these three matrices gives the identity matrix I in solution space, It is obvious thatV is invertible, thus by multiplyingV T and (V T ) −1 from left and right respectively, we obtain the identity matrix in permutation space, where by (3.1) we find the right-hand-side of the second equality is precisely given m scalar The inverse of the KLT matrix has also been discussed in [13] where it was related to the field-theory limit of string disk integrals 5 , following computations in [14,15]. As explicit examples, the inverse was given in [13] for up to seven points, and the result agrees with that of the double-partial amplitudes. It would be interesting to explore the connections further.

Color-Kinematics Duality
In the introduction we illustrated how our formulation relates scalar-, gluon-and gravitonamplitudes by simple transformations (C → E orC →Ẽ or both). This replacement occurs solution by solution of the scattering equations. More explicitly, (4.1) Very interestingly, in [16] Hodges made the observation that using twistor diagrams the statement that "gravity is the squared of Yang-Mills" is replaced by "gravity times φ 3 is the square of Yang-Mills". Indeed, our construction shows that this is precisely true solution by solution of the scattering equations! Also mentioned in the introduction is the fact that the color factor C and the kinematic factor E have very similar properties and this is the reason why they are exchangeable. It is very natural to suspect that this color-kinematics duality must have a connection to the color-kinematics duality introduced by Bern, Carrasco and Johansson (BCJ) [8].
Indeed, in this section we show that these simple transformations are equivalent to the BCJ doubling procedure to obtain gravity amplitudes from Yang-Mills ones.

Expansion of the Pfaffian
Here we show that the Pfaffian can be expanded in a way similar to the expansion of the color factors. We choose to expand the color factor, C (C), in terms of the basis given in (2.4), and a similar formula holds forC. This form of C hints that a similar form for E should exist. More explicitly, there must exist functions, denoted as n, which only depend on kinematic data { µ a , k µ a }, such that Straightforwardly expanding the reduced Pfaffian Pf Ψ( , k, σ) leads to an expression very different from the right hand side. In fact, at first sight it seems difficult to rewrite it into the form of (4.3). Luckily, not only (4.3) holds but the proof is very simple again thanks to the KLT orthogonality [3,6]. In order to prove that the expansion in (4.3) exists it is enough to rewrite (3.13) as (recall the definitions in (3.11),(3.12)) then to multiply it by Pf Ψ(σ J ) and sum over J get Using (3.10) and (3.14) one finds that the left hand side of (4.5) is given by where note that α and β are permutations of labels 2, 3, . . . , n−2. The sum over (n − 3)! solutions appearing in this equation is nothing but the Yang-Mills partial amplitude M (1) n (1, β, n, n − 1). Therefore, on the support of scattering equations, This concludes the proof of (4.3) since we have found an explicit form of the numerators in terms of only external kinematic invariants. More explicitly, (4.8) The form given in this proof (which was also discussed in e.g. [17]) is not unique and not very useful for practical computations as it is tautological in nature. However, it shows that expansions of the form (4.3) do exists if we regard the M n 's in (4.8) as functions of external data, and this is all we need in order to make the connection to the BCJ doubling construction. As we will see shortly, the (n−2)! n's play the role of a basis for BCJ numerators, and in practice it is always possible to derive local expressions for them by carefully using scattering equations when expanding E; explicit expressions for such BCJ numerators can be found in [18,19]. Now we can unify (4.2) and (4.3) by denoting both c and n as e in all three theories. Then if we expand the product of the two summations, each term is of the form e 1α(2)···α(n−1)nẽ1β(2)···β(n−1)n m (0) n (1α(2) · · · α(n − 1)n|1β(2) · · · β(n − 1)n) ≡ e αẽβ m (0) (α|β). where for e,ẽ we have s kinematic numerators and 2−s color numerators, with s = 0, 1, 2.
In this form the color-kinematics correspondence of our formula becomes more transparent: the σ-independent color and kinematic factors c α and n α are on a equal footing, and by exchanging them one relates amplitudes of scalar theory, Yang-Mills and gravity.

Relation to BCJ Color-Kinematics Duality
To see that (4.10) actually gives a representation which respects the BCJ color-kinematics duality, we can start by exchanging the two summations. First note that the union of all sets T (α) ∩ T (β) with α, β ∈ S n−2 is the complete set of trivalent diagrams with n legs, T n . For each trivalent diagram g, we define the set of pairs of orderings 6 that can generate g, P O(g), and it is obvious that (α, β) ∈ P O(g) ⇔ (β, α) ∈ P O(g).
The formula can be written as a sum of all cubic diagrams, each with propagators and two numerators, Note that the signs satisfy a composition identity, which simply follow from its definition in the previous section: for any γ ∈ O(g), (−1) n flip (α|β) = (−1) n flip (α|γ)+n flip (γ|β) , thus the double sum over α, β factorizes into two sums, 12) where we have defined e g andẽ g to be the combinations in the two brackets respectively. Note that this definition depends on the choice of γ that is fixed for the sums. Different choices of γ may leads to different signs, and when {e g } corresponds to color factors, this is the standard sign ambiguity for generic trivalent diagrams (in the special case when e g receives contributions from only one e α , we can choose γ = α).
In this way, up to possible signs, our formula defines numerator factors e g ,ẽ g for every trivalent diagram, and the physical quantity which is the product of the two copies, e gẽg , has no such ambiguities. We now argue that these numerators automatically satisfy Jacobi-like identities. Given three trivalent trees which differ only by a four-particle subdiagram, as illustrated in figure. 3, we need to prove that the numerators e g (similarly forẽ g ) satisfy a Jacobi identity, e gt = ±(e gs − e gu ), (4.13) and here we emphasize that the relative sign of the two terms on the right-hand-side is important. Let us denote the four trees attached to the four legs of the subdiagram as A, B, C and D, and we first consider the simple case where particles 1, n are attached to two different trees, which, without loss of generality, we assume to be A, D respectively (see figure. 3). The idea is that the set of orderings for the complete diagram can be obtained from putting together the four sets corresponding to the four trees. For the tree B(C), we remove the label of the internal leg from each α ∈ O(B) (α ∈ O(C)), and define the new set which only has labels for external particles as O * (B) (O * (C)); we do the same thing for A and D, but in the definition of O * (A) we only keep the orderings with 1 in one end, and in O * (D) only those with n in the other. It is straightforward to find that, where γ t,s,u ∈ O(g t,s,u ) are the elements chosen in the definition of e gt,s,u , and in the second line the two summations give e gs and e gu respectively. The sign for an individual e does not matter, but the relative sign of the two terms, (−1) n flip (γs|γu) , is unambiguous: from our definition of n flip , this sign is determined by looking at the relative orientation for each cubic vertex of g t in the ordering γ s and γ u , and the only difference appears in the cubic vertex with three internal legs, where B and C are exchanged from γ s to γ u , thus the relative sign of e gs and e gu is always −1. Figure 3. Trivalent diagrams g t , g s , g u when particle 1, n are contained in two different trees, e.g. A, D, attached to the four-particle subdiagram. Red and blue regions correspond to O(g s ) and O(g u ) respectively, the union of which gives O(g t ).
The case when particle 1, n are attached to the same tree, e.g. A , see figure. 4, can be argued similarly. To keep the expressions short, we denote the ordering corresponding to tree A by two groups of ellipses which contain 1 and n respectively. The sets of orderings for g s , g t , g u can be written as  where γ ∈ O, and the three sign factors all give −1. Note that the second and the fourth term cancel with each other, and the first and the third term combine into, up to a sign, e gt , thus we again find e gt = ±(e gs − e gu ). Since the above statement holds regardless of e,ẽ being color or kinematic numerators, it guarantees that the representation of M (s) n as a sum of cubic diagrams given by our formula automatically respects color-kinematics duality. Moreover, the simple transformations between C (C) and E (Ẽ) in our formula correspond to the substitution of the set of e factors into a set of e factors which satisfy the same algebraic relations (orẽ intoẽ ), thus they are equivalent to the double-copy procedures that relate scalar, Yang-Mills and gravity amplitudes. Proofs of the existence of dualtiy-respecting numerators and the double-copy relations can also be found in [20] [18]. Given any three graphs, as the ones discussed above, g s , g u , g t , which only differ in a single propagator, i.e., one can freely redefine a kinematic set of numerators, say, {e g } according to e gs → e gs = e gs + ωs, e gt → e gt = e gt + ωt, e gu → e gu = e gu − ωu, (4.19) and the value of the formula remains invariant. This is simply because the set {ẽ g }, which can be either color or kinematic numerators, satisfy Jacobi-like identities.
In general the new numerators {e g } obtained by applying any number of transformations of this form, with completely arbitrary parameters w, will not satisfy Jacobi-like identities. Applying this to Yang-Mills amplitudes one can obtain formulas of the form where {n g } satisfy the Jacobi-like relations. Clearly, one can start with a gravity formula given by e g = n g andẽ g = n g and apply the same transformations that led to (4.20) on theẽ g factors to get This result can be interpreted as saying that the gravity formula (4.22) can be obtained by using a double-copy procedure applied to (4.20) and (4.21) and therefore conclude that the double-copy procedure works even when one of the two Yang-Mills numerators do not satisfy Jacobi-like identities 7 .

Special Kinematics: Generating Catalan Numbers
The Feynman diagrams of a colored cubic scalar theory are all possible trivalent, colorordered planar trees. Each internal edge of the tree is dressed with a propagator factor and its contribution to the amplitude is the product of all propagators in the tree. In spacetime dimensions large enough compared to the number of particles, the set of kinematic invariants At the special kinematic point it is easy to compute that all two particle kinematic invariants are {s i,i+1 = 1, s i,i+2 = −1} for i ∈ {1, . . . , n} (indexes are understood modulo n) and zero otherwise.
The scattering equations simplify dramatically and become The most direct way to solve the equations is to use n−4 equations to express n−4 σ's in terms of the remaining one, and when plugging into the last equation, it becomes a polynomial equation of the remaining σ, which one can solve to get all the solutions. Here we follow a different route which turns out to be more instructive. The equations are SL(2, C) invariant, and one can rewrite them in terms of the cross-ratios: where by definition u i,i±1 = 0 and u i,i = ∞. It is straightforward to see that the a th scattering equation at the special kinematic point becomes, Since the second factor does not vanish (unless σ a−1 = σ a+1 which represents a unwanted, singular kinematic point), we conclude that the scattering equations are equivalent to u 1,3 = u 2,4 = . . . = u n,2 , (5.7) which corresponds to a symmetric configuration with all cross-ratios built from four consecutive points equal to each other. Therefore, in terms of the cross-ratios, the scattering equations become trivial. On the other hand, by definition the cross-ratios satisfy constraints which correspond to the so-called "Y-system" equations [21] ( ), (5.8) and it turns out we can directly solve the combined system of (5.7) and (5.8)! Given the fact u i,i+2 is independent of i, we immediately see that for m = 2, 3, . . . , n−2 (m = 0, 1, n−1, n are trivial), the cross-ratios u i,i+m are independent of i. Thus we can define R m ≡ u i,i+m , and by a simple change of variable R m = Ym 1+Ym the equations in terms of the new variables become precisely the Y-system equations in the high temperature limit considered in [21], with the boundary condition Y 0 = −1, Y 1 = 0. It is well known that there are (n−1)/2 solutions, labeled by i, of the following form, m ) means that the full integrand is a rational function of the Y m variables. For general n and kinematic points, the integrand is a rather complicated function of cross-ratios; on the other hand, with the special solutions (5.10), as we have checked up to n = 20, the summand simplifies significantly, which we conjecture to be true for general n.
Conjecture: At the special kinematic point, the formula as a function of Y whereî = i for n even, andî = i − 1 2 for n odd. From this it is straightforward to see that the amplitude at the special kinematic point indeed counts the number of planar trivalent diagrams, which is the Catalan number.
Corollary: The double-partial amplitude at the special kinematic point is given by

13)
Proof: Consider the Dynkin diagram of A n−1 , which consists of n−1 vertices and n−2 edges connecting them in sequence. By the combinatoric definition of Catalan number, C k is the number of closed paths with length 2k from one end of the diagram to itself. Denote the vertices as 1, 2, . . . , n−1, then the adjacency matrix of the diagram is given by (A n−1 ) i,j = δ i,j−1 + δ i,j+1 . If we start from vertex 1, the number of closed paths with length 2n−4 is given by the (1, 1) component of the matrix A 2n−4 n−1 . It is straightforward to see that A n−1 has eigenvectors v i = 2 n (sin( πi n ), sin( 2πi n ), . . . , sin( (n−1)πi n )), and eigenvalues λ i = 2 cos( πi n ), for i = 1, 2, . . . , n−1, thus where note the summand is symmetric under i ↔ n−i (for n odd, the term with i = (n−1)/2 vanishes), thus one can replace the range of summation by i = 1, . . . , (n−1)/2 and multiply the result by 2. A simple rewriting gives (5.13), which concludes the proof. One of the most interesting future directions is to consider the scattering equations for general kinematics in light of the special kinematic calculation. It is possible that by rewriting the scattering equations in terms of cross-ratios, one may encounter more general Y-systems. It would be fascinating to find a physical interpretation of such system of equations, especially an interpretation of our formula as certain physical quantities, in analog of the minimal area as the free energy of the Y-system. Doing so may allow us to find simple expressions for the integrand (in particular a proof for the conjecture above), and possibly to connect scattering equations to some integrable system.

Comments on the Special Kinematics
As the first step towards more general kinematics, we study the cases where all kinematic invariants of the form s i,i+m with any fixed m are equal, as long as the diagrams contributing to the double-partial amplitude do not diverge, which is the most general kinematics cyclically symmetric with respect to the canonical ordering. Up to a rescaling, we can always set s i,i+1 = 1 for all i, then the simplest cases of this type is s i,i+1 = 1 and s i,i+m = −1 for some fixed m, while all the remaining s a,b are identically zero (when n is even and m = (n − 2)/2, we need to set s i,i+m = −2 due to momentum conservation). The above case correspond to m = 2.
However, since these configurations are in general highly singular with respect to the scattering equations, one may encounter problems if evaluating the formula directly on them. For example, in the case s i,i+1 = 1, s i,i+2 = − 1 3 and s i,i+3 = − 2 3 at 7 points, directly solving the equations shows that apart from three isolated solutions, there is a continuous region of infinite number of solutions, and naively our previous way of evaluating the formula fails. Even for our favorite special kinematic points m = 2, which has only isolated solutions and they give the correct result, the number of solutions seems to contradict the general counting, (n − 3)!. To understand these situations properly, one should regard these kinematics as a limit, i.e. starting by adding to the kinematics small deviations controlled by a scale and evaluate the formula, and then take the → 0 limit. To illustrate this, here we provide two examples.
First we give an interpretation to the known result for the kinematics s i,i+1 = −s i,i+2 = 1 at 6 points from this point of view. There we get two solutions, and if we fix SL(2, C) by setting {σ 1 , σ 2 , σ 3 } = {0, 1, −1} they are 15) and the formula gives − 27 2 on the first solution and − 1 2 on the second, thus adding up to −14 which matches the result calculated from Feynman diagrams. Alternatively, we can start from the configuration where the parameters x, y, z all deviate from −1 by some small values proportional to a scale . This configuration is generic enough to produce all (n − 3)! = 6 solutions, and whatever we choose the summation of the formula evaluating on these solutions gives the correct value matching the diagrams. As gradually approaches zero, we observe that the evaluation of the formula also approaches zero on three of the solutions. In these solutions, some of the σ i − σ i+3 become O( ), and their role is to compensate the infinitesimal s i,i+3 to produce a finite value that keeps the scattering equations satisfied. However, at = 0 terms with s i,i+3 disappear, leaving the equations un-balanced, and so these three sets of σ values are excluded. For those remaining three solutions, for infinitesimal the evaluation of the formula approaches − 27 2 on one of them and − 1 4 on the other two, and the two solutions that gives the same value approaches each other. Although collision of σ i and σ i+3 may still occur in this case, the scattering equations are still balanced at = 0, so all the three solutions remain. Since the multiplicity of the two identified solutions has already been taken into account by the Jacobian of the delta constraints, upon the limit = 0 we only see two distinct solutions and we only need to evaluate the formula once on each to produce the correct result.
For the second example, we look at a more interesting case which gives rise to a continuous region of solutions. This happens first time at 7 points, when we consider a generic cyclic symmetric configuration, the most general form of which is (after normalizing s i,i+1 = 1) parameterized by a single variable t. In this specific example we set t = 1 3 , and we are going to show that from the view of limit the formula still gives correct answer, which is 462 25 . To approach this kinematics from a generic data, we start with the configuration s i,i+1 = 1, ∀i, s 14 = p, s 25 = q, s 36 = u, s 47 = v, s 15 = x, s 26 = y, s 37 = z, (5.18) and obtain all the remaining kinematic invariants by momentum conservation, and we set all the seven parameters p, q, u, v, x, y, z to deviate from − 2 3 by some small values proportional to the scale . This in general gives all the (n − 3)! = 24 solutions. Again we let to approach zero, and we will observe three types of behaviors for the formula: (a) the formula approaches zero on 7 solutions. The collision of σ's again occurs in all these solutions, and at = 0 these solutions are excluded also due to the fact that the equations are no longer balanced. (b) there are 3 solutions upon which the evaluation of the formula remains finite and the solutions are regular. However, the summation of the formula evaluated just on these solutions does not yet add up to the correct result. (c) the evaluation of the formula diverges on all the remaining 14 solutions, and these solutions are also regular. However, they always add up to be a finite number and when combined with the evaluation upon the 3 solutions in case (b) the total summation gives the correct result 462 25 in the → 0 limit! On the other hand, if we work at = 0 right at the start, we will find only 3 isolated solutions, which match exactly with the → 0 limit of the 3 solutions in case (b). The 7 solutions in case (a) are already excluded. Moreover, here we will see a 1-dimensional continuous region of solutions to the scattering equations, and if we pick up the 14 solutions in case (c) and take the → 0 limit, we will see that they ultimately sit within this region. In other words, upon = 0, the original 14 isolated solutions emerge into a continuous region of solutions, which also gives a non-trivial contribution to the final result.
This continues to be true for generic t, where there is always a continuous region of solutions together with 3 isolated ones at = 0. Given a certain SL(2, C) fixing the 3 isolated solutions in case (b) are actually independent of the value of t, and it is possible to work out that the formula evaluated over these 3 solutions add up to On the other hand, the total value determined from the scalar diagrams is which only contains one physical pole corresponding to s i,i+1,i+2 = 0, and so all the poles in (5.19) are spurious. Subtracting the above two expressions, we know that the contribution from the continuous region should be (5.21) In particular, this contribution vanishes when t assumes one of the roots of its numerator. Actually upon such values one can verify that the continuous region of solutions is absent and all the solutions left are the 3 isolated solutions, which by themselves determines the correct result. Among these, t = 1 is what we have mainly studied in this section. If we choose to approach this limit from generic kinematics, we will find that in the neighborhood of the limit the formula remains finit on the 3 solutions in case (b), while it approaches zero on all the remaining solutions, which are ultimately excluded at = 0.
The above explains that although the formula naively fails for such singular kinematics, it gives the correct answer when regarded as a limit from generic kinematic data. Due to the existence of the continuous region of solutions and its non-trivial contribution to the final result as in (5.21), it is interesting to look for a way to determine the contribution from such region directly. A possible solution may be similar to the method proposed in [22].

Consistency Checks
In this section we perform the two standard consistency checks on tree-level amplitudes. The first is the behavior when one of the particles becomes soft. The second is the presence of a simple pole whenever a kinematic invariant of the form ( i∈I k i ) 2 for a set I vanishes with residue equaling to the product of lower-point amplitudes.
Since σ a 's are taken to be complex numbers in this paper, the delta functions imposing the scattering equations are in fact poles and all our integrals are contour integrals. In using the residue theorem, one finds that there is no contribution at infinity; for each term in the sum over i, j, there is a pole at σ n = σ 1 , and only for each of the terms with i = j, we have a pole at σ n = σ i which gives non-zero residue. Upon each pole only one term from a =n sn,a σn,a will contribute to the residue, giving rise to 2) where we have pulled out one structure constant of U (N ) from c α i and one ofŨ (N ) from c β j , and rewritten the remaining factors as full amplitudes M (0) n−1 with the color indices of particle i (and those of particle j for the first sum inside the bracket) being summed over; note that we have suppressed color indices of other particles, which are a l , b l for particle l. Let us consider the combination and by combining the sum over i and that over α ∈ S n−2 we get a sum over α ∈ S n−1 , i.e. permutations of 1, . . . , n−1; for each permutation there are two terms that differ only by a sign, thus the combination vanishes. This means the sum over i = 2, . . . , n−1 in the double sum of (6.2) gives minus the term with i = 1, and similarly for the sum over j = 2, . . . , n−1, thus (6.2) can be simplified as which is the correct soft behavior as one can see from Feynman diagrams. This can be compared with the soft limit of Yang-Mills full amplitude, n−1 (i c ). (6.5) The soft limit of the partial amplitude, M

Factorization
For the purpose of showing factorizations of the formula (1.1) in the scalar theory, it sufficies to have a look at only the kinematics singularity defined by with 2 ≤ n L ≤ n − 2, and we denote L = {1, . . . , n L } and R as its complement set, with n R = n − n L . Upon such a limit, only a subset of (n L − 2)! × (n R − 2)! from all the solutions are singular [6]. To make factorization manifest, we start by choosing a special redefinition where we regard v n−1 as being fixed to a specific value v * n−1 , and leave s as well as the remaining u's and v's as variables to be integrated over. Now the measure transforms to where u denotes the product of all u a with a ∈ L. We choose to fix {u 1 , u 2 , v n } to get rid of the SL(2, C) redundancy. Faddeev-Popov method in this gauge-fixing gives rise to a Jacobian 8 − 2u 1,2 v n−1,n (−s 4 + u 1 u 2 v n−1 v n ) s 2 v n−1 . (6.10) If s becomes infinitesimal in the neighborhood of singular solutions, we are able to approximate (6.10) and we can think of the u 1 in the parentheses above as (u 1 − 0) where the "0" corresponds to the puncture of the internal particle in the factorization (the same for The reason that s is constrained to be infinitesimal near singular solutions can be justified by the behavior of the delta constraints. Here we choose to eliminate the constraints corresponding to particles {1, 2, n}. Then for a ∈ R, we can expand the n R − 1 constraints with respect to s s s a,I R v a + s a,n L +1 v a,n L +1 + · · · + s a,n v a,n + O(s 3 ) = 0. (6.11) However, we only need n R − 2 of them for the right sub-amplitude, and so one delta function needs to be isolated to produce a constraint on the invariant mass k 2 I R , so that when it approaches zero a massless internal particle is emerged. For this purpose, we dress each constraint in (6.11) with a factor vavn,a svn and sum them up, we obtain − k 2 I R + s 2 (F (k, u, v) + O(s)) = 0, (6.12) with F some function independent of s. From this, we explicitly see that upon the factorization limit (6.7) we do have a subset of solutions which gives s → 0. In producing the constraint (6.12) from the original ones we get an additional Jacobian v n−1 v n,n−1 svn . Moreover, for a ∈ L, the n L − 2 delta constraints are expanded to − u 2 a s s a,1 u a,1 + · · · + s a,n L u a,n L + s a,I L u a + O(s) = 0. (6.13) Given (6.12), now we focus only in the regions where s is constrained to be infinitesimal, we are allowed to pick up only the leading terms in (6.10), (6.11), (6.12) and (6.13). Collecting these results together with (6.9), we conclude that the measure and the delta constraints approximates to −ds 2 n L a=3 du a n−2 ). (6.14) Then we go on with the summation of color-dressed Parke-Taylor factors. It is convenient to start with the form as in (1.8) and (1.9) where the color factors are expressed in terms of structure contants. With infinitesimal s, it is not hard to check that the leading terms are contributed solely by terms where both the two Parke-Taylor forms are of the pattern such that all labels in L sit in front of all the labels in R, which behave as 1 σ 1,α(2) · · · σ α(n L ),β(n L +1) · · · σ β(n−1),n σ n,1 −→ (−1) n L s −n L +n R +2 ( u) 2 (u 1 u 1,α(2) · · · u α(n L −1),α(n L ) u α(n L ) )(v β(n L +1) v β(n L +1),β(n L +2) · · · v β(n−1),n v n ) , where α is any permutation within the label set L\{1}, and β any permutation within the label set R\{n}. All the other Parke-Taylor forms are of higher order in s compared to (6.15), and so they are irrelavent in the factorization. For these leading terms, it is clear that each copy of their corresponding color factors also breaks into two parts, with one new fixed end arising in each part, and the two new ends are glued by a Kronecker delta with indices in the adjoint representation of the color groups {c} f a 1 ,a α(2) ,c 1 · · · f c n L −2 ,a α(n L ) ,c n L −1 f c n L −1 ,a β(n L +1) ,cn L · · · f c n−3 ,a β(n−1) ,an = a L ,a R ( {c } f a 1 ,a α(2) ,c 1 · · · f c n L −2 ,a α(n L ) ,a L )δ a L ,a R ( {c } f a R ,a β(n L +1) ,c n L · · · f c n−3 ,a β(n−1) ,an ).
As a consequence, each copy of the summation over Parke-Taylor forms exactly splits into the product of summations on the left part and that on the right part. Combining two copies of this summation together with (6.14), it is easy to see that if we regard the factors {u 1 , u 2 , v n−1 , v n } therein as {u 1 − u I L , u 2 − u I L , v n−1 − v n R , v n − v n R }, with u I L and v I R as punctures for the internal particles on the left part and right part, which are fixed to be zero, then we have an emergent SL(2, C) × SL(2, C) redundancy in these leading terms, acting on L ∪ {u I L } and R ∪ {v I R } respectively. When we integrate s 2 out, we see that M (6.17) We can also apply the above discussion to the double-partial amplitudes m (0) (α|β) with any orderings α and β. By (6.15), we see that m (0) (α|β) will have diverging leading terms in the form m (0) n (α|β) −→ m n R +1 (I R , α R |I R , β R ) (6.18) if and only if the labels contained in k I forms a consecutive subset in both α and β, and since m (0) (α|β) is a function of only kinematic invariants, it vanishes whenever there doesn't exist such a factorization channel. By studying factorizations recursively it is straightforward to observe that up to a sign m (0) (α|β) has the form as presented in (3.4). For the overall sign, we can just pick up any tree diagram g ∈ m (0) (α|β) and do n − 3 consecutive factorizations according to the propagators therein to fully factorize m (0) (α|β) down to n − 2 cubic vertices. On the one hand (6.18) indicates that each propagator gives rise to a minus sign, and on the other hand as in (3.2) it is easy to see that each cubic vertex gives +1 if the α and β orderings of the three labels are the same and −1 if flipped, so the overall sign is (−1) n−3+n flip . Since as argued in Section 3, this value is independent of the diagram g we choose, we have thus verified (3.4).

Discussions
In this paper we provided a unified description of the tree-level S-matrix of a colored massless cubic scalar theory, Yang-Mills and gravity. The new description manifests a relation between factors that contain the color information and factors that contain the polarization information. All amplitudes are written as integrals over the moduli space of an n-punctured where I (s) with s = 0, 1, 2 represent the integrand for scalars, gluons or gravitons respectively. The integrands, in their simplest forms, are given by This shows that solution by solution one has (I (1) ) 2 = I (2) I (0) . This connection between the square of Yang-Mills and the product of gravity with a φ 3 theory has the same structure as that found by Hodges using twistor diagrams in [16]. Very nicely, there is a third way to make this connection between Yang-Mills, gravity and φ 3 explicit. Recall that the KLT construction represents a gravity amplitude as a linear combination of products of partial Yang-Mills amplitudes. The sum is over certain pairs of (n − 3)! permutations and it is given schematically as It would be interesting to explore possible connections among all three descriptions.
Finding a direct proof of the equivalence of the two formulas for m (0) n (α|β) is an important problem which might also give a reason why the formulas for Yang-Mills and gravity amplitudes are strikingly similar to those of string amplitudes in the Gross-Mende or high energy limit [23]. Hints in this direction were already discussed in [6].
Clearly there are other pressing issues which we leave for future work. The first is finding a formula analogous to the reduced Pfaffian, Pf Ψ, which could accommodate the scattering of fermions. In dimensions less than twelve for gravity and less than eleven for Yang-Mills one could try and find a supersymmetric version of the reduced Pfaffian. In four dimensions, our formula has been shown to agree with the Witten-RSV formula for the scattering of gluons for up to eight particles and in all helicity sectors [5]. The agreement happens solution by solution of the scattering equations which means that the formulas are equivalent at the level of the integrand. It is well known that the Witten-RSV formula for gluons is only a set of components of an elegant N = 4 supersymmetric formula [1,2]. An analogous formulation for N = 8 supergravity is also known [3,4]. This suggests that a supersymmetric generalization of the Pfaffian exists in four dimensions and perhaps in higher. Finally, all the formulas obtained so far have been restricted to tree level S-matrices. Extending this formalism to loop level will likely involve higher genus Riemann surfaces and constrains on the theories. It would be fascinating to understand the class of theories whose S-matrices can be recast in this way.