Loop Equation and Exact Soft Anomalous Dimension in N=4 Super Yang-Mills

BPS Wilson loops in supersymmetric gauge theories have been the subjects of active research since they are often amenable to exact computation. So far most of the studies have focused on loops that do not intersect. In this paper, we derive exact results for intersecting 1/8 BPS Wilson loops in N=4 supersymmetric Yang-Mills theory, using a combination of supersymmetric localization and the loop equation in 2d gauge theory. The result is given by a novel matrix-model-like representation which couples multiple contour integrals and a Gaussian matrix model. We evaluate the integral at large N, and make contact with the string worldsheet description at strong coupling. As an application of our results, we compute exactly a small-angle limit (and more generally near-BPS limits) of the cross anomalous dimension which governs the UV divergence of intersecting Wilson lines. The same quantity describes the soft anomalous dimension of scattering amplitudes of W-bosons in the Coulomb branch.


Introduction
The loop equation was proposed initially in [1,2] as an alternative way to formulate, and possibly solve, the gauge theories (see e.g. [3] for a review). It has the conceptual advantage that it directly constrains the most basic observables, namely the Wilson loops. In lowerdimensional theories such as matrix models [4] and two-dimensional Yang-Mills theory [5][6][7][8][9], it has proven to be a powerful tool for solving the theories exactly. Unfortunately, solving the loop equation is much harder in higher dimensions and progress remains to be made.
In this paper we demonstrate the power of the loop equation in higher dimensions when used in conjunction with other non-perturbative techniques. Specifically, we consider intersecting 1/8-BPS Wilson loops in N = 4 supersymmetric Yang-Mills (N = 4 SYM) and compute their expectation values at finite coupling and finite N using a combination of the loop equation and supersymmetric localization [10,11].
The 1/8-BPS Wilson loop is a supersymmetric Wilson loops in N = 4 SYM which can be defined on a arbitrary contour on a two-sphere and preserves four fermionic charges. It was conjectured in [12,13] and later supported by supersymmetric localization [11] that its expectation value (as well as general correlation functions of any number of loops) conicides with that of the standard Wilson loop in two-dimensional Yang-Mills theory (2d YM) in a zero-instanton sector. Based on this, the expectation value of a non-intersecting 1/8-BPS loop was computed in [12,13] generalizing the famous result for the 1/2-BPS loop [10,14,15]. It was also used to derive multi-matrix models for correlators of local operators and nonintersecting BPS loops [16,17]. By now these results have been tested by numerous direct computations [18][19][20][21][22][23] and they provide convincing evidence for the equivalence between the BPS sector of N = 4 SYM and 2d YM.
The goal of this paper is to generalize them to intersecting loops. The strategy is simple: Using the conjecture above, we relate the intersecting 1/8-BPS loops in N = 4 SYM to intersecting Wilson loops in 2d YM on S 2 in the zero-instanton sector. We then solve the loop equation of 2d YM exactly at finite N . The loop equation in 2d YM was first solved for loops on R 2 at large N [5]. The result was generalized to loops on R 2 at finite N in [6] and to loops on S 2 at large N in [7][8][9]. While all these results take simple and compact forms, no such expressions were known for loops on S 2 at finite N : the only known expression in the literature involves a rather complicated sum over the Young diagrams [24]. We show that a simple closed-form expression does exist for loops on S 2 at finite N if the theory is restricted to the zero-instanton sector-the sector relevant for the BPS loops in N = 4 SYM.
The result of our computation is a coupled system of multiple integrals and a Gaussian matrix model. For instance, the expectation value of the figure eight loop with areasĀ 1 and A 2 reads (see figure 1) . (1. 3) The integration contours C 1,2 encircle the eigenvalues of the Gaussian matrix model (1.3) and they are placed far apart from each other (see sections 2 and 3 for more details). The integral can be evaluated explicitly at large N and it gives ,ā i :=Ā i − 2π 2 , gā := g 1 −ā 2 π 2 , ρā := π +ā π −ā , (1.4) and Iā k := I k (4πgā) is the modified Bessel function. The result at strong coupling reproduces the area of the minimal surface as we show in section 3.2.1.
As an important application, we compute a small angle limit of the cross anomalous dimension of intersecting Wilson lines. The cross anomalous dimension controls the mixing of two different configurations of Wilson lines, two intersecting lines and two touching lines (see figure 1), under the renormalization group. It also describes the soft anomalous dimension, which controls how the soft gluons transfer the color degrees of freedom of partons in a scattering process. We generalize the analysis of the Bremsstrahlung function in [25,26] and relate the small angle limit-and more generally the near-BPS limits-of these anomalous dimensions to our localization computations. The results are exact at finite λ and N and reproduce the answers at weak coupling in the literature [27][28][29].
The formalism in this paper will be used in our upcoming paper [30] on the defect CFT correlators on the higher-rank Wilson loops. There, we consider a configuration in which multiple fundamental Wilson loops with different areas are joined together by a projector to a higher-rank representation. We compute its expectation value by taking an appropriate linear combination of multiply intersecting loops.
The rest of the paper is organized as follows. In section 2, we review the 1/8 BPS Wilson loops and their relation to 2d YM. We then present a new integral representation for correlators of non-intersecting Wilson loops which simplifies the analysis of the loop equation in the subsequent sections. In section 3, we review the loop equation in 2d YM and solve it in the zero-instanton sector on S 2 . As a result, we obtain a closed-form integral representation for intersecting loops. We then demonstrate how to evaluate the integral at large N using examples of a figure eight loop and a two-intersection loop and study the weak-and strongcoupling limits. In section 4, we use these results to compute the small-angle and near-BPS limits of the cross anomalous dimension at finite λ and N . We then discuss future directions in section 5. A few appendices are provided to explain technical details.
Here all the indices run from 1 to 3 and x j 's are the coordinates of S 2 , 3 j=1 (x j ) 2 = 1. Throughout this paper, we consider the U(N ) gauge group unless otherwise stated.
It was conjectured through the perturbative computations [12,13] and later supported by the localization [11] that the computation of the 1/8 BPS Wilson loop reduces to that of the standard Wilson loop in 2d YM on S 2 in the zero-instanton sector, where the coupling constants of 2d YM (g 2d ) and N = 4 SYM (g YM ) are related by and the action of 2d YM is The expectation value of a single 1/8-BPS loop was computed by resumming the perturbative series in 2d YM and re-expressing it as a matrix integral 1 , By evaluating this matrix integral, we get where λ := g 2 YM N is the 't Hooft coupling constant, L 1 N −1 is an associated Laguerre polynomial and A is the area of the region encircled by the Wilson loop. The large N limit of this result gives a generalization of the famous result for the 1/2-BPS Wilson loop, with I n being the modified Bessel function.
As another application of the relation (2.2), a multi-matrix model for the correlation functions of 1/8-BPS Wilson loops was derived in [17]. A heuristic way 2 to derive the matrix model is as follows: First we rewrite the action of 2d YM as a deformed BF theory, If one integrates out B, one recovers the standard 2d YM action. Second we use the fact that the theory reduces essentially to the abelian theory and localizes to configurations for which B is piecewise constant with discontinuities across the Wilson loops. In such configurations, the first term of the action (2.8) becomes a boundary term where Σ is a subregion of S 2 whose boundaries are given by the Wilson loops and X := i ∂Σ A is a boundary holonomy. On the other hand the second term simply gives where A Σ is the area of the region Σ.
Performing such rewritings to each region Σ m , we obtain the following multi-matrix model action: Here s (m) j is an orientation factor which takes value +1 when the holonomyX j is oriented in the same direction as the boundary ∂Σ m , and −1 otherwise. From this action, the correlator of the Wilson loops can be computed as follows, where r k is the representation of the k-th Wilson loop W k . As was done in [17], one can integrate outB Σm 's and obtain an action purely in terms ofX j 's. This however is not convenient for analyzing the loop equation since the resulting action depends nonlinearly on the areas A Σm . In this paper, we instead integrate outX j and derive a new integral representation which is more convenient for the application of the loop equation.

