Bootstrapping Massive Quantum Field Theories

We propose a new non-perturbative method for studying UV complete unitary quantum field theories (QFTs) with a mass gap in general number of spacetime dimensions. The method relies on unitarity formulated as positive semi-definiteness of the matrix of inner products between asymptotic states (in and out) and states created by the action of local operators on the vacuum. The corresponding matrix elements involve scattering amplitudes, form factors and spectral densities of local operators. We test this method in two-dimensional QFTs by setting up a linear optimization problem that gives a lower bound on the central charge of the UV CFT associated to a QFT with a given mass spectrum of stable particles (and couplings between them). Some of our numerical bounds are saturated by known form factors in integrable theories like the sine-Gordon, E8 and O(N) models.


Introduction
One can define a quantum field theory (QFT) non-perturbatively as a renormalization group (RG) flow from the UV to the IR fixed point. The fixed points are assumed to have conformal invariance and are described by conformal field theories (CFTs). In this work, we will focus on massive QFTs, i.e. the IR CFT is trivial. To specify a particular QFT it is sufficient to provide the UV CFT and the relevant deformation triggering the RG flow. We would like to determine IR observables like the mass spectrum and scattering amplitudes from the UV data. However, generically, the RG flow is strongly coupled and it is not possible to compute IR observables using perturbation theory around the UV CFT. In these cases, one has to resort to numerical methods (like lattice field theory, Hamiltonian truncation, tensor networks, etc) that require a UV cutoff and a costly extrapolation to the continuum limit. This calls for modern non-perturbative bootstrap methods that can constrain the space of QFTs directly in the continuum. Unfortunately, the present bootstrap methods study the UV and the IR separately. Namely, one can use the conformal bootstrap [1] 1 to study the UV CFT data or the S-matrix bootstrap [3][4][5][6][7][8] to study scattering amplitudes of light particles. Ideally, we would like to connect these two bootstrap approaches. This work is a step in this direction. 2 Our strategy is simple. We consider a set of states that include asymptotic scattering states and states created by local operators acting on the vacuum. In the simplest setting, we consider the following three states 3 |ψ 1 = |p 1 , p 2 in , |ψ 2 = |p 1 , p 2 out , |ψ 3 = dxe i(p 1 +p 2 )·x O(x)|0 . (1.1) Unitary implies positive semi-definiteness of the matrix 3 where S denotes the 2 to 2 scattering amplitude, F 2 denotes the two-particle form factor of the operator O and ρ its spectral density. Notice that the S-matrix bootstrap can be formulated in the same way if we only use scattering states (|ψ 1 and |ψ 2 in this case). The presence of the spectral density ρ ∼ OO in the setup makes it straightforward to establish a connection with the UV CFT: at large energies all the correlation functions should coincide with the ones of the UV CFT. For example in the case of the two-point function of the stress-tensor, conformal invariance in the UV fixes its form uniquely up to a constant known as the C T central charge [9]. Notice that even if we knew the S-matrix for all energies it would not be easy to extract from it information about the UV CFT. 4 In this work we test this strategy in 1+1 dimensional QFTs. In this case, we can write the central charge c of the UV CFT as an integral over the spectral density of the trace of the stress tensor. This allows to address the following question: what is the minimal central charge of a UV CFT that can give rise to a massive QFT with a given set of masses and couplings 5 of stable particles?
In practice, we use analyticity of the amplitude and form factor, to write a general ansatz for the amplitude, the form factor and the spectral density. Then we numerically optimize the parameters of the ansatz such that unitarity is obeyed at all energies and the value of the central charge is as low as possible. In several cases we find that the optimal form factors are given by known integrable theories such as the sine-Gordon, the E 8 and the O(N ) models.
The structure of the paper is as follows. In section 2 we review all the basic ingredients in a consistent manner for generic number of dimensions and provide all the normalization conventions. In section 3 we formulate unitarity as the semipositive definite condition on the three by three matrix ψ a |ψ b and discuss its implications. Then, we illustrate how this works with analytic examples from 2d integrable models like the sine-Gordon, the E 8 and O(N ) models, which we review in section 4. In section 5 we define and set up the numerical linear optimization problem for 2d QFTs. We also present our results and compare them with the analytic formulas for integrable models. We conclude and briefly discuss applications to higher dimensions in section 6. We derive various auxiliary results in appendices A, B and C.

Review of basic ingredients
We work in (1, d − 1) Minkowski space with the mostly plus metric η µν = {−, +, . . . , +}. (2.1) The position and momentum in this d-dimensional space are denoted respectively by where x and p are the position and momentum in the (d − 1) Euclidean subspace. We refer to p 0 as energy and x 0 as time.
We study unitary quantum field theories with restricted Poincaré symmetry group. 6 When working in d = 2 in addition we will also assume parity. In what follows we summarize our conventions and review basic ingredients such as asymptotic states, scattering 4 In d ≥ 3 using holography one can argue that the regime of hard scattering (high energy and fixed angle) should be directly related to the UV CFT [10]. However, the present S-matrix bootstrap methods are not precise enough in this regime. 5 We define couplings from the physical S-matrix. For example, cubic couplings are given by residues of poles of the 2 to 2 scattering amplitudes. 6 It is defined as the Poincaré group without parity and time-reversal discrete subgroups. and partial amplitudes, spectral density and form factors. We will conclude with a discussion of unitarity and its implications.

States
The state of a system described by the unitary QFT is represented by a "state" vector in an infinite dimensional Hilbert space. In this space it is convenient to choose a basis of state vectors (or simply states) in such a way that they are eigenstates of the generators of translations P µ with eigenvalues p µ and transform in the irreducible representation of the Little group SO(d − 1) which leaves invariant the d-vector {p 0 , 0}. We will always work with states which have a strictly positive energy p 0 > 0. We also restrict our attention to traceless symmetric representations of the SO(d − 1) Little group. Then any state will have at least three labels |p, j, µ . (2. 3) The label j is a non-negative integer called spin (j = 0, 1, 2, . . .). The label µ denotes the components of the spin j irreducible representation of SO(d − 1). In the case of d = 4, one can choose µ = −j, . . . , +j to be the helicity, i.e. the projection of spin j on the direction of p. The normalization of the states (2.3) is chosen as p , j , µ |p, j, µ = (2π) d δ (d) (p − p)δ j j δ µ µ (2.4) except for the special case of one particle states that are normalized as in (2.8). We note that from this normalization condition it follows that the dimensionality of the state vector is |p, j, µ = − d 2 . (2.5) For further discussion, see construction of irreducible unitary representations of the restricted Poincaré group [11,12].

Free particle states
As we review in the next section, the asymptotic states of an interacting massive QFT are in one-to-one correspondence with the states of a non-interacting QFT. For this reason, we first consider a free QFT. There is a special set of states, called the one particle states, which describe a single freely propagating particle. The tensor product of one particle states defines multi particle states which describe a system of multiple non-interacting particles. The Hilbert space spanned by all the possible one and multi particle states is called the Fock space.
The states (2.3) which obey the "mass-shell" condition where m is a discrete real non-negative number, called mass, are referred to as one particles states (1PS). We can denote them as |m, p . (2.7) We focus only on scalar particles in this work, thus we omit the spin labels j = µ = 0. The one particle states are normalized as m , p |m, p = 2p 0 δ m m × (2π) d−1 δ (d−1) ( p − p). (2.8) From the above normalizations it is clear that the one particle states have the following mass dimensions |m, p = − d − 2 2 . (2.9) We define the n particle state as |n ≡ |m 1 , p 1 ⊗ . . . ⊗ |m n , p n . (2.10) The n particle state has a well defined total d momentum which reads as p µ = p µ 1 + . . . + p µ n . (2.11) Due to the very definition of the Fock space one can write the completeness relation in this space by summing over n particle states and integrating over their phase space as where the phase space Φ n for n identical particles is defined in (A. 19). Let us finish this section by focusing on two particles states |2 of identical particles with mass m. Writing all the labels explicitly we denote it by |m, p 1 ; m, p 2 ≡ |m, p 1 ⊗ |m, p 2 . (2.13) This state does not transform in the irreducible representation of the restricted Poincaré group (like any other n particle state with n ≥ 2) simply because it is not in the irreducible representation of the SO(d − 1) Little group. We can project it however to irreducible representations. For simplicity we focus on the two particle states in the center of mass frame defined as p 2 = − p 1 and the vector p 1 has an angle θ 1 with the x 1 axis and θ 2 = . . . = θ d−2 = 0. See (A.4) for our conventions for spherical coordinates. The projection is done by integrating over the (d − 1) scalar spherical harmonics, which are the Gegenbauer polynomials, as 7 where γ j is some coefficient fixed by the normalization, which we derive in (2.18). In the left-hand side (2.14) we dropped the label µ because we are considering states with zero spin projection along p 1 and invariant under SO(d − 3) rotations that leave the scattering plane (θ 2 = . . . = θ d−2 = 0) invariant. The normalization of the two particle state (2.13) is fixed by the normalization of one particle states (2.8). One has m, p 1 ; m, p 2 |m, p 1 ; m, p 2 = 4p 0 Notice that the normalization (2.15) reflects explicitly that the system is symmetric under the permutation of particles 1 and 2. In the second line of (2.15) we have performed a change of variables, see (A.17). For identical particles the factor N d , derived in (A.18), reads as In (2.15) the spherical angles Ω and Ω correspond to the (d−1) vectors p 1 and p 1 respectively. The δ-function in spherical coordinates is defined in (A.8). 8 The Mandelstam variable (2.17) defines the square of the total energy for the two particle state in the center of mass frame.
We can now evaluate the value of the constant γ j . Using (2.4), (2.16) and the orthogonality relation (A.11), we get 9 where the coefficient ν is defined in (A.11). Finally, let us invert the projection (2.14) by means of the orthogonality relation (A.12). One finds where the Clebsch-Gordan coefficient C j (cos θ 1 ) reads as

