text:
# On systematic errors due to quenching in 𝜀'/𝜀

## Abstract

Recently, we pointed out that chiral transformation properties of penguin operators change in the transition from unquenched to (partially) quenched QCD. As a consequence, new operators appear in (partially) quenched QCD penguins, introducing ambiguities which should be considered a quenching artifact. Here we discuss more specifically the effects of this phenomenon on the quenched $`\mathrm{\Delta }I=1/2`$ $`K\pi \pi `$ amplitude, and in particular, its potential numerical effect on recent lattice estimates for $`\epsilon ^{}/\epsilon `$.

The strong $`LR`$ penguin operators $`Q_{5,6}`$ mediating $`\mathrm{\Delta }S=1`$, $`\mathrm{\Delta }I=1/2`$ weak transitions are of the form

$$Q_{penguin}^{QCD}(\overline{s}\gamma _\mu (1\gamma _5)d)(\overline{q}\gamma _\mu (1+\gamma _5)q),$$
(1)
with $`q`$ summed over $`u`$, $`d`$, $`s`$, and thus transform trivially under SU(3)R. However, this is no longer true when one quenches the theory. The quenched theory can be defined by introducing a bosonic “ghost” quark $`\stackrel{~}{q}`$ for each of the three valence quarks $`q`$ , and thus the symmetry group is enlarged to a graded group transforming all six quarks into each other, SU(3)$`{}_{L}{}^{}\times `$SU(3)$`{}_{R}{}^{}`$ SU(3$`|`$3)$`{}_{L}{}^{}\times `$SU(3$`|`$3)R , under which $`Q_{penguin}^{QCD}`$ is no longer a singlet. This was observed in ref. , in which also the consequences for $`K0`$ and $`K\pi `$ matrix elements were studied in ChPT, both in the fully and partially quenched cases. Here we extend the discussion to $`K\pi \pi `$ amplitudes, restricting ourselves to the quenched three-flavor theory. The partially quenched case, as well as more general results, are contained in ref. .

The right-handed current in eq. (1) can be split into a singlet ($`S`$) and a non-singlet ($`NS`$) part, $`(\overline{q}q)_R\frac{1}{2}(\overline{\psi }\psi )_R+\frac{1}{2}(\overline{\psi }\widehat{N}\psi )_R`$ (in an obvious shorthand), where $`\psi `$ runs over valence and ghost quarks, and $`\widehat{N}=\mathrm{diag}(1,1,1,1,1,1)`$. Accordingly, $`Q_{penguin}^{QCD}`$ splits up into a singlet and a non-singlet part $`Q^{QS}`$, resp. $`Q^{QNS}`$, which are represented in leading-order ChPT by

$`Q^{QS}\alpha _{q1}^{(8,1)}\mathrm{str}(\mathrm{\Lambda }L_\mu L_\mu )+\alpha _{q2}^{(8,1)}\mathrm{str}(\mathrm{\Lambda }X_+),`$
$`Q^{QNS}f^2\alpha ^{NS}\mathrm{str}(\mathrm{\Lambda }\mathrm{\Sigma }\widehat{N}\mathrm{\Sigma }^{}),`$ (2)
with $`\alpha ^{NS}`$ a new low-energy constant (LEC), which only exists in the quenched theory, and thus should be considered an artifact of the quenched approximation .111See also this reference for notation. While a similar phenomenon also occurs for $`LL`$ penguins, we note that in the $`LR`$ case, the new non-singlet operator is $`O(p^0)`$ in ChPT, thus enhancing its contribution to weak matrix elements. (This is similar to the enhancement of the EM penguins $`Q_{7,8}`$.)

At tree level, the new operator does not contribute to matrix elements with only physical (valence) mesons on the external lines, but it does contribute at one loop, thus competing at $`O(p^2)`$ with the tree-level contributions from $`Q^{QS}`$. For the $`K\pi `$ (with degenerate masses) and $`K0`$ matrix elements used recently by the CP-PACS  and RBC  collaborations to extract the value of $`\alpha _{q1}^{(8,1)}`$ we find

$`\pi ^+|Q_{penguin}^{QCD}|K^+=`$ (3)
$`{\displaystyle \frac{4M^2}{f^2}}\left[(\alpha _{q1}^{(8,1)}\widehat{\beta }_{q1}^{NS}{\displaystyle \frac{1}{2}}\widehat{\beta }_{q2}^{NS})(\alpha _{q2}^{(8,1)}+\widehat{\beta }_{q3}^{NS})\right],`$
$`0|Q_{penguin}^{QCD}|K^0=`$ (4)
$`{\displaystyle \frac{4i}{f}}\left[(M_K^2M_\pi ^2)(\alpha _{q2}^{(8,1)}+\widehat{\beta }_{q3}^{NS})+\alpha ^{NS}\mathrm{logs}\right],`$
and for the $`K\pi \pi `$ matrix element with physical kinematics we find

$`\pi ^+\pi ^{}|Q_{penguin}^{QCD}|K^0={\displaystyle \frac{4i}{f^3}}[`$ (5)
$`(M_K^2M_\pi ^2)(\alpha _{q1}^{(8,1)}\widehat{\beta }_{q1}^{NS}{\displaystyle \frac{1}{2}}\widehat{\beta }_{q2}^{NS})+\alpha ^{NS}\mathrm{logs}].`$
Here $`\widehat{\beta }_{qi}^{NS}(4\pi )^2\beta _{qi}^{NS}`$ are $`O(p^2)`$ LECs which appear at next-to-leading order in the non-singlet sector, and absorb the scale dependence from the chiral logs . Note that the singlet and non-singlet LECs always appear in the same linear combinations. This can be understood from the fact that only physical mesons appear on the external lines (see ref.  for details).

From these results, we conclude that there are at least three different strategies one could follow to estimate the real-world $`\mathrm{\Delta }I=1/2`$ decay rate from a quenched computation:
$``$ Ignore $`\alpha ^{NS}`$, but not the $`\widehat{\beta }_{qi}^{NS}`$. This introduces a scale dependence, but makes sense if $`\alpha ^{NS}\mathrm{logs}`$ is small in the physical matrix element, at a reasonable scale. Thus $`\alpha _{q1}^{(8,1)}\widehat{\beta }_{q1}^{NS}\frac{1}{2}\widehat{\beta }_{q2}^{NS}`$ from $`K\pi `$ and $`K0`$ is taken as the “best” estimate of unquenched $`\alpha _1^{(8,1)}`$, and the physical $`K\pi \pi `$ matrix element is calculated from the unquenched expression $`4i(M_K^2M_\pi ^2)\alpha _1^{(8,1)}/f^3`$. This is the strategy followed by CP-PACS and RBC .
$``$ Drop all the non-singlet operators, following the simple prescription given in ref. . This was investigated for $`Q_6`$ in ref. ; see below.
$``$ Keep everything in eqs. (3-4), and also calculate the physical $`K\pi \pi `$ matrix element using eq. (5), assuming that it gives the best estimate of the real-world matrix element.

Naively, the second strategy might be considered the obvious one, since the aim is to compute the (unquenched) value of $`\alpha _1^{(8,1)}`$, which (at least at tree level in unquenched ChPT) determines the penguin contribution to the physical decay rate. However, we do not know whether the quenched (e.g. $`\alpha _{q1}^{(8,1)}`$) and unquenched (e.g. $`\alpha _1^{(8,1)}`$) LECs have the same values; it might happen that $`\alpha _{q1}^{(8,1)}\widehat{\beta }_{q1}^{NS}\frac{1}{2}\widehat{\beta }_{q2}^{NS}`$ is a better estimate of $`\alpha _1^{(8,1)}`$. For this reason, we believe that the spread in values obtained by following different strategies should be taken as a systematic error due to quenching. The origin of this uncertainty comes from the fact that in the case of electro-weak operators, defining the quenched theory as QCD with the fermion determinant set equal to a constant is not enough: in addition, a prescription has to be given for the “embedding” of electro-weak operators in the quenched theory . In some cases this is straightforward, but in the case at hand it is not. Note that in the case of partially quenched QCD with three dynamical light flavors, the choice of strategy is unambiguous: the second strategy should be followed , since the LECs of this theory are those of unquenched QCD .

We will now comment on the three different strategies listed above. First, we address the issue of ignoring $`\alpha ^{NS}`$. We found that if $`\alpha ^{NS}`$, with the normalization defined in eq. (2), is of the same order as $`\alpha _{q1}^{(8,1)}`$, its numerical contribution to the $`K\pi \pi `$ matrix element in eq. (5) is small ($`10^2`$ of the total), for a reasonable choice of the ChPT scale . If so, this justifies ignoring $`\alpha ^{NS}`$. The reason appears to be purely numerical, and can be understood as a suppression by the typical factor $`1/(4\pi )^2`$ arising at one loop.

It would of course be better to determine $`\alpha ^{NS}`$ from a lattice computation. In principle, it can be determined from existing data for the $`K0`$ matrix element (cf. eq. (4)), but this may be difficult to do in practice, since it is hard to disentangle logs from the linear term in the quark mass at the typical values of quark masses employed in current simulations . In fact, it is easier to determine $`\alpha ^{NS}`$ from a process in which it appears at leading order in ChPT. This can be done as follows. First, rotate $`Q_{penguin}^{QCD}`$ by a SU(3$`|`$3)L rotation into $`(\overline{s}\gamma _\mu (1\gamma _5)\stackrel{~}{d})(\overline{\psi }\widehat{N}\gamma _\mu (1+\gamma _5)\psi )`$, which is in the same representation of SU(3$`|`$3)$`{}_{L}{}^{}\times `$SU(3$`|`$3)R. Next, consider the matrix element between a fermionic kaon $`\stackrel{~}{K}\overline{\stackrel{~}{d}}\gamma _5s`$ and the vacuum. ChPT tells us that

$$0|(\overline{s}\stackrel{~}{d})_L(\overline{\psi }\widehat{N}\psi )_R|\stackrel{~}{K}=2if^2\alpha ^{NS}+O(p^2).$$
(6)
The key observation is that no simulations with ghost quarks are needed. After carrying out all Wick contractions on the left-hand side of eq. (6), one simply uses the fact that ghost and valence propagators are equal, $`\stackrel{~}{d}(x)\overline{\stackrel{~}{d}}(y)=d(x)\overline{d}(y)`$, etc. Estimating $`\alpha ^{NS}`$ is thus as simple as estimating $`\alpha _{2q}^{(8,1)}`$, once the usual valence-quark propagators have been computed. It is important to determine $`\alpha ^{NS}`$ in order to see whether indeed its contribution can be considered a small effect.

The remaining choice of strategy boils down to whether one considers $`\alpha _{1q}^{(8,1)}`$ or $`\alpha _{1q}^{(8,1)}\widehat{\beta }_{q1}^{NS}\frac{1}{2}\widehat{\beta }_{q2}^{NS}`$ to be a better estimate of unquenched $`\alpha _1^{(8,1)}`$. Since the issue only arises for QCD penguins, and is a non-leading effect in ChPT for the $`LL`$ case , it is unlikely to have much influence on the real part of the $`\mathrm{\Delta }I=1/2`$ amplitude. This is not the case for $`\epsilon ^{}/\epsilon `$, where $`Q_6`$ is expected to be one of the largest contributions. In ref.  it was found, following the simple prescription for dropping non-singlet operators given in ref.  (and using staggered fermions, while refs.  used domain-wall fermions), that leaving out the non-singlet operator in $`Q_6`$ enhances $`B_6^{(1/2)}`$ by a factor two at the kaon mass. To leading order in ChPT, this means that $`\alpha _{1q}^{(8,1)}`$ is twice as large as $`\alpha _{1q}^{(8,1)}\widehat{\beta }_{q1}^{NS}\frac{1}{2}\widehat{\beta }_{q2}^{NS}`$, potentially enhancing the contribution of $`Q_6`$ to $`\epsilon ^{}/\epsilon `$ by a factor two relative to the contributions found in refs. .

No estimate of this enhancement with domain-wall fermions is available. In order to get an idea of how large the effect on $`\epsilon ^{}/\epsilon `$ can be, we therefore consider what a factor-two enhancement of $`B_6^{(1/2)}`$ would imply for the central value of $`\epsilon ^{}/\epsilon =4\times 10^4`$ reported by RBC . We find that this value would shift to $`+10\times 10^4`$. Considering the CP-PACS results  the situation is less clear. They typically find a smaller (in absolute value) negative contribution from $`Q_6`$, and there appears to be a stronger dependence on quark mass. When we take their value of $`2\times 10^4`$ at a degenerate Goldstone-boson mass of about 600 MeV, we find that enhancing $`B_6^{(1/2)}`$ by a factor two changes this into $`+2\times 10^4`$.