A new integral representation for Wilson loop correlators
To derive an integral representation, it is convenient to first rescale the matrices as Then, the action and the correlator read with In the last equality, we used the notation for the 't Hooft coupling constant commonly used in the integrability literature, (2.16)

A single fundamental loop
Let us first discuss the expectation value of a single fundamental Wilson loop with no selfintersections. In the presence of a Wilson loop, S 2 is divided into two regions, one with the area 4π − A(=: A 0 ) and the other with the area A(=: A 1 ) (see figure 2), and the action of the matrix model reads To compute the expectation value of the Wilson loop in this matrix model, we use the Harish-Chandra-Itzykson-Zuber identity where A and B are diagonal matrices with entries a j 's and b j 's respectively, dΩ is an integral over the unitary matrices and ∆ is the Vandermonde determinant ∆(a) := i<j (a i − a j ).
Using (2.18), we can reduce the partition function of the matrix model (2.17) to integrals of eigenvalues, where the integrand I is given by . (2.20) We next expand the determinants into a sum over permutations: Owing to the symmetry of the rest of the integrand, all these permutations give the same answer . We thus pick the simplest one (σ j = j and σ j = j) and multiply a factor (N !) 2 . After integrating out x's, we get (2.22) In the second line, we integrated out b (1) 's and denoted the remaining variables b (0) j as b j . We also used A 0 + A 1 = 4π to simplify the exponent.
Similarly, the expectation value of the fundamental Wilson loop can be reduced to the following eigenvalue integral: with I defined in (2.20). To evaluate this integral, we again expand the determinants and integrate out x's. The only modification is that one of the delta functions gets shifted by −i owing to the factor e x k . We then get (2.24) Now the crucial observation is that the sum k can be recast into a contour integral, (2.25) We can then interpret this as an expectation value of a Gaussian matrix model and get (2.26) Here the integration contour C encircles all the eigenvalues b k 's and f A is given by The symbol • M denotes the expectation value of • in a Gaussian matrix model with the action S M := 8π 2 tr (M 2 ) /g 2 YM : The result (2.26) may appear more complicated than the expressions known in the literature [12,13]. However, for the analysis of the loop equation, (2.26) is more convenient since the area-dependence is simple. See section 3 for more details.