Asymptotic states
This section is based on chapter 3.1 of [11]. We work with states in the Heisenberg picture (states do not evolve in time) and describe the entire evolution of the system. They are defined however with an implicit choice of a 8 Given a spherical angle Ω of a d − 1 vector p, we schematically denote by −Ω the spherical angle of a d − 1 vector − p. If the former has the angles (θ1, . . . , θ d−3 , θ d−2 ), the latter has all the angles shifted as (π − θ1, . . . , π − θ d−3 , π + θ d−2 ). This is easy to see from (A.4). 9 Without loss of generality the factor γj is chosen to be purely real in the rest of the paper.
reference frame f . Suppose we have another reference frame f with time t = t + τ . 10 If a state |ψ is seen by an observer in f , the same state will be seen by an observer in f as |ψ . Due to time translation invariance these two states are related as In strongly interacting theories with a mass gap one can define a (complicated) state in the reference frame f , which for an observer in the reference frame f either in the far past (τ → −∞) or in the far future (τ → +∞) however will look like a set of non-interacting (free) particles. We call states with such a property the asymptotic in and out states. In what follows we will make this statement formal.
We assume that the strongly interacting Hamiltonian of our system can be written in the following way where H 0 is a free Hamiltonian (with the mass spectrum including stable composite particles and bound states) and H int is the "interaction part". Note, that H int is not a small perturbation around H 0 and we do not know how to construct it explicitly. 11 As a consequence the expression (2.22) is highly formal. The eigenstates of the free Hamiltonian H 0 are nothing but the n particle states defined in (2.10), in other words where p µ is the total d-momentum of the n-particle state (2.11).
We can now define the in state |n in and the out state |n out via the following conditions 12 Here the asymptotic and the free n particle states are defined in the reference frame f and are required to match in the reference frame f . Given the condition (2.24) one can express the in and out states in terms of the free n particle states as where we have defined the operator known as the Møller operator. For details see [13]. Clearly, the Møller operator is unitary 10 If some event happens at t = 0 in f , the very same event happens at t = τ in f . 11 Not all the systems can be written in such a way. Notable example are system with long-range interactions. 12 The relation below should be understood in a sense of wave packets.
From this it follows that the normalization of the asymptotic states is the same as the one of the n particle states in m|n in = out m|n out = m|n . (2.28) Let us from now assume that all the asymptotic states in the theory span a complete basis of states. Then the completeness relation (2.12) can also be written for the asymptotic states. Multiplying (2.12) by Ω(∓∞)Ω † (∓∞) and using (2.27) and (2.25) we simply get I = n |n in in n| = n |n out out n|. (2.29)

Scattering and partial amplitudes
The scattering process of n free particles in the far past and m free particles in the far future is described by the n → m scattering amplitude defined as follows.
Here p i and p i describe the d-momenta of the one particle states constituing |n and |m , p and p denote the total incoming and outgoing momenta. In (2.30) we have explicitly extracted the overall δ-function. The scattering operator S due to (2.25) reads as It can be split into the trivial and the interacting part as From now on let us focus on the 2 → 2 processes of identical scalar particles. The expression (2.30) then reads as where we have defined the three Mandelstam variables as which obey the standard constraint Notice, that the s variable has already appeared in (2.17). The partial amplitude of the 2 → 2 process is defined as a matrix element of the S operator between the states (2.14), namely (2.36) We would now like to write the relation between the partial amplitude S j (s) defined in (2.36) and the scattering amplitude (2.33). In principle this can be done by simply plugging (2.14) into (2.36), however it is easier to derive this relation in the following way. Take the two particle states in the center of mass frame m, p ; m, − p | and |m, p; m, − p . We align the incoming particle p with the x 1 axis. The outgoing particle p will have an angle θ 1 with the x 1 axis. All the other angles for incoming and outgoing particles are chosen to be zero. In this frame it is very convenient, instead of using (s, t, u) obeying the constraint (2.35), to use the variables (s, cos θ 1 ). The variables t and u can then be written as We consider the scattering amplitude (2.33) and apply the decomposition of states (2.19).
Using the definition (2.36) we can write This relation can be inverted by means of (A.11) and leads to where the coefficient κ j reads as Notice, that for identical particles only the partial amplitudes with even spin j exist. The partial amplitudes (2.39) with odd spin j vanish due to the x ↔ −x symmetry of the scattering amplitude and the antisymmetry of the Gegenbauer polynomial.
To conclude let us address the consequences of (2.32). The scattering amplitude (2.33) can be split into the connected and the disconnected parts according to (2.32). Using (2.15) we can write where we have defined The connected (interacting) part of the amplitude T (s, t, u) should not be confused with time-reversal operator T which we unfortunately denote in the same way. Combining (2.39) and (2.42) we can also write a similar expression for the partial amplitude This matches precisely the expression given in equation (10) of [5].

Spectral density
In this section we will discuss two important instances of two-point correlation functions, namely the two-point Wightman and time-ordered correlation functions. We will define the notion of spectral density and show how both types of two-point functions can be rewritten in terms of the spectral density.
Wightman correlation functions Let us consider a local operator O(x) and study the Wightman two-point correlation function where i are infinitesimal positive numbers. See (A.20) for slightly more details. In what follows we will not display the i 's in order not to complicate the notation. They are however always present and must be taken into account when we deal with Wightman functions. By using translation operators we can write 13 Assuming that the basis of asymptotic states is complete, we can inject the completeness relation (2.29) into (2.46) and using (2.47) we find where θ is the step function and ρ is the spectral density defined via 14 The Fourier transform of the Wightman function (2.46) is related to the spectral density in the following simple way (2.51) When taking the Fourier transform we integrate over coincident points which is potentially dangerous. The presence of i 's in the time components ensures that the integral always converges since it gives a dumping prefactor e −H with > 0. The spectral representation (2.48) can be further rewritten by injecting an additional δ-function and integrating over it. One then has where we have defined the Wightman propagator as which satisfies the Klein-Gordon equation Time-ordered correlation functions Now let us consider the time-ordered correlation function Plugging here the expression (2.52) we obtain the Källén-Lehmann spectral representation of the time-ordered two-point correlation function where the Feynman propagator is defined as Equivalently, we can write the Feynman propagator (2.57) in its standard form 15 From the above expression it is clear that the Feynman propagator satisfies the Klein-Gordon equation with a source Finally, the Fourier transform of the time-ordered two-point function reads as High energy behavior In the UV, due to the presence of conformal symmetry, the Wigthman function (2.46) is completely fixed and reads as 16 where ρ CFT (s) simply follows from (2.61). The precise value of the constant factor is irrelevant for this work. It is found straightforwardly by performing the Fourier transform carefully. For its value see (2.4) in [14] and section 2 in [15].

