A note on the S-matrix bootstrap for the 2d O(N) bosonic model

In this work we apply the S-matrix bootstrap maximization program to the 2d bosonic O(N) integrable model which has N species of scalar particles of mass m and no bound states. Since in previous studies theories were defined by maximizing the coupling between particles and their bound states, the main problem appears to be to find what other functional can be used to define this model. Instead, we argue that the defining property of this integrable model is that it resides at a vertex of the convex space determined by the unitarity and crossing constraints. Thus, the integrable model can be found by maximizing any linear functional whose gradient points in the general direction of the vertex, namely within a cone determined by the normals to the faces intersecting at the vertex. This is a standard problem in applied mathematics, related to semi-definite programming and solvable by fast available numerical algorithms. The information provided by the numerical solution is enough to reproduce the known analytical solution without using integrability, namely the Yang-Baxter equation. This situation seems quite generic so we expect that other theories without continuous parameters can also be found by maximizing linear functionals in the convex space of allowed S-matrices.


The 2d O(N) bosonic model, general properties
The 2d O(N ) non-linear sigma model model was solved by Zamolodchikov and Zamolodchikov in [5]. Its action is given by and ⃗ n is an N -dimensional unit vector, ⃗ n 2 = 1. It is an asymptotically free theory that develops a mass gap [6] and therefore it can be considered as a toy model for QCD. The dynamics is described by N equivalent species of scalar particles with an O(N ) symmetry that relates them. There are no bound states, an interesting property for the S-matrix bootstrap maximization program since in previous cases the main idea was to maximize the coupling between particles and their bound states. In the case of only one spatial dimension, in 2 → 2 particle scattering the scattered particles have to have the same momenta as the initial ones. Since we do not initially assume integrability we should allow for possible particle production. Consider the 2 → 2 scattering matrix for (p 1 , a) + (p 2 , b) → (p 3 , c) + (p 4 , d), where p 1,2 (p 3,4 ) are the initial (final) momenta and a, b, c, d are O(N ) indices. The S-matrix can be written in terms of three functions of s = (p 1 + p 2 ) 2 , namely S T (s), S R (s) and S A (s) representing the transmission, reflection and annihilation amplitudes: S ab→cd = [δ ab δ cd S A (s) + δ ac δ bd S T (s) + δ ad δ bc S R (s)] δ(p 1 − p 3 )δ(p 2 − p 4 ) + (p 3 ↔ p 4 )(c ↔ d), (1.2) or equivalently S ab→cd = 1 N δ ab δ cd S I (s) + 1 2 (δ ac δ bd + δ bc δ ad − 2 N δ ab δ cd )S + (s) with S I = N S A + S T + S R , S ± = S T ± S R . (1.4) The functions S I (s) and S ± (s) represent the scattering amplitudes in the three isospin channels: isoscalar, symmetric and antisymmetric. The functions S a , a = I, +, − can be analytically continued in the s plane, defining the physical sheet. They have a cut for real s and s > 4m 2 (the physical line) and another cut for real s and s < 0 corresponding to the t-channel (in 2d, t = 4m 2 − s). The values of S a right above and below the cut are related by the real analyticity constraint: It is convenient to use the relative rapidity variable θ defined through where we also introduced the variable z since it will be useful later (see Figure  1). The physical line is then θ ∈ R, or z = e iφ , 0 ≤ φ ≤ π and the physical region 0 ≤ Imθ ≤ π, or z ≤ 1. The crossing symmetry relates the functions at θ and iπ − θ (or z and −z) as S T (iπ − θ) = S T (θ), S R (iπ − θ) = S A (θ). In terms of the isospin channels the crossing symmetry reads where from now on the indices a, b, . . . label isospin channels and the matrix C satisfies C 2 = 1 and is given by (1.8) It is natural to diagonalize C by introducing the functions Up to now we have only described general constraints on the S-matrix due to standard properties of the field theory.

