Quantum Wronskian approach to six-point gluon scattering amplitudes at strong coupling

We study the six-point gluon scattering amplitudes in N=4 super Yang-Mills theory at strong coupling based on the twisted Z_4-symmetric integrable model. The lattice regularization allows us to derive the associated thermodynamic Bethe ansatz (TBA) equations as well as the functional relations among the Q-/T-/Y-functions. The quantum Wronskian relation for the Q-/T-functions plays an important role in determining a series of the expansion coefficients of the T-/Y-functions around the UV limit, including the dependence on the twist parameter. Studying the CFT limit of the TBA equations, we derive the leading analytic expansion of the remainder function for the general kinematics around the limit where the dual Wilson loops become regular-polygonal. We also compare the rescaled remainder functions at strong coupling with those at two, three and four loops, and find that they are close to each other along the trajectories parameterized by the scale parameter of the integrable model.


Introduction
The gluon scattering amplitudes in N = 4 super Yang-Mills theory are a subject of great interest in recent years. In the planar limit, they are dual to the nullpolygonal Wilson loops whose segments are light-like and proportional to external gluon momenta [1][2][3][4][5]. The duality implies a conformal symmetry in the dual space [1,[6][7][8]. This dual conformal symmetry strongly constrains the form of the amplitudes.
In particular, the maximal helicity violating (MHV) amplitude is expressed as a sum of the Bern, Dixon and Smirnov (BDS) formula [9] and a finite remainder (remainder function), which is a function of the cross-ratios of the cusp coordinates for the null polygon.
The amplitudes have been studied intensively from both weak-and strong-coupling sides. At weak coupling, recent developments using the mixed motive theory have made it possible to evaluate the remainder function for the six-point amplitudes up to four-loop level [10]. Moreover a method based on the OPE and integrability has been proposed to calculate the scattering amplitudes in this theory, which is expected to be applicable to the intermediate coupling region [11][12][13][14].
At strong coupling, the AdS/CFT correspondence asserts explicit relations between N = 4 super Yang-Mills theory and the superstring theory in AdS 5 ×S 5 . Alday and Maldacena have thereby proposed that the MHV amplitude can be evaluated by the area of the minimal surfaces in AdS with a null-polygonal boundary along the Wilson loop [1]. It turns out later that the remainder function for the nullpolygonal minimal surfaces is calculated with the help of integrability [15]. Namely it is obtained by solving the Y-system or the thermodynamic Bethe ansatz (TBA) system related to certain two-dimensional quantum integrable systems [16][17][18]. The cross-ratios are given by the Y-functions at special values of the spectral parameter and the remainder function is expressed by the free energy of the TBA system and the Y-functions.
For the null-polygonal minimal surfaces in AdS 3 and AdS 4 space-time, the relevant integrable systems are the homogeneous sine-Gordon models [19] with purely imaginary resonance parameters [18], which are the perturbed SU(N) k /U(1) N −1 coset conformal field theory (CFT) at level k = 2 and 4, respectively. Around the limit where the null boundary becomes regular polygonal, corresponding to the UV limit of the two-dimensional systems, the remainder functions are calculated analytically for lower point amplitudes [20][21][22]. There, the free-energy part is evaluated by the standard bulk conformal perturbation theory (CPT). In order to evaluate the Y-functions, the g-function or the boundary entropy is utilized because the Y-function itself is not well incorporated in quantum field theory. The boundary CPT then efficiently yields the analytic expansions of the Y-functions. The resultant remainder functions are observed to be close to the two-loop results after an appropriate normalization/rescaling. Numerically, one also finds that this similarity extends beyond the UV limit. The minimal surfaces in these cases, however, give the amplitudes with some specific kinematic configurations of gluon momenta.
The null-polygonal minimal surface in AdS 5 with six cusps is the simplest nontrivial example that allows the most general kinematic configuration. At strong coupling, the relevant two-dimensional system is the Z 4 -symmetric integrable model [23][24][25] with a boundary twist [26]. The remainder function around the UV limit in this case has been studied in detail in [27]. Although the free-energy part is analytically evaluated by the bulk CPT for the twisted Z 4 -parafermion, the expansion of the Y-functions there is determined by numerical fitting. The difference from the the AdS 3 and AdS 4 cases come from the fact that the TBA equations in the AdS 5 case have a twist parameter, and it is unclear how to construct the g-function with this twist parameter.
Given the analytic results at weak coupling as well as the OPE method for finite coupling, the analytic data at strong coupling would provide pieces of the whole picture of the scattering amplitudes. They would also be useful for a check of the finite-coupling analysis. In this report, we thus decide to devote ourselves to the analytic expansions of the remainder function for the general kinematics.
In order to overcome the problem mentioned above, we take below another route for the UV expansion of the Y-functions, which does not rely on the g-function. Instead, our analysis is based on a seminal work by Bazhanov, Lukyanov and Zamolodchikov [28,29], where the role of quantum monodromy matrix is clarified in the minimal CFT, M 2,2n+3 , perturbed by the Φ 1,3 operator. Most remarkably, a new object in field theories, Baxter's Q operator, is introduced in this work. They noted the importance of the fundamental relations among the T-and Q-functions, the quantum Wronskian relation [29]. The T-and Y-systems can be regarded as colloraries of this. The Q operators and the quantum Wronskian relation have also played important roles in the non-equilibrium current problem [30], in the ODE/IM correspondence [31], in the spectral problem of the AdS 5 /CFT 4 correspondence [32] and so on.
In this report, we provide a yet another application: the quantum Wronskian relation is very efficient in obtaining the analytic expansion of the Y-functions particularly in the CFT limit. Specifically, we apply it, for the first time, to the Z 4symmetric integrable model or its lattice regularization. The lattice regularization adopted here allows one to elucidate the analyticity of the T-/Y-functions numerically. The expansion of the Y-functions around the UV limit is then determined analytically up to and including the terms of order (mass) 4 3 . A series of the higherorder coefficients is also determined recursively.
Combined with the free-energy part, the UV expansion of the Y-functions gives the analytic expansion of the six-point reminder function for the general kinematic configuration. We also compare the strong-coupling results with the perturbative ones, and find that the rescaled remainder functions are close to each other for large ranges of the parameters. This is in accord with the previous observations in the AdS 3 and AdS 4 cases [20][21][22]33] as well as in the perturbative cases [10,34]. This paper is organized as follows: In Section 2, we review the TBA-system for the six-point gluon scattering amplitudes at strong coupling and express the remainder functions using Y-functions. In Section 3, we reconsider the Y-/T-functions based on the lattice model. Taking the scaling limit, we derive the TBA system for the amplitudes. In Section 4, we study the CFT limit of the TBA system. Based on the T-Q relation and the quantum Wronskian relation, we calculate the analytic expansion of the Y-function in the CFT limit. The detailed analysis of the asymptotics of the related spectral determinant based on non-linear integral equations is studied in Appendix A. In Section 5, we apply the analytic expansion of the Y-functions to determine the leading expansion of the remainder function for the six-point amplitudes and compare it with the perturbative calculations.
2. Y-system and TBA for scattering amplitudes at strong coupling Let us begin with a review on the evaluation of the six-point MHV amplitudes at strong coupling using TBA of the twisted Z 4 -symmetric integrable model.