Form factors
Consider the following matrix elements called the form factors where p is the total d-momentum of the in and out asymptotic states. The Fourier transform of the matrix elements (2.66) reads as Let us now discuss the structure of F n (p 1 , . . . , p n ) for n = 0, 1, 2. When n = 0 or n = 1 we simply get F 0 , F 1 = const.
(2.69) since it is impossible in these cases to construct a scalar function out of zero or one d-momenta.
When n = 2 we can form only one scalar object out of two d-momenta p 1 and p 2 which is simply the s Mandelstam variable (2.34). We can write then Analogous statements hold for G n .
CPT invariance As a consequence of the CPT theorem [16] there is always an anti-unitary operator Σ in the theory which acts on scalar local operators and asymptotic states (of scalar neutral particles) as We can use this fact to write the following equality From the definition of the spectral density (2.49) we see that the one particle state according to (A.19) gives the following contribution Analogously the two particle states contribute to the spectral density as where we have performed the change of variables according to (A.17). Combining the above we get the following expression for spectral density ρ(s) = ρ 1 (s) + ρ 2 (s) + . . . , (2.75) where the one and two particle contributions read as and the dots denote the contribution of n ≥ 3 particle form factors.
Crossing symmetry Let us now consider the following matrix element out p 1 |O(x)|p 2 in , where both p 1 and p 2 satisfy the "mass-shell" condition (2.6). We demand that this matrix element satisfies crossing 17 which can be written as Notice, that the right-hand side of (2.77) is not the usual form factor defined in (2.65), because it has a negative energy −p 0 2 . Thus, the expression in the right-hand side of (2.77) is related to the usual two particle form factor by an analytic continuation.
Constraint from the UV Due to the relation of the form factors with the spectral density (2.76) and the UV behavior of the spectral density (2.63) one obtains the following bound on the large s behaviour of the two particle form factor (2.78)

Unitarity constraints
We are now ready to discuss the implications of unitarity for scattering amplitudes, partial amplitudes and form factors. In order to do this we will exploit the unitarity of the S operator which follows from (2.31) and (2.27). Taking into account (2.32) it reads as (2.79)

Appearance of poles
The main goal of this section is to argue that the interacting part of the scattering amplitude and the two particle form factor contain simple poles and show how they are related. 18 Scattering amplitude We focus here on the interacting part of the two to two scattering amplitude. Using (2.79) and the completeness relation (2.12) we can write m, p 1 ; m, p 2 |T |m, p 1 ; m, p 2 − m, p 1 ; m, p 2 |T † |m, p 1 ; m, p 2 = i n m, p 1 ; m, p 2 |T |n n|T † |m, p 1 ; m, p 2 . (2.80) Let us focus on the left-hand side of (2.80) and evaluate it in the center of mass configuration where p 2 = − p 1 and p 2 = − p 1 . Furthermore we use the following spherical angles (0, 0, . . . , 0) 17 The crossing equations for the form factors in 2d are discussed for example in [17] and [18]. In general dimensions they can be derived in the QFT framework using the LSZ procedure. For the derivation of crossing equations in the case of scalar form factors in 4d see chapter 7.2 in [19]. For the derivation of crossing equations in the case of scattering amplitudes in 4d see section 5.3.2 in [20]. 18 Strictly speaking the presence of poles cannot be deduced from the pure S-matrix approach and should be accepted as an additional assumption. One can however trade this assumption for another one, namely the existence of the relation (2.80) for complex values of external momenta, which in turn allows for all the derivations in this section. In order to discuss rigorously the presence of poles one needs to appeal to a higher level framework, e.g. quantum field theory, see section 10.2 in [11].
for the vector p 1 and (θ 1 , 0, . . . , 0) for the vector p 1 . Labeling the states for transparency by the square of the total energy s and the angle θ 1 we have In the last equality we have used rotational invariance. Using the fact that the matrix element depends on the angle via cos θ 1 we conclude that the left-hand side of (2.80) reads as Let us now discuss the right-hand side of (2.80) in a generic frame. We focus on the special case of n = 1 where (2.82) Due to translation invariance we can extract an overall delta function of the matrix elements entering (2.82) as where s = −p 2 = −(p 1 +p 2 ) 2 is the total energy. Notice that (2.83) does not exist for physical (d − 1)-momenta and is defined via an analytic continuation. Plugging (2.83) into (2.82) we perform the integral and get the final form of (2.82) which reads as We can now evaluate (2.84) in the center of mass frame and plug it into (2.80) together with (2.81). Dropping the overall delta function we get the following expression where . . . denote the continuous part due to n ≥ 2 particle states. This corresponds to a pole in the s complex plane, Due to presence of n ≥ 2 particle states the imaginary part of the amplitude (2.85) is non-zero for s ≥ 4m 2 . This implies that the amplitude itself develops a discontinuity or equivalently has a branch cut 19 along the real axis for s ≥ 4m 2 .
Form factors Given the definitions of the form factors, we can use the completeness relation (2.29) to write the following equality Using the definition of the scattering amplitude (2.30) and its splitting into the trivial and interacting part (2.32) we can re write the above relation as In the first term of the right-hand side (2.88) we have used the normalization of multi particle states which removes the sum over m and integration over the phase space. (This normalization follows from (2.8). See the first line of (2.15) for the example of two identical particles). Let us focus on the n = 2 case, using (2.73) this allows to write (2.88) as Plugging here (2.83), analogously to (2.84), we get where F * 1 is a constant as discussed below (2.69). We thus get an analogous expression to (2.85). Assuming analyticity in some region of the complex plane we obtain the pole structure of the form factor Analogously to the discussion below (2.86) the form factor develops a branch cut along the real axis for s ≥ 4m 2 due to the contribution of two particles states (and higher) in (2.90).

Watson's equation
Let us consider the two particle contribution to the completeness relation (2.12). Using the change of variables to spherical coordinates (A.17) and (2.19) we can write In the last equality we have used the explicit expression of the Clebsch-Gordan coefficient (2.20) and the orthogonality of the Gegenbauer polynomial (A.11).
We can now consider matrix elements of the first unitarity condition (2.79) with two particle states of definite spin (2.14). Injecting the completeness relation (2.12) and focusing on the two particle contribution (2.92) we get the standard unitarity condition on the partial amplitudes where the dots denote the n ≥ 3 particle contribution which starts at s ≥ (3m) 2 .
Let us now consider the two particle form factor and use the decomposition of the two particle state into irreducible representations (2.19) to write the following equality In the last equality we have used the rotation invariance and the fact that the local operator O(0) is a scalar which selects j = 0 representations only. In (2.94) C 0 is simply a real constant which follows from (2.20) since the Gegenbauer polynomial is one for j = 0. Consider now the relation (2.87) where instead of the two particle asymptotic state we use the two particle state projected into j = 0 irreducible representation. We also rewrite the two particle contribution in the completeness relation according to (2.92) in order to get where S 0 (s) is the zero spin partial amplitude and the dots denote the n ≥ 3 contribution which starts at s ≥ (3m) 2 . This can be simply rewritten by using (2.94) to get the final expression of interest Using the equality (2.73) we obtain the Watson's equation [21] S 0 (s) = In integrable models, there is no particle production and the Watson's equation is valid to all energies. One can solve the Watson's equation (2.97) and obtain the form factor in terms of the partial amplitude up to an analytic function which is real on the real axis s, see for example [22]. 20,21 20 See also [22] for a nice application of pion form factors in phenomenology. 21 In the S-matrix literature it is common to use the term real analytic function which means an analytic function which takes real values on some interval of the real axis. In mathematics instead the term real analytic function means a function on (an interval of) the real axis which allows analytic extension to its neighborhood. The function itself may be real or complex on the real axis, see for example [23].

