Octagon at finite coupling

We study a special class of four-point correlation functions of infinitely heavy half-BPS operators in planar N=4 SYM which admit factorization into a product of two octagon form factors. We demonstrate that these functions satisfy a system of nonlinear integro-differential equations which are powerful enough to fully determine their dependence on the 't Hooft coupling and two cross ratios. At weak coupling, solution to these equations yields a known series representation of the octagon in terms of ladder integrals. At strong coupling, we develop a systematic expansion of the octagon in the inverse powers of the coupling constant and calculate accompanying expansion coefficients analytically. We examine the strong coupling expansion of the correlation function in various kinematical regions and observe a perfect agreement both with the expected asymptotic behavior dictated by the OPE and with results of numerical evaluation. We find that, surprisingly enough, the strong coupling expansion is Borel summable. Applying the Borel-Pade summation method, we show that the strong coupling expansion correctly describes the correlation function over a wide region of the 't Hooft coupling.

D Profile function in the null limit 41 1 Introduction In this work, we continue the study of four-point correlation functions of heavy singletrace half-BPS operators in maximally supersymmetric N = 4 super Yang-Mills theory initiated in Ref. [1]. It has recently been realized in Refs. [2,3] that, under a judicial choice of polarizations of half-BPS operators, these correlation functions can be constructed at finite 't Hooft coupling λ = g 2 YM N c in terms of a fundamental building block -the octagon O. Taking advantage of integrability of the effective two-dimensional worldsheet in the AdS/CFT correspondence, the octagon can be expressed in terms of the hexagon form factors [4][5][6][7] describing dynamics of magnons propagating on the string worldsheet.
At present, explicit expressions for the octagon are known in planar N = 4 SYM both at weak [2,3] and strong coupling [8,9]. At finite coupling, a concise representation for the octagon as a determinant of a semi-infinite matrix was worked out in Refs. [10,11]. This result laid out the foundation for our analysis in Ref. [1], where the octagon was further recast as a Fredholm determinant of an integral operator defined on a semi-infinite line. The kernel of the operator in question turned out to be closely related to the Bessel kernel that previously appeared in the study of the Laguerre ensemble in random matrix theory [12,13].
At weak coupling, the octagon is given by a multilinear combination of the so-called ladder integrals [3]. The corresponding expansion coefficients can be determined to any loop order by requiring the octagon to have an appropriate asymptotic behavior in different kinematical limits. At strong coupling, for g 2 = λ/(4π) 2 ≫ 1, the octagon possesses a semiclassical asymptotic behavior O ∼ e −gA 0 +O(g 0 ) , where A 0 is the minimal area of a string that ends on four geodesics in AdS [8,9].
In this paper, we determine the octagon at finite 't Hooft coupling. To achieve this goal, we exploit its representation as a Fredholm determinant of a (modified) Bessel kernel. In this formulation, the problem is akin to that encountered in study of two-point functions in integrable low-dimensional theories [14]. Following a general strategy after Its-Izergin-Korepin-Slavnov [15], we derive a system of exact equations obeyed by the octagon and demonstrate that its solution interpolates between weak and strong coupling expansions. Our consideration adds to a scarce list of examples where correlations functions in an interacting field theory can be determined for arbitrary coupling.
The starting point of our analysis is a four-point correlation function of single-trace half-BPS operators, dubbed the simplest in Ref. [2], where the four operators are built out of two complex scalars Z and X and their complex conjugate partners, O 1 = tr(Z K/2X K/2 ) + permutations, O 2 = tr(X K ) and O 3 = tr(Z K ).
Here G(z,z) is a function of the cross ratios and the coupling constant g 2 = g 2 YM N c /(4π) 2 . The correlation function (1.1) is dual to a scattering amplitude of four closed string states. In the limit when the operators become infinitely heavy, for K → ∞, it factorizes into a double copy of an open-string partition function, the octagon [2] G(z,z) (1.3) Our goal is to find O at finite coupling constant and for generic values of the cross ratios. As a first step in this direction, we can consider the null limit x 2 12 , x 2 13 , x 2 24 , x 2 34 → 0, or equivalently z → 0 − andz → ∞, when the four operators in (1.1) are light-like separated in a sequential manner. Introducing convenient kinematical variables ξ = − 1 2 log(zz) and y = − 1 2 log(z/z), one finds that in the null limit y ≫ ξ the octagon takes on a remarkably simple form at weak coupling [3,11] log O = − Γ(g) 2π 2 y 2 + C(g) 8 + g 2 ξ 2 + . . . , (1.4) where the ellipses denote subleading terms which vanish for y → ∞ with fixed g and ξ.
The dependence on the coupling constant enters (1.4) through the two functions Γ(g) and C(g). In our previous study [1], we used the aforementioned determinant representation of the octagon to find their exact expressions Γ(g) = log(cosh(2πg)) , C(g) = − log sinh(4πg) 4πg . (1.5) Taking the relation (1.4) at face value and going to the strong coupling limit results in some puzzling consequences. As alluded to above, log O is expected to have a linear growth in g at most. It is easy to see that the first two terms in the right-hand side of (1.4) do have such a scaling whereas the third term grows as g 2 .
The strong coupling limit of the octagon was studied in Ref. [9] using a clustering technique previously developed in Ref. [16] in application to three-point correlation functions. In the limit g ≫ y ≫ ξ, the octagon takes the form log O = − g π y 2 − gπ + g π ξ 2 (log y + 1 + γ − log(2π)) + O(ξ 4 ) + O(g 0 ) , (1.6) where γ is the Euler-Mascheroni constant. Comparing this relation with (1.4) and (1.5), we observe that the leading O(y 2 ) terms coincide at strong coupling. For the second term in Eq. (1.4), however, the strong coupling limit of C(g)/8 yields a result which is twice smaller than that in (1.6). Finally, the ξ−dependent terms in (1.4) and (1.6) exhibit different dependence on g and y. These observations are troublesome, but a priori not completely unexpected since the two results (1.4) and (1.6) were obtained in two different regions of the parameter space, for g ≪ y and g ≫ y, respectively. This suggests that the difference between (1.4) and (1.6) at strong coupling is due to an order of limits [9]. We show below that this is indeed the case.
To address this question, we develop a systematically expansion of O in powers of 1/g for general values of the cross-ratios (1.2). Doing so, one may attempt to apply the clustering technique [16]. The method is very efficient for extracting the leading large g asymptotics but it generalization to subleading terms suppressed by powers of 1/g proves to be difficult. 1 Presently, we demonstrate that these difficulties can be avoided by solving the abovementioned system of exact equations for the octagon at large g. We work out its strong coupling expansion in different kinematical limits and compare it with numerical results. We observe that the strong coupling expansion is given by a Borel summable series in 1/g. Applying the Borel-Padé summation method, we show that the resulting expression for the octagon agrees with its numerical value over a wide range of the coupling constant.
The paper is organized as follows. In Section 2, we review the representation of the octagon as a determinant of a semi-infinite matrix obtained in Refs. [10,11]. We demonstrate that this matrix can significantly be simplified by an appropriate similarity transformation allowing us to express the octagon as a Fredholm determinant of the Bessel kernel modified by some function of the coupling constant and cross ratios. Applying the method of differential equations [15] to this determinant, we derive a system of nonlinear integro-differential equations for the octagon. Its solution at weak coupling is presented in Section 3. In Section 4, we discuss the properties of the octagon at strong coupling and calculate the leading term of its expansion at large g. In Section 5, we develop a systematic series of the octagon in the inverse coupling and present analytical expressions for accompanying expansion coefficients. We use these results in Section 6 to study properties of the octagon in different kinematical regimes. In Section 7, we compare the obtained expressions with the numerical results. Section 8 contains concluding remarks. Technical details of our analysis are summarized in four appendices.

