N = 2* Phase Transitions and Holography

TL;DR

This paper explores the connection between supergravity probe analysis and large-N solutions in N=2* super-Yang-Mills theory, revealing that quantum phase transitions appear perturbatively at strong coupling.

Contribution

It introduces a supergravity-inspired formalism to solve the localization matrix model at strong coupling, linking phase transitions to perturbative phenomena.

Findings

Quantum phase transitions are visible at third order in strong-coupling expansion.

The formalism systematically solves the matrix model at strong coupling.

Phase transitions are perturbative on the string worldsheet.

Abstract

We clarify the relationship between probe analysis of the supergravity dual and the large-N solution of the localization matrix model for the planar N = 2* super-Yang-Mills theory. A formalism inspired by supergravity allows us to systematically solve the matrix model at strong coupling. Quite surprisingly, we find that quantum phase transitions, known to occur in the N = 2* theory, start to be visible at the third order of the strong-coupling expansion and thus constitute a perturbative phenomenon on the string worldsheet.

Equations

$$\Phi =\mathrm{diag}(a_1-u/N,a_2-u/N,....,a_N-u/N,u),\sum _i{a}_i=0,$$

$$\tau =\frac{4\pi i}{g_{\mathrm{YM}}^2}+\frac{\theta }{2\pi }=\frac{{\partial }^2\mathcal{F}}{\partial u^2},$$

$$\mathcal{F}(u)$$

$$\tau =\mathrm{e}{}^{-\varphi }+C_0,$$

$${\tau }_{\mathrm{sugra}}=\frac{4\pi Ni}{\lambda }\frac{u}{\sqrt{u^2-{\mu }^2}}.$$

$$\tau =\frac{4\pi Ni}{\lambda }+\frac{i}{2\pi }\sum _i\ln \frac{(u-a_i{)}^2}{|(u-a_i{)}^2-M^2|}.$$

$$\tau \approx \frac{4\pi Ni}{\lambda }+\frac{i}{2\pi }\sum _i\frac{M^2}{(u-a_i{)}^2}=\frac{4\pi Ni}{\lambda }+\frac{M^{2}Ni}{2\pi }{\int }_{-\mu }^{\mu }\frac{dx\rho (x)}{(u-x{)}^2}.$$

$$M^2{\rho }^{\prime }(x)\approx -\frac{8\pi }{\lambda }\frac{x}{\sqrt{{\mu }^2-x^2}}.$$

$$\rho (x)\approx \frac{8\pi }{\lambda M^2}\sqrt{{\mu }^2-x^2},$$

$$\mu =\frac{\sqrt{\lambda }M}{2\pi },$$

$${\tau }^{\prime }(u)=\frac{Ni}{2\pi }{\int }_{-\mu }^{\mu }dx\rho (x)\left(\frac{2}{u-x}-\frac{1}{u-x+M}-\frac{1}{u-x-M}\right).$$

$${\tau }_{\mathrm{sugra}}^{\prime }(u)=-\frac{4\pi Ni}{\lambda }\frac{{\mu }^2}{(u^2-{\mu }^2{)}^{3/2}}.$$

$$2\rho (x)-\rho (x+M)-\rho (x-M)=\frac{8\pi }{\lambda }\frac{{\mu }^2}{({\mu }^2-x^2{)}^{3/2}}.$$

$$Z=\int d^{N-1}a\prod _{i<j}\frac{(a_i-a_j{)}^{2}H(a_i-a_j)}{H(a_i-a_j+M)H(a_i-a_j-M)}\mathrm{e}{}^{-\frac{8{\pi }^{2}N}{\lambda }{\sum }_i{a}_i^2}{\left|{\mathcal{Z}}_{\mathrm{inst}}\right|}^2,$$

$$\Phi =\mathrm{diag}(a_1\dots a_N)$$

$$H(x)=\prod _{n=1}^{\infty }{\left(1+\frac{x^2}{n^2}\right)}^n\mathrm{e}{}^{-\frac{x^2}{n}}$$

$${\text{strokedint}}_{-\mu }^{\mu }dy\rho (y)\left(\frac{2}{x-y}-\frac{1}{x-y+M}-\frac{1}{x-y-M}\right)=0,$$

$${\int }_{-\mu }^{\mu }dx\rho (x)\ln \frac{\left|M^2-x^2\right|}{x^2}=\frac{8{\pi }^2}{\lambda }$$

$${\int }_{-\mu }^{\mu }dx\rho (x)=1.$$

$$\rho (x)\approx \frac{2}{\pi {\mu }^2}\sqrt{{\mu }^2-x^2},$$

$$\mu =\frac{\sqrt{\lambda }M}{2\pi }\left(\lambda \to \infty \right).$$

$$n=\left[\frac{2\mu }{M}\right],\Delta =\left\{\frac{2\mu }{M}\right\},2\mu =(n+\Delta )M,$$

$$r(u)=\frac{i}{2\pi }{\int }_{-\mu }^{\mu }dy\rho (y)\left(\frac{2}{u-y}-\frac{1}{u-y+M}-\frac{1}{u-y-M}\right),$$

$$\mathrm{Cont}r(x)=0,x\in (-\mu ,\mu ).$$

$$2\rho (x)-\rho (x-M)-\rho (x+M)=2r(x).$$

$$-M^2{\rho }^{\prime \prime }=2r(M\to 0).$$

$$\rho (x)\stackrel{M\to 0}{=}\frac{1}{\mu M^2}{\int }_{-\mu }^{+\mu }dy\left[\theta (x-y)\left(\mu -x\right)\left(\mu +y\right)+\theta (y-x)\left(\mu +x\right)\left(\mu -y\right)\right]r(y).$$

$$k=\left[\frac{\mu +x}{M}\right]+1,a=\left[\frac{\mu -x}{M}\right]-n+k.$$

$$\rho (x)=\frac{2k}{n+a}\sum _{l=0}^{n+a-k-l}(n+a-k-l)r(x+lM)+\frac{2(n+a-k)}{n+a}\sum _{l=1}^{k-1}(k-l)r(x-lM).$$

$$r(x)\simeq \frac{1}{\sqrt{M\left(\mu \pm x\right)}}\left(x\to \mp \mu \right).$$

$$\left\langle\mathcal{F}(x)\right\rangle=\frac{\left\llangle\mathcal{F}(x)\right\rrangle}{\left\llangle 1\right\rrangle}\,.$$

Full text

NORDITA-2019-003

UUITP-1/19

ICCUB-19-001

$N=2^{*}$** Phase Transitions and Holography **

Jorge G. Russo1,2, Erik Widén3,4 and Konstantin Zarembo3,4†††Also at ITEP, Moscow, Russia

1* Institució Catalana de Recerca i Estudis Avançats (ICREA),

Pg. Lluis Companys, 23, 08010 Barcelona, Spain.

