Quantum Monte Carlo for Gauge Fields and Matter without the Fermion Determinant

Ab-initio Monte Carlo simulations of strongly-interacting fermionic systems are plagued by the fermion sign problem, making the non-perturbative study of many interesting regimes of dense quantum matter, or of theories of odd numbers of fermion flavors, challenging. Moreover, typical fermion algorithms require the computation (or sampling) of the fermion determinant. We focus instead on the meron cluster algorithm, which can solve the fermion sign problem in a class of models without involving the determinant. We develop and benchmark new meron algorithms to simulate fermions coupled to $\mathbb{Z}_2$ and $U(1)$ gauge fields in the presence of appropriate four-fermi interactions. Such algorithms can be used to uncover potential exotic properties of matter, particularly relevant for quantum simulator experiments. We demonstrate the emergence of the Gauss' Law at low temperatures for a $U(1)$ model in $(1+1)-$d.

Introduction.-Microscopicmodels involving fermions that strongly interact with each other, either directly or mediated via gauge fields, are essential ingredients of many theories in condensed matter and particle physics [1][2][3].From the Hubbard model describing the physics of correlated fermions, to the quantum Hall effect and high-temperature superconductivity, fermions subjected to various interactions have been studied both perturbatively and non-perturbatively [4,5].Fermions constitute the matter component of all microscopic theories of particle physics (as leptons in electromagnetic and weak interactions, as quarks in strong interactions) and interact with gauge fields (the photon, the W ± , Z, and the gluons respectively) [6].Gauge fields are also becoming increasingly important to condensed matter systems, from frustrated magnetism to theories of deconfined quantum criticality [7].
While Quantum Monte Carlo (QMC) methods are robust for non-perturbative studies of the aforementioned systems, they are also vulnerable to the sign problem [8].QMC methods work by performing importance sampling of fermion and gauge field configurations that make up the partition function.Since fermions anti-commute, their sign problem can be straightforwardly understood when the configurations considered are worldlines: whenever fermions exchange positions an odd number of times, the configuration weight acquires another negative sign factor, leading to huge cancellations in the summation, accompanied by an exponential increase of noise [9].
A large family of QMC methods deal with the fermion sign problem using determinants: they introduce auxiliary bosonic fields and integrate out the fermions, or expand the partition function, Z, as powers of the Hamiltonian (or parts of the Hamiltonian) and get fermion determinants for the resulting terms [10,11].Because these methods result in weights that are the sums of many worldline configurations, they can be used more generically to simulate the largest classes of sign-problem-free Hamiltonians, with Auxiliary Field QMC as the most widely applicable method [12][13][14][15][16][17][18][19][20].Determinantal methods in general scale with either the spatial lattice volume, or the spacetime lattice volume, which in terms of imaginary time β and spatial lattice dimension N goes either as O(β 3 N 3 ) or O(βN 3 ), depending on the method [21][22][23][24][25][26].While this polynomial scaling is much better than the exponential scaling from a straightforward exact diagonalization (ED), it is much worse than the linear scaling achievable for spin systems, where worldline-based methods can simulate systems orders of magnitude larger than typical simulable fermionic systems [27][28][29][30][31][32].
An alternative approach-well-utilized in the lattice quantum chromodynamics (LQCD) community, is the Hybrid Monte Carlo technique [33][34][35], which computes the fermion determinant stochastically, and theoretically scales linearly rather than cubicly with the spatial volume.For systems of massless fermions and thus zero modes in Z, however, the method can run into complexities and the scaling significantly worsens, closer to the cubic scaling of before [36].More recently, it has been applied to problems in condensed matter with some promising optimizations [37][38][39].
It is however possible to develop worldline-based algorithms for fermionic Hamiltonians [40][41][42][43].Meron cluster methods [44], so named due to the presence of merons (half-instantons) in the first model for which they were developed to simulate (the 2d O(3) sigma model with θ = π), can be used to solve sign-problems in fourfermion Hamiltonians for certain parameter regimes, as well as for free fermions with a chemical potential [45].Because these methods sample worldlines, computing the weights scales linearly with the volume of the system, and negative terms in the partition function are taken care of by avoiding merons -this is what distinguishes them from bosonic simulations.The relative simplicity of these methods, with each weight corresponding to a worldline configuration and the lack of stabilization issues that can arise in determinantal methods [11], as well as the favorable scaling of the weight computations, makes them an attractive choice for simulation when applicable.Correspondingly, exciting opportunities open up when new interesting physical models are found which can be simulated using this method [46][47][48].
Recently, there has been intense experimental development to study the physics of confinement and quantum spin liquids [49][50][51][52][53][54][55] using tools of quantum simulation and computation.The microscopic models used to capture the physics involve fermions interacting with (Abelian) gauge fields.In this Letter, we develop meron cluster algorithms for a class of experimentally relevant models [56,57], enabling a robust elucidation of their phase diagrams.We also introduce new classes of Z 2 and U (1) multi-flavored gauge-fermion theories, which might be realized in cold-atom setups and also be further studied using Monte Carlo techniques.Notably, the U (1) family of these models seems to be one of the few families that falls outside the class of models known to be simulable by auxiliary-field methods, as are [23,58,59].Moreover, the worldline nature of the method makes it easily employed to study the corresponding phases in these theories in higher spatial dimensions, and the resulting physically-relevant configurations are promising inputs for machine learning algorithms.
Models.-We start with the half-filled t-V model-a spinless fermionic Hamiltonian involving only the most local interactions, Here ⟨xy⟩ are the nearest neighbor site pairs, c † , c are the creation and annihilation operators respectively, and the repulsive interaction V is given in terms of the number operator n = c † c.It is simulable by meron clusters for V ≥ 2t [45,48].In this Letter, we extend the meron cluster method to physically interesting Hamiltonians involving gauge fields, which are lower-dimensional versions of quantum electrodynamics (QED) [60,61].The Z 2 -and U (1)-gauge symmetric families are given by: The label g ∈ {U (1), Z 2 } is the gauge symmetry, with The hopping of spinless fermions between the nearest neighbors ⟨xy⟩ are now governed by the presence of gauge fields, represented by spin-1/2 operators, s k xy , on the bond.Fig. 1 illustrates the model degrees of freedom.Then H (g),des ⟨xy⟩,f is a designer term [62] that makes the mod- els simulable by the meron algorithm, These terms serve a similar role as the V = 2t term in the meron cluster algorithm applied to the t-V model [45], and similarly an additional particle-hole symmetric V ≥ 0 term can be added to the models here in a signproblem-free way.In the Z 2 gauge theory, the gauge field s 1 = σ 1 /2 couples to fermions, and the local Z 2 symmetry is manifest via the commutation [Q x , H Z2 ] = 0, where Q x = (−1) f n x,f f, α s 3 x,x+ α,f s 3 x− α,x,f , α are the unit vectors in a d-dimensional square lattice.For the U (1) theory, the unitary operator V U (1) , which commutes with H U (1) is given by V U (1) = x e iθxGx , with x − 1)/2 .In the terminology of gauge fields, our microscopic models are quantum link models [63], which realize the continuous gauge invariance using finite-dimensional quantum degrees of freedom.The identification with usual gauge field operators is . We note that a straightforward application of the meron idea necessitates the introduction of an equivalent flavor index for gauge links as fermion flavors.Naively, the total Gauss law can be expressed through a product (Z 2 ) or sum (U (1)) of the Gauss law of individual flavors degrees of freedom, and the resulting theories have and U (1) ⊗N f gauge symmetry.However, flavored gauge-interactions can also be turned on in the U (1) model (as explained in the Supplementary Material), or through a Hubbard-U interaction for both Z 2 and U (1)-symmetric models [48].These additions directly cause ordering for either the gauge field or the fermions, and the coupling between them leads to the interesting question of how the other particles are affected by this ordering.In similar contexts, interesting simultaneous phase transitions of both the fermions and gauge fields have been found [47,61,[64][65][66][67][68] or conjectured [69,70].
Algorithm.-Thealgorithm is best understood through the worldline configurations for the models defined in Eqs.(2) to (4) in the occupation number basis for the fermions and the electric flux (spin-z) basis for the gauge links.The partition function in (1 + 1)-d is where H = H e +H o , and H e (H o ) consists of Hamiltonian terms that correspond to even (odd) links.This Trotterized approximation, is a sum of terms over discrete timeslices 1, • • • , 2N t , each with locally defined electric flux and fermion occupation numbers.All terms within H e and H o commute with each other (there are straightforward generalizations for higher dimensions) [71].Each of the terms in Eq. ( 6) is a worldline configuration, and the rules for allowable worldline configurations apply consistently to all models within each of the symmetry families.The worldline configurations are a tool to obtain meron cluster configurations by introducing appropriate breakups, which decompose the terms in Eq. ( 6) into further constituents.In considering the allowed worldline configurations given in Fig. 2 for the U (1) theory, for example, each of the active plaquettes in each timeslice (shaded in gray) must be one of the plaquettes given in Table I.The plaquettes in each row share the same weight, computed using 6), where b is a nearest neighbor bond, b = {x, y}.The corresponding breakup cell for each row gives allowable breakups: if all fermion occupations/spins are flipped along any one of the lines, the resulting plaquette also exists in this table.From the table, we see two such breakups are defined, A and D. Although these breakups resemble those from the original meron algorithm, in our case the breakups involve the link variables as well-either as additional lines for the A breakups, or as binding lines extending outward from the horizontal D breakup lines.This is a key difference for the gauge extension of the algorithm.By computing the matrix elements that correspond to the plaquettes in each grouping, we find that for the U (1) theory, the corresponding breakup weights w A and w D must obey: to satisfy detailed balance.Moreover, the choice of the breakups is such that the total sign of a configuration factorizes into a product of the signs of each cluster: Sign where the configuration C has been decomposed into N c clusters.We can thus simulate this system by exploring a configuration space where each configuration is defined according to the combination of worldlines and breakups.By assigning breakups to all active plaquettes, clusters are formed, and then updates involve flipping all fluxes and fermions within a cluster, which generates a new worldline configuration.The algorithm begins with putting the system in a reference configuration, defined by the fermionic worldlines only, where the weight is known to be positive, and it is always possible to reach this configuration by appropriately flipping a subset of clusters in a given configuration.For both U (1) and Z 2 theories, the reference configuration has a staggered fermionic occupation (charge density wave, or CDW), where fermions and fluxes are stationary throughout imaginary time.Fluxes can be in any spatial configuration (because they do not contribute a sign), and the breakups are all A. Fluxes and breakups may be initially attached to the plaquettes in any way allowed by Table I.A QMC sweep is then: change.If it contains a cluster where flipping the fermion occupation causes the fermions to permute in a way that produces a negative sign, then it is a meron.In that case, restore the breakup back to its initial state.Rules for identifying merons generalize [45,72] and are in the Supplementary Material.This describes sampling of the zero-meron sector only, but sectors with other numbers of merons may become relevant depending on the observable [45].We note that the cluster rules implement the Hamiltonian dynamics, but the constraints due to Gauss' law are not included.Like any cluster algorithm, once the detailed balance conditions have been satisfied, the meron algorithm is expected to be efficient in any space-time dimension [46,72].We provide a demonstration of the efficiency in (2 + 1)−d in the Supplementary Material, with an extensive investigation left for future work.Numerical Results.-To illustrate the efficacy of the algorithm, we discuss results obtained by simulating the 2), which is related to the massless quantum-link Schwinger model [57,73] and the PXP model [74,75], where quantum scars were first demonstrated experimentally [50].We simulate the model for different temperatures β = 1/T , without imposing Gauss' law.A filter may then be applied to study the physics in the desired Gauss law sector.The onedimensional nature of the problem forbids the presence of merons, providing a technical simplification.The first non-trivial result is the emergence of two Gauss' law sectors at low temperatures, as shown in Fig. 3.For the Z 2 theory in (1 + 1)-d, this result was also observed in [68].
Generating different Gauss' Law sectors has the benefit that the physics in each sector can be easily studied by applying a filter.At low temperatures, this is an O(1) effort, but becomes exponentially difficult at higher temperatures, since exponentially many sectors will be populated.Hence, we note that the efficiency of this meron algorithm for true gauge theories (where Gauss' law is imposed) is more suited to the study of quantum phase transitions rather than finite-temperature ones.For theories where multiple non-trivial Gauss' Law sectors emerge at low temperature, it is possible to study the physics in all sectors without any extra effort.
Similar to the well-studied Schwinger model [76][77][78][79][80][81][82], our model has the following discrete global symmetries: Z 2 chiral symmetry, charge conjugation, C, and parity, P , [57], whose breaking depends on the strength of the four Fermi coupling.The order parameter sensitive to the P or the C symmetry is the total electric flux, E = 1 Lt x,t s 3 x,x+1 , while the one for Z 2 chiral symmetry is the chiral condensate, ψψ = x (−1) x n x .In Fig. 3 we show the probability distribution for ψψ, which samples the two vacua very well, indicating that at T = 0 the symmetry is spontaneously broken.We use these operators to check the algorithm against exact diagonalization results, as well as explore other features of the phase at low temperatures.We leave these discussions to the Supplementary Material.Here we concentrate on the meron algorithm's performance measured via the autocorrelation function: where O(i) is the measured value at the i-th step of the appropriate operator (whose average is O), and is the running index summed over the MC data, while the autocorrelations are measured τ steps apart.Fig. 3 shows the C O (τ ) for three different operators: E, ψψ, and CDW.
We note that the bosonic E relaxes the slowest, while the fermionic operators relax faster.Even for the slowest relaxing mode, the autocorrelation decreases by more than an order of magnitude within 10 MC steps for the largest lattice at the lowest temperature, demonstrating the efficiency of the algorithm.Finally, we also show the behaviour of the normalized susceptibilities corresponding to E and the CDW operator as a function of temperature for smaller lattices up to L = 22 in the G x = 0 sector.We are able to capture the finite temperature crossover.
Conclusions.-We have generalized the construction of the meron algorithm to cases where staggered fermions are coupled to quantum link gauge fields.This construction of the Monte Carlo algorithm is agnostic to the space-time dimension, and paves the way for ab-initio studies of large scale gauge-fermionic system with odd or even numbers of fermionic flavors, and includes models not simulable using DQMC.While we are able to simulate low temperatures at fixed values of gauge coupling by using two breakups, A and D, it is possible to add different microscopic terms by increasing the allowed ways of bonding the fermions and gauge links.We have also indicated how to include multiple flavors, and multiflavor interactions.Our investigations open up avenues to study quantum link gauge theories coupled to fermions in higher dimensions, which are almost certain to exhibit quantum phase transitions [83].Since the physics of Abelian gauge fields represented by half-integer spins are sometimes related to quantum field theories at θ = π [60], where θ is the topological angle, our numerical method also promises to increase our knowledge of quantum field theories with non-trivial topologies.Possible future extensions include gauge fields with larger spin representation and non-Abelian gauge fields as well.Our methods can be extended to gauge fields with larger spin representation, and hopefully to non-Abelian gauge fields as well, to tackle realistic interacting systems of increasing complexity in particle and condensed matter physics.