We emphasize that these are rather unreliable estimates of the effect of keeping or dropping the non-singlet operators in $`Q_6`$. Many other effects are not under control, such as the fact that all computations were done at only one value of the lattice spacing, and that our estimates rely on leading-order ChPT, to mention but a few. However, these estimates do make it clear that the issue is important, not only theoretically, but also phenomenologically. While the quenched approximation does remarkably well in “simple” strong interaction physics, such as the light hadron spectrum, it is a major obstacle to obtaining phenomenologically relevant estimates of finely-tuned quantities like $`\epsilon ^{}/\epsilon `$. Partially quenched computations with three dynamical light flavors are needed in order to resolve the ambiguities introduced by quenching.

Finally, we recall that similar ambiguities also affect $`Q_5`$ and strong $`LL`$ penguins. However, these are expected to be numerically less important. For $`Q_5`$ the effect occurs at leading order in ChPT, just as for $`Q_6`$, but the corresponding Wilson coefficient in the weak $`\mathrm{\Delta }S=1`$ hamiltonian is much smaller, suppressing the contribution of $`Q_5`$ altogether. For $`LL`$ penguin operators, the effect is non-leading in ChPT, and affects $`K\pi \pi `$ matrix elements only at $`O(p^4)`$.

Acknowledgements. We would like to thank N. Christ, R. Mawhinney, M. Savage and S. Sharpe for discussions. EP thanks the Dept. of Physics and Astronomy of San Francisco State University for hospitality. MG is supported in part by the US Dept. of Energy; EP was supported in part by the Italian MURST under the program Fenomenologia delle Interazioni Fondamentali.


format: text

text:
# Gauge Problem of Monopole Dynamics in SU(2) Lattice Gauge Theory

## 1 Introduction

It is important to understand color confinement mechanism in Quantum Chromodynamics (QCD). Many numerical simulations have been done and they support the dual superconductor scenario of the QCD vacuum as a confinement mechanism . Magnetic monopoles are induced by performing an abelian projection , i.e., a partial gauge-fixing that keeps $`U(1)U(1)`$. It is known that the string tension calculated from the abelian and the monopole parts reproduces well the original one when we perform an abelian projection in Maximally Abelian (MA) gauge where link variables are abelianized as much as possible. In addition to the string tension, many low energy physical properties of QCD are reproduced from the abelian and the monopole degrees of freedom alone. It is called as “abelian and monopole dominance”. These facts suggest that monopoles play an important role for the confinement mechanism. Actually, a low energy effective theory which is described in terms of monopole currents has been derived by Shiba and Suzuki  and an almost perfect monopole action showing the scaling behavior is derived by Chernodub et al . Monopole condensation occurs due to energy-entropy balance . Abelian color-electric flux is squeezed into a string-like shape  by the superconducting monopole current. This squeezed color flux causes a confinement potential between quarks.

We note that we have infinite degrees of freedom when we perform an abelian projection. That is to say, which gauge should be chosen? Recently Laplacian Abelian (LA) gauge was proposed and it appears to have similar good properties . Actually MA and LA gauges are very similar. Are MA and LA gauges exceptional? If such is the case, there must exist any reason to justify it, although it seems very difficult to find such a reason. Another interpretation is that monopole dynamics does not depend on the choice of gauge in the continuum limit, although it seems dependent of the gauge choice at the present stage of lattice study. In other words, MA gauge and LA gauge are considered to have a wider window even at present to see the continuum limit than other gauges have.

Our aim of this paper is to show first that MA gauge is not a special choice of the good gauge for color confinement. We restrict ourselves to pure $`SU(2)`$ QCD for simplicity. Here we discuss two new gauges in addition to LA gauge. They have a different continuum limit but they all can reproduce well the $`SU(2)`$ string tension. The second is to derive an effective monopole action and to study the blockspin transformation of the monopole currents in various abelian projections. If their renormalization group (RG) flows converge onto the same line with a finite number of blockspin transformations, we can expect gauge independence of monopole dynamics in the infrared region. The paper is organized as follows: In Section 2, we present some theoretical and phenomenological arguments which support gauge independence of abelian and monopole dominance. In Section 3, we describe gauge fixing procedures being used. In Section 4, we show that the $`SU(2)`$ string tension is well-reproduced from abelian or monopole degrees of freedom alone in four different abelian projection schemes. In section 5, we present our results from RG flow study of effective monopole actions in various abelian projections. In Section 6, we summarize our conclusions.

## 2 Theoretical and phenomenological background

### 2.1 Gauge fixings and abelian dominance

It is known that the abelian Wilson loop reproduces well the $`SU(2)`$ string tension numerically, if MA gauge or LA gauge is applied . In the case of Polyakov gauge, the string tension which is calculated from abelian Polyakov loop correlators is exactly the same as that of $`SU(2)`$ . Shoji et al. have developed a stochastic gauge fixing method which interpolates between the MA gauge and no gauge fixing . They have found that abelian dominance for the heavy quark potential is realized even in a gauge which is far from the MA gauge. In a finite temperature system, abelian Polyakov loops in various gauges reproduce the phase transition behavior of $`SU(2)`$ Polyakov loop . See, Figure 1.

Abelian dominance is shown also analytically. Abelian Wilson loops constructed without any gauge fixing give the same string tension as that of $`SU(2)`$ Wilson loops in the strong coupling expansion . The same fact for any coupling region has been proved by Ogilvie using the character expansion . An abelian Wilson loop operator is given by

$$W_A\left[C\right]=\frac{1}{2}\text{Tr}\left[\underset{s,\mu C}{}u_\mu (s)\right],$$

where $`u_\mu `$ is an abelian projected $`U(1)`$ link variable. Since $`W_A`$ is not a class function of $`SU(2)`$ group, only the $`SU(2)`$ invariant part extracted from $`W_A`$ is non-vanishing in the expectation value. This can be written as

$$W_A^{inv}=\frac{1}{2}𝒟g\text{Tr}\left[\underset{s,\mu C}{}g(s)u_\mu (s)g^{}(s+\widehat{\mu })\right].$$

Using a character expansion, we get an expression for the expectation value of the abelian Wilson loop in terms of $`SU(2)`$ Wilson loops:

$$W_A^{inv}=\left(\frac{2}{3}\right)^{P(C)}W_{SU(2)}_{1/2}+\text{(half integer higher rep.).}$$

Since the lowest representation is dominant, we can show that the $`SU(2)`$ string tension $`\sigma _{SU(2)}`$ can be reproduced perfectly from the abelian string tension $`\sigma _A`$:

$$\sigma _A=\underset{I,J\mathrm{}}{lim}\mathrm{ln}\frac{W_A(I+1,J+1)W_A(I,J)}{W_A(I+1,J)W_A(I,J+1)}=\sigma _{SU(2)}.$$

Furthermore, Ogilvie has shown that similar arguments hold even with a gauge fixing function

$$S_{gf}=\lambda \text{Tr}\left[U_\mu (s)\sigma _3U_\mu ^{}(s)\sigma _3\right],$$

if the gauge parameter $`\lambda `$ is small enough.

### 2.2 Monopole dominance

There are numerical results supporting monopole dominance. $`SU(2)`$ string tension is well reproduced only from the monopole part of abelian Wilson loops in MA gauge  and LA gauge . We note also that monopole Polyakov loops in various gauges reproduce the phase transition behavior of $`SU(2)`$ Polyakov loop . See, Figure 2.

In addition to these numerical evidence, we can prove analytically gauge independence of monopole dominance if abelian dominance is gauge independent . If abelian dominance is gauge independent, a common abelian effective action $`S_{eff}`$ written in terms of abelian gauge field surely exists in any gauge and works well in the infrared region as in MA gauge. Since $`S_{eff}`$ takes a form of modified compact QED, an effective monopole action can be derived analytically. One can evaluate the contribution of monopoles to the abelian Wilson loop using this effective monopole action.

In MA gauge, it is known numerically that an effective monopole action composed of two-point self+Coulomb+nearest neighbor interactions is a good approximation in the infrared region. The action can be transformed exactly into a modified compact QED action in the generic Villain form:

$$Z=_\pi ^\pi 𝒟\theta \underset{n}{}\mathrm{exp}\left[\frac{1}{4\pi ^2}(d\theta +2\pi n,\mathrm{\Delta }D(d\theta +2\pi n))+i(J,\theta )\right],$$

where $`D\beta \mathrm{\Delta }^1+\alpha +\gamma \mathrm{\Delta }`$. An expectation value of the abelian Wilson loop $`W=e^{i(\theta ,J)}`$ can be estimated using this action, where $`J`$ is the color electric current which takes a value $`\pm 1`$ on a closed loop. When we use the BKT transformation , we get an expectation value of the abelian Wilson loop in terms of monopole currents $`k`$:

$$W=\frac{1}{Z}\underset{k,dk=0}{}\mathrm{exp}\left[(k,Dk)2\pi i(k,\delta \mathrm{\Delta }^1M)\pi ^2(J,(\mathrm{\Delta }^2D)^1J)\right],$$

where $`M`$ takes a value $`\pm 1`$ on a surface whose boundary is $`J`$ ($`J=\delta M`$). Electric-electric current ($`J`$-$`J`$) interactions are of a modified Coulomb interaction and have no line singularity leading to a linear potential. A linear potential of the abelian Wilson loop originates from the second term of the monopole contribution. Gauge independence of monopole dominance is derived from that of abelian dominance. Gauge independence of an order parameter is also observed in Ref. .

### 2.3 An objection to gauge independence

As we have shown in previous subsections, there is encouraging evidence which supports gauge independence of the confinement scenario in terms of monopoles. On the other hand, there is a strong objection to the idea of gauge independence.

Consider a gauge called Polyakov gauge where Polyakov loop operators are diagonalized in the continuum finite-temperature QCD. It is proved  that the singularities of the gauge fixing run only in the time-like direction. This means that there are only time-like monopoles in the system when the Polyakov gauge is employed, if the degeneration points in abelian projection only correspond to monopoles as ’t Hooft argued. Since such time-like monopoles do not contribute to the physical string tension , monopole dominance is violated.

But numerically the above theoretical expectation seems to be inconsistent with numerical data. We show our preliminary result in Figure 3. The spatial and temporal monopole densities are plotted in Figure 3 as a function of lattice spacing $`a`$ in the unit of physical string tension $`\sqrt{\sigma _p}`$. These densities are defined as

$$\rho _s(\beta )=\frac{\frac{1}{3}_s_{i=1,2,3}|k_i(s)|}{(N_sa)^3N_4},\rho _t(\beta )=\frac{_s|k_4(s)|}{(N_sa)^3N_4},$$

respectively. Figure 3 shows that spatial (lattice) monopole density may take non-zero value even in the $`a0`$ limit. This is not compatible with the theoretical expectation above. In the authors’ opinion, the continuum limit of lattice monopoles must contain extra ingredients different from the expected monopoles corresponding to singularities of Polyakov loop operators. We will give a detailed analysis elsewhere.

## 3 Various abelian projections on a lattice

To check gauge (in)dependence of monopole dynamics, we study the abelian projection in various gauges.

1. MA gauge:
The most known is Maximally Abelian (MA) gauge. It is defined by maximizing the following quantity $`R_{MA}`$:

$$R_{MA}=\text{Tr}\underset{s,\mu }{}U_\mu (s)\sigma _3U_\mu ^{}(s)\sigma _3.$$
(1)
This is achieved by diagonalizing an operator

$$X_{MA}(s)=\underset{\mu }{}\left(U_\mu (s)\sigma _3U_\mu ^{}(s)+U_\mu ^{}(s\widehat{\mu })\sigma _3U_\mu (s\widehat{\mu })\right).$$

That is,

$$X_{MA}(s)X_{MA}^{}(s)=V(s)X_{MA}(s)V^{}(s)=\text{diag}\{\lambda _1,\lambda _2\},$$

where $`V(s)`$ is a gauge transformation matrix. The diagonalization corresponds to the condition,

$$\underset{\mu }{}\left(_\mu iA_\mu ^3\right)A_\mu ^\pm =0,$$
(2)
in the continuum limit.
2. LA gauge :
First consider MA gauge again. To maximize $`R_{MA}`$ in Eq.(1) is to minimize the functional

$`S_{MA}`$ $`=`$ $`{\displaystyle \underset{s,\mu }{}}\left\{1{\displaystyle \frac{1}{2}}\text{Tr}[\mathrm{\Phi }(s)U_\mu (s)\mathrm{\Phi }(s+\widehat{\mu })U_\mu ^{}(s)]\right\}`$ (3)
$`=`$ $`{\displaystyle \underset{s,\mu }{}}\left\{\varphi ^a(s)R_\mu ^{ab}(s)\varphi ^b(s+\widehat{\mu })\right\},`$
where $`R_\mu `$ is a gauge field in the adjoint representation,

$$R_\mu ^{ab}(s)=\frac{1}{2}\text{Tr}\left(\sigma _aU_\mu (s)\sigma _bU_\mu ^{}(s)\right).$$

$`\mathrm{\Phi }`$ is parameterized by a spin variable $`\varphi `$ which satisfies a local constraint