Octagon as a Fredholm determinant
In this section, we discuss the determinant representation of the octagon and use it to derive a system of exact equations that it obeys.

Kinematical limits
We will study the octagon in Euclidean and Lorentzian kinematical regimes. Depending on the choice of the region of interest, it is convenient to introduce auxiliary kinematical variables where φ = π + iy. The variables φ and ξ are real in the Euclidean regime, so that z * =z. The corresponding expressions for the cross ratios (1.2) are For future reference, we describe different OPE limits that we examine below in detail.
In this limit, two operators in the correlation function (1.1) collide and the asymptotics of O is controlled by the operators with the lowest scaling dimensions propagating in the corresponding OPE channel. We can distinguish two different OPE limits corresponding to x 1 → x 2 and x 1 → x 4 , respectively.
In the first case, the leading contribution to the octagon comes from single-trace operators with the scaling dimension ∆ = K + γ(g) and it takes the form O ∼ e −ξ(∆−K)/2 = e −ξγ(g)/2 . At strong coupling, the anomalous dimension γ(g) increases with g and the octagon is expected to vanish O → 0.
In the second case, the OPE is dominated by double-trace operators with the scaling dimension ∆ = 2K. Their contribution to the octagon scales as ∼ φ ∆/2−K = φ 0 and the octagon approaches a finite value O → 1.
In what follows we refer to the two limits in (2.3) as single-and double-trace OPE limits, respectively.
Lorentzian null limit x 2 12 , x 2 13 , x 2 24 , x 2 34 → 0 In this limit, the operators in (1.1) are located at the vertices of a null rectangle and we have The leading contribution in each OPE channel x 2 i,i+1 → 0 comes from single-trace operators with the leading twist K and large Lorentz spin. As in the previous case, their anomalous dimension grows indefinitely at strong coupling leading to exponentially small octagon O → 0, see Eqs. (1.4) and (1.6).
In this case, we have z =z * and z → −1 in the Euclidean regime. This kinematical configuration does not correspond to a particular OPE limit. The reason why we consider it is that the octagon simplifies and reveals some interesting properties.

Determinant representation of the octagon
As was mentioned in the Introduction, the octagon at finite coupling admits a representation in terms of a determinant of a semi-infinite matrix [10,11] where λ = −2(cosh y + cosh ξ) is a scalar factor depending on kinematical variables and the matrices C and K are where m, n ≥ 0 and J m is a Bessel function. In the hexagonalization approach [4][5][6], the relation (2.5) comes about as a result of resummation over an arbitrary number of elementary excitations and their bound states propagating on the worldsheet of an open string describing the octagon. At weak coupling, the representation (2.5) can be used efficiently to expand O in powers of g 2 . Since the matrix elements (2.6) scale as K nm = O(g n+m+1 ), to any given order in g 2 we can replace the K−matrix in (2.5) by its finite-dimensional minor and expand O in powers of CK. As was shown in Ref. [11], this leads to the weak-coupling expansion of the octagon in terms of known ladder integrals as was previously bootstrapped in Refs. [2,3].
At finite coupling, the calculation of (2.5) may seem to be problematic due to a rather complicated form of the K−matrix. This matrix depends in a nontrivial way on the coupling constant g and the kinematical variables y and ξ encoding the cross ratios (2.2). As we demonstrated in Ref. [1], the octagon simplifies significantly for ξ = 0. In this case, the matrix (2.6) reads where we inserted the subscript to indicate that this relation holds for ξ = 0. Notice that in distinction to (2.6), the matrix elements (K 0 ) mn vanish for even m − n and, as a consequence, the matrix λCK 0 has a block structure. Taking advantage of this property, we were able to recast the octagon (2.5) at ξ = 0 into the form of a Fredholm determinant of a (modified) Bessel kernel [12,13]. We will see momentarily that the octagon admits a similar representation for arbitrary ξ.

Similarity transformation
Denoting H = λCK, we observe that the octagon O = det (1 − H) is invariant under a similarity transformation H → Ω −1 H Ω. Choosing Ω appropriately we can simplify the form of the semi-infinite matrix H. 2 As a hint, we examine a trace of this matrix, tr H = λ tr(CK) = −2 m≥0 K m,m+1 . Replacing K m,m+1 with its expression (2.6) we get where we applied a summation formula for the Bessel functions and changed the integration variable to z = t 2 − ξ 2 . Notice that, aside from the factor λ = −2(cosh y + cosh ξ), the ξ−dependence of (2.8) only resides in the denominator of the integrand. In particular, the dependence of the integrand in (2.8) on ξ can be restored from its value at ξ = 0 by simply replacing z → z 2 + ξ 2 in the denominator of (2.8).
Applying the same recipe to Eq. (2.7), we define the matrix By construction, it satisfies the relation tr(CK) = tr(CK Ω ). It is a straightforward but tedious exercise to verify that analogous relations holds for powers of the two matrices, tr[(CK) n ] = tr[(CK Ω ) n ] with n = 2, 3, . . . . This suggests that the two matrices H = λCK and H Ω = λCK Ω are related to each other by a similarity transformation We verified that det(1 − H) = det(1 − H Ω ) at weak coupling 3 and checked this relation numerically by truncating the semi-infinite matrices to a finite size and evaluating the determinants for various values of the coupling and kinematical parameters. Although we do not need the matrix Ω for our purposes, it would be interesting to construct it explicitly. It is remarkable that the ξ−dependence of (2.9) is much simpler as compared to that of (2.6). At the same time, the matrix (2.9) has many properties in common with the matrix (2.7) evaluated at ξ = 0. In particular, the matrix elements (H Ω ) mn vanish if indices m and n have different parity. The nonzero entries of this matrix with even and odd indices respectively define two irreducible blocks. They are related to each other by a similarity transformation and, as a consequence, the octagon for ξ = 0 possesses the same blockdiagonal form as for ξ = 0 case (see Ref. [1]) (2.11) Here k − is a semi-infinite matrix which is given by (2.10) with odd indices, (k − ) mn = (H Ω ) 2m+1,2n+1 . It follows from (2.10) and (2.9) that it admits two equivalent representations where the notation was introduced for The function χ(z) approaches 1 at the origin and decreases exponentially fast at large z. It suppresses the contribution from large z to (2.12) and serves as an ultraviolet cut-off. The representation (2.11) and (2.12) is advantageous as compared to (2.5) and (2.6) because the dependence of the octagon on the coupling constant and the kinematical parameters is confined to the cut-off function (2.13). We exploit this property below to formulate a system of equations for the octagon. Moreover, the relation (2.11) can be efficiently used to compute the octagon numerically for arbitrary coupling (see Section 7 below).