Exact solution using the Yang-Baxter equation
Now let us describe the exact S-matrix obtained in [5] by using the Yang-Baxter equation. It can be written in terms of the function (1.13) The exact S-matrix reads From here one can compute S I,+,− using eq. (1.4). In the physical region 0 ≤ Im θ ≤ π, S + has a zero at θ = iλ, both S ± have a zero at θ = iπ and none of them has poles. All other zeros and any poles are outside the physical region. Therefore, in this case, the"minimal" solution, namely the one without poles in the physical region, corresponds directly to the O(N ) sigma model. Multiplying by a CDD factor [7] one can obtain the S-matrix of the Gross-Neveu model [5]. In the rest of this paper we analyze the general constraints on the S-matrix and how the minimal solution, namely the O(N ) sigma model can be reproduced by a maximization procedure without any input from integrability.

The convex space of allowed S-matrices
Given a general S-matrix the unitarity constraint is Consider a subsector D of the space of states and a given state ψ α ⟩ ∈ D then we have Namely, the matrix S D defined as S projected to the subspace D satisfies Consequently, the difference S D S † D − I is negative definite. This is equivalent to the positive semi-definite condition Since the space of positive semi-definite matrices is a convex space, so is the space of allowed S-matrices. In the case discussed in this paper, the subsector D will correspond to two-particle states. Using various trial wave-functions, it is easy to see that these constraints read simply Intuitively, for any given total isospin state and fix energy and center of mass momentum (P CM = 0), there is only one two-particle state. Therefore S a (θ) 2 ≤ 1, θ ∈ R since that channel cannot receive contributions from any initial two-particle state other than itself. Thus, if we parameterize the space of allowed S-matrices by the (complex) values S a (θ = σ), σ ∈ R, namely on the physical line, then each variable is restricted to the unit disk which is a convex space. The space is further restricted by the linear crossing constraint (1.8). Therefore, the space can be described as the intersection of the (infinite) product of unit disks and the linear subspace defined by crossing. In the following section we characterize this space more precisely.