$$\mathrm{\Phi }(s)=V^{}(s)\sigma _3V(s)=\underset{a=1}{\overset{3}{}}\varphi ^a(s)\sigma _a,\underset{a=1}{\overset{3}{}}\left(\varphi ^a(s)\right)^2=1.$$
(4)
Because of the local constraint from the normalization, it is very difficult to find a set of $`\varphi `$ which realizes the absolute minimum of Eq.(3).

The key idea of the LA gauge fixing is to relax this constraint:

$$\underset{a=1}{\overset{3}{}}\left(\varphi ^a(s)\right)^2=1\underset{s}{}\underset{a=1}{\overset{3}{}}\left(\varphi ^a(s)\right)^2=1$$

The functional to minimize becomes

$$S_{LA}=\frac{1}{2}\underset{x,a}{}\underset{y,b}{}\varphi ^a(x)(\mathrm{}_{xy}^{ab})\varphi ^b(y),$$
(5)
where

$$\mathrm{}_{xy}^{ab}=\underset{\mu }{}\left(2\delta _{xy}\delta ^{ab}R_\mu ^{ab}(x)\delta _{y,x+\widehat{\mu }}R_\mu ^{ba}(x\widehat{\mu })\delta _{y,x\widehat{\mu }}\right).$$
(6)
Minimizing Eq.(5) amounts to finding the eigenvector belonging to the lowest eigenvalue of the covariant laplacian operator. This eigenvalue problem can be solved numerically (we used an implicitly restarted Arnoldi method. For example, see Ref. ). The gauge transformation matrix $`V(s)`$ is defined by

$$V^{}(s)\sigma _3V(s)=\underset{a=1}{\overset{3}{}}\widehat{\varphi }^a(s)\sigma _a,$$
(7)
where

$$\varphi ^a(s)=\rho (s)\widehat{\varphi }^a(s),\rho ^2(s)=\underset{a=1}{\overset{3}{}}\left(\varphi ^a(s)\right)^2.$$
(8)
In the continuum limit, LA gauge corresponds to the gauge condition,

$$\underset{\mu }{}\left(_\mu iA_\mu ^3\right)\left(\rho ^2A_\mu ^\pm \right)=0.$$
(9)
3. MAWL gauge :
Maximally Abelian Wilson Loop (MAWL) gauge is a gauge which maximizes a Wilson loop operator written in terms of abelian link variables:

$$W_A=\mathrm{cos}\theta _{\mu \nu }(s),$$
(10)
where $`\theta _{\mu \nu }(s)=\theta _\mu (s)+\theta _\nu (s+\widehat{\mu })\theta _\mu (s+\widehat{\nu })\theta _\nu (s)`$. It is achieved by diagonalizing the following operator

$`X_{MAWL}(s)`$ $`=`$ $`{\displaystyle \underset{\mu \nu }{}}[{\displaystyle \frac{\mathrm{sin}\theta _{\mu \nu }(s)\mathrm{sin}\theta _{\mu \nu }(s\widehat{\nu })}{U_0^2(s,\mu )+U_3^2(s,\mu )}}(U(s,\mu )\sigma _3U^{}(s,\mu ))`$
$`+{\displaystyle \frac{\mathrm{sin}\theta _{\mu \nu }(s\widehat{\mu }\widehat{\nu })\mathrm{sin}\theta _{\mu \nu }(s\widehat{\mu })}{U_0^2(s\widehat{\mu },\mu )+U_3^2(s\widehat{\mu },\mu )}}(U^{}(s\widehat{\mu },\mu )\sigma _3U(s\widehat{\mu },\mu ))].`$
In the continuum limit, we get the following gauge condition:

$$\underset{\mu \nu }{}_\nu f_{\mu \nu }A_\mu ^\pm =0.$$
(11)
4. L-type gauge:
There are infinitely many gauges similar to MA gauge. Here we show one simplest extension called L-type gauge. It is defined by maximizing

$$R_L=\text{Tr}\underset{s,\mu \nu }{}U_\mu (s)\sigma _3U_\nu (s+\widehat{\mu })\sigma _3U_\nu ^{}(s+\widehat{\mu })\sigma _3U_\nu ^{}(s)\sigma _3.$$
(12)
This is given by diagonalizing

$`X_L(s)`$ $`=`$ $`{\displaystyle \underset{\mu \nu }{}}[U_\mu (s)\sigma _3U_\mu ^{}(s)\sigma _3U_\nu (s)\sigma _3U_\nu ^{}(s)`$
$`+U_\mu ^{}(s\widehat{\mu })\sigma _3U_\nu (s\widehat{\mu })\sigma _3U_\nu ^{}(s\widehat{\mu })\sigma _3U_\mu (s\widehat{\mu })].`$
A schematic representation of $`R_L`$ is shown in Figure 4.

In the continuum limit, we get the following gauge condition:

$$\underset{\mu \nu }{}\{(_\mu \pm iagA_\mu ^3)+(_\nu \pm iagA_\nu ^3)\}(A_\mu ^{}+A_\nu ^{})=0.$$
(13)
5. There are various gauges called unitary gauge. Polyakov gauge and F12 gauge are defined with the following operators which are diagonalized:

$`X_{Pol}(s)`$ $`=`$ $`{\displaystyle \underset{i=1}{\overset{N_4}{}}}U_4(s+(i1)\widehat{4}),`$ (14)
$`X_{F_{12}}(s)`$ $`=`$ $`U_1(s)U_2(s+\widehat{1})U_1^{}(s+\widehat{2})U_2^{}(s),`$ (15)
respectively.

In the continuum, the Polyakov gauge is reduced to

$$A_0^\pm (x)=0,$$
(16)
whereas F12 gauge gives

$$F_{12}^\pm (x)=0.$$
(17)
6. We also consider simple abelian components extracted without gauge-fixing, where exact abelian dominance is proved analytically .

## 4 String Tension

As a first step, we measure abelian and monopole contributions to the string tension in various abelian projections. We used 100 configurations of $`32^3\times 16`$ lattice for the measurement. In this case, the critical point lies near $`\beta _c2.7`$. We set the gauge coupling $`\beta `$ to 2.5, so that the system is in the confinement phase. To reduce the statistical errors efficiently, we adopted the hypercubic blocking  to the original configurations.

The value of Polyakov loop correlators corresponds to the static potential between one pair of quark and anti-quark;

$$\text{Tr}P(0)\text{Tr}P^{}(R)=e^{V(R)/T},$$
(18)
where $`P(R)`$ is the Polyakov loop operator Eq.(14). $`V(R)`$ gives the inter-quark potential

$$V(R)=\sigma R\frac{\alpha }{R}+c$$
(19)
and $`T=1/(N_4a)`$ is the temperature of the system.

The abelian Polyakov loop operator is written as

$$P_a=\mathrm{exp}[i\underset{i=0}{\overset{N_41}{}}\theta _4(\stackrel{}{s}+i\widehat{4})].$$
(20)
Eq.(20) can be decomposed to photon and monopole parts  as follows:

$`P_a`$ $`=`$ $`P_pP_m,`$
$`P_p`$ $`=`$ $`\mathrm{exp}[i{\displaystyle \underset{i=0}{\overset{N_41}{}}}{\displaystyle \underset{s^{}}{}}D(\stackrel{}{s}+i\widehat{4}s^{})_\nu ^{}\overline{\mathrm{\Theta }}_{\nu 4}(s^{})],`$
$`P_m`$ $`=`$ $`\mathrm{exp}[2\pi i{\displaystyle \underset{i=0}{\overset{N_41}{}}}{\displaystyle \underset{s^{}}{}}D(\stackrel{}{s}+i\widehat{4}s^{})_\nu ^{}n_{\nu 4}(s^{})],`$
where we use an identity

$$\theta _4(s)=\underset{s^{}}{}D(ss^{})\left[_\nu ^{}\mathrm{\Theta }_{\nu 4}(s^{})+_4(_\nu ^{}\theta _\nu (s^{}))\right].$$

The abelian field strength tensor

$$\mathrm{\Theta }_{\mu \nu }(s)\theta _\mu (s)+\theta _\nu (s+\widehat{\mu })\theta _\mu (s+\widehat{\nu })\theta _\nu (s),(4\pi \theta _{\mu \nu }(s)<4\pi )$$

can be decomposed into two parts:

$$\mathrm{\Theta }_{\mu \nu }(s)\overline{\mathrm{\Theta }}_{\mu \nu }(s)+2\pi n_{\mu \nu }(s),(\pi \overline{\mathrm{\Theta }}_{\mu \nu }(s)<\pi ).$$

Here, $`\overline{\mathrm{\Theta }}_{\mu \nu }(s)`$ is interpreted as the electro-magnetic flux through the plaquette and the integer valued $`n_{\mu \nu }(s)`$ corresponds to the number of Dirac string piercing the plaquette. $`D(ss^{})`$ is the Coulomb propagator on a lattice.

Figures 6, 6, 8 and 8 show the values of $`SU(2)`$, abelian and monopole Polyakov loop correlators in MA, LA, MAWL and L-type gauges, respectively. The values of abelian and monopole Polyakov loop correlators in each gauge almost degenerate. The string tension $`\sigma `$ can be extracted from these values by fitting them to Eq.(18). Fitted lines are also plotted in the same figure. In the case of MA gauge, fitted values are consistent with the results by Bali et al . In the case of other gauges like a unitary gauge, one can not extract the string tension clearly from the abelian and the monopole Polyakov loop correlators due to large statistical errors.

Explicit values of the fitted string tension are shown in Table 1. They almost agree each other, although these four gauges have different gauge fixing condition in the continuum limit.

## 5 RG flows of the effective action in various abelian projections

To clarify what is happening in the monopole dynamics, we study the effective monopole actions in various gauges in this section.

### 5.1 Simulation method

Our method to derive an effective monopole action is the following. We generate $`SU(2)`$ gauge fields $`\{U_\mu (s)\}`$ using the standard $`SU(2)`$ Wilson action. We consider $`48^4`$ hypercubic lattice for $`\beta `$ from 2.1 to 2.5. We took 50 independent configurations after 10,000 thermalization sweeps. Then, we perform the abelian projection in six different gauge fixings to extract abelian gauge fields $`\{u_\mu (s)\}`$ from $`SU(2)`$ gauge fields.

One can define magnetic monopole currents from abelian field strength tensor along the way of DeGrand and Toussaint . We can define the monopole current $`k_\mu (s)`$ as

$$k_\mu (s)=\frac{1}{2}ϵ_{\mu \nu \rho \sigma }_\nu n_{\rho \sigma }(s+\widehat{\mu }).$$
(21)
By definition, it satisfies a current conservation law

$$_\mu ^{}k_\mu (s)=0,$$

where $`_\mu `$ and $`_\mu ^{}`$ denote the forward and the backward differences in $`\mu `$-direction respectively.

We want to get an effective monopole action $`S[k]`$ on the dual lattice, integrating out all degrees of freedom except for monopoles:

$`Z`$ $`=`$ $`{\displaystyle 𝒟Ue^{S[U]}\delta (X^\pm )\mathrm{\Delta }_F(U)}`$
$`=`$ $`{\displaystyle 𝒟u\left[𝒟Ce^{S[U]}\delta (X^\pm )\mathrm{\Delta }(U)\right]}`$
$`=`$ $`{\displaystyle 𝒟ue^{S_{eff}[u]}}`$
$`=`$ $`\left({\displaystyle }\right){\displaystyle 𝒟ue^{S_{eff}[u]}\delta (k,f(u))}`$
$`=`$ $`{\displaystyle \underset{s,\mu }{}}{\displaystyle \underset{k_\mu (s)=\mathrm{}}{\overset{\mathrm{}}{}}}\left({\displaystyle \underset{m,\nu }{}}\delta _{_\nu ^{}k_\nu (m),0}\right)e^{S[k]},`$
where $`U_\mu =C_\mu u_\mu `$ and $`X^\pm `$ is the off-diagonal element of the matrix $`X`$ which is diagonalized in the procedure of abelian projection. $`\mathrm{\Delta }_F(U)`$ is the Faddeev-Popov determinant and $`\delta (k,f(u))`$ gives the definition of the monopole current $`k`$ as a function of abelian gauge field $`u`$.

Above integrations are done numerically. We create vacuum ensembles of monopole currents using the Monte-Carlo method and the definition of the monopole current Eq.(21). Then, we construct the effective monopole action from monopole vacua using the Swendsen’s inverse Monte-Carlo method which was developed originally by Swendsen  and extended by Shiba and Suzuki .

We consider a set of independent and local monopole interactions which are summed up over the whole lattice. We denote each interaction term as $`S_i[k]`$. Then the effective monopole action can be written as a linear combination of these operators:

$$S[k]=\underset{i}{}g_iS_i[k],$$
(22)
where $`g_i`$ denotes the effective coupling constants. Explicit forms of the interaction terms are listed in Table 2 and 3. We determine the set of couplings $`\{g_i\}`$ from the monopole current ensemble $`\{k_\mu (s)\}`$ with the aid of an inverse Monte-Carlo method. Practically, we have to restrict the number of interaction terms. The form of action adopted here is 27 quadratic interactions and 4-point and 6-point interactions .

