arXiv:1302.1622v2 [hep-lat] 22 May 2013
Thirring model at finite density in 0 + 1 dimensions
with stochastic quantization: Crosscheck with an exact solution
Jan M. Pawlowski
1, 2
and Chr istian Zielinski
1, ∗
1
Institut für Theoretische Physik, Universität Heidelberg, Philosophenweg 16, 69120 Heidelberg, Germany
2
ExtreMe Matter Institute EMMI, GSI, Planckstraße 1, D-64291 Darmstadt, Germany
(Dated: July 8, 2018)
We consider a generalized Thirring model in 0 + 1 dimensions at finite density. In order to
deal with the resulting sign problem we employ stochastic quantization, i.e., a complex Langevin
evolution. We investigate the convergence properties of this approach and check in which parameter
regions complex Langevin evolutions are applicable in this setting. To this end we derive numerous
analytical results and compare directly with numerical results. In addition we employ indirect
indicators to check for correctness. Finally we interpret and discuss our findings.
PACS numbers: 05.50.+q, 71.10.Fd
Keywords: Thirring model, finite density field theory, complex Langevin evolution, stochastic quantization
I. INTRODUCTION
Despite all effor ts, one of the outstanding problems of
lattice field theory until this day is the sign problem. The
introduction of a finite chemical potential µ > 0 renders
the path integral measure complex and rapidly oscillating
in many theories of interest, like quantum chromodynam-
ics (QCD) in 3 + 1 dimensions. The oscillatory behavior
significantly increases the numerical costs, in particular
in the continuum limit. A hard sign pr oblem ex ists for
theories where the costs grow more than polynomially
with the volume. This obstacle hinders numerical ab ini-
tio studies of str ongly interacting matter under extreme
conditions and the understanding of the phas e diagram of
QCD. There is no satisfactory solution known to the sign
problem, despite the large number of propose d solutions.
Among the proposed solutions we can find reweighting
techniques, Taylor expansions about µ = 0, ex trapola-
tions from imaginary chemical potential, the introduction
of dual variables and a canonical ensemble approach. For
recent reviews see e.g. [
1, 2]. However, due to the over-
lap problem reweighting techniques are computationally
expensive and can only be used for small µ, while the
numerical determination of Taylor coefficients is no isy
and the expa nsion converges slowly [
3]. Also the contin-
uation from imaginary chemical potential is a nontrivial
task [
4]. The applica tion of dual variables and the canon-
ical ensemble approach is still under active research, see
e.g. [
5, 6] for dua l observables and [7, 8 ] for simulations
with canonical ensembles.
In this paper we employ a different a pproach. Parisi
proposed already in 1983 that stochastic quantization
[
9]—for a review see e.g. [10]—could circumvent the sign
problem in terms of a complex Langevin evolution [1 1].
However, it is well known that the Langevin evolution
may converg e towards unphysical fixed points. It has
∗
Permanent address: Division of Mathematical Sciences, Nanyang
Technological University, Singap ore 637371
bee n successfully applied to the SU(3) spin model [12, 13],
to an effective theo ry of QCD in the strong- c oupling limit
[
14], simple models of quantum chromodynamics [15, 16]
and to the relativistic Bose gas [17, 18] at finite den-
sity. Furthermore it has been applied to quantum fields
in Minkowski time [
19, 20], also in nonequilibrium [21].
Counterexamples are given by the three-dimensional XY
model at finite chemical potential for small β [
22] and in
cases of gauge theories with static charges [
23]. Early in-
vestigations of complex Langevin evolutions can be found
in [
24–26], while for reviews see e.g. [27, 28]. Rec e ntly,
a set of consistency conditions indicating correct conver-
gence could be derived [29–31]. When truncating this
infinite towe r of identities one obtains necessary condi-
tions for correctness.
In this work we apply a complex Langevin evolution
to a generalized Thirring model at finite density. Here it
serves us as a model theor y to check for the applicabil-
ity of this method. Our results extend the studies ca r-
ried out in [
32], which led to ambiguous results. In this
paper we restrict our selves to the case of 0 + 1 dimen-
sions and dea l with the ques tion of whether a complex
Langevin evolution can enable finite density c alculations
in this setting. Further investigations of this approach
in the 2 + 1-dimensional generalized Thirring model are
presented in [
33]. The 2 + 1-dimensional model appears
for example in effective theories of high temperature su-
perconductors and graphene, see e.g. [34] and references
given therein. It is also worth mentioning that in the
case of the three-dimensional massless Thirring model, a
fermion bag approach was successfully applied in [
35].
We organized the paper as follows: In Sec .
II we intro-
duce a generalized Thirring model and its formulation on
the lattice. We discuss the Langev in equation and its nu-
merical implementation. In Sec.
III we present a closed
expression for the partition function of the lattice theory
and derive some obser vables of interest. We also discuss
additional indicators to evaluate the convergence prop-
erties of the complex Langevin e volution. In Sec .
IV we
discuss the results of the numerical part o f this wor k and
aim to answer the question in which parameter regime
2
results are reliable. We end this paper with concluding
remarks in Sec. V.
II. THE GENERALIZED THIRRING MODEL
A. Continuum formulation
We consider a generalization o f the Thirring mo del.
The historical model was intro duced in 1958 by Walter
E. Thirring and is one of the rare examples of an exactly
solvable quantum field theory [
36]. While the original
model describ es self-interaction fermions in 1 + 1 dimen-
sions, we c onsider N
f
fermion flavors at finite density.
We begin with a ge neralization to d dimensions and
then later specialize to the case of 0 + 1 dimensions. The
Euclidean Lagr angian in the continuum reads
L
Ψ
=
N
f
X
i=1
Ψ
i
/
∂ + m
i
+ µ
i
γ
0
Ψ
i
+
g
2
2N
f
N
f
X
i=1
Ψ
i
γ
ν
Ψ
i
2
. (1)
The index i = 1, . . . , N
f
enumerates fermion flavors, m
i
and µ
i
denote the bare mass and bare fer mion chemi-
cal potential of the r e spective flavor a nd g
2
is the bare
coupling s trength. The γ matrices satisfy the Clifford
algebra {γ
µ
, γ
ν
} = 2δ
µν
1
.
The model shows breaking of chiral symmetry at µ = 0
in 2 + 1 dimensions [
37]. For the 1 + 1-dimensional
Thirring model, the e quivalence to the sine-Gordon
model can be shown [
38, 39].
The four-po int intera ction can be resolved w ith the
introduction of an auxiliary field A
ν
. This formulation
reads
L =
X
i
Ψ
i
/
∂ + i
/
A + m
i
+ µ
i
γ
0
Ψ
i
+ N
f
βA
2
ν
. (2)
Here we introduced the inverse coupling β = 1/
2g
2
.
When integrating A
ν
out, we recover (
1). Although the
auxiliary field A
ν
is not a gauge field, the model can be
interpreted as a more ge neral gauge theory after gauge
fixing, see e.g. [
40]. After integrating out the fermionic
degrees of freedom we find
Z =
ˆ
DA
Y
i
det K
i
!
e
−S
A
=
ˆ
DA e
−S
eff
,
S
A
= N
f
β
1/T
ˆ
0
dt
ˆ
d
d−1
x A
2
ν
. (3)
Here we introduced the temperature T and K
i
=
/
∂ +
i
/
A+ m
i
+ µ
i
γ
0
. Including the fermion determinant in the
exponential term yields
S
eff
= S
A
−
X
i
Tr log K
i
. (4)
For the fermion determinant the relation
det K
i
(µ) = [det K
i
(−µ
⋆
)]
⋆
(5)
holds, thus rendering the path integr al measure complex
for µ > 0. At vanishing or purely imaginary chemical
potential, the determinant is real and the theory is free
of a sign problem. If the fermion determina nt is replaced
by its modulus, we refer to this as the phase-quenched
case. Physically this corr e sponds to the introduction of
an isospin chemical potential.
Like quantum chromodynamics, the Thirring model
exhibits Silver Bla z e behavior [
41, 42]. It implies that
at vanishing temperature there is a thre shold µ
c
, so that
observables are independent of the chemical potential µ
for µ < µ
c
. While in the full theory the onset is given by
the physical fermion mass m
phys
, in the phase-quenched
theory we have µ
c
= m
π
/2, where m
π
is the physical
pion mass.
B. Lattice formulation
We consider the case of 0 + 1 dimensions—
corresponding to a quantum mechanical system—with
lattice spacing a and N
t
lattice po ints. We employ stag-
gered fermions [
43–46] and denote the number of lattice
flavors, i.e., the number of staggere d fermion fields, by N.
Furthermore we assume that N
t
is even, as otherwise the
formulation of stagge red fermions is conceptually prob-
lematic and (
5) is violated. In order to introduce a finite
chemical potential µ, we use the prescription by Hasen-
fratz and Ka rsch [
47]. Fo r notational ease we refer to
the one-component auxiliary field as A
t
= A
0
(x = t)
and all dimensionful quantities are scaled dimensionless
by appropriate powers of a. Furthermore we introduce
the hopping parameter κ = 1/ (2m). The tempera tur e
corresponds to the inverse temporal extension T = N
−1
t
.
Using this fo rmulation the lattice partition function
reads
Z =
∞
ˆ
−∞
N
t
Y
t=1
dA
t
Y
i
det K
i
!
e
−S
A
(6)
with S
A
=
1
2
Nβ
P
t
A
2
t
and flavor index i = 1, . . . , N.
The fermion matrix takes the form
K
i
(t, τ) =
1
2
(1 + iA
t
) e
µ
i
δ
t+1,τ
−
1
2
(1 − iA
τ
) e
−µ
i
δ
t−1,τ
+ m
i
δ
tτ
, (7)
where we impose antiperiodic boundary conditions,
cf. [
32, 48]. In our analysis we focus on a few observ-
ables, namely the fermion density and condensate, the
3
energy density and the phase factor o f the fermion de-
terminant. In the following, sums over the flavor index i
are not implied. The fermion density of a given flavor is
given by
hn
i
i =
1
N
t
∂ log Z
∂µ
i
T
=
1
N
t
Tr
∂K
i
∂µ
i
K
−1
i
. (8)
The fermion condensate follows from
h
χ
i
χ
i
i =
1
N
t
∂ log Z
∂m
i
T,µ
i
=
1
N
t
Tr K
−1
i
(9)
and the energy density reads
hε
i
i = −
∂ log Z
∂N
t
µ
i
+ µ
i
hn
i
i, (10)
which we normalize to hε
i
i(µ = 0) = 0.
The phase factor of the determinant is defined by
exp (iφ) = det K/ |det K|. It can be expres sed in terms
of the partition function
Z
N
=
∞
ˆ
−∞
Y
t
dA
t
(det K)
N
e
−S
A
, (11)
for N degenerated flavors and Z
pq
N
for the phase-
quenched case, where the fermion determinant in (
11)
is r e placed by its modulus. The expectation va lue of
exp (iNφ) follows in the N flavor phase-q uenched theory
[
49, 50] as
e
iN φ
pq
N
=
Z
N
Z
pq
N
∈ [0, 1] . (12)
A value close to zero indicates a rapidly oscillating path
integral measure with a severe s ign problem.
C. Complex Langevin evolution
The idea of stochastic quantization is that observables
in a Euclidean quantum field theory can be obtained as
the equilibrium values of a statistical system coupled to
a heat bath [
10]. The problem of quantizing a field the-
ory is then reduced to finding the static solutions of an
associated Langevin equation. If the action is real and
bounded from below, c orrectness of this approach can be
ensured. We can also formally generalize to the case of a
complex action [
11]. This situation naturally arises when
considering field theor ie s at finite density. Until this day
there is a lack of rigor mathematical understa nding re-
garding the validity of this procedure. However, in cases
where it is converging correctly one has a very elegant
solution for the sign problem at hand.
We aim to check for the applicability of complex
Langevin evolutions to the Thirring model. To this end
we have to find the static solution of the Langevin equa-
tion
∂
∂Θ
A
t
(Θ) = −
δS
eff
[A]
δA
t
(Θ)
+
√
2 η
t
(Θ) , (13)
where Θ denotes a fictitious time. The noise term η
t
(Θ)
follows a Gaussian distribution with
hη
t
(Θ)i = 0,
hη
t
(Θ) η
t
′
(Θ
′
)i = δ (t − t
′
) δ (Θ −Θ
′
) . (14)
A simple approach to solve the Langevin equation nu-
merically is a first order integration scheme with fixed
stepsize ǫ
L
. Higher order integration schemes of O(ǫ
3/2
L
)
have been employed in the literature too [
13, 51]. How-
ever, in some models fixe d stepsize integration schemes
fail due to the occurrence of r un-away trajectories, which
can be avoided by the use of an ada ptive stepsize [
52, 53].
Although a constant stepsize proved here to be sufficient
[
32], we employ an adaptive stepsize algorithm due to
better convergence properties. For N deg enerated fla -
vors our discretization of (
13) reads
A
t
(Θ + ǫ
L
) = A
t
(Θ) + ǫ
L
D
t
(Θ) +
√
2ǫ
L
η
t
(Θ) (15)
with drift term
D
t
(Θ) = −NβA
t
(Θ)
+
Ni
2
K
−1
(t + 1, t) e
µ
+ K
−1
(t, t + 1) e
−µ
. (16)
After each integration step the stepsize ǫ
L
will be up-
dated according to
ǫ
L
≡ ǫ
L
(Θ) =
δ
max
t
|D
t
(Θ)|
(17)
with stepsize parameter δ = 10
−3
(compare to [
32]).
It is possible to generalize the real noise term in (15)
to an imagina ry one [
32] via the replacement
η
t
(Θ) →
√
I + 1 Re η
t
(Θ) + i
√
I Im η
t
(Θ) (18)
with I ≥ 0. The noise correlators then read
hRe η
t
(Θ) Re η
t
′
(Θ
′
)i = hIm η
t
(Θ) Im η
t
′
(Θ
′
)i
= δ (t − t
′
) δ (Θ −Θ
′
) (19)
and hRe η
t
(Θ) Im η
t
′
(Θ
′
)i = 0. Assuming correctness of
the complex Langevin evolution and numerical stability,
we expect ex pectation values to be independent of I.
III. ANALYTICAL RESULTS
A. Exact partition function
We begin with the partition function for one staggered
fermion field, i.e., N = 1. We incorpo rate antiperiodic
boundary conditions and for brevity we introduce
B
±
=
1
2 (2κ)
N
t
1 ±
p
B
c
/β
N
t
,
4
(a) Plot of the fermion density hni.
(b) Plot of the fermion condensate hχχi.
Figure 1: In the phase structur e in d = 0 + 1 we find a
condensed phase for large µ.
B
c
= β + 4 (β + 1) κ
2
. (20)
Then the partition function (
6) reads
Z
1
= 2
π
2β
N
t
/2
[B
+
+ B
−
+ cosh (N
t
µ)] . (21)
This can be shown for example by systematic saturation
of the Grassmann integral or the help of the determinant
identities in [
54]. For the fermion density we find
hni =
sinh (µ/T )
B
+
+ B
−
+ cosh (µ/T )
, (22)
while the fermion condensate is given by
h
χχi =
2κ
p
β/B
c
(B
+
− B
−
)
B
+
+ B
−
+ cosh (µ/T )
. (23)
The expression for the energy density is rather lengthy
and we will not quote it here explicitly. Figure
1 shows
the dependence of these observables on both β and µ . For
large µ we find a condensed phase, which is well separated
for large N
t
.
As it turns out, we can ta ke the continuum limit of
the density hni analytically. To this end we recover the
physical units of all dimensionful quantities by r e intro-
ducing the lattice spacing a . We fix the dimensionful
temper ature T
−1
= aN
t
, express the lattice spacing a as
a function of the number of lattice points N
t
and take
the limit N
t
→ ∞. We obtain
hni
cont
=
sinh
µ
T
1
2
exp
κ−β
2T βκ
1 + exp
1
T κ
+ cosh
µ
T
,
(24)
where all units are explicitly dimensionful. In the zero
temper ature limit T → 0 we find hni
cont
= Θ (µ − m
phys
)
with physical fermion mass m
phys
= m+g
2
and bare mass
m.
B. Several flavors
We also cons idered the case of more than one lattice
flavor, but only quote here the simplest case of (
11) for
N
t
= 2 lattice points and N = 2 sta ggered fermion fields
Z
2
=
π
32β
3
κ
4
2β
2
+ 4βκ
2
+ 8β
2
κ
2
+ 5κ
4
+ 12βκ
4
+ 12β
2
κ
4
+ 8β
2
κ
2
cosh (2µ)
+ κ
4
cosh (4µ) + 8βκ
4
cosh (2µ) − 4βκ
4
cosh (4µ)
+16β
2
κ
4
cosh (2µ) + 4β
2
κ
4
cosh (4µ)
. (25)
We observe that for two flavors the density, conden-
sate and energy density have plateaus in the range of
β ≈ 0.3–1.2, see also Fig .
9. They appear between the
onset to the condensed phase and saturation of the den-
sity. The height of the plateau depends on the masses
of the flavors, and in the case of degenerated flavors it
is exa c tly at half of the maximum value of hni or h
χχi,
respectively. If β is very small or large, the platea us
eventually disappear. For large β this can be unders tood
from the weak coupling limit, see Sec.
III C. T he exis-
tence of these plateaus on the lattice is further confirmed
by a heavy dense limit [
33] and the Monte-Carlo studies
in Sec.
IV D. In the general case of N flavors, we can find
up to N − 1 intermediate plateaus in these observables.
We give a natural explanation of these structures in the
Appendix.
5
C. Weak coupling limit
In the limit of β ≫ 1 the path integral measure has a
strong peak at the origin and can be approximated by a
Dirac δ function. For N flavors we find
Z
weak
1
=
2π
Nβ
N
t
2
2
2
N
t
N
×
Y
i
B
+
i
+ B
−
i
+ cosh (N
t
µ)
(26)
for the partition function in the weak coupling limit.
Here we introduced
B
±
i
=
1
2 (2κ
i
)
N
t
1 ±
q
1 + 4κ
2
i
N
t
. (27)
This can be shown again using the identities in [
54]. No te
that in this limit the contribution from different flavors
factorize and the phase-quenched case equals the full the-
ory. For N degenerated flavors, i.e., B
±
i
= B
±
, the total
fermion density is given by
hni =
N sinh (µ/T )
B
+
+ B
−
+ cosh (µ/T )
(28)
and the fermion co ndens ate by
h
χχi =
2κ N (B
+
− B
−
)
√
1 + 4κ
2
(B
+
+ B
−
+ cosh (µ/T ))
. (29)
While approaching the noninteracting limit, the plateau
structures observed for N > 1 disapp e ar.
D. Analyticity in µ
2
An observable O which is even in µ can be interpreted
as a function of µ
2
. Assuming analyticity in µ
2
, we can
analytically continue O to purely imaginary chemical po-
tential, i.e., µ
2
≤ 0. In this cas e the fermion determinant
is real and free of a sign problem due to (
5) and we can
employ a real Langevin evolution. For µ
2
> 0 we employ
a complex Langevin evolution.
Assuming the correctness of the complex Langevin evo-
lution, O sho uld be analytic at µ
2
= 0. Any nonanalytic-
ity would be an indicator for incorrect convergence. This
criterion was previously employed in [
13] and in this work
we apply it to the condensate h
χχi.
E. Consistency conditions
In [
31] the authors derived a set of identities indicating
correct converge nce of expectation values obtained by a
complex Langev in evolution. These consistency condi-
tions state that for all entire holomo rphic observables O
the expectation value hLOi = 0 has to vanish. Here
L
t
=
d
dA
t
+ D
t
d
dA
t
(30)
10
3
10
4
10
5
10
6
0 1 2 3 4 5 6
Propagator |G(∆)|
2
Distance ∆
Full theory - β = 1, m = 1, N
t
= 4
µ = 0.0
µ = 0.1
µ = 0.2
Figure 2: The propagator at different values of µ.
denotes the Langevin operator and D
t
= −dS
eff
/dA
t
the
drift term.
As this defines an uncountable number of identities,
we re strict the analysis to a finite subset. We follow [
31]
and consider here observables O(t, k) = exp (ikA
t
). The
resulting conditions take the form
hL
t
O(t, k)i =
ik [ik + D
t
] e
ikA
t
= 0, (31)
which have to hold for ∀t and ∀k. Without loss of gen-
erality we set t = 1. Besides O(t, k), we also check the
consistency conditions for the fermion density and the
condensate.
F. Propagator at finite density
We define the propagator at finite temperature T =
N
−1
t
and chemical potential µ by
hχ (t
1
)
χ (t
2
)i =
1
Z
1
K
−1
(t
1
, t
2
)
(32)
with partition function Z
1
given by (
21). It is helpful to
introduce the notation
G (∆) =
D
χ (1)
χ
1 +
ˆ
∆
E
(33)
with
ˆ
∆ = ∆ mod N
t
. For small lattices, the inversion of
the fermion matrix and the calculation of the expectation
value c an b e done analytically. A typical example can
be found in Fig. 2. For small distances the pro pagator
G (∆) falls off exponentially, but due to the finite lattice
extension eventually rises again. The introduction of a
finite chemical potential µ shifts the minimum and results
in G (∆) 6= G (−∆).
6
-0.2
0
0.2
0.4
0.6
0.8
1
1.2
0 0.5 1 1.5 2
Fermion density <n>
Chemical potential µ
Full theory - β = 3, m = 1, N
t
= 4
Theory
CLE
MC
Figure 3: Typical errors estimated by the standard
deviation.
IV. NUMERICAL RESULTS
A. Implementation
In the previous section we derived numerous analyt-
ical results, which allow us to benchmark the complex
Langevin evolution and check for its correctness. For
the numerical s imulations we implemented an adaptive
numerical integration scheme of the associated Langevin
equation as described in Sec.
II C. Using a computer alge-
bra system, we can determine the inverse of the fer mion
matrix in (
7) and find an analy tical expression of the
drift term in (
16) for small lattices in order to minimize
numerical errors.
The evaluation of the observables for a given set of
parameters begins with a hot start, i.e., the random ini-
tialization of the auxiliary field, followed by typically 10
4
steps for thermalization. The field configuration is then
sampled about O
10
5
times and used to evalua te the
observables. To reduce potential autocorrelation effects,
two subsequent samples are separated by ten dummy
steps. In all following plots, numerical values obtained by
a complex Langevin evolution will be denoted by “CLE,”
while analytical results are denoted by “Theory.”
A simple method to estimate the e rror is to take the
standard deviation of the different samples of a given
observable in a particular run. However, the resulting
errors are much larger than the empirical observed sta-
tistical fluctuations between different runs and are over-
estimated, see Fig.
3 for a typical example. In order to
obtain more reliable error bounds, we employ a bootstrap
analysis [
55]. Besides the statistical error, we also have to
deal with systematic err ors induced by a finite stepsize.
To this end we checked that the stepsize parameter δ in
(
17) is sufficiently small.
As an additional test we used an adaptive quasi-Monte
Carlo strategy [56] for the direct evaluation of expecta-
tion values, which works well for sufficiently small N
t
and µ. On larger lattices the sign problem is severe
and the algorithm fails. As the path integral measure
falls off rapidly for large field configurations, we replace
the numerically difficult noncompact integration over the
auxiliary field by a compact domain [−Λ, Λ]
N
t
. We
parametrize this ad hoc cutoff Λ by
Λ =
r
2σ
Nβ
. (34)
If the geometric mean of the field configura tion is Λ, then
the path integral measure is dampened by a factor of
exp (−σ). If we choose σ too small, we introduce a large
cutoff effect. If σ is very large, we still have the problem
that the integrand is close to zero in most o f the inte-
gration reg ion. As a compromise here we use σ = 10.
In figures, we refer to numerica l r esults obtained by this
method as “MC.”
B. Comparison of results
We can see the results of a typical run for two different
sets of parameters in Figs.
4 and 5. Because of the high
sample size, the error bars tend to be very small. We see
that the results obtained with an adaptive Monte-Carlo
method are in good agreement with analytical results.
The numerical results obtained with a c omplex
Langevin evolution show some deviations for intermedi-
ate µ. We systematically observe that the ga p widens for
decreasing values of β, i.e., strong er couplings. In a ddi-
tion the numerical evaluation becomes more noisy as the
relative magnitude of the noise term η
t
increases. For
large β the a greement becomes better and the numerics
are very stable.
From the phase factor of the fermion determinant, we
see how the sign problem gets more severe for increasing
lattice sizes. However, our approach is not affected by
this. It is interesting to note that the sign problem seems
to be less pronounced when using complex Langevin evo-
lutions.
In Fig.
5 we can already see how in the limit T → 0
the Silver Blaze b ehavior becomes apparent. In the limit
N
t
→ ∞ the observables will eventually b ecome indepen-
dent of the chemical potential µ up to some threshold µ
c
.
It seems that the positio n of this onset differs slightly
from analytical predictions. However, for large β they
are in good agreement.
That the observed deviations are not just an effect of
a finite stepsize can be seen in Fig. 6. For a typical set
of parameters we calculate the fermion density and co n-
densate for varying stepsiz e factors δ. The ex trapolation
shows that also in the limit δ → 0 there is a discrepancy
between theory and complex Langevin evolutions.
7
-0.1
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
0 0.5 1 1.5 2
Fermion density <n>
Chemical potential µ
Full theory - β = 3, m = 1, N
t
= 4
Theory
CLE
MC
(a) Fermion density hni.
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0 0.5 1 1.5 2
Fermion condensate <χ
-
χ>
Chemical potential µ
Full theory - β = 3, m = 1, N
t
= 4
Theory
CLE
MC
(b) Fermion condensate hχχi.
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
1.1
0 0.5 1 1.5 2
Phase factor Re <e
iφ
>
Chemical potential µ
Full theory - β = 3, m = 1, N
t
= 4
CLE
MC
(c) Phase factor of determinant
e
iφ
.
Figure 4: Benchmarking results for β = 3.
-0.2
0
0.2
0.4
0.6
0.8
1
1.2
0 0.5 1 1.5 2
Fermion density <n>
Chemical potential µ
Full theory - β = 1, m = 1, N
t
= 8
Theory
CLE
MC
(a) Fermion density hni.
-0.1
0
0.1
0.2
0.3
0.4
0.5
0.6
0 0.5 1 1.5 2
Fermion condensate <χ
-
χ>
Chemical potential µ
Full theory - β = 1, m = 1, N
t
= 8
Theory
CLE
MC
(b) Fermion condensate hχχi.
-0.2
0
0.2
0.4
0.6
0.8
1
0 0.5 1 1.5 2
Phase factor Re <e
iφ
>
Chemical potential µ
Full theory - β = 1, m = 1, N
t
= 8
CLE
MC
(c) Phase factor of determinant
e
iφ
.
Figure 5: Benchmarking results for β = 1.
8
0.45
0.5
0.55
0.6
0.65
0.7
0.75
0.8
0.85
0 0.05 0.1 0.15 0.2
Fermion density <n>
Stepsize δ
Full theory - β = 1, m = 1, µ = 1, N
t
= 4
δ = 0, Theory
CLE
Linear Fit
(a) Fermion density hni.
0.1
0.12
0.14
0.16
0.18
0.2
0.22
0.24
0.26
0.28
0.3
0 0.05 0.1 0.15 0.2
Fermion condensate <χ
-
χ>
Stepsize δ
Full theory - β = 1, m = 1, µ = 1, N
t
= 4
δ = 0, Theory
CLE
Linear Fit
(b) Fermion condensate hχχi.
Figure 6: Extrapolation to vanishing stepsize δ → 0.
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
0 0.5 1 1.5 2
Fermion density <n>
Chemical potential µ
Full theory - β = 100, m = 1, N
t
= 4
Weak Coupling Limit
CLE
MC
Figure 7: Density in the weak coupling regime.
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
1.1
0.5 1 1.5 2 2.5 3
Fermion density <n>
Inverse Coupling β
Full theory - m = 1, µ = 1, N
t
= 4
Theory
CLE
MC
(a) Fermion density hni.
-0.05
0
0.05
0.1
0.15
0.2
0.25
0.3
0.5 1 1.5 2 2.5 3
Fermion condensate <χ
-
χ>
Inverse Coupling β
Full theory - m = 1, µ = 1, N
t
= 4
Theory
CLE
MC
(b) Fermion condensate hχχi.
Figure 8: Observables as function of β.
C. Coupling parameter dependence
In the weak coupling limit of Sec.
III C, we observe
very go od agreement of all numerical and analytical re-
sults. For increasing β the numerical results are asymp-
totically approaching analytical re sults. In Figure
7 we
see that for β & O(100) the method is then able to de-
liver reliable results. In the XY model at finite density
[
57] a s imilar behavior was observed. For small β the
authors observed incorre c t convergence, while for large β
they found agreement. However, in the XY model there
is a clear separation between both regimes.
As already pointed out, the applicability of complex
Langevin evolutions depend on the magnitude of β. In
Figure
8 we find the dens ity and condensate as a function
of β. We observe that for weak couplings, i.e., in the
limit β → ∞, it slowly converges towards the ana lytical
results. For β = O(100) and above we then find good
agreement. Contrariwise for strong couplings we obs e rve
a large discrepancy. Moreover for β . 1 the e valuation
of observables is very noisy, resulting in large error bars.
9
-0.2
0
0.2
0.4
0.6
0.8
1
1.2
0 1 2 3 4 5
Fermion density <n>
Chemical potential µ
Full theory - N
f
= 2, β = 0.6, m = 3, N
t
= 4
Theory
CLE
MC
(a) Fermion density hni.
-0.05
0
0.05
0.1
0.15
0.2
0.25
0.3
0.35
0 1 2 3 4 5
Fermion condensate <χ
-
χ>
Chemical potential µ
Full theory - N
f
= 2, β = 0.6, m = 3, N
t
= 4
Theory
CLE
MC
(b) Fermion condensate hχχi.
Figure 9: Observables for N = 2 degener ated flavors.
D. Several flavors
In Sec.
III B we found that in a certain range of β,
the density, the condensate and the energy density show
up to N − 1 intermediate plateaus when considering N
flavors.
We give an example for N = 2 flavors at a coupling
of β = 0 .6 in Fig. 9. For small lattice sizes and small µ
Monte Carlo studies confirm these predictions. Because
of the sign problem we are restricted to µ . 2 with these
particular parameters, before the numerical evaluation
breaks down.
When trying to reproduce the plateaus with a com-
plex Langevin evolution, we notice that for every value
of β the result qualitatively resembles the N = 1 case.
Directly after the onset the fer mion density rise s until
it reaches saturation. Only in the limit of large β the
plateaus disappear in the analytical solution and there is
agreement with numerical results.
0.01
0.1
1
-1 -0.5 0 0.5 1
Fermion condensate <χ
-
χ>
Chemical potential µ
2
Full theory - m = 1, N
t
= 4
β = 0.5, Theory
β = 1, Theory
β = 3, Theory
β = 10, Theory
β = 0.5, Langevin
β = 1, Langevin
β = 3, Langevin
β = 10, Langevin
Figure 10: Logarithmic plot of hχχi over µ
2
.
E. Analyticity at µ
2
= 0
As described in Sec.
III D, we check for possible no n-
analytic behavior of the c ondensate h
χχi at µ
2
= 0. In
Fig. 10 we see that for µ
2
≤ 0 the rea l Lange vin evolu-
tion is, as expected, in agreement with theo ry. Because
of the abse nce of a sign problem numerical re sults are
reliable in this regime. Also for small µ
2
& 0 there is no
statistically significant deviation. Hence the condensate
h
χχi is analytic within our numerical accuracy.
If we incr e ase µ further, we obse rve that at some point
a disagreement becomes apparent. The resulting gap is
more pronounced for small β. It does not appear sud-
denly, but develops smoothly and is visible for all non-
large β.
F. Consistency conditions
We checked the consistency conditions of Sec.
III E for
the dens ity hni, the condensate hχχi, and O(t = 1, k)
for different sets of parameters at δ = 10
−3
. Before we
begin with the interpretation of the consistency condi-
tions, we point out again the problem of using the s tan-
dard deviation to estimate error bounds. The res ult-
ing error is much larger than the actual observed sta-
tistical fluctuations and typically we obtain r e sults like
e.g. Re hLO(1, 1)i = 0.015 2 ± 6.1221 for N
t
= 4, I = 0,
N = β = m = µ = 1. Here N
t
denotes the temporal
extension of the lattice, N the number of dege nerated
flavors, I the imaginary noise, β the inverse coupling
constant, m the mass and µ the chemical potential. The
overestimation of the er ror makes a meaningful interpre-
tation of the conditions impossible. Hence we estimate
the error with a bootstra p analysis.
We begin with the consistency conditions of the density
and the co ndensate. T he evaluation is extremely noisy
and only becomes slightly mor e stable for large values
of β. Without loss of g e nerality we restrict ourselves to
10
100
1000
10000
-3 -2 -1 0 1 2 3
Counts
Re Ln
Full theory - β = 1, m = 1, µ = 1, N
t
= 4
Samples
(a) Histogram of hLni.
100
1000
10000
100000
-3 -2 -1 0 1 2 3
Counts
Re Lχ
-
χ
Full theory - β = 1, m = 1, µ = 1, N
t
= 4
Samples
(b) Histogram of hLχχi.
10
100
1000
10000
-3 -2 -1 0 1 2 3
Counts
Re Le
iA
1
Full theory - β = 1, m = 1, µ = 1, N
t
= 4
Samples
(c) Histogram of hLO (1, 1)i.
Figure 11: Histograms for different consistency
conditions.
100
1000
10000
0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5
Counts
Im A
-
Full theory - β = 1, m = 1, µ = 1, N
t
= 4
Samples
Figure 12: Distribution of
Im A
.
the case of t = 1. As a n example we quote Re hLni =
(143 ± 305) ×10
3
and Re hLχχi = (−52 ± 167) ×10
3
for
N
t
= 4, I = 0, N = β = m = µ = k = 1. The
resulting expressions seem to be numerically ill-defined
and despite large sample sizes, we are unable to draw any
conclusions about the validity of the conditions. In our
case both obser vables proved to be unsuitable to check
for the correctness o f the complex Langevin evolution.
The evaluation of the conditions for the observable
hO(t, k)i on the other hand is stable for all checked pa-
rameter sets. If we take the erro r bounds seriously, we
have to conclude that some conditions are not compat-
ible with a vanishing value and are v iolated. We quote
here Re hLO(t, k)i = −0.1395±0.0110 for N
t
= 4, I = 0,
N = β = m = k = 1, µ = 3 as a typical example of a
violated condition and Re hLO(t, k )i = 0.0061 ± 0.0208
for N
t
= 4, I = 0, N = β = m = µ = k = 1 as one which
is compatible with a vanishing va lue. We als o checked
the conditions for k = 2, 3, which turned out to be more
noisy compared to the k = 1 case. We interpret the vio-
lated conditions as an indicator for incorre c t convergence
of the complex Langevin evolution.
When considering the distributions of the evaluated
consistency conditions, we can gain further insights. In
Fig.
11 we find histograms of Ln, Lχχ, and LO(1, 1).
In general the distributions are asymmetrical a nd non-
Gaussian. I n the case of the density and condensate we
saw that the distribution falls off extremely slowly and we
observed values with an absolute value of up to O
10
9
,
resulting in a very noisy evaluation. In contrast the dis-
tribution of LO(1, 1) falls off rapidly.
It is a lso interesting to check the distribution of the
imaginary part of the auxiliary field itself. In the
derivation of the consistency condition, one ass umes that
boundary terms of the field would vanish. It is then im-
portant to check that the imaginary part falls off rapidly
enough. I n Fig.
12 we find a typical histogram of the
average imaginary part Im
A = N
−1
t
P
t
Im A
t
with 10
5
samples. We see that the resulting distribution appears
11
Figure 13: Eigenvalues of the fermion matrix.
to be compatible with our hypothesis.
G. Eigenvalues of the fermion matrix
In Figure
13 we find a scatter plot of the eigenvalues of
the fermion matrix K defined in (7). We sampled approx-
imately 6 × 10
4
eigenvalues during a complex Langevin
evolution. Eigenvalues come in point-r e flected pairs at
point (m, 0), where m denotes the mass of the fermion.
For a vanishing chemical p otential µ = 0 a ll e igenvalues
lie on the line Re λ = m.
H. Imaginary noise
Assuming correct converg e nce of the complex La ngevin
evolution, observables should turn out indepe ndent of the
imaginary noise term controlled by I. Howe ver, in Fig.
14
we a c tually obser ve such a dependence. This was ob-
served previously in other systems too [3 0]. Attempts
to fine-tune I so that the complex Langevin evolution
is deformed and correctly reproduces analytical results
did not succeed. In almost all cases an imaginary noise
I > 0 ca used a more severe disagreement between nu-
merical and analytical r e sults.
I. Periodic continuation of the ima ginary part
A nece ssary condition for the correctness of the com-
plex Langevin evolution is that the distributions of the
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0 0.5 1 1.5 2
Fermion condensate <χ
-
χ>
Chemical potential µ
Full theory - β = 3, m = 1, N
t
= 4
Theory
I = 0.0, CLE
I = 0.1, CLE
I = 0.3, CLE
I = 1.0, CLE
Figure 14: Condensate hχχi with imaginary noise I.
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0 0.5 1 1.5 2
Fermion condensate <χ
-
χ>
Chemical potential µ
Full theory - β = 3, m = 1, N
t
= 4
Theory
Λ = ∞, CLE
Λ = 1.0, CLE
Λ = 0.2, CLE
Λ = 0.1, CLE
Figure 15: Condensate hχχi with periodic cutoff Λ.
imaginary parts of the fields have to fall off rapidly
enough. Although in the generalized Thirring model this
seems to be the case, we can restrict the imaginary part
of the field to the interva l [−Λ, +Λ] a nd make a perio dic
continuation. A similar a pproach was employed in [
30]
(see also [29]), where a fine-tuning of the cutoff enabled
simple toy models to repr oduce correc t results. However,
in Fig.
15 we see that every finite value of Λ widens the
gap to the correct results. Hence in this simple form this
approach cannot be applied to the generalized Thirring
model.
V. CONCLUSIONS
In this paper we applied a complex Langevin evolution
to a generalized Thirr ing model in 0 + 1 dimensions in
order to deal with the resulting sign problem fo r µ >
0. For intermediate values of the chemical potential we
found a gap between analytical and numerical results,
which size depends on the inverse coupling β. While for
12
small β the discr e pancy is large, we observe agreement
for large β. In particular, for β & O(100) we are usually
able to reproduce ana lytical results with high accuracy.
Furthermore for small µ & 0 we did not observe any
significant deviations to theoretical predictions and the
fermion condensate is analytic at µ
2
= 0.
However, in the case o f more than one flavor we ob-
served a qualitative disagreement. O ur approach seems
to be unable to reproduce the plateaus we found for
certain ranges of β. Another interesting observation is
the violation of se veral consistency conditions, indicating
that the complex Langevin evolution might not converge
correctly in general. Attempts to force correct conver-
gence by a n ad hoc fine-tuning of a periodic cuto ff Λ or
an imaginary noise term I did not succeed.
In a subsequent paper, we will pres ent our findings for
the generalized Thirring model in 2 + 1 dimensions [
33].
Further investigations have to deal with the question of
how to address the aforementioned problems. In par-
ticular, coordinate transforma tions as suggested in [
58]
and gaug e coo ling pr oc e dures like the one employed in
[
59] might allow a stabilization of the complex Langevin
evolution.
ACKNOWLEDGMENTS
We thank I.-O. Stamatescu fo r uncountable discussions
and useful remar ks on the manuscript. We also thank
C. Gattringer for a proo f of (
21) by systematic satura-
tion of the Grassmann integral and V. Kasper for a proof
using the determinant identities in [
54]. Furthermore we
acknowledge E. Seiler’s derivation of the plateau struc-
tures as presented in the Appendix. Finally, we thank
G. Aarts and D. Sexty for discussions. This work is sup-
ported by the Helmholtz Alliance HA216/EMMI and by
ERC-AdG-290623 . C. Z. thanks the German National
Academic Foundation for financial support.
APPENDIX: PLATEAUS
As previously discussed, in the case of N > 1 flavors we
observe intermediate plateaus in the considered observ-
ables. In the nonlinear O(2) sigma model, the authors of
[
60] found similar structures. They explained this finite
size behavior with the crossing of energy levels at finite
chemical potential. However, here we re produce them
with the help of a continuum model, which corresponds
to the 0 + 1 dimensional generalized Thirring model. To
this end we consider the Hamiltonian
H =
N
f
X
i=1
m
i
(n
i
+
n
i
) + g
2
N
f
X
i=1
Q
i
2
(35)
with flavor index i = 1, . . . , N
f
and bare masses m
i
. We
introduced Q
i
≡ n
i
−
n
i
, the particle number operator
0
0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
1.8
2
0 1 2 3 4 5
Fermion density <n>
Chemical potential µ
Toy model - β = 10, g = 1, m
1
= 1.0, m
2
= 0.9
<n
total
>
<n
1
>
Figure 16: Densities in the toy model.
n
i
= a
†
i
a
i
and the antiparticle number operator n
i
= b
†
i
b
i
,
which a re given in terms of ladder operato rs. From the
grand canonical partition function
Z = Tr exp
"
−β
T
H −
X
i
µ
i
Q
i
!#
(36)
with β
T
= T
−1
and temperature T , we c an derive the
fermion densities
hni
i
=
1
β
T
∂ log Z
∂µ
i
, hni
total
=
X
i
hni
i
. (37)
The plateaus a re a res ult of the competition between the
Q
i
and Q
2
i
terms. A typical plot for N
f
= 2 can be found
in Fig .
16. Like in the Thirring mo del they app e ar for
intermediate values of g
2
and eventually disappear in the
limit of small and large values.
[1] P. de Forcrand, PoS LAT2009, 010 (2009),
arXiv:1005.0539 [hep-lat].
[2] M. P. Lombardo, Mod.Phys.Lett. A22, 457 (2007),
arXiv:hep-lat/0509180 [hep-lat].
[3] F. Karsch, B.-J. S chaefer, M. Wagner, and J. Wambach,
Phys.Lett. B698, 256 (2011), arXiv:1009.5211 [hep-ph].
[4] M. P. Lombardo, PoS CPOD2006, 003 (2006),
arXiv:hep-lat/0612017 [hep-lat].
[5] A. Schmidt, Y. D. Mercado, and C. Gattringer, PoS
LATTICE2012, 098 (2012),
arXiv:1211.1573 [hep-lat].
[6] Y. D. Mercado, C. Gattringer, and A. Schmidt, Com-
put.Phys.Commun. 184, 1535 (2013)
, arXiv:1211.3436
13
[hep-lat].
[7] A. Alexandru, M. Faber, I. Horvath, and K.-F. Liu,
Phys.Rev. D72, 114513 (2005), arXiv:hep-lat/0507020
[hep-lat]
.
[8] A. Alexandru and U. Wenger, Phys.Rev. D83, 034502
(2011)
, arXiv:1009.2197 [hep-lat].
[9] G. Parisi and Y. Wu, Sci.Sin. 24, 483 (1981).
[10] P. H. Damgaard and H. Huffel, Phys.Rept. 152, 227
(1987)
.
[11] G. Parisi, Phys.Lett. B131, 393 (1983).
[12] F. Karsch and H. W. Wyld, Phys.Rev.Lett. 55, 2242
(1985)
.
[13] G. Aarts and F. A. James, JHEP 1201, 118 (2012),
arXiv:1112.4655 [hep-lat].
[14] N. Bilic, H. Gausterer, and S. Sanielevici, Phys.Rev.
D37, 3684 (1988)
.
[15] G. Aarts and I.-O. Stamatescu, JHEP 0809, 018 (2008),
arXiv:0807.1597 [hep-lat].
[16] G. Aarts and K. Splittorff, JHEP 1008, 017 (2010),
arXiv:1006.0332 [hep-lat].
[17] G. Aarts, Phys.Rev.Lett. 102, 131601 (2009),
arXiv:0810.2089 [hep-lat].
[18] G. Aarts, JHEP 0905, 052 (2009), arXiv:0902.4686 [hep-
-lat]
.
[19] J. Berges, S. Borsanyi, D. Sexty, and I .- O. Stamatescu,
Phys.Rev. D75, 045007 (2007), arXiv:hep-lat/0609058
[hep-lat]
.
[20] J. Berges and D. Sexty, Nucl.Phys. B799, 306 (2008),
arXiv:0708.0779 [hep-lat].
[21] J. Berges and I.-O. Stamatescu, Phys.Rev.Lett. 95,
202003 (2005)
, arXiv:hep-lat/0508030 [hep-lat].
[22] G. Aarts and F. A. James, JHEP 1008, 020 (2010),
arXiv:1005.3468 [hep-lat].
[23] J. Ambjorn, M. Flensburg, and C. Peterson, Nucl.Phys.
B275, 375 (1986)
.
[24] H. W. Hamber and H. Ren, Phys.Lett. B159, 330 (1985).
[25] J. Flower, S. W. Otto, and S. Callahan, Phys.Rev. D34,
598 (1986)
.
[26] E.-M. Ilgenfritz, Phys.Lett. B181, 327 (1986).
[27] G. Aarts, PoS LATTICE2012, 017 (2012),
arXiv:1302.3028 [hep-lat].
[28] G. Aarts, PoS LAT2009, 024 (2009), arXiv:0910.3772
[hep-lat]
.
[29] G. Aarts, E. Seiler, and I.-O. Stamatescu, Phys.Rev.
D81, 054508 (2010)
, arXiv:0912.3360 [hep-lat].
[30] G. Aarts, F. A. James, E. Seiler, and I.- O. S tamatescu,
Eur.Phys.J. C71, 1756 (2011), arXiv:1101.3270 [hep-lat].
[31] G. Aarts, F. A. James, E. Seiler, and I.- O. S tamatescu,
PoS LATTICE2011, 197 (2011),
arXiv:1110.5749 [hep-
-lat]
.
[32] D. S pielmann, Aspects of confinement in QCD from
lattice simul ations, Ph.D. thesis, Ruprecht–Karls–
Universität Heidelberg (2010).
[33] J. M. Pawlowski and C. Zielinski, Phys.Rev. D87, 094509
(2013)
, arXiv:1302.2249 [hep-lat].
[34] H. Gies and L. Janssen, Phys.Rev. D82, 085018 (2010),
arXiv:1006.3747 [hep-th].
[35] S. Chandrasekharan and A. Li, Phys.Rev.Lett. 108,
140404 (2012)
, arXiv:1111.7204 [h ep-lat].
[36] W. E. Thirring, Annals Phys. 3, 91 (1958).
[37] S. Christofi, S. Hands, and C. Strouthos, Phys.Rev.
D75, 101701 (2007)
, arXiv:hep-lat/0701016 [hep-lat].
[38] S. R. Coleman, Phys.Rev. D11, 2088 (1975).
[39] G. Benfatto, P. Falco, and V. Mastropietro, Com-
mun.Math.Phys. 285, 713 (2009)
, arXiv:0711.5010 [hep-
-th]
.
[40] T. Itoh, Y. Kim, M. Sugiura, and K. Yamawaki,
Prog.Theor.Phys. 93, 417 (1995), arXiv:hep-th/9411201
[hep-th]
.
[41] T. D. Cohen, Phys.Rev.Lett. 91, 222001 (2003),
arXiv:hep-ph/0307089 [hep-ph].
[42] K. S plittorff and J. J. M. Verbaarschot, Phys.Rev.Lett.
98, 031601 (2007)
, arXiv:hep-lat/0609076 [hep-lat].
[43] J. B. Kogut and L. Susskind , Phys.Rev. D11, 395 (1975).
[44] T. Banks, L. Susskind, and J. B. Kogut, Phys.Rev. D13,
1043 (1976)
.
[45] T. Banks et al. (Cornell-Oxford-Tel Aviv-Yeshiva Collab-
oration),
Phys.Rev. D15, 1111 (1977).
[46] L. Susskind, Phys.Rev. D16, 3031 (1977).
[47] P. Hasenfratz and F. K arsch, Phys.Lett. B125, 308
(1983)
.
[48] L. Del Debbio and S . Hands, Phys.Lett. B373, 171
(1996)
, arXiv:hep-lat/9512013 [hep-lat].
[49] J. Han and M. A. Stephanov, Phys.Rev. D78, 054507
(2008)
, arXiv:0805.1939 [h ep-lat].
[50] J. O. Andersen, L . T. Kyllingstad, and K. Splittorff,
JHEP 1001, 055 (2010), arXiv:0909.2771 [hep-lat].
[51] C. C. Chang, Math. Comp. 49, 523 (1987).
[52] J. Ambjorn and S. K. Yang, Phys.Lett. B165, 140
(1985)
.
[53] G. Aarts, F. A. James, E. Seiler, and I.-O. Stamatescu,
Phys.Lett. B687, 154 (2010), arXiv:0912.0617 [hep-lat].
[54] L. G. Molinari, Linear Algebra and its Applications 429,
2221 (2008)
.
[55] T. DeGrand and C. E. Detar, Lattice m ethods f or quan-
tum chromodynamics (World Scientific, 2006).
[56] W. J. Morokoff and R. E. Caflisch , Journal of Computa-
tional Physics 122, 218 (1995)
.
[57] G. Aarts and F. A. James, PoS LATTICE2010, 321
(2010),
arXiv:1009.5838 [hep-lat].
[58] G. Aarts, F. A. James, J. M. Pawlowski, E. Seiler,
D. Sexty, et al.,
JHEP 1303, 073 (2013), arXiv:1212.5231
[hep-lat]
.
[59] E. Seiler, D. Sexty, and I.-O. Stamatescu, (2012),
arXiv:1211.3709 [hep-lat].
[60] D. Banerjee and S. Chandrasekharan, Phys.Rev. D81,
125007 (2010)
, arXiv:1001.3648 [h ep-lat].