2 Departament de Física Cuántica i Astrofísica and ICCUB

Universitat de Barcelona, Martí Franquès, 1, 08028 Barcelona, Spain.

3Nordita, KTH Royal Institute of Technology and Stockholm University, Roslagstullsbacken 23, SE-106 91 Stockholm, Sweden.

4Department of Physics and Astronomy, Uppsala University

SE-751 08 Uppsala, Sweden.

jorge.russo@icrea.cat, erikwid@kth.se, zarembo@nordita.org

Abstract

We clarify the relationship between probe analysis of the supergravity dual and the large-$N$ solution of the localization matrix model for the planar $\mathcal{N}=2^{*}$ super-Yang-Mills theory. A formalism inspired by supergravity allows us to systematically solve the matrix model at strong coupling. Quite surprisingly, we find that quantum phase transitions, known to occur in the $\mathcal{N}=2^{*}$ theory, start to be visible at the third order of the strong-coupling expansion and thus constitute a perturbative phenomenon on the string worldsheet.

1 Introduction

The ${\mathcal{N}}=2^{*}$ theory is obtained by softly breaking the ${\mathcal{N}}=4$ super Yang-Mills theory to ${\mathcal{N}}=2$ supersymmetry by the addition of a mass term. It provides an extremely interesting setup to understand holography in non-conformal systems, as the precise supergravity dual is known [1, 2]. One well known feature of the ${\mathcal{N}}=4$ theory is its smooth dependence on the ’t Hooft coupling $\lambda=g_{{\rm YM}}^{2}N$ all the way from $\lambda=0$ to $\lambda=\infty$. This feature changes dramatically when the theory is deformed. The resulting ${\mathcal{N}}=2^{*}$ theory undergoes an infinite number of phase transitions [3, 4]. The first transition occurs at $\lambda\approx 35.42$ and is followed by a sequence of secondary transitions appearing at regular intervals $\sqrt{\lambda}\approx\pi n$ with integer $n\gg 1$. The critical behavior is caused by rearrangement of the vacuum whereby new massless resonances enter the spectrum [3].

Holography allows one to explore the vacuum structure at strong coupling by placing a D3-brane probe in the dual supergravity background [2, 5]. Quite remarkably, the eigenvalue distribution predicted by the supergravity analysis obeys Wigner semicircle law [2], a hallmark of random matrix theory. This finding led the authors of [6] to speculate that Gaussian random matrices underlie vacuum structure of the $\mathcal{N}=2^{*}$ SYM and other strongly-coupled $\mathcal{N}=2$ gauge theories at large-$N$. We now know that random matrices arise in $\mathcal{N}=2$ theories upon supersymmetric localization, which maps the partition function on $S^{4}$ to a zero-dimensional matrix model [7]. Albeit the matrix model is not Gaussian in the $\mathcal{N}=2^{*}$ case (it is only such for $\mathcal{N}=4$ SYM [8, 9, 7]), its course-grained, averaged eigenvalue distribution indeed takes semicircular Wigner shape at strong coupling, in spectacular confirmation of the holographic predictions [10].

The gauge/gravity duality for the ${\mathcal{N}}=2^{*}$ theory has withstood many other remarkable, rigorous tests, by confronting the exact results obtained from supersymmetric localization with the corresponding gravity dual observables, such as Wilson loops [10, 11, 12, 13] and free energy on $S^{4}$ [14]. Matrix model tools [10, 3, 4, 15, 16, 17, 18], advances in supergravity [19, 14, 13] and deeper understanding of string [20, 10, 12] and D-brane [2, 5, 11] probes have been instrumental in these developments. While localization presents abundant evidence for phase transitions [3, 4, 16, 17, 18, 21, 22], they so far remained in limbo on the gravity side.

A challenge is to understand how the non-analytic features of the quantum phase transitions are reflected in the holographic description. Strictly speaking, the supergravity background describes the $\lambda=\infty$ theory. In this limit, the gauge theory has infinitely many phases and it approaches an average, coarse-grained, smooth description where the non-analytic behavior does not show up in the leading approximation [3, 4]. For example, the vacuum expectation value of the 1/2 BPS circular Wilson loop approaches $\ln\langle W(C)\rangle\approx ML_{C}\sqrt{\lambda}$, which would seem to have a continuous dependence on the coupling. This average description of the strong coupling regime makes particularly difficult and challenging to reveal the existence of phase transitions on the string theory side of the duality.

Keeping this mind, we revisit the supergravity analysis of the vacuum structure [2] in a hope to establish a tighter connection to localization. A closer look at the supergravity calculation reveals striking similarities to the strong-coupling solution of the matrix model [17]. Armed with this observation we extend and streamline large-$N$ analysis of the localization matrix model and push the strong-coupling solution to the third order in $1/\sqrt{\lambda}$. Quite surprisingly, phase transitions become fully visible starting with this order enabling, for example, computation of the critical indices.

The upshot of our analysis is that the phase transitions constitute a perturbative phenomenon from the supergravity/string theory point of view. They were not visible so far for a simple reason that existent string or gravity calculations have explored the LO and NLO of the strong-coupling expansion, while the phase transitions appear at NNLO.

2 Probe analysis and eigenvalue density

The probe analysis of the Pilch-Warner geometry defines an exact distribution of D3-branes on the enhançon [2], which appears to match the density of eigenvalues in the gauge theory, obtained from localization [10, 3, 4]. Here we review the supergravity calculation, slightly generalizing it to include finite-mass effects. This is crucial for potential non-analytic features responsible for the phase transitions. We will later incorporate elements of the supergravity analysis in the solution of the localization matrix model.

In the gauge-theory language, a probe D3-brane placed at position $u$ on the locus of marginal stability corresponds to a point on the moduli space of vacua of the $SU(N+1)$ theory parametrized by

[TABLE]

where $\left\{a_{i}\right\}$ describe the vacuum of the $SU(N)$ theory without the probe, sourced by the geometry itself. The kinetic term of the probe is characterized by the effective coupling

[TABLE]

where $\mathcal{F}$ is the prepotential. The instanton part of the prepotential is irrelevant in the large $N$ limit, so we just need to include the classical and one-loop parts. These are given by

[TABLE]

where we have omitted an additive term that does not depend on $u$ and thus does not contribute to the effective coupling.

Expanding the probe action to quadratic order in the field strength, one can identify the effective coupling with the supergravity axi-dilaton:

[TABLE]

It was verified in [2] that the locus where the probe experiences no force is real two-dimensional (one complex dimensional). Evaluating the axi-dilaton on the locus of marginal stability, [2] found that

[TABLE]

On the other hand, from (2.2) and (2.3):

[TABLE]