Maximization of a linear functional in a convex space
Let us start with a simple example of convex maximization. Consider a convex polyhedron in R n defined by linear constraints that are given by normals where x A ∈ R n and b j are real constants. Any generic linear function F w = ∑ n A=1 w A x A inside such polyhedron attains its maximum at a unique vertex on the boundary (for exceptional values of w A the maximum can be an edge or a face). This property of convex maximization allows for efficient numerical algorithms to obtain the maximum, in particular in this work we used [8].
In this section we argue that the problem at hand is a generalization of this simple problem to infinite dimensions and quadratic constraints. The analogy is made clearer by using a notation where the index A should be thought as (a, σ) or (a, j). Here a = I, +, − labels the different functions in the S-matrix and σ is a continuous index that labels the physical line. Alternatively, as described in the appendix, the index j = 1 . . . M labels a set of points on the physical line used to interpolate the functions and obtain a discretized version of the problem, i.e. M → ∞ recovers the continuum. Also, we do not use the repeated index convention, sums over A are written explicitly and in the case of σ should be understood as integrals.
With this in mind, the unitarity constraints are for every value of the index A. For any analytic function in a region, the real part at the boundary is enough to determine the function inside and, in particular, the imaginary part at the boundary. Thus, if ReS A , is given, namely the real part of the S-matrices on the real line, using crossing symmetry S a (iπ − σ) = ∑ b C ab S b (σ) the real part on the line θ = iπ − σ can be computed. Since now the real part is known in the whole boundary of the physical region one can compute the imaginary part. This results in a linear relation for some kernel K AB . In the θ plane, the explicit form of the kernel K AB is In both formulas we included the term ImS a (θ = iπ 2 ) for completeness. In the case of interest here ImS a (θ = iπ 2 ) = 0 by real analyticity; therefore, we ignore this term in the discretized version used for numerical computations in the appendix. Now, the unitarity constraints read where Since S A 2 ≥ 0, these are positive semi-definite quadratic forms labeled by the index A. Therefore the intersection is a convex space. It is non-empty since the origin ReS A = 0 is in the space, and it is contained in the hypercube ReS A ≤ 1. The (curved) faces of this generalized polyhedron are given by the points that satisfy all constraints but saturate only one (namely for one value of A). The edges are determined by the points saturating all constraints except for one value of A.
The integrable model saturates all the constraints and therefore is at the boundary of this space. Moreover, since the integrable model has no free parameters, we do not expect to have any other near-by point that saturates all constraints and therefore the integrable model should be at a vertex (or corner) of this space. In fact this expectation can be checked numerically. By the discussion at the beginning of the section, to find such a vertex, it is sufficient to maximize a linear functional such that w A points in the general direction of the vertex. Near the vertex that we take to be at ReS A = ReŜ A the space has the shape of a cone that we can obtain by linearizing the quadratic forms using a perturbation ξ A around the maximum: It follows that the space of S-matrices is given, at the linearized level (see fig.2), by the constraints The matrix V C B can be thought as a set of vectors labeled by C and normal to the faces of the cone whose tip is at the vertex. The edges are given A B Figure 2: Schematic representation of the convex space of allowed S-matrices. The non-linear sigma model is at a vertex like A where there are no directions such that at the linearized level satisfy the constraints in both directions. At a point like B there are tangents that are parallel to the boundary and appear as zero modes of the linearized constraints. The dotted lines represent the linearized constraints described by by allowing only one constraint not to be saturated. Arranging the vectors parallel to the edges in a matrix we get the inverse matrix of V . Numerically, once a reasonable direction is found it can be improved iteratively by taking The rationale is that such direction is the sum of the normals to the faces and therefore should point in the direction of the vertex. The numerical calculations agree well with this naive expectation. There are many choices for the initial direction w A but a little experimentation shows that adequate functionals are Here, S 1,2 are the functions defined in eq.(1.11). There are various choices for the point θ 0 . One we choose frequently is θ 0 ≃ i (corresponding to z 0 = 0.3i) but it does not have to be on the imaginary axis. There is a range of N α 5 13 6 7. 5  7  6  8  5  20  3  100 2.4 Table 1: Values for α coefficient in maximization function coefficients α that can be chosen depending on θ 0 and N . We give suggested values in Table 1. As mentioned, the functional, namely the choice of vector w A , can be improved iteratively using eq.(3.10). After a few iterations this leads to a rapid decrease of the error as measured by how close the solution is to saturate the unitarity bounds. After the iterative procedure converges, the maximum of the functional is simply where M is the number of interpolating points on the physical line. The reason is that the integrable modelŜ A saturates all the constraints and therefore Finally, we can compare the maximization results with the exact integrable model as shown in Figures 3-6. It is clear that there is a very good agreement. We emphasize again that this agreement was obtained by maximizing a functional without making any assumptions on the S-matrices such as position of its zeros. Only the constraints of unitarity, crossing and real analyticity were used. Finally, for the integrable model to be at a vertex, the matrix V B A should have no zero modes. The reason is that we cannot move away from a vertex and still saturate all constraints (see Figure 2). If there is a zero mode ξ then a deformation along ξ (0) A will still saturate all constraints. Numerically we checked that the matrixṼ = V T V has no zero eigenvalues in agreement with V having no zero modes.