Hitchin system and Stokes data
Alday and Maldacena proposed a method of computing the gluon scattering amplitudes in N = 4 super Yang-Mills using the AdS/CFT correspondence [1,15].
Consider the the scalar part of the n-point gluon MHV scattering amplitudes in the strong coupling limit, and factor out the contribution of the tree amplitudes. According to [1], the result can be evaluated by computing the area A of the corresponding classical open string solutions in AdS 5 spacetime: 1 where λ denotes the 't Hooft coupling. A string solution represents a minimal surface whose boundary is a polygon located on the boundary of AdS 5 . The polygon consists of n null edges given by the n momenta of incoming gluons. Note that the amplitudes are defined in space-time with signature (3,1), (2,2) or (1,3).
The equations of motion of the string under the Virasoro constraints are rephrased as the SU(4) Hitchin equations for the connections (A z , Az) and the adjoint scalar fields (Φ z , Φz) with the Z 4 automorphism [16]. They are equivalently presented by the linear equations for a four component vector q(z,z; ζ), with appropriate boundary conditions. Here D z and Dz denote the covariant derivatives and ζ stands for the spectral parameter. The explicit forms of D z and Dz are given in [16]. The information of the null polygon is encoded in the asymptotic behavior of Φ z , which is diagonalized at infinity by an appropriate gauge transformation: The polynomial degree of P (z) is n − 4 and its coefficients parameterize the shape of the polygon. When n > 4, the linear equation necessarily possesses irregular singularity at infinity, which implies the Stokes phenomena. Customarily, the whole complex plane is divided into sectors, We denote by s k (z,z; ζ), the most recessive solution as |z| → ∞ in W k . One can consistently choose (s k , s k+1 , s k+2 , s k+3 ) as a linearly independent basis in W k . This implies a linear dependent relation among five neighboring s j 's, where b k+1 and c k+1 are some constants. Note the periodicity b i+3 = b i .
The normalization of s j (z,z; ζ) is fixed such that s j , s j+1 , s j+2 , s j+3 = 1, (2.5) where s i , s j , s k , s l ≡ det(s i s j s k s l ). The Z 4 automorphism results in the following relations for the Stokes data, Below, we confine our argument to the n = 6 case where For W j+6 = W j , we shall impose the boundary condition, The multiplier µ is parametrized as where φ is real for solutions in the (1,3) or in the (3, 1) signature of the fourdimensional space-time whereas it is purely imaginary for the (2, 2) signature. It also appears, e.g., in the relation among b i , (2.11)