The eigenvalue density can be computed by comparing (2.5) and (2.6). This calculation was carried in [2] using the approximation $M^{2}\ll(u-a_{i})^{2}$. Such an approximation is justified at strong coupling, because, as we shall see, $a_{i}\sim u\sim\sqrt{\lambda}\,M$. In this approximation,

[TABLE]

Taking discontinuity across the cut and comparing to (2.5), one finds

[TABLE]

This gives

[TABLE]

the result found in [2]. It is easy to check that not only the discontinuity matches, but the whole function (2.5) is reproduced by the Wigner density. The width of the eigenvalue distribution is fixed by normalization:

[TABLE]

which justifies a posteriori the simplifying assumption that $M^{2}\ll(u-a_{i})^{2}$. It would have been technically cleaner to solve for the density without making this assumption and to expand in $M$ only after establishing that $\mu\gg M$. Let us repeat the calculation without making approximations in intermediate steps.

Differentiating (2.6) with respect to $u$, we obtain

[TABLE]

This is to be compared with the supergravity prediction, obtained by differentiating (2.5):

[TABLE]

Equating, as before, discontinuities of the two functions, we get:

[TABLE]

Taylor expanding in $M$ approximates the left-hand side to $-M^{2}\rho^{\prime\prime}(x)$, resulting in a differential equation solved by the Wigner distribution (2.9). The finite-difference equation can be solved as well, without the small-mass approximation. The technical details are postponed till the next section, but the key qualitative features are apparent just from the form of the equation. The boundary conditions on the density are set at $x=\pm\mu$. In between of those points the evolution goes in steps of $M$. It thus crucially depends on how many steps separate the endpoints of the interval. If one gets from $-\mu$ to $\mu$ in an integer number of steps, namely when $2\mu=nM$, a resonance occurs and the solution undergoes a discontinuous mutation. These are precisely the phase transitions visible on the gauge-theory side [3].

Still, it is fair to say that there is no direct evidence of phase transitions from the pure gravitational description. The gravitational input is (2.12), a nice continuous function of $\mu$. Irregularities appear when the gravitational data is recycled into the eigenvalue density which, by itself, does not have any direct holographic description. A clear signature of the phase transitions would be a non-analytic behavior of the expectation values (moments of the eigenvalue density). Such non-analyticities do not occur in the supergravity approximation, a property which we will be able to quantify by a systematic strong-coupling calculation on the gauge-theory side.

The supergravity analysis suggests that $\tau^{\prime}(u)$ is a more convenient characterization of the eigenvalue distribution than the density. The latter has a very irregular structure at strong coupling [16, 17]. An ansatz, first proposed in [17, 23], expresses it through a regular function bearing a lot of similarity to (2.12). The ansatz can be shown to solve the finite-difference equation (2.13), and in the next section we sharpen and extend this connection between supergravity and gauge-theory quantities.

3 Saddle point equation

The partition function of the $\mathcal{N}=2^{*}$ theory on $S^{4}$ can be computed exactly by supersymmetric localization [7]:

[TABLE]

where $a_{i}$ are the eigenvalues of the scalar from the vector multiplet:

[TABLE]

and

[TABLE]

The instanton contribution can be neglected at large-$N$111The instanton weight $\,{\rm e}\,^{-8\pi^{2}N/\lambda}$ is exponentially small at large-$N$, but the moduli integral can potentially overcome the exponential suppression [24]. It has been explicitly checked that the instanton explosion never happens in the planar $\mathcal{N}=2^{*}$ theory [4]., while the remaining integral is of the saddle-point type and is solved by a configuration of the eigenvalues that minimizes the effective action.

The eigenvalue integral (3.1) is written in terms of the dimensionless variables. The canonical dimensions are recovered by rescaling $a_{i}\rightarrow a_{i}R$, $M\rightarrow MR$, where $R$ is the radius of the four-sphere. In the decompactification limit $R\rightarrow\infty$ the problem boils down to solving a singular integral equation[4]222Which arises upon differentiating the original saddle-point equation twice and using the large-argument asymptotics $H(x)\approx|x|^{-x^{2}}$.

[TABLE]

supplemented by two normalization conditions:

[TABLE]

The eigenvalue density is defined on the interval $(-\mu,\mu)$ and has the inverse square root singularities at its endpoints. The symmetry-breaking mass scale $\mu$ is a function of the ’t Hooft coupling, implicitly determined by the auxiliary normalization conditions. A known exact solution of these equations in terms of theta-functions [3, 4] is valid for sufficiently small $\lambda$ and terminates at a quantum phase transition at $\lambda_{c}\approx 35.42$. Remarkably, the same $\lambda_{c}$ appears to be a special value of the coupling for all $SU(N)$ theories with $N\geqslant 3$, where an exact, analytic solution is also known up to $\lambda_{c}$ [25]. On the other hand, the theory with $SU(2)$ gauge group has a smooth behavior with $\lambda$ all the way from 0 to infinity [21, 25].

At large $\lambda$, as discussed above, the eigenvalue density approaches

[TABLE]

with

[TABLE]

which matches the prediction from the Pilch-Warner supergravity solution. This is, however, a coarse-graining description which hides a complicated structure involving multiple discontinuities.

In general, the structure of the density qualitatively changes when the length of the eigenvalue interval crosses an integer multiple of $M$. A new resonance appears at that point and the system undergoes a phase transition. It is convenient to introduce variables $n$ and $\Delta$ that characterize the size of the interval relative to $M$ [23]:

[TABLE]

where $[x]$, $\{x\}$ denote integer/fractional part of $x$. An integer $n$ enumerates phases of the theory.

The resonance terms in the saddle-point equation (3.4) induce singularities in the eigenvalue distribution, at points $\mu-lM$ and $-\mu+lM$ with $l=1\ldots n$. As a result, the density has a complicated, irregular structure, especially at strong coupling, featuring a growing number of peaks and jumps. But, as we already discussed, the density is not the most convenient characterization of the eigenvalue distribution.

Inspired by the supergravity analysis we introduce a resolvent-type function

[TABLE]

which coincides with $\tau^{\prime}(u)$ up to normalization.

The resolvent has a cut across the eigenvalue interval. The continuous part has to vanish according to the eigenvalue equation, which is equivalent to the condition

[TABLE]

Since its continuous part vanishes, the resolvent equals to $\pm$ its discontinuity on the cut. Defining $r(x)\equiv r(x+i0)$ for definiteness, we have

[TABLE]

This is as a second order finite-difference equation for $\rho(x)$.

Although the equation is the same as (2.13), the actual problem here is different. Before, the function $r(u)$ was given, being extracted from the supergravity probe analysis. Now $r(u)$ has to be self-consistently determined. We will do it by solving the difference equation for the density, substituting the solution back into the saddle-point equation and reformulating the latter as an integral equation for the resolvent.

In the limit $M\rightarrow 0$, the difference equation becomes differential:

[TABLE]