Modified Bessel kernel
Following Ref. [1], we can rewrite the octagon (2.11) as a Fredholm determinant of an integral operator. This can be achieved by expanding Eq. (2.11) in terms of the traces of powers of the matrix k − where we replaced the matrix k − by its expression (2.12) and introduced a notation for 4 This function has previously appeared in the study of level spacing distributions in random matrices and it is known as an integrable Bessel kernel [13]. 4 It also admits a compact integral representation K(x 1 , The relation (2.14) suggests to define an integral operator whose kernel is given by the Bessel kernel (2.15) modified by the cut-off function (2.13) where f (x) is a test function. In this way, we obtain from (2.14) an equivalent representation of the octagon It is remarkably similar to the analogous relation for the octagon at ξ = 0 derived in Ref. [1]. The relation (2.17) holds for arbitrary values of the coupling g and generic kinematical variables y and ξ. The dependence of O on these parameters resides in the cut-off function χ(x). This function plays a central role in our subsequent analysis. To elucidate its meaning, we examine (2.13) for large y and x. In this limit, we find that χ(x) takes the form of the Dirac-Fermi distribution where the temperature T , chemical potential µ and the energy ε are related to the 't Hooft coupling and kinematical variables as Replacing x = (2gξ) 2 sinh 2 θ we note that ε = 2gξ cosh θ coincides with the energy of a relativistic particle with mass m = 2gξ and rapidity θ. Substitution of (2.18) into (2.17) and (2.16) yields an expression for O that resembles a Fredholm determinant representation of two-point correlation functions in integrable models at finite temperature [14]. Together with (2.19) this suggests that the asymptotic behavior of the octagon at weak and strong coupling should be similar to that of two-point correlation functions at low and high temperature, respectively. Indeed, a two-point function of currents in one-dimensional Bose gas model is known to decay exponentially at high temperature in the so-called pre-asymptotic region, G 2 (x) ∼ exp(−T x 2 /2) [17]. The linear temperature dependence in the exponent translates into the linear dependence of log O on the coupling, see Eq. (1.6).

Method of differential equations
A powerful technique for studying Fredholm determinants has been developed in Ref. [15] in application to two-point correlation functions in integrable models. In our previous paper [1] we extended this technique to the octagon (2.17) at ξ = 0. Due to a particular form of the kernel in (2.16), generalization to arbitrary ξ is straightforward.
To start with, we introduce the so-called potential u, in the terminology of [15], Bessel function. Following [1] we can show that it is related to the logarithmic derivative of the octagon Being a function of g, y and ξ, it satisfies differential equations where the cut-off function χ(x) is given by Eq. (2.13) and an auxiliary function Q(x) is defined as The latter, in turn, obeys a differential equation which involves the potential u and its derivative with respect to the coupling constant. This equation should be supplemented with boundary conditions at weak coupling (2.25) They follow from the expansion of (2.20) and (2.23) in powers of K χ . The relations (2.22) and (2.24) define a system of nonlinear equations for the potential. Having solved them, we can determine u for arbitrary g, y and ξ and, then, apply the relation (2.21) to compute the octagon as This equation can be efficiently used at weak coupling whereas at strong coupling it requires knowing the potential at finite coupling. Another relation for log O was found in Ref. [1] where α = {g, y, ξ}. It can be used to determine the octagon from the solution to (2.24).
We can simplify the relations (2.22) by taking into account the property of the cut-off function (2.13) We then deduce from (2.22) Being combined together with (2.21) and (2.24), these equations allow us to determine the dependence of the octagon on the kinematical variables y and ξ for any value of the coupling constant g.

Moments
To analyze the relations (2.29), it is convenient to introduce the moments where the g−dependent normalization factor was introduced for convenience. It follows from the second relation in (2.29) that for ℓ = 1 the moment is related to the derivative of the potential Moreover, it is possible to show using (2.24) that the moments satisfy a differential equation (see Ref. [1]) It relates the moments with sequential indices. In particular, for ℓ = 0 the relations (2.32) and (2.31) lead to the following equation for the moment Q 0 Obviously, Q 0 = 1 is a special solution to (2.33). At weak coupling, substituting u = O(g 2 ) in (2.33) and expanding Q 0 − 1 in powers of g 2 , it is possible to show that corrections to Q 0 vanish to all orders in g 2 . This suggests that for arbitrary values of g, ξ and y the solution to (2.33) is We demonstrated in Ref. [1] that for ξ = 0 this condition leads to the expected result for the octagon (1.4). We argue below the same is true for ξ = 0. We will use the relation (2.34) to develop a systematic expansion of the octagon at weak and strong coupling.

Octagon at weak coupling
To illustrate the power of the method of differential equations described in the previous section, we develop in this section the perturbative expansion of the octagon. For small coupling, the octagon can be easily found from the recurrence relation for the moments (2.32). It is convenient to rewrite (2.32) as Then taking into account that u = O(g 2 ) as g → 0, we observe that the expression in the right-hand side of (3.1) is suppressed by a power of g 2 . This suggests us to solve (3.1) by iterations. Namely, replacing Q ℓ and u in (3.1) by their perturbative expansions we compare the coefficients in front of powers of g 2 on both sides of (3.1) and obtain recurrence relations between the coefficients Q (k) ℓ . Solving them, we can express Q ℓ in terms of Q (0) ℓ and the expansion coefficients of the potential u (k) where a notation was introduced for the moments at zero coupling q ℓ ≡ Q Here in the first relation, we used the definition (2.30) and replaced Q(x) with its expression (2.25) at weak coupling. In the second relation, we changed the integration variable to z = (2g) √ x, applied (2.13) and replaced the Bessel function by its leading asymptotics at small g. Replacing χ(z) in (3.4) with its explicit expression (2.13), we get It is easy to see that q 0 = 1. For ℓ ≥ 1, we can show making use of the integration by parts that q ℓ is a positive definite function of ξ and y. The same is true in the Euclidean regime, for y = i(π − φ) and φ being real. According to Eq. (2.31), the moment Q 1 is related to the derivative of the potential with respect to the coupling. Substituting ℓ = 1 in (3.3) and matching the result to a weak coupling expansion of ∂ g u/(8g), we can compute the coefficients u (k) . 5 The resulting expression for the potential (3.2) is Here the coefficients in front of the powers of g 2 are given by multi-linear combinations of functions (3.5).
Note that all terms in (3.6) except the first one would vanish if the functions q ℓ had the form q ℓ = z ℓ 0 , or equivalently the function χ(z) in (3.4) was given by the step function χ(z) = θ(z 0 − z) with some z 0 . Recalling (2.18), we observe that the cut-off function takes such a form in the limit of zero temperature provided that z 2 0 = y 2 − ξ 2 . In this case, vanishing of the coefficients in (3.6) is in agreement with the known property of the Fredholm determinant of the Bessel kernel at zero temperature [12,13,18].
Substituting (3.6) into (2.26), we can compute the octagon at weak coupling This relation should be compared with the analogous expression for the perturbative octagon derived in Refs. [2,3]. We find that the two expressions coincide after we take into account that the moments (3.5) can be expressed in terms of ladder integrals where z,z are related to y, ξ through (2.1).
In the rest of the section, we examine (3.7) in the kinematical limits described in Section 2.1.
We find from (3.5) that q 0 = 1 and q 1 = 8 log 2. For ℓ ≥ 2, the function (3.5) looks as where ζ is the Euler-Riemann zeta function and η is the Dirichlet eta function. The octagon (3.7) admits the form It is convenient to assign to log 2 and ζ(k) the weight 1 and k, respectively. Then, the first few terms of the expansion have a homogenous weight but starting from order O(g 8 ), the coefficients have an admixture of lower-weight contributions.
Single-trace OPE channel ξ → ∞ and φ = 0 In this case, we substitute y = iπ in (3.5) and simplify the integrand at large ξ to get where K ℓ+ 1 2 (ξ) is the modified Bessel function. A close examination exhibits that q ℓ is a polynomial in ξ of degree ℓ with integer positive coefficients.
The corresponding expression for the octagon then reads Note that the expansion goes in powers of g 2 ξ and 1/ξ.
Double-trace OPE channel ξ = 0 and φ → 0 Replacing y = i(π − φ) in (3.5) and going through the calculation of the integral at small φ, we get q 0 = 1 and In this case, the octagon is given by 14) The expansion coefficients of O depend on log φ. Such terms come from q 1 in (3.13). Since the expression on the right-hand side of (3.7) is linear in q 1 , the expansion coefficients are linear in log φ to any order in g 2 [3].
Notice that the weak-coupling corrections to O vanish as φ → 0. In this limit, the leading contribution to the octagon comes from double-trace half-BPS operators with the scaling dimension K propagating in the OPE channel x 2 14 → 0. It is protected from quantum effects and leads to O → 1 as φ → 0.
Null limit ξ = fixed and y → ∞ At large y, the dominant contribution to (3.5) arises from z = O(y). Replacing the integration variable z → yz, we get from (3.5) in the large-y limit where B 2k+2 (x) is a Bernoulli polynomial. The function q ℓ is a polynomial in y 2 of degree ℓ.
In this case, it is advantageous to consider the logarithm of the octagon (3.7) This relation is in agreement with (1.4) and (1.5). In particular, the expansion coefficients of log O are linear in y 2 and the ξ−dependence cancels to all orders except the lowest one. The last property can be understood using the last relation in (2.29). The first term in its right-hand side is proportional to the moment Q 0 . Taking into account (2.34), we get We can show that the second term in this equation is exponentially small at large y. Indeed, at small g, or equivalently at low temperature, the function χ(x) effectively reduces the integration region in (3.17) we find that the integral in (3.17) scales as 4g 2 y 2 + O(g 4 ). It is accompanied, however, by the factor of sinh ξ/(cosh y + cosh ξ) ∼ e −y and, therefore, provides a vanishing contribution to (3.17) at large y. Thus, ∂ ξ u = −8g 2 ξ at weak coupling leading to ∂ ξ log O = 2g 2 ξ in virtue of (2.26). As was mentioned in Section 2.1, the relation (3.16) describes the asymptotic behavior of the null octagon. It arises as a result of infinite resummation of contributions of leading twist operators in different OPE channels. Exponentially small corrections, generated by the second term in the right-hand side of (3.17) at weak coupling, are induced by high twist operators. We show below that at strong coupling the situation is completely different. Anomalous dimensions of the operators grow with the coupling and their classification with respect to twist (= bare dimension minus spin) becomes redundant. We demonstrate in Section 6 that the second term in (3.17) scales at strong coupling as 8g 2 ξ + O(g), so that u = O(g) in agreement with the expected behavior of the octagon (1.6).