Analytic results from the numerics
In this section we show how one could extend the numerical solution into an analytical one by using some simple inputs from the numerics. This might seem as a moot point since the analytical solution is already known but it is useful to see how the solution arises without using factorization. The numerical solution clearly shows that there is a zero of S + in the imaginary axis of θ and that S ± are zero at (or near) θ = iπ. See Figure 7. Moreover, since the functions are related by crossing, it is clear that they should have most zeros and poles in common. We then expect that the ratios have a small number of zeros and poles. Actually since the S-matrices saturate the constraints ( S a (θ) = 1, θ ∈ R), then R 1 (θ) = R 2 (θ) = 1 on the real axis, the zeros and poles should be situated symmetrically with respect    to the real axis: As discussed, the numerical solution is such that, on the imaginary axis, S I has no zeros, S − has one zero at iπ and S + has two zeros, one at iπ and the other at an intermediate point θ = iτ 0 that can be determined from the plot (e.g. Figure 7 for N = 8). The naive ratios that take these zeros into account fit the numerical results perfectly well as seen in Figure 8. From the crossing conditions: it results that the only solution where R 1,2 are ratios of linear functions is the one described with one zero at iπ, and the other at exactly τ 0 = λ = 2π N −2 , in good agreement with the numerics. It is now convenient to use as a reference the crossing invariant function S 2 in eq.(1.11). Some simple algebra gives 20) Since S 2 (iπ − θ) = S 2 (θ) and the poles and zeros of S I,+,− are symmetric with respect to the real axis, it is clear that, for example, in S + , the zero at iλ implies a pole at −iλ that should be in S 2 . But then S 2 should have a pole at iπ + iλ that appears in S + implying a zero at −iπ − iλ and so on. Thus, the factor (θ − iλ) should be extended as and similarly with the factor (θ − iπ). If we were to do the same in S I the factor (θ + iλ) should be extended without generating a pole at iλ and that is accomplished by having the pole at −iλ, namely the same S 2 as needed in S + . Finally, by normalization, the overall factor N −2 π 2 should be absorbed in S 2 . All in all we get for S + from where we can compute S I and S − satisfying all the required properties and in agreement with the integrable model. It seems that one can find other functions that saturate all constraints by assuming that R 1,2 are ratios of higher order polynomials since that appears to allow more freedom given that we have more roots to choose. However it is easy to see that the number of equations grows faster suggesting that it is not easy to solve the crossing constraints in terms of other than linear polynomials.

Saturating the unitarity constraints
As we saw, the integrable model saturates unitarity, however they are not the only functions to do that. In this subsection we discuss other functions that also saturate those constraints as a way to study the boundary of the space of allowed S-matrices.

The CDD factors
The 2d O(N ) bosonic integrable model admits a CDD type ambiguity [7,5]. The minimal solution can be multiplied by any number of CDD factors: and the result preserves crossing symmetry and unitarity saturation. Notice that a CDD factor will introduce a pair of poles into the physical strip at θ = iα and θ = iπ − iα if α is taken to be real and Re α ∈ [0, π]. In particular one can use such a factor to obtain the Gross-Neveu model S-matrices [5].
Here we want to stay in the space of functions analytic in the physical region and avoid adding any poles. One can then construct S-matrices of the typẽ These are at the boundary of the convex space because they satisfy all the constraints and saturate the unitarity condition. We use the explicit form of such CDD dressed integrable model to construct the w A according to eq.(3.10) and indeed we find them numerically through maximizing the functional eq.(3.7). We also check the rigidity of the CDD dressed integrable model in this space, namely, whether it is a vertex. We find that the matrix has 6 zero modes for each CDD factor which corresponding to 6 independent infinitesimal deformation that preserve the crossing and unitarity condition. We include modes that do not satisfy real analyticity since they can be useful when considering fluctuations of several CDD factors. To understand these 6 zero modes, it is helpful to consider the following simpler example of S-matrix which corresponds to a free theory S-matrix dressed with one CDD factor. Using the change of variables θ → z in eq.(1.6) the CDD factor can be written as where z 0 is the corresponding zero inside the unit disk and f (z) satisfy crossing and saturate unitarity: Now consider a small deformation of the S-matrix S a (θ) → S a (θ)±ξ a (θ) which preserves the crossing symmetry and unitarity saturation for both signs (see fig.2). We then have ξ a (iπ − θ) = C ab ξ b (θ) (4.9) and to first order S * a (σ)ξ a (σ) + S a (σ)ξ * a (σ) = 0. (4.10) A described in Appendix C, finding the analytic functions ξ a (z) is the so called Riemann-Hilbert problem that, in this case can be solved by elementary methods. Indeed, since S * a (σ) = 1 Sa(σ) , we can therefore define the meromorphic functions with the property Re[η a (σ)] = 0. (4.12) In the case of this simple example (4.6), the η's satisfy the same crossing condition as the ξ's (4.9), and have poles where the f has zeros. Therefore, one has In terms of the unit disk in z, this corresponds to meromorphic functions with zero real part on the boundary of the unit disk. One can construct such functions explicitly: The last term in the expression of h 2 (z) is added to ensure that h 2 (z = 0) = 0. We include h 2 (z) for completeness but it does not satisfy the real analyticity constraint so it is a deformation outside the space of allowed S-matrices.
With h 1 and h 2 , we can also write crossing symmetric and antisymmetric combinations: They give rise to the following zero modes for {η I , η + , η − }: plus another three zero modes obtained replacing H 1 → H 2 giving the 6 zero modes of (4.6). It should be noted that onlyη 1,2,3 satisfy real analyticity. Consider now the case whereŜ a indicate the integrable model S-matrix. Numerically one still finds 6 zero modes. Analytically we were able to identify 4 of them: Again, only two of them are real analytic.

