The phase diagram of the $U(2) \times U(2)$ Sigma Model and its Implications for Chiral Hierarchies

Motivated by the issue of whether it is possible to construct phenomenologically viable models where the electroweak symmetry breaking is triggered by new physics at a scale $\Lambda \gg 4\pi v$, where $v$ is the order parameter of the transition ($v\sim 250$ GeV) and $\Lambda$ is the scale of new physics, we have studied the phase diagram of the $U(2) \times U(2)$ model. This is the relevant low energy effective theory for a class of models which will be discussed below. We find that the phase transition in these models is first order in most of parameter space. The order parameter can not be made much smaller than the cut-off and, consequently a large hierarchy does not appear sustainable. In the relatively small region in the space of parameters where the phase transition is very weakly first order or second order the model effectively reduces to the O(8) theory for which the triviality considerations should apply.


Introduction
In the minimal version of the Standard Model of electroweak interactions the same mechanism (a one-doublet complex scalar field) gives masses simultaneously to the W and Z gauge bosons and to the fermionic matter fields (other than the neutrino). This remains so in many extensions of the minimal Standard Model, such as those consisting of more scalar doublet fields or even fields in other representations of SU (2) L .
On the contrary, the mechanism that gives masses to the W and Z bosons and to the matter fields remains somewhat distinct in models of dynamical symmetry breaking (such as technicolor theories [1]). In these models, there are interactions that become strong, typically at the scale 4πv (v = 250 GeV), breaking the global SU (2) L × SU (2) R symmetry to its diagonal subgroup SU (2) V and producing Goldstone bosons which eventually become the longitudinal degrees of freedom of the W ± and Z. Yet, to transmit this symmetry breaking to the ordinary matter fields one usually requires additional interactions, characterized by a different scale M . Generally, it is assumed that M ≫ 4πv. It seems then natural to ask whether it is necessary to have these two very different scales and whether it would not have been possible to arrange things in such a way that the SU (2) L × SU (2) R → SU (2) symmetry breaking takes place at a scale Λ, where 4πv ≪ Λ. Although the scale where the symmetry breaking takes place, Λ, and the scale characterizing the new physics, M , need not be exactly the same, we shall assume that Λ ≃ M and refer only to Λ hereafter.
Two phenomenologically viable paradigms of the above possibility are the strong-coupling extended technicolor (ETC) models [2], and the top-condensate (TopC) models [3], in which the underlying dynamics is, typically, a spontaneously broken gauge theory, characterized by a scale Λ, with Λ ≫ 1 TeV. At sufficiently low energies, the dynamics can be modeled by four-fermi interactions (of either, technifermion doublets in ETC models, or the quarks of the third family in TopC models) which are attractive in the scalar channel. Then, it appears possible that, by tuning the corresponding coupling sufficiently close, but above a critical value, chiral symmetry breaking occurs, but the condensate itself is of the order of the weak scale v, much smaller than its natural value O(Λ). It has been pointed out in ref. [4], that a necessary condition for this to happen is that the low-energy effective theory, which retains the light degrees of freedom below the scale Λ, where chiral symmetry breaking occurs, possesses itself a second order phase transition. It is only then that there can consistently be a hierarchy between the order parameter v and the scale Λ.
If the strongly interacting fermions at scale Λ are electroweak doublets, then the chiral symmetry is U (2) L × U (2) R and the relevant lagrangian at that scale consists in a bunch of four-fermion interactions. The precise form of these four-fermi interactions will not concern us here. It has been argued, using analytical methods [5] as well as lattice simulations [6,7], that four-fermion models are, at low energies, equivalent to an effective theory consisting in a scalar-fermion model with the appropriate symmetry, which in our case it will be U (2) L × U (2) R .
Using an effective theory description also frees us from having to appeal to any particular model. Thus, the analysis may remain valid beyond the models we have just used to motivate the problem. In this effective theory we need to retain only the particles that remain light after chiral symmetry breaking. Thus, it must necessarily contain the Goldstone bosons emerging from the breaking of the global U (2) L × U (2) R symmetry. There may also be some light (compared to the scale Λ) scalars. If present, the U (2) L × U (2) R symmetry can be linearly realized. If not, the symmetry should be realized non-linearly. Of course the presence or not of such additional scalars depends only on the underlying additional sector (or, equivalently, on the four-fermion interactions it leads to). However, the linear and non-linear theories differ, for sufficiently large values of the scalar masses, by terms of O(µ 2 /16π 2 v(µ) 2 ). In the coming paragraphs we just argue that v is generally large. Hence, at low energies these terms are small and thus using a linear realization is really no restriction at all.
We shall thus assume that the low-energy theory is a general linear sigma model with U (N ) L × U (N ) R (which we later take N = 2) symmetry, whose effective action, for general N , is given by where φ(x) is a complex N × N matrix, the order parameter of the high-energy phase transition. The action (1) is invariant under the global symmetry transformation φ → LφR † , where L, R are U (N ) matrices. The electroweak interactions are obtained by gauging an SU (2) × U (1) subgroup of this global symmetry, which after the U (N ) L × U (N ) R symmetry breaking gives masses to the W ± and Z. In some models additional fermions remain in the spectrum below Λ. We have not considered them and thus our results do not apply to such models.
This model, as will be discussed in detail below, possesses for λ 2 = 0 a first order phase transition whose strength varies in the (λ 1 , λ 2 ) space. As one adjusts m 2 past a critical value (m 2 = 0 in mean field theory), the vacuum expectation value v jumps discontinuously from zero in the unbroken phase to some finite nonzero value in the broken phase. If the couplings (λ 1 , λ 2 ), which are obtained by matching with the underlying strong dynamics at the cutoff Λ, belong to a region where the phase transition is weakly first order, then the vacuum expectation value (v.e.v.) v can be small, v/Λ ≪ 1, and the U (2) L × U (2) R model can still be a valid description of the low-energy dynamics. If, on the other hand, the couplings (λ 1 , λ 2 ) belong to a region where the transition is strongly first order, such a hierarchy is not possible and models leading to such (λ 1 , λ 2 ) values should be excluded.
It is our purpose in this paper to investigate in detail the phase structure of the above model, and to conclude whether large chiral hierarchies are sustainable in this context. A number of authors have previously considered this possibility. In [4] such a model was considered as a possible parametrization of the low-energy physics of top-condensate models or strongly coupled extended technicolor models. The authors concluded that, within the framework of a perturbative one-loop analysis, a large hierarchy is unlikely unless λ 2 is close to zero. Later, in [8], it was pointed out that a leading order 1/N c analysis combined with two-loop beta functions can change the conclusions, and that consistent large hierarchies were not disallowed. Unfortunately, all these analysis rely on perturbation theory which is unreliable at strong coupling. To settle the issue we have performed Monte Carlo simulations in terms of the lattice regularized version of the action (1) above. A preliminary investigation along these lines was undertaken in ref. [9].
We confirm the first-order character of the transition for λ 2 = 0. We have obtained a detailed picture of the behaviour of the order parameter for nonperturbative values of the sigma model couplings and semi-quantitative estimates of the correlation length. Where meaningful we compare our results to those obtained via the effective potential. We have found that a large hierarchy is untenable in most of parameter space. v is typically several orders of magnitude too large.
To get the above results we have performed high-statistics runs using a hybrid algorithm (with and without Fourier acceleration) which, to our knowledge, had not been used before for this type of systems. We have also written a traditional Metropolis Monte Carlo code for comparison. Therefore we believe that our results are also of some interest to the lattice expert.
In the case N = 2 the action (1) depends on eight degrees of freedom and the field φ can be conveniently parametrized by where τ a are the Pauli matrices for a = 1, 2, 3, and the identity matrix for a = 0. The action (1) is invariant under the rigid symmetry, φ → LφR † , where L, R ∈ U (2). If we set λ 2 = 0, then the symmetry is enhanced to The pattern of symmetry breaking depends on the sign of the coupling λ 2 . If λ 2 > 0 the vacuum expectation value (v.e.v.) can be rotated to < σ 0 >= v, and the breaking occurs according to with v 2 = −m 2 /(4λ 1 + 2λ 2 ). The π α are then the Goldstone bosons while the masses for the other states are If λ 2 < 0, the symmetry breaks along the τ 0 + τ 3 direction, with v 2 = −m 2 /2(λ 1 + λ 2 ), and the masses are In this case though, the mean field solution is not a real minimum but a saddle point.
Although these are tree-level relations, they are preserved by quantum corrections in terms of renormalized parameters appropriately defined. In the physically interesting region, and barring any unexpected non-trivial fixed point, four dimensional scalar theories such as in (1) are believed to be infrared free and have Landau poles in the ultraviolet. Therefore if we wish to have large masses for all these scalar resonances we must set at least λ 2 (Λ) ≫ 1.
(Recall that the bare coupling and the renormalized coupling at the cut-off scale are the same thing.). Let us keep in mind, though, that even after giving these scalars a large mass, some non-decoupling effects remain, as pertain to a spontaneously broken theory.
The electroweak interactions are included by identifying SU (2) L with SU (2) W and the should remain light due to their small Yukawa couplings, y and, will be ignored henceforth.
Similarly, we have neglected the effects of gauge couplings since they are small at scale Λ.
We have also ignored the possible presence of custodially non-invariant interactions.
As discussed, a non-perturbative method, such as lattice techniques, is called for to determine the possibility of making v of the order of weak scale. After introducing the lattice regulator, one expects that Λ/v ∼ ξ β , with β the appropriate critical exponent, so a hierarchy will only be possible if the correlation length ξ is big, or what is the same, the transition is second order or weakly first order.
High-energy physicists have not devoted to first order phase transitions nearly as much interest as in the case of continuous ones since there is no natural way to define a continuum limit, i.e. to shrink the lattice spacing to zero, because the correlation length never becomes infinite. However, this is not a problem from the point of view of an effective field theory because our continuum theory has most definitely an ultraviolet cut-off Λ, above which it is no longer valid. The lattice cut-off can then be identified with this continuum ultraviolet cut-off, i.e. a = 2πΛ −1 . The relation between the lattice and the continuum cut-off can be unambiguously established, but since we will not work it out here all we can say is that the relation between the lattice cut-off and the physical scale Λ can be defined up to terms of only. The physical parameter controlling the size of the corrections is naturally the correlation length of the system. If the correlation length is relatively large (in lattice units), corrections will be small, non-universal cut-off effects controlable, and the results meaningful.
It turns out that in most of the interesting regions of parameter space, the transition is of first order, but with relatively large correlation lengths. We can thus deposit some confidence in our conclusions.

