<?xml-model href='http://www.tei-c.org/release/xml/tei/custom/schema/relaxng/tei_all.rng' schematypens='http://relaxng.org/ns/structure/1.0'?><TEI xmlns="http://www.tei-c.org/ns/1.0">
	<teiHeader>
		<fileDesc>
			<titleStmt><title level='a'>Sine-Gordon model at finite temperature: The method of random surfaces</title></titleStmt>
			<publicationStmt>
				<publisher>American Physical Society</publisher>
				<date>04/01/2025</date>
			</publicationStmt>
			<sourceDesc>
				<bibl> 
					<idno type="par_id">10643245</idno>
					<idno type="doi">10.1103/PhysRevB.111.155112</idno>
					<title level='j'>Physical Review B</title>
<idno>2469-9950</idno>
<biblScope unit="volume">111</biblScope>
<biblScope unit="issue">15</biblScope>					

					<author>M Tóth</author><author>J H Pixley</author><author>D Szász-Schagrin</author><author>G Takács</author><author>M Kormos</author>
				</bibl>
			</sourceDesc>
		</fileDesc>
		<profileDesc>
			<abstract><ab><![CDATA[Not Available]]></ab></abstract>
		</profileDesc>
	</teiHeader>
	<text><body xmlns="http://www.tei-c.org/ns/1.0" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xmlns:xlink="http://www.w3.org/1999/xlink">