Octagon at strong coupling
In this section, we study the octagon in the limit when g → ∞ with the kinematical variables y and ξ held fixed.
In this limit, the octagon is expected to have the following form where the coefficient functions depend on the kinematical variables y and ξ and relative rational prefactors are introduced for convenience. The sum contains terms of the form The first term in (4.1) was computed in Ref. [9] using the clustering procedure developed in Ref. [16] (see Eq. (4.10) below). It was argued there that A 0 should be related to the minimal area of a string that ends on four geodesics in AdS and rotates on the sphere. The remaining terms in (4.1) remain unknown.
The subleading term A 2 1 log g/2 + B in (4.1) describes quadratic fluctuations around the minimal area. It is enhanced by log g and produces an overall g−dependent factor of the Here O(1/g) terms in the exponent take into account yet higher order quantum fluctuations.
We are going to determine the coefficient functions in (4.1) from the system of equations (2.21), (2.29) and (2.24). Combining together (4.1) and (2.21) we expect the potential to have the following form at strong coupling Comparing this relation with (4.1) we notice that B does not contribute to u. This means that having determined the potential u, we will be able to determine the octagon (4.1) up to the function B(y, ξ).

Leading order
To find the leading term A 0 , we use the relation (2.14) and examine the properties of its iterated integral representation in the limit g → ∞.
It follows from the explicit form of the cut-off function (2.13) that the dominant contribution to (2.14) comes from integration over x i = O(g). This allows us to replace the K−kernel in (2.14) by its leading asymptotic behavior at large x i Here the ellipses denote subleading terms suppressed by powers of 1/x i . We notice that for x i = O(g), the second term inside the square brackets is a rapidly oscillating function of x i . At large g it scales as O(1/g 1/2 ). At the same time, the first term is peaked around x 1 = x 2 and scales as O(g 0 ) for x 1 − x 2 = O(1/g). This suggests that for g → ∞, the integral in (2.14) is localized at To show this, we use the identitŷ where f (x) is a smooth test function. Here in the second relation, we changed the integration variable to z = √ x 2 − √ x 1 and expanded the integral at large x 1 and fixed z.
Subsequently applying (4.5), we get from (2.14) As before, we can replace K(x 1 , x 1 ) with its leading asymptotic behavior at large x 1 . Applying (4.4) for . Then, changing the integration variable to x 1 = (2gz) 2 , we finally obtain in the leading large g limit According to its definition (2.13), the function χ(z) is independent of the coupling constant but carries dependence on the kinematical variables y and ξ. Comparing (4.7) with (4.1), we deduce that As follows from this relation, A 0 is a positive definite function of real y and ξ. The same is true in the Euclidean regime for y = i(π − φ) with φ real.
As was mentioned at the beginning of this section, the function A 0 was also computed in Ref. [9] using a different technique. To compare the two results, we change the integration variable in (4.8) to z = ξ sinh θ and introduce a new function . (4.9) Then, the relation (4.8) takes the form and it coincides with the analogous expression obtained in Ref. [9].