Y-functions, TBA and evaluation of area
The key ingredients in the following discussion are the Y-functions, which are defined explicitly by 2 Y 2 (θ) = s 1 , s 2 , s 3 , s 5 s 2 , s 4 , s 5 , s 6 (e θ+πi/4 ), (2.13) The Y-functions here are identified with those in [17,22] as The cusp coordinates which appear below are also related as Then it was shown in [16] that the Hirota bilinear identities (or Plücker relations), as well as the relations (2.9), (2.6) and (2.7), lead to the following Y-system, The asymptotic behavior of the Y-functions is shown to be where Z is a complex parameter with phase ϕ. This is related to the moduli param- We further introduce The asymptotic behavior (2.17), together with the assumption of the analyticity of ǫ,ǫ in the strip Im θ ∈ (− π 4 , π 4 ), leads to the integral equations, 22) and the symbol * denotes the convolution, . These turn out be identical to the TBA equations for the Z 4 -symmetric integrable model [23][24][25]  In the original setting, the geometric data such as cross-ratios are given first, and then the area of the surfaces should be evaluated. Below, we slightly deform this logic: the TBA equations are given first, then the cross-ratios and the area are evaluated second. Once the TBA equations are solved, the coefficient b k in eq. (2.4) and the cross-ratios of gluon momenta are given by The cusp coordinates x j are related to the external momenta through In the literature, u k := 1/U k are often used as a basis of independent cross-ratios for the six-point case. The number of the independent cross-ratios matches that of the parameters in the TBA system (|Z|, ϕ, µ).
Naively, the values of Y j outside the analytic strip are necessary in order to evalu- Alternatively, one may also use the Y-system (2.15), (2.16) as recurrence relations to generate Y j (kπi/4) for any k ∈ Z from a set of Y j (k ′ πi/4) in the analytic strip. As discussed shortly, the Y-functions have the periodicity Y j (θ + 3πi/2) = Y j (θ), and thus the procedure terminates after a few steps.
We are now in position to write down the area A of the 6-cusp solutions or the scalar magnitude of the gluon scattering amplitudes in the strong coupling limit.
Instead of A itself, we deal with the finite remainder defined by where A BDS is the all-order ansatz for the MHV amplitude proposed by Bern, Dixon and Smirnov [9], including the divergent part. The present formulation then yields and each part reads The minus of the last term F = −A free coincides with the free energy whereas 2|Z| is identified with the mass/scale parameter of the Z 4 -symmetric integrable model.
The overall coupling dependence √ λ has been omitted above.
Within the framework described above, the numerical solutions of (2.20) and (2.21) yield explicit evaluation of the gluon scattering amplitudes [27] . Although the analytic solution to the TBA equations for generic Z and µ is beyond our reach, two limiting cases are accessible [36]. One is the limit where |Z| → ∞. In this case, the integrable model reduces to a free massive theory, and the free-energy part A free and the Y-functions Y j (θ) are expanded by multiple integrals. Via analytic continuation, this limit is also relevant for the amplitudes in the Regge limit [37,38].
Another limit is |Z| → 0, which we are interested in here.
When |Z| is strictly zero, A free and Y j are obtained as the central charge of the Z 4parafermion theory and a solution to the constant Y-system, respectively. Moreover, for small |Z| the free-energy part is expanded by the bulk conformal perturbation theory (CPT). The Y-functions are expanded by the boundary CPT through the relation to the g-function for µ = 1 [22], corresponding to the minimal surfaces in AdS 4 . For generic µ, however, it is still unclear how to incorporate µ in the framework of the boundary CPT.
In the following, we take an approach to the problem, which is different from any of the above, and is based on the integrable field theoretical structure proposed in [29]. This allows us to analytically evaluate the Y-functions for small |Z|, as shown in section 4. .