The solution of the latter is given by the Green’s function with Dirichlet boundary conditions at $x=\pm\mu$:

[TABLE]

When applied to (2.12), this formula reproduces (2.9).

To write down the Green’s function of the finite-difference operator, we introduce integer variables $k$ and $a$, illustrated in fig. 1:

[TABLE]

Then

[TABLE]

Replacing sums by integrals we get back to (3.14) in the continuum limit.

Exactly the same formula was introduced in [17] on phenomenological basis, just as an ansatz motivated by the structure of the saddle-point equation (3.4). We see that this representation naturally arises when eigenvalue distribution is characterized by the resolvent, thus making direct contact to supergravity.

The saddle-point equation (3.4), the finite-difference equation (3.12) and its solution (3.16) are invariant under simultaneous rescaling $r\rightarrow Cr$, $\rho\rightarrow C\rho$ with constant $C$. We can use this freedom to normalize the resolvent such that near the boundaries of the interval it behaves as

[TABLE]

Fixing the edge behavior is in general inconsistent with canonical normalization of the density. But for our purposes edge normalization is more convenient and we prefer to abandon (3.6) in favor of this new condition. The resulting averages, denoted by \left\llangle\ldots\right\rrangle, will be wrongly normalized and will have to be divided by a common factor to get correct expectation values:

[TABLE]

While the resolvent is a continuous function, the discrete map (3.16) induces singularities at integer multiples of $M$ away from the endpoints, resulting in a cuspy profile in the middle of the eigenvalue distribution. The detailed structure of cusps is quite intricate. The map from $r$ to $\rho$ splits the interval $(-\mu,\mu)$ into $2n+1$ subintervals according to the value of $a$, equal to either $1$ or to $2$. The density is a smooth function on the $a=1$ intervals (“b-intervals”, in terminology of [17, 23]) and has inverse square root singularities at the endpoints of the “a-intervals”, those with $a=2$. The overall structure of the map $r(x)\rightarrow\rho(x)$ is illustrated in fig. 2.

3.1 Average formulas

The next step is to express average quantities in terms of the resolvent. In computing averages it is convenient to split the integral over $x$ into summation over $k$ defined in (3.15) and integration over $\xi\in(0,1)$ defined as

[TABLE]

Then

[TABLE]

where, in these variables,

[TABLE]

The formula looks rather ugly and uninformative but, recalling that the expression in square brackets is the Green’s function of the discrete Laplacian, we may anticipate simplifications to occur for functions that are total second derivatives:

[TABLE]

It is indeed easy to show by direct computation that for such functions summation over $l$ localizes to $l=m$, up to boundary terms at $l=0$ and $l=n+a$:

[TABLE]

The second line can be written as an integral of a smooth function over the whole interval from $-\mu$ to $\mu$. Changing summation variable in the last term from $m$ to $n+a-m$, and integration variable from $\xi$ to $2-a+\Delta-\xi$ brings the average formula to a very neat form

[TABLE]

where in the last step we also changed $\xi$ to $1-\xi$, and we have introduced a function:

[TABLE]

The average formula (3.24) has a simple meaning. The first term is the continuum limit that arises upon substituting (3.14) for the density, approximating the finite-difference Laplacian by the second derivative and dropping the boundary terms. The last term and the function $h(\xi)$ encode boundary contributions and are the sole remnants of the discreteness inherent to the map (3.16).

The $h(\xi)$ function will play a prominent rôle in the formalism that we are developing and for future use we list here its salient features. The function is defined on the interval $(0,1)$ and has inverse square root singularities at $\xi=1$ and $\xi=1-\Delta$:

[TABLE]

as follows from (3.17). Notice that when $n\rightarrow\infty$, which corresponds to the strong coupling limit, the singularity at the endpoint is parametrically stronger than the midpoint one. The overall shape of the function $h(\xi)$ is illustrated in fig. 3.

3.2 Regular form of integral equation

The average formula (3.24) applied to the function $\mathcal{F}(y)=1/(x-y)$ transforms the original saddle-point equation (3.4) into an equation for the resolvent:

[TABLE]

This equation, first derived in [17], is much better behaved than the original integral equation for the density. Now the integral operator has the standard Hilbert kernel and the source term is regular on the whole interval from $-\mu$ to $\mu$. General theorems about singular integral equations [26] guarantee existence and uniqueness of the solution with boundary conditions (3.17), which in addition is well-behaved in the interior of the eigenvalue interval. The well-known inversion formula for the Hilbert kernel [26] yields [17]:

[TABLE]

Having the resolvent, we can formulate the problem entirely in terms of the auxiliary function $h(\xi)$.

We start with the average formula. Substitution of the integral representation for $r(x)$ into (3.24) gives

[TABLE]

This result can be further simplified. The integral over $x$ is equivalent to a contour integral, with the contour of integration encircling the cut counterclockwise. Enlarging the contour as shown in fig. 4, from $C^{\prime}$ to $C$, picks the poles at $z=\pm(\mu+M)$. The residues conspire to cancel the last two terms, leaving a nice compact form of the average formula:

[TABLE]

The hat-operator is defined as

[TABLE]

The contour of integration should leave outside all singularities of $\hat{\Delta}^{-1}\mathcal{F}(z)$.

This form of the average formula is ideally suited for the strong-coupling expansion. Inverting the discrete Laplacian is easy for simple functions, and upon computing the contour integral defining $\widehat{\mathcal{F}}(\xi)$ the average reduces to a simple integral with the weight determined by a so far unknown $h(\xi)$.

Let us illustrate this strategy on some simple examples:

[TABLE]

For

[TABLE]

we get:

[TABLE]

As an application, we can reformulate the normalization condition (3.5) in terms of the auxiliary function $h(\xi)$. Since (3.5) is equivalent to

[TABLE]

we find:

[TABLE]

Once $h(\xi)$ is known, this equation would determine $\lambda$ as a function of $\mu/M$.

Another quantity of interest is the vacuum susceptibility:

[TABLE]

related to the derivative of the partition function with respect to the coupling: $\chi\propto\partial\ln Z/\partial\ln\lambda$. The susceptibility can be computed from (3.32), (3.33).

All averages can thus be expressed through a single function $h(\xi)$, obtained by applying a linear finite-difference operator (3.25) to the resolvent. The resolvent, in its turn, admits an integral representation (3.28) in terms of $h(\xi)$. Equations (3.28) and (3.25) can be solved numerically by iteration. These equations determine $r$ and $h$ up to a common normalization factor, which can be fixed by the asymptotic conditions (3.17) or (3.26), but in practice normalization is not really important because in physical averages an overall factor cancels out, and one can use any convenient normalization whatsoever.

It is possible to eliminate the resolvent from the integral equation and to write down a closed equation for the function $h(\xi)$ only. Applying the difference operator from (3.25) to both sides of (3.28), we get a linear integral equation