Z2-Symmetric Algorithm
The breakups for the Z 2 -symmetric model defined Eq. ( 2) and Eq. ( 3) are given in Table II.Their weights must satisfy the same equations as for the U (1) theory, which are consistent and solved.The key difference is that more plaquettes are allowed for the Z 2 theory compared to the U (1) theory, as seen in the Table II.
The algorithm then proceeds the same way described in the main text.Note that there are the same number and types of breakups, and their weights are the same as for a simulation of the t-V model.The key difference for the gauge theory is that the A-breakup involves an additional vertical spin line.This causes clusters rather than loops to be formed, which constrain the types of fermionic worldlines that are allowed compared to the loops.

More than One Layer of Fermions and Spins
As mentioned in the main text of the paper, the models in (Eq.( 2)) are sign-problem-free and the meron cluster method applies to them.The weights of the matrix elements for one layer of the U (1) model were shown explicitly.Here we give the weights for the N f -layer model for the U (1) theory and Z 2 theory explicitly so that it is clear that they also are sign-problem-free and amenable to the meron cluster method.
First, for N f = 1 of the discrete-time version U (1) theory, we can write out the matrix elements for a single bond as (10) These expressions for the diagonal and off-diagonal matrix elements in the 2 × 2 space are the e ϵt cosh(ϵt) and e ϵt sinh(ϵt) from the main text, and because we are able to use the A and D breakups to simulate the system.This mathematical form for the weights is more obviously generalizable to the N f -flavor version of the U (1) theory, where the matrix elements for a bond are where m = 2 N f and n = 2 3N f − 2 N f , performing the exponentiation then yields where we define O U (1) m to be an "off-diagonal" matrix for the U (1) theory, which is m × m and has entries of 1 everywhere except for on the diagonal, where the entries are zero.
As required all weights are positive, and we have allowing us to use A ⊗N f and D ⊗N f breakups to simulate the system.
The Z 2 matrix elements are found from where this time m = 2 2N f and n = 2 3N f − 2 2N f .Performing the exponentiation then yields where the off-diagonal O Z2 m matrix is defined as As discussed in the main text, the U (1) model families have additional diagonal terms for the spins that may be added: Here we discuss why these additions are sign-problemfree.Note that we will focus on the discrete-time version the Meron Cluster algorithm, but the result also straightforwardly carries to continuous time.This addition is sign-problem-free for a similar reason that the addition of the Hubbard-U term is sign-problem-free, as discussed in [48].We apply the discrete time partition function expansion defined in the main text with general d spatial dimensions, which is given for the new Hamiltonian H ′ by This Hamiltonian is Trotterized into T sets of operators where the operators within each set commute with each other, H ′ = H ′ a1 + H ′ a1 + ... + H ′ a T , thus there are T • N t timeslices in this formula.In one spatial dimension, as seen in the main text, T is 2 and the Hamiltonian was divided into H e and H o , for even and odd, respectively.
Writing in terms of the original U (1) Hamiltonian H (as simulated in this paper), we have where for example in (1 + 1)-d we would have a 1 = o, a 2 = e, and Here we note that because these additional terms respect the gauge symmetry, we have for a generic Trotter operator set a, where h s a (s i ) is just a number that depends on the spin configuration at particular timeslice i, which is s i .
Thus, we end up with a partition function very similar to our original one, but with the worldlines weighted slightly differently: Because of these additional weight factors, the merons no longer have zero weight, so the question is instead whether their weight can be guaranteed to be positive.
To figure out the meron weights, we consider a generic meron, a portion of which is illustrated in Fig. 5, noting that merons only appear in two spatial dimensions and higher, so the images in the figure would represent a cross-section of a (2 + 1)-d (or higher dimensional) lattice, with two Trotter time slices shown, and the others squashed since they are irrelevant to this portion of the meron.Here, we assume a staggered reference configuration of fermions that is identical for both spin layers.In the portion highlighted, the weight factors for the highlighted meron portions of the configurations (a) and (b) will have a positive sign because both layers have the same fermion configurations, and thus the sign contribution from each layer will be identical, leading to a positive product.The configurations (c) and (d) will have a negative sign factor coming from the highlighted meron because fermion occupations are flipped between the two layers.We have already established that the weights coming from the fermionic term are identical.The contributing factors for this meron portion that come from the spins are for configurations (a) and (b), and for configurations (c) and (d), because in the U (1) theory spin degrees of freedom are tied to the fermion occupation numbers on either side of the D-breakups, and so if the fermions match in the two layers the spins must also match, and similarly if the fermions do not match in the two layers for the meron the spins also cannot match.Thus, the sum of the four configurations must be a positive number.We note that this also can be done for negative J-coupling, but there the staggered reference configuration changes to be opposite occupation numbers for the two layers.