Beyond leading order
To derive the strong coupling expansion (4.1), we use the system of integro-differential equations (2.21), (2.24) and (2.29). Namely, we will construct the function Q(x) at strong coupling for arbitrary u and, then, apply the second relation in (2.29) to determine the expansion coefficients of the potential in (4.3).
The function Q(x) depends on the coupling constant and obeys Eq. (2.24). It is convenient to replace x = (2gz) 2 and introduce the function (4.11) It follows from (2.24) that it satisfies the following differential equation where a notation was introduced for . . (4.13) and the potential u was replaced with its general expression at strong coupling (4.3). Note that the leading O(g) term in (4.3) proportional to A 0 does not contribute to w(g). Substituting (4.11) into (2.29), we get the system of equations relating the potential and solutions to (4.12) with χ(z) defined earlier in Eq. (2.13). In a similar manner, it follows from (2.34) that Solution to (4.12) yields the function q(z) that depends on u(g) in a nontrivial way. Its substitution into (4.14) and (4.15) yields a complicated system of equations for the potential. Luckily, these equations can be solved at strong coupling.

Strong coupling expansion
In this section, we use the relations (4.12) -(4.15) to systematically calculate the coefficients A 1 , A 2 , A 3 , . . . of the strong coupling expansion of the octagon (4.1).

Next-to-leading order
To begin with, we neglect the O(1/g) correction to (4.13) and examine the differential equation (4.12). In this case, for w = −A 2 1 , a general solution to (4.12) reads where J α and Y α (with α = A 1 ) are Bessel functions of the first and second kind, respectively, and c(z), c ′ (z) are arbitrary functions of z. The last term in the right-hand side of (5.1) denotes corrections due to O(1/g) terms in (4.13). For arbitrary α, the function Y α has singularity at the origin and is multivalued. This suggests to put c ′ (z) = 0.
To find the function c(z), we examine the second relation in (4.14) where in the last relation we replaced u with its expression (4.3). Evaluating the integral in (5.2), we can replace the function q(z) with its leading asymptotic behavior at large g, Finally, substituting this expression in the previous relation, we deduce Here we took into account that the dominant contribution to the integral comes from z = O(g 0 ) and, as a consequence, a rapidly oscillating sinus function does not contribute. Matching the last relation to the right-hand side of (5.2) and replacing A 0 with (4.8), we conclude that c 2 (z) = 1/(1 − χ(z)). Thus, the solution (5.1) looks as (5.5) Here the second term in the numerator denotes corrections due to O(1/g) terms in (4.13). Substituting (5.5) into (4.14) and repeating the same analysis, it is straightforward to verify that the two remaining relations in (4.14) are automatically satisfied.
Having constructed the function q(z), we can examine now the normalization condition (2.34) In distinction to (5.4), the leading contribution to this integral comes from z = O(1/g). In this region, we can replace the cut-off function (2.13) by its leading behavior around the origin ∂ z log(1 − χ(z)) = 2/z + . . . . We can thus continue the previous relation to get We conclude that Thus, the subleading O(log g) correction to the octagon (4.1) is universal, i.e., it does not depend on the kinematical variables y and ξ.
Combining together (5.5) and (5.8), we find that the solution to the differential equation (4.14) looks as We can check this relation by computing the function q(z) numerically at some reference value of the coupling and the kinematical parameters as described in Appendix C. Its comparison with the leading term in (5.9) is shown in Figure 1.
It is important to stress that at strong coupling the integrals in (5.2) and (5.6) receive the leading contributions from two different regions, z = O(g 0 ) and z = O(1/g), respectively. This suggests that, in order to compute the corresponding expressions for ∂ g u and Q 0 , we have to construct the function (5.9) in these two regions. For z = O(g 0 ) this is done in the next subsection and for z = O(1/g) in Appendix B.

Solution at finite z
At large g and z = O(g 0 ), it is convenient to switch from g to x = gz and expand the function (5.9) in powers of 1/x. Since the function f differs from q(z) by a g−independent factor, it satisfies the very same differential equation (4.12) where we used (4.13) and replaced g = z/x.
To leading order in 1/x, the solution is given by (5.9). Using the properties of Bessel functions, we find from (5.9) that f is an oscillating function Going beyond leading order, we look for a general solution to (5.10) in the form where a and b are given by infinite series in 1/x Substituting (5.12) into (5.10) and matching the coefficients in front of powers of 1/x, we can determine the functions a(x, z) and b(x, z) to any order in 1/x As a check, we verify that for z = 0 the functions a(x, 0) and b(x, 0) coincide with the coefficients accompanying trigonometric functions inside the brackets in Eq. (5.11).
Let us now examine (5.2) and replace the function q(z) with its expression (5.9) and (5.12) We recall that the cut-off function (2.13) does not depend on the coupling constant and the dependence on g resides in the first factor in the right-hand side of (5.16). At strong coupling this factor contains rapidly oscillating trigonometric functions. Expanding the integral at large g, we can safely replace them by their averaged values. This leads to ∂ g u = 2 πˆ∞ 0 dz z∂ z log(1 − χ(z)) a 2 (2gz, z) + b 2 (2gz, z) .

(5.17)
It is important to stress that this relation was derived under the assumption that the dominant contribution to the integral comes from z = O(g 0 ). Using (5.14) and (5.15), we find Notice that the expansion coefficients contain only even powers of 1/z. We would like to emphasize that the expansion (5.18) is well-defined for g ≫ 1 and z = O(g 0 ).
Substituting (5.18) into (5.17), we can expand ∂ g u in powers of 1/g with a notation I n introduced for the integral As we will see in a moment, the dependence of the octagon on the kinematical variables y and ξ enters through I n (y, ξ). This is the reason why we will refer to I n (y, ξ) as a profile function. Strictly speaking, the integral in (5.20) diverges at the lower limit for n ≥ 1 and requires a regularization. We address this issue in Section 5.2.2.

Quantization conditions
The leading term in the expansion (5.19) looks as where we integrated by parts and matched the result to Eq. (4.8). We verify that this relation is in agreement with (4.3). Then, we replace u on the left-hand side of (5.19) with its general expression (4.3) and compare the coefficients in front of powers of 1/g. This gives recurrence relations between the coefficients A k These relations allow us to express all the coefficients A k in terms of the functions I n (y, ξ) defined in (5.19). The explicit expressions for the first few of them are where we also included for completeness the values of the lowest coefficients A 0 and A 1 . It is straightforward to extend these relations to arbitrary order in 1/g expansion. We present the explicit expressions for the coefficients A k (with 0 ≤ k ≤ 40) in an ancillary file with our paper. Relations (5.23) combined with (4.1) yield a strong coupling expansion of the octagon.
We notice that A 3 = −A 2 2 . Similar relations hold for all coefficients with odd indices. Namely, the coefficients A 2k+1 can be expressed in terms of the even coefficients A 2m (with m ≤ k), e.g., 15 , The origin of these relations is elucidated in Appendix B.