<div xmlns="http://www.tei-c.org/ns/1.0"><head>I. INTRODUCTION</head><p>Sine-Gordon (sG) quantum field theory is the low-energy (or, equivalently, long-distance) description of a large spectrum of gapped one-dimensional (1D) systems <ref type="bibr">[1]</ref>, with numerous applications ranging from quasi-1D antiferromagnets, carbon nanotubes and organic conductors <ref type="bibr">[2,</ref><ref type="bibr">3]</ref>, and trapped ultracold atoms <ref type="bibr">[4]</ref><ref type="bibr">[5]</ref><ref type="bibr">[6]</ref><ref type="bibr">[7]</ref><ref type="bibr">[8]</ref> to quantum circuits <ref type="bibr">[9]</ref> and coupled spin chains <ref type="bibr">[10]</ref>. It is also a paradigmatic example of quantum field theory with highly nontrivial nonperturbative dynamics and a strong-weak coupling duality <ref type="bibr">[11,</ref><ref type="bibr">12]</ref>. Due to its integrability, many exact results are available, including its exact S matrix <ref type="bibr">[13,</ref><ref type="bibr">14]</ref>, its thermodynamics via the thermodynamic Bethe ansatz <ref type="bibr">[15]</ref> and the nonlinear integral equation <ref type="bibr">[16,</ref><ref type="bibr">17]</ref>, its form factors <ref type="bibr">[18]</ref>, and exact zero-temperature expectation values of local fields <ref type="bibr">[19]</ref>. However, determining finite-temperature expectation values and, more generally, correlation functions presents serious challenges.</p><p>Correlation functions are among the most fundamental quantities characterizing spatial correlations and dynamics and are also. measurable in experiments. Besides their use in reconstructing physical quantities such as susceptibility or conductivity, cold atom experiments can directly access sG correlation functions using atom wave interferometry <ref type="bibr">[7]</ref>. Despite the integrability of the sine-Gordon model, the available theoretical approaches face severe limitations. Perturbative expansions only work for short distances/high temperatures. For zero-temperature correlators, a long-distance expansion can be written in terms of the exact form factors <ref type="bibr">[3]</ref>, which also enables the computation of the low-temperature behavior of expectation values of local fields <ref type="bibr">[20]</ref>. However, extending the form-factor expansion to finite temperature is much more limited in scope <ref type="bibr">[21,</ref><ref type="bibr">22]</ref>.</p><p>An alternative is to consider approaches which are independent of integrability. One possible approach is discretizing the fields, putting the theory on a space-(imaginary) time lattice, and using standard Monte Carlo techniques. The classical 2D sine-Gordon model was studied this way in connection with the BKT roughening transition <ref type="bibr">[23]</ref><ref type="bibr">[24]</ref><ref type="bibr">[25]</ref>. However, from the point of view of the 1D quantum sine-Gordon theory, these studies were carried out at zero temperature (infinite 2D geometry).</p><p>Another possibility is provided by Hamiltonian truncation <ref type="bibr">[26]</ref>, which was extended to the sG theory some time ago <ref type="bibr">[27]</ref>; for recent reviews see Refs. <ref type="bibr">[28,</ref><ref type="bibr">29]</ref>. It was recently utilized to evaluate the time evolution of correlation functions both in and out of equilibrium <ref type="bibr">[30,</ref><ref type="bibr">31]</ref> and provides a fully nonperturbative approach to quantum dynamics. Still, it is limited by the truncation to low-energy densities and by finite-size effects that appear due to the finite volume setting <ref type="bibr">[32]</ref>.</p><p>Here we consider an alternative numerical approach to sine-Gordon quantum field theory that relies on reformulating the problem in terms of random surfaces, which was previously used to compute full distribution functions of interference amplitudes <ref type="bibr">[33]</ref>. The method of random surfaces (MRS) is a type of Monte Carlo method, with the main advantage given by its simplicity and speed, with the necessary computing resources scaling polynomially with the number of modes retained in the computation. It is also free of the well-known sign problem and can access multipoint functions (including finite-temperature ones), which is a great advantage given their experimental availability <ref type="bibr">[7]</ref>. In addition, it is also possible to implement spatially inhomogeneous couplings, which are generally present in experimental realizations due to the traps used to localize the atomic clouds. Previous applications of the method include experimental analysis of variations of interference fringe contrast for pairs of independently created one-dimensional Bose condensates <ref type="bibr">[34]</ref>, and dissipative Landau-Zener transitions <ref type="bibr">[35,</ref><ref type="bibr">36]</ref>.</p><p>The present paper is devoted to developing the method for the sine-Gordon model and demonstrating its validity and efficiency, with results for the free-energy and operator expectation values validated by using the nonlinear integral equation (NLIE) as well as the minisuperspace extended truncated Hamiltonian approach (MSTHA) previously developed to study nonequilibrium dynamics in the model <ref type="bibr">[37]</ref>. Application of the MRS approach to more general correlation functions is relegated to a follow-up paper <ref type="bibr">[38]</ref>.</p><p>The outline of the present work is as follows. Section II introduces the model and the method of random surfaces, the latter in two stages: first for the partition function and then for the expectation values of vertex operators. Section III presents our numerical results and their validation using the NLIE and the MSTHA, while we give our conclusions and outlook in Sec. IV. To keep the main argument uninterrupted, some details were relegated to the Appendixes including technical aspects of the method of random surfaces, a short summary of the NLIE, and some additional numerical results.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>II. MODEL AND METHOD</head><p>In this section, we introduce the sine-Gordon model and provide a detailed derivation of the partition function and onepoint functions using the MRS approach.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>A. The sine-Gordon model</head><p>The sine-Gordon model is a 1 + 1-dimensional relativistic quantum field theory that at finite temperature T is defined by the Euclidean action</p><p>Here R = 1/T is the inverse temperature, 0 &#964; R is the imaginary time, and &#966;(x, &#964; ) is a real scalar field satisfying periodic boundary conditions in the &#964; direction: &#966;(x, &#964; + R) = &#966;(x, &#964; ). We work with units such that both the Boltzmann constant and the speed of light are one (k B = 1, c = 1). The model can be viewed as a perturbed conformal field theory where the free boson CFT with action S 0 is perturbed by the cos(&#946;&#966;) operator, which is relevant for 0 &lt; &#946; 2 &lt; 8&#960; . This work will focus on the so-called attractive regime 0 &lt; &#946; 2 &lt; 4&#960; .</p><p>The partition function of the model at an inverse temperature R can be written as a coherent state bosonic path integral Z = Z 0 Z I , where Z 0 = D[&#966;] e -S 0 is the partition function of the free boson, and</p><p>where we have introduced the average over free bosonic fields . . . 0 = D[&#966;] (. . . )e -S 0 /Z 0 . Using this representation, the partition function can be evaluated by the MRS introduced in Ref. <ref type="bibr">[33]</ref>. Here we develop the MRS further to compute the free energy F = -R ln Z and extend it to the one point functions e i&#945;&#966; . In both cases, we ascertain the accuracy and speed of the MRS when compared to other known available techniques to study the sG model.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>B. Method of random surfaces: Partition function</head><p>The MRS can be derived by considering a perturbative expansion of Z in terms of &#955; 0 ,</p><p>where r i = (x i , &#964; i ) and the subscript 0 refers to the expectation value in the unperturbed free boson theory as defined in Eq. ( <ref type="formula">2</ref>). Expanding the product leads to the sum of n-point functions of the operators e &#177;i&#946; &#966;(r i ) . The multipoint functions appearing in the integrand can be calculated by carrying out the Gaussian path integral, which leads to the multiplicative form of Wick's theorem:</p><p>if the so-called neutrality condition j &#945; j = 0 is satisfied and zero otherwise <ref type="bibr">[39]</ref>. Here</p><p>is the Green's function, i.e., a solution of the Poisson equation on an infinite cylinder of circumference R. Here we introduced a regulator length a, and (x jk , &#964; jk ) = (x jx k , &#964; j -&#964; k ) so the Green's function is translational invariant; moreover, it is periodic in the &#964; direction with period R. Plugging Eq. ( <ref type="formula">6</ref>) into Eq. ( <ref type="formula">5</ref>) we obtain</p><p>The neutrality condition ensures the result is finite in the R &#8594; &#8734; (zero-temperature) limit. It is useful to note that G(r j , r k ) can be shifted by an arbitrary constant without changing the result in Eq. ( <ref type="formula">7</ref>), which we will exploit shortly.</p><p>Substituting Eq. ( <ref type="formula">7</ref>) into Eq. ( <ref type="formula">3</ref>) and collecting the combinatorial factors we arrive at</p><p>where only the n = 2m even terms survive due to the neutrality condition. The R dependence is only through G(r j , r k ) = 4&#960; G(r j , r k ), = &#946; 2 /4&#960; is the scaling dimension of the operators e &#177;i&#946; &#966;(r) , and j = &#177;1, m of them being +1 and m being -1. The explicit a power in the correlation function <ref type="bibr">(7)</ref> has been absorbed into a multiplicative renormalization of the coupling constant,</p><p>If we define the so-called vertex operators as <ref type="bibr">[40]</ref> V&#945; (r) = a -&#945; 2 4&#960; e i&#945; &#966;(r) , <ref type="bibr">(10)</ref> then the perturbation can be written as</p><p>From Eq. ( <ref type="formula">7</ref>) we can read off the two-point function of vertex operators,</p><p>which for |r 12 | &#8594; 0 behaves as V &#945; (r 1 )V -&#945; (r 2 ) &#8764; a -&#945; 2 /(2&#960; ) (or without the regulator a, V &#945; (r 1 )V -&#945; (r 2 ) &#8764; |r 12 | -&#945; 2 /(2&#960; ) ), obeying the standard conformal normalization of scaling operators.</p><p>In the attractive regime &#946; 2 &lt; 4&#960; ( &lt; 1), the multiple integrals in Eq. (8b) are finite as a &#8594; 0 so we are allowed to drop a in the expression of G(r j , r k ).</p><p>The series in Eqs. ( <ref type="formula">8</ref>) can be interpreted as the grandcanonical partition function of a classical 2D system of particles with electric charge j that interact via the 2D Coulomb potential G(r j , r k ). In this picture, is related to the inverse temperature, while &#955;R -is related to the fugacity of the classical system. The Coulomb gas partition function in Eq. ( <ref type="formula">8</ref>) appears in the field-theoretical description of the 2D XY model, 2D superfluids, and thin superconducting films <ref type="bibr">[41]</ref>.</p><p>The first step in the MRS approach is to write G(r 1 , r 2 ) in a projector decomposition form (i.e., diagonalize with Fourier series)</p><p>This can be achieved by expanding the left-hand side in a double Fourier series. To obtain a discrete series, we restrict integration over the x j variables to the interval [-L/2, L/2] and extend the function G(r 12 , 0) on x 12 &#8712; [-L, L] periodically in the x direction with period 2L (see Appendix A for details). As a result, we integrate over a finite cylinder but using the Green's function of the infinite cylinder. This introduces a finite-size effect that must be addressed in the calculation. A necessary condition for L is that it must be much larger than the correlation length of the exponential decay of the two-point function <ref type="bibr">(12)</ref>, which gives</p><p>Since G is a real symmetric function, it can be expanded in a double cosine series, so</p><p>A mn cos(&#960; mx 12 /L) cos(2&#960; n&#964; 12 /R) <ref type="bibr">(14)</ref> which, on using trigonometric subtraction identities, leads to the form of Eq. ( <ref type="formula">13</ref>) where f = (m, n, &#8467;) is a multi-index with integers m, n &gt; 0, 1 &#8467; 4, and the eigenfunctions for m = 0, n = 0 are</p><p>If m = 0, n = 0 or n = 0, m = 0, then there are only two functions (sin and cos), and for m = n = 0 there is only one, the &#968; 0,0,1 constant term.</p><p>The coefficients G (m,n,&#8467;) = A mn can be obtained by numerical integration. For L/R 1, they turn out to be very close to the following simple expressions (cf. Appendix A):</p><p>and</p><p>The number of Fourier modes is truncated in a numerical calculation to a finite value m max , which naturally introduces an ultraviolet cutoff. From now on, we assume that the number of modes is finite and vary m max to ensure that the results have converged at a given L/R ratio and dimensionless coupling constant &#955;R -.</p><p>Using the decomposition in Eq. ( <ref type="formula">13</ref>) in Eq. (8b), we can write Z 2m as</p><p>The last term factorizes into a 2m-fold product in terms of the integration variables. However, due to the specific trigonometric form of the &#968; f functions in Eq. (15a), it is a constant equal to exp(-m m,n A mn ). To factorize the first term, we observe that all G f are positive except G 0,0,1 for which the whole exponent is zero due to the j j = 0 neutrality condition. To avoid numerical instabilities, we set G 0,0,1 to zero, which amounts to a constant shift of G(r 1 , r 2 ), which can be accounted for at the end of the calculation. We can thus apply a Hubbard-Stratonovich decoupling using the identity</p><p>e -y 2 2 &#177;ixy dy <ref type="bibr">(18)</ref> and introducing a Gaussian field t f for each mode, which yields</p><p>Swapping the integrals over the t f with the space-time integrals, the latter factorizes to a product of m + m identical integrals. This leads to</p><p>where</p><p>with a prefactor</p><p>The high-dimensional integral over {t f } can be performed in a Monte Carlo fashion by averaging the random t f variables drawn from a Gaussian distribution with zero mean and unit variance. In this picture, the integrand of g({t f }) is a random function of two variables, i.e., a fluctuating surface, hence the method's name. Plugging Z 2m in Eq. ( <ref type="formula">20</ref>) into the expansion (8a), and swapping the integral with the infinite sum, the latter can be carried out exactly to yield</p><p>where I 0 (x) is the modified Bessel function of the first kind. Remarkably, this sums up the entire perturbative series, exactly. Importantly, as demonstrated by the following numerical calculations presented in Sec. III, the MRS result fully captures all the nonperturbative contributions as well, which is a feature arising (somewhat miraculously) from the above manipulations.</p><p>In the numerical calculations, dimensionful quantities are measured in powers of R, i.e., one sets R = 1. We choose an L/R ratio and determine the coefficients A mn up to a maximal mode number m max . We then draw a set of random {t f } = {t (m,n,&#8467;) } from the Gaussian distribution with zero mean and variance equal to unity. We compute the 2D integral <ref type="bibr">(21)</ref> for g({t f }) by discretizing the mode functions &#968; f on a grid. We study the effect of the finite resolution of the grid in Appendix D 1. Z I is then obtained by averaging over N s different {t f } samples. Note that once the g({t f }) are computed and stored, the partition function can be easily evaluated for any value of the coupling &#955;, showcasing one powerful aspect of the MRS approach.</p><p>In the following, we consider aspect ratios of the twodimensional system in the range 6 L/R 12; this requires taking on the order of 20 m max 60 Fourier modes and the number of Monte Carlo samples 10 4 N s 10 6 . A major benefit of the MRS approach to the sG model is that it is simple to implement (the length of the necessary code in C is of the order of 50 lines), understand, and execute. Moreover, being a Monte Carlo technique, it is amenable to parallelization, which trivially reduces the running time. However, it is far from obvious how the accuracy of the MRS approach depends on the parameters used in it, and therefore, it is imperative to compare against exact results to establish both the validity and power of the technique. As a result, we can understand the convergence of the MRS as we tune the number of Fourier modes, the number of Monte Carlo samples, and the numerical accuracy in the numerical integration of Eq. ( <ref type="formula">21</ref>). Therefore, the present work serves as a benchmark for the MRS and a stepping stone for the MRS to be further applied to compute properties of models where such comparisons are no longer possible.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>C. Method of random surfaces: Vertex operator expectation values</head><p>The MRS method can also be used to compute observables that can be measured experimentally, such as expectation values and correlation functions of local operators. One family of such observables is the finite-temperature expectation values of the vertex operators V&#945; (r) defined in Eq. <ref type="bibr">(10)</ref>. From the theoretical perspective, these are truly nonperturbative objects in an interacting quantum field theory, so they are of fundamental interest. From a more pragmatic standpoint, they are related to physical observables in the various realizations of the sine-Gordon model. For example, in the tunnel-coupled cold atomic condensates system, this expectation value is called the coherence factor and is routinely measured via matter-wave interferometry <ref type="bibr">[42,</ref><ref type="bibr">43]</ref>. Comparing the measured values with theoretical results provides important information about the effective temperature so it can be used for thermometry <ref type="bibr">[44]</ref>.</p><p>For generic exponent &#945;, the expectation value can be written as</p><p>Expanding the numerator results in</p><p>Due to the neutrality condition, the expectation value in the integral is only nonvanishing when &#945; is an integer multiple of &#946;, i.e., &#945; = s&#946; with s &#8712; Z. The nonvanishing neutral terms are those where out of the n terms, we pick the + exponentials from n + = (ns)/2 and theexponentials from n -= (n + s)/2, so the parity of s and n must be the same. This gives</p><p>where the primed sum runs over integers n that have the same parity as s, n + + n -= n, n +n -= s, and</p><p>dr j e -n j&lt;k j k G(r j ,r k ) e -s n j j G(r,r j ) , <ref type="bibr">(27)</ref> where &#945; = &#945; 2 /(4&#960; ) and now n + of the signs are j = +1 and n -of them is j = -1. We proceed as above and decouple the integration variables using the Hubbard-Stratonovich transformation. It turns out to be convenient to perform the decoupling in the G(r, r j ) terms as well, which leads to</p><p>where</p><p>with s = &#945;/&#946; and C s = exp s 2 /2 m,n A mn = C s 2 . Plugging everything back into Eq. ( <ref type="formula">26</ref>), the infinite series can again be summed up, yielding</p><p>Compared to the partition function, the new ingredient is h &#945; ({t f }, r), the only term depending on the position r. In the ideal L &#8594; &#8734; case, the final result should be independent of the position r, but in a numerical calculation, there are deviations near the integration boundaries x = &#177;L/2. As shown in Appendix D 2, these deviations decay exponentially fast with the distance from the boundaries. It is easy to see from Eq. ( <ref type="formula">30</ref>) that V-&#945; (r) = V&#945; (r) * , as it should be. However, the result for V&#945; (r) is real, so</p><p>As a nontrivial consistency check, consider the case &#945; = &#946; and integrate both sides of Eq. ( <ref type="formula">30</ref>) over r. Using</p><p>where in the last step we used I 0 (x) = I 1 (x). Since from Eqs. ( <ref type="formula">9</ref>), <ref type="bibr">(10)</ref>, and (31) &#955; 0 cos &#946; &#966;(r) = &#955; V&#946; (r) , the above result is nothing but the Hellmann-Feynman theorem <ref type="bibr">[45,</ref><ref type="bibr">46]</ref>, demonstrating the validity of the MRS equations.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>III. RESULTS</head><p>This section presents our numerical results obtained with the MRS approach discussed in the previous section. To provide evidence for the strong claim about the nonperturbative nature of the approach and to justify the somewhat dubious steps, such as swapping integrals and summation in the perturbative expansion during the derivations in Secs. II B and II C, we compare our predictions with exact nonperturbative results obtained by exploiting the integrability of the sine-Gordon model via the NLIE approach, as well as with numerical results from the MSTHA.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>A. Free energy</head><p>Let us start with the free-energy density, which in our approach is given by</p><p>where the last term is the free-energy density of the free massless boson at inverse temperature R = 1/T . The exact free energy can be computed from the NLIE, which we review in Appendix B. It can also be determined using the thermodynamic Bethe ansatz, which was derived recently for general (rational) values of the coupling parameter <ref type="bibr">[15]</ref> </p><p>The exact relation between the (renormalized) perturbing parameter &#955; and the soliton mass is <ref type="bibr">[47]</ref> </p><p>where</p><p>We then obtain the free-energy density in the form</p><p>where f (MR) is the free-energy density calculated from the NLIE/TBA and B(&#955;) is the so-called bulk-energy density</p><p>which accounts for the different normalization of the free energy. Alternatively, f /T can also be interpreted as the (zerotemperature) ground-state energy in a finite volume R = 1/T . In the MRS calculation, we can conveniently tune the dimensionless inverse temperature MR through M by changing &#955;. The resulting free-energy density is plotted in Fig. <ref type="figure">1</ref> for four different L/R ratios and the NLIE result (dashed line). The deviation is a finite-size effect which decreases with increasing L/R. This is shown in the inset, where the difference between our results and the NLIE value is shown as a function of L/R on a log-log scale at three different MR values. The symbols in the inset represent three different mode cutoffs, demonstrating that the results have already converged with respect to m max . The dashed lines are linear fits, showing that the deviations follow a power law with an exponent close to 1.</p><p>Free-energy density at &#946; 2 = 8&#960;/25 ( = 2/25) as a function of the dimensionless inverse temperature MR for four different cylinder ratios L/R = 6, 8, 10, 12 for a fixed maximal mode m max = 20. The inset shows the difference between the NLIE and MRS results for three different inverse temperatures MR = 4, 8, 12 and for different mode cutoffs (m max = 20 in red triangle, 25 in blue crosses, 30 in black plus signs). We find the results are smoothly approaching the thermodynamic limit and have already converged as a function of the number of Fourier modes at m max = 20. The MRS results were obtained using &#8764;10 5 Monte Carlo samples, the statistical errors are smaller than the symbol size.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>B. Vertex operator expectation values</head><p>Despite the integrability of the sine-Gordon model, there is no analytic result for its finite-temperature expectation values, only for ratios of expectation values of descendant and primary vertex operators <ref type="bibr">[48]</ref>. The only exact result is for zero temperature <ref type="bibr">[19]</ref>, which we quote in Appendix C. At nonzero temperature, we benchmark the numerical results of the MRS approach with results obtained using the recently developed MSTHA method <ref type="bibr">[37]</ref>. Exploiting the Lorentz invariance of the sine-Gordon theory, expectation values calculated at zero temperature on a ring of circumference R are equal to finitetemperature expectation values at inverse temperature R in infinite volume. The accuracy of the truncated Hamiltonian approach can be well controlled by analyzing the convergence of the results with respect to increasing the energy cutoff. In the following, the truncation parameter is chosen so that the results are accurate up to (at least) four decimal places.</p><p>We show the results of the two methods for the one-point functions V&#946; (r) and V2&#946; (r) , together with the analytic T = 0 result, in Fig. <ref type="figure">2</ref> for &#946; 2 = 8&#960;/25 ( = 2/25) and in Fig. <ref type="figure">3</ref> for &#946; 2 = 8&#960;/15 ( = 2/15). To track the finite-size effects in the MRS method, we used different L/R ratios and calculated the expectation value at the origin (at the center of the cylinder). At small and medium inverse temperatures, the deviations are due to finite-size effects, and by increasing the L/R ratio, the agreement converges as a power law in L/R. This is also demonstrated in the insets of both figures, where we plot the difference between the MSTHA and the MRS results. In the inset of Fig. <ref type="figure">2</ref>, these differences are shown for two different MR values and for three different mode cutoffs (see the caption for details), demonstrating that the results are nicely converged in the Fourier mode cutoff. For large MR (low temperatures), we find deviations at larger L/R due to an insufficient number of Monte Carlo samples that become more severe as we increase L/R. This is natural as larger L/R requires more modes, and more modes require more averaging. We examine this further in Fig. <ref type="figure">3</ref>, where the convergence in L/R is no longer monotonic at sufficiently low temperatures. In the insets of Fig. <ref type="figure">3</ref>, we examine the T L/R L/R L/R L/R L/R L/R FIG. 3. Expectation value of the V&#946; (r) (left) and V2&#946; (r) (right) vertex operator for &#946; 2 = 8&#960;/15 ( = 2/15). Main panels: expectation values as functions of the dimensionless temperature MR. The MRS results obtained with mode cutoff m max = 20 for different L/R ratios are shown in symbols, while the dashed line represents the MSTHA data. The horizontal solid blue lines show the exact T = 0 result. In the MRS N s = 5&#215;10 5 Monte Carlo samples were used, the error bars show the standard deviation of the sampling, increasing whenever MR or the ratio L/R becomes larger. In particular, the observed deviation at low temperatures (large MR) signals the need for more MC samples to get accurate results. Insets: Difference between the MRS and the MSTHA results at fixed values of MR = 3 for different numbers of Monte Carlo samples (magenta diamond and cyan diamond N s = 3&#215;10 5 , green down-facing triangle and up-facing triangle 4&#215;10 5 , black cross and black plus 5&#215;10 5 , for MR = 2.23 and MR = 4.46, respectively).</p><p>dependence on the number of Monte Carlo samples N s of the results, finding that our estimates are nicely converged at small L/R ratios, but this increases the number of required Monte Carlo samples. We show some additional numerical results for other values of &#946; in Appendix E. At small &#946; values much larger L/R ratios are necessary to avoid finite-size effects, while for larger couplings more Monte Carlo samples are necessary to obtain accurate results.</p><p>The expectation value of V&#946; (r) can also be computed from the free energy via the Hellmann-Feynman theorem <ref type="bibr">[45,</ref><ref type="bibr">46]</ref>, making it accessible to the NLIE approach as well; we note that the result obtained from the NLIE agrees with the MSTHA to a high accuracy. This is not the case for the other vertex operators for which the only viable method is the MSTHA, which we use to verify the MRS results in the right panels of Figs. <ref type="figure">2</ref> and <ref type="figure">3</ref>.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>IV. CONCLUSIONS AND OUTLOOK</head><p>We generalized the method of random surfaces to compute one-point functions of exponential operators in the sine-Gordon quantum field theory at finite temperature. We demonstrated the simplicity and effectiveness of the method of random surfaces and validated its convergence by comparing the free energy and the one-point functions to the results of alternative nonperturbative approaches, such as the exact NLIE relying on integrability and the numerical MSTHA truncated Hamiltonian method, finding excellent agreement.</p><p>This work constitutes the first step in developing the method. The next step is the computation of finite-temperature multipoint correlation functions, which will be reported soon in another publication <ref type="bibr">[38]</ref>, where such comparisons are no longer possible. The method is clearly versatile enough to treat additional refinements of the model, such as spatially inhomogeneous couplings and the inclusion of a disorder (notably being able to continue to work when direct quantum Monte Carlo treatments suffer from a sign problem). It is capable of treating finite systems with boundary conditions. Additionally, the method of random surfaces has the potential to be extended to time evolution after a quantum quench, providing a promising alternative to the existing semiclassical approximations <ref type="bibr">[49]</ref><ref type="bibr">[50]</ref><ref type="bibr">[51]</ref><ref type="bibr">[52]</ref><ref type="bibr">[53]</ref> and Hamiltonian truncation methods <ref type="bibr">[30,</ref><ref type="bibr">37,</ref><ref type="bibr">53]</ref>. These features make the MRS an invaluable tool for the description of experimental realizations of sine-Gordon quantum field theory <ref type="bibr">[5]</ref><ref type="bibr">[6]</ref><ref type="bibr">[7]</ref><ref type="bibr">[8]</ref><ref type="bibr">[9]</ref><ref type="bibr">[10]</ref> and we expect that its applicability will be even further reaching than previously anticipated.</p><p>NSF Career Grant No. DMR-1941569. This work was performed in part at the Aspen Center for Physics, which is supported by National Science Foundation Grant No. PHY-2210452 (J.H.P.).</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>APPENDIX A: FOURIER MODES</head><p>In this Appendix, we analyze the Fourier coefficients appearing in Eqs. <ref type="bibr">(16)</ref> in the main text. Let us recall that we consider the translational invariant Green's function of the infinite cylinder</p><p>It is periodic in the &#964; direction with period R, and we periodically extend it in the x direction with period 2L. Therefore, it can be expanded in a double Fourier series:</p><p>B mn cos(&#960;mx/L) cos(2&#960;n&#964; /R).</p><p>(A2)</p><p>We observe that the Fourier coefficients B mn are numerically very close to the simple expressions</p><p>and</p><p>This is shown for L/R = 10 in Fig. <ref type="figure">4</ref> for the coefficients A 0m and A 1m and in Fig. <ref type="figure">5</ref> for A m0 and A m1 . Increasing the accuracy of the numerical integration decreases the difference between the analytical and numerical values, showing that the analytical expressions are indeed correct.</p><p>To understand these coefficients better, we take the Laplacian of the function &#7712;(x, &#964; ) defined by the analytic coefficients. The (n, m) dependence of the coefficients cancels with the factors coming from the derivatives of the cosine functions, and we are left with</p><p>cos(&#960; mx/L) cos(2&#960; n&#964; /R) .</p><p>(A4)</p><p>Applying the Poisson summation formula, this is equal to</p><p>so it is the Green's function of the Laplace equation with periodic boundary condition in the &#964; direction. As L &#8594; &#8734;, it becomes the Green's function on the infinite cylinder. However, this does not specify the function completely, since it also depends on the boundary conditions. We now provide some analytic evidence supporting that &#7712;(x, &#964; ) goes m m FIG. 5. Numerical and analytical values for the Fourier modes A m0 and A m1 for L/R = 10.</p><p>to G(r, 0) for L R. First, consider our function at x = 0,</p><p>For L R, the coefficients are exponentially close to 2/n, and we obtain</p><p>for r = (0, &#964; ). As a second check, let us evaluate the function at x = &#177;L:</p><p>For L R, the coefficients in the sum are exponentially suppressed and we are left with the first two terms. It is easy to check that for r = (&#177;L, &#964; )</p><p>up to exponential precision in 2&#960;L/R, so &#7712;(&#177;L, &#964; ) rapidly approaches G(r, 0). These analytic derivations hold for infinitely many modes. In the numerical calculations, the double Fourier series is truncated at a mode number m max , which causes deviations from the exact Green's function even at large-enough L/R. We show this truncation effect in Fig. <ref type="figure">6</ref> where we compare &#964; = 0 and x = 0 sections of the Green's function with the truncated Fourier series approximation. Naturally, there are significant deviations near the divergence at the origin. Nevertheless, the fact that it is an integrable singularity makes it possible to obtain accurate results. </p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>APPENDIX D: ADDITIONAL NUMERICAL RESULTS</head></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="1.">Numerical integration via the trapezoidal rule</head><p>As discussed in the main text, we compute the integrals</p><p>by discretizing the mode functions &#968; f (u) calculating the Riemann sum. The discretization must be fine enough to get a sufficiently accurate result. In Fig. <ref type="figure">7</ref>, we plot the value of the integral for three different mode cutoffs m max versus the ratio of the number of discretization points n d (used both in the x and &#964; directions) and the mode cutoff m max . A similar accuracy is achieved for the same ratio in the different cases, showing that the necessary resolution of the grid (i.e., n d ) for a desired accuracy is proportional to m max . In our simulations, we used n d /m max between 2.5 and 4. </p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="2.">Position dependence of the one-point function</head><p>Since we integrate the multipoint functions over a finite space-time domain, our numerical results for the vertex operator expectation values are not translationally invariant in the x direction. We plot the x dependence of the expectation value V &#946; (x, 0 in Fig. <ref type="figure">8</ref>. The inset shows that the deviation decays exponentially with the distance from the edges. The associated correlation length depends on the dimensionless inverse temperature MR.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>APPENDIX E: ADDITIONAL NUMERICAL RESULTS FOR SMALLER AND LARGER COUPLINGS</head><p>In this Appendix, we present some additional results for vertex operator expectation values for smaller and larger values of &#946;.</p><p>In Fig. <ref type="figure">9</ref> we compare our results obtained for the small coupling &#946; 2 = 4&#960;/125 ( = 1/125) with two other approaches. T L/R L/R FIG. 9. Expectation value of V&#946; (r) (left) and V2&#946; (r) (right) vertex operator for &#946; 2 = 4&#960;/125 ( = 1/125) as functions of the dimensionless temperature MR. The MRS results, shown in symbols, were obtained with mode cutoff m max = 20 for two different L/R ratios and using 5&#215;10 5 Monte Carlo samples. The MSTHA data is shown in dashed line, the horizontal solid blue line shows the exact T = 0 result, and the red crosses represent the results of Ref. [56] obtained in the classical sine-Gordon model.</p><p>As can be seen from the plot, this &#946; requires a large L/R ratio due to the large correlation length of the conformal correlator. As discussed in the main text, the condition for V&#945; is L/R 1/&#945; 2 , which for &#945; = &#946; yields L/R 10. We compare our results with the MSTHA data and, for &#945; = &#946;, with the result of Ref. <ref type="bibr">[56]</ref> obtained in the classical (&#946; &#8594; 0) limit of the sine-Gordon model. We find excellent agreement between the three methods.</p><p>We present our results for the larger coupling &#946; 2 = 6&#960;/5 ( = 3/10) in Fig. <ref type="figure">10</ref>. Here a larger number of Monte Carlo samples are necessary to achieve an agreement with the MSTHA data similar to the &#946; values shown in the main text. In the case of the vertex operator V&#946; (r), 3&#215;10 7 samples were required to get a nice agreement. For V2&#946; (r), even more samples would be necessary for accurate results at larger values of MR. </p></div></body>
		</text>
</TEI>