The Coleman-Weinberg mechanism
In this section, we analyze the phase transition within (renormalized) perturbation theory.
Although this applies strictly only at weak coupling, that will provide a qualitative feeling about the transition. Moreover, it will be suggestive of the regions in coupling constant space where we should perform Monte Carlo simulations and provide a qualitative understanding of our results.
The model described by (1) is very similar to the Ginzburg-Landau phenomenological model of continuous transitions. However, theoretical considerations [10] and numerical work show that, whenever λ 2 = 0, by tuning the parameter m 2 (m 2 < 0) the system undergoes a first-order transition, which particle physicists know as the Coleman-Weinberg mechanism [11]. Due to quantum fluctuations the system develops a vacuum expectation value at a finite value of the correlation length ξ.
The Coleman-Weinberg mechanism has been given a nice geometrical interpretation in a massless theory, due to Yamagishi [12], in terms of the β functions, and its associated fixed points. The general solution of the renormalization group equation (in a dimensionless regulator, such as dimensional regularization) for the effective potential V (ϕ) where ϕ = (Trφ † φ) 1/2 , when λ 2 > 0 is given by where t = ln(ϕ/µ). Then, the condition for the existence of a local extremum away from the Eq. (9) is referred to as the "stability line".
The corresponding RG equation for λ 2 < 0 is given by and the stability line is described by If there exists a certain value of t where the conditions (9) or (11) are satisfied, in a region where V ′′ > 0 and V < 0 at the minimum (see for the corresponding equations in terms of the appropriate β functions in ref. [9]) then ϕ = v = 0 and the transition is of first order.
The expressions for the β functions can be found at one-loop level in [13,4,9] and are plotted as solid lines in fig. (1). The stability line is indicated as a dotted line. Then, starting from some value (λ 1 , λ 2 ) at the scale Λ and following its RG trajectory, one flows in the infrared either towards the stability line or towards the infrared fixed point at the origin (if λ 2 = 0). If the RG trajectory crosses the stability line, then the transition must necessarily be of first order at that particular value of (λ 1 , λ 2 ) we started. Were it of second order, the correlation length would be divergent and it cannot possibly become finite again after a finite number of renormalization group blockings. For |λ 2 | small the couplings flow towards the region λ 1 , λ 2 ≪ 1 and even if they cross the stability, they do so after very many decades of running; the phase transition in this case is weakly first order, the more so as |λ 2 | → 0.
The corresponding β functions at two-loop level can be found in ref. [14,8], whose solutions are plotted (for zero Yukawa coupling) in fig. (1) with dashed lines. One finds out that the stability is improved by the two-loop corrections [8]. For a bare theory with λ 2 < 0 or one that is close to the stability line the flow is again to the left towards the stability line, however the flow is slower than at one-loop. Furthermore, there exists a region with λ 1 , λ 2 > 0 where the flow is reversed and it appears that it never crosses the stability line. However, this only hints upon the breakdown of perturbation theory and a nonperturbative analysis is called for.
Although strictly valid only with a dimensionless regulator, and hence definitively linked to perturbation theory, the conclusions of the above analysis are expected to remain approximately valid in the lattice regularization where all sorts of irrelevant operators appear in the effective potential. As long as the correlation length is large enough, the continuum physics can be used as a guidance. Checking to what extend these arguments are valid is one of the motivations of the present work.
It is also essential that no other fixed point unreachable in perturbation theory exists.
Should one be present, the RG trajectories would be distorted and there could be regions where the transition is second order. We found no evidence of such a fixed point. Even with only the gaussian fixed point it would still be conceivable that there might exist a region of non-zero measure whose RG trajectories end in the gaussian fixed point at the origin.
For these values the transition would be of second order. For small values of (λ 1 , λ 2 ) an effective potential calculation shows that this happens only if λ 2 = 0, so assuming that the RG trajectories follow a potential flow this possibility appears to be ruled out too.
The above picture suggests that, if the couplings of a given bare theory are located in the region limited by the stability line and the straight lines 2λ 1 (Λ) + λ 2 (Λ) = 0 for λ 2 > 0 and λ 1 (Λ) + λ 2 (Λ) = 0 for λ 2 < 0 (so that the potential is positive definite at large φ), In this case, there is no proper continuum limit and strictly speaking no field theory at all.
However, if the transition is sufficiently weakly first order one can speak of an approximate continuum limit, with lattice artifacts being still relatively small. In contrast to the λ 2 = 0 case though, the correlation length is not tunable (by m 2 ) but is rather determined by the quartic couplings (λ 1 , λ 2 ). Moreover, a large hierarchy, Λ/v is in general not tenable.