Profile function
A close examination shows that the integral in (5.20) is not well-defined for positive integer n and requires a regularization. Indeed, it is easy to see that log(1 − χ(z)) ∼ log z for z → 0 and, therefore, the integral in (5.20) diverges at the lower limit as´dz/z 2n . To understand the origin of this divergence, we turn back to Eq. (5.17). The integral in (5.17) is convergent at small z and the divergences appear in (5.19) after we exchanged the integration and series expansion of the integrand at large g. As mentioned above, the series expansion is well-defined for z = O(g 0 ) and, therefore, performing the integration we should have imposed a lower cut-off on the possible values of z.
The simplest way to do it is by introducing the so-called analytical regularization (see Appendix A). Namely, we modify the definition (5.20) by inserting an addition factor z ǫ that suppresses the contribution of small z For any given n we choose ǫ sufficiently large so that the integral is convergent. Then, integrating by parts, we can reduce the strength of singularity at small z Subsequently integrating by parts n times we arrive at the integral that is well-defined for ǫ → 0. Taking the limit I n = lim ǫ→0 I n (ǫ) we arrive at It is easy to verify using z∂ z log(1 − χ(z)) ∼ 1 + O(z 2 ) that it is well-defined for z → 0. The explicit expressions for the profile function I n (y, ξ) in different kinematical limits are given in the next section.

Properties of strong coupling expansion
At strong coupling, the octagon is given by (4.1) with the expansion coefficients defined in (5.23). In this section, we examine properties of the obtained expressions.

Improved expansion
According to (5.23), the expansion coefficients A n (with n ≥ 2) are given by multi-linear combinations of the profile functions I k A n = c n I n−1 where the nonnegative integers p i satisfy the equation The rational coefficients c n and c p 1 ,...,p ℓ exhibit a remarkable regularity. To make this manifest, we rewrite the octagon (4.1) in the following form where log O q is given by Replacing the coefficients A n with their expressions (5.23), we find Notice that the coefficients involve terms of the form I ℓ 1 , I 2 I ℓ 1 , I 3 I ℓ 1 , . . . . Such terms can be resummed to all orders in ℓ where the dots denote further terms of the infinite series. The relation (6.6) suggests to change the expansion parameter to g ′ = g − I 1 /2. Expanding log O q at large g ′ , we get a remarkably simple expression where the expansion coefficients do not depend on I 1 . The expressions for the octagon (6.3) and (6.7) up to order O(1/g 38 ) can be found in an ancillary file. The situation here is similar to that of the strong coupling expansion of the cusp anomalous dimension. In the latter case, all corrections proportional to log 2 could be absorbed into the redefinition of the coupling constant [19]. In the present case, the shifted coupling constant, g ′ = g − I 1 /2, depends on the kinematical variables y and ξ.

Kinematical limits
In this subsection, we evaluate the profile function (5.27) and study the properties of the octagon at strong coupling in different kinematical regimes.

Symmetric point y = ξ = 0
In this kinematical point, the profile function (5.27) is given by I 0 = π/2 and I 1 = −2 log 2/π. For n ≥ 2, we have This relation should be compared with its counterpart at weak coupling (3.9). Substituting (6.8) into (5.23), we obtain for the octagon (4.1) As explained earlier, all terms in this expression involving log (2) can be absorbed into the redefinition of the coupling g ′ = g + log 2/π. Assigning a weight 1 to π and log(2) and weight k to ζ(k), we observe that the first two terms in (6.9) possess weight 1, whereas all terms suppressed by powers of 1/g have a uniform weight 0. A close examination of (6.9) shows that, starting from the O(1/g) term, the expansion coefficients have alternating signs and grow factorially at higher orders. 7 This suggests that the strong coupling expansion (6.9) is Borel summable.
As was mentioned above, B in (6.9) arises as an integration constant in (2.21) and, therefore, it does not depend on the coupling constant. For y = ξ = 0, we can fix its value by comparing (6.9) with the numerical result for the octagon at some reference value of the coupling 1 < g < 10. In this manner, we deduce that B = 0.960877 .
We observe that in all terms in Eq. (6.11), except the second one, the coupling constant is accompanied by (2πξ) 1/2 , so that the expansion parameter is effectively (2πξ) 1/2 g. Assuming that the O(g 0 ) term in (6.11) has the same property, we can predict the ξ−dependence of the constant term B B = 1 4 log(2πξ) + c . (6.12) We verified that this relation agrees with the numerical values of the octagon at large ξ and extracted the constant c = 0.343754. It is interesting to note that the weak coupling expansion of the octagon (3.12) also runs in (even) powers of (2πξ) 1/2 g. This suggests to introduceg = (2πξ) 1/2 g and compare the dependence of the octagon ong at weak and strong coupling. Neglecting corrections suppressed by powers of 1/ξ, we observe (see Figure 2) that the two curves defined by (3.12) and (6.11)

merge at intermediate values ofg.
Double-trace OPE channel ξ = 0 and φ → 0 As was already emphasized in Section 2.1, the leading contribution to the correlation function in this limit comes from double-trace operators. We expect that for φ → 0 with g kept fixed, the octagon should take a finite value O → 1.
To approach this limit, we replace y = i(π − φ) and ξ = 0 in the profile function (5.27) and take φ → 0. The leading contribution to the integral in (5.27) comes from z = O(φ) whereas integration over z > 1 in (5.27) yields a subleading contribution as φ → 0. Replacing the integration variable to x = z/φ with 0 ≤ x ≤ 1/φ, we can expand the integrand in (5.27) in powers of φ 2 using the relation This leads to the following result for the profile function I n for n ≥ 1 where the last term accounts for the contribution from z = O(φ 0 ). For n = 0 we have Substituting (6.14) into (4.1) and (5.23), we obtain the strong coupling expansion of the octagon where the dots stand for higher-order terms in 1/g as well as terms vanishing as φ → 0. Notice that the expansion in (6.15) runs in powers of 1/(gφ) and, in order to take the limit φ → 0, the series needs to be resummed. It is convenient to switch to α = 8gφ and examine the limit φ → 0 with α held fixed and small. Then, the relation (6.15) takes the form Assuming that log O depends on φ and g through α, we can determine the asymptotic behavior of B at small φ B = 1 2 log(8φ) . (6.17) In general, we can add an arbitrary constant to the right-hand side of this relation. As we show in a moment, the condition for O to vanish for φ → 0 leaves however no room for it.
To get a closed expression for the series (6.16), we examine the expansion coefficients of ∂ α log O. It turns out that they form a sequence that had already appeared in the literature, see Refs. [20,21], leading to , (6.18) where K 1 and I 1 are modified Bessel functions. The integrand decreases exponentially fast at large s and behaves as α/(16(α + s)) at small s. Integrating by parts, we get for small α Imposing the condition that log O has to vanish for α → 0, we finally arrive at We recall that this relation was obtained for φ → 0 and α = 8gφ fixed. The comparison of (6.21) with the numerical value of log O at g = 10 is shown in Figure 3. Notice that the resummed expression (6.21) vanishes for α → 0 faster than the leading contribution at strong coupling. The latter is described by the first term on the right-hand side of (6.16). It is also interesting to note that the leading asymptotic behavior of the octagon at strong coupling log O ∼ 2g 2 φ 2 log(gφ) is remarkably similar to that at weak coupling log O ∼ 2g 2 φ 2 log φ (see Eq. (3.14)) even though the two expressions are valid in different regions of the parameter space.
Null limit y ≫ 1, ξ < 1 In this limit, the profile function (5.27) can be expanded in powers of ξ 2 where the coefficient functions I The leading term is with the notation introduced for L = log y + γ − log(2π). The subleading terms are given by As before, the expansion coefficients in front of powers of 1/g in the expressions for f (i) have a uniform weight 0.
The expressions for f (i) contain terms enhanced by powers of L = log y + γ − log(2π). Such terms can be resummed to all orders using the property of the octagon exhibited in Eqs. (6.3) and (6.7), leading to (6.28) where g ′ = g − I 1 /2 and the profile functions I 0 , I 1 , . . . are given by The first term in (6.28) agrees with (1.6), as anticipated. The expansion (6.28) is well-defined for large g ′ = g − I 1 /2, or equivalently for g ≫ L/(2π). Comparing (6.23) with (1.4), we notice that the dependence of the null octagon on ξ changes as we go from weak to strong coupling. In the former case, it is one-loop exact log O ∼ g 2 ξ 2 whereas in the latter case, log O scales as O(g) and its small ξ−expansion does not truncate.
The transition between the two regimes can be understood using the last relation in (4.14). We demonstrated in Section 3, that at weak coupling the second term in the righthand side of this relation is exponentially small at large y. At strong coupling, we replace the function q(z) in (4.14) with its leading asymptotic behavior (5.9) and take into account (4.15) to get where the ellipses denote corrections suppressed by 1/g. At large g, the leading contribution to the integral comes from z = O(1/g). Changing the integration variable to x = 2gz and expanding the integrand at large g, we find that the second term on the right-hand side of (6.30) behaves at small ξ as (6.31) Substituting this relation into (6.30), we verify that the leading O(g 2 ) term cancels yielding ∂ ξ u = O(gξ) or equivalently log O = O(gξ 2 ) in agreement with (6.23).

