Branch Point Twist Field Form Factors in the sine-Gordon Model I: Breather Fusion and Entanglement Dynamics

The quantum sine-Gordon model is the simplest massive interacting integrable quantum field theory whose two-particle scattering matrix is generally non-diagonal. As such, it is a model that has been extensively studied, especially in the context of the bootstrap programme. In this paper we compute the form factors of a special local field known as the branch point twist field, whose correlation functions are building blocks for measures of entanglement. We consider the attractive regime where the theory posesses a particle spectrum consisting of a soliton, an antisoliton (of opposite $U(1)$ charges) and several (neutral) breathers. In the breather sector we exploit the fusion procedure to compute form factors of heavier breathers from those of lighter ones. We apply our results to the study of the entanglement dynamics after a small mass quench and for short times. We show that in the presence of two or more breathers the von Neumann and R\'enyi entropies display undamped oscillations in time, whose frequencies are proportional to the even breather masses and whose amplitudes are proportional to the breather's one-particle form factor.


Introduction
The quantum relativistic sine-Gordon model is a paradigmatic example of an integrable quantum field theory (IQFT) that is amenable to solution by the bootstrap program. It provides the simplest example of a theory that is interacting and has a non-diagonal S-matrix, famously obtained in [1]. This means that the theory allows for backscattering or, in a different language, the S-matrix is a non-trivial solution of the Yang-Baxter equation. The theory has a rich particle spectrum containing two fundamental particles known as the soliton (s) and the antisoliton (s) and a tower of breathers (b k ) which can be interpreted both as soliton-antisoliton bound states and as bound states of lighter breathers. The number and masses of these breathers depend on the model's coupling constant. Although the theory is non-diagonal in the standard scattering matrix sense, the breather sector is diagonal and this simplifies form factor computations considerably. In addition, in a certain coupling constant regime, the sine-Gordon model can be seen as the continuum limit of another paradigmatic integrable theory, namely the spin- 1 2 XXZ quantum spin chain.
In the context of the bootstrap program for IQFTs, the matrix elements of local operators (e.g. form factors) of the sine-Gordon model have been extensively studied by many authors. Some of the earliest results are due to F.A. Smirnov [2,3], followed by a long series of works by the Berlin group [4][5][6][7][8]. In both cases, the form factors are given in terms of integral representations and expanded on a basis of off-shell Bethe vectors. A different approach known as free field representation was employed in [11,12] and the fermionic structure of the model was exploited in [14,15], where on-shell Bethe vectors are proposed as a basis for the form factor decomposition. Form factors of the order and disorder fields were first obtained in [9] and applied in the study of the Ashkin-Teller model in [9,10]. Of particular interest to us is the work [12] which focused on breather form factors and used the fusion technique in order to obtain form factors of heavier breathers from those of the lightest one (note that this technique had already been employed much earlier in the work of Smirnov [3]). There has also been intense study and a large body of applications of sine-Gordon form factors in various other contexts such as the case of finite temperature one-point functions [16,17], quantum quenches [18][19][20] and boundary field theory [21] and, in particular, in finite volume [22][23][24] where once again fusion techniques can be employed.
Finally, it is important to note that many studies of the breather form factors (particularly those where fusion is used) exploit the relationship between the sine-Gordon and sinh-Gordon theories. At Lagrangian level the two theories are identical up to the complexification of the coupling constant. In addition, the two-particle S-matrix of the first (lightest) breather is mapped to the two-particle scattering matrix of the sinh-Gordon particle under the same transformation. This implies that the form factors of the first breather (and by fusion, also those of higher breathers) can be obtained from those of the sinh-Gordon field by simply changing the coupling constant dependence. Then the sinh-Gordon form factors computed in various papers [7,13,25,26] become the starting point of computations in the sine-Gordon model.
The works we have referred to so far are concerned with "standard" local fields of the sine-Gordon theory, such as the sine-Gordon field ϕ, its powers and, especially, exponential fields of the form e iaϕ which are of particular interest as they are related to the trace of the stress-energy tensor. In the present work our main aim is to generalize these results to branch point twist fields, starting with the branch point twist field and associated form factor program introduced in [27], and then continuing (in part II) with the symmetry resolved branch point twist field recently introduced in [28,29]. Twist field form factors of the sine-Gordon model were first studied in [30] but only in the so-called repulsive regime where no breathers are present. In this paper we extend those results to the situation when several breathers are present focussing on all non-vanishing one-and two-particle form factors. In the breather sector we employ the results of [27] and [31] where the two-and four-particle form factors of the sinh-Gordon field were obtained, respectively. These will constitute our starting point when employing the fusion procedure to obtain lower particle form factors of higher breathers.
The paper is organized as follows: In Section 2 we review some general results for the sine-Gordon model, notably its S-matrix and particle spectrum. In Section 3 we review the definition of the branch point twist field and the main equations satisfied by its form factors. In Section 4 we diagonalize the two-particle form factor equations to compute the two-particle soliton-antisoliton form factor. We put special emphasis on the discussion of its dynamical pole structure. In Section 5 we use fusion to compute one-and two-particle breather form factors, including up to four breathers and carry out some simple consistency checks of our solutions. In Section 6 we evaluate the ∆ sum rule in several coupling regimes, finding very good agreement with the exact value of the branch point twist field conformal dimension for all coupling choices. In Section 7 we discuss one application of our results to the study of the entanglement dynamics following a mass quench. We conclude in Section 8. The more technical details of our work are presented in various Appendices. Appendix A summarizes some useful formulae for the minimal form factors. Appendices B and C give details of the computation of breather form factors for the branch point twist field and the trace of the stress-energy tensor, respectively. In both cases we use the fusion procedure. In Appendix D we analyse in more detail the dynamical pole axiom for the soliton-antisoliton form factors. In Appendix E we present some additional numerical results concerning our evaluation of the the ∆ sum rule.