The phase transition on the lattice
In Table I which corresponds to the susceptibility, wherē The expectation value of the above operator is then proportional to v 2 in the broken phase and zero in the unbroken phase, modulo finite size corrections. We also measured the expectation value of O ′ = Tr(φ † φ) as an alternative order parameter.
We have used lattices of sizes ranging from L = 4 to L = 14. The exact procedure we have used depended somewhat on the region of (λ 1 , λ 2 ) in which the simulations were performed. In general, to obtain information about the order of the transition at each given (λ 1 , λ 2 ), we have searched for hysteresis effects in the measurement of the order parameter: we performed thermal cycles in the relevant parameter, m 2 , across the critical region where the field configuration from the last run was used as input for the next run. Strong hysteresis loops are an indication of a strong first order transition. On smaller lattices, we have also looked at the histogram distribution of Trφφ † , where a double-peak structure across the critical region is an indication for two coexistent minima. This procedure helps us identify the critical point of equivalent minima, as well as the range of m 2 where metastability is observable. Tunneling, of course becomes rare on larger lattices, which we have used to look for coexistence. The measurements we provide in Table I for the order parameter v 2 always refer to the largest lattice used.
Near a first order transition the effective potential develops two minima. One of the minima, say ϕ = 0, is the lowest, and as we increase the relevant parameter (−m 2 ) the nontrivial minimum becomes dominant, and the system acquires a v.e.v. ϕ = v = 0, and eventually the minimum at the origin disappears. This does not imply that as soon as one of minima becomes dominant the system jumps to it; near the transition there is a potential barrier between them whose height determines the strength of the transition and the tunneling rate. The relation tells us that the weaker the transition, the less steep the effective potential will be, and the more difficult it will be to observe metastable states. Of course, if the correlation length near the transition is bigger than, or of the order of, the lattice size itself, one should not expect to see metastable states because the system is not able to see the distinction between the However, for this relation to be true the residues of all particles have to be properly normalized. This is not necessarily so on the lattice and we are forced to distinguish between v phys , the physical value, and v latt , the value we obtain from our simulations. The relation between the two is where Z G is the Goldstone wave-function renormalization defined by the unit residue condition. The value of Z G has been estimated in [15] and found to be always close to one.
Although these results were derived at λ 2 = 0, we expect that the wave function renormalization Z G stays close to but smaller than one and thus taking it into account can only increase the value of v phys .