[TABLE]

where the kernel is given by

[TABLE]

The formulas (3.30), (3.31), (3.39) and (3.40) completely characterize the eigenvalue distribution, without the need to compute the density nor the resolvent. An arbitrary average can be computed by first solving the integral equation for $h(\xi)$ and then evaluating the integral in (3.30). The resolvent and the density can be reconstructed from (3.28) and (3.16) if necessary, but they are not needed for evaluating expectation values. The system of equations for the eigenvalue distribution written in this form is ideally suited for the strong-coupling expansion.

4 Strong-coupling expansion

As anticipated from (3.8) (we will rederive this formula shortly within the formalism described above), $\mu$ grows with $\lambda$ and the strong-coupling regime corresponds to $\mu\gg M$. The relationship between the ’t Hooft coupling and the dimensionless ratio $\mu/M$ follows from (3.37). It takes zero effort to solve this equation at strong coupling. Expanding in $M/\mu$ and taking into account that

[TABLE]

we get

[TABLE]

which confirms (3.8). All the dependence on the function $h(\xi)$ drops out!

We can find the large-$\lambda$ asymptotics of the vacuum susceptibility, likewise without any detailed knowledge of the function $h(\xi)$. Keeping only the leading terms in (3.33), we get:

[TABLE]

Inspecting the average formulas (3.30), (3.31) and the kernel (3.40) of the integral equation (3.39) one can conclude that

• •

Any expectation value has a regular expansion in $M/\mu$ whose coefficients are expressed through the moments333Our definition differs from the one in [17] by an unimportant normalization factor. of the function $h(\xi)$:

[TABLE]

• •

According to (3.8), this expansion translates into a regular expansion in $1/\sqrt{\lambda}$. This structure of the strong-coupling expansion is in accord with expectations from holography, since in the dual description $2\pi/\sqrt{\lambda}$ is identified with the coupling constant of the string sigma-model.

• •

However, the function $h(\xi)$ and consequently its moments depend on $\Delta$ defined in (3.9), which is clear from fig. 3 or from the defining integral equation for $h(\xi)$.

• •

At strong coupling,

[TABLE]

Due to the dependence of $\mathbbm{K}_{n}$ on $\Delta$, expectation values are not really analytic in $\lambda$, and undergo quantum phase transitions each time $\sqrt{\lambda}$ passes an integer multiple of $\pi$:

[TABLE]

This an approximate formula, valid at asymptotically large $\lambda$. An exact equation for the critical points is $2\mu(\lambda_{c}^{(n)})=nM$ [3].

The first strong-coupling correction to $\mu(\lambda)$ was calculated in [16, 17] and does not show any non-analyticity in $\lambda$, it was later confirmed by a direct one-loop calculation in the string sigma-model [12]. We will push the strong-coupling expansion one order further. Our goal here is two-fold: to test the formalism developed in the previous section, and to see if the phase transitions are indeed visible within the strong-coupling expansion.

Expanding (3.30), (3.32) – (3.35) to the third order in $M/\mu$ we get:

[TABLE]

It is convenient, at this point, to introduce the reduced moments, which do not depend on normalization of $h(\xi)$:

[TABLE]

In terms of those, (3.36) becomes

[TABLE]

Inverting this relationship we obtain:

[TABLE]

The vacuum susceptibility, to this order in $M/\mu$ is given by

[TABLE]

Substituting $\lambda$ using (4.8), we get:

[TABLE]

The term linear cancels. In terms of the ’t Hooft coupling,

[TABLE]

Interestingly, the first, NLO correction to the susceptibility vanishes.

In order to compute $\mu$ and $\chi$ to NNLO we need the coefficients $\mathbbm{K}_{0}$, $\mathbbm{K}_{1}$, $\mathbbm{K}_{2}$ to NLO. In fact, the NLO correction enters only for the ratio $\mathbbm{k}_{1}=\mathbbm{K}_{1}/\mathbbm{K}_{0}$. To compute the moments we need to know the function $h(\xi)$, determined by the integral equation (3.39). We solve this equation at strong coupling by expanding in $M/\mu$.

Taking the limit $M\rightarrow 0$, $n\rightarrow\infty$ in the kernel (3.40) we get:

[TABLE]

The solution at this order is [17]

[TABLE]

Normalization here is chosen to comply with (3.26), but in principle it is arbitrary, not fixed by the equation, nor affecting any average quantities. It will be important to keep in mind this renormalization ambiguity. Application to $h_{\infty}$ of the integral operator with the kernel $G_{\infty}$ returns a telescoping sum, that combines back into $h_{\infty}$:

[TABLE]

which proves that (4.14) solves (3.39) at the leading order in the strong-coupling expansion.

The moments, computed with the leading-order solution, are

[TABLE]

In particular,

[TABLE]

From those we can compute the first-order effective string tension444The part of the expansion contributed by the leading order $h_{\infty}$ is given by

$\frac{4\pi^{2}}{\lambda}=\frac{M^{2}}{\mu^{2}}\left[1-\frac{M}{\mu}+\frac{13M^{2}}{16\mu^{2}}-\frac{315M^{3}}{512\mu^{3}}+...\right]\ \longrightarrow\ \mu=\frac{\sqrt{\lambda}\,M}{2\pi}\left(1-\frac{\pi}{\sqrt{\lambda}}+\frac{\pi^{2}}{8\lambda}+\frac{5\pi^{3}}{128\lambda^{\frac{3}{2}}}+...\right).$

$O(1/\lambda)$ and successive terms will be affected by corrections to (4.14). [16, 17]:

[TABLE]

as well as the second-order vacuum susceptibility:

[TABLE]

To compute the effective string tension to the second order, we need to know the auxiliary function to the first order in $M/\mu$:

[TABLE]

Expanding the kernel (3.40) to the same order,

[TABLE]

and substituting into the integral equation (3.39), we get an integral equation for $h_{1}$:

[TABLE]

The first-order corrections to the kernel (3.40) come from two effects: from $\mathcal{O}(M/\mu)$ corrections to the terms in the sum with $m=\mathcal{O}(1)$ and from the terms with $n-m=\mathcal{O}(1)$. The latter contribution can be obtained by changing summation variable $m\rightarrow n+a-m$ and subsequently extending the summation range to infinity. Altogether we get:

[TABLE]

Applying it to $h_{\infty}$ we get:

[TABLE]

Notice that this expression has the same structure of singularities as (3.26). Using the identity

[TABLE]

we get the explicit form of the equation for $h_{1}$:

[TABLE]

There is a freedom of shifts by the asymptotic solution (4.14):

[TABLE]

with arbitrary constant $C$. This ambiguity reflects arbitrary normalization of the exact solution. It is easy to check that the normalized averages do not depend on the constant $C$ to the requisite order in $M/\mu$, and so the constant $C$ can be chosen arbitrarily. Iterative solution corresponds to one possible choice.