Stress-energy tensor
We discuss here the stress-tensor in general dimensions. We show that there is an integral of the two-point function of the stress-tensor which gives the central charge C T of the UV CFT. We will work in Euclidean signature within this section. The final results however are independent of the signature. Consider the stress-energy tensor T µν (x) operator, which satisfies the following constraints The most general form of the two-point function of the stress-tensor [24] with appropriate mass dimensions which respects Lorentz invariance is 22 where h i (r) are scalar dimensionless functions and T i are dimensionless linearly independent tensor structures which read as Notice, that away from fixed points one is required to have dimensionful parameters to construct dimensionless functions h i (x). At fixed points there are no dimensionful parameters and thus all the functions h i (x) are simply constants. Let us also define the following three contracted two-point functions Comparing these expressions with (2.99) we can write . 22 In d = 2 and d = 3 there are additional parity odd tensor structures which we omit here.
Let us apply the conservation equation (2.98) to the two-point (2.99) using (2.104). Setting to zero coefficients of three independent tensor structures we get the following three conditions Here we have used the fact that There are five functions h i (r) with three differential constraints on them. There are thus only two independent functions which define the two-point function of the stress-tensor. We can take various linear combinations of three equations (2.103) to form a single differential equation. Following [24] we can write for example Another expression, more convenient for d ≥ 3 is as follows Integral expressions for the central charges In the presence of conformal symmetry the form of the two-point function (2.99) is severely restricted. According to [9] it reads as 23 Here C T is one of the central charges of the UV CFT. Comparing this form with (2.99) we deduce that in CFT Provided that our QFT is defined as a flow between the UV and IR fixed points (which are reached at r = 0 and r = ∞) we can write the differential conditions (2.106) and (2.108) in an integral form using (2.111). We get two equivalent expressions Notice, that the latter holds only for d ≥ 3. In a massive QFT, the theory in the IR is empty and thus we have Stress-tensor form factor Let us consider the two particle form factor of the stress-tensor. It has the following most generic form which is symmetric in both indices and satisfies the conservation condition (2.98) written as Taking trace of (2.115) we obtain the form factor of the trace of the stress-tensor 2 (s). Normalization of the stress-tensor We can form the following conserved charges which are the generators of translations and Lorentz transformations respectively. In particular the Hamiltonian is H = P 0 as in (2.22). Let us now evaluate the matrix elements of P µ with one particle states. Since they are the eigenstates of P µ we get the following expression where p 0 1 satisfies the "mass-shell" condition (2.6). On the other hand we have We can now compare the matrix elements in (2.121) and (2.115) taking into account the crossing relation (2.77) which effectively makes a replacement p µ 2 → −p µ 2 . This leads to and F 2 (s) ∼ s around s = 0. For more details see appendix G in [25]. The Lorentz generators in (2.118) do not provide any further conditions. Special case of 2d Let us focus now on the specific case of 2d [26]. It is conventional to use complex coordinates defined as (2.123) In these coordinates we can write the components of the stress-tensor as Notice the presence of 2π factors in the definitions (2.124). Conservation implies At the critical point we have where c =c is the standard central charge in parity preserving 2d CFTs. Using (2.124) we can compare this form with (2.109). We conclude that In this convention the central charge of a free boson is c = 1, see (C.12). We can rewrite the integral expression (2.112) using (2.101), (2.124) and (2.127) in the following form Due to reflection positivity of the two-point function of the stress-tensor, we can conclude that c U V > c IR . This is Zamolodchikov's c-theorem [27,28]. No such statement can be made about C T in higher dimensions.
The Euclidean two-point function in (2.128) is time-ordered. We can then use the Euclidean Källén-Lehmann spectral representation (B.5) to relate the central charge c U V with the spectral density ρ of the trace of the stress-tensor. We have Using the integration by parts we arrive at the final expression where in the second equality we have used (2.75) and (2.76).

Unitarity as positive semidefiniteness
We are now ready to present the main idea of this paper. We will construct a hermitian matrix which must be semipositive definite in a unitary theory. This requirement intertwines the partial amplitudes, the form factors and the spectral density and puts constraints on them.

General spacetime dimension
We will work with the simplest case of identical particles with mass m. Let us define the following three states where p µ is the total d momentum The first two states |ψ 1 and |ψ 2 are the in and out two particle states projected to the irreducible spin representation according to (2.14) with the total d-momentum p µ . The third state is the Fourier transform of the state generated by the local operator O(x) acting on the vacuum. It also has p µ total d-momentum. The extra factor m −δ is injected in order to make all three states to be of the same mass dimension This follows from (2.5) and (2.9). Notice, that the parameter δ has already appeared in (2.63). Let us now construct a 3 by 3 matrix out of all possible inner products of the states where a, b = 1, 2, 3 and the total d-momentum of the states |ψ a and ψ b | are p µ and p µ respectively.
Entries of the B-matrix Let us now inspect the entries of the matrix (3.6). The entries 11 and 22 on the diagonal are simply fixed by the normalization condition (2.4) since the Møller operators are unitary, see (2.28), and thus read as Using (2.51) we can write the entry 33 as Let us address now the off-diagonal elements. Since the matrix (3.6) is hermitian we will only need to discuss the elements 12, 13 and 23. The element 12 reads as where we have used the definition of the S operator (2.31) and partial amplitude (2.36). The element 13 reads as where in the second line we have used (2.14), (2.47) and the results of section 2.4. The coefficient ω is defined as Simplifying we get the following compact result (3.14) Analogously, the element 23 can be written as Positivity constraint Plugging all these expression into (3.6) and using (2.73) we recover the final form of the B-matrix The matrix B is hermitian by construction and must be positive semidefinite in unitary theories. This can be easily seen as follows. The matrix B is positive semi-definite if and only if its eigenvalues are non-negative. One can show that the latter is the case by taking a linear combination of states (3.1) -(3.3) for which the B-matrix is diagonal. The elements of this matrix are simply the norms of the news states. Unitarity of the theory requires these norms to be non-negative. Thus, B j (s) 0, ∀s ≥ 4m 2 and ∀j. Consequences of the positivity constraint Let us start with the minor associated to removing the third row and column. The Sylvester's criterion leads to This is the standard unitarity constraint for the partial amplitude already obtained in (2.93). Now consider instead the minor associated to removing the first row and column. The Sylvester's criterion leads then to This inequality also follows straightforwardly from (2.75) for s ≥ 4m 2 . The minor associated to removing the first two rows and columns leads to the following requirement ρ(s) ≥ 0, (3.21) which was already obvious from the definition (2.49). Finally, the determinant of the B matrix must be non-negative, This provides a non-trivial positivity condition which mixes together the amplitudes, the form factors and the spectral density.
Degenerate situation Let us now investigate a very particular situation when only one state out of the three (3.1) -(3.3) is linearly independent. This for instance happens in the energy range 4m 2 ≤ s ≤ 9m 2 . (3.23) We refer to this situation as the absence of "particle production" in the range of energies (3.23). One can imagine even a more extreme case when there is no "particle production" for the whole range of energies s ∈ [4m 2 , +∞). In d ≥ 3 according to the Aks theorem [29] this situation leads to a trivial theory. Theories in d = 2 escape this constraint however and we enter the realm of integrable models.
In what follows we investigate the consequence of having only a single linearly independent state among (3.1) -(3.3) or equivalently the situation when rankB j (s) = 1. (3.24) The characteristic polynomial in λ is then required to have the following form where λ 0 is the only non-zero eigenvalue of the matrix B. Let us now compute the characteristic polynomial for the B matrix (3.17), it gets the required form (3.25) with if the following conditions are fulfilled The latter equation is solved by

Special case of 2d
Let us summarize here the unitarity constraints for the special case of 2d. We will then generalize them to include the O(N ) global symmetry. Let us start with partial amplitudes. Since the Little group is the discrete Z 2 group, effectively we have a single partial wave with spin j = 0. From now on we denote it aŝ S(s) ≡ S 0 (s). (3.30) Moreover, in 2d there is not much difference between the partial and the scattering amplitude. In our conventions they simply differ by a normalization where T (s) as before is the interacting part of the scattering amplitude and the normalization factor is given by (A.18) and reads in 2d as From now on we will measure every dimensional quantity in units of mass m. This is equivalent to setting m = 1. (3.33) The unitarity constraint (3.17) and (3.18) read in 2d as For the future purposes it is also convenient to rewrite this expression in the following way The full scattering amplitude can be decomposed into three independent scattering amplitudes σ i (s), i = 1, 2, 3. In the notation of [30] we have cd m, p 3 ; m, p 4 |S|m, p 1 ; m, p 2 ab =(2π) 2 δ (2) (p 1 + p 2 − p 3 − p 4 )× σ 1 (s)δ ab δ cd + σ 2 (s)δ ac δ bd + σ 3 (s)δ ad δ bc . (3.37) The st-crossing symmetry (under exchanging particles 1 and 3) relates the amplitudes σ i as (3.38) The two-particle state is in the reducible O(N ) representation and can be further decomposed into three irreducible representations as where we have defined where the tensor structures associated to the three irreducible representations are defined as The relation between two sets of amplitudes σ 1 , σ 2 , σ 3 and S • , S S , S A simply reads as The normalization of two particle states in the irreducible representation of the O(N ) group follows from (3.36). We have

