Integrable structure of Quantum Field Theory: Classical flat connections versus quantum stationary states

We establish a correspondence between an infinite set of special solutions of the (classical) modified sinh-Gordon equation and a set of stationary states in the finite-volume Hilbert space of the integrable 2D QFT invented by V.A. Fateev. The modified sinh-Gordon equation arise in this case as a zero-curvature condition for a class of multivalued connections of the punctured Riemann sphere, similarly to Hitchin's self-duality equations. The proposed correspondence between the classical and quantum integrable systems provides a powerful tool for deriving functional and integral equations which determine the full spectrum of local integrals of motion for massive QFT in a finite volume. Potential applications of our results to the problem of non-perturbative quantization of classically integrable non-linear sigma models are briefly discussed.


Introduction and summary
It is difficult to assign a precise mathematical meaning for the concept of integrability in Quantum Field Theory. A naive intuition goes back to Liouville of the 19 th century and suggests an existence of a sufficiently large set of mutually commuting operators whose joint spectra fully specify stationary states of the quantum system. For deeper insights, it is useful to consider 2D Conformal Field Theory (CFT), where significant simplifications occur due to the presence of an infinite dimensional algebra of (extended) conformal symmetry [1]. For a finite-size 2D CFT (with the spatial coordinate compactified on a circle of the circumference R), a mathematically satisfactory construction of an infinite set of mutually commuting local Integrals of Motion (IM) can be given and the simultaneous diagonalization of these operators turns out to be a well-defined problem within the representation theory of the associated conformal algebra.
Different conformal algebras, as well as different sets of mutually commuting local IM yield a variety of integrable structures in CFT. The series of works [2][3][4] was dedicated to the simplest of these structures, associated with the diagonalization of the local IM from the quantum KdV hierarchy [5][6][7][8]. Subsequent studies of this problem culminated in a rather surprising link between the integrable structures of CFT and spectral theory of Ordinary Differential Equations (ODE) [9][10][11]. In particular, in [11] a one-to-one correspondence was conjectured between the joint eigenbasis of the IM from the quantum KdV hierarchy and a certain class of differential operators of the second order −∂ 2 z + V L (z), with singular potentials V L (z) ("monster" potentials in terminology of [11]). Apart from a regular singularity at z = 0 and an irregular singular point at z = ∞, the monster potentials possess L regular singular points {x a } L a=1 . These potentials are not of much immediate interest in quantum mechanics, but arise rather naturally in the context of the theory of isomonodromic deformations. Solutions of the corresponding Schrödinger equations are single valued (monodromy-free) at z = x a and their monodromy properties turn out to be similar to that of the radial wave functions for the three-dimensional isotropic anharmonic oscillator. The monodromy-free condition was formulated in a form of the system of L algebraic equations imposed on the set {x a } L a=1 . 1 The correspondence proposed in [11] precisely relates the set of monster potentials V L (z) and the joint eigenbasis for all quantum KdV integrals of motion in the level L subspace of the highest weight representation of the Virasoro algebra. In particular, this implies that a number of the potentials V L (z) for a given value of L exactly coincides with a number of partitions p 1 (L) of the integer L into parts of one kind.
Since 1998, the link to the spectral theory of ODE have been extended to a large variety of integrable CFT structures (for review, see [14]), so that a natural question has emerged on whether a similar relation exist for massive integrable QFT. This question remained more or less dormant until the work [15], after which the so-called thermodynamic Bethe Ansatz equations have started to appear in different contexts of SUSY gauge theories [16][17][18][19]. These remarkable developments have led to the work [20], which established a link between eigenvalues of IM in the vacuum sector of the massive sine/sinh-Gordon model and some new spectral problem generalizing the one from [9,10].
This work is aimed to extend the results of [11,20] and provide an explicit example of the correspondence between stationary states of massive integrable QFT in a finite volume and singular differential operators of a certain class. At first glance, the best candidate for this purpose should be the sine-Gordon model, which always served as a basis for the development of integrable QFT. However, in spite of some technical complexity, a more general model introduced by Fateev [21] (which contains the sine-Gordon model as a particular case) turned out to be more appropriate for this task. The situation here is analogous to that in the Painlevé theory. Even though the Painlevé VI is the most complicated and general equation in the Painlevé classification, geometric structures behind this equation are much more transparent than those related to its degenerations. From this point of view, the fact that the sine-Gordon model is a certain degeneration of the Fateev model, could be understood as a QFT version of the relationship between the Painlevé VI and a particular case of Painlevé III.
The organization of this paper is as follows. In Section 2 we introduce the notion of Generalized Hypergeometric Opers (GHO's) -a special class of Fuchsian differential operators of the second order with 3 + L regular singular points at z = z 1 , z 2 , z 3 and z = x 1 , . . . x L . The variable z can be regarded as a complex coordinate on the Riemann sphere with 3 + L punctures. Projective transformations of z allows one to send the three points z i to any designated positions. At the same time other parameters of GHO are chosen in such a way that the remaining L regular singular points satisfy the monodromy-free condition. Therefore, the monodromy properties of GHO for L > 0 turn out to be similar to those for L = 0 (i.e. the ordinary hypergeometric differential operator of the second order). The complex numbers {x a } can be thought as local coordinates in the L-dimensional moduli space of GHO's. In Section 3 we consider more general differential operators, which inherit the monodromy-free property of GHO's. We call them the Perturbed Generalized Hypergeometric Opers (PGHO's). These operators have the form where P(z) = (z 3 − z 2 ) a 1 (z 1 − z 3 ) a 2 (z 2 − z 1 ) a 3 (z − z 1 ) 2−a 1 (z − z 2 ) 2−a 2 (z − z 3 ) 2−a 3 (1. 3) and the parameters 0 < a i < 2 satisfy the constraint a 1 + a 2 + a 3 = 2 . (1.4) Due to the last relation the quantity P(z)(dz) 2 transforms as a quadratic differential under PSL(2, C) transformations and the points z 1 , z 2 , z 3 on the Riemann sphere can still be sent to any given positions. In the presence of "perturbation" the monodromy properties of the operators (1.2) are changed dramatically in comparison with λ = 0 case. However, one can still find positions of the punctures x 1 , . . . x L , so that they remain monodromy-free singular points for any values of λ. In this case the coordinates {x i } L i=1 obey a system of L algebraic equations similar to that from [11][12][13]. Therefore, the moduli space of the PGHO's constitute a finite discrete subset A (L) in the moduli space of GHO's. 2 It appears that, for a given L, the cardinality of A (L) coincides with a number of partitions p 3 (L) of the integer L into parts of three kinds. In Sections 4-6 we interpret this fact in the spirit of [11] and present arguments in support of existence of a one-toone correspondence between elements of A (L) and the level-L common eigenbasis of the local IM in the integrable hierarchy introduced by Fateev in [21]. The arguments closely follow the line of [2-4, 10, 11] adapted to the algebra of extended conformal symmetry, which can be regarded as a quantum Hamiltonian reduction of the exceptional affine superalgebraD(2, 1; α) [22] (the "corner-brane" W -algebra, in terminology of [23]).
The generalization of the above constructions to the case of massive QFT is given in Sections 7-9. It is based on the idea from [20], which was inspired by the works [16][17][18][19]. As far as our attention has been confined to the case of CFT, there was no need to separately consider the antiholomorphic PGHO,D(λ) = −∂ 2 z +TL(z) +λ 2P (z), since there is only a nomenclature difference between the holomorphic and antiholomorphic cases. In massive QFT, following [20], one should substitute the pair of PGHO's (D(λ),D(λ)) by a pair of (2 × 2)-matrix valued differential operators where σ 3 , σ ± = (σ 1 ± iσ 2 )/2 are the standard Pauli matrices. In fact, (A z , Az) forms an sl (2) connection whose flatness is a necessary condition for the existence of solution of the linear problem The zero-curvature condition leads to the Modified Sinh-Gordon equation (MShG): ∂ z ∂zη − e 2η + ρ 4 P(z)P(z) e −2η = 0 , ρ 2 = λλ . (1.8) We consider a particular class of singular solutions of this equation, distinguished by special monodromy properties of the associated linear problem (1.7). The set of constraints imposed on these solutions is discussed in Section 7. In summary, e −η should be a smooth, single valued complex function without zeroes on the Riemann sphere with 3 + L +L punctures. Since z = ∞ is assumed to be a regular point on the sphere, then e −η ∼ |z| 2 as |z| → ∞ . (1.9) At the same time, e −η develops a singular behavior at z = z i , and also at z = x a andz =ȳ b e −η ∼z −x a z − x a (a = 1, . . . L) (1.11) The arbitrary parameters m i in the asymptotic formulae (1.10) should be restricted to the domains 3 3 At m i = 1 2 , 1 4 (2 − a i ) the leading asymptotic (1.10) involves logarithms. Here we ignore such subtleties.
whereas positions of the punctures (1.11) are constrained by a certain monodromy-free condition. The latter is now understood as a requirement that e ± 1 2 ησ 3 Ψ (where Ψ is a general solution of the auxiliary linear problem (1.7)) is single-valued in the neighborhood of the punctures z = x a (a = 1, . . . L) andz =ȳ b (b = 1, . . .L). Following the consideration from [13], the monodromy-free condition can be transformed into a set of L +L constraints imposed on the regular part of local expansions of (∂ z η, ∂zη) at the monodromy-free punctures: where γ a = ∂ z log P(z)| z=xa ,γ b = ∂z logP(z)|z =ȳ b . We expect that as far as of solutions satisfying all the above requirements. We can now define the moduli space A m which is the union of such finite sets: ( 1.15) Notice that, to a certain extent, A m can be regarded as a Hitchin moduli space [27]. An essential ingredient of the formal theory of the partial differential equation (1.8) is the existence of an infinite hierarchy of one-forms, which are closed by virtue of the equation (1.8) itself. In the case under consideration this formal property leads to the existence of an infinite set of conserved charges {q 2n−1 ,q 2n−1 } ∞ n=1 , which can be used to characterize the elements of the moduli space A m . The proof of this statement goes along the following lines. It easy to see that the flat connection A = A z dz + Az dz (1.6) associated with an element of A m is not single-valued on the punctured sphere. However, it does return to the original value after a continuation along the Pochhammer loop -the contour γ P depicted in Fig.1. Therefore one can consider the Wilson loop whose definition does not depend on the precise shape of the integration contour. In particular, it is not sensitive to deformations of γ P which sweep through the monodromy-free punctures. By construction the Wilson loop is an entire function of the spectral parameter θ, (1.17) Furthermore, since the shift of the argument θ → θ + iπ does not affect the connection A, the Wilson loop W = W (θ) is a periodic function of the period iπ. The textbook calculation [28] yields the following asymptotic expansions: stand for the constants that set a conventional multiplicative normalization (see Eqs.(7.13)-(7.16) bellow) for the conserved charges {q 2n−1 ,q 2n−1 } ∞ n=1 . The main result of this work is presented in Section 8, where we conjecture a correspondence between elements of the moduli space A m (1.15) and a subset H for the three scalar fields ϕ i = ϕ i (x, t). Here α i are coupling constants satisfying the constraint (1.20) In this work we focus on the case where α 2 i > 0. The parameter µ in the Lagrangian sets the mass scale, µ ∼ [ mass ]. We will consider the theory in a finite-size geometry (where the spatial coordinate x compactified on a circle of circumference R) with the periodic boundary conditions (1.21) Due to the periodicity of the potential term in (1.19) in ϕ i , the space of states H splits on the orthogonal subspaces H k 1 ,k 2 ,k 3 := H k characterized by the three "quasimomenta" k i ∈ [− 1 2 , 1 2 ]: Similar to the quantum mechanical problem of a particle in a periodic potential, the subspaces H k possess the band structure; they split into discrete components labeled by three integers: (1.23) The QFT (1.19) is integrable, in particular, it has an infinite set of commuting local integrals of motion I 2n−1 , with 2n = 2, 4, 6, . . . being the Lorentz spins of the associated local densities [21]: where ∂ ± = 1 2 (∂ x ∓ ∂ t ) and . . . stand for the terms involving higher derivatives of ϕ i , as well as the terms proportional to powers of µ. The constant C (n) ijk was found in [23] where (x) n stands for the Pochhammer symbol. Note, that the displayed terms in (1.24) with C It seems natural to expect that for 0 ≤ k i ≤ 1 2 , the sets of eigenvalues {I In the recent paper [29] it was argued that the vacuum eigenvalues {I and With the normalization conditions described above, the constants d n and E 0 reads explicitly as . (1.30) These relations should be supplemented by the identification of the parameters from the quantum and classical integrable problems: whereas the relation between dimensionless parameter µR and ρ is given by In this work we promote Eqs.(1.27)-(1.32) to a general relations between the joint spectra of the local IM in the subspace H (0) k corresponding the first Brillouin zone and the set of the conserved charges associated with the elements of the moduli space A m . For the values of k i restricted to the segment [ 0, 1 2 ], this gives a remarkable bijection between the joint eigenbasis of the local IM in H (0) k and the elements of A m . In Section 9 we demonstrate that the correspondence between the classical and quantum integrable systems provides a powerful tool for deriving integral equations which determine the full spectrum of local IM in the massive QFT.
We conclude this paper with few remarks concerning the QFT (1.19) in the regime where one of the couplings α i is pure imaginary.

Monodromies of the Fuchsian differential equations
In this preliminary subsection we include some basic concepts and results about the Fuchsian differential equations.
Let z stands for the complex coordinate on CP 1 \{z 1 , z 2 , . . . z n }, the Riemann sphere with n punctures. Consider the second order Fuchsian differential operator −∂ 2 z + T (z), where T (z) is given by The equation is a general second-order differential equation with n regular singular points. We will always regard the parameters δ i as fixed numbers. The positions of the singularities z i and the coefficients c i (which are usually referred to as the "accessory parameters") will be treated as variables. The accessory parameters c i are constrained by the elementary relations ensuring that T (z) has no additional singularity at z = ∞. Thus only n − 3 of these parameters are independent. Also, the projective transformations of the variable z allows one to send three of the points z i , say (z 1 , z 2 , z 3 ), to any designated positions, usually (0, 1, ∞). Therefore, with fixed δ i , the differential equation (2.2) essentially depends on 2 (n − 3) complex parameters. The equation (2.2) generates a monodromy group -a homomorphism of the fundamental group of the sphere with marked points into the group SL(2, C),  Figure 2: The elements Γ 1 , Γ 2 , . . . Γ n−3 of the fundamental group π 1 CP 1 \{z i } . Choosing the accessory parameters in (2.1) according to (2.11) fixes the conjugacy classes of the associated elements of the monodromy group of (2.2) as given in (2.9).
Let (ψ 1 (z), ψ 2 (z)) is a basis of linearly independent solutions of (2.2). Then its continuation along any closed path γ defines the monodromy matrix which depends only on the homotopy class of γ ∈ π 1 (CP 1 \{z i }). Let γ i ∈ π 1 (CP 1 \{z i }), i = 1, 2, . . . n be the elementary paths around the points z i , and the associated elements of the monodromy group of (2.2). The parameters determine the conjugacy classes of M (i) via the equation Let Γ a n−3 a=1 be the system of contours shown in Fig. 2, such that Γ a loops around the punctures z 1 , z 4 , . . . z 3+a only; and let the set ν = (ν 1 , . . . ν n−3 ) parameterize the conjugacy classes of the corresponding monodromy matrices M (Γ a ), For given conjugacy classes (2.9) (i.e., for a given set ν), the accessory parameters are determined in terms of the so-called classical conformal block f ν (X 1 , . . . X n−3 ) corresponding to "haircomb" diagram shown in Fig. 3 (for details, see e.g. [30,31]). Namely, where the shortcut notation F stands for with z ij := z i − z j and the arguments of the conformal block are cross ratios A certain additive normalization of the classical conformal block is usually assumed. For this reason, we reserve the room for an arbitrary coordinate-independent constant F 0 in (2.11).

Definition of GHO
Here we consider a special class of second order Fuchsian differential operators. We will always assume that the three parameters p i defining conjugacy classes of the matrices M (i) in (2.6)-(2.8), associated with the elementary paths around the "fixed" punctures z 1 , z 2 , z 3 , are positive numbers satisfying the following constraints For the remaining L := n − 3 punctures we require that both linearly independent solutions of (2.2) are single-valued (or monodromy-free) in the vicinity of these punctures. It is well known [32] how to reformulate this condition as a set of algebraic relations imposed on the corresponding parameters δ 3+a , C a := c 3+a and x a := z 3+a , (a = 1, . . . L). Namely, suppose that T (z) has a Laurent expansion at z = x a of the form 14) It is easy to see that l a must be an integer. We will focus on the case l a = 1, i.e. δ 3+a = −2 (a = 1, . . . L) . (2.15) To ensure that solutions of Eq.(2.2) are single-valued in the vicinity of the punctures at z = x a (a = 1, . . . L), the expansion coefficients t 1 in (2.14) should be constrained as The prime in T ′ (z) stands for the derivative w.r.t. the variable z. This system of algebraic equations should be supplemented by the three conditions (2.3), specialized to the case (2.15): As far as positions of the punctures z 1 , z 2 , z 3 and corresponding parameters p = (p 1 , p 2 , p 3 ) are fixed, Eqs.(2.17), (2.19) define an algebraic variety which will be denoted by V p , then a system of L locally defined functions C a (x 1 , . . . x L ) satisfy the integrability conditions These relations can be verified by the brute-force calculation using Eqs.(2.17) and (2.19), but, of course, they follows immediately from the general relation (2.10). In this particular case the classical conformal block in Eq.(2.11) is related to the classical limit of the (3 + L)-point correlator involving three generic chiral vertex operators V i with conformal dimensions ∆ i (i = 1, 2, 3) and L degenerate vertices V (3,1) with dimensions ∆ (3,1) : where the parameter b 2 and other conventional notations are inherited from the quantum Liouville theory (see Ref. [33] for details) Due to the well known fusion rule for the degenerate vertex V (3,1) [34], only a discrete set of the parameters ν = (ν 1 , . . . ν L ) is allowed (see Fig. 4): In what follows we will refer to the differential operators (1.1), whose moduli space coincides with the algebraic variety V (L) p as to the Generalized Hypergeometric Opers (GHO's). 4 The marked points x a (a = 1, . . . L) will be called as monodromy-free punctures.
on the sphere with three punctures.

Connection matrices for GHO
We have introduced the concept of GHO because the monodromy group of such opers coincides with the monodromy group of the conventional hypergeometric equation. Let us recall some facts about this group. In the case under consideration there are only three elementary SL(2, C)matrices M (i) , M (j) and M (k) (2.8), corresponding to the contours γ i , γ j and γ k , shown in Fig. 5. Here (i, j, k) is any cyclic permutation of (1, 2, 3). These matrices satisfy an obvious relation because the contour γ i • γ k • γ j is a contractible loop. Further, since 0 < p i < 1 2 , one of these matrices, say, M (i) can always be chosen diagonal (2.25) Here and below we use standard notation for the Pauli matrices σ 1 , σ 2 , σ 3 . Then Eqs.(2.8), (2.24), (2.25) define M (j) and M (k) up to a diagonal similarity transformation. In particular, The quantity ω i is an arbitrary complex number and The prefactor here is chosen to satisfy the normalization condition where W[f, g] = f g ′ − gf ′ stands for the Wronskian. If the constraints (2.13) are imposed, the asymptotic conditions (2.29) define 5 three different bases (for i = 1, 2, 3) in the two-dimensional linear space of solutions of (2.2). Let us combine the solutions (2.29) for given i into the row Then the two sets of basis vectors χ (i) and χ (j) are related through a linear transformation Here S (j,i) stand for SL(2, C)-matrices, satisfying the relations det S (j,i) = 1 , where again (i, j, k) is any cyclic permutation of (1, 2, 3). It is easy to see that one needs six independent complex numbers to parameterize the three matrices S (i+1,i) (i ∼ i + 3) satisfying (2.33). Moreover, these connection matrices are subject to three additional complex constraints. Indeed, the monodromy matrix M (j) (2.26) can be expressed in term of the connection matrix S (j,i) : This relation combined with Eqs.(2.33) leads to with where the matrix elements read explicitly (2.37) 5 It is worth noting, however, that (2.29) define these solutions only up to phase factors of the form ±e 2πipiM (M ∈ Z). and Λ is given by (2.27). Note that Eq.(2.35) can be equivalently rewritten as a formula for the Wronskians: The complex parameters ω i (i = 1, 2, 3), entering the expression (2.35), remain undetermined. These parameters do not affect the conjugacy class of the representation of π 1 CP 1 \{z 1 , z 2 , z 3 } . Nevertheless, they are important characteristics of the GHO itself. In the next section we argue that, in the case of GHO, a coordinate-independent additive normalization of the function F (2.11) can be chosen in such a way that An immediate consequence of this fact is that Finally, note that exp(ω i ) is defined up to the phase factor ±e 2πipM (M ∈ Z). This ambiguity is inherited from the similar ambiguity in the definition of χ

GHO and complex solutions of the Liouville equation
Until now we have discussed holomorphic GHO only. Of course, with minor modifications all the above can applied to the antiholomorphic GHŌ In what follows we assume that the triple , the correspondinḡ δ i = δ i are real and, furthermore, where G (i) is an arbitrary 2 × 2 matrix and superscript T stands for the matrix transposition. We specialize G (i) by the requirement that τ (z,z) is a single-valued function on the punctured sphere. Imposing this condition in the vicinity of the puncture z i , one concludes that G (i) is a diagonal matrix. With the connection formula (2.32), the single-valuedness implies that S (ji) G (i) S (ji) T is also a diagonal matrix. Using the explicit form of connection matrices (2.35)-(2.37), one finds If the undetermined constant is chosen to be ±1, then the complex function η : At the monodromy-free punctures, z = x a (a = 1, . . . L) andz =ȳ b (b = 1, . . .L), e −η becomes singular, where The constants η (reg) i can be regarded as regularized values of the Liouville field at the punctures z i . They are be expressed in terms of the regularized Liouville action [33,36]. To explain this important relation we recall that the Liouville equation and the asymptotic conditions (2.50) follow from the variational principle for the functional (2.52) Since the fields configuration is singular at z → z i , we cut out a small disk of radius ǫ i around the point z i and add the boundary terms with to ensure the behavior (2.50) near z i . To control the large |z|-behavior we regularize the action for large values of z and add the boundary term with In addition, we include some field independent terms such that A Liouv is finite and independent on ǫ i and R when ǫ i → 0, R → ∞. Contributions of the monodromy-free punctures (2.49) to the functional (2.52) are finite 6 and therefore there is no need to include additional regularization terms to the action. It is now easy to show that where A * Liouv stands for the stationary value of the functional (2.52) calculated on the field configuration η defined by Eq.(2.48). More generally, the stationary value of the Liouville functional depends on 6 + L +L variables, and its total deferential is given by [33,36] dA * Liouv = 2 As it follows from (2.10), A * Liouv can be expressed in terms of F and its antiholomorphic counterpart F : The coordinate-independent constant F 0 in Eq.(2.11) has not yet been fixed. Therefore there is no need to add a p i -dependent constant in (2.58), as it can always be absorbed by F 0 andF 0 . The number and positions of the holomorphic and antiholomorphic monodromy-free punctures are fully independent. For instance, we can consider the general holomorphic GHO, whereas the antiholomorphic differential operator (

Definition of PGHO
Consider the universal cover of the Riemann sphere with three marked points z 1 , z 2 , z 3 and P(z) (dz) 2 , where P(z) is given by (1.3) with positive parameters a i . If a i satisfy the constraint (1.4), the quantity P(z) (dz) 2 transforms as a quadratic differential under PSL(2, C) transformations and the punctures z 1 , z 2 , z 3 on the Riemann sphere can still be sent to any desirable positions. Suppose we are also given a GHO, D = −∂ 2 z + T L (z) which has its first three punctures at the branching points of P(z) plus L monodromy-free punctures at z = x a (a = 1, . . . L). Remind, that previously we have required that the parameters p i obey the constraints (2.13). In what follows we will impose somewhat stronger constraints on these parameters. Namely we replace (2.13) by The rôle of this constraint will be explained in Section 5.1 bellow. An immediate object of our interest is an ODE of the form where λ stands for an arbitrary complex parameter. The properties of the differential equation (3.2) are essentially affected by the presence of the λ-dependent term. Nevertheless, one can still find particular values of its parameters to make the marked points z = x a (a = 1, . . . L) to be monodromy-free punctures for arbitrary values of λ. Indeed, the conditions (2.16) can be easily generalized for λ = 0. In this case the system of algebraic system (2 which determine the values of x a (a = 1, . . . L). These algebraic equations have a finite discrete set of solutions. Therefore, for any given L, there only a finite number N L of sets of monodromy-free punctures Notice that Eqs.(2.17)-(2.19), (3.3) are symmetric upon permutations of (x 1 , . . . x L ), therefore we will not distinguish sets, which differ only by a permutation of the positions of the punctures. Let us illustrate the situation on the simplest L = 1 example. As in Appendix A, we set (z 1 , z 2 , z 2 ) = (0, 1, ∞), so that together with (A.1) one has an additional relation This leads to a cubic equation for the position x of the monodromy-free puncture: In many extents these equations are similar to the Bethe Ansatz equations. In particular, using Eq.(2.10), they can be written in a compact Yang-Yang form Once the algebraic system (3.8) is solved, the function T L (z) in (3.2) can be written in the form and We will refer to the differential operator D(λ) of the form (3.2) with P(z) and T L (z) are given by (1.3) and (3.10), respectively, as Perturbed Generalized Hypergeometric Oper (PGHO). The finite set A (L) p (3.4) can be regarded as a moduli space of the PGHO's. It is a finite discrete subset in the moduli space of GHO's. (Notice that we slightly modify the notation used in the introduction by including the subscript p = (p 1 , p 2 , p 3 ).)

Wilson loop for PGHO
As we have just explained, the position of the punctures x a (a = 1, . . . L) can be specially chosen so that solutions of ODE (3.2) still remain single-valued in the vicinity of these points. However, contrary to T L (z), the term λ 2 P(z) is not single-valued on the punctured sphere. Thus, even with the special choice of x a , the monodromy group of the differential operator D(λ) turns out to be essentially different from that in the case λ = 0. Here we begin to explore the monodromy properties of PGHO.

Definition of the Wilson loop
Let us consider the contour γ P depicted in Fig. 1. It is usually called the Pochhammer contour (loop). As an element of the fundamental group π 1 (CP 1 \{z 1 , z 2 , z 3 }), it can be expressed in terms of the elementary loops γ i , γ j and γ k which wind around the punctures z i , z j and z k , respectively: Since the Pochhammer loop winds around each of the three punctures and the relation (1.3) is imposed, the value of the function P(z) does not change upon the analytic continuation along the contour γ P . Therefore the coefficients of PGHO return to their original values and it makes sense to introduce the quantity where M (γ P |λ) is the monodromy matrix for D(λ) corresponding to the Pochhammer loop. A significant advantage of W(λ) is that it does not depend on the precise shape of the integration contour. In particular, it is not sensitive to deformations of γ P which sweeps through the monodromy-free punctures. In what follows we will refer to (3.14) as the Wilson loop (corresponding to PGHO D(λ)). The second order differential operator D(λ) depends analytically on λ 2 and hence W(λ) is an entire function of λ 2 , i.e., the series expansion converges for any complex λ. Its value at λ = 0 can be found using Eqs.(2.25), (2.26) from Section 2.3: Note, that this expression does not depend on the number of the monodromy-free punctures L. Higher expansion coefficients in the series (3.15) can, in principle, be calculated using the standard perturbation theory.

Large-λ asymptotic expansion
The leading large-λ asymptotic of the Wilson loop can be obtained within the WKB approach. It is easy to see that Here the r.h.s. is written as a sum of two WKB exponents. Of course, for different values of arg(λ 2 ) only one term dominates whereas another exponent should be neglected. The quantity P(z) is a multivalued function on CP 1 \{z 1 , z 2 , z 3 } whose phase has not been yet uniquely specified. To resolve this phase ambiguity, we consider the Möbius transformation which sends (z 1 , z 2 , z 2 ) to (0, 1, ∞). With this change of the integration variable, the integral in (3.17) transform to the form The Pochhammer contour now looks as in Fig. 6. Let us choose the base point z * ∈γ P within the real segment [0, 1] and assume that z contour. This convention removes the phase ambiguity of P(z) for z ∈ γ P . The integral which appears in the r.h.s. of (3.18) is well known in the theory of the hypergeometric equation: We can now rewrite Eq.(3.17) in the form .

(3.21)
It is not difficult to extend the above leading asymptotics to a complete asymptotic expansion for large values of λ. For this purposes, we perform the change of variables in ODE (3.2) This transformation brings Eq.(3.2) to the form of an ordinary Schrödinger equation with the potentialT It is well known how to develop the large-λ asymptotic expansion of monodromy coefficients of Eq.(3.23). The procedure leads to the following asymptotic series In the last formula U n [û ] are homogeneous (grade(û) = 2, grade(∂) = 1, grade(U n ) = 2n) differential polynomials inû of the degree n (known as the Gel'fand-Dikii polynomials [37]), and prime stands for the derivative. Thus, where the last line shows the overall normalization of the polynomials. There is no need to do describe the contourγ P in (3.27) explictly, since we now change the integration variable in Eq. (3.27) back to the original coordinate z. In this way one obtains (3.31) (Notice that in the derivation of the second formula we dropped the term − 1 3 ∂ 2 wT L in Eq.(3.27) with n = 2, which do not contribute to the integral.) Of course, it is straightforward to perform the change of variables in Eq.(3.27) for any given n. We do not present explicit formulae for n > 2, but note that where the omitted terms contain derivatives and the lower powers of T L .

Expansion coefficients q
Using the formulae (3.31), (3.32), one can perform some explicit calculations of the coefficients in the asymptotic series (3.25). Let us first consider of the perturbed hypergeometric oper, i.e., PGHO without monodromy-free punctures. Using Eqs.(3.10), (3.31) (specialized for the case L = 0) and the integral (3.19), one can show that and and where the numerical coefficients E i , E ij and H i are given by The indices (i, j, k) represents any permutation of the numbers (1, 2, 3). For n > 2 the calculation of q (0) 2n−1 is straightforward, but rather long. It is much easy to establish the following general structure: where R n stands for n-the degree polynomials in the variables P 2 (the dots represent the sum of monomials of degrees lower than n). One can show that . for arbitrary L. Indeed, the third term in (3.10) do not contribute to the integral (3.31) for q (L) 1 . The contribution of the first term in (3.10) is given by (3.33). The second term in (3.10) gives a contribution proportional to L. The final result reads as The calculation of q (L) 3 is very cumbersome and we do not describe it here. Bellow we quote the result which is expressed in terms of the parameters δ i , a i (i = 1, 2, 3) and the coordinates x a (a = 1, . . . L) of the monodromy-punctures. Also, it is assumed that (z 1 , z 2 , z 3 ) = (0, 1, ∞);

Hidden algebraic structures behind PGHO
We have already mentioned a remarkable property of the algebraic system (2.17)-(2.19), (3.8). Our numerical work shows that for given L the number N L of its solutions (i.e., the cardinality of the set A  .7). In this pattern one can recognize the first values for the number of partitions of the integer L into integer parts of three kinds, which we denote as p 3 (L). This sequence is generated by the series We now interrupt our formal study to discuss remarkable algebraic structures behind PGHO. Introduce the three-component chiral Bose field φ = (φ 1 , φ 2 , φ 3 ), i.e. the operator valued function where Q i , P i and a i (n) are operators satisfying the commutation relations of the Heisenberg algebra Let P s+1 (∂φ, ∂ 2 φ, . . .) be a local field of spin s + 1, which is a local polynomial of ∂φ and its higher derivatives ( ∂ stands for ∂ ∂u here). All such fields are periodic in u, therefore one can introduce the integral, Bellow the shortcut notation I s for I[P s+1 ] is used. Suppose we are given a special infinite sequence of operators I s (corresponding to special infinite sequence of the polynomials P s+1 ) which are mutually commutative operators, We will refer to the operators {I s } as the (chiral) local Integral of Motions (IM). A complete algebraic classification of all possible infinite sets of local IM seems to be a hopeless task. However some non trivial examples are available. Among them there is a two-parameter family discovered by Fateev in [21]. The first two representatives from this set are given by and Numerical coefficients in the last formula depends on three parameters α 1 , α 2 and α 3 obeying the quadratic constraint α 2 1 + α 2 2 + α 2 3 = 1 2 . It turns out that E j , E mj and H j are given by Eqs. (3.36), provided the parameters are identified as whereas An explicit form for the higher spin representatives is not available. However it is known that [23] where the dots stand for the terms, involving higher derivatives of φ i and the constant C (n) ijk is given by (1.25). There are good reasons to expect (see Ref. [23] and Section 6.1 bellow) that I 1 and I 3 are just the first two representatives of an infinite two-parameter family of mutually commuting IM, {I 2n−1 } ∞ n=1 . Let F P with P = (P 1 , P 2 , P 3 ) be the Fock space, i.e., the space generated by the action of a i (n) with n < 0 on the vacuum state | P which satisfies the equations The space F P naturally splits into the sum of finite dimensional "level subspaces" The dimensions of the level subspaces do not depends on P. Obviously, it coincides with the number of integer partitions of L into parts of three kinds, defined in (4.1), (4.14) The grading operator L essentially coincides with I 1 (4.6): Therefore all local IM from the Fateev family act invariantly in the level subspaces F (L) P . The diagonalization of I 2n−1 in a given level subspaces reduces to a finite-dimensional matrix problem which however rapidly becomes very complex for higher levels.
Of course, the highest weight vector of the Fock space (the "vacuum" vector) is an eigenvector for all integrals of motion I 2n−1 . Let I (0) 2n−1 be the corresponding eigenvalues. The results from Section 3.2.3 and Eqs.(4.6), (4.7) imply that for n = 1 and n = 2 the following relation holds where the parameters a i and p i of q 2n−1 are related to α i and the zero mode momentum P i as in Eqs. (4.8) and (3.35), respectively. Moreover, for any value of n both sides of (4.16) are polynomials in the variables (P 1 ) 2 , (P 2 ) 2 , (P 3 ) 2 of the degree n. Comparing Eqs.(3.37), (3.38) with (4.10), (1.25), it is easy to check that all leading n-th degree monomials are exactly the same in the both sides. Thus one can reasonably expect that (4.16), involving the vacuum eigenvalues of the integral of motion and the expansion coefficients of the Wilson loop for the PGHO with L = 0 holds exactly for any value of n ≥ 1.
Actually, we expect that (4.16) can be extended to the relation between the whole spectrum of I 2n−1 in any level subspace F P . We expect that, up to the overall normalization factor, the set q for any fixed n and L. Thus the subscripts α and β can be identified and Eq.(4.16) is generalized as follows:

Connection matrices for PGHO
In the previous sections we have discovered interesting properties of the Wilson loop (defined in (3.14)) by studying its asymptotic expansions at large values of λ, using the WKB approximation. Even though this asymptotic analysis has led to remarkable insights into the algebraic structure of the problem, considered in Section 4, it does not solve the mathematical problem of an exact calculation of the Wilson loop as entire functions of the variable λ 2 . In this section we address this problem. Actually, here we solve a more general problem of an exact calculation of all connection matrices for the PGHO. By doing this we employ and extend ideas and methods previously developed in [2][3][4]10]. The matrix elements of the connection matrices are entire functions of λ 2 . Additional information about their analytic properties, namely, asymptotic distributions of their zeroes, is deduced from the standard WKB analysis. We use various symmetries of the differential operator (3.2) and derive a system of functional relations, which allows one to completely determine all the connection matrices. Interestingly, these functional relations have only a discrete (albeit infinite) set of solutions, which possess the required analytic properties. We conjecture that these solutions precisely correspond to PGHO's with an arbitrary number of monodromy-free punctures. The results are supported by several analytical and numerical checks for PGHO with L = 0.

Functional relations for the connection matrices
A proper definition of the bases of solutions (2.29) for λ = 0 requires some additional considerations. First of all, one needs to take into account that (unlike the λ = 0 case) analytic continuations along infinitesimal loops around the singular points z 1 , z 2 and z 3 affects the PGHO itself. Therefore, in order to define solutions by asymptotic conditions at these points one needs to make suitable brunch cuts. Let us chose an extra point, say z = ∞, and cut the Riemann sphere along the lines, connecting this point with the branching points of P(z). In Fig. 7 these cuts are shown by the dashed lines. Next, the asymptotic conditions (2.29) must be slightly modified since the order of the correction term is changed with respect to that in the λ = 0 case. The above conditions uniquely define the solutions χ (i) σ (z) provided that the parameters p i satisfy an additional constraints 0 < p i < a i 4 , which were already enforced in Eq.(3.1) above. The connection matrices S (j,i) (λ) for λ = 0 can be defined in the same way (2.32) as in the case of unperturbed GHO: They satisfy the same relations (2.33) as for λ = 0: Throughout this section we assume that (i, j, k) is a cyclic permutation of (1, 2, 3). In Fig. 7 the matrices S (j,i) (λ) are associated with the oriented lines connecting the points z i and z j . The main rôle in the following analysis belongs to symmetry transformations which essentially allows one to connect solutions χ (i) on different sheets of the Riemann surface of the PGHO. Let be a transformation, involving a translation of the independent variable z along the contour γ i , accompanied by the substitution λ → q −1 i λ, where q i = e iπa i , q 1 q 2 q 3 = 1 .

(5.5)
It is easy to see, that the substitutions (5.4) leave PGHO unchanged. Therefore they act as linear transformations in the space of solutions. Namely, in the basis χ (i) they read The most fundamental property of the differential operator (3.2) is that a combined transformation Ω k • Ω j • Ω i , where (i, j, k) is a cyclic permutation of (1, 2, 3), is equivalent to the identity transformation in the space of solutions of (3 .2), The proof follows from the relation (5.5) and the fact that γ k • γ j • γ i is a contractible contour, which loops around a regular point (z = ∞) of the PGHO (see Fig. 5). Combining (5.6) and (5.7) with the definition (5.2) one easily obtains Consider now the transformation Ω k • Ω i • Ω j , where the indices i and j are interchanged with respect to (5.7). Repeating the above arguments (again with an account of (5.5)) one can show that this transformation is equivalent to a linear transformation of solutions where M (γ P |λ) can be interpreted as a monodromy matrix of the Pochhammer loop depicted in Fig. 1. Then using (5.6) one obtains, We would like to stress that the above considerations apply to all PGHO's with an arbitrary number of the monodromy-free punctures. This means that the connection matrices will always satisfy the same relations ( In what follows the set of functions W (k) σ ′ σ from Eq.(2.37). As well as the Wilson loop, the connection coefficients are entire functions of the variable λ 2 , i.e. they can be represented by power series in λ 2 with infinite radius of convergence. Our next goal is to describe their characteristic properties.

Large-λ asymptotic
Consider the large λ behavior of the connection coefficients. In the leading order one has where the integrals taken along the oriented links depicted in Fig. 7. Introduce the constants r k > 0 and b k : (5.14) Then Eq.(5.13) can be equivalently written in the form Note that, as it follows from the definition (5.14), the positive constant r k is given by To assign precise meaning to an individual phase factor e −iπb k in Eq.(5.15), one needs to resolve the overall phase ambiguity of P(z). Following the procedure from Section 3.2.2 we send (z 1 , z 2 , z 3 ) to (0, 1, ∞), then In fact, it is not difficult to calculate explicitly the subleading term in the asymptotic formula (5.15). In order to simplify formulae bellow we make use the notation σ ′ σ is given by Eq.(2.37). Then The above formula can be applied for large λ 2 such that | arg(λ 2 )| < π. In the case of real λ 2 < 0, i.e. when λ = i |λ|, the asymptotic is given by As it was discussed at the end of Section 2.3, the combinations (2.40), which appear in the formula (5.22), are functions of the L projective invariants. In the case L = 0, they are given by equation (B.2) from Appendix B. For this reason it is convenient to write the subleading terms in the asymptotic formulae (5.21) and (5.23) as and g (L,α) i (g (0) i = 1) stand for λ-independent constants corresponding to a given set of monodromyfree punctures (3.4). 7 The asymptotic formula (5.21) can be extended to the following systematic asymptotic series Here the quantity B (k) (λ) is a formal power series In the case L = 1, the formula (B.4) from Appendix B leads to Here we have omitted the indices σ ′ , σ in the notation of zeroes. This dependence will be implicitly assumed. We will also assume that the sign of λ (k) n is fixed by the requirement −π/2 < arg λ Due to the cyclic symmetry, it is sufficient to consider one value of k, say k = 3. Let us set (z 1 , z 2 , z 3 ) to (0, 1, ∞) and then make the change of variables (3.22). The transformation w = w(z) is the Schwartz-Christoffel mapping which sends (0, 1, ∞) to (0, r 3 , r 2 e iπ 2 a 1 ), whereas the functionψ satisfies the ordinary Schrödinger equation (3.23) with the potentialT L (w) given by (3.24). Consider a zero λ If the parameters p i are restricted by the condition 0 < 2p i a i < 1 2 (see Eq.(3.1)), the above asymptotic conditions lead to well-defined spectral problems for all σ, σ ′ = ±1. An immediate consequence of this fact is that all the zeroes of A σ ′ σ (λ) are simple. In the simplest case of the perturbed hypergeometric oper (i.e., for L = 0), the potential in the Schrödinger equation (3.23) is real and positive. Therefore all the zeroes λ (3) n are also real and positive. Then the large-λ asymptotic formulae (5.21)-(5.23) imply that the zeroes accumulate at the infinity along the positive real axis and for large integer n ≫ 1 one has (Because of cyclic symmetry, the last formula is valid for any cyclic permutation (i, j, k).) For a general case of PGHO with L > 0 the potential in the equation (3.23), in general, becomes complex-valued for w ∈ [0, r 3 ], so that the zeros λ n also become complex. However, they still 8 An explicit form of the expansion coefficients x (i) σ,n are not known. The sole exclusion is the first coefficient in the case of PGHO with L = 0, which reads explicitly remain simple and accumulate at infinity in the vicinity of the positive real axis. The asymptotic formula (5.31) continues to hold for L > 0. Moreover, we would like to stress, that for large n this formula gives the asymptotics of precisely the n-th zero λ (k) n (in the sense that n coincides with the number of zeroes, whose absolute value is less or equal than |λ (k) n |). Similar considerations apply to all functions A (k) σ ′ σ (λ); they can be written in the form of convergent products σ ′ σ is given by (2.37). At this stage it is convenient to introduce spectral ζ-functions, which capture all information about the distribution of zeroes λ (recall that we assume that −π/2 < arg λ (k) n ≤ π/2). As follows from the asymptotic formula (5.31) the function ζ k (ν) is analytic in the lower half plane ℑm(ν) ≤ 0, except the point ν = −i, where it has a simple pole with the residue −i r k /π. Using these properties the product formula (5.32) can be transformed into an integral representation Closing the integration contour in this formula in upper half plane and comparing the result with the asymptotic expansion (5.26) one concludes that the function ζ k (ν) has zeroes at ν = 2i, 4i, . . . and additional simple poles on the imaginary axis ν with the following residues res ζ k (ν) ν=i(2n−1) = i Γ(n + 1 2 ) q where n = 1, 2, . . . . Moreover it follows from (5.21) that where Λ i (λ) is defined in (5.22).

Bethe Ansatz equations
The nine non-linear functional equations (5.3) and (5.8) involve too many unknown functions (twelve) and, in fact, appear to be rather complicated for a direct analysis. Fortunately, it is possible to reduce these equations to eight sets of rather compact equations of the Bethe Ansatz type, where each set involves only three unknown functions. In principle, this could be done by direct manipulations with the equations (5.3) and (5.8), but here we prefer a more efficient approach involving direct calculations of the Wronskians. It is based on the relation (5.7) and the following simple properties: (i) The solutions χ (i) (z) are simply transformed under the action of Ω i with the same i (see Eq.(5.6a)). This follows from the fact that the asymptotic condition (2.29), defining the solution χ (ii) For any two solutions ψ 1 (z) and ψ 2 (z) of (3.2) one has The proof follows from the definition (5.4) and the fact that the Wronskian in (5.38) does not depend on the point z where it is calculated.
Let (i, j, k) be a cyclic permutation of (1, 2, 3) and σ, σ ′ , σ ′′ = ±1. Together with functions A The same quantity can be calculated in a different way, using the additional relations (5.7) and (5.38). First, from (5.6a) and (5.38) it follows that Similarly, using also the property (5.7), one can easily show that Making simultaneous cyclic permutations of the indices (i, j, k) and the values (σ, σ ′ , σ ′′ ) one obtains another two equations of the same type, which contain the same three functions A σσ ′′ (λ) and A (k) σ ′ σ (λ) as in the equation (5.44). By definition, T (i) σ (λ) is an entire function of λ 2 , therefore the l.h.s. of (5.44) vanishes at all zeroes of A (i) σ ′′ ,σ ′ (λ). Proceeding in this way one obtains a system of three coupled Bethe Ansatz type equations for the position of zeroes where n = 1, 2, . . . and λ σ ′ σ (λ) enters into the two sets of these equations. As an immediate consequence one can derive an "asymptotically exact" Bohr-Sommerfeld quantization condition for the roots λ It is convenient to introduce a new function and another two functions ǫ j (λ) and ǫ k (λ), which are obtained from (5.48) by simultaneous cyclic permutations of the indices (i, j, k) and (σ, σ ′ , σ ′′ ). To simplify the following equations we have omitted the indices σ, σ ′ , σ ′′ in the notation of ǫ-functions since their arrangement is firmly connected to the indices (i, j, k) and so that they can always be restored. We expect that the Bethe Ansatz equations (5.45) combined with the asymptotic formula (5.31) have an infinite number of solutions, corresponding to PGHO's with different configurations of monodromy-free punctures. These solutions are distinguished by different phase assignments in the logarithmic form of the Bethe Ansatz equations (5.45),  Further, although we have previously assumed that the parameters p i , p j , p k obey the constraints (3.1), the resulting Bethe Ansatz equations (5.45) make sense for any complex values of p i . Most importantly, their solutions continuously depend on these parameters. Below we will use this fact to enumerate all solutions of (5.45), following the line of Appendix A of Ref. [11]. Fix the values σ, σ ′ and σ ′′ and assume that and that |p i |, |p j | and |p k | are of the same order of magnitude. Then the asymptotics (5.31) (as well as the numerical analysis of (5.49)) suggests that for sufficiently large values of the parameters (5.51) all roots λ Here We conjecture that the solutions of the Bethe Ansatz equations (5.49), associated with such set of integers {µ (i) } and {μ (i) } correspond to PGHO's with exactly monodromy-free punctures. For a given value of L the number of the integer sets μ (1) , µ (1) ,μ (2) , µ (2) ,μ (3) , µ (3) , satisfying this equation is equal to p 3 (L) (which is the number of partitions of L into integer parts of three kinds, already defined in (4.1)).

Non-linear integral equations for L = 0
The entire function A σ ′ σ (λ) from this data is similar to the one which emerged long ago in the context of the analytic Bethe Ansatz [38][39][40][41]. For the sine-Gordon model the problem was solved by Destri and De Vega [42,43], who have reduced it to a single complex non-linear integral equation. Similar equation was earlier derived in the lattice XXZ-model in Ref. [44]. Here we consider the the non-linear integral equations determining A (k) σ ′ σ (λ) in the simplest case of PGHO without monodromy-free punctures, i.e, L = 0.
Using (5.34) define spectral ζ-functions ζ i (ν), ζ j (ν) and ζ k , associated with A σσ ′′ (λ) and A (k) σ ′ σ (λ), respectively. It is convenient to introduce a new variable θ = log(λ). The Bethe Ansatz equations (5.45) allows one to derive a non-trivial relation between ǫ-and ζ-functions. For the case when all roots lie on the positive real axis, it reads (see [42,43] for details of a similar derivation) The integral (5.56) converges in the half plane ℑm(ν) > 0, but it can be analytically continued to the whole complex plane of ν. In fact, as it was remarked before, the function ζ i (ν) is analytic in the lower plane ℑm(ν) ≤ 0 except a simple pole at ν = −i. Combining the relations (5.34), (5.48) and (5.56) it is easy to show that where we use function which is analytic in the upper half-plane, ℑm(ν) > 0. The function f i (ν) has a simple pole at ν = 0, In a view of Eqs.(5.37) and (5.56), it is easy to see that The equation (5.58) has been solved numerically for various values of the parameters a 1 , a 2 , a 3 and p 1 , p 2 , p 3 . Using the obtained numerical values of ǫ i (θ) we calculated (5.60) for n = 1, 2 and (5.61) for n = 1 and checked that they are in an excellent agreement with Eqs.

Hidden algebraic structures (continuation)
In a view of identification (4.17), the formal power series B (k) (λ) in the asymptotic formula (5.26) can be understood as eigenvalues of the formal operator in the Fock space F P with P i = 2p i √ a i . (Here and below, we always assume that the parameters α i and a i are related as in (4.8).) In fact, all other terms in (5.26) can be also understood as eigenvalues of certain operators commuting with the local IM.

Corner-brane W -algebra and reflection operators
Here we argue that the factor g (L,α) j σ ′ g (L,α) i σ in (5.24) can be identified with an eigenvalue of certain λ-independent operator R (k) σ ′ σ acting in the Fock F P space and commuting the local IM I 2m−1 : The operators R (k) σ ′ σ are similar to the reflection operator from Ref. [33]. The main part in the construction belong to a W -algebra whose rôle is analogous to that of the Virasoro algebra in the quantum Liouville theory. This W -algebra was introduced in Ref. [47] and studied in Ref. [22]. Bellow we closely follow the consideration from Ref. [23], where this W -algebra was called "corner-brane" W -algebra.
Let us introduce four vectors and define the exponential vertex operators is the three-component chiral Bose field (4.2) and the dot product stands for x · y = 3 i=1 x i y i . We now choose the first three vectors α 1 , α 2 and α 3 from the set (6.3) and define the algebra W (1,2,3) as an algebra generated by the holomorphic currents W s of spin s characterized by the condition that they commute with three "screening charges" The integration here is taken over a small contour around the point u. For small s the condition (6.5) can be straightforwardly analyzed. In particular, one can show that spin-1 currents satisfying (6.5) are absent, but there is one (up to an overall multiplier) spin-2 current which generate the Virasoro subalgebra with the central charge Furthermore, there are no non-trivial spin-3 currents since the spin-3 fields satisfying (6.5) turns to be the derivative ∂W 2 . For spin-4 there are three fields -two "descendent" currents ∂ 2 W 2 and (W 2 ) 2 , but also one new current W 4 . Explicit form of W 4 is somewhat cumbersome and can be found in Appendix A of Ref. [23]. For s > 4, the calculations based on definition (6.5) become very complicated. However one can argue (see Ref. [23]) that there is exactly one independent current W 2n at each even spin s = 2n, having the form where the non-derivative term W (sym) 2n (but not V 2n−1 ) is symmetric with respect to all 180 o rotation around the coordinate axes of the (φ 1 , φ 2 , φ 3 ) space: (6.10) The above construction can be repeated for any choice of three vectors α A , α B and α C from the set (6.3) to yield four corner-brane W -algebras which are labeled by are triple integers (A, B, C): 9 W (1,2,3) , W (2,3,4) , W (3,4,1) , W (4,1,2) . (6.11) To simplify formulae, bellow we will use the shortcut notations 14) It is also expected that the whole W (A) -algebra is generated by the spin-4 current W (A) 4 (u), so that relations (6.14) for any s follow from s = 4 case. The operators (6.13) will be referred to bellow as reflection operators.
A rigorous proof of existence of the reflection operators is absent. However, assuming that they are exist, it is not difficult to describe the procedure which allows one to construct them explicitly.
Let us denote the Fourier coefficients of the W These W -algebras are naturally associated with four corners of the pillow-brane from Ref. [23]. The notations X i and α i from Ref. [23] coincides with ours φ i and 2i α i , respectively.
values of A = 1, 2, 3, 4 and therefore, using these bases, one can introduce the set of linear operators according to the rule The matrix of the linear operator R (A,B) in the Heisenberg basis is given by In a view of Eq.(6.9), the operatorsW 2n (u) are elements of all W -algebras W (A) . We introduce special notations for these elements: (u) . As for W (4) -basis, one can use the following three linear independent vectors w 1 =W (4) in (6.19). 10 Having at hand an explicit expression for this matrix, other matrixes (R (A) ) β b can be obtained by means of the formal substitutions Then Eq.(6.20) allows one to construct twelve 3 × 3-matrices R (A,B) β ′ β and then check the commutativity condition (6.22).
Returning to general properties of the reflection operators, it can be easily seen from Eqs.(6.20), (6.22) that they are mutually commute and satisfy the relations For our purposes it is useful to enumerate twelve reflections operators in slightly different manner than in the definition (6.14). Namely we define by means of the following relations R (k) Then, Eqs.(6.22), (6.27) and (6.28) imply that eigenvalues of the operators R where S(p i |q) is given by (5.25) and P i related to p i as in (3.35), i.e., P i = 2p i √ a i .

Large-λ asymptotic expansion and dual non-local IM
Here we discuss the formal asymptotic series X (i) σ (λ) (5.28) which appears in the large-λ asymptotic expansion of the connection coefficients A (k) σ ′ σ (iλ) (5.26). 10 We are grateful to A.V. Litvinov for writing a computer code for this calculation.
In the previous section we have described the characteristic property of the local IM -they are integrals over the local densities W 2n−1 are some local fields. In fact, there exists another set of vertex operators satisfying similar conditions [22,47]. Namely, consider six vertex operators (1, 2, 3) , then using explicit formulae for the first two W -currents, it is straightforward to check that for m = 1 and m = 2 We expect, that both Eqs.(6.32) and (6.34) hold for any m = 1, 2 . . . ∞. Let us introduce the following notations for the integrals of the vertex operators (6.33): Repeating the calculations from Ref. [4], one can show that x satisfy the Serre relations for the quantum Kac-Moody algebra Uq i sl(2) : and the conventional notation [n] q = (q n − q −n )/(q − q −1 ) is applied. We may now employ the whole machinery developed in the work [4], to construct families of mutual commuting operators which are also commute with the local IM. In particular, let us introduce the operators Here ρ ± are representations of the so-called q-oscillator algebra generated by the elements H (i) , E (i) 39) and such that the traces exist and do not vanish for complex P i belonging to the lower half plane ℑm(P i ) < 0. The operator X (i) ± (λ i ) can be understood as the series expansion in (λ i ) 2 41) and, as it follows from the result of Ref. [4], each of the expansion coefficient commute with the local IM In fact, the definition (6.38) and the series expansion (6.41) can be applied literally only within the domain In this case all the matrix elements of X (i) σ,n are represented by convergent 2n-fold P-ordered integrals. Furthermore, all the matrix elements are entire functions of (λ i ) 2 in this case. Following the approach developed in Ref. [4], one can show that P-ordered integrals X (i) σ,n can be always rewritten as contour integrals. The such representation allows one to define the operators X (i) σ,n outside the domain (6.43) through the analytical continuation. In particular, the action of X (i) σ,n can be defined in the domain of our current interest, i.e. for 0 < a i < 2. Following the terminology of Ref. [3] we will referred to X (i) σ,n in the domain 0 < a i < 2 as dual non-local integrals of motion. Notice that, the possibility of analytical continuation of the coefficients in the expansion (6.41) does not necessarily imply the convergency of the series. As 0 < a i < 2, Eq.(6.41) should be understood as a formal series expansion with zero radius of convergence. 11 The analytical calculation of the spectrum of dual non-local IM is a complicated unsolved problem. An explicit result can be obtained only for the vacuum eigenvalue of the first IM, X (i) σ,1 . It suggests that, in all likelihood, the formal asymptotic series X (i) σ (λ) (5.28) coincides with the eigenvalue of the formal operator (6.41) provided the following relation between the expansion parameters holds 6.3 Relation to quantum superalgebra U q D(2, 1; α) The appearance of the exponential fields V A (u) and the relation (6.32) suggests a strong connection of our problem with the quantized exceptional affine superalgebra U q D(2, 1; α) . This algebra is generated by twelve elements h 1 , h 2 , h 3 , h 4 , x 1 , x 2 , x 3 , x 4 and y 1 , y 2 , y 3 , y 4 . An unusual feature of this superalgebra is that its Cartan matrix contains an arbitrary (complex) parameter α. So, that together with the "deformation" parameter q the algebra U q D(2, 1; α) has two continuous parameters. For our purposes it is convenient to connect these parameters to our constants a 1 , a 2 and a 3 , used before in (5.5), 46) and introduce additional notations The generating elements of the algebra satisfy the following commutation relations [22,45] [ and also the Serre relations and for any triple (A, B, C) with A, B, C all different 12 Remarkably, as shown in [22], the integrals of the vertex operators satisfy the above Serre relations (6.49), (6.50). Therefore, one could again execute the program of the work [4] now based on the quantum affine algebra U q ( D(2, 1; α)) and define families of commuting transfer matrices, which are entire functions of the variable λ 2 and act directly in the Fock spaces discussed in Section 4. Moreover, in view of the relation (6.32), these transfer matrices will commute with all local integrals of motion I 2n−1 . This direction, however, would requires a lot of additional work, since, to our knowledge, the representation theory of U q ( D(2, 1; α)) is not sufficiently studied. Nevertheless, we expect that the values of the Wilson loop W(λ) (3.14) for various PGHO's will be given by eigenvalues of the transfer matrix obtained as a trace over some finite-dimensional representation of U q ( D(2, 1; α)). Moreover, we expect that the corresponding values of the connection coefficients A (k) σ ′ σ (λ) (5.20) will be given by eigenvalues of appropriate analogs of Baxter Q-operators, obtained as traces over some special oscillator-type representations, first introduced in [3] in the context of U q ( sl(2)).
As a justification of the above picture, consider the value W 0 , given by (3.16). This expression looks like a diagonal character Tr exp( A β A h A ) of an 18-dimensional representation (indeed, it contains 18 exponential terms, where some exponents vanish). Quite excitingly, Zengo Tsuboi has pointed out that his calculations with analytic Bethe Ansatz suggest [46] that the algebra U q ( D(2, 1; α)) does indeed have an 18-dimensional representation. Note, that the corresponding non-affine algebra has only 17-dimensional representation, but there is no an "evaluation map", so that dimensions of representations of the affine and non-affine algebras should not necessarily coincide in this case [46]. We hope to return to this interesting question in the future, as well as to other question relevant to an algebraic construction of commuting transfer matrices in this case.

MShG equation and the auxiliary linear problem 7.1 Complex solutions of MShG
We now turn to further development of the concept of PGHO, where the central rôle is played by the modified sinh-Gordon (MShG) equation where P(z) is still given by (1.3), i.e., P(z) stands for complex conjugate of P(z) and ρ is an arbitrary constant. For ρ = 0, this partial differential equation reduces to the Liouville equation, ∂ z ∂zη − e 2η = 0. In what follows, the field η is understood as a solution of the MShG equation, rather than the Liouville equation.
The subject of our interest are solutions of (7.1), which can be thought as "ρ-deformation" of the complex solutions of the Liouville equation from Section 2.4. To describe their properties it is still convenient to employ the function e −η (see Eq.(2.46)). As before, we assume that e −η is a smooth, single valued complex function without zeroes on the sphere with 3 + L +L punctures. Since z = ∞ is a regular point on the Riemann sphere, e −η satisfy the asymptotic condition e −η ∼ |z| 2 as z → ∞. The asymptotic behavior at the punctures z = z 1 , z 2 , z 3 are given by the same formulae as (2.50), i.e., e −η ∼ |z − z i | −2m i as |z − z i | → 0. Notice that as m i < − 1 4 (2 − a i ) the first term in the r.h.s. of (7.1) dominates as |z − z i | → 0. Therefore, the term ∝ e −2η can be neglected for sufficiently small |z − z i | and we return to the Liouville equation. From the other hand, it is easy to see that in the case of the Liouville equation, the parameters m i should be bounded from below. For this reason we assume that the constraints − 1 2 < m i < − 1 4 (2 − a i ) are enforced.
In the case of the Liouville field the behavior at the punctures z = x a (a = 1, . . . L) andz =ȳ b (b = 1, . . .L) are given by (2.49). This singular behavior consistent with the Liouville dynamics, however the term ∝ e −2η in the MShG equation (7.1) can not be treated as a small perturbation in the vicinity z = x a andz =ȳ b -it essentially modifies the singular behavior at these points. A brief analysis of (7.1) suggests to replace (2.49) by e −η ∼z −xa z−xa and e −η ∼ z−y b z−ȳ b , i.e. the asymptotic formulae (1.11) from the introduction.
We also impose certain "monodromy-free" constraints on positions of the punctures (1.11). For this purpose, let us recall that the MShG equation constitute the flatness condition for sl(2) connection (1.6). Suppose Ψ is a general solution of the associated linear problem (∂ z − A z )Ψ = (∂z − Az)Ψ = 0. The monodromy-free constraints mean that e ± 1 2 ησ 3 Ψ is single-valued in the neighborhood of the points z = x a (a = 1, . . . L) andz =ȳ b (b = 1, . . .L). (Notice that, since e −η is a single valued function, the factor e ± 1 2 ησ 3 does not essentially affect the monodromy properties -its rôle is clarified by the forthcoming consideration.) It is not difficult to reformulate the monodromy-free constraints as local conditions at the punctures imposed on the MShG field η [13]. For this purpose, it is useful to rewrite the matrix differential operators in (1.5) in the form where Let us focus on the differential operator D in the vicinity of the puncture z = x a , where It can be easily seen that as z → x a Here we use the notation γ(z) = ∂ z log P(z) , (7.6) and the gauge transformation is performed by the singular, but single-valued at z = x a , matrix Hence for D is gauge equivalent to a nonsingular at z = x a differential operator. Similarly in the case one can consider the gauge transformation Therefore, asḡ a = 0 , (7.12) D is gauge equivalent to a regular atz =x a differential operator. An immediate consequence of our analysis is that Eqs.(7.8) and (7.12) constitute the single-valuedness condition for e ± 1 2 ησ 3 Ψ at z = x a . Of course, similar consideration can be done for the punctures atz =ȳ b .
The partial differential equation ( where {ω 2n−1 ,ω 2n−1 } ∞ n=1 constitute an infinite hierarchy of one-forms, which are closed by virtue of the MShG equation only, dω 2n = dω 2n = 0 . (7.14) Explicit formulae for {ω 2n−1 ,ω 2n−1 } ∞ n=1 are not particular important for us. They can be found in Ref. [29] (Here we closely follow notations from this paper.) It is useful to mention that ω 2n are usually normalized by the condition (7.15) where dots in the first bracket involves terms with higher derivatives of ∂ z η and/or P(z). Similarlȳ The one-forms are not single-valued on the punctured sphere due to the presence of the multivalued functions P(z),P(z). However, the restriction of ω 2n to the Pochhammer loop depicted in Fig.1 are single-valued and the integrals (7.13) are not sensitive to continuous deformations of the contour. (The second integral in (7.13) is taken over the contourγ P which is complex conjugate γ P .)

Relation to PGHO
Here we describe a relation between the linear problem associated with the complex solutions from the set A (L,L) m and the Perturbed Generalized Hypergeometric Opers.
It is well known, the matrix linear problem (∂ z − A z )Ψ = (∂z − Az)Ψ = 0 can be reduced to second order linear differential operators. The relations (7.2), (7.3) allows one to write its general solution as , (7.17) where ψ andψ solve the equations In the vicinity of the monodromy-free puncture z = x a whereas the monodromy-free conditions (1.13) imply that where u (0) a and u (1) a are defined through the expansion Notice that u(z,z) remains finite at the monodromy-free punctures at z = y b (b = 1, . . .L).
Similarly the fieldū(z,z) is nonsingular at z = x a (a = 1, . . . L), whereas We may now consider the limit ρ → 0. Contrary to the MShG field, the composite fields u(z,z) andū(z,z) admit small-ρ perturbative expansion even in the vicinity of the monodromyfree punctures. It can be easily seen that, with the identification

Local IM versus MShG conserved charges
In the case without the monodromy-free punctures (i.e. L =L = 0), the values of conserved charges q 2n−1 andq 2n−1 coincide: In the recent work [29], a relation between these quantities and vacuum eigenvalues of local IM for the Fateev model was proposed. In this section we discuss a natural generalization of that relation to the excited states spectrum. The definition and some basic properties of the Fateev model [21] have been already described in the introduction. The following clarifying remark on the decomposition (1.23) is in order at this stage. A brief inspection of the Lagrangian (1.19) reveals that the model can be understood within the Conformal Perturbation Theory -the Gaussian theory of three-component Bose field perturbed by the relevant operator. As µ = 0 the general solution of the equation of motion can be written in the form, where we use the "right-moving" chiral Bose field (4.2). The "left-moving" chiral Bose fieldφ(ū) is defined by similar formulae, in particular, As it follows from the periodic boundary condition the zero-modes eigenvalues satisfy the condition P+P = 0, and hence the Hamiltonian of the Gaussian theory acts irreducibly in the tensor product Here F P is the Fock space defined in Section 4 andF −P is similar space generated by the "leftmoving" chiral Bose field. The Hamiltonian of the model (1.19) with µ = 0 acts in the space and sum is taken over the vectors of the form The component H where P = (P 1 , P 2 , P 3 ) related to k = (k 1 , k 2 , k 3 ) as in Eq. (8.6). In this work we do not consider the stationary states corresponding to the higher Brillouin zones, i.e., the subspaces H (n 1 ,n 2 ,n 3 ) k of (8.5) whose CFT limit is described in terms of the tensor product F P+α ⊗ F −P−α with a given vector α = 0 of the form (8.7). Also we do not consider "charged" sectors of the model associated with quasiperiodic boundary conditions The QFT (1.19) possesses an infinite set of commuting local integrals of motion (1.24). In the CFT limit the operators I where P = (P 1 , P 2 , P 3 ) related to (k 1 , k 2 , k 3 ) as in Eq.(8.6). In Ref. [29], the k-vacuum eigenvalues were considered, They are correspond to the vacuum states from H . The vacuum energy is composed of an extensive part proportional to the length of the system, where E 0 stands for the specific bulk energy (1.30) [21]. The main observation of Ref. [29] is that the vacuum eigenvalues (8.14) can be expressed in terms of the classical conserved charges (8.1). The relations are described by Eqs.  The usual approach for studying off-shell physics is the Thermodynamic Bethe Ansatz (TBA). Its key input is factorizable scattering theory underlying integrable QFT. In principle, TBA is a mathematically well defined method for evaluating thermodynamic quantities by solving a set of coupled integral equations. However, in the case of a non-diagonal scattering this method requires many ad hoc assumptions (such as "string hypotheses") and, therefore, is not very practical for complicated theories. The model described by the Lagrangian (1.19) in the regime where all the couplings α i are real, seems to be a good illustration of this statement. Even though the corresponding factorizable scattering theory has been proposed quite a while ago [48], the derivation TBA equations for this complicated theory (to the best of our knowledge) has never appeared in the literature.
The purpose of this section is to demonstrate that the correspondence between classical and quantum integrable systems proposed in the previous section, provides an alternative powerful tool for deriving functional and Bethe Ansatz type equations which determine the full spectrum of local IM in the massive QFT. For the vacuum sector of the Fateev model, we convert our functional equations into the non-linear integral equations and numerically study their solutions, extending the similar analysis from Section 5.5. Note, that a system of integral equations corresponding to the ground state was independently proposed by Fateev [47].

Connection matrices for MShG linear problem
Let us consider to the axillary linear problem (1.7) associated with some element of the finite set A (L,L) m . We introduce three matrix solutions For given i, Ψ (i) solves the linear problem and satisfies the following asymptotic condition Here p i = m i + 1 2 whereas β i stands for arbitrary constant which will be fixed later (see Eq.(9.17) bellow). The connection matrices are defined as follows: At this stage we will treat them as matrix functions of the spectral parameter θ (1.17); that is indicated explicitly in (9.3). They are entire functions of θ satisfying the equations identical to (5.3): The axillary linear problem is invariant with respect to the symmetries analogous to (5.4) 1, 2, 3) .
Similar to derivations of Eqs.(5.8) and (5.10), the symmetry transformations (9.6) allow one to obtain the relation and express the Wilson loop (1.16) in terms of the connection matrices Another easily established symmetry of the axillary linear problem (1.7) involves the operation Using this symmetry it is easy to show that S (j,i) are quasiperiodic matrix functions of the spectral parameter θ: To describe properties of S (j,i) (θ) it is convenient to use the matrices Q (k) (θ) defined through the relation where b i stand for the constants given by Eq. (5.19). In terms of the matrix Q (k) (θ) the quasiperiodicity condition (9.10) looks somewhat simpler: The shift of the argument in the definition (9.11) makes simpler the relation between the matrix elements Q (k) σ ′ σ (θ), (9.14) .

Reconstruction of the connection matrices for L =L = 0
The quasiperiodic entire function Q (k) σ ′ σ (θ) is completely determined by its zeroes in the strip |ℑm(θ)| ≤ π 2 and the leading asymptotic behavior given by (9.18). On the other hand, positions of the zeroes are restricted by the Bethe Ansatz equations similar to (5.45). In fact, the problem of reconstruction of the connection matrices S (ij) (θ) is almost identical to that is studied in Section 5. Thus, we will not repeat the analysis in detail but quote the non-linear integral equation determining the connection matrices in the case without monodromy-free punctures. In this case, using the arguments similar to those given in Section 3.2 from Ref. [20], one can argue that all the roots of Q (k) σ ′ σ (θ) are simple and located at the lines ℑm(θ) = i n + 1 2 π (n ∈ Z). After that the derivation becomes straightforward and yields the system of non-linear integral equations which differs from (5.58) in the source terms only: Once the numerical data for ǫ i (θ) are available, Q (k) σ ′ σ (θ) can be computed by means of the relation sin(νθ) , (9.22) and Φ kl (ν) are defined by (5.57). Notice that, in the case L =L = 0 the conserved charges q 2n−1 andq 2n−1 have the same value which is given by Also the subleading term in the asymptotic (9.18) is given by We solved the integral equation (9.20) numerically for various values of the parameters a 1 , a 2 , a 3 and p 1 , p 2 , p 3 and then calculated the values the conserved charge q 1 from the formula (9.23). The results are in an excellent agreement with the expression for the vacuum energy given by Eqs.(3.25)-(3.28) of Ref. [29]: . Note that the third term in this expression involves the solution of the MShG equation without monodromy-free punctures (its contribution is essential for large values of ρ). The MShG equation has been solved numerically to find the value of the integral, entering (9.26), and calculate the constants η (reg) i , entering (9.19). Using the latter, we have verified the numerical agreement between (9.19) and (9.24).

k-vacuum eigenvalues of local IM in the Fateev model
First of all let us recall some facts concerning the factorizable scattering theory associated with QFT (1.19). All the details can be found in Appendix F in Ref. [48].
The spectrum consists of three quadruplets of fundamental particles with the masses Γ a i 2 (9.28) and their bound states. (Here the relation a i = 4α 2 i is assumed to hold.) The Zamolodchikov-Faddeev commutation relations for the fundamental particles read where (i, j, k) = cyclic perm(1, 2, 3) and Also S a (θ) stands for the conventional S-matrix in the quantum sine-Gordon theory [49] with the renormalized coupling constant a, related to the Coleman coupling β 2 C [50] as follows In particular, the sine-Gordon soliton-soliton scattering amplitude reads explicitly as In a view of our previous discussion it is expected that the non-linear integral equations (9.20) solves the problem of calculation of the k-vacuum eigenvalues (8.14) in the Fateev model in the finite volume. To make the link more explicit let us note the kernels (5.59) in the equations (9.20) are simply related to the amplitude (9.32) and Namely, it is easy to see that These formulae confirm an empirical rule that the kernels of the non-linear integral equations coincide with the logarithmic derivative of some diagonal elements of the S-matrix, which has been previously observed for some other models (e.g., the sine-Gordon model). The µ − ρ relation (1.32) combined with the definition of r i (5.16), implies that the combination 4ρr i , which appears in the source term of the integral equation (9.20), is simply expressed in terms of the particle mass M i (9.28): Also, as it follows from (1.31), the parameters 2 p j /a j should be identified with the quasimomentum magnitudes |k j |. Finally the k-vacuum eigenvalues (8.14) can be calculated using Eqs. (1.27), (1.28) and (9.23). Unfortunately, at the moment there is no independent derivation of the results of Section 9.2 from the field theory side -the Fateev model does not have any known lattice analog, and neither it has any known coordinate or algebraic Bethe Ansatz solutions. 13 10 Concluding remarks In this paper we have described the relation between the MShG equation, on one hand and the Fateev model on the other. We belive that the outlined results open a new general way of approaching integrable QFT.
As an immediate (but perhaps not entirely straightforward) application one could consider various Toda QFT's. This would involve differential operators of higher orders (the g-opers [35]) and the classical modified Toda equations. Some basic ingredients, required for this development have already been revealed. Among them a classification of third order differential operators with monodromy-free singular points, corresponding to stationary states in CFT's with the extended W 3 -symmetry [52], and the relation between the vacuum sector in theÂ (2) 2 Toda QFT and the modified Bullough-Dodd equation [53]. However, perhaps the most important potential outcome of our approach is related to the problem of non-perturbative quantization of classically integrable non-linear sigma models. Here, we are motivated by the following consideration.
This work has been focused on the "symmetric" regime of the Fateev model where all the couplings α i in (1.19) are real, so that the Lagrangian is completely symmetric under simultaneous permutations of the real fields ϕ i and the real couplings α i . The theory is apparently non-unitary in this case. In the most interesting regime one of the couplings, say α 3 , is pure imaginary α 2 1 > 0 , α 2 2 > 0 , α 2 3 := −b 2 < 0 , (10.1) and the theory is governed by the real Lagrangian − 2µ e bϕ 3 cos(α 1 ϕ 1 + α 2 ϕ 2 ) + e −bϕ 3 cos(α 1 ϕ 1 − α 2 ϕ 2 ) , where The physical content in the unitary regime is different from the symmetric one. However, assuming the same periodic boundary conditions for each field ϕ i (i = 1, 2, 3), we can use the same symbols H and H k to denote the spaces of states and their certain linear subspaces in the both cases. Just remember that because of the lack of periodicity in ϕ 3 -direction in the unitary regime, Eq. k and A m in the symmetric regime, let us consider the MShG equation in the regime a 1 , a 2 > 0, a 3 < 0 (the constraint a 1 + a 2 + a 3 = 2 is still assumed). A brief inspection shows that set of requirements (1.9)-(1.14) imposed on the MShG field looks quite meaningful in this case. Only the formulae (1.10) and (1.12) which describe the behavior of the solution in the vicinity of the third puncture z 3 , call for a special attention. As a i > 0 we had a freedom to control the asymptotic behavior of η as z → z i , with the free parameter m i . If a 3 < 0, the situation is different -the leading asymptotic behavior of the solution at z = z 3 is fixed by the MShG equation itself [20]: (10.5) With this modification, we expect Eq.(1.15) to remain a meaningful definition of the moduli space A m , provided m is understood now as a pair (m 1 , m 2 ). Using the intuition gained from the study of the sine and sinh-Gordon models [20], we expect that the relations (1.27)-(1.32) remain valid for the case a 3 = −b 2 /4 < 0. The only exemption is the second formula in Eq.(1.31) for i = 3evidently it cannot be applied literally. Notice also that the definition of the set of the conserved charges {q 2n−1 ,q 2n−1 } ∞ n=1 remains unchanged. We expect that, with these simple modifications the relation between the subspace H (0) k and the moduli space A m holds for the case a 3 = −b 2 /4 < 0. The Fateev model in the unitary regime admits a dual description in terms of the action d 2 x G µν (X) ∂ a X µ ∂ a X ν , where G µν is a certain two-parameter families of metric on the topological three-sphere which possesses two U(1) Killing vector fields [21]. The sigma-model description is especially useful in the strong coupling limit (α 2 i , b 2 → ∞ with α 2 i /b 2 kept fixed), which can be regarded as the classical limit. Notice that the classical integrability of the theory was established only recently in Ref. [54]. while the corresponding accessary parameter will be denoted as C. Eqs. (2.19) allows one to express the accessory parameters c 1 , c 2 , c 3 in terms of C and x, and, thereby, to reduce Eq.(2.17) to a single algebraic equation. The later can be brought to the form For generic values of δ i , (A.1) considered as a cubic equation for C, has three different roots. We will label them by an integer ǫ = 0, ±1. For small x the roots admit Laurent expansions, which can be related to the series expansions for the classical conformal blocks depicted in Fig. 8: To simplify the notations, we shall denote by f (x) the classical conformal block associated to the "principle" branch C (−) (A.5c): f (x) = (1 + 2 p 1 ) log(x) + This defines f (x) unambiguously in the neighborhood of x = 0. In general, the integral depends on a integration contour connecting x to the origin, and the classical conformal blocks corresponding to ǫ = 0 and ǫ = +1 are just different branches of the multivalued function (A.6).
To describe global properties of the multivalued functions C(x) and f (x) we use the fact that any nondegenerate cubic is homeomorphic to an elliptic curve. In the case under consideration the corresponding elliptic modulus (denoted by k bellow) can be chosen as . (A.7) (Recall that k 2 is defined up to modular transformations k 2 → 1 − k 2 , 1/k 2 ). Bellow we use the nome q which is related to the elliptic modulus as k 2 = ϑ 4 2 (0, q) ϑ 4 3 (0, q) .

(A.8)
For the purpose of uniformization of the cubic (A.1), it is useful to introduce three parameters u 1 , u 2 , u 3 such that Here ρ(u, q) stands for the double periodic function ρ(u, q) = ϑ 3 (u, q)ϑ 4 (u, q) ϑ 1 (u, q)ϑ 2 (u, q) = ϑ 4 (2u, q 2 ) ϑ 1 (2u, q 2 ) , (A. 10) and (i, j, k) is an arbitrary permutation of (1, 2, 3). Notice that the above relations define u i up to the overall shift u i → u i + const. For real p i restricted as in Eq.(2.13), the elliptic nome is real and 0 < q < 1, whereas the parameters u i can be chosen in the form In terms of an uniformizing variable u : u ∼ u + Nπ + i M log q (N, M ∈ Z) , (A.12) the cubic (A.1) is described as follows x = − ϑ 3 (u − u 3 − u 2 + u 1 , q)ϑ 1 (u − u 1 , q)ϑ 2 (u − u 1 , q) ϑ 3 (u + u 3 − u 2 − u 1 , q)ϑ 1 (u − u 3 , q)ϑ 2 (u − u 3 , q) ϑ 1 (u 3 − u 2 , q)ϑ 2 (u 3 − u 2 , q) ϑ 1 (u 2 − u 1 , q)ϑ 2 (u 2 − u 1 , q) y = ρ(u − u 3 , q) ρ(u 2 − u 1 , q) , (A. 13) whereas the accessory parameter is given by . (A.14) These equations imply that C is a single-valued doubly-periodic function with simple poles located at u ∈ u i , u i + 1 2 π , u j + u k − u i + 1 2 (π + i log q) (A. 15) corresponding to x = 0, 1 and ∞ at the three sheets of the Riemann surface. Since the classical conformal block (A.6) has the logarithmic branching at these points, it is not a single-valued function on the two-torus. Note that the residues of C at u = u j + u k − u i + 1 2 (π + i log q) do not depend on p i , whereas all the residues of ∂ p i C equals to ±2 (see (A.5)). Using this observation one can show that exp 1 2 ∂f ∂p i = ξ i ϑ 1 (u − u i , q) ϑ 2 (u − u i , q) , (A. 16) i.e., it is a double periodic function as well as the accessary parameter itself. The constant ξ i depends on the normalization prescription for the classical conformal block. For our assignment (A.6), it reads explicitly as ξ i = ϑ 3 (u ji + u ki , q) ϑ 1 (u kj , q) ϑ 2 2 (0, q) ϑ 1 (u ji , q) ϑ 2 (u ji , q) ϑ 1 (u ki , q)ϑ 2 (u ki , q) (u ji = u j − u i ) . (A.17) where Ω(p, p ′ ) = Γ( 1 2 + p − p ′ ) Γ( 1 2 + p + p ′ ) In the case L = 1, one can show that where the notations are inherited from Appendix A. The derivation is based on the general facts (2.39) and (2.11) specialized for GHO with L = 1. Combining these relations with the result (A.16) from Appendix A, one obtains ω i up to an additive coordinate-independent constant ∂ p i F 0 . We may now consider the limit when the monodromy-free puncture approaches to z i . At this limit the differential equations (2.2) becomes the hypergeometric one, and the limiting behavior of ω i can be analyzed explicitly. This fix the value of ∂ p i F 0 for L = 1 and yields formula (B.4).