We perform a blockspin transformation in terms of the monopole currents on the dual lattice to study the RG flow. The $`n`$-step blocked current is defined by

$$K_\mu (s^{(n)})=\underset{i,j,l=0}{\overset{n1}{}}k_\mu (ns^{(n)}+(n1)\widehat{\mu }+i\widehat{\nu }+j\widehat{\rho }+l\widehat{\sigma }).$$
(23)
The blocked lattice spacing $`b`$ is given as $`b=na(\beta )`$ and the continuum limit is taken as the limit $`n\mathrm{}`$ for a fixed physical scale $`b`$. We determine the effective monopole action from the blocked monopole current ensemble $`\{K_\mu (s^{(n)})\}`$. Then one can obtain the RG flow in the 29-dimensional coupling constant space.

### 5.2 Numerical Results

The effective monopole action is determined successfully. All coupling constants which are contained in the effective monopole action are obtained with relatively small errors. We use the jackknife method for the error estimation. These effective monopole actions except in MA gauge are determined for the first time in this paper. Moreover, these effective monopole actions are determined from the blocked monopole configurations, too. The results are summarized as follows:

1. Only the quadratic interaction subspace seems sufficient in the coupling space for the low-energy region of QCD. Figure 10 and 10 show coupling constants 111Effective coupling constants for the blocking factor $`n=1`$ are omitted in Figures 10, 10 and 15-23. for 4-point and 6-point interaction terms versus physical scale $`b`$. In the case of MA, LA, MAWL and L-type gauges, these coupling constants take relatively larger absolute values for small $`b`$ region. They become negligibly small for large $`b`$ region. In the case of Polyakov, F12, no gauge fixings, coupling constants for 4-point and 6-point interaction terms take the values very close to zero in the whole region of $`b`$.
2. Typical case of the coupling constants for quadratic interaction terms versus squared distances in the lattice unit are shown in Figure 11. We see that coupling constants for the self interaction term $`g_1`$ and the nearest-neighbor interactions $`g_2`$ and $`g_3`$ are dominant, and $`g_2g_3`$. Other couplings decrease exponentially as distance between the two monopole currents grows. Such a behavior does not depend on a gauge coupling constant $`\beta `$. Therefore, we concentrate our analysis on the coupling constants of quadratic interaction terms, especially $`g_1`$ and $`g_2`$.
3. We used a standard iterative gauge fixing procedure for MA, MAWL and L-type gauges. In such a case, gauge fixing sweeps may be stuck for some local minima of a gauge fixing functional. Different local minima give rise to different gauge transformations, but they can not be distinguished from the view point of the iterative gauge fixing procedure. These are the lattice Gribov copies. Indeed, Bali et al. showed that the effect of such copies to the abelian string tension is not so small . To check the effect of copies to the effective couplings, we generate 100 of $`SU(2)`$ configurations on $`24^4`$ lattice at $`\beta =2.5`$. Then, we generate 7 of gauge equivalent configurations (i.e., copies) via a random gauge transformation. Using these gauge copies, we constructed effective monopole actions and compared their effective couplings. Figure 13 shows $`g_1`$ in the case of MA gauge. $`g_1`$ for the different blocking factors are described in different symbols. We see some fluctuations in $`g_1`$ for MA gauge. This is nothing but the effect of lattice Gribov copies. The effect of the copies, however, are negligibly small. Therefore, qualitative analyses which are given later will not be affected. In principle, LA gauge does not have such a copy . Indeed, we confirmed that effective couplings for LA gauge are not affected by Gribov copies (Figure 13).
4. Figure 15 and Figure 15 show the most dominant quadratic self coupling constant $`g_1`$ and quadratic nearest-neighbor coupling constant $`g_2`$ versus physical scale $`b`$ in the case of MA, LA, MAWL and L-type gauges, respectively. In these gauges, effective coupling constants take large values in small $`b`$ region and the scaling behavior (i.e., a unique curve for different blocking factor $`n`$) is seen even in small $`b`$ region. The effective actions which are obtained here appear to be a good approximation of the action on the renormalized trajectory corresponding to the continuum limit. In addition to this, coupling constants for these four gauges are very close to each other, although these gauges have a completely different form in the continuum limit.
5. However, in the case of Polyakov, F12 and no gauge fixings, coupling constants are different from those in the above four gauges (See, Figure 17 and Figure 17). In these gauges, coupling constants take smaller values and the scaling behavior is not seen especially in small $`b`$ region. To clarify the scaling properties of these coupling constants, we show the figures as showing a distinction between the different blocking factors $`n`$ in two typical gauges. In the case of Polyakov gauge (Figure 19), the coupling constants depend on the blocking factor $`n`$ strongly in small $`b`$ regions. On the other hand, in the case of LA gauge (Figure 19), renormalized coupling constants are lying on the unique curve.
6. Once the effective actions are fixed, we can see from energy-entropy balance of the system whether monopole condensation occurs or not. If the entropy of a monopole loop exceeds the energy, the monopole loop condenses in the QCD vacuum. In four-dimensional lattice theory, the entropy of a monopole loop can be estimated as $`\mathrm{ln}7`$ per unit loop length. It is determined by the random walk without backward tracking. The action can be approximated by the self interaction term $`g_1`$ alone since the interactions with two separate currents are almost canceled . The free energy per unit monopole length is approximated by

$$Fg_1\mathrm{ln}7,$$

since $`g_1`$ can be regarded as the self energy per unit monopole loop length. If $`g_1<\mathrm{ln}7`$, the entropy dominates over the energy, that is, monopole condensation occurs. In Figure 15 and Figure 17, we see that the entropy of the system dominates over the energy in the large $`b`$ region for all gauges. In other words, monopole condensation occurs  in the large $`b`$ region for all gauges.
7. Figures 21,21,23 and 23 show the RG flows projected onto $`g_1`$-$`g_2`$, $`g_1`$-$`g_5`$, $`g_1`$-$`g_7`$ and $`g_1`$-$`g_{10}`$ coupling planes, respectively. The effective coupling constants for all gauges seem to converge to the identical line for large $`b`$ region. This may show gauge independence of the monopole condensation in the low energy region. Although all coupling constants become very small in the large $`b`$ region, it is important that the slopes of the renormalization flows seem to converge in all gauges.

## 6 Summary

We have measured first the abelian and the monopole contributions to the string tension in four types of abelian projections, i.e., MA, LA, MAWL and L-type gauges. They show a good agreement with each other. Monopole string tension are extracted in the same manner as abelian string tension, and they agree also with each other. MA and LA gauges are not unique good gauges.

Next, we have determined the effective monopole actions in various gauges from monopole vacua using the modified Swendsen’s method. In the case of MA gauge, an effective monopole action has already been obtained in Ref. . In addition to this action, the effective monopole actions in Polyakov gauge, F12 gauge, LA gauge, MAWL gauge, L-type gauge and no gauge fixing are also determined for the first time in this paper. Moreover, these effective actions are determined on the blocked monopole vacua, too. In these effective actions, two point interactions are dominant, whereas 4-point and 6-point effective coupling constants are negligibly small in the infrared region. The RG flows seem to converge to the identical line when repeating the blockspin transformation. It is important that the slopes of renormalization flows in all gauges seem to converge. The data are compatible with the assumption of gauge independence of monopole dynamics in the continuum limit. Energy-entropy balance also tells us the monopole condensation occurs in the large $`b`$ region for all gauges.

###### Acknowledgments.

The authors would like to thank Fumiyoshi Shoji at Hiroshima University for fruitful discussions. This work is supported by the Supercomputer Project of the Institute of Physical and Chemical Research (RIKEN). A part of our numerical simulations have been done using NEC SX-5 at Research Center for Nuclear Physics (RCNP) of Osaka University.

## Appendix A Maximally Abelian Wilson Loop (MAWL) gauge

$`SU(2)`$ gauge field $`U_\mu (s)`$ can be parameterized by its isospin components. In this section, we denote each isospin component of $`U_\mu (s)`$ as $`U_0(s,\mu )`$, $`U_1(s,\mu )`$ and so on, for simplicity. This gauge is realized with maximizing the abelian Wilson loop of $`1\times 1`$ size:

$$R=\underset{s,\mu \nu }{}\mathrm{cos}\mathrm{\Theta }_{\mu \nu }(s),$$
(24)
where the abelian link field is extracted as

$$\theta (s,\mu )=\mathrm{arctan}(U_3(s,\mu )/U_0(s,\mu )).$$
(25)
Let us consider an infinitesimal gauge transformation of $`U`$,

$$U^{}(s,\mu )=(1+i\alpha _i(s)\sigma _i)(U_0(s,\mu )I+iU_j(s,\mu )\sigma _j)(1i\alpha _k(s+\widehat{\mu })\sigma _k).$$
(26)
This gives

$`\delta U_0(s,\mu )`$ $`=`$ $`(\alpha _i(s)\alpha _i(s+\widehat{\mu }))U_i(s,\mu ),`$ (27)
$`\delta U_k(s,\mu )`$ $`=`$ $`(\alpha _k(s)\alpha _k(s+\widehat{\mu }))U_0(s,\mu )ϵ_{ijk}(\alpha _i(s)+\alpha _i(s+\widehat{\mu }))U_j(s,\mu ).`$ (28)
Then $`R`$ changes as

$$\delta R=\underset{s,\mu \nu }{}\mathrm{sin}\mathrm{\Theta }_{\mu \nu }(s)(\delta \theta (s,\mu )+\delta \theta (s+\widehat{\mu },\nu )\delta \theta (s+\widehat{\mu },\mu )\delta \theta (s,\nu )),$$
(29)
where

$$\delta \theta (s,\mu )=\frac{U_0(s,\mu )\delta U_3(s,\mu )U_3(s,\mu )\delta U_0(s,\mu )}{U_0^2(s,\mu )+U_3^2(s,\mu )}.$$
(30)
One can check that $`R`$ is invariant under the $`U(1)`$ transformation. Hence we do not need to consider the $`\alpha _3(s)`$ part. First, let us consider the $`\alpha _1`$ part. Since there is the sum over whole lattice sites $`s`$, one can shift the site variable, for example, $`s`$ to $`s\widehat{\mu }`$. Also one can use the (anti)symmetric property with respect to $`\mu `$ and $`\nu `$ directions. Finally one gets,

$`{\displaystyle \frac{\delta R}{2}}`$ $`=`$ $`{\displaystyle \underset{s,\mu \nu }{}}(\alpha _1(s)X_1(s,\mu ,\nu )+\alpha _2(s)X_2(s,\mu ,\nu )),`$
$`X_1(s,\mu ,\nu )`$ $`=`$ $`W_1{\displaystyle \frac{U_1(s,\mu )U_3(s,\mu )U_0(s,\mu )U_2(s,\mu )}{U_0^2(s,\mu )+U_3^2(s,\mu )}}`$
$`+W_2{\displaystyle \frac{U_1(s\widehat{\mu },\mu )U_3(s\widehat{\mu },\mu )+U_0(s\widehat{\mu },\mu )U_2(s\widehat{\mu },\mu )}{U_0^2(s\widehat{\mu },\mu )+U_3^2(s\widehat{\mu },\mu )}},`$
$`X_2(s,\mu ,\nu )`$ $`=`$ $`W_1{\displaystyle \frac{U_2(s,\mu )U_3(s,\mu )+U_0(s,\mu )U_1(s,\mu )}{U_0^2(s,\mu )+U_3^2(s,\mu )}}`$
$`+W_2{\displaystyle \frac{U_2(s\widehat{\mu },\mu )U_3(s\widehat{\mu },\mu )U_0(s\widehat{\mu },\mu )U_1(s\widehat{\mu },\mu )}{U_0^2(s\widehat{\mu },\mu )+U_3^2(s\widehat{\mu },\mu )}},`$
$`W_1`$ $`=`$ $`\mathrm{sin}\mathrm{\Theta }_{\mu \nu }(s)\mathrm{sin}\mathrm{\Theta }_{\mu \nu }(s\widehat{\nu }),`$
$`W_2`$ $`=`$ $`\mathrm{sin}\mathrm{\Theta }_{\mu \nu }(s\widehat{\mu }\widehat{\nu })\mathrm{sin}\mathrm{\Theta }_{\mu \nu }(s\widehat{\mu }).`$
When we write $`X^\pm =X_1\pm iX_2`$, it is easy to see $`X^\pm `$ transforms covariantly under the residual $`U(1)`$.

Finally, one gets the matrix which is diagonalized in this gauge,

