S-duality invariant perturbation theory improved by holography

We study anomalous dimensions of unprotected low twist operators in the four-dimensional $SU(N)$ $\mathcal{N}=4$ supersymmetric Yang-Mills theory. We construct a class of interpolating functions to approximate the dimensions of the leading twist operators for arbitrary gauge coupling $\tau$. The interpolating functions are consistent with previous results on the perturbation theory, holographic computation and full S-duality. We use our interpolating functions to test a recent conjecture by the $\mathcal{N}=4$ superconformal bootstrap that upper bounds on the dimensions are saturated at one of the duality-invariant points $\tau =i$ and $\tau =e^{i\pi /3}$. It turns out that our interpolating functions have maximum at $\tau =e^{i\pi /3}$, which are close to the conjectural values by the conformal bootstrap. In terms of the interpolating functions, we draw the image of conformal manifold in the space of the dimensions. We find that the image is almost a line despite the conformal manifold is two-dimensional. We also construct interpolating functions for the subleading twist operator and study level crossing phenomenon between the leading and subleading twist operators. Finally we study the dimension of the Konishi operator in the planar limit. We find that our interpolating functions match with numerical result obtained by Thermodynamic Bethe Ansatz very well. It turns out that analytic properties of the interpolating functions reflect an expectation on a radius of convergence of the perturbation theory.


Introduction
In the last two decades, there has been huge progress in understanding the four-dimensional N = 4 supersymmetric Yang-Mills theory (SYM). The N = 4 SYM is supposed to have many special properties such as superconformal symmetry [1,2,3,4,5], UV finiteness [5,6,7], S-duality [8,9,10,11,12], integrability [13,14], dual conformal symmetry [15] and so on. The N = 4 SYM also provides the canonical example of the AdS/CFT correspondence [16,17,18], where the N = 4 SYM is dual to type IIB superstring on AdS 5 × S 5 . While the AdS/CFT correspondence has stimulated our understanding on the N = 4 SYM and vice versa, most of the progress is based on weak coupling perturbation theory, planar limit or protected observables 1 .
Recently it has turned out that the conformal bootstrap approach [19,20,21,22] is very powerful tool also for the N = 4 SYM [23] (see also [28,29,30,31]). Indeed the N = 4 superconformal bootstrap [23] gives strong constraints on dimensions of unprotected leading twist 2 operators with various spins, which are SU(4) R singlets and non-BPS primary operators belonging to long representation of the N = 4 superconformal group. The leading twist operators in the N = 4 SYM at classical level are the so-called twist-two operators defined by where φ I is the adjoint scalar. The N = 4 superconformal bootstrap finds rigorous bounds on the dimensions within numerical errors, which are fully non-perturbative even for finite N and independent of the complex gauge coupling τ = θ 2π + 4πi The authors in [23] have also conjectured that the upper bounds on the dimensions are saturated at either τ = i or τ = e iπ/3 , which are the duality invariant points under S-transformation and (T · S)transformation, respectively.
Main purpose of this paper is to study the dimensions of the leading twist operators in the SU(N) N = 4 SYM by a completely different approach. We find approximate formula of the dimensions by resumming the perturbative data with help of the S-duality [8,9,10,11,12] and AdS/CFT correspondence [16,17,18]. We approximate the dimensions through interpolating functions which interpolates between two perturbative expansions around two different points in parameter space. The standard approach is to apply (two-point) Padé approximation, which is a rational function encoding the two expansions up to some orders. Recently Sen constructed another type of interpolating function, which has the form of a Fractional Power of Polynomial (FPP) [32]. A more general form of the interpolating function with the form of Fractional Powers of Rational function (FPR) has been constructed out by one of the current authors [33]. It has turned out that these interpolating functions usually provide better approximations than each perturbative expansion in intermediate regime of the parameter. See [32,33,34,35,36,37,38] for various applications 3 .
In this paper, we construct a class of interpolating functions to approximate the anomalous dimensions, which are consistent with known results on the perturbation theory, holographic computation and full S-duality SL (2, Z). In other words, our modular invariant interpolating functions reproduce the correct weak coupling expansion and large 't Hooft coupling limit in the planar limit. Such attempt was initiated in [41], which has constructed interpolating functions invariant under one particular element of SL(2, Z) either S or (T · S)transformations. After a while, Alday and Bissi constructed a class of interpolating functions [37], which is similar to FPP [32] but invariant under all elements of SL(2, Z). Here we construct a new class of interpolating functions by combining the ideas of FPR [33] and Alday-Bissi [37], and further imposing a consistency with the previous holographic results.
Our interpolating functions give predictions for arbitrary values of N and the complex gauge coupling τ . We compare our result with the recent results by the N = 4 superconformal bootstrap and test the conjecture that the upper bounds on the dimensions are saturated at one of the duality-invariant points τ = i and τ = e iπ/3 . As a conclusion we find that when we expect reasonable approximation by the interpolating functions, the interpolating functions have their maximal values at τ = e iπ/3 , close to the conjectural values of the N = 4 superconformal bootstrap.
In terms of the interpolating functions, we also study an image of conformal manifold in the space of the dimensions of the leading twist operators. We find that despite varying the coupling τ in the (real) two-dimensional region, the image is a very narrow line, which is almost one-dimensional. The narrow line is almost straight for N = 2 as in the result by different interpolating functions [41] while it is somewhat curved for N ≥ 3 contrary to [41].
We also construct interpolating functions for the subleading twist operator and study level crossing phenomenon between the leading and subleading twist operators. We use the terminology "level crossing" in the following two senses. As we increase the coupling, the dimensions of the leading and subleading operators approach each other with two possibilities.
Firstly, there is no operator mixing due to additional protected symmetries and their dimensions cross over. Secondly, operator mixing occurs and the dimensions of new eigenstates repel. We refer to the both as level crossing. We check that the interpolating functions for the leading and subleading twist operators with spin-0 do not cross each other for finite N. Namely the level crossing in the first sense does not occur for finite N. For large but finite N, we observe that the dimension of the leading twist operator becomes very close to the subleading one around τ = i and τ = e πi/3 . This implies that the level crossing in the second sense occurs for large but finite N.
We also study the dimension of the Konishi operator 4 in the planar limit. We construct interpolating functions consistent with the weak coupling expansion and holographic computation. We compare our result with numerical data obtained by Thermodynamic Bethe Ansatz (TBA) and find that our interpolating functions match very well with the TBA result. We also discuss that analytic property of the interpolating function reflects expectations on radius of convergence from the weak coupling perturbation theory.
This paper is organized as follows. In section 2 we briefly explain the previous results obtained by the weak coupling perturbation theory, holographic computation, and superconformal bootstrap. In section 3 we begin with introducing interpolating functions constructed previously. Then we construct our interpolating functions for the anomalous dimensions, which are consistent with the known results on the weak coupling perturbation theory, holographic results and full S-duality. Finally we discuss which of our interpolating functions would give the best approximation. Section 4 is the main section of this paper. We compare our results with the recent results by the N = 4 superconformal bootstrap. We also draw the image of the conformal manifold in the space of the dimensions of the leading twist operators and study the level crossing phenomenon between the leading and subleading twist operators for finite N. In section 5 we study the dimension of the Konishi operator in the planar limit. Section 6 is devoted to conclusion and discussions.