Merons and Autocorrelation in the Gauge Models
For purely fermionic models, as shown in [45], a loop is a meron if the quantity n w + n h /2 is even, where n w is the number of temporal windings and n h is the number of fermionic hops in the loop.For the gauge extension of this algorithm, we now have the possibility of binding two or more loops to each other through the gauge, and thus we have the following criterion for a cluster to be a meron

Meron if
n w + n h /2 odd, even # of loops n w + n h /2 even, odd # of loops It can be seen immediately that (26) reduces to the original definition in the case of one loop, and to understand the extension we note that in the case of two loops bound together (or more generally an even number of loops bound together), the result will be a meron in the case of one loop (or an odd number of loops) being a "loop meron" by the original loop meron definition, and the rest of the loops being non-"loop merons."This results in the total n w + n h /2 being an odd number of even numbers plus an odd number of odd numbers, which is an odd number.A similar argument follows for the clusters that consist of an odd number of loops.
As mentioned in the main text, in one dimension there should be no merons.Figure 7 plots the n h and n w for 100 equilibrated clusters at low temperature, and indeed we see that by the criterion, none of these clusters are merons.
In two dimensions, we indeed get merons, and we can compute observables by sampling the zero-and twomeron sectors according to the procedure in [45].Figure 6 gives a similar autocorrelation calculation to the one given in the main text for the E and CDW susceptibilities on a 10 × 10 lattice in two dimensions, and we see that while the autocorrelation times are marginally longer than those in one dimension, they are still quite short.