Other Vertices
By choosing other linear functionals to maximize, it is possible to find other vertices in the boundary of the space. We give an example in Figure 9. This is purely a numerical result, that is, we obtain functions that saturate the constraints and such that the corresponding linearized constraints had no zero modes. It will be interesting if the existence of these vertices can be confirmed by finding their analytic expressions. Since the Yang-Baxter equation leads to the integrable model [5] these vertices should not correspond to integrable theories, at least in the S-matrix factorization sense.

Conclusions
It was recently argued that the S-matrix of certain theories can be computed by maximizing a coupling between particles and their bound states under the assumption of a fixed spectrum [3,2,4]. It was not clear how the program proceeded if the particles formed no bound states. For that reason, in this paper we considered the 2d O(N ) bosonic integrable model that is exactly solvable and has no bound states. What we found was that there is an infinite number of functionals that give rise to this model by maximization. We emphasize that the maximization procedure reproduces the S-matrix of the integrable model using only the assumptions of crossing, real analyticity and the unitarity bounds. We do not use any extra information about zeros, poles or other properties of the S-matrix. The reason is that this theory has no parameters other than N which appears in the crossing matrix and therefore the S-matrix should be fixed without any extra information.
Given that there is an infinite number of functionals that one can maximize, instead of concentrating on the functional, we proposed that the defining property of the model is that it is at a vertex of the convex space of allowed S-matrices. Then, it is clear that any linear functional whose gradient points in the general direction of the vertex is maximized precisely at the vertex. In fact, maximizing a linear functional in a convex space is a standard numerical problem related to semi-definite programming and solvable by fast algorithms. There has to be some initial trial and error in choosing an appropriate functional but once found, it can be iteratively improved by linearizing the convex space near the vertex. We expect this can lead to a general method to define field theories without continuous parameters 1 . Such theories should lie at particular points of the boundary of the convex space of allowed S-matrices and should be found by maximizing linear functionals in such space. Of ciurse it will be of great interest if this idea applies to QCD, see the recent work [4,9] for interesting ideas in higher dimensions.

Acknowledgments
We are very grateful to Lucía Gómez Córdova and Pedro Vieira for discussions on their related research. We are also grateful to Oleg Lunin and Anatoly Dymarsky for comments and suggestions. This work was supported in part by DOE through grant de-sc0007884.
Note: While this work was being completed we learned of [10] that also studies the 2d O(N ) non-linear sigma model using the S-matrix bootstrap. After this work was submitted [11] appeared with related results.