Previous results on leading twist operators
In this paper we mainly study the dimensions of the leading twist operators in the 4d SU(N) N = 4 SYM by using the technique of interpolating functions. Before introducing the interpolating functions, we review some relevant previous results know in the literature.
The leading twist operators under consideration are SU(4) R singlet and non-BPS primary operators belonging to long representation of the N = 4 superconformal group. At classical level, these operators are so-called twist-two operators: where φ I are the adjoint scalars in the N = 4 SYM and I is an index in the 6 of SU(4) R .
The leading twist operator has the dimension (2 + M) classically but receives quantum corrections: depending on the complex gauge coupling .
It is known that its weak coupling perturbative expansion is independent of θ while nonperturbative corrections depend on θ generically 5 . It is expected that the N = 4 SYM possesses the S-duality [8,9,10,11,12] described by the SL(2, Z) transformation which is a combination of S-and T-transformations: Note that there are two special values of τ : which are invariant under S-transformation and (T · S)-transformation, respectively. In this paper we assume SL(2, Z) invariance of the dimensions of the leading-twist operators: and construct the interpolating functions based on this assumption.

Weak coupling expansion
In perturbative regime, the leading twist operator is the twist-two operator O M (1.1), which is the Konishi operator especially for M = 0 (see e.g. [42,44]). The anomalous dimension of the twist-two operator has been computed up to four-loop for M = 0, 2 and three-loop for M = 4 in the weak coupling perturbation theory 6 [48,49,50,51,52,53,54]: There is a statement that the two-point function of Konishi operator does not receive instanton corrections [42]. Recently it is stated in [43] that the instanton correction starts at O(g 4 ) . 6 If it was limited to the planar limit, there are higher order computations [45,46,47].
These data will be used to construct our interpolating functions in subsequent sections.

Supergravity limit
The N = 4 SYM is expected to be dual to type IIB supergravity on AdS 5 × S 5 at large 't Hooft coupling λ = gN [16,17,18] in the planar limit. In the planar limit, the anomalous dimension of the twist-two operator O M (1.1) typically grows as ∼ λ 1/4 [17] and the twist-two operator is no longer leading twist operator for large 't Hooft coupling. The leading twist operator in the supergravity limit is a double trace operator which has a schematic form where tr(φ (i φ j) ) is symmetric traceless part of tr(φ i φ j ) and chiral primary operator belonging to 20 ′ representation 7 of SU(4) R . The double trace operator (2.8) is not protected in general but the large-N factorization implies that the dimension becomes the sum of the protected single trace operators and therefore protected in the planar limit. One can compute the anomalous dimension of the double trace operator (2.8) by the supergravity [55,56,57] and then the one of the leading twist operator in the supergravity limit is (2.9) Note that the first terms are easily understood by the large-N factorization.

The N = 4 superconformal bootstrap
In the past few years the N = 4 superconformal bootstrap approach [23] has obtained a relatively satisfying upper bounds on the dimensions of the unprotected leading twist operators by studying the four-point function (see also [28,29,30]) where O I 20 ′ is a superconformal primary scalar operator of dimension two in energy-momentum tensor multiplets transforming as 20 ′ representation in SU(4) R . The N = 4 superconformal symmetry allows us to describe the four-point function in terms of the N = 4 superconformal block [58,55,59,60].
In [23] the upper bounds on the dimensions of the leading twist operators with spin-0, 2 and 4, were obtained which are rigorous within numerical errors. Exclusion plots on the anomalous dimensions (γ 0 , γ 2 , γ 4 ) are presented in fig.1 of [23]. While the shape of not-excluded region roughly looks like a cube, its precise shape is complicated function of (γ 0 , γ 2 , γ 4 ). The "bound" values are listed in table 1 for each maximal value of (γ 0 , γ 2 , γ 4 ) in the not-excluded region.
The "bound" values in table 1 are somewhat conservative. This is because if actual values of (γ 0 , γ 2 ) were not equal to the bound values in table 1 (namely smaller than the bound values) for example, then possible value of γ 4 would generically be more strongly constrained. Hence we might have better estimates from fig.1 of [23] than the bound values in table 1. The authors in [23] have conjectured that this better estimate is given by the value at the corner of the cube-like region and that this is saturated by values of γ M at one of duality invariants points τ = τ S or τ = τ T S . This conjecture essentially claims the following two things: 1. The corner value obtained by the conformal bootstrap is saturated by the maximal value of γ M (τ, N).
2. The maximal value of γ M (τ, N) in the physical region of τ is given by τ = τ S or τ = τ T S .
The first point of this conjecture is closely related to whether the constraints from the conformal bootstrap is sufficiently strong or not. Namely, the upper bound of the bootstrap is greater than all the possible values of the anomalous dimension in general and may have a gap from the maximal value in principle. However, if the upper bound is maximally strong, then there is no such gap and the upper bound is the same as the maximal value though it is currently unclear if this is true. Regarding the second point, we do not know a priori which value of the coupling realizes the maximal value but it is natural to expect that such special thing happens in some special values in the τ -space, which are only the duality invariant points τ = τ S and τ = τ T S to our knowledge. Main purpose of this paper is to test the conjecture by using the interpolating functions.

Interpolating functions
In this section we introduce some classes of interpolating functions constructed in the literature [32,41,37,33] and then in the remaining part of this section we construct new class of interpolating functions for the anomalous dimensions of the leading twist operators. We impose the following conditions to the interpolating functions: 1. Real for Imτ ≥ 0.
2. Small-g expansion agrees with the weak coupling expansion of γ M (g) up to certain order.
4. Reproduce the holographic result in the planar limit at large 't Hooft coupling 8 .

Interpolating functions without S-duality (FPR)
Before considering the S-duality invariant interpolating functions, we introduce usual interpolating functions, which can be applied to problems without S-duality. Suppose that we would like to approximate a function F (g), which has the small-g expansion around g = 0 and large-g expansion around g = ∞ taking the forms F (g) = g a (s 0 + s 1 g + s 2 g 2 + · · · ) = g b (l 0 + l 1 g −1 + l 2 g −2 + · · · ). (3.1) The author in [33] constructed the following type of interpolating function for the function F (g): Here the coefficients c k and d k are determined such that power series expansions around g = 0 and g = ∞ agree with the ones of F (g) up to O(g a+m+1 ) and O(g b−n−1 ), respectively. By construction, the interpolating function reproduces both the small-g and large-g expansions of F (g). Since this interpolating function is described by Fractional Power of Rational function, we call this FPR. Note that we need

4)
8 Note that S-duality does not automatically imply the holographic matching and this condition is not redundant since the S-duality acts on g rather than λ = gN . To see this explicitly, let us consider a S-duality invariant quantity f (g, N ), with the 't Hooft expansion f (g, N ) = ∞ k=0 f 2k (λ)/N 2k . Then the S-duality implies . In the leading planar limit, the LHS has a contribution only from genus-0 while the RHS may receive all genus corrections. Thus the matching of the small-λ expansion does not imply the holographic matching in general. for m + n : odd , with ℓ ∈ Z. (3.5) If we take 2ℓ + 1 = a − b for a − b ∈ Z and m + n to be even, then this becomes the Padé approximant: while taking 2ℓ + 1 = m + n + 1 (2ℓ = m + n + 1) for even (odd) m + n gives the Fractional Power of Polynomial (FPP): recently constructed in [32]. In next subsection we will introduce interpolating functions invariant under the full S-duality inspired by the FPR. In sec. 5 we will use the FPR to study the dimension of the Konishi operator in the planar limit.