$`X(s)`$ $`=`$ $`{\displaystyle \underset{\mu \nu }{}}[{\displaystyle \frac{\mathrm{sin}\mathrm{\Theta }_{\mu \nu }(s)\mathrm{sin}\mathrm{\Theta }_{\mu \nu }(s\widehat{\nu })}{U_0^2(s,\mu )+U_3^2(s,\mu )}}(U(s,\mu )\sigma _3U^{}(s,\mu ))`$
$`+{\displaystyle \frac{\mathrm{sin}\mathrm{\Theta }_{\mu \nu }(s\widehat{\mu }\widehat{\nu })\mathrm{sin}\mathrm{\Theta }_{\mu \nu }(s\widehat{\mu })}{U_0^2(s\widehat{\mu },\mu )+U_3^2(s\widehat{\mu },\mu )}}(U^{}(s\widehat{\mu },\mu )\sigma _3U(s\widehat{\mu },\mu ))].`$
Because of the non-locality of the gauge condition, one can not calculate the gauge transformation matrix which diagonalizes $`X(s)`$ in a simple way. Therefore, we employed an iterative updation procedure to satisfy the gauge condition.

1. Make a trial gauge transformation, adopting $`\alpha _1`$ and $`\alpha _2`$ as follows:
$`\alpha _1(s)=\kappa X_1(s)`$,
$`\alpha _2(s)=\kappa X_2(s)`$.
2. Measure $`R`$. If $`R`$ becomes larger than before, accept this trial and repeat step 1. If $`R`$ becomes smaller than before, take $`\kappa _{new}=\kappa _{old}/2`$ and adopt the gauge transformation using this $`\kappa _{new}`$ with respect to the configuration before trial, and then repeat step 1.
3. If the off-diagonal element of $`X(s)`$ becomes smaller than a suitable threshold (we set this to $`1.0`$), one can regard the gauge fixing procedure as having been completed.

We set an initial value of $`\kappa `$ to $`0.1`$. $`R`$ can be maximized as long as we take $`\kappa >0`$. We apply the MA gauge fixing as a preconditioning for MAWL gauge fixing and then we perform the above procedure to the MA fixed configuration. This preconditioning is required to improve a convergence property of the MAWL gauge fixing. We have to note that the configurations which obtained via the above procedure are not perfectly gauge fixed because the off-diagonal element of $`X(s)`$ still remain not very small.


text:
# Measuring the gap in ARPES experiments

## I Introduction

The question of symmetry of the superconducting gap is crucial to understand the nature of superconductivity in cuprates.TsueiRMP00  Although not being phase sensitive, the angle-resolved photoemission spectroscopy (ARPES) has been considered as the only experimental tool from which both the superconducting and pseudo-gap anisotropy could be quantitatively derived.ShenPRL93 ; DingPRL95 ; MarshallPRL96 ; DingPRL97 ; MesotPRL99  The binding energy of the leading edge of the photoemission spectrum, usually called the “leading edge gap” (LEG),MarshallPRL96 ; DingPRL97  is the model-independent quantity which can be measured in the modern ARPES experiments with the highest, up to 1 meV accuracy,BorisXXX02  and is, actually, the best quantity from which we can judge about the real gap values—in absence of a widely accepted model for the gap formation in these compounds, there is no direct way to extract the real gap value from the ARPES spectra.

The LEG is determined as a lowest binding energy at which the energy distribution curve (EDC) reaches half of its maximum (here and thereafter we keep the same notation “LEG” for either gapped or non-gapped spectra). Therefore it is understandable, and usually admitted, that the LEG should depend on any parameters that determine the EDC line shape and, unless the relation between the LEG and the real gap is known, can be considered only as a qualitative representation of the real gap. Moreover, the mentioned parameters (and EDC line shape) are substantially changing during any experiment, e.g. with temperature, along the Fermi surface (FS), etc. So, analysing the experimental data it is very important to distinguish between the artificial variations of the LEG and changes caused by the real gap in the electronic density of states. Despite the big importance of this question, there is no systematic study of the LEG available today. Hence, even the qualitative relation between the LEG and the real gap can be doubtful.

In this paper, by means of a numerical simulation we present a systematic investigation of the LEG as a function of a number of physical and experimental parameters which govern the ARPES spectra of the Bi-2212 cuprates.

Starting with the non-gapped case, we demonstrate that the dependencies of the LEG on the temperature and slight deviations in the momentum space (even within the momentum resolution window) from the FS ($`\text{k}_F`$) result in strong variations of the LEG values which are comparable with the values of the superconducting and pseudo-gaps. We also show that at the low temperatures, choosing an appropriate criterion for the $`k_F`$ determination—the “maximum intensity” or “minimum gap locus”—the absolute value of the LEG is not very sensitive to the “physical” parameters of the model spectral function and therefore can be easily obtained. As an example we consider the experimental data of the LEG for an underdoped Bi-2212 in the superconducting state,BorisXXX02  and prove that the nodal direction there is really the node—the gap is zero within the experimental accuracy. Finally, switching to the gapped case, we show that the absence of a “cusp” of the LEG vs FS angle dependencies (so called “U-shape”),BorisXXX02  cannot be explained neither by the interplay between the temperature and Fermi-function nor by other “artificial” parameters which may vary over the FS: momentum resolution, self-energy and band-structure. As an intermediate result, we present the tight-binding parameters for the considered compounds.

## II Model description

We start with a model spectral function

$`A(\omega ,ϵ,T){\displaystyle \frac{|\mathrm{\Sigma }^{\prime \prime }(\omega ,T)|}{(\omega ϵ\mathrm{\Sigma }^{}(\omega ,T))^2+\mathrm{\Sigma }^{\prime \prime }(\omega ,T)^2}},`$ (1)
which is an essential part of the photocurrent measured in the experiments

$`I(\text{k},\omega ,T,h\nu )[(M(\text{k},h\nu )A(\text{k},\omega ,T)f(\omega ,T)]`$
$`R_\text{k}R_\omega +B(\omega ,T).`$ (2)
Here $`M(\text{k},h\nu )`$ represents the dependence of the squared one-electron matrix element on in-plane electron momentum k and excitation energy $`h\nu `$, $`f(\omega ,T)=1/[\mathrm{exp}(\omega /T)+1]`$ is the Fermi function. $`R_\text{k}`$ and $`R_\omega `$ represent the functions of the momentum and energy resolutions respectively. We include the energy resolution as a convolution with a Gaussian of the $`R_\omega `$ full width at half maximum (FWHM) and the momentum resolution through a simple summation over a part of the detector area,

$`AR_\text{k}={\displaystyle \underset{\delta \text{k}}{}}A,`$ (3)
where $`\delta \text{k}`$, in general case, is a 2D area with the $`R_k\times R_k`$ dimensions (momentum resolutions parallel and perpendicular to the entrance slit of an analyser respectively) whose orientation in k-space depends on k and the experimental geometry but, in the case of the big FS of superconducting cuprates, one can consider $`\delta \text{k}`$ as an one-dimensional cut $`\pm R_k(\text{k})/2`$ perpendicular to the FS.

For the momentum independent background we assume an empirical relation: $`B(\omega ,T)(1+b\omega ^2)f(\omega \mathrm{\Delta }_b,T+T_b)`$ with $`b`$ = 1 eV-2, $`\mathrm{\Delta }_b`$ = 5 meV and $`T_b`$ = 90 K. The parameters $`\mathrm{\Delta }_b`$ and $`T_b`$ depend on doping and $`\mathrm{\Delta }_b`$ also depends on temperature KordPRL  but we neglect this here because, as it is shown below, the typical signal-to-background ratios are large and the influence of the background on the LEG value is rather weak.

For the imaginary part of the self energy we also use an empirical formula Self-energy  $`\mathrm{\Sigma }^{\prime \prime }(\omega ,T)=\sqrt{(\alpha _0\omega )^2+(\beta _0T)^2}`$ which fits best the experimental data (see Refs. BorisPRB01, ; KordPRL, ). In the vicinity of the Fermi-level ($`\pm `$ 50 meV) the real part of the self-energy can be well approximated by its linear term $`\mathrm{\Sigma }^{}(\omega ,T)\lambda (T)\omega `$ ($`\lambda >0`$) and the spectral function (1) can be rewritten in a renormalized form

$`A(\omega ,\text{k},T){\displaystyle \frac{1}{1+\lambda }}{\displaystyle \frac{|\mathrm{\Sigma }^{\prime \prime }(\omega ,T)|}{(\omega \epsilon (\text{k}))^2+\mathrm{\Sigma }^{\prime \prime }(\omega ,T)^2}},`$ (4)
with the renormalized imaginary part of the self-energy $`\mathrm{\Sigma }^{\prime \prime }=\sqrt{(\alpha \omega )^2+(\beta T)^2}`$, $`\alpha =\alpha _0/(1+\lambda )`$, $`\beta =\beta _0/(1+\lambda )`$, and the renormalized dispersion $`\epsilon =ϵ/(1+\lambda )`$. The last one, in a small vicinity of the Fermi energy $`E_F`$ on the path perpendicular to FS, can be written using the renormalized Fermi velocity, $`\epsilon =v_Fk`$, which in the same way relates to the bare one: $`v_F=u_F/(1+\lambda )`$.

At last, when we take into account the bilayer splitting we include it as a simple superposition of the photocurrent from bonding (“b”) and antibonding (“a”) bands:

$`I[M_aA(\epsilon _a)+M_bA(\epsilon _b)]fR_\text{k}R_\omega +B.`$ (5)
## III Gapless case

First, we examine the LEG behaviour assuming that there is no real gap at all. It is natural to expect that the leading edge position of an EDC should depend on every physical parameter which forms the spectral function and the background as well as on energy and momentum resolutions. Among others, the temperature seems to be the most crucial here (see Ref. BorisPRB01, ).

### III.1 Temperature dependence

Fig. 1 shows the temperature dependence of the model photoemission spectra (a–c) defined by Eqs. II and 4 and the LEG values (d–f) derived therefrom. Here we use the following parameters: $`\alpha `$ = 1, $`\beta `$ = 2, $`v_F`$ = 2 eVÅ, $`R_\omega `$ = 20 meV, and $`R_k`$ = 0.07 Å-1 that corresponds to 0.2 angular resolution for the nodal point at 21 eV excitation energy. Panel (a) represents the photocurrent (II) at $`k=k_F`$ ($`k_F`$-EDCs) at different temperatures from 10 to 300 K (shown in color from blue to red with 10 K step). Panel (b) shows two extreme spectra but normalized to their maximum in order to illustrate the motion of the leading edge with temperature. The LEG($`T`$) dependencies are shown on panel (d) and (e) for different momentum and energy resolutions.

The signal-to-background ratio (which can be estimated as a ratio of peak-to-tail photocurrent values) for the EDCs presented on panel (a) is chosen to be similar to typical photoemission data, and panel (e) also shows that the influence of the background in this case is negligible—the LEG($`T`$) for the spectra without the background are shown as dashed lines.

It can be seen that, although the LEG strongly depends on temperature and moves about 15 meV going from 10 to 300 K, which should be taken into account when the temperature dependence of the real gap is studied, this dependence, for reasonable energy and momentum resolutions ($`R_\omega `$ = 20 meV, $`R_k`$ = 0.07 Å-1: shown as a bold red curve in Fig. 1, d–f), is quite monotonic (i.e. has no maxima on both itself and its derivative) and cannot imitate a rapid gap opening at a certain temperature.

The situation is different if EDCs are taken from a $`k`$-point which does not exactly coincide with $`k_F`$. Panels (c) and (f) of Fig. 1 illustrate this case. Panel (c) shows EDCs for $`k=k_F`$ \+ 0.01 Å-1 for the same temperature range (10–300 K) and panel (f) shows the LEG($`T`$) dependencies for different $`kk_F`$ (from -0.01 to 0.01 Å-1). One can see that the unoccupied part of the Brillouin zone ($`k>k_F`$) is the most dangerous in this sense. If one steps away from $`k_F`$ in the unoccupied direction in 0.01 Å-1 (about 0.3 for the nodal point at 21 eV excitation energy), the LEG changes about 30 meV going from 10 to 300 K but, what is most important, a kink appears at about 80 K on the LEG vs $`T`$ dependence. This kink can be easily misinterpreted as a gap opening, so we may conclude that the uncertainty in $`k_F`$ determination (see discussion in Ref. BorisPRB01,  and below) is the most crucial parameter for the correct LEG evaluation.

Regarding Fig. 1 (a) we also would like to note that due to the temperature dependent self-energy the area under each EDC is not constant with temperature, as one could expect according to the corresponding sum rule (e.g., see Ref. RanderiaPRL95, ).

### III.2 Resolution

The dependence of the LEG on resolutions can be estimated from Fig. 1 (d) and (e). The dependence on energy resolution (Fig. 1 (e)) is stronger but it is not a problem because we believe that normally it remains constant during the experiment and the absolute value can be easily taken into account. The momentum resolution can change during the experiment. For example, for the experimental geometry described in Ref. BorisPRB01, , the actual momentum resolution changes going from the $`\mathrm{\Gamma }`$$`(\pi ,\pi )`$ to $`(\pi ,\pi )`$$`(\pi ,0)`$ FS crossings: it is defined by the angular resolution along the entrance slit of the analyser for the former and by the angular resolution perpendicular to the slit for the latter. As a result of this, the related momentum resolution can change by a factor of two but, as one can see from Fig. 1 (d), this will move the leading edge about 1 meV only (even for the ”worst” case of 10 K).

