text:
# 1 Introduction

## 1 Introduction

The 1 GeV energy region is a particularly challenging domain. On one side, it is far below from the perturbative QCD regime and, on the other, Chiral Perturbation Theory (ChPT) is not expected to make reliable predictions at these energy values where resonance effects are known to be present. Among the latter, those proceeding by the exchange of the scalar resonances $`f_0`$(980) and $`a_0`$(980) should dominate their respective channels. The controversial nature of these resonances  and the poor knowledge of their properties  adds then further complexity (and interest) to this 1 GeV energy region. Indeed, several proposals have been suggested along the years concerning the constitution of these scalars as complex $`qq\overline{q}\overline{q}`$ states , $`K\overline{K}`$ molecules  or ordinary $`q\overline{q}`$ mesons .

The Novosibirsk CMD-2 and SND Collaborations have reported very recently, among others, the branching ratio and the $`\pi ^0\eta `$ invariant mass distribution for the $`\varphi \pi ^0\eta \gamma `$ decay. For the branching ratio, the CMD-2 Collaboration reports $`B(\varphi \pi ^0\eta \gamma )=(0.90\pm 0.24\pm 0.10)\times 10^4`$ , while the SND result is, consistently, $`B(\varphi \pi ^0\eta \gamma )=(0.88\pm 0.14\pm 0.09)\times 10^4`$ . The observed invariant mass distribution shows a significant enhancement at large $`\pi ^0\eta `$ invariant mass that, according to Refs. , could be interpreted as a manifestation of a sizeable contribution of the $`a_0(980)\gamma `$ intermediate state. This and other radiative $`\varphi `$ decays are also expected to be intensively investigated at the Frascati $`\varphi `$-factory DA$`\mathrm{\Phi }`$NE .

On the theoretical side, the $`VP^0P^0\gamma `$ decays have been considered by a number of authors . In particular, it has been shown that the intermediate vector meson contributions to $`\varphi \pi ^0\eta \gamma `$ lead to a small $`B(\varphi \pi ^0\eta \gamma )_{\mathrm{VMD}}=5.4\times 10^6`$ , whereas a chiral loop model closely linked to standard ChPT predicts $`B(\varphi \pi ^0\eta \gamma )_\chi =3.0\times 10^5`$ . Needless to say, the scalar resonance effects and, in particular, the resonance pole associated to the $`a_0`$(980) were not contemplated in these two approaches. The recent experimental data from Novosibirsk —for both the branching ratio and the $`\pi ^0\eta `$ invariant mass spectrum showing an enhancement around the $`a_0`$(980) mass— seem therefore to disfavour these predictions based on vector meson exchange and/or a simple extrapolation of ChPT ideas.

If we rely on the resonance picture, it is clear that the $`a_0(980)`$ scalar meson —lying just below the $`\varphi `$ mass and having the appropriate quantum numbers— should play an important rôle in the $`\varphi \pi ^0\eta \gamma `$ decay. Several theoretical attempts to describe the effects of scalars in $`\varphi `$ radiative decays have appeared so far. Among others, we would like to refer to the “no structure” model , to the $`K^+K^{}`$ model , where the $`\varphi a_0\gamma `$ amplitude is generated through a loop of charged kaons, and to the chiral unitary approach $`(U\chi PT)`$ , where the decay $`\varphi \pi ^0\eta \gamma `$ occurs through a loop of charged kaons that subsequently annihilate into $`\pi ^0\eta \gamma `$. In the two former cases the scalar resonances are included ad hoc while in the latter they are generated dynamically by unitarizing the one-loop amplitudes.

In this letter, we are mainly interested in incorporating scalar resonances and their pole effects into a ChPT inspired context . While vector and axial-vector resonances can be included in a transparent and successful way, offering some theoretical basis to conventional vector meson dominance (VMD) ideas , the incorporation of scalar resonances has been more ambiguous and less successful up to now . In order to take explicitly into account scalar resonances and their pole effects, we propose to use the linear sigma model (L$`\sigma `$M). This will allow us to take advantage of the common origin of ChPT and the L$`\sigma `$M to improve the chiral loop predictions for $`VP^0P^0\gamma `$ exploiting the complementarity of both approaches for these specific processes. On one side, ChPT is the established theory of the pseudoscalar interactions at low energy. However, it is not reliable at energies of a typical vector meson mass and, as just stated, scalar resonance poles are not explicitly included. As a consequence, ChPT inspired loop models can give rough estimates for $`B(VP^0P^0\gamma )`$ but will hardly be able to reproduce the observed enhancements in the invariant mass spectra. On the other side, the L$`\sigma `$M is a much simpler model dealing similarly with pseudoscalar interactions but incorporating scalar resonances in a systematic and definite way. Thanks to this, the L$`\sigma `$M should be able to reproduce the resonance peaks in the spectra and, although it does not provide a systematical framework for the pseudoscalar meson physics, this model could be of relevance in describing the scalar resonances when linked to a well established ChPT context. In order to show in detail the proposed framework, we will focus our attention on the $`\varphi \pi ^0\eta \gamma `$ decay mode. Other decay modes are somewhat more involved and will be analyzed in forthcoming work.

## 2 $`\varphi \pi ^0\eta \gamma `$ and chiral loops

The vector meson initiated $`VP^0P^0\gamma `$ decays cannot be treated in strict Chiral Perturbation Theory (ChPT). This theory has to be extended to incorporate on-shell vector meson fields. At lowest order, this may be easily achieved by means of the $`𝒪(p^2)`$ ChPT Lagrangian:

$$_2=\frac{f^2}{4}D_\mu U^{}D^\mu U+U^{}\chi +\chi ^{}U,$$
(1)
where $`f=f_\pi =92.4`$ MeV at this order, $`U=\mathrm{exp}(i\sqrt{2}P/f)`$ with $`P`$ being the usual pseudoscalar nonet matrix, and $`\chi =2B_0`$ with $`=\text{diag}(m_u,m_d,m_s)`$. The covariant derivative, now enlarged to include vector mesons, is defined as $`D_\mu U=_\mu UieA_\mu [Q,U]ig[V_\mu ,U]`$, with $`Q=\text{diag}(2/3,1/3,1/3)`$ being the quark charge matrix and $`V_\mu `$ the additional matrix containing the nonet of ideally mixed vector meson fields. The diagonal elements of $`V`$ are $`(\rho ^0+\omega )/\sqrt{2},(\rho ^0+\omega )/\sqrt{2}`$ and $`\varphi `$, thus following the same conventional normalization as for the pseudoscalar nonet matrix $`P`$.

There is no tree-level contribution from this Lagrangian to the $`\varphi \pi ^0\eta \gamma `$ amplitude and at the one-loop level one needs to compute the set of diagrams shown in Ref. . A straightforward calculation leads to the following finite amplitude for $`\varphi (q^{},ϵ^{})\pi ^0(p)\eta (p^{})\gamma (q,ϵ)`$ (see Ref.  for further details):

$$\begin{array}{cc}\hfill 𝒜(\varphi \pi ^0\eta \gamma )_\chi & =\frac{eg}{2\pi ^2m_{K^+}^2}(ϵ^{}ϵq^{}qϵ^{}qϵq^{})L(m_{\pi ^0\eta }^2)\hfill \\ & \times 𝒜(K^+K^{}\pi ^0\eta )_\chi ,\hfill \end{array}$$
(2)
where $`ϵ^{}ϵq^{}qϵ^{}qϵq^{}`$ makes the amplitude Lorentz- and gauge-invariant, $`m_{\pi ^0\eta }^2s(p+p^{})^2=(q^{}q)^2`$ is the invariant mass of the final pseudoscalar system and $`L(m_{\pi ^0\eta }^2)`$ is the loop integral function defined as

$$\begin{array}{cc}\hfill L(m_{\pi ^0\eta }^2)& =\frac{1}{2(ab)}\frac{2}{(ab)^2}\left[f\left(\frac{1}{b}\right)f\left(\frac{1}{a}\right)\right]\hfill \\ & +\frac{a}{(ab)^2}\left[g\left(\frac{1}{b}\right)g\left(\frac{1}{a}\right)\right],\hfill \end{array}$$
(3)
where