Modular invariant interpolating functions
Here we introduce interpolating functions, which are consistent with the weak coupling expansion (2.7) and full S-duality. Such an attempt was initiated in [41], where the author constructed interpolating functions invariant under one specific element of SL(2, Z) such as S-and (T · S)-transformations. Then Alday-Bissi constructed a class of interpolating functions, which are similar to FPP but invariant under all the elements of SL(2, Z), namely modular invariant interpolating functions [37]. Here we would like to have a new class of modular invariant interpolating functions, whose form is similar to FPR.

Alday-Bissi's interpolating function
Alday and Bissi constructed the following type of interpolating function [37] where the coefficient c k is determined such that expansion ofF 1 |m + nτ | 2s (Imτ ) s .
(3.9) 9 Note that s can be non-integer and E s (τ ) has a pole at s = 1. Hence we take s > 1.
where f np s (q) is the non-perturbative contribution containing powers of q = e 2πiτ (see app. A for details). Hence, we easily find that the expression inside of the bracket of F c (3.11) Thus an appropriate choice of c k correctly gives the weak coupling expansion of γ M (τ ). Since the interpolating function is similar to FPP, it is natural to consider FPR-like duality invariant interpolating functions as in next subsection.

FPR-like duality invariant interpolating function
We propose FPR-like generalization 10 of the Alday-Bissi's interpolating function: where we determine the coefficients c k and d k such that expansion ofF (s,α) m around g = 0 agrees 11 with the one of γ M (τ ) up to O(g m+1 ). Matching at O(g) leads us to Since the interpolating function is invariant under c k , d k → λc k , λd k , we can take without loss of generality. Imposing matching at other orders leads and hence we find (3.16) 10 We can also construct FPR-like generalization of the interpolating functions of [41], which is invariant under the particular elements of SL(2, Z), but we do not use it here. 11 Note that m should be m ≥ 2 since we need two coefficients at least for this interpolating function.
We also require p, q ∈ Z ≥1 , which implies for m : odd 1 2ℓ+1 for m : even , with ℓ ∈ Z. (3.17) Note that although the interpolating function (3.12) is inspired by FPR, this does not include the Alday-Bissi's interpolating function (3.8) as some special case. In Appendix B we also construct another type of FPR-like interpolating function invariant under the S-duality, which includes the Alday-Bissi's interpolating function as a special case. In next subsection we will further improve the interpolating functions of the type (3.12) by holography.

Further improvement by holographic computation
In previous subsection we have introduced the FPR-like interpolating functions consistent with the weak coupling expansion and full S-duality but not necessarily with the holographic result (2.9). Here we impose further consistency with the holographic computation. Let us consider which is formally the same as (3.12). However, we determine the coefficients c k and d k except d 1 such that expansion of F Without loss of generality, we can again take d q = 1. The remaining coefficient d 1 is determined as follows. Let us consider 't Hooft expansion of the interpolating function 12 : (3.20) Then we determine d 1 to satisfy where γ SUGRA M is the result in the supergravity limit given by (2.9). Imposing matching of other orders leads us to p + q − 2 = m, (3.22) and therefore we get (3.23) 12 Since we do not know f 4 (λ), we take f 4 (λ) = 0 for simplicity.
We also require p, q ∈ Z ≥1 , which constrains α as α = 1 2ℓ+1 for m : odd 1 2ℓ for m : even , with ℓ ∈ Z. (3.24) In this paper we apply the interpolating function (3.18) to approximate the dimensions of the leading twist operators. By construction, the interpolating functions should give good approximations around g = 0 for any (θ, N), its SL(2, Z) transformations and the supergravity limit. It is a priori unclear how nice the approximations are beyond these regimes.
In general this depends on details of the interpolating functions, which are specified by the parameters (m, s, α). Since we know information on the weak coupling expansions up to three or four loops, the numbers of possible (m, α) are finite but we have still infinite choices of s, which provide infinite choices of interpolating functions as well. We would like to know which (m, s, α) gives the best approximation or reduce the number of candidates. In next subsection we will discuss which interpolating function should give the best approximation by imposing some physical consistencies.
As we argued, we impose the constant behavior (2.9) to the interpolating functions in the large (λ, N) limit. One might wonder whether one can construct another modular invariant interpolating functions, which have the same weak coupling expansions but the different behaviours ∼ λ 1/4 in that regime as in the Konishi operator (1.1). This may not make sense physically since the dimension of (1.1) would not be modular invariant 13 but this may be useful in future for constructing interpolating functions for other quantities, which are modular invariants and have different behaviours in the classical string regime. In Appendix E, we try to construct a class of modular interpolating functions, with the same weak coupling expansion and λ 1/4 behaviour in the classical string limit.

Further constraints on interpolating function
In the previous subsection we have seen that we can construct enormous number of interpolating functions (3.18), which are consistent with the weak coupling expansions, S-duality and holographic results. This situation leads to "landscape problem of interpolating functions" as pointed out in [33]. Namely, it is a priori unclear which interpolating function gives the best approximation. In this subsection we discuss which value of (m, s, α) would give the best approximation. As a result, we will effectively find the best value of (m, s, α) for every spin. Because this subsection is not necessary to understand the main results of this paper, you can skip this subsection if you are interested only in the results.

Choice of m
By definition, our interpolating function F (s,α) m (τ ) reproduces the correct weak coupling expansion up to m-loop correction. In general the best value of m depends on details of problems and other parameters of interpolating functions. Probably most important point on this is convergence property of the weak coupling expansion. Namely, if the weak coupling expansion was convergent, then we should take m as large as possible, while if asymptotic, then we should be more careful.
Let us gain some intuitions from experiences on one-point Padé approximation of small parameter expansion. It is known that one-point Padé approximation including more terms often give more precise approximation even if the small parameter expansion is asymptotic. For instance, such behaviour appears in the series n (−1) n n!g n . A sufficient condition for convergence to exact result has been found in [62]. Similar results are obtained in two-points approximation by FPR analysis in 0d φ 4 theory, average plaquette in 4d pure SU(3) YM on lattice, and so on. When we do not know about properties of expansions sufficiently, we should conservatively choose m to be close to the optimized value 14 m * (g) in a range which we would like to approximate. So, independent of problems, when we would like to have better approximation in the range g ∈ [0, g * ] we expect that larger m gives more precise approximation until m ≃ m * (g * ).
In our problem, we expect that the weak coupling expansion is asymptotic and behaves as ∼ m! at m-loop for large m as in typical of field theory 15 . We do not know whether we should take m to be as large as possible or not. If this is the case, then we should take m to be our maximal value, namely m = 4 for spin-0, 2 and m = 3 for spin-4. If not, then we should think of optimization for the weak coupling expansion in the range g ∈ [0, 1] since S-duality relates this region to the other region. Ideally, we would like to know the optimized value m * (g) at g = 1 but our current information is not sufficient to estimate the optimized value. However, from many examples with factorial behaviour, we expect that m * (g = 1) is larger than 4. Thus we shall take m = 4 for spin-0, 2 and m = 3 for spin-4.

Constraints from weak coupling perturbation theory
The anomalous dimension γ M (g) has the small-g expansion where only positive integer powers of g appear. On the other hand, the perturbative part of the interpolating function takes the form whose small-g expansion contains fractional powers of g for general s. In order to guarantee absence of such fractional powers, we should take 2s ∈ Z. (3.27)