47)
A [cd] m, p 3 ; m, p 4 |m, p 1 ; m, p 2 The crossing equations (3.38) in terms of the amplitudes (3.45) read as  Let us now consider the unitarity constraints. We have three states transforming in irreducible representations of O(N ). They cannot mix with each other, in other words nonzero inner products can be formed only between the states in the same representation. Let us start with the trivial representation where the hatted amplitudes are defined according to (3.31) and the form factor is defined as For the symmetric and antisymmetric representations we consider only the in and out states because they do not overlap with the state created by the O(N ) invariant local operator. 24 The unitarity conditions then simply read as

Analytic examples in 2d
In this section we provide a uniform summary of the exact analytic expressions of the partial amplitudes and form factors in several 2d integrable models, namely the sine-Gordon, the E 8 model (also known as the 2d Ising model with magnetic deformation) and the O(N ) σ-model with N ≥ 3.
θ variable In 2d instead of the Mandelstam variable s it is convenient to use the rapidity variable θ. Given a particle with the 2-momentum p µ i and the mass m i we can define  The Yang-Baxter equations together with unitarity and crossing (4.4) allow to obtain exact analytic expressions for partial amplitudes up to a CDD ambiguity [31]. The latter states that given the solution to all the above constraints, one can obtain another solution by multiplying it with any number of CDD factors (and their inverses) 25 defined as Here α ∈ (0, 1) is a real parameter. 26 Notice the factors of π in (4.5) compared to the standard definition. The CDD factor (4.5) satisfies automatically both constraints (4.4). It contains a pole at θ = iπα and thus encodes the contribution of a given asymptotic state to the amplitude. The correct choice of the CDD factors is usually postulated and then gets checked in perturbation theory for some range of parameters in the model when it is applicable.
Form factors In 2d the form factors satisfy the following equations The former is the familiar Watson's equation and the latter encodes crossing symmetry. Given the analytic expression for the partial amplitude, the equations (4.6) can be solved analytically [32]. The solution reads as where R(θ) is an arbitrary rational function of cosh(θ) since cosh(θ) automatically satisfies the second condition in (4.6). More precisely it can be written as where we have defined 2 cos 2 (πα/2) cos(πα) − cosh θ . (4.9) In (4.8) the parameters α i define the positions of poles and parameters β i (or equivalently A, B, C, etc) define the positions of zeros. The constant factors in the numerator of (4.9) are introduced for convenience, they allow for the following normalization K α (iπ) = 1. The "minimal" form factor F min (θ) in (4.7), also known as the Omnès solution in higher dimensions [33], is defined as a function without poles or zeros. According to [32] in 2d due to (4.6) it can be expressed in terms of the partial amplitude as follows (4.10) 25 An inverse of the CDD factor (4.5) introduces zeros in the amplitude. 26 One can consider the CDD factors with negative or even complex values of the parameter α. For a discussion see page 12 of [4].
The functions (4.9) can be thought of as analogues of the CDD factors (4.5) for the form factors. The choice (4.8) is usually postulated first and then gets checked with perturbation theory when applicable. Finally the form factor (4.7) for a given operator must obey the bound (2.78) which reads in 2d as 27 (4.11)

sine-Gordon model
The quantum sine-Gordon model is defined as the renormalization group (RG) flow triggered by the deformation of the free scalar UV CFT by the vertex operator in the following way Here φ(x) is the real scalar field, ∆ V β is the UV scaling dimension of the vertex operator 28 , m 0 is a mass-like parameter and β is a real coupling constant. The sine-Gordon model possesses several remarkable properties. First, the model is dual to the Thirring model [35]. 29 Second, it possesses the O(2) topological symmetry [30] and thus can also be regarded as the O(2) σ-model. The mass spectrum of the sine-Gordon model was first found with semi-classical methods [37,38] and later argued to be exact [39,40]. It consists of a soliton and an antisoliton with mass m and a number of breathers (soliton -antisoliton pairs) denoted by b n with the masses , n = 1, 2, . . . , 8π γ . (4.14) Here x denotes the greatest integer not larger than x. The breathers exist only in the following range of parameters It is interesting to study the sine-Gordon model in the regime when at least two breathers exist. Then one can define a ratio of masses for the first two lightest breathers (4.16) 27 In 2d this bound was first derived in [34], see formulas (3.33) and (3.34). 28 The scaling dimension of the vertex operator can be straightforwardly deduced from the Euclidean two point function computed in the free massless theory (which posses conformal invariance for m = 0) In the second step we have used the propagator 0|φ(xE)φ(0)|0 = − 1 4π log x 2 E . 29 See also chapter 6 of [36].

The range of parameters which allow for this is
Let us now discuss partial amplitudes for the scattering of asymptotic states in the sine-Gordon model. The soliton -(anti)soliton scattering was computed in [30]. The (anti)soliton -breather and breather -breather scattering was computed in [41]. The uniformed treatment for all these cases was done in [42]. In this work we are concerned only with the lightest breather -breather scattering in the parameter range (4.17) which reads aŝ It is given by the single CDD factor (4.5) and thus contains a single pole at θ = iγ/8 or equivalently at √ s = 2m 1 cos(γ/16). From (4.16) we see that this pole is simply at the mass of the second breather √ s = m 2 .
Let us now address the form factors of a scalar operator. The scalar soliton -(anti)soliton form factors were computed in [43]. The scalar breather -breather form factors were found in [17]. 30 The latter form factor corresponds to the partial amplitude (4.18) and reads as where A is a normalization constant, K is given by (4.9) and F SG min is the minimal form factor of the sine-Gordon model where we have defined The form factor (4.19) for the vertex operator (4.12) is the most general solution (with a single pole due to b 2 ) which satisfies the bound (4.11) for the whole range of parameters (4.17) since lim and ∆ V ∈ [0, 2/3]. The form factor of the trace of the stress-tensor is proportional to the UV deforming operator (vertex operator in our case). It is thus also given by (4.19). The value of the constant A follows from the normalization convention (2.122) and reads as A = −2m 2 1 .
(4.23) 30 For some more recent work see [44]. 31 In order to show this, one can make a variable redefinition x → x ≡ xθ. Keeping x fixed, we can then consider only the leading behavior of the integrand at large θ and perform the integration.  Here we wrote explicitly only the pole and denoted by . . . the finite part at s = m 2 2 . Similarly for the form factor we have In case of the trace of the stress-tensor we can use (2.76), (2.130) and the explicit expressions for the form factor (4.19) to estimate the contribution to the total central charge of the one particle state of the second breather c b 2 and of the two particle state of the first breather c b 1 b 1 . The total central charge then reads as For concreteness, on figure 1 we provide the numerical value of c b 2 as a function of m 2 2 . All the contributions in (4.27) should sum up to c = 1 which is the central charge of a free boson, see (C.12). For more detailed investigation of the sine-Gorden model see [45], in particular figures 6 and 7.