The first-order correction to the normalized moment

[TABLE]

is expressed through the solution of the integral equation as

[TABLE]

The shift symmetry (4.27) leaves this expression invariant. Substituting this result, along with the unperturbed moments (4.17), into (4.9) gives the effective string tension at NNLO of the strong-coupling expansion:

[TABLE]

It is important to note that at this order $\mu$ starts depending on $\Delta$, which is a fractional part of $\sqrt{\lambda}$ in units of $\pi$. This implies that the dependence on $\lambda$ is not analytic, even though $\mu$ can be expanded in regular power series in the inverse coupling. The non-analytic behavior occurs at $\Delta=0,1$ and will be discussed in the following section.

The combination of moments that appears in the strong-coupling expansion of susceptibility (4.12) is expressed through another function

[TABLE]

We find:

[TABLE]

Because the susceptibility does not receive corrections at $\mathcal{O}(1/\sqrt{\lambda})$, the dependence on $\Delta$ and the ensuing non-analyticity is postponed by one order in the strong-coupling expansion. We have:

[TABLE]

where $C$ is some numerical constant that we are not going to compute.

We have calculated the functions $u(\Delta)$ and $v(\Delta)$ numerically (fig. 5), solving the integral equation (4.26) by two different methods, either expanding $h_{1}(\xi)$ in the moments (appendix A) or, more directly, by Galerkin method (appendix B). The results agree within error bars (see appendix B) and demonstrate that both functions depend non-trivially on $\Delta$.

5 Phase transitions

When (4.8) was inverted to arrive at (4.9) the moments $\mathbbm{k}_{n}$ were assumed to be constant, but in fact they are also functions of $\mu$ because of their dependence on $\Delta$. Fortunately, this dependence starts at a rather high order in the strong-coupling expansion and can be taken into account perturbatively. But still, the right-hand side of (4.8), taken at face value, has different functional dependence on $\mu$ for $\mu=M(n-\epsilon)/2$ and $\mu=M(n+\epsilon)/2$. In the former case, $\Delta=1-\epsilon$, while in the latter, $\Delta=\epsilon$. In practice this means that we should use $u(1-\epsilon)$ in (4.28) slightly below $\mu=Mn/2$ and $u(\epsilon)$ slightly above. In other words, near the critical point functional dependence of the coupling constant on $\mu=M(n+\epsilon)/2$ is expressed through $u(\epsilon)$ for $\epsilon>0$ and through $u(1+\epsilon)$ for $\epsilon<0$. Is a function so defined continuous? We are going to argue that it is, but that a milder non-analyticity occurs at $\epsilon=0$, leading to critical behavior at $\mu=Mn/2$ or for $\lambda=\lambda_{c}^{(n)}$.

Continuity in $\epsilon$ trivially follows from the integral equation (4.24), because there $\Delta$ appears only in combination $\Delta-\theta(\Delta-1+\eta)$, which equals zero for both $\Delta=0$ and $\Delta=1$, independently of $\eta$. The solution therefore is the same for $\Delta=0,1$ and so are all the moments of the function $h_{1}$. The question is what happens slightly below and slightly above the critical point, where

[TABLE]

The analysis of the integral equation (4.26) in the critical region is rather intricate, and we carry out the detailed calculation in the appendix C. Here we present a simplified qualitative argument that illustrates all the salient features.

5.1 Edge singularity

We are interested in averages of the form

[TABLE]

where $f(\xi)$ is a polynomial. To get a rough idea of what happens for small $\Delta$ or $1-\Delta$, we keep only the source term in the integral equation (4.26):

[TABLE]

This can be regarded as the zeroth-order approximation in the iterative solution of the integral equation (4.26), albeit there is no small parameter that would justify neglecting the integral term.

It is convenient to use the identity

[TABLE]

Taking into account that $\zeta_{\frac{1}{2}}(\Delta+\eta+1)$ is regular in the full intervals $\eta\in[0,1]$, $\Delta\in[0,1]$, we can write:

[TABLE]

For $\Delta\ll 1$, the first term blows up at $\eta\rightarrow 0$ and the second term blows up at $\eta\rightarrow 1$. We can thus set the argument of $f(\xi)$ to zero in the first term and to one in the second term, if we are only interested in the non-analytic behavior at $\Delta\rightarrow 0$:

[TABLE]

where $\stackrel{{\scriptstyle!}}{{=}}$ denotes equality up to an analytic function of $\Delta$.

When $\Delta\rightarrow 1$, the first term in (5.5) is regular everywhere, while the second term blows up at $\eta\rightarrow 0$. This justifies setting $f(\xi)\simeq f(0)$ and yields:

[TABLE]

Next iterations change these results in two ways. The terms containing $\sqrt{\Delta}$ cancel exactly, as we show in appendix C. The logarithmic terms remain, but their coefficients get modified. The exact non-analytic contributions take a neat form in terms of the variable $\epsilon$ introduced in (5.1):

[TABLE]

where the contour of integration encircles the interval $(0,1)$ clockwise, leaving singularities of $f(z)$ outside555The previous simplified calculation neglecting the integral term gives the same formula with $f(z)-f(0)$ replaced by $-f(0)$. Iterations of the integral equation generate the $f(z)$ term.. For a polynomial this gives:

[TABLE]

A derivation of these results is given in appendix C.

Applying these findings to $u(\Delta)$ and $v(\Delta)$ from (4.29), (4.31), we get:

[TABLE]

We have solved for the complete $u(\Delta)$, $v(\Delta)$ functions numerically and, in particular, checked that coefficients of the logarithmic terms at $\Delta\rightarrow 0,1$ match with the analytic predictions.

The numerical results are shown in figures 5 and 5. Figures 6, 6 display critical behavior in the vicinity of $\Delta=0$ and $\Delta=1$. The numerical data accurately fits the expected $\Delta\log\Delta$ behavior. We find

[TABLE]

in good agreement with the analytic result for the coefficient of $\Delta\log\Delta$, $2/\pi\cong 0.637$.

The function $v(\Delta)$ is less singular. Numerics in this case can be well fit by

[TABLE]

The absence of the log-linear term is in agreement with the analytic predictions, and actually follows from the symmetry under $\Delta\rightarrow 1/2-\Delta$ reflection. The function $v(\Delta)$ is even under reflection, while $u(\Delta)$ is odd upon a shift by a constant.

The coefficient of the singular term in $v(\Delta)$ is numerically close to $6/\pi\cong 1.910$, and we give the following prediction for the singular part of the function $v(\Delta)$:

[TABLE]

It would be interesting to derive this result analytically by pushing the perturbative calculation in the appendix C one order further.

5.2 Critical indices

The system undergoes phase transitions each time $\Delta$ jumps from $1$ to [math], which happens at