Main Features of the Model
The sine-Gordon model is characterized by the following euclidean action where g and µ are coupling constants and ϕ is a scalar field. As anticipated in the introduction, this action becomes that of another theory, know as sinh-Gordon model under the mapping g Þ Ñ ig with g P R.
For generic values of the coupling, the theory has a rich particle spectrum consisting of a soliton (s) and anti-soliton (s) of opposite U p1q charge and a family of bound states known as breathers. Defining the new coupling we have that the masses of the breathers take the form where m is the mass of the soliton and the anti-soliton and pξq " 1 ξ´1 if 1 ξ P Z and r 1 ξ s otherwise, where r¨s denotes the integer part. The mass m is related to the couplings µ and g through the mass-coupling relation first found in [32]. There are various interesting regimes: • For ξ ą 1 there are no bound states and the full spectrum consist only of the soliton and the antisoliton. This is called the repulsive regime. In this regime, the theory is equivalent to the massive Thirring model, a perturbation of the massive Dirac theory that preserves the U p1q symmetry. We studied the entanglement entropy in this particular regime in [30].
• The point ξ " 1 is special as can be seen more precisely from the S-matrices given below. From (5) we have that S ss ss pθq " Sss ss pθq "´1 and also S ss ss pθq "´1 and Ss s ss pθq " 0. At this point the theory becomes a Dirac free fermion.
• For ξ ă 1 the model is in the attractive regime were bound states (breathers) are formed with the masses (3).
• In particular, whenever 1 ξ " n, with n P Z`the non-diagonal scattering amplitude Ss s ss pθq " 0 is vanishing and the theory becomes diagonal. In fact, it reduces to the D n -minimal Toda field theory.
where S ss ss pθq and Ss s ss pθq are the off-diagonal amplitudes. Useful linear combinations are S`pθq " S ss ss pθq`S ss ss pθq S´pθq " S ss ss pθq´S ss ss pθq .
The remaining S-matrices are diagonal and can be expressed in terms of the standard blocks: for k ă p and, finally A good summary of all the S-matrices, and of how to derive Gamma-function representations from integral representations can be found for instance in [4].