Numerical checks
In this section, we compare numerical values of the octagon in different kinematical limits with the corresponding analytical expressions obtained in the previous section. The numerical results presented in this section were obtained in collaboration with Riccardo Guida.
To compute O numerically for finite 't Hooft coupling, we apply the relations (2.11) and (2.12) and truncate the size of the matrix k − to N max = 100. In addition, we also compute a logarithmic derivative of the octagon with respect to the coupling constant, ∂ g log O. According to (2.21), it is related to the potential, u = −2g∂ g log O, which can be found from (C.6) by replacing the matrix k − with its finite-dimensional minor.
To study its asymptotic properties, it is advantageous to consider the logarithmic derivative ∂ g log O. It does not contain the constant term B and its expansion involves only inverse powers of the coupling.
To begin with, we consider the octagon at the symmetric point y = ξ = 0. According to (6.9), it is given by a sign-alternating series in 1/g with factorially growing coefficients. Following the standard procedure [22], we can improve the strong coupling expansion (7.1) by applying the Borel transformation and replacing a partial sum B(t) by its Padé approximant [n/m] = P n (t)/Q m (t).
For n = 2 and m = 3, the resulting expression for ∂ g log O is shown in Figure 4. It agrees with the numerical values of the octagon up to g = 0.1. We encounter similar situation in the single-trace regime, for φ = 0 and ξ ≫ 1. Strong coupling expansion of the octagon (6.11) is Borel summable and it can be improved using the Borel-Padé method. We verified that for ξ = 4 the resulting strong coupling expansion agrees with numerical values of ∂ g log O up to g = 0.1.

Order of limits phenomenon
As mentioned in the Introduction, in the null limit, for y ≫ 1, the octagon has different asymptotic behavior depending on the hierarchy between y and g. At weak coupling, for y ≫ g, the octagon scales as (1.4). At strong coupling, for g ≫ y ≫ 1, its leading behavior is given by (1.6) and subleading corrections are defined in Eqs. (6.23) and (6.28).
We recall that at strong coupling the relations (1.4) and (1.6) differ by the terms proportional to gπ and they exhibit different dependence on ξ. We will verify them in two steps. First, we put ξ = 0 and examine the dependence of ∂ g log O(0) on the coupling constant at y = 10. Second, in order to check the ξ−dependence, we consider the logarithmic derivative of the ratio of octagons ∂ g log(O(0)/O(ξ)) and examine its g−dependence for ξ = 1/10. We observe (see Figure 5) that the relations (1.4) and (6.28) are in agreement with numerical values at lower and higher values of g, respectively. The transition between the two regimes occurs for small g ′ = g − L/(2π), or equivalently for g ∼ log y/(2π).