E 8 model
The 2d Ising model is a 2d conformal field theory with a Z 2 symmetry. It contains only two relevant operators σ and with the scaling dimensions ∆ σ = 1/8 and ∆ = 1 respectively. The former is Z 2 odd and the latter is Z 2 even. In the lattice formulation of the 2d Ising model the operator σ couples to the magnetic field. We consider here a QFT obtained by deforming the 2d Ising model with the operator σ. This QFT defines an integrable model [46] which we refer to as the E 8 model. It contains eight asymptotic states m 1 , . . . , m 8 . Given the value of m 1 , the spectrum in the E 8 model reads as m 2 = 2m 1 cos(π/5), m 3 = 2m 1 cos(π/30), (4.28) where we have ignored particles with masses m 4 , . . . , m 8 since their masses lie above the two particle threshold 2m 1 and are thus invisible to the techniques of section 5. The partial amplitude for the scattering of the lightest asymptotic state reads aŝ where t are the CDD factors (4.5). The form factor for a scalar relevant operator with the lightest asymptotic states was computed in [34], it reads as where A and B are independent parameters and the minimal form factor for the E 8 model reads as It has thus the most general form which obeys the bound (4.11) for both σ and operators. According to [34] the expression (4.30) provides the two particle form factor for the σ and operators, given the following ratios of parameters The overall normalization of the form factor depends as usually on the chosen normalization of the operators σ and and is not important in our work. We set for convenience Consider now the form factor of the trace of the the stress-tensor Θ. It is proportional to the form factor of the deforming operator which is σ in our case. Thus, the coefficients A Θ and B Θ for the trace of the stress-tensor Θ are completely fixed by the following conditions where the second equation follows from the normalization condition (2.122). Let us discuss now the interacting part of the scattering amplitude in s variable which is related to the partial amplitude via (4.24). Given (4.29), we can write its pole structure as The values of the trilinear couplings read as Similarly for the form factor of the trace of the stress-tensor we can write where the one particle form factors read as In case of the trace of the stress-tensor, using (2.76), (2.130) and the explicit expressions for the form factor (4.30), (4.35) we can estimate the central charge contribution of the one particle states of the first three lightest asymptotic states c i and the two particle contribution of the very first asymptotic state c 11 . The total central charge reads