T-functions from the lattice regularization and their scaling limit
In this section we embed the Y-system into another tractable object in integrable systems, the T-system. This enable us to apply the machinery of the latter to evaluate the Y-functions for small |Z| in the next section.
There are several ways to introduce the T-system which is equivalent to the Ysystem in (2.15) and (2.16). Here we start with the lattice regularization [39]. One advantage of this choice is that analyticity assumptions, necessary to derive the TBA equations, can be checked numerically.
As is well known, the Z 4 parafermion model is related to the spin-1 2 XXZ model. Its Hamiltonian and spectrum can be studied from the transfer matrix. In order to define the latter, we introduce R(v), the U q ( sl 2 ) R matrix of spin 1 2 representation: . Then one can construct an inhomogeneous transfer matrix T 1 (x) acting on 2N sites by, where suffix 0 denotes the auxiliary space. The spectral parameter x is set via v = ix for later convenience. We have also introduced the diagonal twist matrix ] under the trace. The quantity µ = e i 3 2 φ is to be identified with the multiplier in (2.10).
To diagonalize T 1 , Baxter [40] ingeniously introduced an operator Q which commutes with T 1 (x). They satisfy Baxter's TQ relation, Below we shall consider T 1 and Q on their common eigenspace, thus we do not distinguish operators from their eigenvalues. The eigenvalue of Q is explicitly written with a set of Bethe roots {x j (φ)} with twist φ, This, together with (3.2), parameterizes the eigenvalue of the transfer matrix.
The fusion method generates a series of vertex models such that the auxiliary with a suitable normalization. 3 By construction {T j (x)} constitute a commutative family and they satisfy the T-system of SU(2) type, Note the periodicity of T j (x) in the present normalization is There is an additional relation when q is at a root of unity, which is crucial in obtaining a closed set of functional relations. From now on we fix Then the desired relation is The equation for j = 2 in (3.3) can be thus rewritten as thereby yielding a closed functional relations among T 1 and T 2 .
Below we will show that (3.3) for j = 1 and (3.7) can be transformed into TBA equations. Before doing this, we elucidate the analytic properties of T j deduced from numerics, as a merit in the lattice regularization. By definition, T j has 2N zeros and has no poles in complex x plane. Led by numerical observations we conjecture that all zeros of T 1 (x) are on the Im x = π 2 line while those of T 2 (x) are on the real axis in the ground state. Below, quantities which have no zeros and poles in the strip including the real axis will play an important role. We thus define Then the closed functional relations now read By adopting the change of variables, it is easily checked that Y 1 and Y 2 satisfy the same Y-system as (2.15) . This motivates us to define the pole-free functions, One then derives the functional equations, where both sides do not have any zeros or poles in the strip. Thanks to the analyticity, one arrives at Now consider the following scaling limit, In this limit, the driving terms become We denote the Y-functions in the scaling limit by Y sc j (x). If changing the variables as θ = 3 2 x, we recover (2.20) and (2.21) by the identification, We also define the T-functions in the scaling limit by where especially T sc 0 (x) = 1. In the scaling limit, the relations between the Yfunctions and the T-functions are drastically simplified, For later use, we shall also discuss the scaling limit of Q. The Bethe ansatz roots are roughly classified into two clusters, x ℓ j ∼ −Λ and x r j ∼ Λ. We thus adopt parameterizations,x Then the scaling limit of Q sc reads where we assumed the numbers of roots in the left and the right clusters are identically equal to N 2 . The prefactor stands for This Q-function plays an important role for studying analytical properties of the Y-functions.