Constraints from 1/N expansion and holography
In the 't Hooft limit λ = gN = fixed, N ≫ 1, γ M has the following 1/N-expansion up to instanton corrections. We have the following two expectations for this expansion.
1. Since the N = 4 SYM has only adjoint fields, we do not have O(1/N 2ℓ+1 ) corrections.
2. Since the leading twist operators are dual to multi particle states appearing in the supergravity with G N ∼ 1/N 2 , large-λ expansion of a ℓ (λ) can be regarded as α ′expansion, where α ′ ∼ 1/ √ λ. Hence, we expect that the α ′ -expansion of a ℓ (λ) begins with some non-negative integer powers, namely O(α ′0 ) at lowest 16 .
These points can be used for constraining interpolating functions because the interpolating functions may not satisfy these conditions in general.
Indeed we find that the interpolating functions with s ∈ Z have odd powers of 1/N in the large N expansion while the interpolating functions with half-odd s do not have this problem. For example, the interpolating functions for the spin-0 operator with (m, α) = (4, 1/4) and s ∈ Z has O(1/N 2s+1 ) corrections. This means that we should take s to be as large as possible to make "wrong 1/N-corrections" as small as possible.
Regarding the second point, we find that interpolating functions with half-odd s have strange α ′ corrections while those with integer s are completely fine. For example, the interpolating function for the spin-0 operator with (m, α, s) = (4, 1/4, (2ℓ − 1)/2) has O(λ 2ℓ−1 ) in large-λ expansion of a ℓ (λ). Thus we should take s to be large as possible for 2s ∈ Z.
functions for all the cases become constant for large-s regime 17 . These behaviours are not only for this particular interpolating function but also for all other interpolating functions as long as we use interpolating functions of the type (3.18). See Appendix C for similar results on the other interpolating functions. Furthermore we can analytically show the saturation for large-s in weak coupling regime and at the duality invariant points τ = τ S and τ = τ T S . For details, see Appendix D. Thus we should pick s from the region having the plateau behaviour. These plots indicate that we can regard s = 30 as sufficiently large s.

Choice of α
The parameter α determines the type of branch cuts of the interpolating functions. In [38] it was discussed for the standard FPR (3.2) that correct values of α would be related to analytic properties of exact results. For example, if observables under consideration had square type of branch cuts, then interpolating functions with α = 1/2 would tend to be better approximations 18 . Since we do not know analytic properties of the anomalous dimensions, we do not know what should be the correct value of α from this viewpoint.
However, we now see that the upper bounds obtained by the bootstrap are useful to find "wrong" choices of α. In the large-N limit 19 , the upper bounds on the anomalous dimensions are [23,30] lim If a planar limit of an interpolating function breaks this bound considerably, then we can regard this interpolating function as the wrong choice. Let us consider the spin-0 and spin-2 cases. For these cases, the maximum of m is 20 m = 4. By the above arguments in this subsection, we expect that the best approximation among our interpolating function is either F (τ ) with sufficiently large s. As we discussed, we can regard s = 30 as sufficiently large s. Therefore F would give the best approximation (Their explicit forms are written in appendix F). In fig. 2 we plot the planar limits of the interpolating functions F (τ ), whose 18 Also note that interpolating functions with many poles may describe different type of branch cuts. For example, it is known that Padé approximant often describes branch cuts by bunch of poles. 19 The large-N limit taken in the context of the bootstrap so far seems g = fixed, N → ∞ rather than the planar limit. We expect that the results include the planar limit because the limit g = fixed, N → ∞ would be equivalent to very strong 't Hooft coupling limit for this case [68,69]. 20 For the spin-4 case, α of the interpolating function with m = 3 is uniquely determined as α = 1/3. (τ ) gives the best approximation and uses the interpolating function F 4 Results on the leading twist operators for finite N In this section we present our result on the leading twist operators and compare this with the N = 4 superconformal bootstrap. We also discuss the image of the conformal manifold in the space of the dimensions. Finally We study the dimension of the sub-leading twist operator and the level crossing phenomenon with the leading twist operator.

Comparison with the N = 4 superconformal bootstrap
In this subsection we compare our interpolating functions with the conjecture [23] by the conformal bootstrap that the upper bounds on the dimensions of the leading twist operators are saturated at one of the duality invariant points τ = τ S = i and τ = τ T S = e iπ/3 .

Spin-0
We begin with the spin-0 leading twist operator. By the arguments in sec. 3.4, we expect that the best approximation among our interpolating functions is F (τ ) for comparison with the conformal bootstrap, whose explicit forms are 21 Note that the planar limit is described only by the perturbative part. In order to compare our interpolating function with the N = 4 superconformal bootstrap, we shall ask where the interpolating function takes its maximal value as a function of τ . We expect that the maximal value is given at either of the duality invariant points τ = τ S or τ = τ T S if the interpolating function reasonably approximates the dimension. Note that the duality invariant points are also quite special for our interpolating functions. By construction our interpolating functions always have local extremum at τ = τ S and τ = τ T S because the building block E s (τ ) of the interpolating functions has local minimum 22 at these points for arbitrary s. What is nontrivial here is whether one of the extremum of the interpolating function at τ = τ S , τ T S is global maximum or not. We will see soon that the global maximum is given by 23 In fig. 3 we plot coupling dependence of the interpolating function F for θ = 0. 22 Global minimum of E s (τ ) is given by τ = τ T S . 23 Interpolating functions used in [37] also have global maximum at τ = τ T S while those used in [41] have global maximum at τ = τ S . (τ = 1/2+i/g) for the spin-0 operator is plotted as a function of g together with the 4-loop result. Note that g = 2/ √ 3 corresponds to τ = τ T S = e iπ/3 . Note 24 that tan −1 g = π/4 corresponds to the duality invariant point τ = τ S = i under the S-transformation. First we easily see that the interpolating function is consistent with the upper bounds for all the values of N. Next we observe that the interpolating function has the peak at g = 1, namely, τ = τ S . This indicates that τ = τ S gives the local maximum of the interpolating function. In fig. 4, we give similar plots for θ = π as fig. 3, whose right end g = 2/ √ 3 corresponds to the duality invariant point τ = τ T S . We again see that the interpolating functions have the local maximum at τ = τ T S .
Which of the duality invariant points does give the global maximum? In table 2, we explicitly write down the values of the interpolating functions at τ = τ S and τ = τ T S for various N. The table tells us that the interpolating function has the larger values at τ = τ T S than the one at τ = τ S . Actually we have checked that this is true for many other values of N. Thus we conclude that the interpolating function has the global maximum at τ = τ T S .
In fig. 5, we plot the values of our interpolating function at the duality invariant points and the data of the N = 4 superconformal bootstrap. The horizontal axis 25 is √ a, where a 24 One might wonder that the interpolating function has a cusp at τ = τ S . But we can prove analytically that the interpolating function is differentiable at τ = τ S . 25 Note that we have constructed the interpolating functions for SU (N ) gauge group and therefore the   (τ T S ), respectively. The red dashed line shows 2 − 16/N 2 , which is obtained by numerical fitting of the corner values in the large-a regime [23].
is the central charge 26 given by 27 Note that information on the gauge group in the N = 4 superconformal bootstrap is packaged into the central charge a.
As discussed above the interpolating function F (τ ) has the greater values at τ = τ T S than the one at τ = τ S . Thus "prediction" for the maximum of the dimension from our smallest value of a, which we can compare with the bootstrap, is a = 3/4 corresponding to the SU (2) case. 26 a is defined as From fig. 5 we see that our result is close to the corner values in the small-a regime and largea regime but there are about 10% discrepancies in the intermediate regime.
To interpret this, note that the accuracy of the interpolating function should depend on a (or equivalently N) for the following two reasons. First of all, we have imposed matching with the holographic computation in the supergravity limit: λ ≫ 1, N ≫ 1 and consistency with the upper bound (3.29) in the planar limit. Since the maximal value of the interpolating functions satisfying these conditions is two, our interpolating function matches with the corner values in the large-a regime almost 29 by construction. Second, recalling that the effective coupling constant is gN rather than g, we know that the weak coupling expansion is more precise for the small-a regime. Thus the intermediate region is harder to be approximated by our interpolating functions compared to the other regime. Hence, the interpolating function would give relatively worst approximation in the intermediate regime.
Thus we interpret the discrepancies as the lack of the accuracy of the interpolating function.
As a result, when we expect good approximations by the interpolating function, we have found that the interpolating function at τ = τ T S is very close to the corner values. This supports the conjecture of the N = 4 superconformal bootstrap [23] that the upper bounds on the dimensions are saturated at one of the duality-invariant points τ = τ S and τ = τ T S . We certainly expect that if we could include more higher order terms of the perturbation theory, then the interpolating function would have a better approach to the corner values. In the rest of this subsection we will see that similar results hold also for the spin-2 and spin-4 operators.

