Simulating lattice field theories on multiple thimbles

Simulating thimble regularization of lattice field theory can be tricky when more than one thimble is to be taken into account. A couple of years ago we proposed a solution for this problem. More recently this solution proved to be effective in the case of 0+1 dimensional QCD. A few lessons we can learnt, including the role of symmetries and general hints on algorithmic solutions.


Introduction
Thimble regularization of lattice field theories was put forward as a possible solution to the sign problem [1,2]. The solution is very elegant and in principle it is a fundamental one: the domain of integration of an integral featuring a (possibly widely) oscillating phase is deformed in such a way that the latter is turned into a stationary phase. There are actually (at least) a couple of caveat attached to the former statement. First of all, by changing the domain of integration one is left with a new phase [3] (aka the residual phase) which roughly speaking comes as a consequence of having changed integration variables. This is not a really serious problem, or at least for every problem which has been worked out till now in this framework it was shown that this phase can be effectively taken into account by reweighting. A second matter of concern is more serious: an integral is in general turned into a sum of integrals. One speaks of thimble decomposition as the original domain of integration is turned into the sum of many thimbles. While there is an argument suggesting that in the thermodynamic limit one single (dominant) thimble provides the relevant contribution one is interested in [1], in recent times examples were provided showing that in finite systems collecting the contributions of different thimbles can be a tricky business [4][5][6] and this was the main motivation for the modified approach known as the holomorphic flow [4]. The thimble regularization of QCD in 0 + 1 dimensions was presented at the last year Lattice conference [7]. Here we improve on that work, showing a better way to take into account the three contributions which are expected in the thimble decomposition of the problem at hand. Results are better than the previous ones due to two improvements: first of all, a symmetry argument can reduce the number of contributions that we have to sum to solve the theory (it turns out that there are essentially two distinct contributions); in the second place, we make use of a better Monte Carlo strategy. These are indeed valuable improvements. Together with a semiclassical argument pointing out in which regions of the parameter space it is essential to collect the results from all the three thimbles, they lead us to the full solution of the problem. All in all, QCD in 0 + 1 dimensions provides a nice example of how a theory can be simulated on multiple thimbles. e-mail: francesco.direnzo@unipr.it 2 Thimble regularization in a nutshell