Expansions of T-/Y-functions around CFT limit
In this section, we consider the expansions of the T-and Y-functions around ℓ = 2|Z| = 0. As mentioned in the previous section, such expansions for the minimal surface in AdS 3 or AdS 4 were studied in detail in [20][21][22] through the bulk and boundary conformal perturbation theory. Also, in [27], the expansions of the Yfunctions for the six-point case in AdS 5 were considered, but there remained an unknown function Y

General argument
Let us start with a general argument on the expansions of the T-and Y-functions.
Solutions to the Y-system have a periodicity as conjectured first in [41], and it played an important role in the analysis of the perturbed CFT (see, e.g. [36]). Here we have an apparent periodicity 4 inherited from the underlying lattice model (see (3.4)), and it motivates the expansions of Y j (θ) and T j (θ), ) as a function of θ. The conformal perturbative argument suggests that the coefficients are also expanded around ℓ = 0, 5 The central issue here is to determine these coefficients.

Quantum Wronskian and CFT limit
We remind that Baxter's TQ-relation (3.2) is a second order difference equation, and there are two linearly independent solutions. Let {x j (−φ)} be a set of BAE roots for negative twist −φ, then the second solution readsQ(x) = e − x γ φ j 2 sinh(x−x j (−φ)). Its scaling limit is obtained from (3.13) by φ → −φ.
There exist remarkable relations (the quantum Wronskian relations) between such two independent solutions of Baxter's TQ relation and the fusion transfer matrices [29], The same relation naturally arises in the context of the massive generalization of the ODE/IM correspondence [46]. A little bit different view from the lattice model is remarked in [47].
Below we will show that coefficients Y (p,0) j (j = 1, 2) are determined by applying the above relations. We start with the scaling limit of (4.3), Note q = e 2π 3 i . We further take the CFT limit ℓ → 0 with the shift λ → ℓ 2 − 2 3 λ. Let us define scaled functions, Then the quantum Wronskian relation takes the form, where we used C(φ)C(−φ) = 1, as {x j (φ)} = {−x j (−φ)} resulting form the Bethe ansatz equations. The same limit of the Y-functions is denoted by Y CFT j (λ 2 ) = lim ℓ→0 Y j (θ − log ℓ 2 ), which has an obvious expansion from (4.1) and (4.2), In the literature, the limit, ℓ → 0 without rescaling in the spectral parameter is also referred to as the CFT limit. In this paper, we call it the UV limit in order to avoid confusion. In the CFT limit, the mass terms in the TBA equations become chiral, and we are left with the massless TBA system.
Thus once the "moments "a n andā n are known, Y Since p n (−a) is of the form −a n + · · · with the ellipses being a polynomial in a 1 , · · · , a n−1 , the above relation reduces to where γ is defined in (3.5). The explicit forms of the first few r n (φ) read + a 2ā1 cos 1 6 (π − 3φ) − a 1ā2 cos 1 6 (3φ + π) .
In general r n is a polynomial of a j andā j for 1 ≤ j ≤ n − 1.
Remarkably, they are sufficient to determine a n (φ) [30,31]. To be precise, we need two assumptions to accomplish this. First, we presume that a n (ā n ) is analytic . This is consistent with the observation from BAE: as φ → − 2 3 π one of the roots λ r → 0 thus a n diverges. This implies the existence of overlapping of the analytic strips near the origin of φ for both a n andā n . Second, the following asymptotic behavior of a n as φ → ∞ is postulated, a n ∼ α n φ 2π where the coefficient reads explicitly, This is deduced from the analysis on the Bethe ansatz equations in the large φ limit, as shortly discussed in Appendix A. The first member, α 1 can also be fixed by the expansion of the Y-functions for the AdS 4 case corresponding to φ = 0 [22].
Once these assumptions are taken for granted, we can successively determine a n .
The first order equation is simply solvable and one finds For Re φ > 0, the second moment is given explicitly, For Re φ < 0, a 2 (φ) is evaluated by the analytic continuation of the above expression.
A double integral formula for the third moment a 3 is derived similarly. For the evaluation up to Y (3,0) j , however, the explicit forms of a 3 andā 3 are dispensable.
This is due to the fact that only the sum, a 3 +ā 3 , appears in Y (3,0) j and this sum also appears in the condition (4.17) for n = 3, in the special case γ = 2π 3 , In this way, one can determine Y (p,0) j successively. In particular, the first nontrivial coefficients are especially simple, to be given explicitly by As a check, one can compare the above result with the relations of Y By recovering the phase ϕ by Y j (θ) → Y j (θ − iϕ), we finally obtain the expansion of the massive Y-functions,  Figure 1 (a) shows a plot of Y 2 (0) as ℓ varies. As an example, the phase is fixed to be ϕ = −π/20 and the chemical potential to be µ = 10. We have chosen a real µ (imaginary φ) so that we can compare our data against the three-loop result in term of the multiple polylogarithms [34] , which is discussed in the next section.
Our expansions are valid also for real µ. We find a good agreement between the numerical results and our analytic expansion. A fit by a function  Figure 1 (b) shows a plot of the expansions of Y 2 (0) obtained from the analytic data in the CFT limit, i.e., Y (p,0) 2 , up to p = 1, 2, 3, respectively. ϕ and µ are the same as in (a). Up to ℓ ∼ 1, the CFT data approximate Y 2 relatively well. Combining them with the Y-system yields better approximation. To evaluate the amplitudes, we first recall that the remainder function in (2.26) consists of three terms. One of the terms denoted by A free is nothing but the free energy of the Z 4 -symmetric integrable model twisted by φ. Since this Z 4 -symmetric integrable model reduces, in the UV limit, to the twisted Z 4 -parafermion CFT, the free energy A free is expanded by the bulk conformal perturbation theory [27]:

Application to six-point amplitudes at strong coupling
where and γ(z) = Γ(z)/Γ(1 − z). The second term cancels A period .
The expansion of the remaining term ∆A BDS is derived by utilizing the results in (2.23), (2.24) and (2.26),  Here, U 0 is the value of U k in the UV limit, common to all k = 1, 2, 3. We note that ∆A BDS is expanded in powers of Z 4/3 = |Z| 4/3 e i 4 3 ϕ and its complex conjugate, but the ϕ-dependence remains only at O(|Z| 8n/3 ) with n ∈ 3Z. This is a consequence of the Z 6 -symmetry of the six-point amplitudes θ → θ + πi/4 or ϕ → ϕ + π/4, which corresponds to cyclically renaming the cusp points x a → x a+1 . This symmetry also explains the absence of the O(|Z| 4/3 ) term.
Combining the above results, we obtain the UV expansion of the remainder function, where Inverting this relation, one can express the parameters in the TBA equations by the cross-ratios [27], which is valid up to O(ℓ   u k obtained from (5.7). Figure 3 (a) shows the remainder function at strong coupling R for the cross-ratios u k given in Figure 2. The points represent the numerical results, whereas the solid line represents our expansion R = r (0) + r (1) ℓ 8/3 . Figure 3 (b) is the same plot for small ℓ. From Figure 2 and 3, we find again a good agreement between the numerical results and our analytic expansions for small ℓ.

Comparison with perturbative results
In [20][21][22]33], the remainder functions of the strong-coupling amplitudes corresponding to the minimal surfaces in AdS 3 and AdS 4 were compared with the two-loop results. It was found there that after an appropriate normalization/rescaling they are close to each other. In this section, we compare the six-point remainder function at strong with those at two [48,49], three [34] and four [10] loops.
Since A free → 0 and A period cancels the leading term from ∆A BDS we are left with R strong IR = π 2 /12. On the perturbative side, R UV are obtained from the cross-ratios U 0 in the UV limit, whereas the IR values just vanish R (L) IR = 0. We evaluate these rescaled remainder functions for the cross-ratios u k (= 1/U k ) given in Figure 2, which are parametrized by ℓ with ϕ = −π/20, µ = 10. At strong coupling, it is readily read off from the results in the previous section. At two loops, we use the simple analytic expression given in [49], whereas at three loops we use the expression given in [34] and evaluate it by using the C++ library GiNaC [50,51].
The direct three-loop expression in terms of the multiple polylogarithms is given for the parameter region with real µ, which corresponds to the (2,2) signature of the four-dimensional space-time. In order to use this expression, we have chosen a real µ. At four loops, it can be evaluated from the analytic expression in [10]. 6 Figure 4 (a) is a plot of the rescaled remainder functions at two, three and four loops and at strong coupling. Figure 4 (b) is a plot of the ratios of these rescaled remainder functions;R (3) /R (2) ,R (4) /R (3) andR (2) /R strong . As ℓ increases, it becomes harder to evaluate R (3) and Figure 4 includes the data for ℓ ≤ 7. It also includes the  Table 1. From these figures and the table, we find that the rescaled remainder functions stay close to each other as ℓ varies, but the perturbative results gradually move away from those at strong coupling as the number of loops increases. Their ratios changes slowly along ℓ, and accumulate to 1 for large ℓ as assured by definition. The ratios of the perturbative results,R (3) /R (2) , , are very similar. These are in accord with the observations in [20][21][22]33] mentioned above and those in [10,34] that the ratios of the remainder functions at two, three and four loops are relatively constant for large ranges of the cross-ratios.
As µ increases, we have observed thatR (2) 6 ,R 6 andR strong 6 change in a similar manner: they tend to start decreasing for larger ℓ. Although their differences increase, they are still kept relatively close to each other. The dependence on µ seems very weak. For example, even for µ = 10 6 , the ratios are still of O(1). At the order of the expansion in (5.5) and (5.6), the ϕ-dependence does not appear. Although it does appear at higher orders, the behavior of the rescaled remainder functions is still qualitatively similar.