Multiple fundamental loops
We now generalize the result (2.26) to correlators of multiple Wilson loops.
To see how it works, let us first consider the correlator of two fundamental Wilson loops with the same orientation. As depicted in figure 3, we denote the areas inside outer and inner loops by A 1 and A 2 respectively and the area of the complement by A 0 . Then the matrix model action is given by withĀ 0 := A 0 ,Ā 1 := A 1 − A 2 andĀ 2 := A 2 . By reducing the matrix integral to the eigenvalues using the identity (2.18), we get with Expanding a product of determinants, we obtain a sum of permutations Again all these permutations turn out to give the same result owing to the symmetry of the integrand. We can therefore replace the integrand I 2 with Plugging this expression into (2.30), we find that the integrals of x (k) j give the delta functions j ). Performing the integrals of b (1) and b (2) , we get (2.34) Here we again denoted b (0) j by b j and used s A s = 4π. To compute the correlator of the Wilson loops, we insert to the partition function (2.30). We then obtain (1,2) j 's and recast the sums into integrals. The first term in (2.37) gives (2.39) Again (2.38) can be interpreted as an expectation value in the Gaussian matrix model (2.40) Similarly, the second term in (2.37) gives 3 Thus the sum of the two contributions (2.40) and (2.41) gives One can further simplify the expression by combining the two terms in (2.42): In the first term of (2.42), the integration contours of u 1,2 are on top of each other and encircle all the eigenvalues of the matrix model (see figure 4). If one deforms the contour of u 2 so that the two contours are far separated, the integral picks up a contribution from poles in the interaction term∆ (2.43) whose residues are The residue at u 2 = u 1 + i turns out to vanish since the product which is nonsingular inside the integration contour C. On the other hand, the residue at u 2 = u 1 − i precisely cancels the second term in (2.42). We therefore arrive at the following simple expression for the correlator of two fundamental loops: (2.46) The notation C 1 ≺ C 2 means that the contour C 1 is inside the contour C 2 and they are far separated from each other.
Carrying out the same analysis for correlators of more than two Wilson loops of the same orientations (see figure 5 for details of the setup), we obtain a simple generalization of (2.46), Interestingly this integral resembles multiparticle integrals [31][32][33][34][35][36][37][38][39][40][41] in the hexagon approach to the correlation functions and we will use this connection in our upcoming paper [30] to analyze the defect CFT correlators on the higher-rank Wilson loop.

Loop Equation and Intersecting Wilson Loops
We now discuss intersecting 1/8 BPS Wilson loops. In section 3.1, we review the derivation of the loop equation in 2d YM in [5]. We then apply it to the results in the previous section and derive integral representations for the intersecting BPS Wilson loops in section 3.2.
Let us make a remark before we proceed: Below we consider the Wilson loops without the normalization factor 1/N since they are often more convenient for analyzing the loop equation. To distinguish them from a more standard definition (2.2), we will put a tilde as Figure 6: (a) Deformation of the contour δC x denoted in red. We pick a point x close to the Wilson loop, draw a small circle and connect it to the original contour. (b) One-loop correction to Sym coming from the gluon exchange. Here we used the double-line notation to make manifest the contraction of indices.