Thimble decomposition
A sign problem is in place when one has to solve a theory, i.e. one has to compute the relevant functional integrals (we make use of a one-dimensional notation: results are to be thought of in and the action is complex: Thimble regularization is built out of three steps: • We complexify the degrees of freedom: x → z = x + iy and as a consequence we consider S (z) = S R (x, y) + S I (x, y).
• We look for (stationary) critical points satisfying ∂ z S = 0.
• To each critical point p σ a thimble J σ is attached and the original function integral is now decomposed into thimble contributions the sum formally extending to all thimbles, even though some n σ can be zero (thus, not all the critical points do contribute).
We still have to define what a thimble is, but we can immediately read interesting features. A positive measure e −S R is in place and constant phases e −i S I (p σ ) have been factored out of the integrals: the imaginary part of the action stays constant on each thimble. A so-called residual phase e iω is there: it accounts for the relative orientation between the canonical complex volume form and the real volume form, characterizing the tangent space of the thimble. So, what is a thimble? It is the union of all steepest ascent (SA) paths attached to a given critical point (we now explicitly write in terms of the components of z) It is easy to see that in this notation the critical point is associated to t = −∞.
A crucial point is that a thimble has the same real dimension of the original domain of integration (on which the functional integral was defined in the first instance). This can be explicitly checked by solving the Takagi problem for the Hessian of the action computed at the critical point Takagi vectors v (i) σ provides a basis for the tangent space at the critical point, while Takagi values λ (σ) i fix the rate at which the real part of action increases along the steepest ascent paths while they leave the critical point. This is not the end of the story, since we have only found out the tangent space at the critical point. In order to get a basis at a generic point on the thimble we have to parallel-transport the basis for the tangent space at the critical point along the flow. Now, any flow (i.e. any SA path) leaves the critical point along one possible direction on the tangent space. If we impose a normalization condition n i=1 n 2 i = R all those directions are mapped to vectors n i=1 n i v (i) . It is thus quite natural to single out any given point on a thimble by the correspondence with S n−1 R the (n − 1)-sphere of radius √ R and where we denote by t the time coordinate parametrizing the flow along the SA path. In order to have a basis V (i) σ (n, t) for the tangent space at the (generic) point associated to directionn and flow time t one has to solve the (associated) flow equations

A crude Monte Carlo on thimbles
By changing variables of integration one can now rephrase the thimble decomposition as where we have defined The formula has to be understood in the following way. The λ (σ)  (3). One then assembles the V (i) σ into the matrix V σ , the modulus of whose determinant enters the definition of S eff . At the same time, the phase of det V σ (n, t) provides the residual phase e iω(n,t) . The notation σ shows that (2) holds for each critical point and at each critical point one has to solve a different Takagi problem, resulting in different ∆ (σ) n (t).
We can now devise a simple crude Monte Carlo scheme for simulating on thimbles: • We pick up randomly (with flat distribution) a directionn.
• Since we want to compute the contribution coming from the SA leaving the critical point p σ alonĝ n, we prepare convenient initial conditions both for the field and for the tangent space basis vectors for such a SA. We can do this, since near the critical point solutions of the flow equations are know as 1 i t which we can compute for t = t 0 0.
• We then integrate the SA equations for the field and the equations for transporting the basis vectors all the way up till we reconstruct the (dt) integrals appearing in (4) and while ascending we compute both the integral in the numerator and the one in the denominator.
3 Thimble regularization for 0 + 1 dimensional QCD 3.1 QCD in 0 + 1 dimensions QCD in 0 + 1 dimensions is given in terms of staggered fermions on a one-dimensional lattice with (even) N t sites in the temporal direction (the temperature being given by aN t = 1/T , where a is the lattice spacing). A genuine sign problem is there as in real QCD, due to the presence of a (quark) chemical potential. The partition function of the theory for N f degenerate quark flavours of mass m is where D is the lattice staggered Dirac operator andδ ii is the anti-periodic Kronecker delta. In a convenient gauge one actually has while A = 2 cosh(µ c /T ) and µ c = sinh −1 (m) (from now on, we set a = 1). We computed three main observables. The chiral condensate is the first one The other two are the Polyakov loop Tr U and the anti-Polyakov loop Tr U † = Tr U µ→−µ . The latter two can be related to the quark number density n ≡ T ∂ ∂µ log Z by a relation which takes quite different forms for different values of N f . There are known numbers to compare to, since analytical results for 0 + 1 QCD are available (see for example [9]).

0 + 1 QCD in thimble regularization
To solve in thimble regularization we need two ingredients: • We need to complexify the degrees of freedom, i.e.

SU (N) U = e ix a T a → e iz a T a = e i(x a +iy a )T a ∈ SL (N, C) .
Notice that SU (N) U † = e −ix a T a → e −iz a T a = e −i(x a +iy a )T a = U −1 ∈ SL (N, C) .
The equations for the Hessian and the flow equations for the vectors can also be easily written down 2 . Notice that we find three critical points {U k = e 2πik/3 1} with k = 0, 1, 2 and they all contribute to the thimble decomposition.
At last year conference we presented results [7], which were obtained using the crude Monte Carlo scheme presented above. We took into account all the three critical point, even if one can notice that there are regions of the parameter space where the thimble attached to the identity captures virtually the complete result. This agrees with semiclassical estimates. The complete partition function is built out of the three contributions attached to the critical points Z = Z 0 + Z 1 + Z 2 and one can define and compute (in the semiclassical approximation)

Reflection symmetry
Simmetries always play an important role: not surprisingly, we have a symmetry at hand which lets us check our results and which also makes our life easier [11]. The QCD action (in our case, we only have the Dirac determinant) is invariant under charge conjugation C defined by Together with the generalization of γ 5 -hermiticity at finite chemical potential We thus expect thimbles to appear in conjugate pairs. Consider the {U k }. U 0 = 1 is real and therefore self-conjugate (i.e. computations on the associated thimble yield real results). Being e 4πi/3 = e −2πi/3 , we immediately see that U 2 = U 1 : U 1 and U 2 are a conjugate pair of critical points and results of integration on U 2 should be the complex conjugate of those on U 1 , yielding an overall real contribution to the partition function and to the expectation value of observables as well. All this is well evident in numerical simulations and in all the numerical results that we present in the following we take advantage of this symmetry: the real part of results from thimble 1 and 2 have been averaged.

Importance sampling
If a single thimble contributed 3 , the computation of (1) would simply amount to 4 In terms of the representation (2), we can now rephrase where Now, fn looks like a functional integral along a single complete flow line. On the other side, (7) is nothing but the average of the fn (i.e. the average of the contributions that a given observable takes from complete flow lines) and the weight Zn/Z σ represents the fraction of the partition function which is provided by a single complete flow line. Zn/Z σ provides a natural setting for importance sampling: directionsn have to be extracted according to the probability P(n) = Zn/Z σ . We thus proceed as follow. Sitting on the current configuration (associated to a directionn), we propose a new one (associated to a directionn ).n is identical ton apart from two randomly chosen components, say (n i , n j ) with i j. Given the normalization n = √ R and the values of all {n k } k i, j we define C as and we observe that there is a coordinate system in which the current values of (n i , n j ) are parametrized as . Now we change φ → φ by extracting φ − φ flat in a given (tunable) range. This results in (n i , n j ) → (n i , n j ), while for all the other components (k i, j) n k = n k . We finally accept the proposed configuration with the standard Metropolis test P acc n n = min 1, Zn Zn .
In our case three contributions should be in principle taken into account. Actually, due to the symmetry of Section 4.1, only two distinct contributions are in place and we have to compute n 0 e −i S I0 Z 0 e i ω 0 + n 12 e −i S I12 Z 12 e i ω 12 (9) in which subscript notations should be evident. In [8] we observed that equation (9) can be rewritten where we defined α ≡ n 12 e −i S I12 Z 12 n 0 e −i S I0 Z 0 The idea is now to take a given observable as a normalization point, thus determining the value of α.
All the other observables of the theory can then be computed using this input. This could indeed be successfully done. Figure 2 confirms the effectiveness of the procedure. For N f = 1, m = 0.1 the point µ/T = 2.0 was completely out of reach for flat, crude Monte Carlo, while results are successfully computed with the improved method. At these values of parameters the tiny value of the chiral condensate (of order 10 −2 ) results from a delicate cancelation of the contributions coming from the different thimbles. This has nothing to do with the original sign problem (nor e.g. with the residual phase), but it is simply a numerical accident occurring for a given observable at a given value of parameters.

Conclusions
There is still quite a long way to go before thimble regularization can tackle the real goal, i.e. QCD at finite density. This work nevertheless shows that simulations on multiple thimbles are viable.