Checks Against ED
For this paper we implemented the Meron Cluster method for the U (1) theory in 1 + 1-d, and for small lattices (L = 6 and L = 8) we have checked the following observables against exact diagonalization (ED) calcula-     tions: Table III and Fig. 8 summarize these results, which are all for the observables in the Gauss law sector G x = 0.

Large volume physics with the Meron Algorithm
Although the fermion gauge coupling in the model we simulate with the meron algorithm is identical to the quantum link Schwinger model with spin-1 2 , we have additional four fermion terms, which behave differently from the staggered mass term.While this interaction also induces an effective mass term, the mapping to a massive Schwinger model with renormalized parameters is non-trivial.This makes it hard to connect some of the large volume quantities to those usually studied in the literature [76][77][78][79][80][81][82].Nevertheless, we present an empirical study of our measured observables.
First, in Fig. 11 we show the behaviour of the normalized susceptibilities corresponding to E and the CDW operator as a function of temperature for smaller lattices up to L = 22, filtered to include only the G x = 0 sector (an exponential effort for high temperatures).While they converge to their zero temperature values relatively quickly as a function of the volume, we are also able to capture the finite temperature crossover.
In the main text we have shown the probability distribution of ψψ which shows a double peaked structure.We emphasize that each of the peaks correspond to one of the two Gauss' Law sectors which emerge at low temperatures, indicating broken chiral symmetry.However, a corresponding histogram for the E, the order parameter for the CP symmetry (see Fig. 9  ric distribution.This indicates constraining dynamics between the fermions and gauge fields: with the fermions constrained to occupy certain lattice sites at low temperatures, the gauge fields can fluctuate smoothly.This in turn means that the susceptibility corresponding to the electric flux (normalized to unit volume) cancel the cross terms (where x ̸ = y) and converges to the trivial value 0.125, contributed from the contact terms (s 3 x,x+1 ) 2 .We can characterize the fermionic order by measuring the divergence of the susceptibility of the charge density wave as a function of spatial volume, L. As we show in Fig. 10, the susceptibility diverges as L 1.94 , for both βt = 4 and βt = 8.At high temperatures, on the other hand, the fermions can hop much more, which due to Gauss' Law implies that the gauge fields are only allowed certain values.This is manifested in Fig. 9 (right) which shows the prominence of certain electric fluxes while others are highly suppressed.Moreover, this also explains the behaviour of the electric flux susceptibility in Fig. 4 of the main text: at low temperatures there is cancellation between the smooth electric fluxes, while at high temperatures the finite number of E values can multiply coherently to give rise to a large value of χ E .
It is also interesting to note that the presence of the large CDW order also indicates that the model develops a large chiral condensate ψψ.For the model with staggered fermions, it is known that Gauss' Law explicitly breaks the Z 2 chiral symmetry of the model, which consists of a single spatial lattice translation [57].The chiral condensate is thus large at small temperatures and as one increases the temperature, the condensate smoothly decreases and is expected to vanish at large temperatures.Since this is a crossover, there is no sharp behaviour in the condensate.Moreover, due to this explicit breaking of the Z 2 symmetry by the Gauss Law, the chiral condensate undergoes additive renormalization in the Schwinger model.A similar phenomenology also happens in our model, but since we have a four-fermi coupling instead of a mass term, it is non-trivial to match the two models quantitatively without doing an extensive program involving non-perturbative renormalization.Nevertheless, on comparing the bare chiral condensate at low temperatures on large lattices in the thermodynamic limit, with the corresponding values quoted for the Schwinger model in Fig 15 of [79], we estimate that our model is equivalent to the Kogut-Susskind Schwinger model on the lattice with a bare mass of am ∼ 0.38 − 0.40.In Fig. 12 (top), we show how the chiral condensate approaches a smooth thermodynamic limit.In the bottom panel of Fig. 12, we also show the measurement of the chiral condensate with increasing temperatures until the condensate vanishes (and even becomes slightly negative at infinite temperatures), indicating a restoration of chiral symmetry at high temperatures.While it is qualitatively the same scenario in the Schwinger model [82], we refrain from getting into a technical discussion about the technical aspects of the renormalized chiral condensates in this paper.Moreover, unlike the Schwinger model CP symmetry is not spontaneously broken in the phase with V = 2t.This leaves open the possibility of accessing the quantum phase transition where CP gets spontaneously broken by tuning other bare parameters in the Hamiltonian.