Review of loop equation in 2d YM
We consider an infinitesimal deformation of the contour C → C + δC x . Specifically, we pick a point x close to but not on top of the Wilson loop. We then add a small circle around x and connect it the original contour as depicted in figure 6. This defines the area derivative Here δσ µν is an area tensor of δC x , which is a product of the area |δσ| and the orientation tensor n µν : The orientation tensor is unit-normalized and parameterizes the orientation of δC x , We then expand the path-ordered exponential to get To proceed, we use the Stokes theorem to the first term and write For the second term, we split it into the symmetric and the antisymmetric parts Here y 1 > y 2 means that the integral is path-ordered and y 2 is always behind y 1 . Since δC x is an infinitesimal contour around x, the antisymmetric part can be approximated as The sum of (3.6) and (3.8) gives On the other hand, the symmetric part at one loop reads (see figure 6) Here 1 is the identity matrix for the color factors and G µν (x 1 − x 2 ) is a propagator of gauge fields without color factors. In the axial gauge ( Plugging this into (3.9), we get A couple of comments are in order: First, the result (3.9) receives higher loop corrections. However since g 2d has mass dimension 1, it follows from dimensional analysis that the loop corrections come with higher powers of δσ µν and therefore can be neglected in the small area limit δσ µν → 0. Second, although (3.11) was computed in a specific gauge, the result is actually gauge-invariant. This is because the tree-level 5 gauge transformation only changes (3.9) by a term that integrates to zero: Third, the result (3.9) is for the U(N ) gauge group and the answer will be different for other gauge groups. It would be interesting to derive the loop equation for general gauge groups. Now, combining the symmetric and antisymmetric parts, we arrive at the following expression for the first-order variation: (3.14) 4 Note that our normalization for the coupling constant in 2d YM is different from [5]: (g [5] 2d ) 2 = (g ours 2d ) 2 /2. 5 As stated above, the higher-loop corrections can be neglected in the small area limit.
Let us make some clarification on the second term in (3.14): First, it contains the orientation tensor n µν and changes discontinuously across the contour. This property plays a crucial role in the derivation below. Second, it may appear to contradict a familar statement that the small deformation of the Wilson loop is simply given by the insertion of F µν . This apparent contradiction comes from the difference of the regularizations: Normally we insert F µν on top of the Wilson loop and regularize divergent diagrams which contract F µν and the loop by the principal value prescription. On the other hand, F µν in (3.14) is inserted slightly above or below the loop, which gives a different regularization. This alternative regularization produces a different answer, but one can show by the explicit computation 6 that the sum of the two terms in (3.14) reproduces the result in the principal value prescription. We thus have which is consistent with the familiar statement.
Now we derive the loop equation. First we differentiate (3.14) and get To evaluate ∂ µ n µν that appears in the second term, it is useful to go back to the definition Here we take x to be a point on the Wilson loop and µ to be a vector of length ε in the µ direction; ( µ ) ρ := εδ ρ µ . Since the points x ± µ /2 sit on different sides of the loop, we have n µν (x ± µ /2) = ±1. Thus the derivative ∂ µ n µν picks up a contribution from a delta function where the contour I x is an infinitesimal one-dimensional segment which passes through x. Second we use the fact that ∇ µ F µν is proportional to the equation of motion 7 We can then replace it inside the expectation value as (3.21) 6 In the holomorphic gauge Az = 0, the contraction of F µν and the straight-line loop reads (see e.g. [23]) In the principal value prescription, this integral vanishes while if we shift x slightly above or below by adding ±i , it gives ±g 2 2d N/(4π), which are precisely canceled by the second term in (3.14). 7 S 2d is the Euclidean action for 2d YM given by (2.4). Here we wrote the gauge indices a and b explicitly for convenience. As a result we get (see also figure 7) The integral (3.22) is nonzero even for nonintersecting Wilson loops since it receives a contribution from a "self contraction", namely a contribution from a trivial coincidence point x = x on the loop. It can be evaluated as (3.23) As can be seen from (3.19), this turns out to cancel the second term in (3.17). Thus, when the two terms in (3.17) are combined, we are left with contributions only from intersecting points and obtain the loop equation where the symbol − C means that we remove an infinitesimal segment around x when we perform the integral. For the actual application, it is useful to choose ν to be the direction of the Wilson loop near the intersection and µ to be the direction orthogonal to it. In addition, we integrate (3.24) along a small path in the µ direction. For the left hand side, this is basically equivalent to going back to the definition (3.18) and removing 1/ε, which can be pictorially represented as 8 We then getL where C 1,2 are the two contours obtained by reconnecting the loop C at the point x. Here we assumed that only two lines intersect at the point x. When more than two lines intersect, the right hand side of (3.26) is replaced by a sum of all possible reconnections.
Using the fact that the Wilson loops in 2d YM only depend on the area, we can convert (3.26) to area derivatives (see figure 8), This formula can be applied straightforwardly when the loop is on R 2 . For loops on S 2 , not all the areas (S i 's) are independent. In such a case, we go back to the pictorial definition of L x given in (3.25) and translate it into area derivatives as we see in the next subsection.

Intersecting BPS loops in N = 4 SYM
We now apply the loop equation (3.27) and compute the expectation value of intersecting 1/8-BPS loops in N = 4 SYM.

Figure eight loop
The simplest intersecting loop is the "figure eight" loop. Since the figure eight loop on S 2 is equivalent to a one-intersection loop (see figure 9), we start with the latter, and the result for the figure-eight can be simply obtained by relating the areas as in figure 9.
Under the action ofL x , the one-intersection loopW A 1 ,A 2 reconnects to a product of two disconnected loops W 1,2 with areas A 1,2 . Using (3.25), this can be translated into area derivatives as follows 9 : This can be solved by replacing the right hand side with the integrals (2.42). The result reads Alternatively, one can use (2.46) to write Translating this to the figure eight loop using the relation (see figure 9), To derive this, one simply needs to note that the first term (3.25) decreases A 2 while the second term increases A 1 when applied to the one-intersection loop. 10 For a k-wound fundamental loop, it is easy to see from the Gaussian matrix model (2.5) that Large N In the large N limit, one can evaluate the integral (3.31) explicitly. To do so, we first replace the expectation values of f A 's with their large N results, where f A is given by and G(u) is the large N resolvent We then get Here we replaced∆/(u 1 − u 2 ) with its large N limit, 1/(u 1 − u 2 ). This coincides with the expression obtained in [8] for the figure eight loop at large N (after the redefinition of the areas (3.32)). Below we evaluate this integral explicitly in terms of the Bessel functions.
To proceed, we introduce the Zhukowski variable 11 and express f A as Substituting this into (3.36), we get . (3.39) We then expand 1/(x 1 − x 2 )(1 + 1/x 1 x 2 ), close the contour C 1 at the origin, and then close the contour C 2 at infinity using the formula [26,43] ρ k a I k (4πg a ) = dx 2πix k+1 e 2πg(x+1/x) e 2ga(x−1/x) , (3.40) with As a result, we arrive at the following expression for the one-intersection loop, However, there is a subtlety in the strong-coupling expansion: The asymptotic expansion of the modified Bessel function reads and it gives a polynomial in k at each order. This leads to a divergence of the infinite sum when ρā 1 > ρā 2 (i.e.Ā 1 >Ā 2 ). A more reliable approach at strong coupling is to perform the saddle-point analysis to the integral (3.39) as we see below.
Weak coupling expansion Here we present the weak-coupling expansion (up to two loops) of the figure eight loop obtained from (3.44): As expected, the result is symmetric under the exchangeĀ 1 ↔Ā 2 and reduces to the one for a single Wilson loop whenĀ 2 = 0. The result for the one-intersection loop can be obtained by the replacementĀ 1 = 4π − A 1 andĀ 2 = A 2 .
Strong coupling expansion To perform the strong coupling expansion of the figure eight loop, we do the saddle point analysis to the integral (3.36) with the replacement (3.32). For definiteness, we assume 12Ā 1 <Ā 2 below. The saddle points are determined purely by the exponent of f A 's and we find two saddle points for each variable Owing to the geometrical constraintĀ 1 +Ā 2 ≤ 4π and the assumptionĀ 1 <Ā 2 , we have |ρā 1 | < |1/ρā 2 |. Thus, the saddle points can be reached by deforming the original contours C 1 ≺ C 2 without making the contours pass through 13 each other. Among those four saddle points, the ones with plus signs are always dominant. Expanding the integral around the saddle point using the coordinates x 1 = it 1 + ρā 1 and x 2 = it 2 + 1/ρā 2 , we get Performing similar analysis to other saddle points, we get a sum gā 1ā 2 + gā 2ā 1 − e 4π(−gā 1 +gā 2 ) − e 4π(gā 1 −gā 2 ) gā 1ā 2 − gā 2ā 1 .

(3.50)
To understand (3.50) from the string worldsheet, it is useful to recall the connected correlator of two Wilson loops W A 1 W A 2 analyzed in [18]. In that case, the BPS equation of the string worldsheet does not allow any connected surface and the only allowed surface is a degenerate cylinder made of two disks connected by a zero-area tube (see figure 10). Physically the zero-area tube corresponds to a propagator of a graviton and gives a factor of 1/N 2 . Combined with a factor g 2 which comes from integrating over the endpoints of the propagator, we obtain the correct g-dependence at strong coupling, (3.51) 12 Since the result (3.44) is symmetric under the exchange ofĀ 1,2 , the result for the other caseĀ 1 ≥Ā 2 can be obtained by a simple replacementĀ 1 ↔Ā 2 .
13 When the contours pass through each other, the saddle point approximation can potentially fail owing to the term 1/(x 1 − x 2 )(1 + 1/x 1 x 2 ) in the integrand. A similar argument seems to hold for (3.50): The idea is to consider a degenerate disk made of two disconnected worldsheet with disk topology (ending on the two loops of the figure-eight), joined together by a zero-area strip (see figure 10). The zero-area strip can be viewed as a propagator of an open string and gives a factor of 1/g. Combined with other factors, it gives which is the correct g-dependence in (3.50). The other three saddle points can be interpreted similarly as contributions from stable/unstable disk solutions [18,44] joined by a zero-area strip. It would be interesting to perform more detailed analysis and reproduce the full oneloop answer including the numerical coefficients. (See [45,46] for recent progress on the one-loop computation on the worldsheet for a single Wilson loop.)

Two-intersection loop
We now generalize the analysis to the two-intersection loopW A 1 ,A 2 ,A 3 depicted in figure 11. Applying the loop equation to the intersection shown in figure 12, we obtain is a one-intersection loop with areas A 1 and A 2 −A 3 . Solving this equation using (3.31) and normalizing it as 14 we get 14 The normalization is chosen such that the expectation value is O(1) in the large N limit.  where F (A 1 − A 3 , A 2 ) is the integration constant. To determine F , we consider the limit A 3 → 0, in which the two-intersection Wilson loop reduces to two disconnected Wilson loops with areas A 1 and A 2 respectively. We then get the constraint which allows us to compute F . As a result, we get

(3.56)
Large N Let us next discuss the large N expansion of (3.56). For later purposes, we compute it up to the first subleading order in 1/N .
First we consider the first term of (3.56) (to be denoted by W 1 ). It coincides with the integral representation of the two disconnected Wilson loops with areas A 1 − A 3 and A 2 (see (2.46)). Thus, using the result in the literature, we get the large N expansion with a ij := (A i − A j − 2π)/2. Here the first two terms come from the disconnected part, which is a product of the expectation values of a single loop 15 [15], while the last term is the connected part computed in [18].
The second line of (3.56) (to be denoted by W 2 ) can be evaluated in a manner similar to the figure eight loop: Namely we replace f A and∆/(u 1 − u 2 ) 2 with f A and 1/(u 1 − u 2 ) 2 respectively, rewrite it in terms of the Zhukowski variables and perform the integrals by closing the contour C 1 at the origin and C 2 at infinity. The result reads Adding the two contributions, we get the large N answer (ρ a 13 ρ a 2 ) k .

(3.59)
Weak coupling expansion From (3.59), the weak coupling expansion can be obtained straightforwardly by expanding each term in powers of g 2 . The result up to one loop reads One can check that the result reproduces the one for the two disconnected loops [15,18] in the limit A 3 → 0.
Strong coupling expansion We now evaluate the strong-coupling limit of (3.56). The first term, W 1 , is easy to evaluate since it coincides with the correlator of disconnected 15 The expectation value of a single loop up to O(1/N 2 ) is given by [15] W A = I 1 (4πg a ) 2πg a + π 2 g 2 a 3N 2 I 2 (4πg a ) + · · · . (3.58) Wilson loops with areas A 1 − A 3 and A 2 . Using the result in the literature [15,18], we get where we kept only the leading exponential. W 1,disc and W 1,conn are the contributions from the disconnected part and the connected part respectively. The second term in (3.56), W 2 , can be evaluated by the saddle point analysis. The result reads (g a 1 a 23 + g a 23 a 1 ) 2 − √ g a 13 g a 2 e 4π(ga 13 +ga 2 ) (g a 13 a 2 + g a 2 a 13 ) 2 . (3.62) Let us discuss the worldsheet interpretation. W 1,disc is simply a product of the contributions from two disconnected surfaces whose worldsheet interpretation is already discussed in [15]. The rest (W 1,disc and W 2 ) scales as W 1.disc , W 2 ∼ e 4πg /(gN 2 ), which coincides with (3.51). This shows that the relevant worldsheet configurations are again two disconnected surfaces connected by a zero-area tube. The only complication here is that there are two different disconnected surfaces ending on the Wilson loop, corresponding to two different exponentials in (3.62). The first one ends on the closed loops with areas A 1 and A 2 − A 3 while the second one ends on the closed loops with areas A 1 − A 3 and A 2 (see also figure 11). However owing to the geometrical constraint A 1 > A 2 > A 3 , we have g a 13 + g a 2 ≥ g a 1 + g a 23 .
Thus the leading strong coupling answer is always given by ∼ e 4π(ga 13 +ga 2 ) .

Cross Anomalous Dimension at Small Angle
We now apply the results in the previous sections to the computation of the cross anomalous dimension. The cross anomalous dimension is a quantity which governs the UV divergence associated to an intersection of two Wilson lines. In some respects, it is similar to the cusp anomalous dimension, which governs the UV divergence associated to a cusp of the Wilson line. However, one notable difference is that the cross anomalous dimension is a 2 × 2 matrix since the intersecting Wilson lines mix with the "touching" Wilson lines (depicted in figure  13-(a)) under the renormalization group flow.
Note that there is another cross anomalous dimension which governs the mixing of two different touching Wilson lines depicted in figure 13-(b). As we show in Appendix B, our formalism can be applied to this quantity as well.  Figure 13: (a) Intersecting lines (denoted by i) and "touching" lines (denoted by t). They are characterized by a geometrical angle φ and mix under the renormalization. In N = 4 SYM, we can consider the generalization of these lines which couple also to scalars. In that case, the black lines couple to n 1 · Φ and the red lines couple to n 2 · Φ. This introduces an additional angle cos θ = n 1 · n 2 . (b) Two touching configurations (1 and 2) whose cross anomalous dimension will be computed in Appendix B.

Cross anomalous dimension in
tersection [28,29]: Here µ is the RG scale and W i and W t denote the intersecting and the touching configurations respectively. The superscript R signifies the fact that the Wilson lines are renormalized and are related to the bare Wilson lines W A by the multiplicative renormalziation, with being the UV cut-off. As we discuss in more detail later, in conformal field theories one can compute the cross anomalous dimension more directly from the expectation value of the bare Wilson lines by reading off the coefficient of log , W ∼ e Γcross log .
The cross anomalous dimension is a function of the angle φ between the two intersecting lines. When the angle is analytically continued as φ → iϕ, it gives the so-called soft anomalous dimension. The soft anomalous dimension controls the IR divergence of the scattering amplitude of two massive quarks in the Regge kinematics, s, m 2 −t Λ 2 QCD , and describes how the soft gluons transfer the color degrees of freedom of the quarks. In that context, ϕ is the boost angle between the two quarks defined by with p 1,2 being the four-momenta of the quarks. See [28,29] for more details.
For the application to QCD, the limit ϕ → ∞ is of particular interest since it describes the high energy scattering of light partons. On the Wilson line side, this corresponds to an intersection of two light-like lines, see for instance [47][48][49]. The limit was studied also in N = 4 SYM: In [50], the self-crossing lightlike loop was analyzed up to nine loops. The analysis was pushed further in [51] in which the anomalous dimensions relevant for the limit were determined exactly in the large N limit using the pentagon OPE decomposition [52]. In this paper, we focus on different limits, namely φ ∼ ϕ ∼ 0 and the near BPS limit, and compute the anomalous dimension exactly at finite N .
Generalization in N = 4 SYM In N = 4 SYM, one can consider a generalization of the cross anomalous dimension Γ cross (φ, θ) which depends on another angle θ. To do so, we consider an intersection of the supersymmetric Wilson lines, where the R-symmetry polarization n is a six-component unit vector which dictates the coupling to the scalars Φ = (Φ 1 , . . . , Φ 6 ). The intersection of such lines can be characterized by the geometric angle φ and the R-symmetry angle θ which is defined by where n 1,2 are the R-symmetry polarizations of the two lines at the intersection point. In this phase, we have two kinds of W -bosons, one coming from the (1, k) (or (k, 1)) component and the other coming from the (2, k) (or (k, 2)) component of the gauge field. In the limit m 1,2 1, they can be treated as classical probe particles and their coupling to the unbroken U(N ) degrees of freedom is given precisely by (4.4). Thus Γ cross (φ, θ) gives a natural generalization of the cross anomalous dimension in QCD and it characterizes the IR divergence of the scattering of massive W -bosons after the analytic continuation φ → iϕ.
When φ = θ, the whole configuration becomes BPS and the anomalous dimension Γ cross (φ, θ) vanishes. Expanding Γ cross (φ, θ) around this limit, we obtain 16 In what follows, we compute the first coefficient γ cross (θ) exactly as a function of λ and N .

Two-point function of intersections
Let us explain in more detail how to extract the cross anomalous dimension from the expectation values of the bare Wilson lines. The key idea is to regard them as the two-point functions of intersections. To be concrete, consider the intersection of the following two Wilson lines on R 4 ,  In addition to the obvious intersection at the origin, the lines intersect also at infinity. This is easier to see if one maps the configuration to S 2 using the conformal transformation (4.10) Here X's are the embedding coordinates of S 2 while x's are the coordinates on R 2 inside R 4 . The two intersections are mapped to the north and the south poles of S 2 , see figure 14.
To extract the cross anomalous dimension, one has to consider the operator mixing: For each intersection, we can either let the lines intersect (to be denoted by i) or resolve the intersection and make the lines touching (to be denoted by t). Since there are two intersections, we have in total four different configurations of the Wilson lines which we denote by W ii , W it , W ti and W tt . They can be regarded as two-point functions of "operators", labeled by i and t, sitting at the intersections.
These four choices of two-point functions can be naturally organized into a 2 × 2 matrix Here i, t| and |i, t denote the intersections at the origin and infinity respectively and D is the dilatation operator. UV and r IR are the UV and the IR cutoffs. W can be expressed alternatively in terms of the cross anomalous dimension Γ cross and the overlap η as W = e Γcross log( UV /r IR ) · η , (4.12) where Γ cross and η are 2 × 2 matrices defined by i| , t| Γ cross := i|D , t|D , η := i|i i|t t|i t|t . (4.13) In the near BPS limit θ ∼ φ, both η and Γ cross can be expanded in powers of (φ − θ), (4.14) We then obtain the following expansion of the two-point function matrix W:

(4.15)
Thus γ cross is given by the coefficient of the logarithm in W 1 , multiplied by (W 0 ) −1 , where W 0 is nothing but the expectation value in the BPS limit: (4.17)

Cross anomalous dimension from localization
We now relate γ cross given in (4.16) to the localization computation. This can be done by following the arguments in [23,25,26], which we briefly review below 17 .
The first step is to start with the BPS limit θ = φ of (4.8) and deform θ slightly. This amounts to inserting a scalar on the second Wilson line with Φ (τ ) = − sin θΦ 1 + cos θΦ 2 | x µ =(τ cos θ,τ sin θ,0,0,0,0) . Now, using the invariance under the dilatation around the origin, we can determine the position dependence of the scalar insertion as where we denoted all the other parts (the Wilson lines W 1,2 and possible resolutions of the intersections etc.) by · · · . From this, one can see that the integral of τ produces the logarithmic divergence and its coefficient is given by the expectation value of the Wilson loops with the Φ insertion. More explicitly, we have the relation where W AB is the Wilson line W AB with Φ inserted at (cos θ, sin θ, 0, 0, 0, 0). The factor of 2 comes from summing up contributions from ∞ 0 dτ and 0 −∞ dτ . The second step is to map the whole configuration to S 2 using (4.10). We then get two great circles whose contour are given by with t ∈ (−π, π). They couple to Φ 1 and cos θΦ 1 + sin θΦ 2 respectively and satisfy the 1/8 BPS condition (2.1). Under this map, the insertion Φ (τ = 1) is mapped to the insertion at a point t e where W 2 intersects with the equator of S 2 (see figure 14).
The third step is to replace Φ at t e with Φ − iΦ 4 =: −Φ. This replacement does not affect the expectation values since the correlator with Φ 4 vanishes owing to the charge conservation. We then use the fact that the insertion ofΦ corresponds to the insertion of the field strength in 2d YM [11,16]. Therefore its expectation value can be computed by taking the area derivative. In our setup there are two regions with area 2(π − θ), and changing θ by δθ leads to a total change (decrease) of the area by 4δθ. We thus have the relation with W 0 is the expectation value in the BPS limit (4.17). Combined with (4.21), it gives Using (4.16), we get the formula relating γ cross to the BPS Wilson loops Before we proceed, let us comment on the normalization. In what follows, we normalize W 0 by dividing the path-ordered exponentials by a common factor N 2 . This is a natural normalization for discussing the renormalization of open intersecting lines (see [27][28][29]) and makes Γ cross more symmetric. However it does not coincide with the normalization used in some of the literature in which they discuss the renormalization of closed loops. We will later translate the final result to such a normalization.
The last step is to compute W 0 . Let us first consider W tt . Since both intersections are resolved, it coincides with a correlator of two disconnected loops. Thus, setting A 1 = 2(π +θ) and A 2 = 2(π − θ) in (2.46), we get (4.26) Using the result in the literature, we can compute its large N limit as 18 (4.27) Second, W ti (= W it ) is a single-intersection loop with areas A 1 = 2(π + θ) and A 2 = 2(π − θ) (or equivalently, a figure-eight loop with areasĀ 1 =Ā 2 = 2(π − θ)). We thus get from (3.31) Here the extra factor 1/N comes from the normalization that we adopted. The large N limit can be computed from (3.42) as Finally, W ii is a two-intersection loop which we computed in (3.56). Setting A 1 = 2(π + θ), A 2 = 2π and A 3 = 2θ, we get (4.30) The large N limit can be computed from (3.59), and the result reads Note that these expectation values satisfy the following relation (at finite N ) They are simply the loop equations written in terms of θ-derivatives, but one can also verify them directly from (4.26), (4.28) and (4.30).
These expressions, together with the relation (4.25), give the exact cross anomalous dimension in the near BPS limit of the U(N ) theory. Using the relations (4.32), we get Note in particular that the entries in the first row are 0 and −4πg 2 /N at all orders in λ and N . Epanding the result at large N , we get with

(4.35)
Here h 0 is related to the Bremsstrahlung function B(λ) in [25] by The eigenvalues (γ ± ) of γ cross are given (up to O(1/N 2 )) by Since the U(1) theory is free, the computation is rather straightforward. In addition, as the gauge group is Abelian, the path-ordering is unnecessary and all the four entries of W become identical. The result can be read off from the Appendix A of [17], or from the two-matrix model in [18] specialized to U (1). For the two disconnected loops with areas defined as in figure 3, this gives (4.39) Setting A 1 = 2(π + θ) and A 2 = 2(π − θ), we get We then obtain where γ cross is the result for the U(N ) theory (4.25) and 1 is the 2 × 2 identity matrix.
Weak coupling Expanding (4.25) and (4.41) at weak coupling, we get the following result for γ cross up to three loops: (4.42)

Closed-loop normalization
. 19 In [27], Γ cross and Γ closed cross were denoted as Γ cross and Γ cross respectively. γ closed (1,2) are in perfect agreement 20 with the near-BPS limit of the two-loop results in [27]. It would be interesting to perform a direct three-loop computation and reproduce γ closed (3) . Strong coupling The strong coupling limit of γ cross can be computed from the results for intersecting loops in section 3.2. The result in the planar limit reads where c 0 and c 1 are θ-independent while c 2 -c 4 are θ-dependent prefactors, among which c 3 and c 4 are relevant for the analysis: We then get the following result in the planar limit: The leading strong-coupling answer for the lower-diagonal component −4π∂ θ g θ reproduces the prediction made in [27] using the classical worldsheet. By contrast, the leading strongcoupling answer for the upper-right component is given by and does not match with the one in [27], which predicts the same answer as the lowerdiagonal component, 4π∂ θ g θ . This however does not immediately imply contradiction: As is well-known, the individual matrix elements of the anomalous dimension depend on the choice of the basis of operators. It is likely that the choice we made here for supersymmetric localization is different from the choice implicitly made in the analysis of the classical worldsheet. To avoid such ambiguities, we should compare the eigenvalues of the anomalous dimension matrix, which in fact agree with the ones in [27]. It would be important to understand this point further and also perform a comparison at the nonplanar level.

Conclusion
In this paper, we computed the expectation values of intersecting 1/8 BPS Wilson loops in N = 4 SYM at finite λ and N using supersymmetric localization and the loop equation. The results are given by a coupled system of the Gaussian matrix model and multiple contour integrals, which in the planar limit give an infinite sum of products of modified Bessel functions. Applying the formalism to near-BPS limits of the cross anomalous dimension, we reproduced the perturbative data in [27].
The main message of this paper is that the loop equation provides a powerful computational tool in N = 4 SYM when combined with localization. 21 It would be interesting to explore the connection with other nonperturbative techniques such as integrability and the conformal bootstrap 22 : The intersecting lightlike Wilson lines were studied in the planar limit [51] from integrability by the pentagon OPE [52]. However the relation to the loop equation was not explored. Studying them through the lens of the loop equation may lead to stronger results, or at least would lead to a deeper understanding of the pentagon OPE.
The 1/2 BPS Wilson loop in N = 4 SYM is an example of a conformal defect [58][59][60]. The insertion of F µν -which plays the central role in the derivation of the loop equation-is a displacement operator, which is present in any conformal defect. Rephrasing the loop equation in the language of the defect CFT may allow us to use it as a dynamical input 23 for the conformal bootstrap. This would be perhaps useful for the nonsupersymmetric Wilson line discussed in [61][62][63]. Of course, in the absence of supersymmetric localization, one would need to study in this case the loop equation in the 4d gauge theory. Another direction is to analyze intersecting conformal defects in general CFTs. For the case of two intersecting 1d defects, one should be able to interpret them as a conformal two-point function of intersections as discussed in section 4.2.
Regarding the cross anomalous dimension, the simplest next step would be to generalize our computation to multiple lines intersecting at a point. This would shed light on the structure of the soft anomalous dimension of multileg amplitudes, studied for instance in [64]. Of course, our analysis only applies to a small angle (or near-BPS) limit but it might be possible to combine it with the bootstrap approach in [65] and constrain the full answer.
Another interesting direction is to perform the computation in different setups. For instance, the exact Bremsstrahlung function in N = 2 SCFT was studied in [66][67][68]. Generalizing it to the small angle limit of the cross anomalous dimension is an important problem. It would also be interesting to study the ladder limit of the Wilson loop in N = 4 SYM, in which the R-symmetry angle θ is sent to i∞ while the combinationλ := λe −iθ is held fixed. This limit selects the ladder diagrams which can be resummed analytically [69][70][71][72][73][74]. Last but not least, it is important to further study the cross anomalous dimension in perturbation theory. In particular, it would be desirable to generalize the result for the supersymmetric Wilson lines in [27] to nonsupersymmetric Wilson lines. and Ian Moult for discussions on the cross anomalous dimension. We thank CERN for hospitality during completion of this work. The work of SG is supported in part by the US NSF under Grants No. PHY-1620542 and PHY-1914860. The work of SK is supported by DOE grant number DE-SC0009988.

A Infinite Sum of Modified Bessel Functions
In this appendix, we derive identities for the infinite sum of modified Bessel functions and apply it to (3.44) to rewrite it in a more symmetric form. The starting point is the integral representation We now factorize the exponential into two pieces with z = z 2 1 + z 2 2 + β + 1 β z 1 z 2 , α = z 1 + z 2 β z 1 + z 2 β −1 .

B Cross Anomalous Dimension of Two Touching Lines
In this appendix, we compute the cross anomalous dimension of two touching Wilson lines. The basic strategy is the same as in section 4: We map it to a sphere, view it as two-point functions of intersections and differentiate it with respect to the angle θ.
In this case, the analogue of (4.17) is given bȳ where W ij denotes a Wilson loop whose intersections at the north and the south poles are resolved into configurations i and j in figure 13-(b). The formula (4.25) applies also to this case and the computation boils down to computing the BPS loops W ij | φ=θ . The main difference from section 4 is that all the relevant loops are non-intersecting and one can simply use the results in the literature.