[TABLE]

We now discuss how non-analyticities that arise at these points affect various quantities of interest. Consider first the dynamical scale $\mu$. To find it as a function of the ’t Hooft coupling we need to invert (4.8). The relevant terms near the critical point are

[TABLE]

where we used (4.28) and omitted unimportant analytic terms. Introducing scaling variables

[TABLE]

expanding (5.16) in $\delta$ and $\epsilon$, and using (5.10) for $u(\Delta)$, we get:

[TABLE]

This formula was derived at large $n\sim\sqrt{\lambda}/\pi$, and in the course of the derivation we assumed that the second term in the brackets is a small correction. But near the phase transition this term is logarithmically enhanced and when $\ln|\epsilon|\sim n$ it is of the same order as the nominally larger leading-order term. Strictly speaking, we cannot infer the behavior of $\mu$ arbitrarily close to the critical point without resumming large logs in the strong-coupling expansion. However, in almost all problems where perturbation theory is logarithmically enhanced, logarithms exponentiate and log behavior seen in perturbation theory signals power-like scaling with a non-trivial scaling exponent. Taking exponentiation as a plausible assumption, we can write, to the same degree of accuracy:

[TABLE]

More generally we can posit power-law scaling of $\mu$ near the critical points:

[TABLE]

This is consistent with the scaling just derived. Another justification comes from the exact solution in the weak-coupling phase [3] which also exhibits scaling behavior of this form with the critical exponent [4]:

[TABLE]

From (5.18) we infer that at strong coupling:

[TABLE]

Interestingly, the strong-coupling approximation is not far off the exact result even for $n=1$.

One can introduce other critical exponents. For example, the vacuum susceptibility is expected to scale as

[TABLE]

The weak-coupling solution predicts [4]666The critical exponent is integer, but it does not mean that susceptibility is analytic. The coefficients in front of the cubic term are different on the two sides of the transition, see [4] for more details.:

[TABLE]

The log-enhancement is lacking in the strong-coupling expansion of the susceptibility, in virtue of (5.11), which means that at this order $\chi\sim\epsilon$, and the scaling behavior is governed by the same critical exponent as in (5.18):

[TABLE]

Finally, there is a critical exponent associated with the end-point behavior of the density right at the critical point:

[TABLE]

This critical exponent has been calculated exactly [17]:

[TABLE]

At the first phase transition [4],

[TABLE]

While there is no breaking of symmetries at each phase transition, one can define a sequence of “order parameters” ${\cal O}_{n}$ of mass dimension one, that signal the onset of each phase transition. They are defined as ${\cal O}_{n}=1/\ell_{n}$, where

[TABLE]

The order parameter $\ell_{n}$ represents a correlation length that diverges at the $n$th-transition point whenever new resonances occur. These quantities generalize the similar parameter introduced in section 3.5.4 of [4] for the first ($n=1$) phase transition. The near-critical behavior can be computed by using (5.25). The integral has a leading contribution given by

[TABLE]

Using now (5.19), we find

[TABLE]

We see that the critical exponent $\kappa_{n}$ is not independent, but derived from $\alpha_{n}$, i.e. it is dictated by the behavior of the eigenvalue density $\rho$ near the endpoint. Thus

[TABLE]

6 Conclusions

The strong-coupling expansion of the effective string tension (4.30) has an expected structure of perturbative series in the string sigma-model. The first term is reproduced by the classical area law [10]. The second term, corresponding to the one-loop sigma-model correction, was computed by a semi-analytic calculation and it is also in perfect agreement with the NLO of the localization result [12]. The NNLO term then corresponds to the two-loop correction in the sigma-model, and it is there that we expect to see signatures of the phase transitions.

It is not entirely clear why two loops are sensitive to the phase transitions, while the first two orders are not. Non-analytic dependence on the coupling constant is not that uncommon in quantum field theory and typically arises as a consequence of IR divergences, which entail resummation of perturbative series. Perhaps a similar mechanism is at work here. If so, it must depend on the structure of interactions on the string worldsheet. This is consistent with the observation that at the two lowest orders non-analyticities do not arise, because the classical area law and the one-loop correction due to string fluctuations are not really sensitive to how string modes interact with each other.

Our methods, at least in principle, allow one to develop the strong-coupling expansion in the planar $\mathcal{N}=2^{*}$ theory to any desired order. It would be interesting to carry out higher-order calculations explicitly, in particular to compute the $1/n^{2}$ correction to the critical indices. Even more interesting would be to push perturbative expansion of the string sigma-model [12] beyond the one-loop order. Our results predict that quantum phase transitions should be visible on the string side of the holographic duality, once these two-loop corrections are properly taken into account. It is quite remarkable that such a dramatic effect appears to have a perturbative origin in string theory.

Finally, it would be extremely interesting to understand the phase structure of the $SU(N)$ ${\cal N}=2^{*}$ theory for any given $N$. At finite $N$, the exact analysis of the different phases including instantons has an elegant description in terms of the Seiberg-Witten curve [21, 22, 25]. In the decompactification limit $R\to\infty$ , Pestun’s partition function is computed by saddle-points corresponding to extrema of the action, given by $-R^{2}{\rm Re}\left(4\pi i{\cal F}\right)$. Thus, the partition function is computed by critical points of the prepotential ${\cal F}(a_{i})$ where the ${\rm Im}(a_{Di})={\rm Im}\frac{\partial{\cal F}}{\partial a_{i}}$ are required to vanish on the integration domain [21, 25]. Such singular points describe massless dyon singularities. Phase transitions may occur as the coupling $g_{\rm YM}$ is gradually increased due to the existence of many different competing saddle-points. A strong indication that there might be similar phase transitions at any finite $N$, with $N\geqslant 3$, was found in [25], where the saddle-point corresponding to maximal degeneration was found to exist only for $\lambda<\lambda_{c}$, $\lambda_{c}\approx 35.42$. Surprisingly, this critical coupling is the same for any $N$ and corresponds to the first critical point of the large $N$ phase transitions discussed here. An open problem is to identify the singular points that dominate the partition function integral at $\lambda>\lambda_{c}$ and the total number of different phases for a given $N$. The structure of degenerate points becomes increasingly more difficult as $N$ is increased, but perhaps for low-rank groups such as $SU(3)$ the different phases can be identified (a discussion of the singularity structure can be found in [27]).

Acknowledgments

We would like to thank B. Assel, I. Kostov, J. Penedones, A. Sever, J. Troost and A. Zhiboedov for discussions. The work of K.Z. and E.W. was supported by the ERC advanced grant No 341222, by the Swedish Research Council (VR) grant 2013-4329, by the grant ”Exact Results in Gauge and String Theories” from the Knut and Alice Wallenberg foundation, and by RFBR grant 18-01-00460 A. J.G.R. acknowledges financial support from projects 2017-SGR-929, MINECO grant FPA2016-76005-C.