Conclusions
In this paper, we demonstrated that the correlation function of four infinitely heavy half-BPS operators defined in (1.1) satisfies a system of nonlinear integro-differential equations in planar N = 4 SYM. These equations are powerful enough to determine the correlation function for arbitrary values of the 't Hooft coupling and for generic values of the cross ratios.
The starting point of our consideration was a representation of the correlation function as a determinant of a semi-infinite matrix previously derived in Refs. [10,11] using integrabilitybased hexagonalization framework [4][5][6]. The matrix in question describes magnons propagating on a worldsheet of the octagon in the dual string theory and it has interesting properties. We found that it can be brought by an appropriate similarity transformation to a block-diagonal form. This allowed us to express the correlation function as a Fredholm determinant of an integral operator on a half-line with a well-known Bessel kernel modified by a cut-off Dirac-Fermi-like function. The resulting Fredholm determinant representation of the four-point correlation function has a striking similarity to two-point functions in integrable low-dimensional integrable models. Taking advantage of this fact, we applied the method of differential equations developed in Ref. [15] and derived a system of nonlinear equations for the octagon.
At weak coupling, a solution to these equations yields the known expansion of the octagon over ladder integrals. At strong coupling, we developed a systematic expansion of the octagon in the inverse powers of the coupling constant and derived a representation of the corresponding coefficients as integrals of the cutoff function. We examined the resulting strong coupling expansion of the correlation function in various kinematical regions and observed a perfect agreement both with its expected asymptotic behavior dictated by the OPE as well as outcomes of numerical evaluation of the octagon. We found that, surprisingly enough, the strong coupling expansion is Borel summable. Applying the Borel-Padé method, we demonstrated that the improved strong coupling expansion correctly describes the correlation function over a wide region of the coupling constant.
There are several avenues for further studies. Analyzing the strong coupling expansion of the octagon, we focused on perturbative corrections in 1/g and systematically discarded exponentially small nonperturbative corrections suppressed by powers of e −πg . The latter are ubiquitous to strong coupling analyses in AdS/CFT and it would be interesting to elucidate their origin. Recall that in the wellstudied case of the cusp anomalous dimension, nonperturbative effects manifest themselves through Borel singularities of the strong coupling expansion in 1/g. They are associated with a formation of a mass gap in the nonlinear O(6) sigma model describing massless excitation of the string worldsheet [23,24]. Particularly, it would be important to reconcile exponentially small corrections to the octagon with the apparent fact of Borel summability of the corresponding perturbative series in 1/g. According to AdS/CFT, the leading term A 0 in the strong coupling expansion of the octagon (4.1) should be given by a minimal area of a string worldsheet residing on four AdS geodesics. It would be very interesting to compute A 0 directly from string theory. Quantum fluctuations of the string generate corrections to (4.1) suppressed by powers of 1/g. It would be even more exciting to reproduce the relations (5.23) and, in particular, unravel the origin of universality of log g enhanced term in (4.1) coming from quadratic fluctuations. The method of differential equations allowed us to compute all the coefficients in (4.1) except the constant, g−independent term B which appears as an integration constant in (2.21).
One can use the relation (2.27) to find its dependence on the kinematical variables but the calculation turns out to be surprisingly complicated. The main reason is that the integration in (2.27) does not commute with strong coupling expansion of the function Q(x). Namely, all terms in 1/g expansion of this function contribute to B and the strong coupling series needs to be resummed to all orders in 1/g. This question deserves a thorough investigation.
Intriguingly, the anomalous dimension Γ(g) describing the leading asymptotic behaviour of the null octagon (1.4) also controls a double-logarithmic behavior of the six-gluon MHV amplitude in a kinematical limit when three adjacent pairs of gluon momenta become collinear simultaneously [25]. We also observed in Ref. [1] that this anomalous dimension admits a complimentary description in terms of a flux-tube-like equation, whose origin remains obscure to us. These two facts hint towards the existence of a nontrivial connection between tessellation of amplitudes in terms of pentagons, on the one side, and hexagonalization of correlators, on the other.
It would be interesting to extend the formalism developed in this paper to more complicated correlation functions. The most immediate example are the octagons with nonzero internal bridges [2]. They admit a Fredholm determinant representation similar to the one studied above [10,11] and determine the so-called asymptotic correlation functions [2,26]. We plan to address these questions in the future.
Let us now reproduce the same relations using an analytical regulation. Assuming that the dominant contribution to (A.1) comes from z = O(g 0 ), we replace the Bessel function in (A.1) with its asymptotic behavior at large g The explicit expressions for a 0 and b 0 can be read from (5.11). Replacing rapidly oscillating functions in the expression for J 2 1 (2gz) with their mean values, we get from (A.1) (A.9) Expanding (A.8) in powers of 1/g, we encounter the same integral as in (5.20). We use the analytical regularization to define it according to (5.27). This leads to f n (g) = 1 2g It is straightforward to verify that this relation coincides with (A.6) for n ≥ 1. However, for n = 0 the relation (A.10) does not capture the leading, O(g 0 ) term in the right-hand side of (A.5). The missing contribution comes from the integration over the region z = O(1/g). Indeed, restricting the integration in (A.1) to 0 ≤ z < c/g, we get f 0 =ˆc /g 0 dz ∂ z log(1 − χ(z))J 2 1 (2gz) + · · · = 2ˆc 0 dx x J 2 1 (2x) + · · · = 1 + . . . , (A.11) where the dots denote contributions from z = O(g 0 ). Here in the second relation, we changed the integration variable to x = gz, took the limit g → ∞ together with c ≫ 1.
The above analysis demonstrates that the analytical regularization captures the contribution to (A.1) coming from integration over z = O(g 0 ) and, therefore, suppressed by powers of 1/g. The remaining O(g 0 ) contribution comes from z = O(1/g) and it is controlled by the behavior of the function q(z) at small z.

B Quantization condition from zeroth moment
We demonstrated in Section 5 that the relation (5.2), combined with the strong coupling expansion of the function (5.9) for z = O(g 0 ), can be used to compute the octagon (4.1) at strong coupling.
In this appendix, we present another approach to computing the octagon. It relies on the relation (4.15) and makes use of the strong coupling expansion of the function (5.9) for z = O(1/g). As we demonstrate below, it yields the same result for the expansion coefficients (5.23) and allows us to establish the relations (5.24) between the coefficients with odd and even indices.
We introduce x = 2gz and look for the solution to (4.12) in the form (5.9) with f = f 0 (x) + 1 g f 1 (x) + 1 g 2 f 2 (x) + . . . , (B.1) where x = O(g 0 ) and g ≫ 1. Substituting this ansatz into (4.12) and matching the coefficients in front of powers of 1/g we obtain the system of differential equations We verify that, in accord with (5.9), the solution to the first equation is Solving (B.2) we assume that the expansion (B.1) is uniform at small x, i.e. f 1 , f 2 , . . . do not modify the leading small x behavior and scale as f n ∼ x for x → 0.
At small x, we substitute f 0 (x) ∼ x/2 into the second relation in (B.2) to find that f 1 (x) ∼ A 2 x/2. In a similar manner, replacing f 2 (x) ∼ c x in the last relation of (B.2), we notice that the terms proportional to c drop out resulting to Equating to zero the coefficient in front of x on the left-hand side we arrive at A 3 + A 2 2 = 0, in agreement with (5.24).
The same mechanism is at work for high-order functions in (B.1). Replacing f ℓ (x) = k≥1 f ℓ,k x k in (B.1) we can obtain from (4.12) the system of linear equations for the coefficients f ℓ,k . For ℓ = 4, 6, . . . this system turn out to be overdetermined. For the solution to exist the coefficients with odd indices A 2k+1 have to satisfy consistency conditions. They take the form (5.24) and allow us to express A 2k+1 in terms of the coefficients with even indices. To determine even coefficients A 2 , A 4 , . . . , we use the relation for zeroth moment (4.15).

(B.6)
We show below that these functions are given by series in 1/g with the coefficients depending on the expansion coefficients A k of the octagon (4.1). Equating to zero the coefficients in front of powers of 1/g on the left-hand side of (B.5), we can determine A k . The leading function Q 0,0 can be computed using the relations (A.1) and (A.10) where I 1 is given by (5.27). According to (B.5), the O(1/g) contribution to Q 0,0 should cancel against O(g 0 ) contribution from Q 0,1 on the left-hand side of (B.5).
To compute Q 0,1 we need an expression for the functions f 1 (x). Solving the second equation in (B.2) we get where J 1 and Y 1 are Bessel functions. Here the integration contours are chosen in such a way that the two integrals are convergent. In general, f 1 is defined up to a solution to the homogenous equation, f 1 → f 1 + c 1 xJ 1 (x) + c ′ 1 xY 1 (x). The value of c ′ 1 = 0 is fixed from the requirement that f 1 = O(x) at small x. To fix c 1 , we require that the contribution of J 1 (x) to f 1 (x) should not modify asymptotic behavior of the leading function (B.3) at large x. This leads to c 1 = 0.
Moving on to the next order, the solution to (B.2) for f 2 is where x ′ = 2gz and the function χ(z) is defined in (2.13).
For the potential we find in a similar manner from (2.20) It follows from the second relation in (C.5) that χ m = O(g 2m+2 ) at weak coupling. As a consequence, the weak coupling expansion of (C.4) and (C.6) contains a finite number of terms at any order in the coupling. At finite coupling, we can approximate the exact expressions for Q(x) and u by retaining a sufficiently large number of terms in (C.4) and (C.6). To this end, we replace the semiinfinite matrix (k − ) nm in (C.4) and (C.6) by its finite-dimensional minor, n, m ≤ N max ∼ 10 2 . This procedure proves to be efficient in numerical analysis of the differential equation (2.24).