Conclusions
In this paper we studied six-point gluon scattering amplitudes in N = 4 super Yang-Mills theory at strong coupling by using the AdS/CFT correspondence. The area of the corresponding null-polygonal minimal surface in AdS 5 is evaluated by solving  Table 1: Samples of the numerical values of the rescaled remainder functions plotted in Figure 4.
the T-/Y-system for the Z 4 -symmetric integrable model with a twist parameter.
The leading expansion to the remainder function is determined around the UV limit explicitly. We compared this result with the recent perturbative calculations, to find that the rescaled remainder functions are close to each other along the trajectories parametrized by the scale parameter.
Our results at the leading order relied on the quantum Wronskian relation, which determines the expansion of the T-/Y-functions around the CFT (massless) limit.
For higher order terms, one needs to study the massive TBA system intrinsically. As in [20][21][22], the boundary CPT would be a way toward this direction based on the relation between the T-/Y-functions and the g-function [28,52,53]. Through the T-/Y-functions, one could also study the quantum Wronskian relation from a different perspective by using the boundary CPT. However, it is yet to be figured out how to incorporate the twist of the perturbed Z 4 -parafermion theory in this framework.
Alternatively, the quantum Wronskian relation exists even in the massive case [46]. The more involved analyticity, however, defies the analytical determination of the moments of A,Ā. Along [28,29], the massive TBA systems have also been analyzed in [54,55]. Hopefully one can detour the difficulties, and obtain the systematic higher expansions.
Finally, given the recent developments on the finite-coupling amplitudes around the collinear limit [12][13][14], it would be worthwhile to explore the extrapolation to the finite coupling also around the regular-polygonal limit. We hope to come back to these issues in near future.
Scientific Research (C) No. 23540290, 24540248 and 24540399. The work of K. I. and Y. S. is also supported in part by JSPS Japan-Hungary Research Cooperative Program.

Appendix A. Non-Linear Integral Equation
We supplement a treatment on the Bethe ansatz equations based on suitably chosen auxiliary functions. The approach utilizes Non Linear Integral Equations satisfied by them [26,56], and thus is referred to as the NLIE approach. 7 Here we choose auxiliary functions as , ,Ā(x) = 1 +ā(x).
The Bethe ansatz equations are equivalent to a(x j ) = −1.
Numerically, one observes the following properties of the finite size system in the ground state: 1. |a(x)| < 1 (|ā(x)| < 1) in the upper (lower) half plane.
2. In a narrow strip including the real axis, the zeros of A(x) (Ā(x)) are located only on the real axis and they coincide with the Bethe ansatz roots {x j }.
3. The extended branch cut function, 1 i ln a(x), is an increasing function of real x.
These properties are enough in deriving the NLIE. After the scaling and the conformal limit, for Im θ positive small, it is explicitly given by ln a(θ) = d a (θ) + ∞+i0 −∞+i0 where we use the same symbols a etc for functions of θ(= 3x 2 ) by abuse of notation.
By comparing the above expansion with (4.10), making use of the explicit form of B in (A.3), we have the following asymptotic behavior of a n (φ), The explicit form of the coefficient α n in (4.18) is then determined as in (4.19).