Spin-2
Next we consider the spin-2 leading twist operator. According to sec. 3.4, we expect that the best approximation is given by F (30,1/4) 4 (τ ) as in the spin-0 case. In fig. 6 and 7, we plot coupling dependence of the interpolating function for θ = 0 and π as in fig. 3 and 4, respectively. From these figures, we see that F   (τ ) for τ = i/g is plotted as the function of g.

SU(2)
SU ( Table 3: The interpolating function for spin-2 at the duality invariant points and data from the N = 4 superconformal bootstrap. In fig. 8, we compare the interpolating function at the duality invariant points with the N = 4 superconformal bootstrap. We easily see that the result of F (τ = τ T S ) is very close to the corner values in the whole region. This situation is different from the spin-0 case, where we have about 10% discrepancies in the intermediate regime. We interpret this as large-a regime being effectively broader for larger spin case. Indeed O(1/N 2 ) correction to the holographic result for spin-M behaves as −96/(M + 1)(M + 6) in the supergravity limit [55] and therefore the anomalous dimension with larger-M converges to γ = 2 faster with increasing N. Since our interpolating function correctly approximates the large-a regime, we expect that the interpolating function for the spin-2 case gives reliable approximation (1/2 + i/g) is plotted as a function of g.
in broader range of a compared to the spin-0 case. Thus, we expect that our interpolating function provides better approximation in the whole range for the spin-2 case. Since the interpolating function F (30,1/4) 4 (τ T S ) is close to the corner values, our result supports the conjecture in [23] also for the spin-2 case.

Spin-4
Finally let us consider the spin-4 case. For this case we know the weak coupling expansion only up to three loop. Hence according to sec. 3.4 we shall consider interpolating functions with m = 3. Contrary to the previous case, now α is uniquely determined as α = 1/3. Thus we expect that the best approximation is given by the interpolating function F (30,1/3) 3 (τ ). As in the spin-0 and spin-2 cases, we find from fig. 9 and 10 that the interpolating function F (τ ) has local maximum at the duality invariant points τ = τ S and τ = τ T S . From table 4, we see that the interpolating function at τ = τ T S is greater than the one at τ = τ S . Thus our prediction for the maximum anomalous dimension is given by 31      (τ ) is less than the spin-0 and spin-2 cases in the small-a region. This is because the result of the weak coupling expansion is available only up to three loop for this case 32 . Secondly, as discussed in the spin-2 case, large-a regime is effectively broader for the spin-4 case. We expect that the interpolating function is not reliable in the intermediate regime but reliable in the small-a and large-a regime.
As a conclusion of this subsection, we have seen that when we expect reasonable approximation by the interpolating functions, the maximal values of the interpolating functions are close to the corner values of the bootstrap. Thus we conclude that our interpolating function approach strongly supports the conjecture by the N = 4 superconformal bootstrap, which 32 It is worthwhile to note that the weak coupling expansion up to three loop is the same as the planar limit and the fourth loop is the first order to deviate from the planar limit. Thus we expect including the four-loop to be important. states that the upper bounds on the dimensions are saturated at one of the duality-invariant points τ = τ S or τ = τ T S . Obviously, if higher orders of the weak coupling perturbative series become available, then we can obtain more precise interpolating functions. It would be nice if one can obtain the higher order results.

Comments on other gauge groups
A priori our interpolating functions are valid only for the SU(N) gauge group since we have used the weak coupling expansion and holographic computation for the SU(N) case. However, both the data of the bootstrap and the interpolating functions looks continuous for a ≥ 3/4, whereas the bootstrap data has a cusp at a = 3/4 (SU(2) case) which is not reproduced by the interpolating functions. Thus we expect that our interpolating functions would reasonably approximate the dimension of the leading twist operator for other gauge groups as long as a ≥ 3/4. For a < 3/4, the only possible gauge group is U(1), which is Abelian. The bootstrap results have cusps at a = 3/4 and this implies a kind of transition from Abelian theory to non-Abelian theory, or from free theory to interacting theory. Contrary to the bootstrap, our interpolating function is smooth across a = 3/4. Presumably the difference comes from the fact that naive continuation of our interpolating function to general a is a continuation suited to interacting theory and therefore we cannot apply the (1/2 + i/g) for the spin-4 operator is plotted as a function of g.
interpolating function for the U(1) case.

Comments on instanton corrections
The small-g expansions of our interpolating functions contain exponentially suppressed corrections, whose weights are the same as the instanton actions. This feature technically comes from imposing the full S-duality. Recently Alday-Korchemsky [43] computed instanton corrections to the dimension of the Konishi operator at O(g 2 ) by expanding around the instanton configuration. They have found that the instanton corrections start from O(g 2 ) and the one-instanton correction for the SU(2) case is given by They also computed n-instanton correction in the large-N limit 33 as − 27g 2 10π 5/2 n 3/2 N 3/2 e 2πinτ + e −2πinτ d|n 1 d 2 , (4.6) Figure 11: A similar plot for the spin-4 case as fig. 5 and 8. The red dashed line shows 2 − 48/25N 2 , which is numerical fitting of the corner values in the large-a regime [23].
by using the technique of [70]. On the other hand, the instanton corrections of our interpolating functions with general parameters start from 34 O(g s+1+min(p,q) ). Since we would like to take s to be sufficiently large as discussed in sec. 3.4.2, the interpolating function F However, it is worthwhile to note that there is a subtlety in [43]. We can also compute instanton corrections to the circular Wilson loop by using the same technique as in [71]. However, the result of [71] states that there are non-trivial instanton corrections to the circular Wilson loop. This does not agrees with the results obtained by summing ladder diagrams [72,73] and the localization method [74], where instanton corrections are trivial. Thus we should be careful on this point.

Comments on higher order corrections in the planar limit
Although we have used the four-loop result of the weak coupling expansion to construct the spin-0 interpolating functions, there is a seven-loop result in the planar limit [45,46,47], whose explicit form is given by (5.2). We did not use the seven-loop result because we also need non-planar corrections to completely fix the coefficients in interpolating functions and the interpolating functions strongly depends on the values of the non-planar higher order corrections for small N. Here we just compare the higher order correction in the planar limit with the ones of the interpolating function. The coefficients of the higher order small-λ expansion of the interpolating function F

Image of conformal manifold
The complex coupling τ is the exactly marginal parameter and hence the coordinate of the conformal manifold in the N = 4 SYM. At every point τ on the conformal manifold, we have a set of dimensions of the leading twist operators: (∆ 0 , ∆ 2 , ∆ 4 , · · · ). Since we have constructed the approximations of (∆ 0 , ∆ 2 , ∆ 4 ) by the interpolating functions, we can draw an image of the conformal manifold projected to the (∆ 0 , ∆ 2 , ∆ 4 )-space.

SU(2) case
In fig. 12 [Left-Top] we plot the image of the conformal manifold for the SU(2) case projected 35 to the (γ 0 , γ 2 , γ 4 )-space. We have sampled the interpolating functions for θ = πn/5 with n = 0, 1, · · · , 5 and many values of g. The figure shows that the image is almost one continuous straight line and this is different from the following naive expectation. Since we are considering the image of the two dimensional conformal manifold, one naively expects that the image was almost a straight line around the origin but it starts spreading as the dimension increases. This is because the anomalous dimensions in the weak coupling regime are almost independent of θ but the θ-dependence becomes important in the strongly coupled regime, where the anomalous dimensions become large. Thus, if the above naive expectation was correct, then we should observe six distinguishable lines as we consider the six different values of θ.
This property has been already observed for the (γ 0 , γ 2 )-plane in [41] by using the different interpolating functions partially invariant under the SL(2, Z) transformations. In [41] the authors found that the slope of the almost straight line in the (γ 0 , γ 2 )-plane is very close to 25/18 ≃ 1.38889, which is the same as the ratio between the one-loop anomalous dimensions: 18 . whose slopes are the same as We find that the images projected to the (γ 0 , γ 4 ) and (γ 2 , γ 4 )-planes are very close to the straight lines as well. This implies that the observation in [41] is also true for the (γ 0 , γ 4 ) and (γ 2 , γ 4 )-planes. In the rest of this subsection we see that the situation is different for the higher N cases.

Higher N
Let us consider the higher N cases. In fig. 13 we give the similar plots as fig. 12 for the SU(3) case. We now find both similarity and dissimilarity from the SU(2) case. The similarity is that the image is still very narrow and therefore looks like one-dimensional. The dissimilarity is that the image is no longer straight line, namely it is curved in the strong coupled regime.
To see this in more detail, let us see fig. 13  the matching with the straight line holds even in the regime, where we cannot naively trust the one-loop approximation for the dimensions. Fig. 13 [Left-Bottom] and 13 [Right-Bottom] imply that similar results hold also on the (γ 0 , γ 4 ) and (γ 2 , γ 4 )-planes, respectively. Fig. 14 shows the similar plots as fig. 12 and 13 for the SU(4) case. We easily see that the result for N = 4 is qualitatively the same as the N = 3 case though the deviations from the straight lines are slightly larger. We have checked that similar results hold for other values of N. The above results for N ≥ 3 are different from those in [41], where the images were still the straight lines even for N = 3 and N = 4 while the common feature is that the images are very narrow. As a conclusion, all the interpolating functions constructed so far give very narrow lines for the image of the conformal manifold and we expect that this property is probably true also for the exact results. However, we do not have a definite conclusion on whether the narrow lines are straight or curved though it would be natural to be curved in the strongly coupled regime.

Level crossing
In this subsection we compare the dimension of the leading twist operator with the one of the subleading leading operator and study level crossing phenomenon between the leading and subleading twist operators. We use the word "level crossing" in the following senses. As we increase the coupling, the dimensions of the leading and subleading operators approach each other and the following two things may occur: 1. Operator mixing between them do not occur by some additional protected symmetries and their dimensions cross over.

Operator mixing occurs and the dimensions of new eigenstates repel.
We refer to both of the above as level crossing.
In the N = 4 SYM, the operator Trφ I D M φ I in (1.1) has the dimension (2 + M) at the classical level while this becomes very large for large 't Hooft coupling in the planar limit. This is because this operator is dual to a massive string state and therefore its dimension behaves as ∼ 1/ √ α ′ ∼ λ 1/4 [17]. This implies that the operator (1.1) is no longer the leading twist operator in the large-λ regime of the planar limit. Indeed we have a family of double trace operators with the same spins, whose dimensions are protected in the planar limit. The double trace operators consist of protected chiral primary operators and their dimensions are independent of λ in the planar limit because of large-N factorization. Thus we expect the level crossing in the first sense between the leading and subleading twist operators in the planar limit.
To interpret this, let us recall the Wigner-von Neumann non-crossing rule known in quantum mechanics, which states that levels of states with the same symmetry cannot cross each other. Since the dilatation operator in the N = 4 SYM on R 4 corresponds to the Hamiltonian on R × S 3 , we expect that the dimensions obey the Wigner-von Neumann noncrossing rule. Recently it was discussed [75] that in the problem of the N = 4 SYM , 1/N plays roles as "interaction energy" of two level system in quantum mechanics. The actual crossing between the Konishi and double trace operators in the planar limit is consistent with the Wigner-von Neumann non-crossing rule if we have additional symmetry in the planar limit. Most promising candidate for such symmetry is the one associated with the integrability, which is supposed to appear in the planar limit. Since we do not expect such additional symmetry beyond the planar limit, we expect that the level crossing in the first sense does not occur for finite N but the one in the second sense occurs for large but finite N.
Here we approach the level crossing problem by using our interpolating function for the spin-0 case 36 . The operator (1.1) for M = 0 is nothing but the Konishi operator Trφ I φ I , whose dimension is two at the classical level and 2(4πλ) 1/4 + (corrections) for large-λ in the planar limit. Next let us also consider the following operators which have naively the same symmetry as the Konishi operator. Particular linear combinations of these operators are eigenvectors of the dilatation operator, which are dimension 4 in the weak coupling limit. Therefore the lowest dimension among those is the dimension of the subleading twist operators in the weak coupling regime. The dimensions of the operators (4.11) at one loop is given by [76,77,37] ∆ sub (τ ) = 4 + Nw(N) 2π g + O(g 2 ), (4.12) where w(N) is roots of the equation The dimension of the subleading twist operator is described by the smallest root w − (N) of this equation, which is always negative 37 . Now we construct interpolating functions for the subleading twist operator by imposing a match with the weak coupling expansion (4.12) and the SL(2, Z) duality. Since we know the perturbative expansion up to only one-loop, we have only one coefficient to be tuned in the S-duality invariant interpolating functions. In this situation, we cannot construct the FPR-like interpolating function (3.12) since we needs m ≥ 2. Hence we use the Alday-Bissi type interpolating function (3.8) for the dimension of the subleading twist operator. For general s, we can easily construct the Alday-Bissi type interpolating function as (4.14) We do have infinitely many s-dependent S-duality invariant interpolating functions. Which value of s is most appropriate? For this purpose, let us consider the small-g expansion of the interpolating function: which is true up to non-perturbative corrections. From this expression we can easily see that after the one-loop correction, the next term is a O(g 2s ) correction and this fact is useful in constraining s. First of all, we need 2s ∈ Z to get integer powers of g in the weak coupling expansion. Next when s is too large, we have large jumps of powers in the weak coupling expansion and therefore want a small value of s as possible. Since the interpolating function is well-defined for s > 1, we conclude that the most appropriate value of s is s = 3/2. In fig. 15 we compare 38 the interpolating function ∆  for θ = 0. We easily see from the figure that the interpolating functions do not cross each other for all the values of N. We have checked that this is true for other values of N, which are not present here. In fig. 16 we give similar plots for θ = π as in fig. 15 in order to test that the above is true also for different values of θ. Again the interpolating functions do not cross each other for all N. For SU (12), we observe that the dimension of the leading twist operator becomes very close to the subleading one around τ = τ S and τ = τ T S . Indeed we have found similar results for larger N. Thus we conclude that the level crossing in the first sense (actual crossing) does not occur for general finite N but the one in the second sense (small operator mixing) occurs for large but finite N.
Interestingly the interpolating function for the subleading twist operator has minimum at the duality invariant points. It would be interesting to find any physical interpretations for that. Finally we have not studied the problem in this subsection for non-zero spin cases. This is because there are no available results for the weak coupling expansions of the subleading twist operators but the one-loop computations for non-zero spins should not be hard. It is nice if one can perform the one-loop computations and the same analysis as in this subsection.

Results on Konishi operator in the planar limit
In this section we analyze the dimension of the Konishi operator in the planar limit. While the Konishi operator is the leading twist operator in the weak coupling regime, this has very large dimension in the classical string regime. We approximate the dimension of the Konishi operator in the planar limit in terms of the standard FPR described in sec. 3.1.

Previous results
First we briefly review previous results on the dimension of the Konishi operator in the planar limit: In this limit there is a 7-loop computation in the weak coupling expansion [45]:  There are also some holographic computation of the Konishi operator 39 : There are some numerical computations using the Thermodynamic Bethe Ansatz (TBA) [78,79] and we will compare this with our results.

Comparison with Thermodynamic Bethe Ansatz
By using the weak and strong coupling expansions, we can construct the interpolating functions. Here we do not use the modular invariant interpolating function (3.11) but use the standard FPR (3.2), which is consistent with the weak coupling expansion (5.2) and the holographic result (5.3) since the dimension of the Konishi operator is not expected to be modular invariant. For this purpose, note that the weak coupling expansion of ∆ Konishi (λ) is the power series expansion of λ while the holographic result is the one of λ −1/4 . However, if we work with ∆ Konishi (λ) + 2, then the holographic result becomes the power series expansion 40 of λ −1/2 . Thus we consider the FPR-type interpolating functions for ∆ Konishi (λ) + 2 rather than ∆ Konishi (λ). More precisely we rewrite the power series expansions of ∆ Konishi (λ) + 2 in terms of x = √ λ instead of λ and consider the interpolating functions F (α) m,n (x) for the expansions 41 . Then the interpolating function approximates the dimension ∆ Konishi (λ) by One of subtleties here is that we can construct enormous number of interpolating functions as in the leading twist operators and appropriate choice of (m, n, α) is a priori unclear. However, we have expectations on the appropriate choice at least for (m, n). First for m, since the weak coupling expansion in the planar limit is expected to be convergent [80], we expect that interpolating functions with larger m give better approximations. For n, we have very few choices, namely n = 0, 1, 2. Although the large-λ expansion would be asymptotic and its optimized order is unclear, it is natural to expect that the optimized order is larger than n = 2. Thus we should take (m, n) to be large as possible. It is still unclear what is an appropriate value of α. However, this is not so problematic because the interpolating functions with (m, n) = (15,2) are weakly dependent on α as we will see shortly.
In fig. 17 [Left] we compare our result with numerical result obtained by Thermodynamic Bethe Ansatz (TBA) [78]. We also draw the one-point Padé approximant of the weak coupling expansion defined as where the coefficients are determined by the correct reproduction of the small-g expansion up to O(g m+n ). The one-point Padé approximant in fig. 17 [Left] is the so-called diagonal Padé with m = n = 3. We easily see that our interpolating function agrees with the TBA result in whole region of λ. The one-point Padé approximation is good up to around λ = 5 but deviates from the TBA result in stronger coupling region. In fig. 17  . 40 We expect this property also for higher orders since α ′ ∼ λ −1/2 . 41 Note that because of this parametrization, the interpolating functions may have half-odd power of λ in their small-λ expansions. relative errors of the approximations by the interpolating functions from the TBA result to study precision of their approximations. We find that all the interpolating functions have errors less than 1% in the whole region.

Analytic property
Here we study analytic properties of the interpolating functions following the spirit of [38]. It has been proposed in [38] that when interpolating functions are to give good approximations along an axis in parameter space, then their branch cuts have the following interpretations: 1. The branch cuts are particular to the FPR and the artifact of the approximation by the interpolating function. Namely, this type of branch cut is not helpful in extracting any physical information.
2. The physical quantity, which we try to approximate by the FPR, has an actual branch cut near the branch cut of the FPR. Namely, the branch cut of the FPR approximates the "correct" branch cut of the physical quantity.
3. There is an anti-Stokes line near the branch cut across which perturbation series of the physical quantity changes its dominant part.
Here we expect that some zeros and poles of the interpolating functions approximate analytic properties of the dimension of the Konishi operator.
In fig. 18 we summarize the analytic property 42 of the FPR with F 1/32 15,2 (x) as a complex function of λ. The symbols denote ends of branch cuts of the interpolating function, which are characterized by zeros and poles 43 of the associated rational function (F (α) m,n ) 1/α . The green dashed line denotes |λ| = π/4, which is the expected radius of convergence in the weak coupling expansion from previous works [81,82,47]. From this figure, we observe that many ends of branch cuts are located around the circle |λ| = π/4. This is indeed reasonable because one expects that the dimension of Konishi operator has some singularities around radius of convergence of the weak coupling perturbative expansion and so does the interpolating function if the interpolating function approximates the proper analytic properties.

Conclusion and discussions
In this paper we have mainly studied the dimensions of the unprotected leading twist operators in the 4d SU(N) N = 4 SYM. We have constructed the class of interpolating functions (3.18) to approximate the dimensions. The interpolating functions are consistent with the previous results on the perturbation theory (2.7), holographic computation (2.9) and full S-duality. and give the predictions for arbitrary value of N and the complex gauge coupling τ in the fundamental region, which particularly includes the duality-invariant points τ = i and τ = e iπ/3 . We have used our interpolating functions to test the recent conjecture by the N = 4 superconformal bootstrap [23], which states that the upper bounds on the dimensions are saturated at either one of the duality-invariant points τ = i or τ = e iπ/3 . It has turned out that our interpolating functions have the maximum at τ = e iπ/3 . In the regime where we expect reasonable approximations by the interpolating functions, the maximal values are close to the conjectural values in [23]. Thus we conclude that our interpolating function approach strongly supports the conjecture of [23] as well as [41,37]. In order to construct the interpolating functions we have used the available four-loop or three-loop results of the weak coupling expansions. Obviously, if higher orders of the weak coupling perturbative series become available, then we can obtain more precise interpolating functions which are expected to better approach the corner value at τ = e iπ/3 . It would be nice if one can obtain the higher order results and repeat our analysis in this paper to construct better interpolating functions.
In terms of the interpolating functions, we have drawn the image of the conformal manifold in the space of the dimensions (∆ 0 , ∆ 2 , ∆ 4 ). We have found that the image is almost a line as in [41] despite the conformal manifold being two-dimensional. For the SU(2) case, the line is almost straight, whose slope is the ratio of the one-loop anomalous dimensions as in [41]. For higher N, we have seen that the line is curved contrast to [41].
We have also constructed interpolating functions for the subleading twist operator and 42 We define locations of branch cuts as (F (α) m,n ) 1/α ∈ (−∞, 0). The ends of branch cuts are given by infinity, zeros or poles of (F (α) m,n ) 1/α . 43 The poles of (F studied the level crossing phenomenon between the leading and subleading twist operators. We have checked that the interpolating functions for the both operators do not have actual crossing with each other for finite N. For large but finite N, we have found that the dimension of the leading twist operator becomes very close to the subleading one around τ = i and τ = e πi/3 . This implies the small mixing between the two operators. To construct the interpolating function for the sub-leading twist operator, we have used only the one-loop result. Computing two-loop order would give more insights to the level crossing problem. In this paper we have considered only the spin-0 case as there are no one-loop computations for the subleading twist operators with non-zero spins. It would be nice if one can perform the one-loop computations and repeat the same analysis for non-zero spin cases.
We have also studied the dimension of Konishi operator in the planar limit. We have found that our interpolating functions match with the numerical result obtained by Thermodynamic Bethe Ansatz very well. Furthermore we have discussed the analytic property of the relatively best interpolating function in the spirit of [38]. It has turned out that analytic property of the interpolating function reflects the expectations on radius of convergence from the weak coupling perturbation theory.
The key to our interpolating functions is their modular invariance. It would be illuminating if we study other modular invariant observables by our interpolating functions. More challenging direction is to construct interpolating functions for modular forms, which is not modular invariant but have particular transformation properties under SL(2, Z) transformations. It would be also interesting to consider other theories, which enjoy SL(2, Z) duality. Indeed many theories with the SL(2, Z) duality were recently found by torus compactifications of 6d (1, 0) theories have the S-duality [83]. e.g. sec.5.3 of [84]) where σ s (k) is the divisor function In terms of (g, θ), this is written as This representation is suitable for numerical computation. When the summation does not converge well, we practically compute the summation at another point connected by the SL(2, Z) symmetry.

B Another FPR-like modular invariant interpolating function including Alday-Bissi's one
We can also construct the following interpolating function, which is inspired by the FPR and a generalization of the Alday-Bissi's interpolating function but a different form: We determine the coefficients c k and d k such that expansion ofF (s,α) m around g = 0 agrees with the one of γ M (g) up to O(g m+1 ). For q = 0, this is nothing but the Alday-Bissi's interpolating function and we consider q = 0 case below.
Matching at O(g) leads us to Imposing matching of other orders leads Therefore we get We also require p, q ∈ Z ≥1 , (B

D Saturation of interpolating function for weak coupling and at duality invariant points for large-s
The coefficients in the interpolating function are determined only by the perturbative part of the Eisenstein series: In particular for large s, the coefficients are determined by the first term. Since ζ(2s) is almost unity for sufficiently large s, the interpolating functions for large s is approximately given by and the coefficients c k and d k are almost independent of s for very large s. Because of this, the interpolating function in weak coupling regime is almost independent of s for large-s. For values of the interpolating functions at the duality invariant points and for large s, we can further solidify our independent of s claim. For this purpose, we should know values of E s (i) and E s (e iπ 3 ) for large s. By using ζ a (s) ≃ a −s with a ≥ 0 for large s, we find (τ ) for the spin-0 leading twist operator at randomly chosen five points τ = (r 1 , · · · , r 5 ). (r 1 : blue circle, r 2 : red square, r 3 : green triangle, r 4 : black inverse triangle, r 5 : purple asterisk) [Left-Top] SU(2) case, (r 1 , r 2 , r 3 , r 4 , r 5 ) ≃ (0.3377 + 0.4745i, 0.8231 + 0.9956i, 0.5698 + 0.5929i, 0.8193 + 0.7029i, 0.7449 + 0.34278i). [Right-Top] SU(3) case, (r 1 , r 2 , r 3 , r 4 , r 5 ) ≃ (0.1185 + 0.8770i, 0.1236 + 0.3748i, 0.09718 + 0.6059i, 0.6656 + 0.8062i, 0.015689 + 0.9403i).
Therefore the interpolating functions at the duality invariant points are given by Thus the interpolating function at the duality fixed points is independent of s for large-s. E S-duality interpolating functions with λ 1/4 in the classical string limit In the main text, we have approximated the dimensions of the leading twist operators by the interpolating functions, which are consistent with the weak coupling expansions, holographic results and full S-duality. As the holographic results, we have used the results (2.9) for the double trace operators, which are the subleading twist operators in the weak coupling regime. This is because the operator (1.1) is dual to the massive string state and acquires the very large dimension in the classical string regime.
However, one may wonder if one can construct another modular invariant interpolating functions, which have the same weak coupling expansions but different behaviours say as ∼ λ 1/4 in the classical string regime. This may not make sense physically since the dimension of (1.1) would not be modular invariant but this may be useful in future for constructing interpolating functions for other modular invariants with different behaviours in the classical string regime. In this appendix, we try to construct a class of modular interpolating functions, with the same weak coupling expansion and λ 1/4 behaviour in the classical string limit.
To be specific, let us consider the Konishi operator. The dimensions of the Konishi operator behaves as (2.7) for weak coupling and as (5.3) for large-λ in the planar limit. Now we would like to construct interpolating functions, which are consistent with (2.7), (5.3) and the full S-duality. As a conclusion, we failed to construct a single interpolating function satisfying these properties but we find that a linear combination of multiple interpolating functions which are slight modifications of (3.18), satisfies the properties.
First we discuss that the interpolating functions of the type (3.18), which have been used for the leading twist operators in the main text, cannot satisfy the above properties. Indeed we have explicitly checked this for various cases and we can also show this for large-s in the following way. Recall that the coefficients c k and d k in (3.18) are determined in terms of only perturbative part of the Eisenstein series E s (τ ). Although the perturbative part of E s (τ ) has O(g −s ) and O(g s−1 ) parts, only the O(g −s ) part is relevant to determine the coefficients (E.1) While c k and d k are nontrivial function of N in general, we know that their planar limits behave as O(N q−k ) 44 since the anomalous dimension is O(1) in the planar limit. Hence in the planar limit, the interpolating function becomes Since this function becomes O(1) for large-λ, the interpolating function (3.18) for large-s cannot have the λ 1/4 -law in the classical string limit. Thus we shall consider different types of interpolating functions. Alternatively let us consider the following type of interpolating functions where c 1 = O(N q−t−1 ) and the only difference from (3.18) is the presence of the new parameter t in the first term of the numerator. Note that t should be integer to get weak coupling expansion with only integer powers of g. By a similar argument as above, we find that the planar limit of this interpolating function for large-s is given by (τ ) with large-s, then the λ 1/4 -law uniquely 45 requires (α, t) = (1/4, −1) and therefore the subleading order is O(λ −3/4 ) rather than O(1). Thus we cannot single interpolating function consistent with (2.7), (5.3) and the full S-duality by the type (E.4). Alternatively we find that the following linear combination of I (τ ) respectively matching to the O(λ 1/4 ), O(1) and O(λ −1/4 ) coefficients. We could also work with variations of the above scheme with more terms from the large-λ expansion in the planar limit and weak coupling expansion though we do not explicitly write their constructions.

F Explicit forms of interpolating functions
In this appendix we present explicit forms of the interpolating functions used in the main text. It should also be noted that we have to solve linear equations to solve for the unknown coefficients and thus it does not involve numerical approximation. Although we often write their coefficients with 6 digits of precision to avoid too long expressions, we practically use analytic expressions or infinite digits in Mathematica files.