Branch Point Twist Fields in a Nutshell
It has been known for some time that several entanglement measures, including the Rényi entropies, can be expressed in terms of correlation functions of a special class of local fields T which have been termed branch point twist fields in [27]. Branch point twist fields are, on the one hand, twist fields in the broader sense, that is, fields associated with an internal symmetry of the theory under consideration [27], and on the other hand related to branch points of multi-sheeted Riemann surfaces [35]. They are twist fields associated to the cyclic permutation symmetry of a model composed of n copies or "replicas" of a given theory, characterized by the exchange relations where O i pyq is any local field on copy number i, and with O n`1 pyq " O 1 pyq. The idea of quantum fields associated with branch points of Riemann surfaces in the context of entanglement appeared first in [35]. The general picture of branch point twist fields as symmetry fields associated to cyclic permutation symmetry of the n Riemann surface's sheets, as per (17), was given in [27], where they were studied in massive IQFT. This description is however independent of integrability, and it was first used in massive QFT outside of integrability in [36].
Cyclic permutation symmetry is not naturally present in most IQFTs, but can be "manufactured" by considering a replica model, composed of n copies of the original QFT (e.g. the sine-Gordon model). The connection to replica theories and multi-sheeted Riemann surfaces arises from the explicit formulae for entanglement measures, which generally depend on the quantity Tr A pρ n A q where ρ A is the reduced density matrix associated to a particular region A of the system. It is possible to show that the quantity Tr A pρ n A q is proportional to a correlation function of branch point twist fields involving as many twist field insertions as boundary points between the region A and the rest of the system. We will see an application of these ideas in Section 7 where we discuss the application of our results to the computation of the entanglement dynamics.

Form Factors and Form Factor Equations
Starting with the exchange relations (17), in IQFT one can formulate twist field form factor equations which generalize the standard form factor program for local fields [3,37]. These equations were first given in [27] for diagonal theories and then in [30] for non-diagonal ones. They have been generalized to symmetry resolved branch point twist fields in [28,29]. We will not review all these equations and their properties here but only those relations that are repeateadly used in the current paper, in particular the equations for one-and two-particle form factors. Let us start by defining F a 1 ...a k pθ 1 ,¨¨¨, θ k ; ξ, nq :" n x0|T p0q|θ 1 ,¨¨¨, θ k y a 1 ...a k ;n , to be a k-particle form factor, that is, a matrix element of the field between the vacuum state and a k-particle state. Here n x0| represents the vacuum state and |θ 1 ,¨¨¨, θ k y a 1 ...a k ;n represents an in-state of k particles with rapidities θ 1 , . . . , θ k and quantum numbers a 1 . . . a k , both in the replica model. These quantum numbers generally contain two indices, one for the particle type and one for the copy number. However, in our computations we will generally restrict ourselves to a single copy and will therefore drop the copy index. This is because form factors of other copies can be obtained from these solutions by repeated use of the form factor equations. The branch point twist field is a neutral field in relation to the sine-Gordon U p1q-symmetry that exchanges soliton and anti-soliton. This implies the vanishing of any twist-field form factors involving a different number of solitons and anti-solitons. At the one and two-particle level this means that F ss pθ; ξ, nq " Fsspθ; ξ, nq " Fs b k pθ, ξ; nq " F sb k pθ, ξ; nq " F s pξ, nq " Fspξ, nq " 0 , @ k P Z( 20) Here, we have used relativistic invariance and spinlessness of the twist field, which imply that the two-particle form factor depends on a single rapidity variable (the rapidity difference of the particles) and the one-particle form factor is rapidity independent. In addition, because of Z 2 symmetry we also have Under these considerations, Watson's equations for non-vanishing two-particle form factors and particles in the same copy can be summarized as F ss pθ; ξ, nq " S`pθqF ss p´θ; ξ, nq " F ss p2πin´θ; ξ, nq, whereas the kinematic residue equations aré iRes θ"iπ F ss pθ; ξ, nq "´iRes where xT y is the vacuum expectation value of the branch point twist field in the ground state of the replica theory. Finally, the bound state residue equations aré where c is any particle that is formed as a bound state of s`s for rapidity difference θ " iπu c ss . In the breather sector we will use the bound state residue equation extensively and repeatedly to obtain lower particle form factors of heavier breathers in a process known as "fusion". For this reason it is convenient to write the more general equatioń iRes θ"θ 0 where a 1 , . . . , a k are any particle combination for which the form factor is non-vanishing and u`ũ " u i`j ij where θ " iπu i`j ij is the pole of the scattering matrix S b i b j pθq corresponding to the formation of breather b i`j . Similarly, u andũ are related to the poles of S b j b i`j pθq and S b i b i`j pθq.

Soliton-Antisoliton Form Factors
In the following we summarise the necessary formulas for the two-particle soliton-antisoliton form factors of the branch point twist field. Although these quantities were already derived in [30], the formulas were strictly speaking only justified in the repulsive regime of the sine-Gordon model. As we show below they are, nevertheless, valid in the attractive regime as well once a proper analytic continuation in the parameter ξ is considered. Let us first discuss the minimal form factor of these objects, which we denote by Gpθ; ξ, nq. This is the "minimal solution" to using Eq. (22) which can be constructed in the manner shown in [27], which itself generalizes a standard method in the context of the form factor program (see e.g. [38]). This method takes as starting point the S-matrix involved in the middle identity (S`pθq in the present case) of (22) and assumes that it admits a representation of the type for some function gptq. If such a representation exists, then a minimal solution the equation (22) is given by where N is a normalization constant. To obtain the minimal form factor Gpθ; ξ, nq of interest, we therefore need to write S`in the form (27). This is straightforward since with spθq :" sin π´iθ The function S ss ss pθq already has an exponential representation (5) and one can easily write a similar representation for the function spθq as well An important remark is that the integral (5) is convergent for any 0 ă ξ ă 1. Nevertheless (31) is only convergent for 1 2 ă ξ ă 1. To be precise, for other values of ξ an alternative representation of the function above has to be used given by with p P Z`. Thus, we have two different representations of the minimal form factor Gpθ; ξ, nq depending on whether or not ξ´1 2 ą 0 or ξ´1 2 ď 0 which we denote by G˘pθ; ξ, nq, respectively. Interestingly, the value ξ " 1 2 is precisely the threshold for the formation of breathers and this is no coincidence. From the symmetry arguments presented in subsection 3.1 we know that the branch point twist field has vanishing one-particle form factors for odd-indexed breathers. However, the presence of non-zero one-particle breather form factors for even indices is allowed as we show later. This means that the two-particle soliton-antisoliton form factor of the branchpoint twist field must have bound state poles at imaginary rapidity values θ " iπp1´2kξq for k " 1, . . . , r 1 2ξ s. Equivalently, we can formulate this statement as the dynamical pole axiom (25) which we now specialise to even breather bound stateś i Res Notice that each new representation of spθq from (32), corresponds to a new breather with an even index entering the spectrum of the theory. Thus, when writing down the minimal part of F ss pθ; ξ, nq we have two alternative representations: if we employ the S-matrix representation (32) together with (5) and apply the standard machinery (28) to obtain the minimal form factor, the result possesses no breather bound state poles. This feature, is generally what is meant by "minimal solution". In this case the dynamical pole equation (33) can only be satisfied by multiplying the minimal form factor with another function which incorporates the required poles, similarly as for kinematic poles [27]. On the other hand, if one uses the analytically continued solution (31) instead of (32), the dynamical pole axiom (33) is automatically satisfied by the minimal form factor. In other words, this form factor is no-longer "minimal" in the standard sense, but includes also poles in the physical sheet corresponding to bound states.
Let us now continue our derivation for the minimal form factor, where the above discussed features can be explicitly demonstrated. The minimal form factor can be written as where the function Φpθ; ξ, nq follows from the integral representation of S ss ss pθq and can be written as or, alternatively, as an infinite product of Gamma functions: However, from a numerical viewpoint, the most useful representation is mixed, combining both a finite product of Gamma-fuctions and an integral part. This representation (86) is given in Appendix A. This kind of mixed representation was first used in [25] and is very rapidly convergent.
The function ϕpθ; ξ, nq in (34) follows from either the representation (31) valid for ξ ą 1 or from (32), which is instead valid for 1 2p ě ξ ą 1 2p`2 and p P Zφ´p As before, we can also write a mixed representations (see Eq. (87) and (88)). Similar to the discussion following (31)-(31), the minimal form factors G˘pθ; ξ, nq " ϕ˘pθ; ξ, nqΦpθ; ξ, nq , are two representations both satisfy Eq. (22), but whereas G`pθ; ξ, nq includes bound state poles at θ " iπp1´2kξq for k " 1, . . . , r 1 2ξ s, G´pθ; ξ, nq does not. Instead the necessary bound state poles can be introduced by simply dividing G´pθ; ξ, nq by standard CDD factors of the type A rigorous demonstration of this fact is presented in Appendix D. In this Appendix, the fulfilment of (25) with our soliton and breather form factors is proven as well, and we also derive some identities involving fractions of the minimal soliton-antisoliton form factors G˘pθ; ξ, nq and breather form factor Rpθ; ξ, nq (derived in the next section) based on (25). Now that we have found a minimal form factor that incorporates also the bound state poles, we just need to introduce the kinematic pole that ensures our solution satisfies (24). This kinematic pole can be introduced by multiplying with a function already presented in [27]. The final formulae for particles on the same copy are cos π n´c os πp1´2kξq n cosh θ n´c os πp1´2kξq n fi ffi fl G´pθ; ξ, nq G´piπ; ξ, nq .
We stress again that the two formulas are completely identical on the physical sheet and that the first line is the same expression derived for the repulsive regime in [27].