Appendix A Computation of $u(\Delta)$ and $v(\Delta)$

Our starting point is $G_{\infty}(\eta,\xi)$,

[TABLE]

We have separated the term $k=1$. The terms with $k\geqslant 2$ can be expanded in powers of $\xi$ and the expansion converges for all $\eta\in[0,1]$. The remaining sum from $k=2$ to $\infty$ can then be expressed in terms of Hurwitz $\zeta$ functions. We obtain

[TABLE]

Now we consider the equation (4.26) for $h_{1}(\eta)$. Using the above expansion of $G_{\infty}(\eta,\xi)$, we find

[TABLE]

where

[TABLE]

This equation can be converted into a linear algebraic equation for the moments by multiplying by $\eta^{n+1/2}$ and integrating over $\eta$. All the integrals are convergent. The first integral may be computed by residues. The result is given by the following formula:

[TABLE]

where

[TABLE]

Expanding in powers of $\xi$, the integral over $\xi$ can be computed in terms of moments. We find

[TABLE]

with

[TABLE]

Note that the first term contains moments

[TABLE]

A second equation for the moments $\tilde{k}_{r-1/2}$ can be obtained by performing the integration $\int d\eta\ \eta^{n}$ in (A.3). We need the formula:

[TABLE]

For integer $n$, the hypergeometric function reduces to an arcsin function combined with square roots.

We thus obtain the following linear system of equations for moments $\{\tilde{k}_{r},\tilde{k}_{r+1/2}\}$

[TABLE]

where

[TABLE]

[TABLE]

Equations (LABEL:eqal) can be solved as an algebraic system of linear equations by truncating the infinite sum to some maximum value $n_{\rm max}$. Some of the coefficients involve $\Gamma$ functions and become large for large values of $r,n$. These coefficients are however multiplied by the moments $\tilde{k}_{r}$ or $\tilde{k}_{r-\frac{1}{2}}$, which for large $r$ are very small. It is in fact more convenient to solve the system by iterations, beginning with vanishing values for $\{\tilde{k}_{-1/2},\tilde{k}_{0},....,\tilde{k}_{N-1/2},\tilde{k}_{N}\}$, since the direct solution of the linear algebraic equations has to deal with a matrix with huge coefficients due to the above $\Gamma$ functions. The moments converge very rapidly after a few iterations.

Appendix B Galerkin Method for $u(\Delta)$ and $v(\Delta)$

As a consistency check of the numerical results for $u(\Delta)$ and $v(\Delta)$, we also made a numerical approximation of the integral equation (4.26) for $h_{1}$ by use of the Galerkin method. The equation was thus projected to the finite dimensional function space spanned by the following $6+8$ basis functions $f_{i}$, piecewisely defined on two subintervals:

[TABLE]

These basis functions were orthonormalized with respect to the scalar product $\langle\cdot,\cdot\rangle_{w}$ with weight function

[TABLE]

The Galerkin method then reduces the integral equation to a linear equation system for the best approximating coefficients in $h_{1}(\eta)\approx\sum_{i}c_{i}f_{i}(\eta)$:

[TABLE]

The integrals were computed numerically using interpolation of the kernel $G_{\infty}(\eta,\xi)$. The results are shown in figure 7, together with those of the numerical approximation of the moments presented in the main text.

Appendix C Living at the edge

Here we complement qualitative arguments of section 5 by a more rigorous analysis of the integral equation (4.26). It follows from the discussion in the main text that the non-analytic terms at $\Delta\rightarrow 0,1$ arise either from parametrically small $\xi$ (the $\ln\Delta$ terms) or from $\xi$ very close to one (the $\sqrt{\Delta}$ term). The expected singularities were estimated from the source term in the integral equation, disregarding higher-order iterations. This approximation cannot be really justified, and below we determine asymptotics of the exact solution near $\eta=0$ and $\eta=1$ without making any uncontrollable approximations. The end-point singularities of the source are captured by (5.5), but there is also an integral term in (4.26) which was neglected in the simplified calculation of section 5.

If we set $\eta=0$ in the kernel of the integral equation, nothing dramatic happens. Explicitly, from (4.13) we get:

[TABLE]

which is a nice, regular function on the whole interval $(0,1)$. The singularity at $\eta=0$, therefore, is entirely determined by the source. Introducing the variable $\epsilon$ from (5.1) we can recast (5.5) in the form valid for both $\epsilon>0$ ($\Delta\rightarrow 0$) and $\epsilon<0$ ($\Delta\rightarrow 1$):

[TABLE]

On the contrary, setting $\eta=1$ in (4.13) yields a function with a singularity at $\xi=0$:

[TABLE]

The integral term in the equation will thus get a singular contribution from $\xi$ very close to zero. But we already know the behavior of $h_{1}(\xi)$ at the left end of the interval, it is given by (C.2). Adding this integral contribution to the singular part of the source term we get777We need to impose a cutoff on the $\xi$-integration, because the approximate formula (C.2) is only valid for $\eta\ll 1$. We assume that $|\epsilon|\ll\ell\ll 1$.:

[TABLE]

The integral evaluates to elementary functions, but for our purposes the integral representation is more convenient.

Equipped with the asymptotic form of the solution at the two extremities of the interval, we can now compute the boundary contribution to the average (5.2). The source term at the right boundary gives $\sqrt{\epsilon}$ for $\epsilon>0$, as shown in (5.6). Having the integral term at our disposal we add the two contributions together. Assuming $\epsilon>0$, we get from (C.4):

[TABLE]

The square root has completely cancelled! For $\epsilon<0$ the source term is not singular and the integral term does not induce any new singularities either:

[TABLE]

The singular contribution from the other end was already computed in (5.6), (5.7):

[TABLE]

The integral term does not contribute, but this is not the end of the story, because iterations of the integral equation induce a small but non-analytic piece in the bulk:

[TABLE]

Indeed, the integral equation (4.26) can be written as

[TABLE]

Taking into account that

[TABLE]

and applying (C.7), we find that $\widehat{h}$ should satisfy the following integral equation:

[TABLE]

The equation can be solved with the help of the identity

[TABLE]

that can be proved by the same of chain of arguments as (4.15). The convolution integral here should be understood in the analytic sense as a contour integral around the cut, which removes an apparent divergence at $\xi=1$. The divergence arises because near $\xi=1$ the bulk solution is no longer accurate and has to be replaced with the asymptotic expression (C.4).

Substituting (C.12) into (C.11) we find:

[TABLE]

The total non-analytic piece of the average combines the contribution from the boundary (C.7) with the contribution from the bulk (C.8):

[TABLE]

The nominally divergent integral should again be understood in the analytic sense. In the main text we use an equivalent contour-intergal representation.