Weak coupling
In the weak coupling region (λ 1 , λ 2 < 1), lattice perturbation theory does apply and can be used to compare with the numerical data. At tree level (equivalent to assuming the mean field approximation) the transition is always of second order. The symmetry breaking pattern and the tree-level relations have been described in section 2. Radiative corrections change this behaviour. The bare one-loop effective potential was computed in [9]. If λ 2 > 0, the result is where p 2 = µ 2 − 2cos(p µ ) is the lattice propagator. The quantity m 2 ≡ m 2 − m 2 c , where m 2 c is the value at which V ′′ vanishes at the origin, resums some two-loop corrections into the mass [18]. On the other hand, if λ 2 < 0, one obtains The effective potential for the (λ 1 , We now look at the numerical results in the weak-coupling regime and compare them to the one-loop potential results. The predictions from (17) hold, on average, at the 10 − 30% level for the values of (λ 1 , λ 2 ) we have analyzed in the weak coupling region, and should be more accurate away from the phase transition region. It is indeed natural to expect that deviations are indeed larger near the phase transition surface, at least when the transition is weakly first order (as exemplified e.g. by fig. (3)), since the precise location of the minimum of the potential is in this case unstable against small corrections in its shape originating from two-loop corrections and beyond.
In fig. (4) we plot the results at (λ 1 = −0.22, λ 2 = 0.5) where the correlation length was estimated to be ξ ∼ 3 from the effective potential. The smallest lattice size where metastability was observed was for L = 6. The transition is a relatively strong first order one.
Comparing the evolution of v as a function of m 2 against the effective potential prediction we see that the agreement is good. Following our general argument we expect the transition at (0,0.5) to be weaker since we are away from the stability line. The one-loop effective potential calculation gives ξ ∼ 40, so it is unlikely that we can see metastability, even on our larger lattices. Also, we expect the effective potential calculation to become less reliable.

Strong coupling
The strong coupling region must be studied numerically. The strategy we employed was the following. We studied the smaller lattices, L = 4, 6, 8, using the hybrid algorithm usually, accummulating about 10 5 configurations. We searched for two minima in the histograms corresponding to the expectation value of the operator Tr(φφ † ). We then moved to bigger lattices L = 12, 14 to look for coexistence. Generally, coexistence was found on larger lattices for values of m 2 slightly more negative than on smaller lattices, due to finite size corrections.
Along the process of increasing the lattice size we eventually begin to see metastability at some size L * . We estimate then the correlation length to be ξ ∼ L * /2. Crude as this procedure may seem, it is physically meaningful and it agrees, where comparison is possible, with the effective potential. There is evidence that the transition is stronger in the second case as the two signal minima can be seen in a L = 10 lattice. The transition gets indeed stronger with increasing λ 2 .
For those points deep in the λ 1 > 0, λ 2 > 0 region that we analyzed, we were able to observe coexistence of phases, but only in lattices of L = 14. The transition is always clearly first order, but characterized by correlation lengths much larger than those obtained close to the stability line (this is of course as it should be, given the form of the RG trajectories).
In figure (9) a plot is shown for the point (8,8), where the system eventually tunnels to the right minima. The symmetric phase is in that case a relatively short-lived metastable state.
In fig. (10) we plot the Monte Carlo time evolution of the operator Trφφ † , starting with ordered-disordered initial conditions for the point (8,16).
For points close to the λ 2 = 0 axis, it is very difficult to differentiate a weak first order from a second order transition. More detailed methods with very high-statistics would be needed [17] complemented with finite-size scaling.
We have summarized the knowledge we have gained about the value of the order parameter at the transition and the corresponding correlation length in Table (1). From these results we see that in most of parameter space (in the region where the symmetry breaks the way we are interested in for phenomenological reasons) the vacuum expectation value v (at scale Λ, v(Λ)) jumps to a value which is typically only one order of magnitude smaller than the cut-off. The physically relevant v.e.v. v(v) (that, after gauging, gives a mass to the W ± and

Conclusions
In this paper we report an extensive Monte Carlo simulation of the U (2) × U (2) model. We have found no evidence of the existence of any fixed point other than the gaussian one at the origin of the (λ 1 , λ 2 ) plane.
We have investigated many points in this plane using a variety of numerical and analytical We found just one region satisfying the requirement that the phase transition is weak enough. This is λ 2 → 0, the limit where the model approaches the O(N ) linear sigma model.
For small values of λ 1 , the tunning in λ 2 must be extraordinarily accurate, probably at the 10 −3 precision level or more. This is evidenced by the effective potential calculations. For larger values of λ 1 this is somewhat relaxed, as the jump in the order parameter seems to increase more slowly as we depart from the λ 2 = 0 line for a fixed value of λ 1 . Phenomenologically viable models must then lead to values for the effective couplings which, at the cut-off scale satisfy the above requirements.
All our data conforms perfectly with the standard picture of first order phase transitions with runaway trajectories deduced from the Coleman-Weinberg analysis. We have some evidence that the running is in some cases very fast.
We have left for this appendix all the technical details of the numerical simulations. We have mostly employed the hybrid algorithm since it allows for a better control of the autocorrelation times and the rejection percentage.
We consider the generalized hamiltonian where H is the hamiltonian of the physical problem at hand (in our case H is just the euclidean action of the U (2) × U (2) model), and π some fake momenta conjugate to each variable φ.
We especify some initial values for the momenta π according to a gaussian distribution, and then numerically integrate the Hamilton equations for the (φ, π) dynamical system. Any algorithm can be used provided that is time-reversal and preserves the area of phase space [19].
A convenient way of satisfying both requirements is to use the leap-frog algorithm where F = − ∂H ∂φ , and A is some (arbitrary) t-independent matrix. In the above expressions we use a vector notation for φ, π and F , the vector index running over all lattices sites. After a number of leap-frog steps, the resulting configuration is subject to a standard Metropolis test. It can be either accepted or rejected, and in the latter case we start anew. Using just one leap-frog step the hybrid algorithm would be strictly equivalent to the Langevin one (the fake conjugate momenta playing the role of the gaussian Langevin noise), except that here we must pass the Metropolis test, which makes the algorithm exact. In general it will be convenient to use several leap-frog steps before attempting the Metropolis test.
We have tried two different choices for the matrix A: the identity, A = I (standard hybrid algorithm (SH)), and where ε(p) = (p 2 + m 2 ) −1 is the free lattice propagator. The latter correspond to the Fourier accelerated hybrid algorithm (FA) [19], and it allows for an update of all modes with similar efficiency. This, combined with the decorrelation induced by the numerical integration, makes for a very robust algorithm as far as beating critical slowing down goes. However, due to the need of performing fast Fourier transforms in four dimensions, FA is intrinsically much slower than SH. The gains in beating critical slowing down are only apparent for large correlation length. However, at (λ 1 , λ 2 ) = (−0.22, 0, 5) in a 8 4 lattice, the autocorrelation time at m 2 = −0.895 (transition point) is about four times bigger for the SH than for FA, making FA useful but not really necessary. This is perhaps not too surprising since the correlation length is in much of the parameter space relatively small, even close to the transition.
Two parameters have to be adjusted in the hybrid algorithm, namely the number of leapfrog steps and the step size δt. They are the equivalent to the fudge parameter one uses in a standard Metropolis to adjust the acceptance rate. If we use a relatively large step size δt, successive configurations soon become more uncorrelated. However a large step will decrease the acceptance rate, so a compromise must be reached. We have done extensive tests in the simple case λ 1 = 0, λ 2 = 0, which can of course be solved analytically. The best situation seems to be to take δt such that the acceptance rate is about 90%. On lattices ranging from 8 4 to 12 4 this corresponds to taking δt ∼ 0.2 (depending somehow on the values of (λ 1 , λ 2 ). The other freedom concerns the number of leap-frogs before the rejection Monte Carlo is performed. The larger the number of leap-frogs, the smaller the autocorrelation time, but the required computer time increases too and, at some point, nothing is gained by decorrelating even less our observables (providing the rejection rate remains low). In our case, the optimal choice for the 8 4 , 12 4 lattices were between 5 and 7 leap-frog steps.
After taking all these precautions, the hybrid algorithm works remarkably well.  The results correspond to a L = 4 lattice with 10 5 configurations after thermalization.