Non-linear sigma model
The O(N ) non-linear sigma model (NLSM) with N ≥ 3 is defined in the UV via the Lagrangian density where n i (x) is a O(N ) vector of real scalar fields and g 0 is a dimensionless coupling. This model can be seen as a marginally relevant deformation of a theory of N −1 free massless scalar fields. 32 The NLSM is asymptotically free in the UV and is gapped in the IR. Away from the 32 One way to see this is to solve the constraint on the scalar fields in (4.42), and write The operator O(x) is marginal since its UV scaling dimension is ∆O = 2.
UV fixed point, its spectrum consists of a single asymptotic state of mass m transforming in the vector representation of the O(N ) group. The scattering of asymptotic states is described according to (3.37) by three amplitudes σ 1 , σ 2 and σ 3 or equivalently byŜ • ,Ŝ S andŜ A according to (3.43). The relation between two sets of amplitudes is given in (3.45). The analytic expressions for σ 1 , σ 2 and σ 3 in the NLSM were found in [30], they read as where σ 2 (θ) is given by the "plus" part of (3.17) and (3.18) in [30]. The results of [30] can be rewritten in a compact integral form [17] aŝ together withŜ To characterize the strength of the interaction in the NLSM one can evaluate the partial amplitudes at the crossing symmetric point s = 2 (which corresponds to θ = iπ/2). We have the following values of the partial amplitudes then together withŜ (4.48) We notice that crossing equations have the simplest form for σ i (θ) partial amplitudes. At the crossing symmetric point they lead to the equality σ 1 (iπ/2) = σ 3 (iπ/2). The form factor associated to the antisymmetric partial amplitude was computed in [17]. It reads as The form factor associated to the scalar operator was reported in [47] and reads as 33 where A is the normalization constant. One can estimate the asymptotic behavior of expressions (4.49) and (4.50) at large energies. 34 One has From these asymptotics one sees for instance that the form factor (4.50) is the most general expression which satisfies the bound (4.11) for a relevant scalar operator since −( λ 2π + 1) < 0 for N ≥ 3.
The form factor for the trace of the stress-tensor takes the same form (one can say that it is proportional to the operator O in (4.43)). It is thus given by the expression (4.50) with the following value of normalization constant which follows from the normalization condition (2.122) and the definition of • states given in (3.40). Notice that there are no poles in the scattering amplitudes (3.45) or in the form form factor (4.50). Also the one particle form factor is zero, F Θ 1 = 0. These follow from the O(N ) symmetry (some matrix elements simply cannot be constructed). One can use the form factor (4.50) to compute the two particle contribution to the central charge c 2 via (2.76) and (2. where the dots represent four and higher particle contributions (notice that the odd number of particles in the majority of cases does not contribute due to O(N ) symmetry). In order to obtain (4.53) we have evaluated (4.50) numerically for multiple values of θ, we have interpolated the results to obtain a continuous function and integrated it numerically to get the value of the central charge for different values of N . 35 The contribution coming from two particles and all the multi particle states should some up to the central charge of N − 1 free bosons which is c = N − 1.

Numerical bootstrap in 2d
We are now in position to formulate the numerical bootstrap problem which allows to obtain the partial amplitude, the form factor and the spectral density of a UV complete massive 34 At large θ the integrand in (4.49) has a highly oscillating piece. In order to study the asymptotics of such integrals one needs to rewrite them as a (generalized) Fourier integral. In our case we have Riemann-Lebesgue lemma states then that such integral vanishes at large θ. Its leading behavior can be estimated by using integration by parts, where the leading behavior comes from the boundary term. One has then d dθ ln F A (θ) ∼ − λ 2π 1 θ . 35 It is important to change the integration variable to θ in (2.130) in order to perform the numerical integration. The reason for that is the very slow convergence of the integral for large values of s. unitary QFT. We will focus on two dimensions in this section for two reasons: to avoid technical complications due to spin in higher dimensions and to be able to compare our results with analytic results for 2d integrable models discussed in section 4. Given a QFT which has at least one asymptotic state with a non-zero mass, we can consider the following three functionŝ S(s), F Θ 2 (s), ρ Θ (s), (5.1) which are the partial amplitude for 2 to 2 scattering of the lightest particle, the form factor and the spectral density of the trace of the stress-tensor respectively. The form factor of the trace of the stress-tensor is normalized according to (2.122). The spectral density of the trace of the stress-tensor is related to the UV central charge according to (2.130). We can write the most general ansatz for the functions (5.1) with real coefficients entering linearly. In case of theŜ(s) function, these coefficients are further restricted to satisfy crossing. We can then solve the following problem: determine the parameters of the ansatze leading to the minimal possible UV central charge such that the functions (5.1) obey the unitarity condition (3.34).
In section 5.1 we will provide the details of the numerical setup. In section 5.2 we will present the numerical results. We will consider three different cases: partial amplitude with a single pole, partial amplitude with three poles and a partial amplitude with no poles but with a global O(N ) symmetry. We will see that these cases will reproduce numerically the known results in the sine-Gordon, E 8 and O(N ) integrable models respectively.

Setting up the optimization problem
It is convenient to introduce the r variable defined as It maps an s complex plane with one branch cut s ∈ [4, ∞) into the unit disc (with the cut mapped to the boundary). The point s 0 < 4 is a free parameter which is mapped to the center of the disc. 36 Another useful variable is φ(s 0 ) defined as .
In what follows we will often use the φ variable defined as which maps a ray into an interval, more precisely Let us discuss now the situation when our QFT has k asymptotic states below the twoparticle threshold According to the discussion of section 2.5.1 these asymptotic states will appear as simple poles in the interacting part of the amplitude and the form factor. We can then write the following ansatze [5], a n × r(s; 2) n + (s ↔ 4 − s), where we have defined for convenience The physical domain is defined via s + i with > 0. We can thus rotate the cuts using the identity These ansatze depend on the set of real parameters a n , b n , c n and d n which enter linearly. The form factor of the trace of the stress-tensor obeys the normalization (2.122). This leads to the linear constraint for the unknown coefficient Taking it into account we can write the final ansatz for the form factor as There is in principle no bound on how big the spectral density ρ Θ (s) can be. However there is certainly at least a trivial lower bound ρ Θ (s) ≥ 0 due to (3.21) which implies c U V ≥ 0. The expression (5.16) has the unknown constants F 1,i entering in a quadratic way. Thus, we cannot directly apply methods of linear programming to minimize (5.16). We can however use a simple trick to rewrite (5.16) in a linear way. Consider the following inequality where we have introduced new non-negative real parameters u i which obey the following constraints The latter inequality is equivalent to ∀i : Now instead of minimizing c U V we can minimize c bound U V defined in (5.17) given that the condition (5.19) is satisfied. At the minimum we will simply get u i = |F 1,i | 2 and c U V = c bound Central charge minimization We are finally ready to formulate the numerical bootstrap problem: given a set of asymptotic states and their masses (5.6), determine the linear coefficients g 2 i , λ i , u i , a n , b n , c n , d n , in the ansazte (5.7) -(5.9) and (5.17) such that the semipositive conditions (3.35), (5.19) are satisfied and the central charge c bound U V in (5.17) has the minimal possible value. Sometimes we will also be fixing the values of g 2 i in order to single out known integrable models. To perform the numerics we use the semipositive program solver SDPB [48,49].

Numerical results
We now solve the optimization problem of the central charge minimization defined in section 5.1 in three different cases. First, in the presence of a single pole. We find a special point on the central charge bound which corresponds to the sine-Gordon model. We will recover its partial amplitude and the two particle form factor. We will then investigate the dependence of the central charge on the parameter of the sine-Gordon model. Second, in the presence of three poles. Injecting the values of the masses and the residue of the lightest asymptotic state in the E 8 we recover numerically the partial amplitude of the E 8 model. We also obtain the form factor consistent with the analytic results. Third, we address the case of no poles in the presence of O(N ) global symmetry. We will minimize the central charge by scanning over different values of the partial amplitudes at the crossing symmetric point.

One pole
We assume that the system is described at least by two asymptotic states with masses This parameter range is chosen to mimic the sine-Gordon behavior, see (4.16) and (4.17). We consider the scattering of the m 1 asymptotic state and assume that there is only one simple pole in the scattering amplitude due to the second asymptotic state. 37 We can now look for a minimum of the UV central charge fixing the value of the trilinear coupling g = g 2 . We also set the mass of the second asymptotic state to be m 2 = √ 3. The numerical results are presented on figure 3. On the plot there appears a special value of the trilinear coupling for which the optimization problem becomes unfeasible. This critical value is g ≈ 4.55901. (5.21) The value (5.21) is in a perfect agreement with the results of [4] where it was found that there is an upper bound on the trilinear coupling, see figure 4. It was also found that this value corresponds to the b 1 b 1 b 2 trilinear coupling of the sine-Gordon model (4.25), where b 1 and b 2 stand for the first and the second breathers. At the critical value (5.21) we also recover the partial amplitude of the first breather b 1 and the form factor of the trace of the stress-tensor. They are presented in figures 4 and 5 respectively. They match precisely the exact analytic expressions (4.18) and (4.19). The numerical expression for the spectral density is given on figure 6. It is completely saturated by the two particle contribution and thus it should be regarded as the two particle part of the spectral density. This is a general feature of our numerical results.
The numerical procedure allows to determine the UV central charges of the sine-Gordon . Therefore, the two particle contribution of the first breather accounts for the whole central charge.

Three poles
We would like to study the E 8 model also known as the 2d Ising model with magnetic deformation. We assume that the system is described by three asymptotic states with masses m 1 = 1, m 2 = 2 cos(π/5), m 3 = 2 cos(π/30). (5.25) Notice that the E 8 model has actually eight asymptotic states but only three of them are below the two particle threshold 4m 2 1 . We consider the scattering of particle m 1 and allow all three poles due to particles m 1 , m 2 and m 3 .
We minimize the central charge in this setup not specifying the values of the trilinear couplings first. Unfortunately this turns out not to be enough to single out the E 8 model. We further specify the value of the very first trilinear couplings g 1 given in (4.37). The values of g 2 and g 3 are obtained during the central charge minimization procedure and match precisely the ones in (4.37). As a result of our numerical procedure we also obtain the partial amplitude shown on figure 8 and the two particle form factor of the trace of the stress-tensor shown figure 9. The partial amplitude perfectly matches the analytic expression (4.29). This was expected since in [4], it was shown that this is the unique amplitude with maximal trilinear coupling g 1 given in (4.37). The form factor however does not match the analytic expression (4.30) with the coefficients (4.35). It matches however the following linear combination of the σ and form factors F Θ 2,fake (s) ≈ −0.79 F σ 11 (s) + 4.06 F 11 (s), (5.26) given by (4.30) with the coefficients (4.33) and (4.34). We refer to this as the fake trace of the stress-tensor form factor. The appearance of such an object is due to a peculiar situation when the form factors of different scalar operators cannot be easily distinguished (because they have the same large s-behaviour). In order to distinguish them one needs a more complicated setup which includes scattering of at least two different asymptotic states. Our results are however consistent with the analytical ones due to (5.26). It is not a surprise then that the central charge value of the fake form factor c ≈ 0.04945 does not correspond to the expected values (4.41) which follows from the analytic results.

Zero poles and O(N ) global symmetry
We consider here a single asymptotic state with mass m = 1 which transforms in the vector representation of the O(N ) global symmetry group. We further assume that there are no poles in the scattering amplitude of such states. (For previous works on 2d QFTs with O(N ) symmetry see [50][51][52][53]). In section 3.2 we have shown that this amplitude decomposes into three amplitudes in the trivial, symmetric traceless and antisymmetric representations. We now write the following ansatz for the partial amplitudeŝ where u n and v n are some constants. Notice, that contrary to (5.7) we parametrize here the entire partial amplitude and not only its interacting part. We plug this ansatz into the system of crossing equations (3.49). It becomes a system of linear algebraic equation on the linear coefficients u n and v n . It can be used for example to express u A 0 and v S n , u A n , v A n for n ≥ 1 in terms of the unknown linear coefficients u • 0 , u S 0 and u • n , v • n , u S n for n ≥ 1. Plugging this solution back into (5.27) we obtain an automatically crossing symmetric ansatz. We demand then that the ansatz obeys the unitarity constraints (3.51) and (3.53). The form factor of the trace of the stress-tensor is defined in (3.52) provided it obeys the following normalization which follows from (3.46) and (2.122). Following [51][52][53] we can minimize the central charge fixing the values of partial amplitudes at the crossing symmetric point. Let us define We remind that the crossing symmetry requires σ * 3 = σ * 1 . Fixing the values (5.29) is equivalent to fixing the values u • 0 , u S 0 and u A 0 in (5.27) due to (3.45). We can now scan for example over σ * 1 and σ * 2 and minimize the central charge to obtain the 3d plot. The allowed values of σ * 1 and σ * 2 form the bounded domain shown in figure 10. 38 To decrease the amount of numerical computations we will focus here only on two sections of the plane σ * 1 and σ * 2 , namely The boundary values from left to right correspond to the "minus" NLSM and NLSM respectively. For instance the right boundary value matches the analytic results (4.47) and (4.48). At the right boundary we have reconstructed the partial amplitudes and the form factor. We have observed that they match very well the analytic results (summarized in section 4.3) for φ ∈ [0, 0.8π] but differ for φ ∈ [0.8π, π]. One reason for that is the almost linear growth of the form factor with s. According to (4.51) we have The ansatz (5.8) however does not reproduce such behavior well for any finite value of N max . As a result the central charge differs from the one expected in the NLSM. Further analysis is required in order to tame the NLSM numerically. For section 2 the optimization problem is feasible for σ * 1 ∈ [−0.328, 0.328]. The boundary values correspond to the periodic Yang-Baxter (pYB) solutions [52,54]. We observe numerically that the value of the central charge exhibits a divergence like behavior when 38 This is figure 7 of [53] which we reproduce here for the reader's convenience. Figure 7: Allowed space of S-matrices in the plane 1 (s ⇤ = 2) vs 2 (s ⇤ = 2). The coloring at its boundary matches the convention in table 1. We have also marked the points corresponding to known integrable S-matrices and the constant solution in (7).
It is worth emphasizing that these were by no means imposed and rather come out as a mysterious outcome. It is amusing to think that had Yang-Baxter not been discovered before and these nice integrable solutions not unveiled decades ago, we could have discovered them here in these numerical explorations.
• Another interesting point is the yellow point between free theory and NLSM in figure 7.
The S-matrix there is a simple constant solution to crossing and unitarity but does not obey Yang-Baxter equations. Notice that in the symmetric channel unitarity is not saturated. To our knowledge this is the first analytic solution to the S-matrix bootstrap problem where unitarity is not saturated. We call it the yellow point.
If we look for constant solutions to the bootstrap problem it is actually easy to derive (7) analytically. First, because of crossing, all possible constant solutions lie on the same plane as the slate (i.e. must be eigenvectors of the crossing matrix). The unitarity inequalities then define a polygon on this plane which is nothing but the innermost curve in figure 4. Such polygon is simply given by S a = C ab S b , |S a |  1 with S a constant. The vertices of this polygon are precisely (±) free theory and the yellow point. These are the only points that touch the boundary of the slate. (No other points could touch it since the slate is a convex space.) • As we move along the boundary we observe that all S-matrices saturate the unitarity condition at all values of energy except for the yellow point discussed above. Unitarity approaching the boundary. This might be a sign that the pYB solutions are unphysical. One has to however take into consideration poor convergence of the ansatz, thus further analysis is also required to make a definite statement. To conclude, let us also minimize the central charge for different values of N without fixing σ * 1 or σ * 2 . The result is presented on figure 12. The bound on the left part of figure 12 is almost linear and can be approximated well by c ≈ 0.644 + 0.334N . The values of σ * 1 and σ * 2 which realize the minimum of the central charge lie on the section 2 (5.31). The values of the optimal σ * 1 as a function of N are presented on the right part of figure 12.

Conclusions
In this paper we have extended the S-matrix bootstrap program to include states created by local operators. This gives rise to a bootstrap setup that mixes scattering amplitudes, form factors and spectral densities of local operators. The latter allows to extract direct information about the UV fixed point which was not possible so far in the pure S-matrix bootstrap approach. We have established the groundwork for future explorations and limited ourselves to testing the approach in two dimensional QFTs. Our main result is the derivation of a lower bound for the central charge c of the UV CFT that can flow to a massive phase with a given particle spectrum (and interactions). presence of a single stable particle transforming in the vector representation, we found the universal lower bound c ≥ c min (N ) ≈ 0.6 + 0.3N , see figure 12(a).
We should however be careful with the meaning of "lower bound". As in all the recent S-matrix bootstrap works [5], the "lower bound" decreases when increasing the number of parameters in the ansatz, i.e. N max in equations (5.7) -(5.9). Strictly speaking our result is an upper bound for the lower bound. In all the plots presented in this paper we have taken the value of N max sufficiently large. However, we have not performed a careful convergence analysis. In order to be rigorous, it is important to generalize the functional method of [53] to our setup and obtain rigorous lower bounds on the central charge.
Another direction worth further exploration, is the inclusion of several states created by different local operators. It seems natural to consider the full set of relevant operators of a given CFT. For example, in the 2d Ising model, it would be interesting to consider both σ and . This setup would include form factors for both operators and a 2 × 2 matrix of Wightman two-point functions.
In two spacetime dimensions, it would be interesting to further explore the connection with integrable models. In our setup, we have observed that the numerical optimization problem tends to saturate the conditions (3.27) -(3.29), which include Watson's equation and absence of particle production. Therefore, similarly to the pure S-matrix bootstrap, we found that the optimal solutions often correspond to integrable theories. In this work, we encountered amplitudes and form factors of the sine-Gordon, E 8 and O(N ) models.
In the presence of continuous global symmetries, it is natural to study states created by the conserved currents. Notice that form factors of conserved currents also have a natural normalization (at s = 0) following from the conserved charges. This should have interesting applications both in d = 2 and in higher dimensions. In d = 2, it seems clear that a detailed study of the O(N ) model with our approach will benefit from the inclusion of states created by the non-abelian currents. Moreover, it would be useful to obtain the 3D plot of the central charge lower bound above the allowed region in figure 10.
Let us conclude by discussing our new bootstrap method in higher dimensions. In d ≥ 3 one has to consider the Wightman two-point function of the full stress-tensor T µν and not only its trace Θ = T µ µ . Such two-point function can be decomposed into two spectral densities ρ Θ (s), ρ 2 (s), (6.1) where the first spectral density represents the trace of the stress-tensor exactly as in d = 2 and the second spectral density is the new object special to d ≥ 3. In all dimensions we have the following asymptotic behavior at large energies 39 In d ≥ 3 there is no known analogous integral expression and the value of the central charge is hidden instead in the asymptotics of the second spectral density in (6.1). More precisely lim s→+∞ ρ 2 (s) = const × C U V T s d/2 . (6.4) Here const is a numerical factor which depends on the precise definition of ρ 2 and is irrelevant for the present discussion. Reiterating, even though the spectral densities (6.1) allow to access the values of the central charge in any number of dimensions, the information about it is encoded differently in them in d = 2 and in d ≥ 3. Interestingly enough, the integral formula (2.112) captures both (6.3) and (6.4) in a uniform way. For more details see [25]. Our numerical bootstrap approach can easily bound integrals like (6.3) but it cannot bound coefficients in the asymptotic behaviour of spectral densities like in (6.4). For this reason, we cannot use our bootstrap method to put non-trivial bounds on the central charge C U V T in d ≥ 3. What we can do instead in higher dimensions is to put lower bounds on the following dimensionless quantity ∞ 0 ds s d/2+1 ρ Θ (s), (6.5) which is a simple generalization of (6.3) to higher dimension. We stress however that contrary to (6.3) this quantity is not directly related to a property of the UV CFT and is thus less interesting than (6.3), see [55] for further discussion. Focusing to d = 4 there is a more interesting quantity we can bound which is the a-anomaly. This requires however some modifications in the present formalism. We plan to address this question in the near future. Spherical coordinates We will need to evaluate n-dimensional integrals in Euclidean signature. It is best done in spherical coordinates which we introduce here. The n-dimensional spherical coordinates consist of the radius r and a set of n − 1 angles with the following ranges It is then straightforward to evaluate the spherical angle Ω n 40 Ω n = n π n/2 Γ(n/2 + 1) . (A.7) The n = 1 case is special. We do not have angles, however we have two (equivalent) directions with x ≥ 0 and x < 0. This fact is already contained in (A.7) since Ω 1 = 2. Finally, the spherical δ-function is given by δ (n−1) (Ω) = δ(θ 1 ) . . . δ(θ n−2 ) J(1; θ 1 , . . . , θ n−2 ) × δ(θ n−1 ). (A.8) 40 Notice that often the spherical angle is denote by Ωn−1 in n-dimensions.
Spherical harmonics Any normalizable function on the sphere S n−1 ⊂ R n can be decomposed in the basis of spherical harmonics. If the function is invariant under SO(n − 1) rotations (therefore only depends on θ 1 ) then it can be expanded in Gegenbauer polynomials C (n−2)/2 j (x), x ≡ cos θ 1 , (A.9) where j = 0, 1, 2, . . . is a non-negative integer which is often referred to as spin. In the n = 3 case, the Gegenbauer polynomials coincide with the Legendre polynomials P j (x) = C The Gegenbauer polynomials satisfy the following orthogonality property Defining the total d-momentum p µ ≡ p µ 1 + p µ 2 , (A. 15) we can write Equivalently we have Here all the energies p 0 i satisfy the mass-shell condition (2.6). Notice the presence of the 1/n! factor which removes overcounting of indistinguishable (identical) particles.
Wightman functions In Euclidean signature only the time-ordered correlation functions make sense. In Lorentzian signature we also have non time-ordered correlators known as the Wightman functions. 42 They are defined as follows 0|φ(t 1 , x 1 ) . . . φ(t n , x n )|0 ,t j ≡ t j − i j , The integration goes along the real q 0 values and thus the i 's are needed to avoid the presence of poles on the line of integration. We can now rotate the line of integration by +π/2. This is done in such a way that we do not cross the poles. It is represented by the change of variables which is known as the Wick rotation. 43 The subscript E stands for Euclidean. Applying the change of coordinates (B.2) to (2.58) we obtain the Euclidean propagator 42 See for example appendix B in [56]. 43 More precisely the rotation of the energy p 0 to purely imaginary values is done by p 0 → p 0 e iφ → ip 0 ≡ p 0 E , where the angle φ changes from 0 to +π/2. The rotation of the time is determined by the condition to keep the scalar product (p · x) = p 0 x0 + p · x real, which is important to require in order not to introduce any divergences in the Fourier transform. This leads to the following x0 → x0e −iφ → −ix0 ≡ xE0 = x 0 E .
where we have

C Free scalar theory
Let us consider the free field theory with a single real scalar field φ(x). It is defined via the Lagrangian density In the remainder of this appendix let us focus on the the case of d = 2. The trace of the stress-tensor then reads Θ(x) = −m 2 : φ(x) 2 : . (C.8) Let us compute the one and two particle form factors for the operator (C.8). In the former case we have F Θ 1 = −m 2 0|a p 1 : φ(0) 2 : |0 = 0. (C.9) In the latter case we have Notice that the two particle states reproduce the entire value of the central charge.