### III.3 Band structure

In contrast to the temperature, the other parameters that enter the $`k_F`$-EDCs are expected to depend on the in-plane momentum of the electrons. One can divide these parameters into “structural”, which come from the bare electronic structure, and “physical”, which come from interactions and reveal themselves as a self-energy.

The structural parameters of interest here are the bare Fermi velocity and the bilayer splitting. Fig. 2 represents the result of the tight-binding fit of the experimental data published in Refs. BorisXXX02,  and KordPRB02, . The precisely determined FS of Bi-2212 (at least its bonding sheet, see Ref. KordPRB02, ) at room temperature (300 K) can be well described within the four-band model introduced in Ref. AndersenJPCS95, . Within this model, the simplified bare dispersion

$`ϵ_{a,b}`$ $`=`$ $`\mathrm{\Delta }ϵ2t(\mathrm{cos}k_x+\mathrm{cos}k_y)`$ (6)
$`+`$ $`4t^{}\mathrm{cos}k_x\mathrm{cos}k_y2t^{\prime \prime }(\mathrm{cos}2k_x+\mathrm{cos}2k_y)`$
$`\pm `$ $`t_{}(\mathrm{cos}k_x\mathrm{cos}k_y)^2/4,`$
where second and third nearest-neighbor intra-plane hopping integrals ($`t^{}`$ and $`t^{\prime \prime }`$) are provided by the ”Cu $`s`$” orbital and interlayer hopping is described by the $`t_{}`$ (the bilayer splitting, the energy distance between bonding and antibonding CuO-bilayer bands, is $`2t_{}`$).AndersenJPCS95 ; LiechtensteinPRB96

To fit the data we use the following procedure. First, the relative parameters (relative to $`t`$) are determined from the room temperature FS maps of Ref. KordPRB02, , then the bilayer splitting is taken into account \[the value $`t_{}/t`$ is determined from the $`(\pi ,\pi )(\pi ,0)`$ FS crossing\]. For example, for the overdoped sample (OD, $`T_c`$ = 69 K), using the assumption that $`t^{\prime \prime }t^{}/2`$,PavariniPRL01  these parameters are: $`t^{}/t`$ = 0.23, $`t^{\prime \prime }/t`$ = 0.11, $`t_{}/t`$ = 0.21, $`\mathrm{\Delta }ϵ/t`$ = 1.08. From here one can estimate the variation of the bare Fermi velocity $`u_F`$ around the FS as a ratio of its value in the antinodal ’A’ point (FS crossing along $`(\pi ,\pi )(\pi ,0)`$ direction) to the nodal ($`(0,0)(\pi ,\pi )`$ crossing) one. This ratio, $`u_F`$(A)/$`u_F`$(N), is 0.71 for the bonding band and 0.25 for the antibondig band (for OD 69 K). We stress that these values are determined by the FS shape, as presented on the left panel of Fig. 2, and do not depend on the energy scale.

In order to estimate the “scale” $`t`$ one needs to know some absolute value of the bare band. We can estimate the ’bare’ Fermi velocity as 4.0 eVÅ, subtracting the real part of the self-energy, which has been determined using the Kramers-Krönig transform of its imaginary part, from the renormalized experimental dispersion. The tight-binding parameters which have been calculated for two samples with this value of $`u_F`$, assuming the constancy of the bilayer splitting,ChuangXXX01  are presented in Table 1. The right panel of Fig. 2 shows the obtained bare dispersion along the $`\mathrm{\Gamma }`$-Y-M-$`\mathrm{\Gamma }`$ path (see the left panel). The inset in Fig. 2 zooms in the Y-M-X region (exactly the $`(\pi ,\pi /6)(\pi ,\pi /6)`$ path). This gives the bilayer splitting in ($`\pi `$,0)-point of 160 meV and the absolute bare Fermi velocities: $`u_F`$(A, bonding) = 2.8 eVÅ, $`u_F`$(A, antibonding) = 1.0 eVÅ, while $`u_F`$(N) = 4.0 eVÅ.

One can expect that the renormalized Fermi velocity $`v_F=u_F/(1+\lambda )`$ varies even more strongly than the bare one going from the node to the antinode due to a change of the ’coupling strength’ $`\lambda `$.lambda ; GromkoXXX02  So, going back to the simulation procedure, we examine an effect of $`v_F`$ variation on the leading edge position in the range from 0.5 to 5 eVÅ. The results are represented in Fig. 3 (a,b).

From Fig. 3 (a) one can see that the influence of $`v_F`$ on LEG increases with the decreasing temperature but remains rather weak – about 4 meV change at 10 K over the whole examined range (0.5–5 eVÅ) or in the more reasonable range from 1.5 to 3.5 eVÅ (shadowed on panel b) inferred from the experiment,KordPRB02  this variation is only 2 meV for $`k=k_F`$. Fig. 3 (b) also shows that if one follows the minimum gap locus (MGL) criteria DingPRL97 —choosing the location in $`k`$ for which the LEG($`k`$) along a FS cut has a minimum—this change is even less but, if $`k`$ deviates more ($`|kk_F|>R_k/2`$ = 0.035 Å-1), the change rapidly increases.

The influence of the bilayer splitting on the leading edge is shown in Fig. 3 (c,d). Here $`\mathrm{\Delta }\epsilon =\epsilon _a\epsilon _b`$ is the value of the splitting, and $`\mu =M_a/M_b`$ (see Eq. 5) and $`k=k_F`$ of the bonding band. One can see that at 120 K the influence of other FS sheet (antibonding in this case) can produce some minor variations (within about 2 meV) in the LEG around the FS but at low temperature (30 K) these variations are negligible.

### III.4 Self-energy

The influence of the self-energy parameters on the leading edge position is presented in Fig. 4. Although, due to the presence of the bilayer splitting, there is no systematic information currently available about how these parameters change over the FS, one can expect that in the antinodal point the coupling strength is larger than in the nodal one and, therefore, the self-energy parameters $`\alpha `$ and $`\beta `$ (see above) can be a few times bigger—two times seems to be a reasonable estimation.BorisPRB01 ; KordPRL

Fig. 4 (a) shows that the variation of $`\alpha `$ even in a much wider range (from 0.2 to 4) moves the leading edge within about 1 meV only, if $`k=k_F`$. Small deviation from the $`k_F`$ noticeably increases the amplitude of the LEG($`\alpha `$) dependencies (Fig. 4 (b)).

The dependencies of the LEG on $`\beta `$ within the same range, being more temperature dependent, are similar in amplitude (Fig. 4 (c,d)).

## IV Gapped case

### IV.1 Minimum gap locus

Among the $`k_F`$ determination criteria, the “minimum gap locus” (MGL) intuitively seems to be the most suitable in a gapped state. On the other hand, as it is shown above, the MGL, if applied to the $`k_F`$ determination, gives overestimated values. What is important here is that the $`k`$-locations which come from the MGL-criterion, $`k_{MGL}`$, are in the “dangerous” region, where a strong dependence of the LEG on temperature and other parameters is expected. We apply the simulation to check at which conditions the MGL method can be used.

Fig. 5 shows the dependencies of the leading edge position on $`k`$ along a cut perpendicular to the FS for temperatures from 10 to 300 K (from top to bottom) with 10 K steps. We use parameters typical for the $`\mathrm{\Gamma }`$-X crossing: $`\alpha `$ = 1, $`\beta `$ = 2, $`v_F`$ = 2 eVÅ, $`R_\omega `$ = 20 meV, $`R_k`$ = 0.07 Å-1. Fig. 5 (a) represents the gapless case and Fig. 5 (b) shows the same but taking into account a gap which opens at $`T_c`$ = 90 K. We used the BCS modified spectral function

$`A{\displaystyle \frac{u_\text{k}^2\mathrm{\Sigma }^{\prime \prime }}{(\omega +E_\text{k})^2+\mathrm{\Sigma }_{}^{\prime \prime }{}_{}{}^{2}}}+{\displaystyle \frac{(1u_\text{k}^2)\mathrm{\Sigma }^{\prime \prime }}{(\omega E_\text{k})^2+\mathrm{\Sigma }_{}^{\prime \prime }{}_{}{}^{2}}},`$ (7)
where $`E_\text{k}=\sqrt{\epsilon _\text{k}^2+\mathrm{\Delta }(T)^2}`$, $`u_\text{k}^2=\frac{1}{2}(1+\epsilon _\text{k}/E_\text{k})`$, and the gap function for $`T<T_c`$ is approximated as

$`\mathrm{\Delta }(T)=\mathrm{\Delta }_0\left[1+\left(1.74\sqrt{1{\displaystyle \frac{T}{T_c}}}\right)^{2n}\right]^{\frac{1}{2n}}`$ (8)
with $`\mathrm{\Delta }_0=1.76T_c`$ (temperature in energy units) and $`n=4`$.

The gapped region is shadowed in Fig. 5 (b) with the blue color, and the gapped LEG($`k`$, 10 K) dependency is also shown in Fig. 5 (a) as a dashed curve to compare with the gapless case.

From the results presented in Fig. 5 one can make the following conclusions. (1) The value of $`k`$ at which the leading edge position reaches a minimum, $`k_{MGL}`$, is strongly temperature dependent—thus it is important to realize that the MGL spectra measured at different temperatures correspond to rather different $`k`$-points. (2) The difference between $`k_{MGL}`$ and $`k_F`$ is improperly big at high temperatures (above $`T_c`$) and remains quite big at lower temperatures for both gapped and non-gapped cases. Using a criterion $`|k_{MGL}k_F|<R_k/2`$ one can conclude that the MGL method cannot be used for precise $`k_F`$ determination for $`T>`$ 40 K, at least in the Bi-2212 case. (3) If one stays at a certain $`k`$-point, measuring the LEG as a function of temperature, then the region where the real gap opening can be distinguished from the artificial one can be estimated as $`k<k_F+R_k/2`$. In this case the maximum intensity method,AebiPRL94  which underestimates the $`k_F`$,BorisPRB01  seems to be the most suitable. (4) At low temperatures ($`T40K`$) the MGL method is quite good as a method of the LEG determination—the experimental uncertainty in $`k_{MGL}`$ is big here because of very flat LEG($`k`$) dependence but, as another consequence of this, the LEG uncertainty should be small.

### IV.2 LEG in the case of real gaps

Now we consider how the real gap can influence the LEG value. Fig. 6 represents the results. Panels (a) and (d) correspond to $`k_F`$-spectra and panels (b) and (e) show the LEG of the spectra from $`k=k_{MGL}`$. The $`k_{MGL}(\mathrm{\Delta })`$ dependencies at different temperatures are shown in panel (c). Note, that all values are shown as function of $`\mathrm{\Delta }(T)`$ rather than $`\mathrm{\Delta }_0`$.

If we accept the use of Eq. 7 for the gapped spectral function the question how the imaginary part of the self-energy depends on the gap is unclear—in other words, the $`\mathrm{\Sigma }^{\prime \prime }(\mathrm{\Delta })`$ function is considered to be “more model dependent” than the spectral function $`A(\mathrm{\Delta })`$. The LEG($`\mathrm{\Delta }`$) dependencies presented in Fig. 6 (a,b) have been obtained with a $`\mathrm{\Delta }`$-independent self-energy. In Fig. 6 (d,e) we model the influence of the gap on the imaginary part of the self energy as $`\mathrm{\Sigma }^{\prime \prime }(\omega ,\mathrm{\Delta },T)=\sqrt{(\alpha _g\omega )^2+(\beta T)^2}`$, where

$`\alpha _g(\omega ,\mathrm{\Delta })=\alpha \{\begin{array}{cc}\omega ^2\mathrm{\Delta }^2\hfill & \text{if }|\omega |<\mathrm{\Delta },\hfill \\ 1\hfill & \text{if }|\omega |>\mathrm{\Delta },\hfill \end{array}`$ (9)
which gives the BCS asymptotic $`\mathrm{\Sigma }^{\prime \prime }(\omega )\omega ^3`$ at $`\omega 0`$ and zero temperature.