$$\begin{array}{c}f(z)=\{\begin{array}{cc}\left[\mathrm{arcsin}\left(\frac{1}{2\sqrt{z}}\right)\right]^2\hfill & z>\frac{1}{4}\hfill \\ \frac{1}{4}\left(\mathrm{log}\frac{\eta _+}{\eta _{}}i\pi \right)^2\hfill & z<\frac{1}{4}\hfill \end{array}\hfill \\ g(z)=\{\begin{array}{cc}\sqrt{4z1}\mathrm{arcsin}\left(\frac{1}{2\sqrt{z}}\right)\hfill & z>\frac{1}{4}\hfill \\ \frac{1}{2}\sqrt{14z}\left(\mathrm{log}\frac{\eta _+}{\eta _{}}i\pi \right)\hfill & z<\frac{1}{4}\hfill \end{array}\hfill \end{array}$$
(4)
and $`\eta _\pm =\frac{1}{2}(1\pm \sqrt{14z})`$, $`a=\frac{m_\varphi ^2}{m_{K^+}^2}`$ and $`b=\frac{m_{\pi ^0\eta }^2}{m_{K^+}^2}`$. The coupling constant $`g`$ comes from the strong amplitude $`𝒜(\varphi K^+K^{})=gϵ^{}(p_+p_{})`$ with $`|g|=4.59`$ to agree with $`\mathrm{\Gamma }(\varphi K^+K^{})_{\mathrm{exp}}=2.19`$ MeV. The latter is the part beyond standard ChPT which we have fixed phenomenologically. The four-pseudoscalar amplitude is instead a standard ChPT amplitude111 $`𝒜(K^+K^{}\pi ^0\eta _8)_\chi =\frac{\sqrt{3}}{4f_\pi ^2}\left(m_{\pi ^0\eta }^2\frac{4}{3}m_K^2\right)`$ if only the $`\eta _8`$ contribution is taken into account as in Ref. . which is found to depend linearly on the variable $`s=m_{\pi ^0\eta }^2`$:

$$𝒜(K^+K^{}\pi ^0\eta )_\chi =\frac{1}{\sqrt{6}f_\pi ^2}\left(m_{\pi ^0\eta }^2\frac{10}{9}m_K^2+\frac{1}{9}m_\pi ^2\right).$$
(5)
In the calculation of the decay amplitudes (2) and (5) we have introduced $`\eta `$-$`\eta ^{}`$ mixing effects. As it is well known, a rigorous and general extension of $`SU(3)`$ ChPT to include the ninth pseudoscalar meson $`\eta _0`$ is not straightforward and requires the introduction of new terms in the chiral Lagrangian . However, if one relies on classical arguments based on nonet symmetry, a phenomenologically successful description of the $`\eta `$-$`\eta ^{}`$ system is achieved . The $`\eta `$-$`\eta ^{}`$ mixing angle is then found to be compatible with $`\theta _P=\mathrm{arcsin}(1/3)19.5^{}`$, quite in agreement with recent phenomenological estimates . In Sect. 3, it will be shown that this choice for the $`\eta `$-$`\eta ^{}`$ mixing angle greatly simplifies the calculation of the $`\varphi \pi ^0\eta \gamma `$ amplitude in the L$`\sigma `$M, and reduces up to a minimum the number of free-parameters.

The invariant mass distribution for the $`\varphi \pi ^0\eta \gamma `$ decay is predicted to be given by the following spectrum (see Fig. 1):

$$\begin{array}{cc}\hfill \frac{d\mathrm{\Gamma }(\varphi \pi ^0\eta \gamma )_\chi }{dm_{\pi ^0\eta }}& =\frac{\alpha }{192\pi ^5}\frac{g^2}{4\pi }\frac{m_\varphi ^4}{m_{K^+}^4}\frac{m_{\pi ^0\eta }}{m_\varphi }\left(1\frac{m_{\pi ^0\eta }^2}{m_\varphi ^2}\right)^3\sqrt{12\frac{m_{\pi ^0}^2+m_\eta ^2}{m_{\pi ^0\eta }^2}+\left(\frac{m_\eta ^2m_{\pi ^0}^2}{m_{\pi ^0\eta }^2}\right)^2}\hfill \\ & \times |L(m_{\pi ^0\eta }^2)|^2|𝒜(K^+K^{}\pi ^0\eta )_\chi |^2.\hfill \end{array}$$
(6)
Integrating Eq. (6) over the whole physical region one obtains for the branching ratio:

$$B(\varphi \pi ^0\eta \gamma )_\chi =0.47\times 10^4.$$
(7)
As expected, Fig. 1 shows that our chiral loop approach gives a reasonable prediction for the lower part of the spectrum but fails to reproduce the observed enhancement in its higher part, where $`a_0(980)`$-resonance effects (ignored up to this point of our approach) should manifest. As a consequence, the predicted branching ratio turns out to be below the experimental value by about a factor of 2.

## 3 Improved approach to $`\varphi \pi ^0\eta \gamma `$

To analyze the scalar resonance effects in the $`VP^0P^0\gamma `$ decay amplitudes, the linear sigma model (L$`\sigma `$M)  will be shown to be particularly appropriate. It is a well-defined $`U(3)\times U(3)`$ chiral model which incorporates ab initio both the nonet of pseudoscalar mesons together with its chiral partner, the scalar mesons nonet. Recently, the model has been resurrected as a framework to study the implications of chiral symmetry for the controversial scalar sector of QCD, and some variations of the basic L$`\sigma `$M Lagrangian have shown to be phenomenologically rather successful .

In this context, the $`VP^0P^0\gamma `$ decays proceed through a loop of charged pseudoscalar mesons emitted by the initial vector. Due to the additional emission of a photon, these charged pseudoscalar pairs with the initial $`J^{PC}=1^{}`$ quantum numbers can rescatter into $`J^{PC}=0^{++}`$ pairs of charged or neutral pseudoscalars. The scalar resonances are expected to play an essential rôle in this rescattering process and the L$`\sigma `$M seems mostly appropriate to fix the corresponding amplitudes.

Several simplifications happen when one considers the $`\varphi \pi ^0\eta \gamma `$ decay mode. As in the analysis of Sect. 2, contributions from charged pions in the loops are highly suppressed because they involve the isospin violating and OZI–rule forbidden $`\varphi \pi \pi `$ coupling; hence, the dominant contributions arise exclusively from loops of charged kaons. The subsequent rescattering of these charged kaon pairs into the final $`J^{PC}=0^{++}`$ $`\pi ^0\eta `$ state is then quite simple. Indeed, the L$`\sigma `$M amplitude for $`K^+K^{}\pi ^0\eta `$ contains a contact term, a term with an $`a_0`$ exchanged in the $`s`$-channel and two terms with a $`\kappa `$ (i.e. the strange $`I=1/2`$ scalar resonance) exchanged in the $`t`$\- and $`u`$-channels. However, the latter $`\kappa `$-exchange contributions are absent for an $`\eta `$-$`\eta ^{}`$ mixing angle $`\theta _P=\mathrm{arcsin}(1/3)19.5^{}`$ since the $`g_{\kappa K\eta }`$ coupling constant appearing in one of the $`\kappa `$ vertices vanishes. The calculation is then reduced to the diagrams shown in Fig. 2, and it is thus much simpler for this particular (but phenomenologically acceptable ) choice of the $`\eta `$-$`\eta ^{}`$ mixing angle. Moderate departures from this value translate into weak $`g_{\kappa K\eta }`$ couplings222The coupling constant $`g_{\kappa K\eta }`$ is proportional to $`\mathrm{sin}\left[\theta _P\mathrm{arcsin}(1/3)\right]`$ and its dependence on the $`\eta `$-$`\eta ^{}`$ mixing angle for values around $`\theta _P=\mathrm{arcsin}(1/3)19.5^{}`$ is soft. appearing in the $`\varphi \pi ^0\eta \gamma `$ amplitude. Their effects seem to be small (see below) for the present and expected levels of experimental accuracy and for this kind of processes governed by poorly known scalar resonances. More importantly, the absence of $`\kappa `$ contributions makes our predictions more definite and solid since we avoid one of the major uncertainties affecting the scalar nonet dynamics, namely, the mass of its strange members. Indeed, recent analysis by various authors require a light $`\kappa `$(900) , while other authors deny the existence of a low mass pole  and identify the PDG $`K_0^{}(1430)`$ state with the strange member of the lowest lying scalar nonet .

A straightforward calculation of the $`\varphi \pi ^0\eta \gamma `$ decay amplitude leads to an expression identical to that in Eq. (2) but with the four-pseudoscalar amplitude now computed in a L$`\sigma `$M context. In this case, the amplitude is just

$$𝒜(K^+K^{}\pi ^0\eta )_{\text{L}\sigma \text{M}}=g_{K^+K^{}\pi ^0\eta }\frac{g_{a_0K^+K^{}}g_{a_0\pi ^0\eta }}{D_{a_0}(m_{\pi ^0\eta }^2)},$$
(8)
where $`D_{a_0}(m_{\pi ^0\eta }^2)`$ is the $`a_0`$ propagator and the coupling constants are333See Ref.  for a detailed calculation of these coupling constants in the L$`\sigma `$M.

$$\begin{array}{c}g_{K^+K^{}\pi ^0\eta }=\frac{g_{a_0\pi ^0\eta }}{2f_K}=\sqrt{\frac{2}{3}}\frac{m_{a_0}^2m_\eta ^2}{2f_Kf_\pi },\hfill \\ g_{a_0K^+K^{}}=\frac{m_{a_0}^2m_K^2}{2f_K}.\hfill \end{array}$$
(9)
This amplitude can then be rewritten as

$$𝒜(K^+K^{}\pi ^0\eta )_{\text{L}\sigma \text{M}}=\frac{1}{\sqrt{6}f_Kf_\pi }(m_{\pi ^0\eta }^2m_K^2)\times \frac{m_\eta ^2m_{a_0}^2}{D_{a_0}(m_{\pi ^0\eta }^2)}.$$
(10)
We would like to make a few remarks on the four-pseudoscalar amplitude in Eq. (10) and compare it with the ChPT amplitude in Eq. (5):

* For $`m_{a_0}\mathrm{}`$ and ignoring $`SU(3)`$-breaking in the pseudoscalar masses and decay constants, the L$`\sigma `$M amplitude (10) reduces to the ChPT one (5). The former consists of a constant four-pseudoscalar vertex plus a second term whose $`s`$-dependence is generated by the $`sm_{\pi ^0\eta }^2`$ piece in the $`a_0`$ propagator $`D_{a_0}(s=m_{\pi ^0\eta }^2)`$, as shown in Eq. (8). Their sum (see Eq. (10)) in the good $`SU(3)`$ and $`m_{a_0}\mathrm{}`$ limits ends up with an amplitude which is linear in $`s=m_{\pi ^0\eta }^2`$ and mimics perfectly the effects of the derivative and massive terms in the ChPT Lagrangian (1) leading respectively to the two terms in the ChPT amplitude (5). This, we believe, is the main virtue of our approach and makes the use of the L$`\sigma `$M reliable at least for amplitudes like ours where $`s`$-channel exchange plays the main rôle and $`(t,u)`$-channel exchange can be ignored.
* The L$`\sigma `$M and ChPT yield slightly different amplitudes in the $`m_{a_0}\mathrm{}`$ limit because of the way $`SU(3)`$-symmetry is broken in the two approaches. In the case of the L$`\sigma `$M , a non $`SU(3)`$ symmetric choice of the vacuum expectation values makes simultaneously $`m_\pi ^2m_K^2`$ and $`f_\pi f_K`$, whereas in ChPT $`m_\pi ^2m_K^2`$ is already present in the lowest order Lagrangian while $`f_\pi f_K`$ is only achieved at higher orders.
* In the $`\varphi \pi ^0\eta \gamma `$ decay, the threshold for $`\pi ^0\eta `$ production is not far from the mass of the $`a_0`$(980) and that makes crucial the incorporation of the $`a_0`$ in an explicit way. Due to the presence of the full propagator $`D_{a_0}(s)`$, as in Eq. (10), such an amplitude —closely linked to that from ChPT and thus expected to be able to account for the lowest part of the $`\pi ^0\eta `$ mass spectrum— should also be able to reproduce the effects of the $`a_0`$ pole at higher $`\pi ^0\eta `$ invariant mass values.
* The need for the $`a_0`$ propagator introduces, however, some uncertainties in our treatment. Indeed, the opening of the $`K\overline{K}`$ channel near the $`a_0`$(980) mass has motivated the use of different expressions for $`D_{a_0}(s)`$. A first possibility consists in using a Breit-Wigner propagator with an energy dependent width (to incorporate the known kinematic corrections):

$$D_{a_0}(s)=sm_{a_0}^2+i\sqrt{s}\mathrm{\Gamma }_{a_0}(s),$$
(11)
where

$$\begin{array}{cc}\hfill \mathrm{\Gamma }_{a_0}(s)& =\frac{g_{a_0\pi ^0\eta }^2}{16\pi \sqrt{s}}\sqrt{\left[1\frac{(m_{\pi ^0}+m_\eta )^2}{s}\right]\left[1\frac{(m_{\pi ^0}m_\eta )^2}{s}\right]}\theta (\sqrt{s}(m_{\pi ^0}+m_\eta ))\hfill \\ & +\frac{g_{a_0K^+K^{}}^2}{16\pi \sqrt{s}}\sqrt{1\frac{4m_{K^+}^2}{s}}\theta (\sqrt{s}2m_{K^+})\hfill \\ & +\frac{g_{a_0K^0\overline{K}^0}^2}{16\pi \sqrt{s}}\sqrt{1\frac{4m_{K^0}^2}{s}}\theta (\sqrt{s}2m_{K^0}).\hfill \end{array}$$
(12)
Another interesting and widely accepted option was proposed by Flatté time ago specifically for the two-channel $`a_0`$ resonance . The relative narrowness of the observed $`\pi \eta `$ peak around 980 MeV is then explained by the action of unitarity and analyticity at the $`K\overline{K}`$ threshold. This amounts to extend the preceding formulae below the $`K\overline{K}`$ threshold to include the now purely imaginary kaon contributions.

Due to these distinct possibilities to deal with the $`a_0`$ propagator, as well as to other differences introduced when implementing and fitting the basic L$`\sigma `$M Lagrangian by several authors, a set of predictions can be obtained for the four-pseudoscalar amplitude (10). In turn, these various amplitudes have to substitute the four-pseudoscalar ChPT amplitude in Eq. (6) to finally obtain the corresponding invariant mass distributions of the $`\varphi \pi ^0\eta \gamma `$ decay mode444Since all the four-pseudoscalar amplitudes we are considering depend only on the variable $`s`$, they factorize out of the loop integration and the structure of Eq. (6) is fully preserved.. Our purpose hereafter is to briefly discuss a few of these treatments in order to show that the observed properties for this specific decay can be accommodated in our ChPT- and L$`\sigma `$M-inspired approach.

We start our discussion along the lines of Ref.  taking for the $`a_0`$ propagator the simple Breit-Wigner prescription in Eq. (11). The use of this propagator for the L$`\sigma `$M amplitude Eq. (10) and its insertion in Eq. (6) predicts the $`m_{\pi ^0\eta }`$ invariant mass spectrum shown by the dotted line in Fig. 3. Integrating over the whole physical region leads to the branching ratio

$$B(\varphi \pi ^0\eta \gamma )_{\text{L}\sigma \text{M}\text{[22]}}=0.80\times 10^4,$$
(13)
in agreement with the experimental data . However, since the simple expression used for the $`a_0`$ propagator implies a large $`a_0`$-width ($`\mathrm{\Gamma }_{a_0\pi \eta }460`$ MeV ), the desired enhancement in the invariant mass spectrum appears in its central part rather than around the $`a_0`$ peak.

This unpleasant feature is easily corrected when turning to the proposal by Törnqvist . Indeed, a Gaussian form factor related to the finite size of physical mesons and depending on the final CM-momentum is introduced to describe the decays of scalar resonances in this approach. As a result, the decay width of $`a_0`$(980) into $`\pi ^0\eta `$ is reduced ($`\mathrm{\Gamma }_{a_0\pi \eta }273`$ MeV ) without affecting that of $`a_0`$(980) into $`K\overline{K}`$. This fact produces an enhancement in the spectrum for the higher values of the $`m_{\pi ^0\eta }`$ invariant mass, as shown by the dashed line in Fig. 3. The integrated branching ratio is then predicted to be

$$B(\varphi \pi ^0\eta \gamma )_{\text{L}\sigma \text{M}\text{[23]}}=0.90\times 10^4,$$
(14)
in good agreement with the experimental data . A possible difficulty of this approach at the phenomenological level is that it predicts an $`\eta `$-$`\eta ^{}`$ mixing angle of $`\theta _P5.0^{}`$ , considerably less negative than the usually accepted value ($`\theta _P`$ between $`20^{}`$ and $`15^{}`$ ) and the value required in our simplified analysis $`\theta _P19.5^{}`$. This allows for an estimate of the typical errors introduced when neglecting $`\kappa `$-exchange in the $`t`$\- and $`u`$-channels as compared to $`a_0`$-exchange in $`s`$-channel. The significant factor is the ratio of coupling constants $`r(g_{\kappa ^\pm K^\pm \eta }g_{\kappa ^\pm K^\pm \pi ^0})/(g_{a_0\pi ^0\eta }g_{a_0K^+K^{}})`$, which vanishes in the ideal situation where $`\theta _P19.5^{}`$ but takes the value $`r1/3`$ for $`\theta _P5.0^{}`$ if one ignores $`SU(3)`$-breaking corrections. To the smallness of $`r1/3`$ one has to add the fact that the amplitude for the almost on-shell $`a_0`$-exchange in the $`s`$-channel is mainly imaginary and does not interfere with the almost real amplitude for off-shell $`\kappa `$-exchange in the $`t`$\- and $`u`$-channels. As a result, the error in Eq. (14) introduced by neglecting this latter contribution can be estimated to be below some 10% even for such unusual values of $`\theta _P5^{}`$.

None of these drawbacks are encountered when turning to the treatment proposed by Shabalin . In the fitting procedure adopted by this author no attempt is made to fix the $`\eta `$-$`\eta ^{}`$ mixing angle within the model. Thanks to this, one minimizes the uncertainties associated with the incorporation of the ninth pseudoscalar meson $`\eta _0`$ via the axial anomaly term. The value of the $`\eta `$-$`\eta ^{}`$ mixing angle is then fixed outside the model to its phenomenologically preferred value $`\theta _P19.5^{}`$. Another relevant feature of Shabalin’s approach is the introduction of the well-known Flatté corrections to the $`a_0(980)`$ propagator. The $`a_0`$-width is then drastically reduced from the uncorrected value $`\mathrm{\Gamma }_{a_0}304`$ MeV to a more acceptable visible width of $`\mathrm{\Gamma }_{a_0}65`$ MeV. With all this information taken from Ref. , our approach predicts the $`m_{\pi ^0\eta }`$ invariant mass spectrum shown by the solid line in Fig. 3 and the integrated branching ratio

$$B(\varphi \pi ^0\eta \gamma )_{\text{L}\sigma \text{M}\text{[24]}}=0.93\times 10^4.$$
(15)
Both the spectrum and the branching ratio are in nice agreement with the experimental data . The fact that Shabalin’s model incorporates the Flatté corrections to the $`a_0`$ resonant shape  has played a substantial rôle in this achievement.

## 4 Conclusions

The main aim of the present letter has been to propose and discuss an amplitude for the radiative $`\varphi \pi ^0\eta \gamma `$ decay exploiting the complementarity between ChPT and L$`\sigma `$M ideas. Thanks to the latter, our amplitude contains the full propagator of the $`a_0`$(980) scalar resonance which dominates the higher part of the $`\pi ^0\eta `$ invariant mass spectrum. In the low invariant mass region, where ChPT is expected to work quite reliably, the proposed amplitude is shown to coincide with that coming from a chiral loop calculation. This, we believe, makes reliable our approach to the $`VP^0P^0\gamma `$ dynamics and, in particular, to the $`\varphi \pi ^0\eta \gamma `$ decay mode for which some simplifying conditions hold and lead to a simple and well-defined amplitude. Then our predictions depend only on a reduced number of parameters which, in principle, can be extracted from independent data. Some of these data refer to scalar meson properties which are not well established and thus affect the accuracy of our predictions although by no means in a drastic way.

We can safely conclude that all the reported properties for the $`\varphi \pi ^0\eta \gamma `$ decay mode can be accommodated in our approach. The branching ratio is predicted to be in the range $`BR(\varphi \pi ^0\eta \gamma )=(0.75`$$`0.95)\times 10^4`$, compatible with the present available data. Similarly, the measured $`\pi ^0\eta `$ invariant mass spectrum is reproduced by our amplitude in a reasonable way. The uncertainties affecting these predictions suggest that further tests and more refined analyses are needed, particularly when the higher accuracy data from ongoing experiments will be available. This should contribute to clarify one of the most controversial aspects of hadron physics: the scalar states around 1 GeV.

## Acknowledgements

We would like to thank A. Farilla for helpful comments and clarifying discussions. Work partly supported by the EEC, TMR-CT98-0169, EURODAPHNE network. Work also partly supported by the CONACyT (project I27604-E). J.L. Lucio M. acknowledges partial financial support from CONACyT and CONCyTEG.


format: text

text:
# Quantum Number Fractionization in Cuprates: U(1) RVB Gauge Theory and Senthil Fisher Constructions

## Abstract

We critically look at Senthil Fisher constructions of quantum number fractionization in cuprates. The first construction, while mathematically correct is of limited relevance to cuprates. The second one has a missing U(1) symmetry. Once this aspect is repaired, it reveals itself as the U(1) RVB gauge theory construction of quantum number fractionization. An approximate local $`Z_2`$ symmetry arises in the spin sector for a reason very different from Senthil and Fisher’s. The charge sector continues to have the local U(1) symmetry and the metallic spin gap phase of the 2d cuprates may, in principle, carry (irrational) fractional charge $`e^{}e`$ obeying Haldane’s fractional exclusion statistics like the 1d models.

In the RVB theory of cupratesthe (s)slave-particle formalism and the U(1) gauge theory development suggested a natural mechanism of quantum number fractionization in undoped and doped cuprates. The spin charge separation simply corresponded to the deconfinement of s-particles carrying a new U(1) charge called RVB charge. As the theory involves matter (a fermi sea of spinons and bosonic holons) interacting with dynamically generated U(1) gauge fields, spinon deconfinement turned out to be a natural (the only ?) possibility even in the 1 and 2 space dimensional t-J models (in the absence of matter the U(1) gauge theory in 1 and 2 space dimensions have only a confined phase) . Soon after the formulation of the RVB gauge theory, a reduction of the local U(1) symmetry to a local $`Z_2`$ symmetry was found leading to a theory of $`Z_2`$ gauge fields coupled to spinon and holon matter utilizing an important identity due to Marston. This approach has remained unexplored. Several authors, Patrick Lee and collaborators in particular have pursued the gauge theory approach in the last one decade.

Recently Senthil and Fisher (SF) have advocated a construction of quantum number fractionization (QNF) $`2ee+e`$ in two different ways. In the most recent approach, SF study the t-J model in the s-particle formalism and split a neutral d-wave spinon pair into two and use it to construct a spinon and a ‘chargon’. We critically analyze this construction and find a missing local U(1) gauge invariance. This invariance is necessary to ensure the crucial double occupancy constraint. We have found a way of correcting this serious difficulty. Once this is done, SF’s construction becomes the same as the U(1) RVB gauge theory construction of QNF . In the earlier approach, a charge ‘$`e`$’ excitation of the charge ‘2e’ fluid in the metallic spin gap phase of the cuprate is invoked. While this in principle is an allowed possibility in the cuprates we will argue that this is only a particular case allowed by the U(1) RVB gauge theory of the t-J model; a general possibility allowed by the U(1) gauge theory is more consistent with the ‘charge gaplessness’ of the metallic spin (pseudo)gap phase.

In the process we also bring out a $`Z_2`$ aspect of the RVB gauge theory that was discussed earlier in the game . A fundamental difference of our $`Z_2`$ gauge symmetry, which is restricted to the spin sector of the doped Mott insulator, from that of SF is also brought out. Having established the connection of SF’s construction to the old U(1) RVB mechanism, we also point out how the local U(1) symmetry of the charge sector (in the metallic spin gap phase) provides a more general QNF and the possibility of fractional charge for holons $`e^{}e`$. In this way we also suggest that the physics of the charge sector is fairly dimension independent. We also point out that in a (gapless) metallic 2d state such as the spin gap phase with extended (in general having power law form factor) spinons and holons it is more meaningful to talk in terms of the Haldane’s exclusion statistics rather than limit oneself to ideas derived from exchange statistics.

First we elaborate how the U(1) RVB gauge theory allows for a unique ‘U(1) string’ representation of the physical spinon and holon, allowing for a general possibility of charge $`e^{}e`$ for the holon. Then show how the recent construction of SF, when corrected, becomes the RVB theory’s U(1) string representation.

In the (s)slave-particle representation of the the t-J model an electron is decomposed as

$$c_{i\sigma }^{}e_is_{i\sigma }^{}$$
(1)
in terms of the ‘bare spinon’ and ‘bare holon’ operators ($`s_{i\sigma }^{}`$,$`e_i^{}`$). The double occupancy constraint takes the form $`_\sigma `$ $`s_{i\sigma }^{}`$$`s_{i\sigma }`$\+ $`e_i^{}e_i=1`$. All physical objects are invariant under the local U(1) gauge transformation applied to the s-particles:

$$(s_{i\sigma }^{},e_i^{})e^{i\theta _i}(s_{i\sigma }^{},e_i^{})$$
(2)
In view of the above transformation property, we say that a bare spinon or a holon carry an unit of U(1) ‘charge’ which we called RVB charge. These constituent s-particles are not physical particles. Physical particles must remain invariant under the local U(1) gauge transformation. Gauge theory suggests that we attach electric field lines to ‘electric charges’ to make them gauge invariant. In our system the field lines that go to infinity (boundary of the sample) can be created using spin singlet ‘U(1)string’ operators such as:

$$U_{li_1}U_{i_1i_2}U_{i_2i_3}\mathrm{}U_l(C)$$
(3)
The link variables $`U_{ij}`$ are chosen from the three types of spin singlet operators: $`U_{ij}=b_{ij}^{}=\frac{1}{\sqrt{2}}(s_i^{}s_j^{}s_i^{}s_j^{})`$, or $`\tau _{ij}=_\sigma `$$`s_{i\sigma }^{}`$$`s_{ȷ\sigma }`$ or$`e_i^{}e_j`$ and their complex conjugates. A simple example of the string operator used in reffor the Mott insulator is $`b_{li_1}b_{i_1i_2}^{}b_{i_2i_3}\mathrm{}..`$ . By construction, under the local U(1) transformation a U(1) string operator transforms as $`U_{li_1}U_{i_1i_2}U_{i_2i_3}\mathrm{}e^{i\theta _l}U_{li_1}U_{i_1i_2}U_{i_2i_3}\mathrm{}`$. The field line operators represent the ‘spin singlet’ disturbance and charge disturbance created in the spin liquid state in the process of QNF , of transporting a partner spin-$`\frac{1}{2}`$ moment to the boundary of the system along the sequence C (not necessarily formed by nearest neighbor lattice sites). In a continuum approximation of the U(1) RVB gauge theory the string operators correspond to: $`e^{i_𝐫^{\mathrm{}}𝐀_{rvb}𝑑𝐥}`$. Here $`𝐀_{rvb}`$ is the dynamically generated U(1) RVB gauge field. That is, schematically a spin one operator or an electron operator undergo the following change during quantum number fractionization:

$`S_l^+`$ $`=`$ $`s_l^{}s_ls_l^{}U_{li_1}U_{i_1i_2}U_{i_2i_3}\mathrm{}U_{i_n\mathrm{}}s_{\mathrm{}}`$ (4)
$`c_{l\sigma }^{}`$ $`=`$ $`e_ls_{l\sigma }^{}e_lU_{li_1}U_{i_1i_2}U_{i_2i_3}\mathrm{}U_{i_n\mathrm{}}s_\mathrm{}\sigma ^{}`$ (5)
The field line operators in combination with the bare s-particles have the full local U(1) gauge invariance and provides the most general localized physical spinon or holon operator characterized by a string C and site l:

$$\zeta _{l\sigma }^s=s_{l\sigma }^{}U_l(C)\text{a}nd\zeta _l^h=e_{l\sigma }^{}U_l(C)$$
(6)
A propagating low energy spinon or holon eigen state of a spin liquid state exhibiting quantum number fractionization (deconfinement) can be created by the operators:

$`\zeta _\sigma ^s(𝐪)`$ $``$ $`{\displaystyle \underset{l,C}{}}e^{i𝐪𝐫_l}\varphi ^s(C,l;𝐪)s_{l\sigma }^{}U_l(C)`$ (7)
$``$ $`{\displaystyle \underset{l}{}}e^{i𝐪𝐫_l}s_{l\sigma }^{}\eta _l^s(𝐪)\text{a}nd`$ (8)
$`\zeta ^h(𝐪)`$ $``$ $`{\displaystyle \underset{l,C}{}}e^{i𝐪𝐫_l}\varphi ^h(C,l;𝐪)e_l^{}U_l(C)`$ (9)
$``$ $`{\displaystyle \underset{l}{}}e^{i𝐪𝐫_l}e_l^{}\eta _l^h(𝐪)`$ (10)
The string amplitudes $`\varphi ^s(C,l;𝐪)`$ and $`\varphi ^h(C,l;𝐪)`$ depends on the form of the Hamiltonian and in general may contain some further quantum numbers such as ‘chirality’ for the spinon and holon.

The operator $`\eta _l^{s,h}=_C\varphi ^{s,h}(C,l;𝐪)U_l(C)`$ creates a quantum delocalized U(1) string excitation starting from a site $`l`$; in a classical picture it is a radially flowing U(1) field lines emanating from a site $`l`$. We call $`\eta _l^s`$ and $`\eta _l^h`$ as the ‘radial field operators’ (RFO) of the deconfined phase for the spinon and holon. By construction the RFO’s transform as $`\eta _le^{i\theta _l}\eta _l`$ under the local U(1) gauge transformation.

In a deconfined state the RFO represents the ‘radial’ fashion in which QNF takes place leading to the possibility of a finite energy spinon or holon. The local U(1) gauge invariant physical spinon operator at a site $`l`$, $`s_{l\sigma }^{}\eta _l^s`$ and the holon operator $`e_l^{}\eta _l^h`$ are the analogue of the finite energy Laughlin quasi hole operator. They create a spin half moment or a missing / excess charge in an otherwise singlet environment, by pushing another spin half moment to the boundary.

In a gapless deconfined state the spin and charge disturbance created (by the RFO), in the process of QNF , will extend to infinity in a power law fashion. We will call this as a power law form factor of the ‘extended’ spinon or holon. The form factors will in general be different for a spinon and a holon. In the 1d Heisenberg and Hubbard models there is a domain wall of spin disturbance of infinite size, around the unpaired spin, as the energy of the spinon tends to zero.

Also the spin SU(2) symmetry (present for zero magnetic field) allows only half integer value of the spin after the quantum number fractionization. In terms of RFO the virtual excitations they create are always of singlet type. Once we have an external magnetic field the spin rotation symmetry is reduced to U(1) symmetry and the spin of the spinon can be different from half, as is also known to be the case for the 1d Heisenberg model. In the presence of the magnetic field the RFO will contain some triplet components, which may change the effective spin of the spinon from half.

The case of holon is more complex in the presence of gapless charge excitation, as in the spin gap metallic phase. The U(1) symmetry associated with the electric charge in general allows any fraction $`e^{}e`$ for the holon excitations. This is formally achieved by the appropriate choice of the ‘pure’ charge fluctuation terms $`e_{i_s}^{}e_{i_{s+1}}`$ in the U(1) string, which can push or attract a fractional average charge away from the site $`l`$. This is indeed what happens in the 1d Hubbard model. In general, in view of the charge gaplessness of the metallic phase of the 2d cuprates, the charge of the holon should continuously changes as a function of doping $`x`$ from $`e^{}=e`$ for the Mott insulator. I do not find a basic principle or a topological reason in the metallic phase of the cuprates which restricts the charge $`e^{}`$ to $`e`$. The power law accumulation of charges and the associated pushing of states to the boundary also leads to the natural possibility of Haldane’s fractional exclusion statistics.

All our discussion was for an infinite system. If we have a finite and large system, a spinon or a holon can not be created in isolation, they have to be created in pairs. Our formal construction can be carried out without any difficulty in this case as well. The above discussion does not prove the existence of deconfined phase in some 2d models, but gives a unique construction for objects with fractional quantum numbers in the s-particle formalism. It should be remarked that the RVB construction is dimension independent, suggesting the likely possibility of a fairly universal character of the quantum number fractionization in the t-J model.

Having discussed the essential aspects in the U(1) gauge mechanism of quantum number fractionization in some detail we go to the construction of SF. Central to their construction, for the t-J model in the s-particle representation is a neutral d -wave spinon pair operator $`d_i`$ and ‘half’ of the pair $`b_i`$ defined in terms of the nearest neighbor valence bond operator $`b_{ij}^{}=\frac{1}{\sqrt{2}}(s_i^{}s_j^{}s_i^{}s_j^{})`$,

$$d_i^{}\underset{j}{}\alpha _{ij}b_{ij}^{}\text{a}nd(b_i^{})^2d_i^{}$$
(11)
And $`\alpha _{ij}=1\text{a}nd1`$ for the horizontal and vertical nearest neighbor bonds respectively and zero otherwise. The half of spinon pair $`b_i`$ is used to define the the ‘spinon’ and ‘chargon’ operators $`(f_{i\sigma },h_i)`$ through the following relations:

$$e_ib_ih_i\text{a}nds_{i\sigma }b_if_{i\sigma }$$
(12)
If $`bi`$ transforms as $`b_ie^{i\theta _i}b_i`$, then $`(h_i,f_{i\sigma })`$ are gauge invariant physical observables. Let us check how $`b_i`$ transforms. On applying the U(1) gauge transformation (equation 2) the $`d_i`$ operator transforms as:

$$d_i\underset{j}{}\alpha _{ij}e^{i\theta _i}b_{ij}e^{i\theta _j}=e^{i\theta _i}\underset{j}{}\alpha _{ij}b_{ij}e^{i\theta _j}$$
(13)
The local phase factors of the U(1) transformation do not completely factor out. That is, the d-symmetric spinon pair, as opposed to a valence bond pair, does not carry a definite U(1) charge of the local symmetry. Consequently, the half of spinon pair $`b_i`$ also does not carry unit U(1) charge. Thus the spinon and chargon defined through equation (4) are not neutral under the local U(1) gauge transformation, making them unqualified as physical excitations. One may invoke a coarse gaining argument and suggest that if the local U(1) gauge transformation is smooth over the scale of the lattice parameter, $`d_i`$ transforms as $`d_ie^{2i\theta _i}`$. However, a coarse graining of a local symmetry is not allowed, as the double occupancy constraint on individual sites is ensured only by the local U(1) invariance.

One may be tempted to split the bond singlet operator as $`b_{ii+\mathrm{\Delta }}b_ib_{i+\mathrm{\Delta }}`$, to extract an RFO. This is dangerous and such a split representation for the bond operator fails to reproduce the important commutation relation:

$$[b_{ij}^{}b_{ij},b_{jk}^{}b_{jk}]=i𝐒_j(𝐒_i\times 𝐒_k)$$
(14)
Is there a way of rectifying the U(1) gauge invariance problem ? The answer is obvious once we recognize that the ‘half d-wave spinon pair’ operator $`b_i`$ is invented to do the same job that the RFO $`\eta _i^s`$ and $`\eta _i^h`$ did in the RVB construction. The only objects that carry a U(1) charge (apart from the s-particles themselves) are the open U(1) strings or suitable linear combinations. Open U(1) strings alone can neturalize the U(1) RVB charge of $`s_{i\sigma }^{}`$ or $`e_i^{}`$. The $`Z_2`$ string of SF, for the t-J model is a consequence of their erroneous construction, and it does not solve our U(1) problem.

It is also clear that once we extract $`\eta _i`$ in a gauge invariant fashion as detailed above, there is no new $`Z_2`$ symmetry - we have the full U(1) symmetry. However, for a reason very different from Senthil and Fisher’s, we will show below an approximate reduction of the original U(1) to a $`Z_2`$ symmetry arising from the s-particle construction. This is related to the ‘RVB flux’ or ‘RVB magnetic field’ that came up in a concrete fashion from the discovery of Affleck-marston’s $`\pi `$-flux phase of the Heisenberg antiferromagnet. Later works have also introduced a related $`Z_2`$ symmetry in frustrated spin systems. This $`Z_2`$ symmetry is also related to the topological order that Kivelson and collabarators were trying to bring out through their study of quantum dimer models and the phase string in Sutherland’s analysis of the short range RVB state.

To understand this in detail in our context, let us go back to the early work by Yu Lu, Tosatti and the present author where the local U(1) gauge theory was reduced to a local $`Z_2`$ gauge theory using an important observation of Marston. Consider a closed string made of $`\tau _{ij}=_\sigma `$$`s_{i\sigma }^{}`$$`s_{ȷ\sigma }`$, forming an elementary pleaquette:

$$W_{ijkl}\tau _{ij}\tau _{jk}\tau _{kl}\tau _{li}e^{i{\scriptscriptstyle 𝐀_{rvb}𝑑𝐥}}$$
(15)
Marston considered the Heisenberg model and found that this plaquette operator satisfies the following important identity in the physical Hilbert space:

$$W_{ijkl}^3=W_{ijkl}$$
(16)
As the plaquette loop $`W_{ijkl}`$ is the analogue of the Bohm-Aharanov phase factor of the RVB magnetic field, Marston concluded that the important RVB fluxes enclosed by the plaquette operator are $`Z_2`$ in character ( $`0`$ or $`\pi `$ flux, modulo $`2\pi `$) and he implemented it through a Chern-Simon term in the U(1) gauge theory. On the other hand the authors of reference used Marston identity in an approximate and physical way and reduced the gauge theory from U(1) to $`Z_2`$. Since Marston identity involves only the spinon s-particles, the reduction of U(1) to $`Z_2`$ is restricted to the spin sector. The charge sector continues to have the U(1) symmetry, but is influenced by the condensation of the $`\pi `$ flux in the spin sector as shown in ref.

We discuss an important limitation in Fisher et al.’s approach in the context of cuprates. They are strongly influenced by the notion of the charge incompressible quantum liquid state of an interacting charge 2e bose fluid in 2d, where the insulating phase can be thought of as a vortex pair condensed state. This could exhibit QNF with a fractional charge $`e^{}(=\frac{2e}{2}`$). This is a non trivial quantum fluid, where the powerful charge-vortex duality can be invoked for further understanding. In view of the gap, the solitonic objects have induced charges that are ‘localized’ - in other words they have a finite size. In this situation the exchange statistics makes sense. And one can use relations such as $`\mathrm{\Phi }_c\mathrm{\Phi }_v=1`$ (analogue of Dirac’s quantization condition derived using spatial exchanges) without any problem and this restricts allowed fractional charges to rational values (because the vorticity is always an integer).

Experiments in cuprates show a metallic spin gap state crossing over to a (spin)gapless metallic non fermi liquid state as we increase doping. Hence a physically meaningful zero temperature reference state is a metallic state rather than an insulating state. This reference non fermi liquid state undergoes a superconducting transition at a $`T_c`$. We have to discuss the quantum number fractionization issue in this reference metallic state rather than a charge insulator at finite doping.

The ‘extended’ character of holon with a power law form factor (discussed earlier) makes the exchange statistics lose its precise meaning and also the corresponding relation $`\mathrm{\Phi }_c\mathrm{\Phi }_v=1`$. This was shown by Haldane and Wu in the context of a suggestion of Chao et al. for fractional exchange statistics for quantized vortices in superfluid 4He film. Even though exchange statistics looses its meaning, a new possibility arises, namely Haldane’s exclusion statistics. This notion is meaningful in the gapless system and is dimension independent. In some approximate studies it has been proposed that the spinons and holons have non trivial exclusion statistics in the 2d t-J model.

In view of the above discussion we conclude that the QNF of $`e^{}=\frac{2e}{2}`$ is an artificial restriction in the spin gap phase of the cuprates. The $`Z_2`$ character of the spin sector leads to the interesting possibility of $`\pi `$ flux condensation (i.e., a cross over to the spin gap phase) from the uniform RVB state. It is very likely that the charge of the holon (and its exclusion statistics) will continuously vary from $`e^{}=e`$ as a function of doping. This is also a fundamental difference of cuprates from the well gapped quantum Hall system or other charge insulators.

In conclusion, if quantum number fractionization occurs in cuprates leading to the experimentally seen non-fermi liquid liquid phase and the d-wave superconducting state, the only viable mechanism seems to be through the RVB mechanism of spin charge separation, which comes out in a natural fashion in the RVB U(1) gauge theory. The variational/RVB mean field theory method combined with Gutzwiller projection, or the tomographic Luttinger liquid approach are different ways of solving the same complex problem.

I wish to thank R. Shankar (Madras) for discussions.


text:
# Quantum Monte Carlo Analysis of Exchange and Correlation in the Strongly Inhomogeneous Electron Gas

\[

## Abstract

> We use variational quantum Monte Carlo to calculate the density-functional exchange-correlation hole $`n_{xc}`$, the exchange-correlation energy density $`e_{xc}`$, and the total exchange-correlation energy $`E_{xc}`$, of several electron gas systems in which strong density inhomogeneities are induced by a cosine-wave potential. We compare our results with the local density approximation and the generalized gradient approximation. It is found that the nonlocal contributions to $`e_{xc}`$ contain an energetically significant component, the magnitude, shape, and sign of which are controlled by the Laplacian of the electron density.

\]

Kohn-Sham density-functional theory (DFT)  shows that it is possible to calculate the ground-state properties of interacting many-electron systems by solving only one-electron Schrödinger-like equations. The results are exact in principle, but in practice it is necessary to approximate the unknown exchange and correlation (XC) energy functional, $`E_{xc}[n(𝐫)]`$, which expresses the many-body effects in terms of the electron density $`n(𝐫)`$. The current popularity of density-functional methods in condensed matter physics, quantum chemistry, and materials science reflects the remarkable success of fairly simple approximate XC energy functionals.

In the local density approximation (LDA), the XC hole $`n_{xc}`$ about an electron at $`𝐫`$ is assumed to be the same as in a uniform electron gas of fixed density $`n=n(𝐫)`$. This works surprisingly well, even though real solids and molecules have strongly inhomogeneous electron densities. Almost all attempts to go beyond the LDA have been based on studies of XC in the weakly inhomogeneous electron gas. The widely used generalized gradient approximation (GGA)  incorporates information about first-order gradient corrections to the LDA, as obtained in the limit of slowly varying electron densities. This is done in a controlled way such that, for example, $`n_{xc}`$ satisfies known sum rules , but without any guidance from the behavior of $`n_{xc}`$ in the strongly inhomogeneous regime. This may explain why current GGAs, although better than the LDA, are not consistently able to deliver the very high accuracy of $``$0.1 eV required to study many chemical reactions. To provide a new direction in the search for more accurate and reliable approximate functionals, we have used variational quantum Monte Carlo (VMC)  to explore exchange and correlation in the strongly inhomogeneous electron gas.

In this letter we calculate the central quantities in DFT: the XC hole $`n_{xc}`$, the XC energy density $`e_{xc}`$, and the total XC energy $`E_{xc}`$. We study an electron gas of average density $`n^0`$ corresponding to $`r_s`$$`=`$2 and induce strong inhomogeneities by applying a cosine-wave potential $`V_q\mathrm{cos}(𝐪.𝐫)`$. Fixing $`V_q`$ at $`2.08ϵ_F^0`$, we consider several values of $`q2.17k_F^0`$, where $`ϵ_F^0`$ and $`k_F^0`$ are the Fermi energy and wave vector of a uniform electron gas with $`r_s`$$`=`$2.

The electron density in our systems varies strongly and rapidly on the scale of the inverse local Fermi wavevector $`k_F(𝐫)^1=(3\pi ^2n(𝐫))^{1/3}`$, producing a strikingly nonlocal behavior of $`n_{xc}`$ that cannot be described by semilocal corrections to the LDA. We show, however, that the resulting LDA errors in $`e_{xc}`$ have a dominant and energetically significant component, the magnitude, shape, and sign of which are controlled by $`^2n(𝐫)`$. The GGA is unable to correct the LDA errors in $`E_{xc}`$ resulting from this component in an adequate way, a deficiency that severely restricts its accuracy. The relevance of Laplacian terms has been pointed out previously , but our calculations provide the first quantitative evidence of their importance in strongly inhomogeneous systems and for $`e_{xc}`$.

The starting point of our calculations is the adiabatic connection formula . In Hartree atomic units ($`e=m=4\pi ϵ_0=\mathrm{}=1`$), this expresses $`E_{xc}[n]`$ as the volume integral of an XC energy density $`e_{xc}`$ defined by:

$$e_{xc}(𝐫,[n])=\frac{1}{2}𝑑𝐫^{}\frac{n(𝐫)n_{xc}(𝐫,𝐫^{})}{|𝐫𝐫^{}|},$$
(1)
where $`n_{xc}`$ is the coupling-constant-averaged XC hole. The definition of $`e_{xc}`$ used in the construction of the GGA differs from ours by an integration by parts that does not affect the integrated $`E_{xc}`$. We favor Eq. (1), however, because of its clear physical interpretation and because it aids comparison with the LDA XC energy density, which is constructed from an approximate hole.

The XC hole $`n_{xc}`$ is obtained via a coupling-constant integration,

$`n(𝐫)n(𝐫^{})+n(𝐫)n_{xc}(𝐫,𝐫^{})`$ (2)
$`=`$ $`{\displaystyle _0^1}𝑑\lambda \mathrm{\Psi }_\lambda |{\displaystyle \underset{i}{}}{\displaystyle \underset{j(i)}{}}\delta (𝐫𝐫_i)\delta (𝐫^{}𝐫_j)|\mathrm{\Psi }_\lambda ,`$ (3)
where $`\mathrm{\Psi }_\lambda `$ is the antisymmetric ground state of the Hamiltonian $`\widehat{H}^\lambda =\widehat{T}+\lambda \widehat{V}_{ee}+\widehat{V}^\lambda `$ associated with coupling constant $`\lambda `$. Here $`\widehat{T}`$ and $`\widehat{V}_{ee}`$ are the operators for the kinetic and electron-electron interaction energies, and $`\widehat{V}^\lambda =\mathrm{\Sigma }V^\lambda (𝐫_i)`$ is the one-electron potential needed to hold the electron density $`n^\lambda (𝐫)`$ associated with $`\mathrm{\Psi }_\lambda `$ equal to $`n(𝐫)`$ for all values of $`\lambda `$ between 0 and 1. At full coupling ($`\lambda `$$`=`$$`1`$), $`V^\lambda (𝐫)`$ coincides with the external potential of the system. For other $`\lambda `$ values $`V^\lambda `$ must be determined.

Our VMC method  for calculating $`n_{xc}`$ and $`e_{xc}`$ from Eqs. (1-2) is a generalization of the scheme used by Hood et al. . It amounts to treating both $`\mathrm{\Psi }^\lambda `$ and $`V^\lambda `$ variationally and determining the variational parameters by simultaneously minimizing the variance of the local energy  and the deviation of $`n^\lambda `$ from $`n`$ . This is achieved by adding a term proportional to $`|n^\lambda `$-$`n|`$ to the variance of the energy and minimizing the resulting penalty function, which attains its absolute minimum of zero only when $`\mathrm{\Psi }^\lambda `$ satisfies the Schrödinger equation associated with $`\lambda `$ *and* $`n(𝐫)=n^\lambda (𝐫)`$.

The VMC calculations are performed for a finite spin-unpolarized electron gas in a FCC simulation cell subject to periodic boundary conditions. We generate the ground-state density by applying the cosine-wave potential to the $`\lambda `$$`=`$$`0`$ system and solving the self-consistent Kohn-Sham equations within the LDA to obtain the single-particle orbitals $`\varphi _i`$. We then define the resulting density to be the ground-state density of the fully interacting system . We consider potentials with $`q=1.11k_F^0`$, $`1.55k_F^0`$, and $`2.17k_F^0`$, and cells containing $`64`$, $`78`$, and $`69`$ electrons, respectively. The large amplitude of the cosine-wave potential ensures that these systems are in the strongly inhomogeneous regime, $`|n(𝐫)n^0|/n^01`$.

The adiabatic calculations were performed using $`6`$ equidistant values of $`\lambda `$ in the range $`[0,1]`$, and with the following Slater-Jastrow ansatz for $`\mathrm{\Psi }^\lambda `$:

$$\mathrm{\Psi }^\lambda =D^{}D^{}\mathrm{exp}\left[\underset{i>j}{}u_{\sigma _i,\sigma _j}^\lambda (r_{ij})+\underset{i}{}\chi ^\lambda (𝐫_i)\right],$$
(4)
where $`r_{ij}=|𝐫_i𝐫_j|`$ and $`D^{}`$ and $`D^{}`$ are spin-up and spin-down Slater determinants constructed using the exact Kohn-Sham orbitals $`\varphi _i`$. The two-body term $`u_{\sigma _i,\sigma _j}^\lambda `$ contained both a fixed part and a variational part as discussed in . The one-body term $`\chi ^\lambda `$, the potential $`V^\lambda `$, and the electron density $`n^\lambda `$ were all expanded in plane waves. At each $`\lambda `$ we used a total of $`20`$ variational parameters in $`u^\lambda `$ and $`\chi ^\lambda `$ and up to $`7`$ coefficients in the plane-wave expansions of $`n^\lambda `$ and $`v^\lambda `$. The optimization of the parameters in $`\mathrm{\Psi }^\lambda `$ and $`V^\lambda `$ was performed using $`96000`$ statistically uncorrelated electron configurations. This was sufficient to reduce the root mean square deviation of $`n^\lambda (𝐫)`$ from $`n(𝐫)`$ to less than $`0.5\%`$ of $`n(𝐫)`$ for all values of $`\lambda `$ and all systems. After optimization, expectation values were calculated using the Monte Carlo Metropolis method  with $`10^6`$ independent configurations of all electrons. Throughout we used the modified electron-electron interaction described in , which virtually eliminates the finite-size errors arising from the long range of the Coulomb potential. The statistical errors were negligible except in the tails of $`n_{xc}`$ in low-density regions, where they were less than $`3\%`$; this is much smaller than the differences between the VMC and LDA XC holes.

The largest systematic errors are caused by the finite size of the system and the approximate nature of $`\mathrm{\Psi }^\lambda `$. These errors combine such that, even in a homogeneous electron gas, $`e_{xc}^{VMC}e_{xc}^{LDA}`$. To circumvent this problem we performed additional VMC calculations for finite homogeneous electron gases with $`N`$$`=`$$`64`$ and $`r_s=0.8`$, 1, 2, 3, 4, 5, 8, and 10. The results enabled us to construct a Perdew-Zunger (PZ) parameterization  of the VMC XC energy per electron of a finite uniform electron gas with $`N`$$`=`$$`64`$. This parameterization was used to calculate $`e_{xc}^{LDA}`$ and $`e_{xc}^{GGA}`$ in all systems studied, ensuring that these quantities were exactly equal to $`e_{xc}^{VMC}`$ in any finite homogeneous system with $`N`$$`=`$$`64`$. This procedure largely eliminates the systematic errors in the calculated differences between $`e_{xc}^{VMC}`$ and $`e_{xc}^{LDA}`$ reported below.

In Fig. 1 we show snapshots of the deformation of $`n_{xc}^{VMC}`$ around an electron moving in the $`q`$$`=`$$`1.55k_F`$ system along a line parallel to $`𝐪`$, the direction of maximum inhomogeneity, from a density maximum towards the tail of $`n(𝐫)`$. The XC hole is plotted as a function of $`𝐫^{}`$ around a fixed electron at $`𝐫`$, with $`𝐫^{}`$ ranging in a plane parallel to $`𝐪`$. Also shown is the corresponding LDA hole $`n_{xc}^{LDA}`$ . At the density maximum (not shown) both $`n_{xc}^{VMC}`$ and $`n_{xc}^{LDA}`$ are centered around the electron. However, unlike $`n_{xc}^{LDA}`$, which is always spherically symmetric, $`n_{xc}^{VMC}`$ is contracted in the direction of inhomogeneity and less compact. As the electron moves away from the density maximum to a point on the slope (top panel), the nonlocal nature of $`n_{xc}^{VMC}`$ becomes manifest. While $`n_{xc}^{LDA}`$ is still centered around the electron and is rather diffuse, $`n_{xc}^{VMC}`$ lags behind near the density maximum and is much more compact. The nonlocal behavior of $`n_{xc}^{VMC}`$ becomes remarkable at the density minimum. Here $`n_{xc}^{VMC}`$ has two large nonlocal minima, each centered around a density maximum $``$$`2.80`$ a.u. away from the electron. The LDA hole at this point is much more long ranged than $`n_{xc}^{VMC}`$ and is spread over the whole system in order to satisfy the sum rule : $`𝑑𝐫^{}n_{xc}^{LDA}(𝐫,𝐫^{})=1`$.

This striking nonlocality of $`n_{xc}^{VMC}`$ also occurs in the other two systems we considered. Clearly, semilocal corrections are unable to significantly improve the LDA description of the XC hole in our systems, and fully nonlocal approximations are required. We found, however, that despite the strong nonlocality of $`n_{xc}^{VMC}`$, the LDA errors in $`e_{xc}^{VMC}`$ *can* be described in terms of a semilocal quantity, the Laplacian of the electron density.

In Fig. 2 we show $`e_{xc}^{LDA}e_{xc}^{VMC}`$ for two of the strongly inhomogeneous systems studied, where $`e_{xc}^{LDA}`$ is calculated using the exact ground-state density $`n(𝐫)`$. The results are plotted along a line parallel to $`𝐪`$ (we call this direction $`y`$). Also shown are $`n(𝐫)`$ and $`^2n(𝐫)`$, plotted along the same line. It is apparent that the shape, magnitude, and sign of the LDA errors in $`e_{xc}`$ closely follow the shape, magnitude, and sign of $`^2n(𝐫)`$. The LDA errors in $`e_{xc}`$ are large and negative in regions where $`^2n(𝐫)`$ is large and negative (around density maxima), and large and positive in regions where $`^2n(𝐫)`$ is large and positive. The GGA XC energy density is not defined via Eq. (1) and so GGA errors are not shown.

The VMC values of the integrated $`E_{xc}`$ are shown in Table I, along with the differences $`\mathrm{\Delta }E_{xc}^{LDA}=E_{xc}^{LDA}E_{xc}^{VMC}`$ and $`\mathrm{\Delta }E_{xc}^{GGA}=E_{xc}^{GGA}E_{xc}^{VMC}`$ (The version of the GGA used here is due to Perdew, Burke, and Ernzerhof ). The LDA errors in $`E_{xc}`$ reflect the profound effect of the Laplacian errors in $`e_{xc}`$ and change sign from positive (for the $`q`$$`=`$$`1.11k_F^0`$ system ) to negative (for the two other systems) as $`q`$ increases and the negative contributions to $`\mathrm{\Delta }e_{xc}`$, which occur where $`^2n(𝐫)<0`$, become dominant. The GGA corrections are by construction always negative ; they improve $`E_{xc}^{LDA}`$ for the $`q`$$`=`$$`1.11k_F`$ system but worsen it for the two other systems.

Since any smooth density $`n(𝐫)`$ may be expanded as a power series, the XC energy density functional may always be written as

$$e_{xc}(𝐫,[n])=e_{xc}(𝐫,n(𝐫),_in(𝐫),_i_jn(𝐫),\mathrm{}).$$
(5)
The GGA may be viewed as an attempt to approximate the right-hand side of Eq. (5) as a nonlinear function of $`n(𝐫)`$ and $`_in(𝐫)`$. In the case of the exchange energy, a uniform scaling argument  shows that

$$e_x=F_x(p,l,\mathrm{})e_x^{LDA},$$
(6)
where $`F_x`$ is an enhancement factor, $`e_x^{LDA}`$ is the LDA exchange energy density, and $`p=|n|/(2k_F(𝐫)n(𝐫))`$ and $`l=^2n/(4k_F^2(𝐫)n(𝐫))`$ are a dimensionless gradient and a dimensionless Laplacian, respectively. We have calculated $`F_x`$ exactly for our systems and investigated its dependence on $`p`$ and $`l`$. Figs. 3a and 3b are scatter plots showing values of $`F_x`$ plotted against values of $`p`$ and $`l`$, respectively.

The two-valued nature of Fig. 3a arises because the mapping from density gradient to position is two valued in our systems. For example, the density gradient is zero both at the minimum and the maximum of the density, where the required corrections to the LDA XC hole are completely different, as can be inferred from Fig. 1. In sharp contrast, $`F_x`$ is a simple and almost unique function of $`l`$. This suggests that the inclusion of Laplacian terms may allow the construction of simpler and more accurate approximate functionals. In most gradient expansions the Laplacian terms allowed by symmetry are transformed into $`|n(𝐫)|^2`$ terms via an integration by parts. This is only possible when the dependence on $`l`$ is linear, however, which is not the case here as can be seen from Fig. 3b. Furthermore, the integration by parts destroys the physical interpretation in terms of the XC hole and so hinders further progress.

We note that there is a one-to-one relationship between position and $`l`$ in our structures, and so an enhancement factor of the form $`F_x(l)`$ might not be as universal as Fig. 3b suggests. However, the energetic significance of the Laplacian terms, and the strong similarity between the form of the Laplacian and the LDA errors in $`e_{xc}`$ (see Fig. 2), give us confidence in the physical importance of the Laplacian in describing inhomogeneity corrections to the LDA. Very recently, Becke  and Perdew et al.  have constructed beyond-GGA functionals that contain the reduced Laplacian as a parameter. We consider this a step in the right direction and believe that the detailed *local* information on $`n_{xc}`$ and $`e_{xc}`$ in the strongly inhomogeneous regime provided by our VMC calculations can be used to improve such functionals.

We thank R.Q. Hood for useful discussions, G. Rajagopal and A.J. Williamson for help with the VMC codes, and L. Smith for help with visualization. M.N. was supported by an EU Marie Curie Fellowship under Grant ERBFMBICT961736. Our calculations were performed on the CRAY-T3E at EPCC.


text:
# A study of center vortices in 𝑆𝑈(2) and 𝑆𝑈(3) gauge theories aafootnote apresented by M. Pepe.

## 1 Introduction

There are several approaches addressing the problem of understanding the mechanism of color confinement in non-Abelian gauge theories. The most popular share the idea that only a subset of the degrees of freedom are relevant for confinement. In the Abelian projection approach$`^{\mathrm{?},\mathrm{?}}`$, one takes into account the maximal Abelian subgroup $`U(1)^{N1}`$ of the gauge group $`SU(N)`$ and monopoles are the effective degrees of freedom under study. In the Center projection approach$`^{\mathrm{?},\mathrm{?}}`$, it is the center group $`Z_N`$ which is considered and center vortices are the effective degrees of freedom under study.

These two schemes are commonly believed to give alternative descriptions of confinement. Many analytical and numerical studies have been and are being performed using these two approaches. They give clear evidence that both Abelian and center degrees of freedom play a relevant role.

The reduction of the gauge symmetry and the selection of the effective degrees of freedom is usually carried out by a partial fixing of the initial $`SU(N)`$ gauge freedom. We show that, in the Laplacian Center gauge, monopoles and center vortices are closely related in a unified description and that the latter are the effective degrees of freedom relevant for confinement. Indeed center vortices are related to a more reduced gauge symmetry than monopoles. A fundamental issue of this study is the use of the Laplacian Center gauge, which is free from the problem of the lattice Gribov copies affecting the widely used Maximal Center gauge.

## 2 The center degrees of freedom

Let us consider an $`SU(2)`$ gauge field defined in a plane and expressed in terms of the radial and angular components $`(A_r,A_\phi )`$. Suppose that the radial component is vanishing, $`A_r=0`$, and the angular one is given by $`A_\phi =\frac{1}{2r}\sigma _3`$. This gauge field configuration describes a magnetic flux tube crossing the plane at the origin. If we consider a Wilson loop encircling the origin, we see that it has a non–trivial value $`\text{e}^{i\pi \sigma _3}=1`$ with respect to the center group $`Z_2`$ of $`SU(2)`$. Conversely a Wilson loop non encircling the origin has a trivial value $`+1`$ with respect to $`Z_2`$. This is an example of what a center vortex is.

Consider now an $`SU(2)`$ gauge theory on the lattice and decompose the link variable $`U_\mu (x)`$ as the product of two parts

$$U_\mu (x)=Z_\mu (x)U_\mu ^{}(x)$$
(1)
where $`Z_\mu (x)`$ lives in the center group $`Z_2`$ and $`U_\mu ^{}(x)`$ is the coset part belonging to $`SU(2)/Z_2`$. For instance, this splitting can be carried out defining $`U_\mu ^{}(x)`$ as having positive trace. If $`W(C)`$ is a Wilson loop along the closed path $`C`$, making use of the decomposition (1), we can write

$$W(C)=\sigma (C)W^{}(C)=\left[\underset{p\mathrm{\Sigma }}{}\sigma (p)\right]W^{}(C)$$
(2)
$`W^{}(C)`$ and $`\sigma (C)_{p\mathrm{\Sigma }}\sigma (p)`$ are respectively the Wilson loops evaluated with the coset links $`U_\mu ^{}(x)`$ and with the center links $`Z_\mu (x)`$. $`_{p\mathrm{\Sigma }}`$ is the product over all the plaquettes $`p`$ belonging to a surface $`\mathrm{\Sigma }`$ supported by $`C`$; the value of $`\sigma (C)`$ does not depend on the choice of $`\mathrm{\Sigma }`$ and, for simplicity, we can choose it as the planar surface bounded by $`C`$. If we fix a gauge where $`U_\mu ^{}`$ is smooth, then $`W(C)`$ has a non–trivial value with respect to $`Z_2`$ if $`\sigma (C)`$ does. Moreover, considering (2) for a single plaquette, $`\sigma (p)=1`$ is a “signal” for a center vortex. In 4–dimensional space–time, center vortices form closed surfaces in the dual lattice.

It is possible to give a handwaving and qualitative argument showing how the described center vortices can give rise to an area law behavior for the Wilson loops. According to the previous argument, we are interested in the plaquettes $`p`$ having $`\sigma (p)=1`$. Suppose that $`t`$ is the probability for a plaquette to have $`\sigma (p)=1`$, consequently $`(1t)`$ is the probability to have $`\sigma (p)=+1`$. Consider a time–slice of a lattice $`SU(2)`$ gauge configuration, then center vortices form 1–dimensional closed strings in the dual lattice. The following figure displays an example of such a time–slice

The drawn plaquettes have $`\sigma (p)=1`$ and the broken line is the closed string in the dual lattice representing the center vortex. Consider a Wilson loop $`W(C)`$ which, in the figure, is represented by the continuous line; the plaquette in bold line belongs to the planar surface $`\mathrm{\Sigma }`$ spanned by $`C`$. Let $`L`$ be the linear extension of $`C`$ and $`\xi `$ the average linear size of the center vortex strings. If $`\xi L`$, the value $`\sigma (p)`$ of each plaquette $`p`$ in $`\mathrm{\Sigma }`$ is independent of the others; then, with respect to the center degrees of freedom, one can write

$$<W(C)><\sigma (p)>^A=\text{e}^{A\mathrm{log}(12t)}$$
(3)
where $`A`$ is the area of $`\mathrm{\Sigma }`$. With this qualitative argument we have obtained on one hand that the center degrees of freedom can give rise to an area law behavior for the Wilson loop and, on the other hand, that the string tension is about twice the probability $`t`$ for a plaquette to have $`\sigma (p)=1`$. We have assumed that $`\xi L`$: in order for this inequality to be satisfied for arbitrarily large Wilson loops, $`\xi `$ must be divergent, that is center vortex strings must percolate through the lattice. Conversely, suppose that the center vortices do not percolate, so that $`\xi `$ has a finite value, then for $`L`$ sufficiently large, one must have $`\xi L`$. The next figure displays such a case

The values $`\sigma (p)`$ of the plaquettes in $`\mathrm{\Sigma }`$ well inside the boundary $`C`$ are strongly correlated by pairs; moreover these pairs do not give a net contribution to $`\sigma (C)`$. So only the plaquettes belonging to a ring of thickness $`\xi `$ around the boundary of $`\mathrm{\Sigma }`$ can give a non–trivial contribution. According to this reasoning and, with respect to the center degrees of freedom, we can write

$$<W(C)><\sigma (p)>^{\xi P}=\text{e}^{\xi P\mathrm{log}(12t)}$$
(4)
where $`P`$ is the perimeter of $`\mathrm{\Sigma }`$. Thus, we have obtained a perimeter law behavior for the Wilson loop. The conclusion from this handwaving argument is that the center degrees of freedom can give rise to an area law for the Wilson loop and that the area/perimeter behavior can be recast in terms of percolation/non-percolation of center vortices.

## 3 Laplacian gauge fixing

Many lattice simulations – starting with the initial results by Greensite and collaborators$`^\mathrm{?}`$ – have been performed in order to investigate the role of the center degrees of freedom in the non perturbative features of the non-Abelian gauge theories. In these studies the selection of the center degrees of freedom is carried out by the numerical partial gauge fixing of an ensemble of configurations. The most widely used gauge is the Maximal Center gauge. In this gauge, it has been shown$`^{\mathrm{?},\mathrm{?}}`$ that the removal of the center vortices leads to the loss of the area law for the Wilson loop, to the restoration of chiral symmetry and to the disappearance of non–trivial topological features. The numerical implementation of the Maximal Center gauge fixing is performed by a local iterative procedure and, as there are many local extrema, lattice Gribov copies are present. This is a serious problem$`^{\mathrm{?},\mathrm{?}}`$ that can lead to a complete loss of meaningful information in the $`Z_2`$ projected model. Thus, it is important to study the role of the center degrees of freedom in the confinement mechanism considering a smooth gauge not affected by this problem of the lattice Gribov copies. In the Laplacian gauge the reduction of the gauge degrees of freedom is carried out in an unambiguous way. This gauge was proposed by Vink and Wiese$`^\mathrm{?}`$: they suggested to use the eigenvectors of the covariant Laplacian operator to fix the gauge. Since we are interested in reducing unambiguously the symmetry of the gauge group from $`SU(N)`$ to its center $`Z_N`$, it is useful to consider the Laplacian operator in the adjoint representation. In fact, as the adjoint representation is invariant under gauge transformations in $`Z_N`$, the adjoint Laplacian procedure fixes unambiguously the gauge up to the center symmetry.

Consider the 4–dimensional lattice $`SU(N)`$ gauge theory. The adjoint covariant Laplacian operator $`\mathrm{\Delta }_{yx}^{ab}(\dot{U})`$ is given by

$$\mathrm{\Delta }_{yx}^{ab}(\dot{U})=\underset{\mu }{}\left(2\delta _{y,x}\delta ^{ab}\dot{U}_\mu ^{ab}(x\widehat{\mu })\delta _{y,x\widehat{\mu }}\dot{U}_\mu ^{ba}(x)\delta _{y,x+\widehat{\mu }}\right)$$
(5)
where $`a`$,$`b=1,\mathrm{},(N^21)`$ are color indices and $`x`$,$`y`$ are space–time lattice coordinates. The dotted $`\dot{U}_\mu (x)`$ are the link variables in the adjoint representation and are related to the links $`U_\mu (x)`$ in the fundamental by

$$\dot{U}_\mu ^{ab}(x)=\frac{1}{2}\text{Tr}\left(\lambda _aU_\mu (x)\lambda _bU_\mu ^{}(x)\right)$$
(6)
$`\lambda _i`$, $`i=1,\mathrm{},(N^21)`$ being the generators of $`SU(N)`$ with the normalization $`\text{Tr}(\lambda _a\lambda _b)=2\delta _{ab}`$. If $`V`$ is the volume of the lattice, $`\mathrm{\Delta }(\dot{U})`$ is a $`[(N^21)V]\times [(N^21)V]`$ real symmetric matrix which depends on the gauge field. The eigenvalues $`\mu _j`$ of $`\mathrm{\Delta }`$ are real and the eigenvector equation is

$$\mathrm{\Delta }_{yx}^{ab}(\dot{U})\varphi _b^{(j)}(x)=\mu _j\varphi _a^{(j)}(y)$$
(7)
where $`\varphi ^{(j)}`$, $`j=1,\mathrm{},[(N^21)V]`$ are the (real) eigenvectors. So we can associate $`(N^21)`$–dimensional real vectors $`\varphi ^{(j)}(x)`$ to every lattice site $`x`$. Consider now a gauge transformation on the links in the fundamental representation $`U_\mu ^{}(x)=\mathrm{\Omega }(x)U_\mu (x)\mathrm{\Omega }^{}(x+\widehat{\mu })`$; the eigenvector equation (7) becomes

$$\dot{\mathrm{\Omega }}^{}{}_{}{}^{ai}(y)\mathrm{\Delta }_{yx}^{ik}(\dot{U}^{})\dot{\mathrm{\Omega }}^{kb}(x)\varphi _b^{(j)}(x)=\mu _j\varphi _a^{(j)}(y)$$
(8)
where $`\dot{\mathrm{\Omega }}`$ is the gauge transformation in the adjoint representation. This relation shows that the eigenvalues are gauge invariant and the eigenvectors transform according to $`\dot{\mathrm{\Omega }}^{ab}(x)\varphi _b^{(j)}(x)=\varphi _a^{(j)}(x)^{}`$. This transformation law can be rewritten as follows

$$\mathrm{\Omega }(x)\mathrm{\Phi }^{(j)}(x)\mathrm{\Omega }^{}(x)=\mathrm{\Phi }^{(j)}(x)^{}$$
(9)
where we have defined the $`su(N)`$ matrices – i.e. in the $`SU(N)`$ algebra – $`\mathrm{\Phi }^{(j)}(x)=_{a=1}^{N^21}\varphi _a^{(j)}(x)\lambda _a`$ and $`\mathrm{\Phi }^{(j)}(x)^{}=_{a=1}^{N^21}\varphi _a^{(j)}(x)^{}\lambda _a`$. Gauge transformations rotate the vectors $`\varphi ^{(j)}(x)`$ in color space and we can fix the gauge by requiring a conventional arbitrary orientation for the $`\varphi ^{(j)}(x)`$. To perform the reduction of the gauge symmetry from $`SU(N)`$ to $`Z_N`$ we only need to fix the orientation of two eigenvectors of the Laplacian operator. As we are interested in the non–perturbative features and in fixing a smooth gauge, we consider the two lowest-lying eigenmodes $`\varphi ^{(1)}`$ and $`\varphi ^{(2)}`$.

The gauge fixing procedure can be split into two steps. In the first, one rotates $`\mathrm{\Phi }^{(1)}(x)`$ at every $`x`$ so that $`\mathrm{\Phi }^{(1)}(x)^{}`$ is diagonal. This leaves a residual symmetry corresponding to gauge transformations in the Cartan subgroup $`U(1)^{N1}`$. To further reduce the gauge freedom we have to consider a second step where the second eigenvector $`\varphi ^{(2)}`$ is taken into account. The gauge transformation that has rotated $`\mathrm{\Phi }^{(1)}(x)`$ to the Cartan subalgebra, maps $`\mathrm{\Phi }^{(2)}(x)`$ to $`\mathrm{\Phi }^{(2)}(x)^{}`$. The remaining $`U(1)^{N1}`$ symmetry can be fixed to $`Z_N`$ by requiring that some conventionally chosen color components of the twice rotated matrix $`\mathrm{\Phi }^{(2)}(x)\mathrm{"}`$ vanish. Now we describe explicitly how to perform the presented two–step program for the $`SU(2)`$ gauge theory$`^\mathrm{?}`$, where, at each $`x`$, we consider the 3–dimensional real vectors $`\varphi ^{(1)}(x)`$ and $`\varphi ^{(2)}(x)`$.

Step 1. We consider equation (9) for the first eigenvector and we conventionally define $`\mathrm{\Omega }(x)`$ to be the gauge transformation that rotates $`\varphi ^{(1)}(x)`$ along the direction 3 in color space: $`\dot{\mathrm{\Omega }}(x)\varphi ^{(1)}(x)=\varphi ^{(1)}(x)^{}(0,0,1)`$. $`\mathrm{\Omega }(x)`$ is unambiguously defined up to gauge transformations $`V(x)`$ in $`U(1)`$. So we have reduced the gauge symmetry to the Cartan subgroup of $`SU(2)`$. This is the Laplacian Abelian gauge$`^{\mathrm{?},\mathrm{?}}`$.

Step 2. We apply the gauge transformation $`\mathrm{\Omega }(x)`$ found in Step 1 to $`\varphi ^{(2)}(x)`$: $`\varphi ^{(2)}(x)^{}=\dot{\mathrm{\Omega }}(x)\varphi ^{(2)}(x)`$. $`\varphi ^{(2)}(x)^{}`$ is not invariant under the rotations $`V(x)U(1)`$ that leave $`\varphi ^{(1)}(x)^{}`$ unchanged; then we can fix this symmetry by requiring, for example, that $`\varphi ^{(2)}(x)\mathrm{"}=\dot{V}(x)\varphi ^{(2)}(x)^{}`$ lie in the 1–3 color plane in the positive direction.

We have completely fixed the gauge symmetry in the adjoint representation and, as it is center–blind, we are left with the center symmetry $`Z_2`$. The described two–step procedure for $`SU(2)`$ can be extended$`^\mathrm{?}`$ to $`SU(N)`$ with the same pattern of gauge symmetry fixing: $`SU(N)U(1)^{N1}`$ in the first step and $`U(1)^{N1}Z_N`$ in the second one.

Two last, important remarks concern the accidental degeneracy of $`\mu _1`$ and $`\mu _2`$ in (7) and the arbitrariness in the eigenvectors $`\varphi ^{(j)}`$. If it happens that $`\mu _1`$ or $`\mu _2`$ is degenerate, the gauge fixing can not be carried out unambiguously and one has a (global) Gribov ambiguity. This case is really exceptional and never occurs in the numerical simulations. The second point is about the scale and sign arbitrariness of the eigenvectors $`\varphi ^{(j)}`$. Rescaling can not give rise to any ambiguity in the procedure while the freedom in the choice of the sign does. This global freedom can be eliminated with a conventional choice on $`\varphi ^{(j)}`$.

## 4 Local gauge ambiguities

The adjoint Laplacian gauge fixing procedure has local defects. Now we discuss how these defects can show up and how they can be identified with monopoles and center vortices in $`SU(2)`$. This discussion can be generalized to $`SU(N)`$$`^\mathrm{?}`$.

Step 2 ill–defined: the second step of the Laplacian gauge fixing is not defined at the points $`x`$ where $`\varphi ^{(1)}(x)//\varphi ^{(2)}(x)`$. In such a case also $`\varphi ^{(2)}(x)^{}`$ is invariant under rotations $`V(x)U(1)`$. The condition $`\varphi ^{(1)}(x)//\varphi ^{(2)}(x)`$ sets two constraints and so these points $`x`$ – where the gauge symmetry is promoted from $`Z_2`$ to $`U(1)`$ – form 2–dimensional surfaces in 4–dimensional space–time.

Step 1 ill–defined: the gauge fixing procedure can not be even started at the points $`x`$ where $`\varphi ^{(1)}(x)=(0,0,0)`$. These defects constitute 1–dimensional strings in the 4–dimensional space–time since 3 constraints must be satisfied. At these points the symmetry is not fixed and the gauge freedom is $`SU(2)`$.

The 2–dimensional surfaces of Step 2 ill–defined can be identified with center vortices. Suppose that at a point $`x_0`$ it happens that $`\varphi ^{(1)}(x_0)//\varphi ^{(2)}(x_0)`$, then moving along a small loop around the singularity point $`x_0`$, $`\varphi ^{(2)}`$ describes a $`2\pi `$ rotation in color space. As a $`2\pi `$ phase in the adjoint representation of $`SU(2)`$ corresponds to a $`\pi `$ phase in the fundamental one, it follows that a Wilson loop encircling $`x_0`$ has a non trivial value with respect to the center $`Z_2`$.

The 1–dimensional strings of Step 1 ill–defined can be identified with monopole world–lines. In fact, in analogy with the Georgi–Glashow model, the points $`x`$ where the Higgs field vanishes and the gauge symmetry can not be reduced from $`SU(2)`$ to $`U(1)`$, correspond to the monopole world–lines.

In the adjoint Laplacian gauge, monopoles and center vortices turn out to be closely related in a unified description. Consider the 2–dimensional surface $`𝒮`$ of the center vortices. At every $`x𝒮`$, $`\varphi ^{(1)}(x)`$ and $`\varphi ^{(2)}(x)`$ are parallel or anti–parallel. So $`𝒮`$ is divided in parts where $`\varphi ^{(1)}(x)`$ is parallel to $`\varphi ^{(2)}(x)`$ and parts where it is anti–parallel. By continuity, these patches must be separated by 1–dimensional strings where $`\varphi ^{(1)}(x)=0`$ or $`\varphi ^{(2)}(x)=0`$. Moreover, in the neighbourhood of a monopole, $`\varphi ^{(1)}(x)`$ has a hedgehog–like shape and there will be a direction where $`\varphi ^{(1)}//\varphi ^{(2)}`$. Thus monopole world–lines are embedded within the 2–dimensional surfaces of center vortices.

## 5 Numerical results and their interpretation

We have performed numerical simulations to investigate the role of the center degrees of freedom in the $`SU(2)`$ and $`SU(3)`$ lattice gauge theories. For $`SU(2)`$ we have collected 1000 configurations at $`\beta =2.3`$, $`2.4`$, $`2.5`$ on a $`16^4`$ lattice; for $`SU(3)`$ we have generated 500 configurations on a $`16^4`$ lattice at $`\beta =6.0`$. The following figure shows the measurement of the Creutz ratios $`\chi (R)=\mathrm{ln}(W(R,R)W(R1,R1)/W(R,R1)^2)`$ for $`SU(2)`$ at $`\beta =2.4`$ ($`W(R,T)`$ is the $`R\times T`$ Wilson loop).

Crosses refer to $`SU(2)`$, circles to center projection after Laplacian gauge fixing and stars to the coset part. The continuous band is the value in the literature$`^{\mathrm{?},\mathrm{?}}`$ for the $`SU(2)`$ string tension at the considered set of parameters. The results show, on one hand, the flattening of the Creutz ratios in the $`Z_2`$ sector and, on the other hand, the vanishing of the Creutz ratios computed with the coset links. We have obtained a similar behaviour for $`\beta =2.3`$ and $`2.5`$. In the case of $`SU(3)`$ also, the following figure shows the Creutz ratios in the $`Z_3`$ sector after Laplacian gauge fixing.

The continuous band is the value in literature$`^\mathrm{?}`$ for the $`SU(3)`$ string tension at the chosen set of parameters. Also in this case, one can clearly see flattening to a non vanishing value for the Creutz ratios evaluated with center projected links.

The good agreement with the values in the literature for the string tension in $`SU(2)`$ and $`SU(3)`$ should not be over–estimated. Numerical simulations are performed at finite lattice spacing and lattice artifacts can give non-negligible effects in the center projected theory. Our conjecture is that even if, at finite lattice spacing, the flattening value of the Creutz ratios in the center sector changes with the particular lattice Laplacian used to fix the gauge, this dependence vanishes in the continuum limit. The observation of such a behaviour would be a robust confirmation of the relevance of the center degrees of freedom in the confinement mechanism. To investigate this issue, we have considered three different lattice Laplacians differing by irrelevant operators: they have been built using smeared links in (6). Every set of 1000 configurations at $`\beta =2.3`$, $`2.4`$, $`2.5`$, has been fixed in each one of the three gauges and the Creutz ratios have been measured after center projection. The table summarizes our results:

$`R_i=\sqrt{\sigma _i/\sigma _{SU(2)}}`$ where $`\sigma _i`$ is the string tension measured in the $`Z_2`$ sector; $`i=0,1,2`$ is an index for the three lattice Laplacians and $`\sigma _{SU(2)}`$ is the value in literature for the $`SU(2)`$ string tension at the three values of $`\beta `$. Thus, it is in the continuum limit ($`\beta \mathrm{}`$) that the string tension measured after center projection correctly reproduces the value of the full gauge theory.

## Acknowledgments

We thank C. Alexandrou, M. D’Elia, S. Dürr and J. Fröhlich for useful discussions and G. Bali for providing us with a code for data analysis.

## References