Breather Form Factors
In this section we focus on the breather sector of the theory, where the S-matrices are diagonal. The form factors can be easily obtained from known results for the sinh-Gordon model under the replacement B "´2ξ. With this identification one can then take the form factor solutions found in [27,31] and employ fusion to construct the chains of form factors and A nice example of this approach was given in Appendix A of [39] for the form factors of exponential fields.

Minimal Form Factor and Form Factors of b 1
Although we take the sinh-Gordon solutions as starting point, it is still useful to say a few words about the basic structure of those solutions, specially the minimal form factor. This function provides a minimal solution to the equations (23) for i " j " 1 and two breathers in the same copy. It can be easily adapted from the solutions presented in various papers [4,12,[23][24][25] and the techniques for the computation of minimal form factors introduced in [27]. The generalization to branch point twist fields of the representation given in [12] takes the form This function has the useful properties: Rpθ; ξ, nq " 1 and Rp0; ξ, nq " 0 .
A similar discussion as presented in the previous section also applies to this solution. First, although Rpθ; ξ, nq is constructed from the sinh-Gordon minimal form factor, it has very different analytic properties. Indeed, once more Rpθ; ξ, nq is not minimal, in the strictest sense of having no poles in the physical sheet. Rpθ; ξ, nq does have poles in the physical sheet, when the coupling allows for the the first breather to form higher breather bound states. Therefore, the solution (47) is valid for all values of the coupling ξ, with the function Rpθ; ξ, nq introducing bound state poles as needed. Second, the formula is once more only convergent for ξ ą 1 2 and this can be numerically addressed by employing the mixed representation (89). The full two-particle form factor is then given by The four-particle form factor can be read off from [31] and takes the form with Hpξ, nq " xT y 4ω 6 sin 2 π n n 2 Rpiπ; ξ, nq 2 , and Here σ i are the elementary symmetric polynomials on variables tx 1 , x 2 , x 3 , x 4 u and the coefficients q i pξ, nq where given in the Appendix of [31] (which unfortunately contains a typo). Calling cpaq :" cos πa 2n .
they can be rewritten as

Fusion Procedure
In this section we present the results of the fusion procedure as described schematically in (43). The simplest form factor to be obtained from the bootstrap approach outlined before is F b 2 pξ, nq.
The breather b 2 is a bound state of two b 1 breathers corresponding to the simple pole of S b 1 b 2 pθq at θ " iπξ. The bound state residue equation simply tells us that i Res We also know that the minimal form factor Rpθ; ξ, nq satisfies the equation and so, at the pole we have that i Res θ"iπξ Rpθ; ξ, nq "´i Res Putting all factors together, this gives the formula For n Ñ 1 the form factor vanishes as expected (since the twist field becomes the identity if the replica number is 1). However the limit is non-zero. This limit plays a role in computations of the von Neumann entropy. Note that the breather b 2 is only present for ξ ă 1 2 . Fig. 1 shows the function (56) for several choices of ξ and n. . Right: The one-particle form factor F b 2 pξ, nq as a function of ξ for n " 2 (red), 5 (blue), 12 (green), 50 (magenta).

Higher Breather Form Factors
Let us now consider a more involved fusion-based computation, namely that giving the form factor F b 2 b 1 b 1 pθ 1 , θ 2 , θ 3 ; ξ, nq from the four-particle form factor (48). The key equation in this case is´i Considering the formula (48) we see once more that the pole will originate from of the R-factors in the product, giving the contribution (55). More precisely, we obtain a solution of the form where β " e´i πξ n and Q 211 px 1 , x 2 , x 3 ; ξ, nq is obtained from evaluating Qpx 1 β´1 2 , x 1 β 1 2 , x 2 , x 3 q which simplifies with part of the denominator of (48) giving and As for the constant, we obtain Rp´iπξ; ξ, nq Rpiπ; ξ, nq 2 " Having now seen two applications of the fusion procedure it is easy to proceed for other form factors. We present more details of those computations in Appendix B. Here we just summarize the main formulae: with H 31 pξ, nq and Q 31 px 1 , x 2 ; ξ, nq given in (93), (92), respectively.
with H 22 pξ, nq and Q 22 px 1 , x 2 ; ξ, nq given by (104) and (102) and, finally F b 4 pξ, nq " xT y sin π n sin π 2n p1`2 cos πξ n q cos πp1´ξq Rp´3πiξ; ξ, nqRp´2πiξ; ξ, nq 2 Rp´iπξ; ξ, nq 3 Rpiπ; ξ, nq 2 , which, as we see in Appendix B can be obtained from either fusing b 3 and b 1 or b 2 with itself, giving identical results. Before ending this section, it is worth noting that Watson's equations and the bound state residue equation for form factors can be repeatedly used to obtain the form factors of breather b k starting with a form factor involving k breathers of type b 1 in a more systematic manner. This technique is described for instance in equation (A.3) of Appendix A in [40]. This method would allow us for instance to reduce (48) to (64) by simultaneously fusing all particles. The result is the same as presented here.

Some Consistency Checks
Apart from the ∆ sum rule that we will discuss later, there are a few properties that the form factors must satisfy and which help us make sure these formulae are correct. One of the strongest tests is the clustering decomposition property which states that in the absence of internal symmetries, form factors factorize into products of lower particle number form factors if a subset of the rapidities is sent to infinity. More precisely, for the form factors above we expect that and lim θ 1 Ñ8 (66) These identities can be easily checked thanks to the first property in (46). The first property in (65) follows from observing that for θ 1 Ñ 8 the denominator of the form factor (62) scales with x 4 1 whereas the numerator (that is, the function Q 31 px 1 , x 2 , x 2 ; ξ, nqq scales as x 3 1 . A similar argument applies to the second equality in (66). The second identity in (65) follows from lim θ 1 Ñ8 and lim θ 1 Ñ8 Rpθ 12 ; ξ, nq 2 Rpθ 12`i πξ; ξ, nqRpθ 12´i πξ; ξ, nq px 1´ω x 2 qpx 2´ω x 1 qpx 1´ω βx 2 qpx 2´ω βx 1 q " 1 ω 2 βx 4 1 (68) together with the formula (104). The first identity in (66) follows from Comparing with (47) and (56) we find that the clustering property is exactly reproduced. We may also check that the solution F b 2 b 2 pθ; ξ, nq satisfies the kinematic residue equation (24) which indeed it does. This can be shown by employing the non-trivial identity Rp´iπξ; ξ, nq 2 Rpiπp1´ξq; ξ, nqRpiπp1`ξq; ξ, nq " n tan πξ 2n sin πp1`ξq 2n sin πpξ´1q 2n 2ω 4 sin π 2n sin πp1`2ξq which can be established with the help of the Γ-function representation given in Appendix A and along similar lines as the proofs presented in Appendix D. We conclude by noting that similar consistency checks for the form factors of local fields in the sine-Gordon model were performed also in [10].

Consistency Checks by ∆ Sum Rule
The ∆ sum rule [41] is one of the most useful and common methods for testing form factor solutions. It gives a relationship between the conformal dimension of a local field T and a certain integral involving the two point function n x0|T p0qΘprq|0y n where Θ is the trace of the stress-energy tensor and |0y n is again the vacuum state in the replica theory. In its integrated form given for instance in [42] and after generalizing to branch point twist fields, the rule can be expressed as follows: where ∆ T " c 24 pn´1 n q is the conformal dimension of the branch point twist field [27,35,43] and we have now added a superindex to the form factors to indicate the quantum field they correspond to. The second sum is over all possible choices of particle types a p with masses m p .
As usual with this type of expansion, convergence of the sum is expected to be quick, and the main contributions come from the one-and two-particle form factors. Hence, if we can show such near saturation we can be confident that our form factors solutions are correct.
Let ∆ p pξqq T be the conformal dimension of the branch point twist field as given by (72) in the regime where pξq breathers are present. Although the exact value of ∆ T is independent of ξ the number and contribution of the terms in the sum changes substantially depending on the coupling. In what follows we present numerical results for the sum above for pξq " 1, 2, 3 and 4. For this we need first to obtain the one-and two-particle breather form factors of the stress-energy tensor in the sine-Gordon model. This can be done in a similar fashion as for the branch point twist fields, namely starting from the sinh-Gordon solutions presented in [25] and carrying out the fusion procedure. The results are presented in appendix E.  Table 1: The contributions to the sum (73) from the soliton-antisoliton term (ss) and the breather-breather term (b 1 b 1 ) for ξ " 0.62734 (left) and ξ " 0.82734 (right). The first column shows the exact values of ∆ T and the last column the sum of ss and b 1 b 1 contributions. As expected, the main contribution comes from the ss term. This contribution gets larger as we approach the threshold value ξ " 1, while the breather contribution is reduced.
Let us consider the regimes when there are one, two, three or four breathers present we have that the expansion above can be approximated as follows: • For ξ ą 1 we are in the repulsive regime where no breathers are present. The main contribution to the ∆ sum rule comes from the soliton-antisoliton form factor and was computed in [30].
• For 1 2 ď ξ ă 1 we have a single breather b 1 present and the main contributions are The sum for two values of ξ is presented in Table 1.
• For 1 3 ď ξ ă 1 2 we have two breathers b 1 , b 2 present and the main contributions are Numerical values of the sum (74) and of individual contributions to it are presented in Table 2 of Appendix E.
• For 1 4 ď ξ ă 1 3 we have three breathers b 1 , b 2 , b 3 present and the main contributions are • Finally, for 1 5 ď ξ ă 1 4 we have four breathers b 1 , b 2 , b 3 , b 4 present and the main contributions are (76) Table 3 gives an example of the evaluation of the sum (76), albeit without including the b 2 b 4 , b 3 b 3 and b 4 b 4 contributions, which we have not evaluated in this paper. Even so, the sum rule is approximately 95% saturated.
In conclusion, our numerical evaluation of the ∆ sum rule in various regions of the attractive regime shows near saturation upon inclusion of all relevant one-particle and two-particle form factors and therefore provides strong backing for our analytical results. It is interesting to note that the deeper we go into the attractive regime (i.e. the smaller ξ is) the more significant breather contributions are, so that for instance, in Table 2(d) the soliton-antisoliton contribution represents only about 20% of the total value of the dimension.

Application: Entanglement Oscillations after a Mass Quench
An interesting application of our results is to the study of the entanglement dynamics of the sine-Gordon model after a global mass quench [44,45]. That is, we want to study the timedependence of a certain measure of entanglement when the mass scale m is abruptly changed at time zero. Then, if the original hamiltonian of the system was Hpmq and m was the pre-quench soliton mass, at times t ą 0 the system will time-evolve with a new Hamiltonian Hpmq, wherê m is the post-quench soliton mass. In such a situation, the reduced density matrix may be formally written as: where A and B are two complementary regions and |0y is the pre-quench ground state. In terms of ρ A the Rényi and von Neumann entropies are defined in the usual form: and if A is a semi-infinite region, these expressions are equivalent to: and its n Ñ 1 limit, where ε is a non-universal UV cut-off which can be eliminated by considering instead the quantities ∆S n ptq :" S n ptq´S n p0q .
and |0y n is the pre-quench ground state in the replica theory. Note that S n p0q is a function of the vacuum expectation value n x0|T p0, 0q|0y n which we have abbreviated as xT y in our form factor formulae. With these definitions, the situation we want to consider here is entirely analogous to the studies performed in [46,47]. In fact, the present model has two key common features with the minimal E 8 Toda field theory studied in [47]. They are the presence of non-vanishing one-particle form factors and a mass spectrum where all masses are proportional to a fundamental scale m (the mass of the soliton/antisoliton). Carrying out the quench perturbation theory proposed in [48], non-vanishing one-particle form factors inevitably lead to entanglement oscillations at first order in perturbation theory. As observed in [46,49] the dynamics of entanglement is closely tied to the dynamics of the one-point function of the order parameter. Indeed, oscillations of the one-point function of the order parameter in the sine-Gordon model, following a mass quench were found in [50] employing perturbation theory.
The formulae involved are almost identical to those presented in [47], specially in the supplementary material. We must just highlight that the field associated with the mass quench in this case is the perturbing field in the sine-Gordon theory, namely the field Ψ " 2 cos gϕ where g is the coupling we first encountered in the action (1). This field is, as usual, proportional to the trace of the stress-energy tensor, hence its form factors are identical to those of Θ up to a proportionality constant (essentially, we need to replace xΘy " 2πm 2 1 with xΨy). Let us consider a perturbation where the original coupling µ in the action (1) is changed by a small amount δ µ , that is µ Þ Ñ µ`δ µ with δµ µ ! 1. Then, the Opδ µ q contribution to the Rényi entropies may be expressed as a series in form factors of T and Ψ, where the leading contributions to the increment of the Rényi entropies, come from one-and two-particle form factors. After various simplifications, the series takes the form ∆S n ptq " 1 1´n Re " rF Ψ ss p2θ; ξqs˚F T ss p2θ; ξ, nqe´2 imt cosh θ ı 2 cosh 2 θ r k "m k m are the scaled post-quench breather masses. The "prime" symbol in the last sum indicates the additional restriction that only terms where k and p are either both even or both odd will be non-vanishing. The "hatted" form factors are scaled versions of the usual form factors where the expectation values of the associated fields have been factored out. This dependency can then be absolved into the ratio of couplings δ µ {µ and the constant C Ψ . The conformal dimension ∆ Ψ " g 2 " ξ 1`ξ and the constant These are the standard scaling laws for vacuum expectation values and the mass-coupling relation. A relationship between the constant A Ψ and κ can be read off from the paper [51] where the expectation values of exponential fields in the sine-Gordon model were obtained. From this formula it follows that It is important to note that this formula is only convergent for g 2 ă 1 2 , which excludes the repulsive regime [51]. The mass-coupling relation was given earlier in (4). This allows us to fix the ratio above to Despite the messy nature of the formula (81) (a very similar formula can be written for the von Neumann entropy) the main features of entanglement are rather clear: for small quenches, there will be undamped oscillations whenever any one-particle form factors are non-vanishing, confirming the general ideas observed in [47,48]. In addition, there will be additional oscillatory terms coming from higher particle form factors which will be suppressed by a power of t that depends on the leading behaviour of the form factors near zero rapidity (this can be analysed further by using a saddle-point approximation). This means that the dynamics of entanglement following a mass quench is rather different in the regime ξ ď 1 2 (undamped oscillations with at least two breathers present) and for for ξ ą 1 2 (damped oscillations with at most one breather present).
We demonstrate these qualitative differences in the entanglement evolution by evaluating (81) numerically for various values of n and two particular values of ξ. In Figure 2 ∆S n ptq is displayed for n " 2, 3, 4, 5 and for ξ " 0.810361 and ξ " 0.420712. Clearly, above the second breather threshold (ξ " 0.810361) no undamped oscillations can be seen, which are, very clearly present when the second breather joins the spectrum (ξ " 0.420712). Note that although ∆S n p0q " 0 by the definition (80), it is not exactly zero numerically (although it is rather small). This is because S n ptq is evaluated at first order in perturbation theory and therefore its value at zero is only an approximation of the exact analytic value S n p0q that is subtracted in (80). Note that these results are only expected to hold for small quenches and times t ă µ´1, as explained in [48] and also discussed in [52,53]. In addition, we know that for large times the leading feature of entanglement (in any regime) is linear growth [44,45,54]. As observed in other studies, this feature is not recovered using first-order perturbation theory as it is a second order effect [46,47]. In addition, in some cases, like for E 8 Toda field theory, linear growth is very slow so that it only becomes apparent for large times in numerical simulations [47]. The same phenomenology is also observed for the von Neumann entropy for the same reasons.
It is worth considering whether or not these behaviours will persist for larger times and quenches. In this regard, arguments have been put forward as to why the undamped oscillations found at first order should be damped when including higher order terms [19,52]. At the same time, we know of at least one theory, E 8 Toda field theory, where this eventual damping is not observed numerically even for large quenches and times [47]. This suggests that this phenomenology still needs to be better understood. Similar behaviours have been observed in [19,20,50,52] for the expectation value of the field Ψ and σ and in [18] for two-point functions of the field ϕ.

Conclusion
In this paper we have carried out an in-depth study of low particle-number form factors of the branch point twist field in the sine-Gordon model in the attractive regime. We have considered up to four breathers in the spectrum and focused on the one-and two-particle form factors. Our work extends results for the repulsive regime that were presented in [30].
Although computations are generally tedious, great simplification comes from the presence of U p1q symmetry in the soliton-antisoliton sector and Z 2 symmetry in the breather sector. The latter imply the vanishing of a large number of form factors so that only form factors containing the same number of solitons and antisolitons as well as an even number of odd breathers are non-vanishing. The two-particle form factors of the soliton-antisoliton sector can be computed by diagonalizing the form factor equations (as the theory is non-diagonal) and incorporating the correct structure of bound state and kinematic poles, as discussed in great detail in Section 4 and Appendix D. For the breather sector the combination of the fusion mechanism with the analytic continuation from sinh-Gordon, provide an effective way of constructing the form factors of heavier breathers from those of lighter ones leading to the results of Secion 5 and Appendices B and C.
Our form factors can now be employed to compute correlation functions of branch point twist fields, hence a number of entanglement measures. In this paper we have highlighted just one such application, namely to the study of the entanglement dynamics following a mass quench in the sine-Gordon model. As observed in a similar study [47] we find that at least for small quenches, undamped oscillations of frequencies proportional to the even breather masses are present and constitute the leading behaviour of Rényi and von Neumann entropies. This is analogous to results found in [19,20,50,52] for the expectation value of the field Ψ and σ and more generally in [18] for two-point functions of the field ϕ.
As anticipated in the introduction, our immediate goal now is to extend these results to symmetry resolved twist fields [28,29].
Acknowledgment: We are grateful to Benjamin Doyon, Gábor Takács and Jacopo Viti for useful discussions. We specially thank Pasquale Calabrese for discussions and for his early stage involvement in this project. DXH acknowledges support from ERC under Consolidator grant number 771536 (NEMO).

A Minimal Form Factors Mixed Representations
As mentioned in Sections 4 and 5 the most useful representation for the minimal form factors is neither exponential nor based entirely on Gamma-function products, but a mixture of the two. This idea was employed first in [25] and can be implemented in a similar way for any minimal form factor of the type described in this paper. The function (36) can be written as Whereas the functions ϕ˘pθ; ξ, nq can be written as This quantity is independent of the choice of N only when ξ ą 1{2, but can be made convergent for any positive real values of ξ if the minimal allowed value for N is suitably chosen. In fact (87) gives a physically motivated (as seen in Appendix D) and correct analytic continuation. Alternatively, for 1 2p ě ξ ą 1 2p`2 and p P Zφ´p Finally, the breather-breather minimal form factor also admits the representation Using fusion again we have that This gives a solution of the form F b 3 b 1 pθ 12 ; ξ, nq " H 31 pξ, nqQ 31 px 1 , x 2 ; ξ, nq Rpθ 12 ; ξ, nqRpθ 12`i πξ; ξ, nqRpθ 12´i πξ; ξ, nq px 1´x2 ωβqpx 2´x1 ωβqpx 1 ω´x 2 qpx 2 ω´x 1 βq .(91) The polynomial Q 31 px 1 , x 2 ; ξ, nq follows from the reduction of Q 211 px 1 β´1 2 , x 1 β, x 2 ; ξ, nq and can be written as if we also identify H 31 pξ, nq "´xT y 2ωβ sin π 2n sin π n cos πp1´ξq 2n p1`2 cos πξ Rp´2πiξ; ξ, nqRp´iπξ; ξ, nq 2 Rpiπ; ξ, nq 2 .(93) which gives F b 4 pξ, nq " xT y sin π n sin π 2n p1`2 cos πξ n q cos πp1´ξq Rp´3πiξ; ξ, nqRp´2πiξ; ξ, nq 2 Rp´iπξ; ξ, nq 3 Rpiπ; ξ, nq 2 , which is plotted in Fig. 3 as a function of ξ and n.  Figure 3: Left: The one-particle form factor F b 4 pξ, nq as a function of n for ξ " 0.24 (red), 0.21 (blue), 0.2 (green), 0.1 (magenta) and 0.05 (brown). Right: The one-particle form factor F b 4 pξ, nq as a function of ξ for n " 2 (red), 5 (blue), 12 (green), 50 (magenta).
We have also that Note that the breather b 4 is only present for ξ ă 1 4 .
This is identical to the result we obtained from fusing F b 3 b 1 pθ; ξ, nq with the identification

C Form Factors of the Stress-Energy Tensor from Fusion
The form factors of the trace of the stress-energy tensor in the sinh-Gordon model were first computed in [25], where closed formulae for special values of the coupling were obtained. We are interested in the case of generic coupling B for which solutions up to 14 particles were given. These solutions will be the building blocks for our fusion procedure. We are particularly interested in the one-particle form factors of the second and fourth breather which requires the two-and four-particle form factors of the stress energy tensor in sinh-Gordon. Replacing B "´2ξ these form factors are given by where m 1 is the mass of the first breather as given in (3) and with x i " e θ i and σ i the elementary symmetric polynomial on variables tx 1 , x 2 , x 3 , x 4 u .

C.1 Computation of F Θ b 2 pξq
Applying the fusion procedure to the two-particle form factor we have that i Res and we get simply which is plotted in Fig. 4.

C.4 Computation of F Θ b 4 pξq
By computing the residue´i Res which gives A plot of F Θ b 4 pξq as a function of ξ is presented in Fig. 5.

C.5 Computation of F
The last two-particle form factor that we can obtain starting with (109) is F Θ b 1 b 3 pθ; ξq, resulting from the fusion process:´i

D Dynamical Poles of the Soliton-Antisoliton Form Factors
In this appendix we first show that the two representations G˘pθ; ξ, nq of the minimal solitonantisoliton form factor are indeed proportional to each other when the proper CDD-factors accounting for the bound state poles are introduced. We also demonstrate the precise working of the dynamical pole axiom (33).
This function is compared to for any integer n ě 1 and 1 2 ě ξ. To prove the above formula, we use the infinite product representations of Rpθ; ξ, nq and G´pθ; ξ, nq. First, one can easily show via (36) the following identity and so multiplying these two formula, we obtain a simple expression for the ratio of G-functions above. For the other ratio of R-functions, we can use (89) after sending N to infinity. This cancels the integral and leaves an infinite product equivalent to (45) but more suitable for our ongoing studies. We have, therefore, the following expression Rp´iπξ; ξ, nq Rpiπ; ξ, nq " (134) Then, putting these results together and carrying out some Gamma-function cancellations, the r.h.s. of (131) becomes Rp´iπξ; ξ, nqG´piπ; ξ, nq Rpiπ; ξ, nqG´piπp1´2ξq; ξ, nq " sin πξ 2n sin πp1´2ξq n Γ´1`ξ 2n¯Γ´1´ξ 2n¯Γ´ξ 2n¯Γ´´ξ 2nΓ´1`ξ´1 Using the standard identity ΓpxqΓp1´xq " π csc πx the first set of Gamma-functions can be simplified to Γ´1`ξ 2n¯Γ´1´ξ 2n¯Γ´ξ 2n¯Γ´´ξ 2nΓ´1`ξ´1 2n¯Γ´1´ξ`1 2n¯Γ´ξ`1 2n¯Γ´1´ξ 2n¯" while the infinite product can be easily rewritten using the standard identity Γpx`1q " xΓpxq yielding lim N Ñ8 where the l.h.s. can be evaluated analytically. For fixed N , the product can be expressed in terms of Gamma functions and Pochhammer's symbols and then the limit can be performed straightforwardly. Putting everything together, we find (131), as expected.

E ∆ Sum Rule Evaluation
In this Appendix we summarize our numerical results for the sum (72) and several distinct values of ξ and n. As discussed in Section 6 we include one-and two-particle contributions. In the regime 1 ě ξ ą 1 3 all the non-vanishing two-and one-particle contributions are taken into account as described in equations (73)-(74).
As we can see in Tables 1 and 2, the contribution from the b 1 b 1 and b 2 b 2 terms is very small compared to those of ss and b 2 . Assuming this tendency to hold for the contributions b 3 b 3 , b 2 b 4 and b 4 b 4 , in the interaction regimes 1 3 ě ξ ą 1 5 we have neglected the corresponding terms and still found good saturation of the rule (see Table 3).  Table 2: One-and two-particle contributions to the ∆ sum rule for four values of ξ P p 1 3 , 1 2 ). It is interesting to observe how the breather contributions become larger as ξ is decreased, sending the theory deeper into the attractive regime. For instance, in Table ( Table 3: One-and two-particle contributions to the ∆ sum rule for ξ " 0.22108. For this value of the coupling the first four breathers can be formed and approximately half the value of ∆ T comes from breather contributions. Even after neglecting the terms b 3 b 3 , b 2 b 4 and b 4 b 4 the rule is 95% saturated.