A Numerical approach
The physical region is defined as 0 ≤ Imθ ≤ π in the θ plane and as z ≤ 1 in the z plane. The relation between the two is The physical line (E > 2m) is the real axis in θ coordinates and z = 1, Imz ≥ 0 in the z plane. By assumption, there are no bound states so the scattering matrices are analytic in the physical region (no poles). Given a general, analytic function in the unit disk, the real and imaginary parts are related by This equation becomes eq.(3.4) after the change of variables σ = − atanh(cos ϕ) which maps the upper half circle in z to the real line of θ. The first step is to find a numerical analog to this relation, that is, given the values of Reφ(z j = e iϕ j ), ϕ j = jπ M for j = 1 . . . 2M we would like to find the approximate values of the imaginary part at those same points, the approximation being better as M → ∞. Before going into the derivation let us mention that the result is a straight-forward discretization of the integration kernel in the previous integral: The derivation also gives us the analytic function inside the disk and proceeds as follows. Given the values of Reφ(z j ) we approximate the real part at the boundary of the disk as andŜ(z) a holomorphic function in the unit disk such thatŜ(z) =Ŝ(−z) by crossing and ImŜ(iξ) = 0 for −1 ≤ ξ ≤ 1 by real analyticity. Consider now the CDD factor f 0 (z) = z 2 + a 2 1 + a 2 z 2 (B.4) such that f 0 (e iφ ) = 1, f 0 (z) = f 0 (−z) and f 0 is analytic in the unit disk with zeros at z ± = ±ia. Define now It is clear that S 1 (z) is analytic in the unit disk, and satisfies S 1 (e iφ ) ≤ 1, crossing and real analyticity. Moreover which is real as required by real analyticity. In [2,3,4] it was proposed to maximizeg or, as we see here, equivalently Re(S 1 (ia)). Since the maximum of the modulus of S 1 (z) should be at the boundary of the disk, and is bounded by S 1 (e iφ ) ≤ 1, the maximum is max[Re(S 1 (ia))] = 1 attained by the constant function S 1 (z) = 1. This gives which agrees with the S (1,1) (θ) scattering matrix in the sine-Gordon model [5] as already explained in [3]. Here we want to point out that we could equally well maximize F max = Re(S 1 (ib)) (B.8) for any −1 < b < 1 and the maximum will be the same S 1 (z) = 1, i.e. there is an infinite number of functionals that we can maximize. Since we can equivalently write the functional F max in terms of the original Smatrix S(z): F max = Re(S(ib)), for b 2 < a 2 , or F max = −Re(S(ib)), for b 2 > a 2 . (B.10) For the purpose of this paper it is useful to check if this S-matrix is at a vertex of the allowed space. Consider a fluctuation S 1 (z) = 1 ± ξ(z) , (B.11) with ξ(z) analytic in the unit disk and satisfying crossing and real analyticity: ξ(z) = ξ(−z), (B.12) Im(ξ(iη)) = 0, −1 ≤ η ≤ 1. (B.13) Since we require S 1 (e iφ ) ≤ 1, at first order for both signs ±ξ (see fig.2) we should have Re(ξ(e iφ )) = 0 . (B.14) This implies that ξ(z) = 0 and therefore there are no zero modes. The integrable model is once again at a vertex of the space of allowed S-matrices.

C Relation to the Riemann-Hilbert problem
In section 4.1 we discussed that the problem of finding fluctuations to the maximum is equivalent to the RH problem. As discussed in eq.(4.10), what we want to find is a fluctuation ξ a (z) analytic and crossing symmetric such that S a (t)ξ * a (t) + S * a (t)ξ a (t) = 0, t = e iφ , 0 < φ < π (C.1) Due to crossing symmetry we also have S a (−t) b C ab ξ * b (t) + S * a (−t) b C ab ξ b (t) = 0, t = e iφ , π < φ < 2π (C.2) To obtain the standard formulation of the RH problem [12], define functions analytic inside the disk ξ If we also define the matrix C ab (t) = Ŝ * a (t) δ ab 0 < φ < π S * a (−t) C ab π < φ < 2π (C.6) then the problem requires finding two analytic functions, ξ (±) a , one outside the disk ξ If the RH problem has no solutions then the theory is at a vertex of the space of S-matrices, if not, there is a continuum of theories satisfying the constraints in the vicinity of the maximum. Unfortunately, the vector RH problem has no analytical solutions, numerical solutions are equivalent to checking the zero modes of the matrix V as previously discussed around eq.(4.5). In the case of the simple S-matrix (4.6) the RH problem is diagonal and can be solved by standard methods [12]. The result is After imposing crossing we obtain ξ a = (1 + z 4 )(α 1 ψ 1a + α 2 ψ 2a ) + z 2 (α 4 ψ 1a + α 5 ψ 2a ) + (α 3 z + α * 3 z 3 )ψ 3a (z 2 − 1 which agrees with the simpler calculation in eq.(4.20). Here the constant vectors ψ 1,2,3 are the eigenvectors of the crossing matrix C: with eigenvalues 1, 1, −1 respectively.