Fig. 2 (
a)-(c) give examples of such configurations for the t-V model as simulated by meron clusters.In the Z 2 case, for each time-slice a fermion has the option of hopping to an unoccupied nearest neighbor site of the same flavor.The hop flips the flux on the bond between the sites of the same flavour index-this is the result of the s 1 xy operator.Fig. 2(d)-(f) gives example configurations for the N f = 1 version of this model.Due to the trace condition, odd winding numbers are ruled out because these would cause mismatch between the spins in the initial and the final state.The possible worldline configurations for the U (1) case are even more restrictive than the Z 2 case.The s + xy and s − xy operators allow the hopping for a given flavor of fermion only in one direction or the other for each bond, depending on the orientation of the same flavored flux on the bond.Fig. 2(g)-(i) illustrates an example configuration and restrictions for the single flavor version of the U (1) model.In (1 + 1)-d, it is clear that all allowed configurations must have zero winding number.
Plaquettes and breakups for the U (1)symmetric Hamiltonian.The middle cluster lines in the A-breakups and binding lines in the D-breakups distinguishes them from the original meron cluster breakups.

1 . 8 FIG. 3 :
FIG.3: Clockwise from top left: (i) and (ii) Number of configurations versus Gauss law sector indexx [G x + 2] • 4 x (not all indices correspond to actual sectors) for 50000 equilibrated configurations.Two sectors emerge at large β: G x = 0 and G x = (−1) x .(iii) The probability distribution of ψψ, with peaks from the two emergent Gauss' Law sectors, indicating that the algorithm efficiently samples all sectors.(iv) The autocorrelation functions for different operators.