One can see that there is some “rounding” (deviation from linear) of the LEG($`\mathrm{\Delta }`$) curves at low $`\mathrm{\Delta }`$ which is explained by cutting finite width EDCs by the Fermi-function. These “rounding” increases with temperature but it is possible to conclude that for low temperatures ($`T<`$ 30 K) it remains rather small (less than 2 meV) for all $`k`$-points and $`\mathrm{\Sigma }^{\prime \prime }(\omega )`$ models considered. This means, for example, that the U-shape of the LEG($`\phi `$phi  observed in Ref. BorisXXX02,  cannot be explained by such a rounding of a sharp cusp expected in a simplest $`d`$-wave gap function.

The comparison of Fig. 6 (a,b) and (d,e) respectively shows that at higher $`\mathrm{\Delta }`$ the LEG($`\mathrm{\Delta }`$) dependencies become model dependent. Then, from the “more physical” model of $`\mathrm{\Sigma }^{\prime \prime }(\omega )`$, determined by Eq. 9, we can get a rough estimation of the coefficient between the LEG and $`\mathrm{\Delta }`$: $`d`$LEG($`\mathrm{\Delta }`$)$`/d\mathrm{\Delta }0.5`$.

## V Discussion

With the above simulations we mostly focused on the gapless case—we examined the leading edge position of the photoemission spectra in case when there is no gap present in the electron density of states.

It has appeared that the leading edge position of such a gapless spectrum depends on temperature, resolutions, band structure and self-energy parameters, which is not surprising of course. A surprising result is that under certain circumstances these dependencies can be quite complicated (exhibit non-monotonic rate) that can be misinterpreted as a physical transition.

These “circumstances” are a wrong position in momentum space, i.e. when $`k`$-location of the spectrum does not coincide with $`k_F`$ but is uncertain. If the parameters mentioned above are well determined experimentally and can be taken into account with the described simulation, the uncertainty in momentum cannot be dealt with in the same way, by definition. More precisely, this uncertainty is the uncertainty in $`kk_F`$ and, therefore, consists of the two: $`k`$-uncertainty, which is determined by the $`k`$ step of the measurements, and $`k_F`$-uncertainty which is related to a $`k_F`$ determination procedure. The first one can be reduced by refining the $`k`$ step of the measurements, which should be significantly less that the momentum resolution. The second problem is well known and heavily discussed in the past but we would like to address this question a bit more.

One way to define an uncertainty of the method is to model how far the directly determined value of $`k_F`$ deviates from the real one—for example, in the “maximum intensity” method, how the value of $`k`$, at which the photocurrent intensity exhibits a maximum, deviates from the real $`k_F`$. It has been shown that for high temperatures the maximum intensity method, especially in combination with the normalization to “highest binding energy”, is much more precise than others, such as “$`n(\text{k})`$”, “minimum gap locus” or “$`\mathrm{\Delta }T`$” methods,BorisPRB01 ; KordPRB02  but modelling this, one can find that it is difficult to suggest such a universal procedure of the $`k_F`$ determination which always gives the best result. All of the methods give some deviations and then one can reason that it is possible to determine the real $`k_F`$ taking these deviations into account turning the uncertainty to zero. Then we come to another uncertainty: how much these deviations depend on the model of the spectral function we use. This is an interesting and complicated question which can be subject of a separate study. Here we apply the word “uncertainty” of the method in the first sense and, in this respect, the maximum intensity method appeares to be the most accurate.

What is the real benefit one can extract from the results obtained or what do these results contribute to the gap problem? There are two kinds of gap measurements which are considered to be “model independent”: the gap as a function of temperature, LEG($`T`$), and momentum, LEG(k), which is usually presented as a function of FS angle, LEG($`\phi `$). In both cases the obtained “gap” is relative: in the first case it is relative to the LEG at the highest temperature, in the second case it is relative to the LEG in nodal direction.

The presented simulation shows that even when a real gap is zero the LEG is strongy temperature dependent (possibly in a complex way), and the amplitude of this dependency is comparable with the estimated real gap values in the superconducting cuprates. On the other hand, the momentum dependence of such a gapless LEG at low temperatures is shown to be rather weak for all possible variations of physical, structural and experimental parameters. But what is most important is that knowing these parameters one can estimate the absolute leading edge position and, therefore, answer the question whether there is any gap in a given spectrum or not.

For example, the LEG at 40 K and with the other parameters used in Ref. BorisXXX02,  for the nodal direction is equal to $`7`$ meV (see Fig. 1, df) which precisely coincides with the value obtained in the cited paper. This fact allowed us to make a firm conclusion that the nodal direction in an underdoped Bi-2212 in superconducting state is really the node—the gap is zero—which is in favor of $`d`$-symmetry of both superconducting and pseudo-gaps.

Another conclusion comes from Fig. 6: the “U-shape” of the LEG($`\phi `$) dependencies, which have been recenly observed for underdoped samples, BorisXXX02  cannot be explained neither by the interplay between the temperature and Fermi-function nor by other “artificial” parameters which may vary over the FS: momentum resolution, Fermi velocity, self-energy parameters. Considering the absence of the gap in the node as an evidence for $`d`$-wave symmetry of the order parameter in Bi-2212, this U-shape can be a consequence of higher harmonics of the gap due to the long range of the pairing interaction, as it was suggested earlier,MesotPRL99  or due to an interplay between the superconducting and pseudo-gaps. We want to note that this conclusion is valid for low temperatures only (superconducting state for Bi-2212), for higher temperatures (e.g. for pseudo-gap state) the influence of the artificial parameters can be essential.

## VI Summary

In conclusion, a systematic study of the LEG is presented. Even in the absence of a real gap, the leading edge position of the ARPES spectrum depends on many factors: energy and momentum resolutions, Fermi-velocity, the self-energy, but mostly on temperature. At some circumstances, when the spectra are measured away from $`k_F`$ (if $`kk_F>`$ 0.005 Å-1 for typical parameters of an experiment on Bi-2212), the temperature dependence of the LEG can be quite complicated and exhibit non-monotonic rate that can be misinterpreted as a physical transition. At low temperatures the dependence of the LEG on other parameters is rather weak and cannot be responsible for observed in experiments the LEG variation over the FS.

The absolute gap values which have been derived from the presented simulation prove that the nodal direction in the underdoped Bi-2212 in the superconducting state is really the node (the gap is zero) but also that the “U-shape” of the LEG vs Fermi surface angle, observed for underdoped samples, is not artificial. This implies the $`d`$-wave symmetry of the superconducting and pseudo-gaps but with essential contributions from the higher harmonics of one or both gaps which can be a result of the long range pairing interaction.

###### Acknowledgements.

We acknowledge the stimulating discussions with O. K. Andersen, M. Eschrig, E. Laptev, and C. Varma.


text:
# Chiral Limit of Staggered Fermions at Strong Couplings: A Loop Representation

## 1 INTRODUCTION

Understanding the chiral limit of lattice gauge theories with fermions still presents a formidable challenge for numerical work. Conventional algorithms like the hybrid Monte-Carlo or multi-boson techniques encounter difficulties in this limit. In this article we discuss how a new loop representation helps one to solve the chiral limit of strongly coupled $`U(N)`$ gauge theories interacting with staggered fermions in two dimensions.

The chiral limit of strongly coupled gauge theories involving staggered fermions in two dimensions contain some unresolved questions. It is well known that staggered fermions preserve a $`U(1)`$ chiral symmetry at finite lattice spacings. Using large $`N`$ (color)  and large $`d`$ (dimension)  approximations, it is possible to show that this chiral symmetry breaks spontaneously at strong couplings. There are also some rigorous results in four dimensions . As a consequence, at least in $`d3`$ the model is expected to contain a single Goldstone boson (pion) which is exactly massless in the chiral limit. What happens in two dimensions? Numerical simulations suggest the existence of a massless pion . This is somewhat surprising since in $`d=2`$ the Mermin-Wagner theorem forbids the breaking of a continuous symmetry and the pion cannot be a Goldstone boson. On the other hand in two dimensions we know that a $`U(1)`$ symmetry is special. Kosterlitz and Thouless have shown that a theory with this symmetry can contain both a massive phase and a massless phase separated by a transition. So an interesting question is whether the strongly coupled $`U(N)`$ lattice gauge theory with staggered fermions contains a massless pion in the chiral limit or not and if it does is the long distance physics in the same universality of the 2d X-Y model.

## 2 LOOP REPRESENTATION

About two decades ago, it was shown that a $`U(N)`$ gauge theory involving a single flavor of staggered fermions in the strong coupling limit in any dimension is equivalent to a statistical mechanics of monomers and dimers . The partition function is given by

$$Z=\underset{[n,b]}{}\underset{x,\mu }{}\frac{(Nb_{x,\mu })!}{b_{x,\mu }!N!}\underset{x}{}\frac{N!}{n_x!}m^{n_x}$$
(1)
where $`x`$ represents a lattice point in the $`d`$ dimensional lattice, $`\mu =1,2,\mathrm{}d`$ the direction, $`\widehat{\mu }`$ the unit vector in the corresponding direction, $`n_x=0,1,\mathrm{},N`$ is the number of monomers located on the site $`x`$ and $`b_{x,\mu }=0,1,2,\mathrm{},N`$ is the number of dimers located on the bond connecting the site $`x`$ and the nearest neighbor site $`x+\widehat{\mu }`$ in the $`\mu `$ direction. We will assume $`b_{x,\mu }b_{x+\widehat{\mu },\mu }`$ for simplicity. The configurations $`[n,b]`$ that contribute to the partition function are constrained by the relation

$$n_x+\underset{\mu }{}b_{x,\mu }+b_{x,\mu }=N$$
(2)
on each site $`x`$.

In the chiral limit $`n_x`$ is forced to be zero and it is difficult to design a local algorithm for the partition function that obeys the constraint. Perhaps for this reason the chiral limit remained unexplored. In order to avoid this problem we rewrite the partition function given in (1) in terms of new “loop” variables. For the case of $`N=1`$ this was recently explained in . We first shade plaquettes of the two dimensional lattice as in a chess board. We then introduce new (dashed) dimers $`d_{x,\mu }=0,1,2\mathrm{}N`$ on each bond such that each corner $`x`$ of any shaded plaquette satisfies

$$b_{x,\mu }+d_{x,\mu }+b_{x,\nu }+d_{x,\nu }=N$$
(3)
where $`x+\widehat{\mu }`$ and $`x+\widehat{\nu }`$ are the two nearest neighbor sites of $`x`$ also belonging to the same shaded plaquette. In (3) $`\mu `$ and $`\nu `$ can be negative if the direction of the nearest neighbor site is backward. Of course there are many ways to introduce the $`d_{x,\mu }`$ bonds that satisfy (3). It is easy to check that for a shaded plaquette with bonds $`b1,b2,b3,b4`$ in a cyclic order, the total number of ways one can introduce the dashed bonds obeying the constraint (3) is $`(N+1\text{max}(b1,b3)\text{max}(b2,b4))`$. Thus in order to keep the partition function the same we re-weight each extended configuration with $`1/(N+1\text{max}(b1,b3)\text{max}(b2,b4))`$. In this construction it is easy to see that each site touches $`Nn_x`$ solid (original) dimers and $`N+n_x`$ dashed dimers. When $`n_x=0`$ one can at random connect every solid dimer to a unique dashed dimers in $`N!`$ ways. In the chiral limit this results in loops made up of an alternate sequence of solid and dashed dimers. In figure 1 an $`N=2`$ configuration of solid dimers on a $`2\times 2`$ lattice is extended to a loop configuration using the above rules.

An update of the system involves flipping of a loop which means converting all solid dimers into dashed dimers within the loop and vice versa. Typically, a point is selected at random and the loop connected to it is flipped using a Metropolis decision. Since there are $`(b+d)!/(d!b!)`$ number of equivalent dimer configurations with $`b`$ solid and $`d`$ dashed dimers on a given bond, it is necessary to take this factor into account in the Metropolis decision to ensure detailed balance.

## 3 CHIRAL SUSCEPTIBILITY

To address the question of the existence of a massless pion in the chiral limit one can measure the chiral susceptibility $`\chi `$ as a function of the lattice size $`L`$. Since it is the first derivative of the chiral condensate with respect to the fermion mass it should diverge as a power of $`L`$ if there is a massless pion in the theory. For example if the chiral symmetry is spontaneously broken we must have $`\chi L^d`$ where $`d`$ is the dimensions of the lattice. In the case of two dimensions where we expect chiral symmetry to remain unbroken, a massless pion leads to a divergence of the form $`\chi L^(2\eta )`$ where the KT prediction is that $`0\eta 0.25`$. Another example is the continuum two flavor Schwinger model where $`\eta =1`$.

To measure the chiral susceptibility, we introduce two monomers into each loop configuration in all possible ways. Since the loops are made up of an alternate sequence of solid and dashed bonds, introducing a monomer introduces a defect in the pattern. However, a single defect is impossible in a loop. Thus, in a two monomer configuration both monomers must be on the same loop. It is easy to see that a loop with two monomers can be generated by partially flipping a loop. Using this observation we find the chiral susceptibility by summing the weights of all the new configurations generated during the loop update. We have checked that our measurement of the susceptibility produces accurate answers on small lattices. Table 1, compares our results to exact calculations.

When we apply this algorithm to larger lattices we find that small clusters often get updated and occasionally even loops of large sizes get flipped. Our results for the chiral susceptibility for various lattice sizes is plotted in figure 2.

We can draw three conclusions from the plot. First, as a function of $`L`$ there is no evidence of a divergence for any $`N`$. Thus the pion is always massive. Second, since the value of $`L`$ where the susceptibility begins to level gives a rough indication of the inverse of the pion mass, the pion mass increases with $`N`$. Finally, the value of the chiral susceptibility diverges with $`N`$, signaling the formation of a condensate. However, our results suggest that this is not driven by long distance fluctuations.

Since the talk was given we have continued to study our algorithm more carefully. It now appears to us that there are configurations which we sample exponentially rarely as the volume increases but whose contribution is still finite to the final result. This problem begins to show up mildly around $`L=32`$ for $`N=1`$ through very large fluctuations. Of course we had collected very very large statistics to be able to control our error bars. However, we think the error analysis is not reliable since we may have uncontrolled auto-correlation times. For this reason we believe the error bars in figure 2 and in the graphs presented in  are incorrect. This also invalidates the conclusions made above and in . Thus, although our representation is a useful first step and perhaps gives a qualitative picture for small volumes more work is needed before we can answer questions we raised in the first section.


text:
# The transfer matrix with Kogut-Susskind fermions

## 1 Introduction

In the investigation of some features of QCD at finite baryonic density  we must be able to write the baryon charge operator in the transfer matrix formalism. Therefore we must express the transfer matrix in terms of quark-antiquark creation-annihilation operators. Such a formulation has been constructed by by Lüscher  for Wilson fermions, but we have not found in the literature anything similar for Kogut-Susskind fermions. Actually we know only one paper on the subject , but in this work the transfer matrix is constructed in the spin-diagonal basis, where it is found necessary to introduce auxiliary fields whose relation to quark-antiquarks is not transparent. Another feature, possibly related to the above, might be a difficulty in the way of the applications we are interested in: the transfer matrix is linear in the creation-annihilation operators of quark-antiquarks, rather than exponential, so that the hamiltonian contains the logarithm of these operators.

We have therefore tackled the problem in the flavour basis, and we have extended the formalism developed by Lüscher for Wilson fermions to the case of Kogut-Susskind under the condition that the elementary gauge variables be attached to the links of the blocks rather than of the lattice. Let us explain the reason of this restriction. In the Conclusion we will point out its consequence.

Lüscher prooved the reflection positivity  of the transfer matrix in the gauge $`U_0=11`$. This gauge fixing is too strong, but Menotti and Pelissetto  have extended the proof to the gauge $`U_011`$, where $`U_0=11`$ with the exception of a single time slice. If for Kogut-Susskind fermions we define the gauge fields on the links of the blocks, nothing changes in the Menotti-Pelissetto’ proof which we have used in our construction. But if we associate the gauge fields to the links of the lattice, their gauge fixing becomes $`U_0(n_0)U_0(n_0+1)11`$ (in standard notation which is however spelled out below). This has a nontrivial Faddeev-Popov determinant which has not been evaluated but must be taken into account to establish the reflection positivity of the transfer matrix.

In this paper we will ignore the pure gauge field part of the action, the partition function and the transfer matrix. Its role in fact is the same as with Wilson fermions, and it can be included in the same way as in Lüscher’s paper.

## 2 The transfer matrix

The quark field $`q`$ in the flavour basis is defined on hypecubic blocks of coordinates $`n_\mu `$ defined in terms of the sites coordinates $`x_\mu `$

$$x_\mu =2n_\mu +\eta _\mu ,x_\mu =0,\mathrm{},2N_\mu 1,n_\mu =0,\mathrm{},N_\mu 1,\eta _\mu =0,1.$$
(1)
It carries Dirac indices $`\alpha `$, flavour indices $`f`$ and color indices $`c`$, sometimes comprehensively denoted in the sequel by $`I`$.

In the transfer matrix formalism most often one has to do with quantities at a given (euclidean) time $`n_0`$. For this reason we use a summation convention over spatial coordinates $`𝐧_j`$ and intrinsic indices $`I`$ at fixed time. So for instance we will write

$$\overline{q}(m_0)Q(m_0,n_0)q(n_0)=\underset{𝐦,𝐧,I,J}{}\overline{q}_{𝐦,I}(m_0)Q_{𝐦,I;𝐧,J}(m_0,n_0)q_{𝐧,J}(n_0),$$
(2)
where $`Q`$ is the quark matrix. It is a function of time dependent link operators $`U_\mu (n_0)`$ which have the standard Wilson variables $`U_\mu (m_0,𝐦)`$ as spatial matrix elements

$$\left(U_\mu (m_0)\right)_{𝐦,𝐧}=\delta _{𝐦,𝐧}U_\mu (m_0,𝐦).$$
(3)
In this notation the quark action can be written

$$S_q=\underset{m_0,n_0}{}\overline{q}(m_0)Q(m_0,n_0)q(n_0).$$
(4)
Its explicit expression  is

$$S_q=\overline{q}\left\{M1111+2K\left[\gamma _\mu 11_\mu \gamma _5t_5t_\mu \delta _\mu \right]\right\}q,$$
(5)
where the lattice spacing is equal to 1, $`M`$ and $`K`$ are the mass and the hopping parameters and $`_\mu `$ and $`\delta _\mu `$ are the derivatives

$`\left(_\mu f\right)(n)`$ $`=`$ $`{\displaystyle \frac{1}{4}}\left(f(n+2\mu )f(n2\mu )\right)`$
$`\left(\delta _\mu f\right)(n)`$ $`=`$ $`{\displaystyle \frac{1}{4}}\left(f(n+2\mu )+f(n2\mu )2f(n)\right).`$ (6)
In the tensor product the $`\gamma `$-matrices act on Dirac indices, while the $`t`$-matrices act on flavor indices and are given by

$$t_\mu =\gamma _\mu ^T.$$
(7)
The minimal gauge interaction (that where the link variables join the quark fields along the shortest path) can be introduced in the standard way by means of the block link variables $`U_\mu `$.

Now we want to identify the components of the quark field which can be related to quarks and antiquarks. They propagate respectively forwards and bakwards in time and are interchanged by charge conjugation. For this purpose we define the projection operators

$$P_\mu ^{(\pm )}=\frac{1}{2}\left[1111\pm \gamma _\mu \gamma _5t_5t_\mu \right].$$
(8)
We can thus rewrite the quark matrix singling out the terms which are not diagonal in time in the form

$`Q(m_0,n_0)`$ $`=`$ $`K\left[\gamma _0P_0^{()}U_0(m_0)\delta _{n_0,m_0+1}\gamma _0P_0^{(+)}U_0^{(+)}(n_0)\delta _{n_0,m_01}\right]`$ (9)
$`+N(n_0)\delta _{m_0,n_0},`$
where

$`N(n_0)`$ $`=`$ $`M\mathrm{\hspace{0.17em}1}111+K\gamma _5t_5t_0+K{\displaystyle \underset{j=1}{\overset{3}{}}}\{\gamma _5t_5t_j`$ (10)
$`+\gamma _j[P_j^{()}U_j(n_0)T_j^{(+)}P_j^{(+)}T_j^{()}U_j^+(n_0)]\}.`$
The translation operators $`T_j^{(\pm )}`$ have matrix elements

$$\left(T_j^{(\pm )}\right)_{𝐦,𝐧}=\delta _{𝐧,𝐦\pm e_j}$$
(11)
with

$$\left(e_j\right)_k=\delta _{j,k}.$$
(12)
We see that the components

$$q^{(\pm )}=P_0^{(\pm )}q$$
(13)
of the quark field propagate backwards/forwards respectively. It remains to check that they are properly related by charge conjugation. A possible definition of this operator, with Lüscher’s convention for the $`\gamma `$-matrices is

$$𝒞=\gamma _0\gamma _2t_0t_2$$
(14)
whose action on the quark field is

$`q`$ $`=`$ $`𝒞^1\overline{q}^{}`$
$`\overline{q}`$ $`=`$ $`𝒞^Tq^{}.`$ (15)
Therefore charge conjugation

$$q^{(+)}=\gamma _2t_0t_2\left(q^{()}\right)^{}$$
(16)
amounts to interchange $`q^{(\pm )}`$ with one another and to rotate the flavours. This action is accettable insofar as the flavour is physically unobservable.

Finally following the Lüscher’s procedure we define the transfer matrix $`\widehat{𝒯}_q`$ in terms of creation and annihilation operators related to $`q^{(\pm )}`$ and show that

$$Z_q=Tr\widehat{𝒯}_q=[d\overline{q}dq]\mathrm{exp}(16S_q).$$
(17)
The factor of 16 in front of $`S_q`$ is the volume of the block. We construct $`\widehat{𝒯}_q`$ in terms of the auxiliary operator

$$\widehat{T}_q(n_0)=\mathrm{exp}\left(\widehat{y}N(n_0)\widehat{x}\right)$$
(18)
according to the ordered product

$$\widehat{𝒯}_q=𝒥\underset{n_0=0}{\overset{N_01}{}}\left(\widehat{T}_q(n_0)\right)^+\widehat{T}_q(n_01).$$
(19)
In the above equations $`𝒥`$ is the jacobian of a transformation which will be defined later, and the operator $`\widehat{T}_q(n_0)`$ depends on the time $`n_0`$ only throu the dependence on it of the gauge fields. Because $`N`$ is hermitean, $`\widehat{𝒯}_q`$ is positive definite. The above expression of $`\widehat{T}_q(n_0)`$ is valid in the gauge $`U_0=11`$, assumed by Lüscher to prove reflection positivity. As already said this gauge fixing is too strong, but to lighten the formalism we will nevertheless put $`U_0=11`$ and we will reinstate $`U_0`$ at the end. The reader can check that keeping $`U_0`$ in the intermediate steps one arrives at the same result.

To prove Eq. (17) we must transform the trace into a Berezin integral. To this end we introduce between the factors in Eq. (19) the identity

$$11=[dx^+dxdy^+dy]\mathrm{exp}(x^+xy^+y)|xy><xy|,$$
(20)
where the basis vectors

$$|xy>=|\mathrm{exp}(x\widehat{x}^+y\widehat{y}^+)>$$
(21)
are coherent states and the $`x,y`$ are Grassmann variables. We then use the following equations. First

$$<x|x^{}>=\mathrm{exp}\left(x^+x^{}\right).$$
(22)
Second, if $`B=B(\widehat{x}^+)`$ and $`C=C(\widehat{x})`$ are operators which depend on $`\widehat{x}^+,\widehat{x}`$ only

$$<x|BC|x^{}>=B(x^+)<x|x^{}>C(x^{}).$$
(23)
By means of these equations we can evaluate the kernel of the transfer matrix

$`<x(n_0),y(n_0)|\widehat{𝒯}_q(n_0)|x(n_01),y(n_01)>=𝒥^1\mathrm{exp}\left(x^+(n_0)N(n_0)y^+(n_0)\right)`$
$`\mathrm{exp}(x^+(n_0)x(n_01)+(y^+(n_0)y(n_01))`$
$`\mathrm{exp}\left(y(n_01)N(n_01)x(n_01)\right).`$ (24)
Collecting all the pieces we arrive at the euclidean path integral form of the partition function

$$Z_q=𝒥^1[dx^+dxdy^+dy]\mathrm{exp}S^{}$$
(25)
with the action

$`S^{}`$ $`=`$ $`{\displaystyle \underset{n_0}{}}x^+(n_0)x(n_01)y^+(n_0)y(n_01)+x^+(n_0)x(n_0)+y^+(n_0)y(n_0)`$ (26)
$`x^+(n_0)N(n_0)y^+(n_0)y(n_01)N(n_01)x(n_01).`$
The change of variables

$$x=4\sqrt{K}q^{(+)},y^+=4\sqrt{K}q^{()},$$
(27)
whose jacobian is the function $`𝒥`$ introduced in Eq.(19) transforms $`S^{}`$ into $`S_q`$ thus proving Eq.(17).

## 3 Conclusion

We have extended the Lusher’s formalism for Wilson fermions to the case of Kogut-Susskind in the flavour basis under the condition that the gauge fields be associated to the links of the blocks. This restriction is due to the use of the Menotti-Pelissetto proof of reflection positivity which we used in the construction.

This result can be sufficient for our purposes  but leaves unsolved the case of the spin-diagonal basis. In fact we cannot get the spin-diagonal formulation by a simple change of basis . There is a formal difficulty and a practical complication in the way. To this end one should associate the gauge fields to the links of the lattice and this requires an extension of the proof of reflection positivity in the presence of the Faddeev-Popov determinant in the gauge $`U_0(n_0)U_0(n_0+1)11`$. Assuming this to be done, however, we would get an awkward result because a minimal gauge coupling in one basis becomes cumbersome in the other . One should then start from the minimal gauge coupling in the spin-diagonal basis, find the coupling resulting in the flavour basis, and show reflection positivity with such coupling. Both shortcomings are hopefully surmountable, but the problem requires further study and we think that deserves a wider interest than our specific motivation.

The alternative is to construct the transfer matrix directly in the spin-diagonal basis in a form convenient for the use at finite baryon density.