2 .
Identify the new clusters formed by the breakups in the new configuration.For each cluster, flip all fermions and fluxes with probability 1/2.

FIG. 4 :
FIG.4: Finite temperature data for U (1) theory in 1 + 1d.The dotted lines show the χ EE , which is the susceptibility corresponding to E. This value rapidly converges to 0.125.On the other hand, the dashed lines trace the χ CDW which display more finite size effects.Thermal behaviour of both observables indicate that the transition from low to high temperature is a smooth crossover.
Plaquettes and breakups for the Z 2symmetric Hamiltonian.

FIG. 5 :
FIG. 5: The spin factors of the weights for the portion of these terms that are shown (assuming the burgundy cluster is a meron) are: e τ1J/4 e τ2J/4 for (a) and (b), and e −τ1J/4 e −τ2J/4 for (c) and (d).Because the merons in (a) and (b) are positive and the merons in (c) and (d) are negative, the fully summed weight of the meron configuration will be positive.This assumes a staggered reference configuration where both layers are identical.

FIG. 12 :
FIG. 12: Bare chiral condensates for U (1) theory in 1 + 1-d at low temperature and large lattices (top), while the temperature dependence of the condensate is bottom part of the figure.
FIG. 7:The n h and n w for 100 equilibrated clusters plotted for the U (1) theory in one dimension.
The dotted lines show the χ EE , which is the susceptibility corresponding to E. This value rapidly to 0.125.On the other hand, the dashed lines trace the χ CDW which display somewhat more finite size effects.The thermal behaviour of both observables indicate that the transition from low to high temperature is associated with a smooth crossover, consistent with the literature.