text:
# Genuine Correlations of Like-Sign Particles in Hadronic Z0 Decays

## 1 Introduction

Correlations in momentum space between hadrons produced in high energy interactions have been extensively studied over many decades in different contexts . Being a measure of event-to-event fluctuations of the number of hadrons in a phase space domain of size $`\mathrm{\Delta }`$, correlations provide detailed information on the hadronisation dynamics, complementary to that derived from inclusive single-particle distributions and global event-shape characteristics.

The suggestion in  that multiparticle dynamics might possess (multi-)fractal properties or be “intermittent”, emphasized the importance of studying correlations as a function of the size of domains in momentum space. A key ingredient for such studies is the normalized factorial moment and factorial cumulant technique (see Sect. 2), which allows statistically meaningful results to be obtained even for very small phase space cells.

Unlike factorial moments, cumulants of order $`q`$ are a direct measure of the stochastic interdependence among groups of exactly $`q`$ particles emitted in the same phase space cell . Therefore, they are well suited for the study of true or “genuine” correlations between hadrons.

Whereas earlier work dealt mainly with correlations between pairs i.e. second-order correlations, the use of factorial cumulants in high-statistics experiments has established the presence of genuine correlations among groups comprising three or more hadrons, hereafter referred to as higher-order correlations.

Experimental results on hadron correlations are reviewed in . Except for heavy ion collisions, where correlations beyond second order are found to be small, as can be understood from a superposition of independent particle-production sources , in all other types of reactions, significant positive higher-order correlations are seen. In two- or three-dimensional (typically rapidity, azimuthal angle, transverse momentum) phase space cells, they rapidly increase as the cell-size, $`\mathrm{\Delta }`$, becomes smaller. Further studies of the dependence on particle charge confirm that, as conjectured in , correlations between hadrons with the same charge play an increasingly important role as $`\mathrm{\Delta }`$ decreases, thus pointing to the influence of Bose-Einstein (BE) interference effects ; for a recent review see . In contrast, correlations in multiplets composed of particles with different charges, which are more sensitive to multiparticle resonance decays than like-sign ones, tend to saturate in small phase space domains .

Two-particle Bose-Einstein correlations (BEC) have been observed in a wide range of multihadronic processes . Such correlations were extensively studied at LEP . Evidence for BEC among groups of more than two identical particles has also been reported . The subject has acquired particular importance in connection with high-precision measurements of the $`W`$-boson mass at LEP-II . For these, better knowledge of correlations in general is needed, as well as realistic Monte Carlo modelling of BEC .

The OPAL collaboration recently reported an analysis of the domain-size dependence of factorial cumulants in hadronic Z0 decays, using much larger statistics than in any previous experiment . In that study, no distinction was made between multiplets composed of like-charge particles and those of mixed charge. Clear evidence was seen for large positive genuine correlations up to fifth order. Hard jet production was found to contribute significantly to the observed particle fluctuation patterns. However, Monte Carlo models based on parton showers and string or cluster fragmentation, gave only a qualitative description of the $`\mathrm{\Delta }`$-dependence of the cumulants. Quantitatively, the model studied, which did not explicitly include BE-type correlation effects, underestimated significantly correlations between hadrons produced very close together in momentum space.

In the present paper, the high statistics OPAL data collected at and near the Z0 centre-of-mass energy are used to measure cumulants for multiplets of particles with the same charge, hereafter referred to as “like-sign cumulants”. They are compared to “all-charge” cumulants, corresponding to multiplets comprising particles of any (positive or negative) charge. Using the factorial cumulant technique, results of much higher precision than previously available are obtained on like-sign correlations up to fourth order. The role of Bose-Einstein-type effects is studied, using recently proposed BEC algorithms111A variety of methods which try to simulate BE-effects in Monte Carlo generators have been developed. A detailed exposition of the basic physics issues involved and a review of the presently most popular models and algorithms can be found in . in the Monte Carlo event generator Pythia for $`e^+e^{}`$ annihilation . Proceeding beyond the usual analyses of two-particle correlations, we show that, at least within the framework of this model, a good description can be achieved of the factorial cumulants up to fourth order in one-, two- and three-dimensional phase space domains.

The paper is organized as follows. Section 2 explains the normalized factorial cumulant method which allows genuine particle correlations to be measured. The OPAL detector, the event selection and the track selection criteria are detailed in Sect. 3. The data on factorial cumulants are presented in Sect. 4. They are compared with Pythia predictions and, in particular, with several variants of a model to simulate the Bose-Einstein effect. The results are summarized in Sect. 5. A brief overview of some of the BEC algorithms implemented in Pythia is given in the appendix.

## 2 Factorial cumulant method

To measure genuine multiparticle correlations in multi-dimensional phase space cells, we use the technique of normalized factorial cumulant moments, $`K_q`$, or “cumulants” for brevity, as proposed in .

In the present paper, the cumulants are computed as in a previous OPAL analysis . A $`D`$-dimensional region of phase space (defined in Sect. 3) is partitioned into $`M^D`$ cells of equal size $`\mathrm{\Delta }`$. From the number of particles counted in each cell, $`n_m`$ ($`m=1,\mathrm{},M^D`$), event-averaged unnormalized factorial moments, $`n_m^{[q]}`$, and unnormalized cumulants, $`k_q^{(m)}`$, are derived, using the relations given e.g. in . For $`q=2,3,4`$, one has

$`k_2^{(m)}`$ $`=`$ $`n_m^{[2]}n_m^2,`$ (1)
$`k_3^{(m)}`$ $`=`$ $`n_m^{[3]}3n_m^{[2]}n_m+2n_m^3`$ (2)
$`k_4^{(m)}`$ $`=`$ $`n_m^{[4]}4n_m^{[3]}n_m3n_m^{[2]}^2+12n_m^{[2]}n_m^26n_m^4.`$ (3)
Here, $`n^{[q]}=n(n1)\mathrm{}(nq+1)`$ and the brackets $``$ indicate that the average over all events is taken.

Normalized cumulants are calculated using the expression

$$K_q=(𝒩)^q\overline{k}_q^{(m)}/\overline{N_m^{[q]}}.$$
(4)
As proposed in , this form is used to correct for statistical bias and non-uniformity of the single-particle spectra. Here, $`N_m`$ is the number of particles in the $`m`$th cell summed over all $`𝒩`$ events in the sample, $`N_m=_{j=1}^𝒩(n_m)_j`$. The horizontal bar indicates averaging over the $`M^D`$ cells in each event, $`(1/M^D)_{m=1}^{M^D}`$.

The factorial moment of order $`q`$ of the multiplicity distribution of particles of the same species (e.g. all charged, negatives only, $`\mathrm{}`$) in a phase space domain $`\mathrm{\Delta }`$ is equal to the integral of the $`q`$-particle inclusive density, $`\rho _q`$, over that domain . As in the cluster expansion in statistical mechanics, $`\rho _q`$ can be decomposed into the sum of contributions from “accidental” coincidences of particles in $`\mathrm{\Delta }`$ and the true or “genuine” correlations. The latter are denoted here by the unnormalised cumulants, $`k_q`$, being the bin-averaged factorial cumulant functions, or “correlation functions”, for short .

Whereas $`n^{[q]}`$ depends on all correlation functions of order $`1pq`$, $`k_q`$ is a direct measure of stochastic dependence in multiplets of exactly $`q`$ particles. By construction, $`k_q`$ vanishes whenever a particle within the $`q`$-tuple is statistically independent of one of the others. For Poissonian multiplicity fluctuations, the cumulants of all orders $`q>1`$ vanish identically. Non-zero cumulants therefore signal the presence of correlations.

In the following, data are presented for “all-charge” and for “like-sign” multiplets. For the former, the cell-counts $`n_m`$ are determined using all charged particles in an event, irrespective of their charge. For the latter, the number of positive particles and the number of negative particles in a cell are counted separately. The corresponding cumulants are then averaged to obtain those for like-sign multiplets. It is to be noted that the cell-counting technique does not allow correlations to be measured directly among groups composed of particles of different charge. For these, other methods, such as correlation integrals (see e.g. Sect. 4.8 in ) have to be used. Nevertheless, cumulants for unlike-sign pairs can be indirectly derived using the equation

$$K_2^{\text{all}}=\frac{1}{2}K_2^{\text{ls}}+\frac{1}{2}K_2^{\text{us}},$$
(5)
which relates second-order cumulants for all-charge pairs to those of like-sign (ls) and unlike-sign (us) pairs.

## 3 Experimental details

The present analysis uses a sample of approximately $`4.1`$$`\times `$$`10^6`$ hadronic Z0 decays collected from 1991 through 1995. About 91% of this sample was taken at the Z0; the remaining part has a centre-of-mass energy, $`\sqrt{s}`$, within $`\pm 3`$ GeV of the Z0 peak.

The OPAL detector has been described in detail in . The results presented here are mainly based on the information from the central tracking chambers, which consist of a silicon microvertex detector, a vertex chamber, a jet chamber with 24 sectors each containing 159 axial anode wires, and outer $`z`$-chambers to improve the $`z`$ coordinate resolution222 The right-handed OPAL coordinate system is defined with the $`z`$ axis pointing in the direction of the e- beam and the $`x`$ axis pointing towards the centre of the LEP ring. $`r`$ is the coordinate normal to the beam axis, $`\phi `$ the azimuthal angle with respect to the $`x`$ axis, $`\theta `$ the polar angle with respect to the $`z`$-axis.. These detectors are located in a 0.435 T axial magnetic field and measure $`p_t`$, the track momentum transverse to the beam axis, with a precision of $`(\sigma _{p_t}/p_t)=\sqrt{(0.02)^2+(0.0015p_t)^2}`$ ($`p_t`$ in GeV/$`c`$) for $`|\mathrm{cos}\theta |<0.73`$. The jet chamber also measures the specific energy loss of charged particles, $`\mathrm{dE}/\mathrm{dx}`$, with a resolution $`\sigma (\mathrm{dE}/\mathrm{dx})/(\mathrm{dE}/\mathrm{dx})\mathrm{\hspace{0.17em}3.5}\%`$ for a track having 159 hits in the jet chamber. The energy loss is used to identify charged particles .

A sample of over $`2`$ million events was generated with Jetset 7.4/Pythia6.1 , including a full simulation of the detector . The model parameters were previously tuned to OPAL data  but Bose-Einstein effects were not explicitly incorporated. These events were used to determine the efficiencies of track and event selection and for correction purposes. In addition, for the evaluation of systematic errors, over $`1.1`$ million events were simulated with Pythia including BEC with the algorithm333We used algorithm BE32 (see  and the appendix) in subroutine PYBOEI with parameters MSTJ(51)=2, MSTJ(52)=9, PARJ(92)=1.0, PARJ(93)=0.5 GeV. PARJ(92) is equal to $`\lambda `$ in Eqs. (8) and (9) of the appendix; the parameter $`R`$ is given by $`\mathrm{}c/\text{PARJ(93)}`$. BE32.

The event selection criteria are based on the multihadronic event selection algorithms described in . It was required that tracks have at least 20 hits in the jet chamber, a first measured point at a maximum radial distance from the interaction point of 70 cm, a minimum transverse momentum with respect to the beam direction, $`p_t`$, of 0.15 GeV/$`c`$, a measured momentum $`p`$ smaller than 10 GeV/$`c`$, a measured polar angle satisfying $`|\mathrm{cos}\theta |<0.93`$, and a measured distance of closest approach to the origin of less than 5 cm in the $`r\phi `$ plane, and less than 40 cm in the $`z`$ direction.

The mean energy loss, $`\mathrm{dE}/\mathrm{dx}`$, of a track was required to be less than 9 keV/cm, thereby rejecting electrons and positrons. About 97% of photon conversion pairs were rejected, reducing the conversion background in the sample to less than 0.1% of the number of tracks . The fraction of pions in the track sample was estimated from Monte Carlo simulation to be about 80%.

Selected multihadron events were required to have at least 5 good tracks, a momentum imbalance (the magnitude of the vector sum of the momenta of all charged particles) of less than $`0.4\sqrt{s}`$ and the sum of the energies of all tracks (assumed to be pions) greater than 0.2 $`\sqrt{s}`$. These requirements provide rejection of background from non-hadronic Z0 decays, two-photon events, beam-wall and beam-gas interactions. In addition, the polar angle of the event sphericity axis, calculated using tracks that passed the above cuts, as well as unassociated electromagnetic and hadronic calorimeter clusters, had to satisfy $`|\mathrm{cos}\theta _{\text{sph}}|<0.7`$ in order to accept only events well contained in the detector. A total of about $`2.3`$$`\times `$$`10^6`$ events were finally selected for further analysis.

The cumulant analysis is performed in the kinematic variables rapidity, $`y`$, azimuthal angle, $`\mathrm{\Phi }`$, and the transverse momentum variable, $`\mathrm{ln}p_T`$, all calculated with respect to the sphericity axis.

* Rapidity is defined as $`y=0.5\mathrm{ln}[(E+p_{})/(Ep_{})]`$, with $`E`$ and $`p_{}`$ the energy (assuming the pion mass) and longitudinal momentum of the particle, respectively. Only particles within the central rapidity region $`2.0y2.0`$ were retained.
* In transverse momentum subspace, the logarithm of $`p_T`$, instead of $`p_T`$ itself, is used to eliminate as much as possible the strong dependence of the cumulants on cell-size arising from the nearly exponential shape of the $`p_T^2`$-distribution. Only particles within the range $`2.4\mathrm{ln}(p_T)0.7`$ ($`p_T`$ in GeV/c) were used.
* The azimuthal angle, $`\mathrm{\Phi }`$, is calculated with respect to the eigenvector of the momentum tensor having the smallest eigenvalue in the plane perpendicular to the sphericity axis. The angle $`\mathrm{\Phi }`$ spans the interval $`0\mathrm{\Phi }<2\pi `$.

The phase space is partitioned into $`M`$ bins of equal size for each of the three variables. $`M`$ varies between $`1`$ (full interval) and 400, corresponding to a smallest bin size in one dimension of $`\delta y_{\mathrm{min}}=0.01`$, $`\delta \mathrm{\Phi }_{\mathrm{min}}=0.9^{}`$ and $`\delta (\mathrm{ln}p_T)_{\mathrm{min}}=0.008`$, each significantly larger than the experimental resolution of the OPAL detector.

The cumulants for like-sign and all-charge multiplets, measured as a function of $`M`$, have been corrected for geometrical acceptance, kinematic cuts, initial-state radiation, resolution, secondary interactions and decays in the detector, using correction factors, $`U_q(M)`$, calculated for all-charge particle combinations and evaluated as in  using the Jetset /Pythia Monte Carlo without BEC.

As systematic uncertainties, we include the following contributions:

* The statistical error of the correction factors $`U_q(M)`$. Statistical errors due to the finite statistics of the Monte Carlo samples are comparable to those of the data.
* Track and event selection criteria variations as in . To this end, the cumulants have been computed changing in turn the following selection criteria: the first measured point was required to be closer than 40 cm to the beam, the requirement of the transverse momentum with respect to the beam axis was removed, the momentum was required to be less than 40 GeV/$`c`$, the track polar angle acceptance was changed to $`|\mathrm{cos}\theta |<0.7`$, and the requirement on the mean energy loss was removed. These changes modify the results by no more than a few percent in the smallest cells, and do not affect the conclusions.
* The difference between cumulants corrected with the factors $`U_q(M)`$ derived from Monte-Carlo calculations with and without Bose-Einstein simulation. The correction factors in these two cases differ by at most 5% in the smallest bins.
* The difference between cumulants corrected with $`U_q(M)`$-factors, calculated for all-charge particle combinations and those calculated for like-sign combinations. The correction factors coincide within 1%.

The total errors have been calculated by adding the systematic and statistical uncertainties in quadrature. It was further verified that our conclusions remain unchanged when events taken at energies off the Z0 peak are excluded from the analysis.

## 4 Results

### 4.1 Like-sign and all-charge cumulants

The fully corrected normalized cumulants $`K_q`$ ($`q=2,3,4`$) for all-charge444The data points for all-charge particle combinations are the same as those published in  but include different systematic uncertainties. and like-sign particle multiplets, calculated in one-dimensional ($`y`$ and $`\mathrm{\Phi }`$) (1D), two-dimensional $`y\times \mathrm{\Phi }`$ (2D) and three-dimensional $`y\times \mathrm{\Phi }\times \mathrm{ln}p_T`$ (3D) phase space cells, are displayed in Fig. 1 and Fig. 2.

From Fig. 1 it is seen that, even in 1D, positive genuine correlations among groups of two, three and four particles are present: $`K_q>0`$. They are substantially stronger in rapidity than in azimuthal angle. The cumulants increase rapidly with increasing $`M`$, the inverse of the bin size, for relatively large domains but saturate rather quickly. For $`K_2`$ this behaviour follows from the shape of the second-order correlation function which is known to be approximately Gaussian  in the two-particle rapidity difference $`\mathrm{\Delta }\delta y`$. The rapid rise and subsequent saturation can be understood from hard gluon jet emission. For small $`M`$, the cumulants are sensitive to the large-scale structure of the event, where any jet non-collinear with the event axis produces a strong fluctuation in the particle density. With $`M`$ increasing, the structure inside a single jet is progressively probed.

In contrast to 1D cumulants, those in 2D and 3D (Fig. 2) continue to increase towards small phase space cells. Moreover, the 2D and 3D cumulants are of similar magnitude at fixed $`M`$, indicating that the contribution from correlations in transverse momentum is small. This can be understood from the importance of multi-jet production in $`e^+e^{}`$ annihilation, which is most prominently observed in $`y\times \mathrm{\Phi }`$ space . Indeed, the 1D cumulants in $`p_T`$ are found to be close to zero and therefore not shown.

It is to be noted that oscillations of the cumulants, clearly seen in the data as well as in the Monte Carlo predictions (see e.g. $`K_q(\mathrm{\Phi })`$ in Fig. 1) are not due to statistical fluctuations. They arise from the jet structure of the events in rapidity-azimuthal angle subspace, and from the phase space partitioning technique used to calculate the cumulants as a function of $`M`$.

The 1D cumulants of all-charge and of like-sign multiplets (Fig. 1) show a similar dependence on $`M`$. The latter, however, are significantly smaller, implying that, for all $`M`$, correlations among particles of opposite charge are important in one-dimensional phase space projections. Using Eq. (5), it is found e.g. that the one-dimensional cumulants for unlike-charge pairs are larger by a factor of $`1.61.7`$ than those of like-sign pairs for $`M5`$. This can be expected in general from local charge conservation and in particular from resonance decays.

In 2D and 3D (Fig. 2), like-sign cumulants increase faster and approach the all-charge ones at large $`M`$. From Eq. (5), it can be inferred that $`K_2`$ for unlike-charge pairs remains essentially constant for $`M`$ larger than about 6. Consequently, as the cell-size becomes smaller, the rise of all-charge correlations is increasingly driven by that of like-sign multiplets.

Similar features, but based on measurements of factorial moments, were observed in other experiments  and, mainly on qualitative grounds, taken as evidence for the large influence of the Bose-Einstein effect on all-charge correlations in small phase space domains555For a recent critical review we refer to .. In the next section, we add quantitative support to those earlier observations, showing that algorithms which simulate BEC indeed allow a quite successful description of the measurements.

### 4.2 Model comparison

In this section, we compare the cumulant data with predictions of the Pythia Monte Carlo event generator (version 6.158) without and with Bose-Einstein effects. Samples of about $`10^6`$ multihadronic events were generated at the Z0 energy. The model parameters, not related to BEC, were set at values obtained from a previous tune to OPAL data on event-shape and single-particle inclusive distributions . In this tuning, BE-effects were not included.

To assess the importance of BE-type short-range correlations between identical particles, and their influence on all-charge cumulants, we concentrate on the algorithm BE32, described in the appendix, using parameter values $`\text{PARJ(93)}=0.26`$ GeV ($`R=0.76`$ fm) and $`\text{PARJ(92)}\lambda =1.5`$. These values were determined by varying independently PARJ(93) and $`\lambda `$ within the range $`0.2`$$`0.5`$ GeV and $`0.5`$$`2.2`$, respectively, in steps of $`0.05`$ GeV and $`0.1`$, leaving all other model-parameters unchanged, until satisfactory agreement with the measured cumulants $`K_2`$ for like-sign pairs was reached666Non-BEC related model-parameters were set at the following values: PARJ(21)=0.4 GeV, PARJ(42)=0.52 GeV-2, PARJ(81)=0.25 GeV, PARJ(82)=1.9 GeV.. Whereas a detailed multi-dimensional best-fit tuning procedure was considered to be outside the scope of the paper, we find that calculations with PARJ(93) in the range $`0.20.3`$ GeV, and the corresponding $`\lambda `$ in the range $`1.71.3`$, still provide an acceptable description of the second-order like-sign cumulants.

Additional studies have revealed some sensitivity of single-particle spectra and event shapes to the inclusion of BEC in Pythia. Nevertheless, it was found that minor variations of the QCD and fragmentation parameters are sufficient to restore agreement with these data, while the changes in the predicted cumulants remain well within the systematic uncertainties of the measurement. These parameters were, therefore, not changed from their values listed in footnote 6.

The dashed lines in Figs. 1 and 2 show Pythia predictions for like-sign multiplets for the model without BEC. Model and data agree for small $`M`$ (large phase space domains), indicating that the multiplicity distribution in those regions is well modelled. However, for larger $`M`$, the predicted cumulants are too small, the largest deviations occuring in 2D and 3D. The model predicts negative values for $`K_4(\mathrm{\Phi })`$ which are not shown.

The solid curves in Figs. 1 and 2 show predictions for like-sign multiplets using the BE32 algorithm. Inclusion of BEC leads to a very significant improvement of the data description. Not only two-particle but also higher order correlations in 1D rapidity space are well accounted for. In $`\mathrm{\Phi }`$-space (Fig. 1), $`K_3`$ and especially (the very small) $`K_4`$ are less well reproduced. Figure 2 also shows that the predicted 2D and 3D cumulants agree well with the data.

For more clarity and later reference, the 1D, 2D and 3D cumulants for particle pairs with the same charge are repeated in Fig. 3. Since BEC occur only when two identical mesons are close-by in all three phase space dimensions, projection onto lower-dimensional subspaces, such as rapidity and azimuthal angle, leads to considerable weakening of the effect. Nevertheless, the high precision of the data in Fig. 3 allows to demonstrate clear sensitivity to the presence or absence of BEC in the model, which was much less evident in earlier measurements .

Whereas the BE-algorithm used implements pair-wise BEC only, it is noteworthy that the procedure also induces like-sign higher-order correlations of approximately correct magnitude. This seems to indicate that high-order cumulants are, to a large extent, determined by the second-order one (see further Sect. 4.3). It is not clear, however, whether the agreement is accidental or implies that the physics of $`n`$-boson ($`n>2`$) BE effects is indeed correctly simulated777For a recent theoretical discussion of multi-boson BEC, see ..

The influence of like-sign BE-type correlations on the correlations in all-charge multiplets becomes clear from Figs. 4 and 5 where all-charge cumulants are compared with Pythia predictions without BEC, and for various BE algorithms.

Large discrepancies, already discussed in  and for factorial moments in , are seen for the model without BEC. Inclusion of BE-effects using the BE32 algorithm leads to considerably better agreement, in particular in domains of $`y`$, $`y\times \mathrm{\Phi }`$ and $`y\times \mathrm{\Phi }\times \mathrm{ln}p_T`$. In $`\mathrm{\Phi }`$, disagreements similar to those for like-sign multiplets are also seen here. In addition, $`K_2(\mathrm{\Phi })`$ may be somewhat overestimated for large $`M`$.

To assess the sensitivity of the cumulants to variations in the BEC algorithms available in Pythia, we have further considered the algorithms BEλ and BE0 (see appendix).

Using the same parameter values as for BE32, we observe that BEλ slightly underestimates $`K_2(y)`$ and overestimates $`K_2(\mathrm{\Phi })`$ for like-sign pairs (Fig. 3), whereas the results coincide with those from BE32 in 2D and 3D. For all-charge multiplets (Figs. 4 and 5), the predicted cumulants generally fall below those for BE32, except for $`K_3`$ and $`K_4`$ in 2D and 3D, where the differences are small. We note, in particular, that the predicted $`K_2(y)`$ (Fig. 4) essentially coincides with the Pythia results without BEC. The differences with respect to BE32 are related to the different pair-correlation functions $`g_2(Q)`$ (Eqs. 8 and 9 of the appendix) used in the algorithms. Although a different choice of the parameters $`R`$ and $`\lambda `$ may improve the agreement with the data, we have not attempted such fine-tuning.

In $`W`$-mass studies using the fully hadronic decay channel of the reaction $`e^+e^{}W^+W^{}`$, and assuming that pairs of pions from different $`W`$’s are fully affected by BEC, algorithm BE32 is found to introduce a negative mass-shift, whereas BEλ generates a positive shift . The mass-shift occurs since the BE-induced momentum changes increase the likelihood that soft particles in an event are assigned to the wrong jet. This, in turn, depends on the strength of the particle correlations in 3D momentum space. Since BE32 and BEλ provide a reasonable description of the domain-size dependence of all-charge cumulants in 3D, both algorithms should be considered in $`WW`$ studies.

Finally, we consider the model predictions based on the algorithm BE0 (dash-dotted curves in the figures) for the same parameter values as quoted above. For like-sign pairs (Fig. 3), $`K_2(y)`$ and especially $`K_2(\mathrm{\Phi })`$ are overestimated. This is also the case for $`K_2(\mathrm{\Phi })`$ for all-charge pairs shown in Fig. 4. In contrast, all-charge higher-order cumulants differ little from those obtained with BE32. It should be noted that the BE0 algorithm, contrary to BE32 and BEλ, enforces energy conservation by a global rescaling of all final-state hadron momenta. This procedure affects the full hadronic final state and induces a large artificial shift in the $`W`$-mass when applied to the reaction $`e^+e^{}W^+W^{}\text{hadrons}`$. It is, for that reason, disfavoured by the authors of .

Traditional studies of the Bose-Einstein effect most often consider only the second-order correlation function for like-sign particle pairs, which is measured as a function of $`Q^2`$, the square of the difference in four-momenta $`p_i`$ ($`i=1,2`$) of particles in the pair. The variable $`Q^2`$ is related to the 3D cell-size used in this paper, via the approximate relation :

$$Q^2m^2\left[\beta ^2\delta \mathrm{\Phi }^2+(1+\beta ^2)\delta y^2+\frac{\beta ^4\delta x^2}{(1+\beta )^2}\right],$$
(6)
valid for small domains in $`y`$, $`\mathrm{\Phi }`$ and $`\mathrm{ln}p_T`$. Here $`\beta =p_T/m`$, $`m`$ is the particle mass and $`\delta x=\mathrm{ln}(p_{T1}/p_{T2})`$. It follows that a small cell-size in 3D corresponds to small $`Q^2`$. The reverse is not true. Equation (6) implies that a measurement of $`K_2`$ versus $`M`$ in 3D is roughly equivalent to that of the inclusive two-particle $`Q`$-distribution.

As a check of the consistency of our results, we have compared predictions for the inclusive two-particle $`Q`$-distribution of like-sign pairs (single-particle rapidity restricted to the interval $`2y2`$) from the BE32 model to that extracted from the data and fully corrected for detector effects. For the values of the BEC parameters $`R`$ and $`\lambda `$ quoted above, satisfactory agreement is obtained (not shown).

In recent work, several LEP experiments have analyzed the second-order like-sign correlation function in terms of individual components of the four-vector $`Q`$ . If interpreted as a measurement of the space-time extent of the particle emitting source, such results show that the emission volume has an elongated shape with respect to the event axis. The BE algorithms discussed here, which treat all components of $`Q`$ symmetrically, are unable to reproduce these measurements . It is therefore possible, although to be verified, that the measured anisotropy of the two-particle correlation function may be responsible for some of the Pythia model discrepancies, e.g. for 1D cumulants in $`\mathrm{\Phi }`$ domains, observed in the present analysis.

In summary, a comparison with Pythia predictions shows that short-range correlations of the BE-type are needed, at least in this model, to reproduce the magnitude and the $`\mathrm{\Delta }`$-dependence of the cumulants for like-sign multiplets. This further leads to a much improved description of not only two-particle but also of higher-order correlations in all-charge multiplets. Since Bose-Einstein correlations are a well-established phenomenon in multiparticle production, it is likely that the above conclusion has wider validity than the model from which it was derived.

### 4.3 The Ochs-Wosiek relation for cumulants

The success of the Pythia model with BEC in predicting both the magnitude and domain-size dependence of cumulants, has led us to consider the inter-dependence of these quantities in more detail.

In Fig. 6 we plot $`K_3`$ and $`K_4`$ in 2D and 3D, for each value of $`M`$ as a function of $`K_2`$. We observe that the 2D and 3D data for all-charge, as well as for like-sign multiplets follow approximately, within errors, the same functional dependence. The solid lines is a simple fit to the function

$$\mathrm{ln}K_q=a_q+r_q\mathrm{ln}K_2.$$
(7)
The fitted slope values are $`r_3=2.3`$ and $`r_4=3.8`$. This is evidence that the slope $`r_q`$ increases with the order of the cumulant. The $`q`$-dependence of the slopes is of particular interest in multiplicative cascade models of hadronisation, and indicative of the mechanism causing scale-invariant fluctuations .

Figure 6 suggests that the cumulants of different orders obey simple so-called “hierarchical” relations, analogous to the Ochs-Wosiek relation, first established for factorial moments . Interestingly, all-charge as well as like-sign multiplets are seen to follow, within errors, the same functional dependence. Hierarchical relations of similar type are commonly encountered, or conjectured, in various branches of many-body physics (see e.g. ), but a satisfactory explanation within particle production phenomenology or QCD remains to be found.

Simple relations among the cumulants of different orders exist for certain probability distributions, such as the Negative Binomial distribution, which often parameterizes successfully the multiplicity distribution of hadrons in restricted phase space regions . For this distribution, one has $`K_q=(q1)!K_2^{q1}`$ ($`q=3,4,\mathrm{}`$), showing that the cumulants are here solely determined by $`K_2`$. This relation is shown in Fig. 6 (dashed line). Comparing to the data, we conclude that the multiplicity distribution of all charged particles, as well as that of like-sign particles, deviates strongly from a Negative Binomial in small phase space domains. For further discussion, in the present context, of this and other much studied multiplicity distributions we refer to .

The Ochs-Wosiek type of relation exhibited by the data in Fig. 6 may explain why the BE algorithms in Pythia generate higher-order correlations of (approximately) the correct magnitude. Assuming that the hadronization dynamics is such that higher-order correlation functions can be constructed from second-order correlations only, methods that are designed to ensure agreement with the two-particle correlation function, could then automatically generate higher-order ones of the correct magnitude.

## 5 Summary and conclusions

In this paper we have presented a comparative study of like-sign and all-charge genuine correlations between two and more hadrons produced in $`e^+e^{}`$ annihilation at the $`Z^0`$ energy. The high-statistics data on hadronic $`\mathrm{Z}^0`$ decays recorded with the OPAL detector from 1991 through 1995 were used to measure normalized factorial cumulants as a function of the domain size, $`\mathrm{\Delta }`$, in $`D`$-dimensional domains ($`D=1,2,3`$) in rapidity, azimuthal angle and (the logarithm of) transverse momentum, defined in the event sphericity frame.

Both all-charge and like-sign multiplets show strong positive genuine correlations up to fourth order. They are stronger in rapidity than in azimuthal angle. One-dimensional cumulants initially increase rapidly with decreasing size of the phase space cells but saturate rather quickly. In contrast, 2D and especially 3D cumulants continue to increase and exhibit intermittency-like behaviour.

Comparing all-charge and like-sign multiplets in 2D and 3D phase space cells, we observe that the rise of the cumulants for all-charge multiplets is increasingly driven by that of like-sign multiplets as $`\mathrm{\Delta }`$ becomes smaller. This points to the likely influence of Bose-Einstein correlations.

The 2D and 3D cumulants $`K_3`$ and $`K_4`$, considered as a function of $`K_2`$, follow approximately a linear relation of the Ochs-Wosiek type: $`\mathrm{ln}K_q\mathrm{ln}K_2`$, independent of $`D`$ and the same for all-charge and for like-sign particle groups. This suggests that, for a given domain $`\mathrm{\Delta }`$, correlation functions of different orders are not independent but determined, to a large extent, by two-particle correlations.

The data have been compared with predictions from the Monte Carlo event generator Pythia, previously tuned to single-particle and event-shape OPAL data. The model describes well dynamical fluctuations in large phase space domains, e.g. caused by jet production, and shorter-range correlations attributable to resonance decays. However, the results of this paper, together with earlier less precise data, show that these ingredients alone are insufficient to explain the magnitude and domain-size dependence of the factorial cumulants. To achieve a more satisfactory data description, short-range correlations of the Bose-Einstein type between identical particles need to be included.

The importance of BE-type effects has been studied using algorithms implemented in Pythia. We find that the model BE32 is able to simultaneously account for the magnitude and $`\mathrm{\Delta }`$-dependence of like-sign as well as of all-charge cumulants. Other models, BE0 and BEλ, when using the same parameters as for BE32, provide also a reasonable description of the data. Although the algorithms implement pair-wise BEC only, surprisingly good agreement with the measured third- and fourth-order cumulants is observed. This could be a consequence of the Ochs-Wosiek type of inter-relationship between cumulants of different orders, exhibited by the data.

Acknowledgements:

We particularly wish to thank the SL Division for the efficient operation of the LEP accelerator at all energies and for their close cooperation with our experimental group. We thank our colleagues from CEA, DAPNIA/SPP, CE-Saclay for their efforts over the years on the time-of-flight and trigger systems which we continue to use. In addition to the support staff at our own institutions we are pleased to acknowledge the
Department of Energy, USA,
National Science Foundation, USA,
Particle Physics and Astronomy Research Council, UK,
Natural Sciences and Engineering Research Council, Canada,
Israel Science Foundation, administered by the Israel Academy of Science and Humanities,
Minerva Gesellschaft,
Benoziyo Center for High Energy Physics,
Japanese Ministry of Education, Science and Culture (the Monbusho) and a grant under the Monbusho International Science Research Program,
Japanese Society for the Promotion of Science (JSPS),
German Israeli Bi-national Science Foundation (GIF),
Bundesministerium für Bildung und Forschung, Germany,
National Research Council of Canada,
Research Corporation, USA,
Hungarian Foundation for Scientific Research, OTKA T-029328, T023793 and OTKA F-023259.

## Appendix: Modelling BEC in PYTHIA

As emphasized in , no rigorous method exists to date which would allow to account for the Bose-Einstein symmetrization of the production amplitudes in Monte Carlo event generators. The algorithms implemented in the Monte Carlo event generator Pythia  which simulate the BE effect, are all based on introducing BEC as local shifts of final-state particle momenta among pairs of identical particles. They differ only in the way global energy and momentum conservation is ensured. Other, so-called “global-weight” methods, have been proposed wherein a single weight-factor is assigned to each event, calculated on the basis of a model for the identical-pair two-particle correlation function .

The Pythia algorithms take the hadrons produced by the string fragmentation, where no BE effects are present, and shift the momenta of mesons $`i`$ and $`j`$ such that the inclusive distribution of the relative separation $`Q`$ of identical pairs is enhanced by a factor $`g_2(Q)1`$. The latter, as used here, is parameterized with the phenomenological form888In naive BE models, based on the analogy with optics and assuming a static incoherent emission source, $`R`$ is related to the size of the source; $`\lambda `$ quantifies the strength of the BE effect.

$$g_2(Q)=1+\lambda e^{R^2Q^2},$$
(8)
where $`Q`$ is the difference in four-momenta of the pair, $`Q^2=(p_1p_2)^2`$.

Short-lived resonances like the $`\rho `$ and $`K^{}`$ are allowed to decay before the BE procedure is applied, while decay products of longer-lived ones (width $`\mathrm{\Gamma }<20`$ MeV/$`c^2`$) are not affected. The procedure also influences groups of non-identical particles causing e.g. a shift of the $`\rho ^0`$ mass peak .

In this paper we compare the data to simulations based on the Pythia algorithms BE0, BE32 and BEλ, following the nomenclature in .

In the BE0 algorithm, momentum is conserved exactly, but energy conservation is explicitly broken in the treatment of individual particle pairs. It is restored only by a global rescaling of all final-state hadron momenta, thus affecting the full hadronic final state.

The other algorithms avoid global energy rescaling by introducing additional momentum shifts for some pairs of particles, ($`k`$, $`l`$), not necessarily identical, in a local region around each identical pair ($`i`$, $`j`$).

The algorithm $`\text{BE}_{32}`$, which is applied to identical pairs only, is based on the ansatz

$$g_2(Q)=\left\{1+\lambda \mathrm{exp}(Q^2R^2)\right\}\left\{1+\alpha \lambda \mathrm{exp}(Q^2R^2/9)\left(1\mathrm{exp}(Q^2R^2/4)\right)\right\},$$
(9)
($`\alpha `$ is an adjustable parameter) which attempts to mimic effects due to oscillations below and above unity of the pair-weights, as found in some models .

In the algorithm $`\text{BE}_\lambda `$, the form (8) of $`g_2(Q)`$ is retained. For each pair of identical particles ($`i,j`$), a pair of non-identical particles, ($`k,l`$), neither identical to $`i`$ or $`j`$, is found close to ($`i,j`$). For each momentum shift among particles $`i`$ and $`j`$, a corresponding shift among the particles $`k`$ and $`l`$ is found so that the total energy and momentum in the $`i,j,k,l`$ system is conserved. As a measure of “closeness”, $`\text{BE}_\lambda `$ uses the so-called $`\lambda `$-measure , related to the string length in the Lund string fragmentation framework.

A given particle is likely to belong to several identical pairs in an event. The net shift in particle momenta and energies therefore depends on the complete configuration of all identical particles. This introduces complex effects, not only among identical pairs but also among unlike-sign pairs and higher multiplicity multiplets of nearby particles. As a result, genuine higher-order correlations emerge. For a full description of the BE algorithms mentioned above, ref.  should be consulted.


format: text

text:
# De Broglie Wavelength Reduction for a Multi-photon Wave Packet

## ACKNOWLEDGMENTS

I wish to thank Janne Ruostekoski for lively discussions.


text:
# Lattice QCD at finite temperature: Evidence for calorons from the eigenvectors of the Dirac operatortalk given by W. Söldner at Lattice 2001, Berlin

## 1 INTRODUCTION

The chiral phase transition is a main feature of QCD. In the instanton picture  the chirally broken phase is represented by an ensemble of weakly interacting instantons and anti-instantons. Due to the fermionic quasi zero modes of the (anti-) instantons, a finite density of eigenvalues, $`\rho _{\mathrm{Dirac}}`$, is generated near zero in a (large) volume $`V`$, and according to the Banks-Casher formula

$$\overline{q}q=\pi \rho _{\mathrm{Dirac}}(0)/V$$
(1)
chiral symmetry is broken. In the chirally symmetric phase instantons and anti-instantons form strongly bound molecules, the eigenvalues of the Dirac operator move towards the bulk of the spectrum and the eigenvalue density vanishes near the origin. Of course one would like to verify the essential parts of this picture in QCD. Therefore we study the localization properties of the eigenmodes of the lattice Dirac operator in the neighborhood of the chiral phase transition. In a previous paper  we have investigated this question with the Wilson gauge action on smaller lattices. Here we present results obtained on larger lattices with the Lüscher-Weisz gauge action.

## 2 PROPERTIES OF CALORONS

The temporal extension $`L_t`$ of the lattice is related to the temperature $`T`$ by $`aT=1/L_t`$, where $`a`$ is the lattice spacing. Due to the periodic boundary condition in the time direction, an instanton becomes a chain of instantons which one calls caloron . The properties of a caloron are similar to that of an instanton: a caloron has topological charge one and a localized left-handed fermionic zero mode. However, the localization pattern of the zero mode is different. Whereas in the case of an instanton the zero mode is localized in space and time, the zero mode of a caloron is localized in space but delocalized in time.

Let us now compare the localization of the zero mode in the different sectors of the Polyakov loop $`P`$. We know that in the deconfined phase $`P`$ clusters around the values $`\mathrm{e}^{\mathrm{i}\theta _P}`$, with $`\theta _P=0`$ (real sector) and $`\theta _P=\pm 2\pi /3`$ (complex sector) . Because the fermionic action does not share the $`Z_3`$ symmetry of the gauge action, the localization of the eigenmodes of the Dirac operator can depend on the Polyakov loop sectors. Indeed one finds that the zero mode $`\psi _0`$ of the caloron is more localized in the real sector than in the complex sector, as can be seen from

$$\left|\psi _0\right|^2\mathrm{exp}[2(\pi |\theta _P|)rT]/r^2.$$
(2)
Next we will search for caloron-like configurations concentrating on the fermionic properties of the caloron described above.

## 3 LATTICE SETUP

We work in the quenched approximation with the Lüscher-Weisz gauge action

$`S_g[U]`$ $`=`$ $`\beta _1_{pl}\frac{1}{3}\mathrm{Re}\mathrm{Tr}[1U_{pl}]+`$ (3)
$`\beta _2_{rt}\frac{1}{3}\mathrm{Re}\mathrm{Tr}[1U_{rt}]+`$
$`\beta _3_{pg}\frac{1}{3}\mathrm{Re}\mathrm{Tr}[1U_{pg}],`$
where $`_{pl}`$, $`_{rt}`$ and $`_{pg}`$ means a summation over all plaquettes, $`2\times 1`$ rectangles and parallelograms, respectively. The coefficients $`\beta _2`$ and $`\beta _3`$ are computed from $`\beta _1`$ via tadpole improved pertubation theory. Further we use the staggered Dirac operator

$$D=\underset{\mu =1}{\overset{4}{}}\frac{1}{2a}\alpha _\mu (x)\left[\delta _{y,x+\widehat{\mu }}U_\mu (x)\delta _{y,x\widehat{\mu }}U_\mu ^{}(y)\right],$$

where $`\alpha _\mu (x)=(1)^{x_1+\mathrm{}+x_{\mu 1}}`$. The eigenvalues of $`D`$ come in pairs of $`\pm \mathrm{i}\lambda `$ with $`\lambda `$ real, so we can restrict ourselves to positive $`\lambda `$ in the following. In the continuum limit this action corresponds to four quark flavors. From now on we set $`a`$ to 1.

We have calculated the lowest eigenvalues and eigenvectors for lattice sizes $`12^3\times 6`$, $`16^3\times 6`$ and $`20^3\times 6`$ below the deconfinement phase transition at $`\beta _1=8.10`$ ($`\beta _{\mathrm{Wilson}}5.8`$ ) and above at $`\beta _1=8.45`$ ($`\beta _{\mathrm{Wilson}}6.0`$ ) with the Arnoldi method . We will present data for the largest lattice at $`\beta =8.45`$.

## 4 CALORONS ON THE LATTICE

When searching for calorons on the lattice with the help of their fermionic properties, the zero mode of the Dirac operator plays the central role. Unfortunately the staggered Dirac operator has no exact zero modes, but, as we can see clearly in the left plot of Fig. 1, there are for both Polyakov loop sectors low-lying eigenmodes with a high value of $`\left|\gamma _5\right|`$ and eigenmodes in the bulk of the spectrum with values of $`\gamma _5`$ around zero. We find that the low-lying eigenvalues are quasi 4-fold degenerate in the sense that the number of low-lying eigenmodes is almost always 4, 8, 12, …(In Fig. 1 this is illustrated by the encircled eigenmodes, which belong to the same gauge field configuration. Note that the total number of low-lying eigenmodes of that configuration is 8.) So one would expect that those low-lying eigenmodes correspond to exact zero modes with definite handedness in the continuum.

Let us now look at the localization properties of these approximate zero modes. We define a gauge-invariant measure of the localization of a quark eigenmode $`\psi _\lambda ^\alpha (x)`$ ($`\lambda `$ is the Dirac eigenvalue, $`\alpha `$ a color index), the inverse participation ratio

$$I_2V\frac{_xp_\lambda (x)^2}{\left[_xp_\lambda (x)\right]^2},$$
(4)
where $`V`$ is the number of lattice sites and $`p_\lambda (x)`$ is the gauge-invariant probability density $`p_\lambda (x)=_{\alpha =1}^{N_c}\left|\psi _\lambda ^\alpha (x)\right|^2`$. For a completely delocalized state (all $`p_\lambda (x)`$ the same) one finds $`I_2=1`$, whereas a state localized on a single lattice site (only one non-zero $`p_\lambda (x)`$) would have $`I_2=V`$. In the right part of Fig. 1 we plot $`I_2`$ vs. $`\left|\gamma _5\right|`$. We see that the modes in the real sector are, in general, more localized than those in the complex sector. In both Polyakov loop sectors the approximate zero modes have a greater value of $`I_2`$ than the modes in the bulk spectrum, which means that they are more localized. These are exactly the properties we expect for caloron configurations.

Finally we consider the localization pattern. In Fig. 2 isosurfaces of $`p_\lambda `$ and the local chirality $`p_\lambda ^5(x)=_{\alpha =1}^{N_c}\psi _\lambda ^\alpha (x)^{}\gamma _5\psi _\lambda ^\alpha (x)`$ are plotted at fixed $`t`$ and at fixed $`z`$ for one of the encircled eigenmodes. This (typical) mode is clearly localized in space but not in time and the isosurfaces for $`p_\lambda `$ and $`p_\lambda ^5`$ are very similar. This is again as is expected for caloron-like configurations.

## 5 CONCLUSIONS

Studying the localization and chirality properties of low-lying quark eigenstates in quenched lattice QCD, we could characterize semiclassical properties of the gauge field configurations without any cooling. For temperatures above $`T_c`$ we found isolated modes with definite handedness. They are localized in space but not in time with characteristic differences between the Polyakov loop sectors. Thus they show essential properties of quark states associated with calorons.
Acknowledgements: We thank the LRZ in Munich for computer time on the Hitachi SR8000. Support from BMBF and DFG is acknowledged.


text:
# Electronic and magnetic structure of CsV2O5

## I Introduction

Low-dimensional spin-gapped quantum systems are of current interest since they show interesting ground states and a variety of unconventional low-lying excitations. Examples of such systems include spin-1 Haldane chains , spin-$`\frac{1}{2}`$ even-leg ladders and spin-$`\frac{1}{2}`$ alternating chains . The discovery of appropriate compounds like Y2BaNiO5 (Haldane chain), Sr14Cu24O41  or SrCu2O3  (ladder systems) has brought new insight into the study of these systems. There is a long list of various spin-$`\frac{1}{2}`$ alternating chain systems which are being intensively studied in connection to their magnetic excitations, to mention a few of them KCuCl3 , TlCuCl3 , Cu(NO3)2.$`\frac{5}{2}`$D2O . If frustration is also present new features appear in the magnetization spectrum of the low-dimensional quantum systems as in CuGeO3 , (VO)2P2O7 (VOPO) , SrCu2(BO3)2  or in the recently synthesized Cu2Te2O5Br2 .

An important family of low-dimensional compounds are the layered vanadates AV2O5 (A=Ca, Mg, Na, Li, Cs). While CaV2O5 and MgV2O5 contain only magnetic V4+ ions and behave like spin-$`\frac{1}{2}`$ ladders with spin-gaps of the order of 600K  and 20K  respectively, $`\alpha ^{}`$-NaV2O5, $`\gamma `$-LiV2O5 and CsV2O5 are mixed-valence systems (V4.5+ on average) with important charge and spin fluctuations. This last property makes this family specially interesting since it allows to study the influence of charge-ordering and the corresponding crystallographic distortions on the magnetic interactions. The quarter-filled ladder compound $`\alpha ^{}`$-NaV2O5  is a highly discussed material and has become a model substance for the study of spin, charge and orbital coupling. It behaves like a spin-$`\frac{1}{2}`$ Heisenberg system at high temperatures and at 34K undergoes a charge ordering (2V4.5+ $``$ V4+ \+ V5+) transition accompanied with a lattice distortion  and the opening of a spin gap. $`\gamma `$-LiV2O5 is a charge ordered system with no spin gap opening at low temperatures . Recent investigations based on an ab initio analysis of the electronic structure show that charge fluctuations in this system are still significant and the experimental observations  can be well described by a spin-$`\frac{1}{2}`$ asymmetric quarter-filled ladder model. Such proposal has also found support by a numerical analysis of optical conductivity experiments. The third member of this vanadate family, CsV2O5, has been given much less attention. Existing susceptibility measurements of Isobe et al. could be rather well explained by considering this material as a system of isolated dimers with a spin-gap of $`12meV`$. Measurements other than that are unknown to us except for older susceptibility measurements  which somewhat disagree with Isobe et al.’s data and different values of the exchange integral $`J`$ and the gyromagnetic ratio $`g`$ are predicted.

There are various reasons why we believe that this system deserves attention: (i) the question of how the electronic and magnetic properties of this material compare with the properties of the other two members of the alkali-vanadate family, i.e. $`\alpha ^{}`$-NaV2O5 and $`\gamma `$-LiV2O5 (ii) this system crystallizes in a similar monoclinic structure as the most discussed VOPO and many other materials e.g. KCuCl3, TlCuCl3 showing an alternating chain behavior of the magnetic interactions. It is then worthwhile to find out whether CsV2O5 has also a similar behavior of the magnetic interactions as its isostructural members. (iii) as will be shown below, a re-examination of the susceptibility measurements reveals that a dimer model is not the only model that can provide reasonable fitting of the susceptibility data.

In view of the above, we have performed a microscopic study of this system by considering Density-Functional-Theory-based (DFT) calculations. The microscopic study of the electronic behavior of this compound reveals that the inter-dimer coupling in this material in addition to the intra-dimer interaction is also significant and shouldn’t be neglected. Further, an analysis of the susceptibility data shows that getting estimates of the coupling constants and the gyromagnetic factor based merely on a fitting procedure is very subtle and with such an analysis one can only test the consistency with the assumed model, but cannot prove the uniqueness of that model. Results on that will be presented in the next sections. We shall also discuss the similarities and differences to the other members of this vanadate family, $`\alpha ^{}`$-NaV2O5 and $`\gamma `$-LiV2O5 as well as to isostructural compounds like VOPO.

## II Crystal Structure

The crystal structure of CsV2O5 is somewhat different from that of $`\alpha ^{}`$-NaV2O5 and $`\gamma `$-LiV2O5. $`\alpha ^{}`$-NaV2O5 and $`\gamma `$-LiV2O5 crystallize in the orthorhombic space group $`Pmmn`$ ($`D_{2h}^{13}`$) and $`Pnma`$ ($`D_{2h}^{16}`$) respectively. They consist of a double chain-structure of edge-sharing distorted square VO5 pyramids running along the orthorhombic $`y`$ axis, which are linked together via common corners to form layers. These layers are stacked upon each other along $`z`$. Na/Li atoms are located between these layers. CsV2O5 crystallizes, on the other hand, in the monoclinic space group $`P2_1/c`$ ($`C_{2h}^5`$). CsV2O5 has also a layered structure with Cs ions in between the layers. These layers, which are somewhat tilted ( i.e. not strictly in yz plane), and are stacked upon each other along $`x`$, are made of two types of vanadium atoms V1 and V2. The V1 atoms are in a distorted square pyramidal environment of oxygens and form pairs by sharing the edge of the pyramids. These pairs are separated from one another by the V2 vanadiums in a tetrahedral co-ordination of oxygens as shown in Fig. 1. Within a pair of square pyramids, the two apical O atoms are pointing in opposite directions and the remaining O atoms of the pyramids are shared with tetrahedra. The V1-V1 distance across the shared edge of the square-pyramids-pair is 3.073 $`\AA `$ and the V1-V2 distance across the shared corner of a pyramid and a tetrahedron is 3.352 $`\AA `$. Since the average V2-O distance within the tetrahedron d=1.718 $`\AA `$ is shorter than the average V1-O distance within the square pyramid d=1.882 $`\AA `$, one can conclude from a bond-valence analysis  that the V1 sites are in the oxidation V4+ and the V2 sites are in the oxidation V5+. As the intra-pair V4+-V4+ distance is 3.073 $`\AA `$ while the inter-pair distance along the $`z`$-axis is 5.501 $`\AA `$, this compound has been assumed to be a system of isolated dimers.

## III Band Structure

In order to analyze the electronic behavior of CsV2O5, we carried out DFT calculations within the Local Density Approximation (LDA) by employing the full-potential linearized augmented plane wave code WIEN97  (LAPW) and by linear-muffin-tin orbital (LMTO) based on the Stuttgart TBLMTO-47 code. Both calculations agree well within the allowed error-bars of the various different approximations involved in two different methods. In both approaches we treated the exchange-correlation part by using the generalized gradient approximation .

In Fig. 2 the energy bands for CsV2O5 along the path $`\mathrm{\Gamma }BDZ\mathrm{\Gamma }YE`$ are shown. Since the unit cell of CsV2O5 contains four formula units, i.e. four V4+ and four V5+, there are four bands around the Fermi level which correspond to V4+-3d orbitals. These bands are separated by an energy gap of about 2 eV from the lower O-p valence bands and an energy gap of about 0.5 eV from the upper conduction bands. The conduction bands up to 4eV are of V-3d nature. There are a few points to be noted here: (i) The energy range on the CsV2O5 bands is of the same order of magnitude as it is in $`\alpha ^{}`$-NaV2O5 and $`\gamma `$-LiV2O5. (ii) The contribution of V5+-3d is mostly in the upper conduction bands, starting from about 0.75 eV above the Fermi level and have a very small contribution at the Fermi level indicating that, contrary to $`\alpha ^{}`$-NaV2O5 and $`\gamma `$-LiV2O5, here the system is close to complete charge ordering. (iii) There is no dispersion along the $`x`$ axis while along the $`z`$ axis (chain direction) the dispersion is non-negligible and of about 300meV. Along the $`y`$ axis the LDA-bands show a small dispersion which corresponds to inter-chain interactions. In order to investigate the dimension-dependent behavior of CsV2O5, we show in Fig. 3(a) the partial density of states (DOS) of the V4+ 3d-orbitals for CsV2O5. For comparison we have also plotted in Fig. 3(b) the partial density of states of V1-3d orbitals for $`\gamma `$-LiV2O5 and in Fig. 3(c) that of V-3d orbitals for $`\alpha ^{}`$-NaV2O5. Note that, while the $`\gamma `$-LiV2O5 and $`\alpha ^{}`$-NaV2O5 DOS show the characteristic quasi-1D van Hove singularities near the band edges around the Fermi level, this is not the case for the Cs-compound. The DOS for CsV2O5 shows more a 2D-behavior than 1D. (iv) In the global frame of reference, the character of the bands around the Fermi level turn out to be that of d$`_{z^2}`$ with small admixture from d$`_{x^2y^2}`$ character. Rotation to local co-ordinate system with the $`z`$-axis aligned along the (V-apical O) direction for individual pyramids transform this to dxy character, in agreement with the crystal field analysis of a V4+ ion in the pyramidal environment.

The band-structure analysis leads us to the conclusion that inter-dimer coupling is non-negligible. If the system would have been a pure dimer system we would have expected no dispersion in any direction. The existence of a non-negligible dispersion tells us that a model of unlinked dimers is too simple and there must be some inter-dimer coupling to explain the band picture. From the DOS analysis we learn that inter-dimer interactions in both $`z`$ and $`y`$ direction are important for the electronic properties. The possible role of the V5+O4 tetrahedra in this context will be investigated in the next section.

## IV Microscopic model (Effective model)

In order to understand the nature of the ground state in CsV2O5 we need to determine the appropriate microscopic model which explains the low-energy physics of this compound. We have therefore applied the so-called downfolding method offered by the new and generalized version of the LMTO  method, together with a tight-binding (TB) analysis on the band-structure results. The downfolding technique consists of integrating out the high energy degrees of freedom so as to describe the details of LDA energy bands close to the Fermi energy in terms of few-orbital effective Hamiltonians. From this analysis, we can extract the effective hopping matrix elements $`t_{ij}`$ between vanadium ions in CsV2O5 (by Fourier Transform of the downfolded Hamiltonian $`H(k)H(R)`$ ) which reproduce the behavior of the LDA bands close to the Fermi level. Since, as seen in the previous section, the contribution of the V5+-3d orbitals around the Fermi level is very small, and the biggest contribution comes from the V4+-3$`d_{xy}`$ orbital, a one-orbital effective tight-binding Hamiltonian should give a good description of the low energy physics. The tight-binding Hamiltonian therefore, can be written as:

$`H_{TB}={\displaystyle \underset{<i,j>}{}}t_{ij}(c_j^+c_i+c_i^+c_j)`$ (1)
where $`i`$ and $`j`$ denote a pair of V4+ ions labeled each one by $`(n,s)`$ where $`n`$ denotes the unit cell and $`s=1,2,3,4`$ V4+ ions in the unit cell. In Fig. 4 all the considered hoppings $`t_1`$, $`t_2`$, $`t_3`$ and t5 have been drawn. The subscript $`1`$, $`2`$, $`3`$ and $`5`$ indicate respectively that the hopping integral is between first nearest neighbors (n.n.), second n.n., third n.n and fifth n.n. The hopping term $`t_4`$, which connects two V4+ ions belonging to two different layers along $`x`$, has been neglected due to its vanishingly small value. Note here that no hopping terms have been considered along the $`x`$ direction, since as we have learnt from the band picture, there is no dispersion along this direction and that means that the interlayer coupling along $`x`$ must be negligible. In Fig. 5 we show a comparison of the downfolding-TB bands to the LDA bands. We observe that not only the hopping corresponding to the intra-dimer interaction $`t_1`$ is important but also the inter-dimer matrix element $`t_3`$ along $`z`$ is considerable. A model with only $`t_1`$ and $`t_3`$, i.e. alternating chain, shows a dispersion relation of the form:

$`E=\pm \sqrt{t_1^2+t_3^2+2t_1t_3cosk_z}`$ (2)
which is symmetric with respect to E=0 (i.e. to $`E_F`$) and doesn’t account for the band-splittings in the paths $`\mathrm{\Gamma }D`$ and $`\mathrm{\Gamma }Y`$. These features are described by the inter-chain hopping parameters along the $`y`$ direction, $`t_2`$ and $`t_5`$, though they are much smaller than $`t_1`$ and $`t_3`$. Note that the TB-results show an artificial band-crossing along the path $`BD`$ and $`Z\mathrm{\Gamma }`$. We have analyzed this result and seen that the inclusion of longer-ranged interactions lifts this artificial crossing.

Concerning the nature of the hopping integrals, we observe that the values of the dominant hopping parameter $`t_1=0.12eV`$ (which agrees with the overall band-width of the V4+-dxy derived low-energy bands in CsV2O5) is smaller than the dominant hopping parameter for $`\alpha ^{}`$-NaV2O5 $`t_a=0.37eV`$ and $`\gamma `$-LiV2O5 $`t_a=0.35eV`$. This difference can be explained by the nature of the path. While $`t_1`$ in CsV2O5 corresponds to a $`VOV`$ path in between two edge-sharing square pyramids, $`t_a`$ in $`\alpha ^{}`$-NaV2O5 and $`\gamma `$-LiV2O5 corresponds to a $`VOV`$ path between two corner-sharing square pyramids. This feature has important implications in the strength of the exchange interaction as has been pointed out in Ref. . Also note that the same kind of edge-sharing pyramids are also present in $`\alpha ^{}`$-NaV2O5 and $`\gamma `$-LiV2O5 with effective hopping values smaller than $`t_a`$ and similar to $`t_1`$ in CsV2O5.

About the nature of the rest of the hopping integrals in CsV2O5, we observe that $`t_2`$, $`t_3`$ and $`t_5`$ correspond to V4+-V4+ paths through the bridging V5+O4 tetrahedra with different lengths and angles. Beltran et al. studied, based on geometrical considerations, the various super-exchange paths for the case of vanadyl phosphates. In particular they considered the VOPO compound. VOPO and CsV2O5 have a similar monoclinic structure if one identifies the P5+O4 tetrahedra and the V4+O5 square pyramids in VOPO with the V5+O4 tetrahedra and the V4+O5 square pyramids in CsV2O5 respectively. As argued by Beltran et al., based on the formation of coherent molecular electron orbitals, the coupling between two V4+ through bridging oxygens for separated pyramids is weaker than those via the pathways involving V5+ ions. Furthermore, double bridging modes, involving two V5+ ions, as in $`t_2`$ and $`t_3`$ (see Fig. 4) can be stronger than single bridging modes involving one V5+ ion as in $`t_5`$. These interactions can be as strong as paths of the type $`t_1`$, i.e. V4+-O-V4+ where the square pyramids are linked to each other by sharing an edge. Geometry considerations then lead to a relative order of importance among the various bridging modes, which for CsV2O5 is observed to be $`t_3>t_5>t_2`$. We note here that while in CsV2O5 we have seen that the dimers are formed by the edge-sharing pyramids ($`t_1`$) with the predominant inter-dimer interaction formed by V4+-V4+ interaction bridged by V5+O4 ($`t_3`$), in VOPO these roles are interchanged in the sense that dimers are formed by V4+O5 square pyramids bridged by P5+O4 groups ($`t_3`$) and the inter-dimer interaction is formed by neighboring V4+O5 pyramids ($`t_1`$). This is due to the different type of distortion in the square pyramids as well as different angle and distance relations between these two compounds.

Using the relation $`J_1=4t_1^2/U`$, $`U`$ being the on-site Coulomb electron-electron interaction and considering U$``$ 2.8 eV  the DFT calculation gives an estimate of the intra-dimer exchange coupling, $`J_1225K`$.

To summarize, the conclusions that we draw from the above analysis is that ab initio calculations support an alternating chain behavior for CsV2O5 along $`z`$ with weaker inter-chain interactions along $`y`$.

## V Susceptibility

In order to test the results obtained from the ab initio calculation we only have available two sets of susceptibility data  as a function of temperature which don’t completely agree quantitatively with each other. We consider here the data by Isobe et al. These authors proposed that their experimental data are well fitted by the susceptibility of a spin-$`\frac{1}{2}`$ dimer system with $`H=J𝐒_1𝐒_2`$ given by:

$`\chi _{raw}=\chi _{CsV_2O_5}+\chi _{imp}`$ (3)
$`={\displaystyle \frac{Ng^2\mu _B^2}{k_BT}}{\displaystyle \frac{1}{3+\mathrm{exp}(J/k_BT)}}+\chi _0+\chi _{imp}`$ (4)
with a g-factor of 1.8, $`J=146K`$ and $`\chi _0=8\times 10^5emu/mol`$. We have re-analyzed the data by considering the temperature dependence of the susceptibility for a spin-$`\frac{1}{2}`$ alternating Heisenberg chain with:

$`H={\displaystyle \underset{n}{}}(J_1𝐒_{2n}𝐒_{2n+1}+J_2𝐒_{2n+1}𝐒_{2n+2}).`$ (5)
with $`J_1>0`$ and $`J_2>0`$ (antiferromagnetic couplings) obtained by the stochastic Quantum Monte Carlo method . We define the parameter $`\alpha =J_2/J_1`$ which measures the ratio of the inter and intra-dimer coupling. $`\alpha =0`$ corresponds to the dimer model. We introduce also an additional parameter $`\delta `$ so that $`J_1=J(1+\delta )`$, $`J_2=J(1\delta )`$ and therefore $`\alpha =\frac{1\delta }{1+\delta }`$.

In Fig. 6(a) we show the comparison of the experimental data to the dimer-model susceptibility given by Eq. 3. We find, in fact, that the best fit corresponds to $`J=146K`$, $`\chi _0=8\times 10^5emu/mol`$ and a value of $`g=1.77`$, somewhat smaller than the one proposed in Ref. . Figs. 6 (b) and (c) show the comparison of the experimental data to the spin-$`\frac{1}{2}`$ alternating Heisenberg chain model for $`\delta =0.8,\alpha =0.11`$ and $`\delta =0.6,\alpha =0.25`$ respectively. We observe that values of $`\alpha `$ up to $`0.25`$ give still very good fits to the susceptibility data by choosing appropriate $`g`$ and $`J`$ parameters. This analysis makes plausible our ab initio results, namely the fact that the inter-dimer interactions can be also significant.

The important conclusion that we gain from this susceptibility analysis is that it is a very hard task to determine the behavior of a compound only from susceptibility fits. Extra information is needed so as to minimize the possible options. For instance, electron spin -resonance (ESR) experiments could complement the susceptibility measurements by delivering the value of the gyromagnetic factor. There are in the literature other examples (e.g. VOPO) where the only examination of the susceptibility led to a wrong interpretation of the behavior of the system.

## VI Conclusions

We have presented ab initio band-structure calculations for CsV2O5, a low-dimensional spin-gapped system and member of the alkali-vanadate family AV2O5 (A=Na, Li, Cs) which has a similar monoclinic structure as the most studied VOPO. By means of a downfolding-tight-binding analysis, we observe that the behavior of this system can be best described by a spin-$`\frac{1}{2}`$ alternating chain model with weak interchain interactions. This proposal has been complemented by an analysis of available susceptibility data for this compound. A re-examination of the temperature dependence of the experimental susceptibility leads us to the conclusion that an alternating chain model can also provide a good fitting by choosing appropriate $`g`$ and $`J`$ values. We have discussed the structural, electronic and magnetic similarities and differences between CsV2O5 and $`\alpha ^{}`$-NaV2O5, $`\gamma `$-LiV2O5 and VOPO.

With this work we would like to draw the attention on this compound which shows interesting features, especially since it has a lot of points of comparison to the most studied VOPO. We think that new experimental studies like ESR, INS, Raman scattering and the study of a possible field-induced magnetic ordering under a very strong magnetic field would be helpful to unveil the behavior of this material.


text:
# Chirally improved Dirac operators: Studying the sensitivity to topological excitations for zero and finite temperatureBased on the contributions of C. Gattringer, C.B. Lang and S. Schaefer at Lattice 2001.

## 1 Motivation and introduction

During the last three years it has been realized that the Ginsparg-Wilson relation  is the key to chirally symmetric fermions on the lattice. However, the overlap operator  which exactly solves the Ginsparg-Wilson equation is numerically very expensive. For many applications an alternative to an exact solution is an approximate Ginsparg-Wilson Dirac operator, such as a finite parametrization of the perfect action  or the recently proposed chirally improved Dirac operator . The latter results from a systematic expansion of a solution of the Ginsparg-Wilson equation.

In this contribution we report on our results for the eigenvectors and eigenvalues of the chirally improved Dirac operator, using smooth instanton backgrounds as well as thermalized (quenched) gauge configurations on both sides of the QCD phase transition. These studies serve two purposes: Firstly, they allow to obtain a sound understanding of the properties of the chirally improved Dirac operator. In particular we also compare our results to calculations done with the overlap operator in order to understand the effects of small violations of chirality. Secondly, the analysis of the eigenmodes provides insights into the nature of the QCD vacuum, since the eigenvectors couple to localized objects in the underlying gauge field such as instantons or instanton anti-instanton molecules.

## 2 Remarks on the construction of the lattice Dirac operators

Chirally improved operator: In  we suggested an expansion of the lattice Dirac operator $`D`$ in terms of a systematic series of simple hopping terms, connecting sites along paths of arbitrary length and allowing for all elements of the Clifford algebra. For the ansatz one respects all standard translational and rotational symmetries of the lattice Dirac operator as well as invariance under C and P and $`\gamma _5`$-hermiticity. We insert that formal series into

$$EDD^{}+D^{}D.$$
(1)
Finding a solution to the Ginsparg-Wilson equation  amounts to solving for $`E=0`$. The product $`D^{}D`$ is again constructed algebraically, leading to terms involving Clifford algebra elements and paths corresponding to ordered products of gauge link variables.

$`E`$ is an infinite series of terms that are linear and quadratic in the coupling coefficients of the Dirac operator multiplying the path terms. Independent terms have to vanish and thus we obtain a system of coupled equations quadratic in the (unknown) coupling constants. We truncate the ansatz for $`D`$ and the series for $`E`$ at paths of length 4. Whereas the complete series would lead to an exact solution of the Ginsparg-Wilson equation, the solution now will only be approximate and it remains to be studied, how good this approximation is.

In solving the set of equations there is much freedom. In order to partially remedy the approximation we introduce two more conditions (and parameters) such that the boundary condition for a massless Dirac operator

$$\widehat{D}(p)i\mathit{}+𝒪(p^2),$$
(2)
is satisfied at a given value of the gauge coupling. The solution we study here has 19 terms and values of the coefficients are listed in the appendix of . It is ultralocal with non-vanishing couplings only for $`|xy|\sqrt{5}`$ decreasing exponentially in size.

Exact solutions to the Ginsparg-Wilson equation ($`E=0`$) have eigenvalues on a unit circle with center 1 in the complex plane. It turned out that our solution, defining a chirally improved Dirac operator, has eigenvalues, that are very close to the unit circle in particular for the small eigenvalue region, which is most relevant for the infrared behavior.

Overlap operator: Below we compare the properties of the sector with small eigenvalues (complex or real) of the introduced chirally improved operator with those of an overlap operator  which has the form

$$D_{ov}=1Z\text{with}Z\gamma _5sign(H),$$
(3)
where $`H`$ is related to the hermitian Dirac operator. This is constructed from an arbitrary Dirac operator $`D_0`$ (e.g. the Wilson operator),

$$H=\gamma _5(sD_0),$$
(4)
and $`s`$ is a parameter which may be adjusted in order to minimize the probability for zero modes of $`H`$. In case $`D_0`$ is already an overlap operator one reproduces $`D_{ov}=D_0`$ for $`s=1`$ because $`sign\left(sign(H)\right)=sign(H)`$.

The sign-function may be defined through the spectral representation

$$signH=\underset{i}{}sign\lambda _i|ii|,$$
(5)
($`|i`$ denote the eigenvectors) but in practical computations this definition cannot be used for realistic QCD Dirac operators, which for $`L^4`$ lattices have dimension $`𝒪(10^510^6)`$ ($`n_{color}n_{Dirac}L^4`$). Note that in the subsequent applications the operator has to be applied many times, since it may be entering a diagonalization or a conjugate gradient inversion tool. One therefore relies on the relation

$$sign(H)=\frac{H}{\sqrt{H^2}}$$
(6)
and approximates the inverse square by some method. In our computations  we follow the methods discussed in . One approximates the inverse square root by a Chebychev polynomial, which has exponential convergence in $`[ϵ,1]`$, where $`ϵ`$ (and thus the order of the polynomial) depends on the ratio of smallest to largest eigenvalue of $`H^2`$. Clenshaw’s recurrence formula provides further numerical stability.

Depending on the input Dirac operator $`D_0`$ the rate of convergence may be quite unfavorable. Although one may adjust $`s`$ it turns out that particular simple operators like the Wilson Dirac operator eventually will give rise to several small eigenvalues of $`H^2`$. Technically, one then proceeds as follows: One determines the subspace of e.g. the lowest 20 eigenmodes of $`H^2`$ and computes the inverse square separately for the subspace (using the spectral representation) and the reduced operator (with polynomial approximation). For this approach it is important to determine the subspace with high accuracy (we request 12 digits).

In our study we need to diagonalize not only the hermitian matrix $`H`$ but, for the later analysis, also $`D_{ov}`$ and the chirally improved operator. All diagonalizations have been done with the Arnoldi method .

CPU-time comparison: The central element of the diagonalization (and conjugate gradient methods) is the multiplication of the Dirac operator with a vector. For the consideration of the computational effort we disregard the startup time for initialization of the operator. Comparing the CPU-time in units of the time used for the Wilson Dirac operator we find that the chirally improved operator needs 24 units and the overlap operator (with Wilson input) typically needs $`300\pm 100`$ units (always with 12 digits accuracy for the matrix norm $`H^21`$). These are the results for lattice size $`12^4`$, but they are not significantly different for $`8^4`$ or $`16^4`$.

As mentioned, whenever $`(1D_0)`$ is already an overlap operator, one finds $`D_{ov}=D_0`$. We are interested in the case when $`(1D_0)`$ is close to, but not quite an overlap operator. Such an operator might accelerate the convergence in the construction of the overlap operator, since the convergence of the polynomial approximation depends crucially on the smallest eigenvalue of $`H^2`$. It is therefore of high practical importance, whether $`H`$ will have eigenvalues close to zero. A Dirac operator $`D_0`$ used as input to the construction of $`H`$ may produce such small eigenvalues. Physically they are related to the appearance or disappearance of an instanton. Different Dirac operators will differ in their sensitivity to “identify” such structures in the gauge configuration.

When using the chirally improved operator – in comparison to using the Wilson Dirac operator – as starting point we need typically a factor 2 less matrix-vector multiplications ($`75\pm 25`$ terms of the Chebychev series for $`12^4`$ or $`16^4`$ lattices). One reason that this improvement factor is not higher may lie in the removal of the small eigenmode subspace, which is relatively more efficient for the “bad” Wilson operator used as input to the overlap construction. Another reason may be that the series approximation forces all eigenvalues to the unit circle; since the chirally improved operator has its large eigenvalues not as close to the circle as the small ones, one has to pay the price with longer polynomials.

## 3 Small eigenvalues and instantons

For an instanton in the continuum it is known that the Dirac operator has a zero mode. On the lattice only an exact solution of the Ginsparg-Wilson equation will display an eigenvalue which vanishes exactly. An approximate Ginsparg-Wilson Dirac operator will then have a small real eigenvalue. Typically the size of the real part will increase as one shrinks the radius of the underlying instanton. The rate at which the real eigenvalue moves away from the origin is a measure for the quality of the approximate solution of the Ginsparg-Wilson equation. A good approximation will keep the eigenvalue near zero also for relatively small instantons. In order to implement this test of the would-be zero modes we use instantons with several values of the radius $`R`$ discretized on the lattice following the method described in .

Fig. 1 shows the position $`x`$ of the real eigenvalue (the would-be zero mode) of different Dirac operators as a function of the radius $`R`$ of the underlying instanton. We show the results computed on $`16^4`$ lattices for the chirally improved operator (circles), the Wilson operator (squares), and the overlap operator (diamonds). The latter, since it is an exact solution of the Ginsparg-Wilson equation, has an exact zero mode. It is, however, interesting to see that the overlap operator can trace the instanton only down to a radius $`R=1.5`$. The Wilson operator has an eigenvalue which very quickly departs from zero as $`R`$ is decreased. The chirally improved operator keeps the real eigenvalue fairly close to 0 even for relatively small instantons and only for $`R2`$ the deviation becomes significant.

For practical applications of the chirally improved Dirac operator, the question is whether also in a gauge configuration with quantum fluctuations the real eigenmodes lie close to the origin. In order to test this property, we computed histograms for the probability distribution of the real eigenvalues. The gauge fields were generated in a quenched simulation, using the tadpole improved Lüscher-Weisz action . The histograms were computed on $`16^4`$ lattices at values of $`\beta _1=`$ 8.10, 8.30 and 8.45 using 200 configurations for each $`\beta _1`$. The corresponding values of the couplings $`\beta _2`$ and $`\beta _3`$ of the gauge action can be found in . There we obtained for the lattice spacing $`a`$ the values 0.127 fm, 0.107 fm, and 0.1 fm for $`\beta _1=`$ 8.10, 8.30, and 8.45, respectively.

From Fig. 2 it is obvious that the eigenvalues are concentrated near the origin. The tail of the distribution decays very quickly and extends almost exclusively towards values $`x>0`$. Both these features are highly welcome: The narrow width of the distribution lets one hope for only a very small additive renormalization of the quark mass. The fact that the distribution does not extend to values $`x<0`$ implies furthermore that exceptional configurations (which spoil the computation of the propagator) are suppressed effectively.

## 4 Localization properties of zero modes

We now analyze the eigenvectors with real eigenvalues, i.e. the modes which correspond to the zero modes in the continuum. For a continuum instanton the zero modes are known to be localized at the same position as the instanton. This property also holds for thermalized configurations on the lattice . The eigenmodes provide an efficient filter for the excitations of the QCD vacuum as seen by the Dirac operator.

To study the eigenvectors we define the scalar density $`p_0(x)`$. Let $`\psi (x)_c`$ be an eigenvector of the Dirac operator with $`x`$ denoting the lattice point while $`c`$ is the color index. A gauge invariant density $`p_\sigma (x)`$ is defined as

$$p_\sigma (x)=\underset{c}{}\psi (x)_c^{}\mathrm{\Gamma }_\sigma \psi (x)_c.$$
(7)
(Here $`\mathrm{\Gamma }_0=\mathrm{𝟏}`$, $`\mathrm{\Gamma }_5=\gamma _5`$; the sum over Dirac indices is implied.) In the last section we have seen that the continuum zero modes correspond to eigenvectors of the lattice Dirac operator with eigenvalue zero (overlap operator) or small real eigenvalue (chirally improved operator, Wilson operator). An interesting question is, to which extent the eigenvectors are different for different lattice Dirac operators.

In Fig. 3 we compare iso-surfaces of the scalar density for the zero-mode of different Dirac operators on the same gauge configuration ($`12^4`$ lattice at $`\beta _1=8.45`$). All three operators see the same localized instanton. Only the details vary slightly as can be seen by the different fluctuations around the core.

An interesting effect was observed for small instantons in  where we compared the localization of the zero mode for the smooth instantons discussed above for different lattice Dirac operators. It was found that for instantons with radius smaller than 3 in lattice units, the overlap operator has a zero mode which is less localized than in the case of a continuum instanton. Since the other two operators (chirally improved, Wilson) do not display such a blowing up of zero modes we attribute this effect to the rather big size of the overlap operator. As opposed to the other two operators the overlap operator is not ultra-local, i.e. $`|D_{ov}(x,y)|`$ decays exponentially with $`|xy|`$ but does not vanish exactly. Objects smaller than the typical size of the overlap operator thus can no longer be resolved. We estimate  that at typical lattice spacings of $`a`$ 0.1 fm this effect sets in at about 0.3 fm.

For further study we introduce the so-called inverse participation ratio, a convenient measure of the localization widely used in solid state physics. Due to the normalization of the eigenmodes we have $`_xp_0(x)=1`$. The inverse participation ratio is then defined as

$$I=V\underset{x}{}p_0(x)^2,$$
(8)
where $`V`$ is the volume of the lattice. An alternative observable for localization, based on the self correlation of the scalar density $`p_0(x)`$ was analyzed in .

In Fig. 4 we show a scatter plot of the inverse participation ratio for zero modes of the overlap operator $`I_{ov}`$ plotted as a function of the inverse participation ratio of the corresponding modes of the chirally improved operator $`I_{ch.imp.}`$. The data were computed on $`12^4`$ lattices using 200 configurations at $`\beta _1=8.10`$ (circles), and $`\beta _1=8.45`$ (filled triangles). The full curve with filled circles is the result for the smooth artificial instantons. The dashed line is a straight line with slope 1. If both operators would see the same localization for the zero modes, all data points would be located on this line. However, obviously this is not the case. Except for a single data point, all data are below this line. This means that the overlap operator systematically produces zero modes with smaller localization when compared to the chirally improved operator; this trend becomes stronger for more localized states.

## 5 Analysis of the near-zero modes

So far we have concentrated on analyzing the eigenvectors with real eigenvalues, i.e. the would-be zero modes. However, also the eigenvectors with small but complex eigenvalues, the so-called near-zero modes, play an important role in QCD phenomenology. In particular the eigenvalue density $`\rho (\lambda )`$ near at the origin is related to the chiral condensate through the Banks-Casher formula ,

$$\overline{\psi }\psi =\pi \underset{\lambda 0}{lim}\underset{V\mathrm{}}{lim}\rho (\lambda )V^1.$$
(9)
According to the instanton liquid picture of chiral symmetry breaking , the near-zero modes come from weakly interacting instantons and anti-instantons. The corresponding eigenvectors are expected to be still quite localized, since they correspond to interacting “lumps” which still resemble instantons or anti-instantons. These predictions may be tested in an ab-initio calculation on the lattice. Again we use the inverse participation ratio defined in Eq. (8) to study the localization of these modes. An interesting question is how the localization of eigenvectors varies as a function of the eigenvalue.

In Fig. 5 we show a histogram of the inverse participation ratio for eigenvectors as a function of the imaginary part of the corresponding eigenvalues. The data were computed on $`16^4`$ lattices at values of $`\beta _1=8.10,8.30`$ and 8.45 using 200 configurations for each $`\beta _1`$. It is obvious that the most localized states have eigenvalues near the origin. Such large values are expected for weakly interacting instantons and anti-instantons. As $`\text{Im }\lambda `$ is increased, the localization decreases. The corresponding modes do no longer feel distinguished topological lumps but become dominated by quantum fluctuations.

An interesting observable was suggested by Horvath et al. in . Their local chirality variable was subsequently studied by several other groups in . The underlying idea is to test whether the near-zero modes are locally chiral, i.e. left-handed for lattice points near the instanton peak of the density and right-handed near the anti-instanton peak. The construction starts with defining left- and right-handed densities $`p_\pm (x)`$,

$$p_\pm (x)=\frac{1}{2}\left(p_0(x)\pm p_5(x)\right),$$
(10)
i.e. the density (7) projected on positive and negative chirality. A mode dominated by instantons and anti-instantons is expected to have $`p_+(x)>0`$ and $`p_{}(x)=0`$ for lattice points $`x`$ near the instanton peak and vice versa for $`x`$ near an anti-instanton peak. Thus the ratio $`p_+(x)/p_{}(x)`$ should be large for $`x`$ close to an instanton peak and small for all $`x`$ close to an anti-instanton peak. In a final step Horvath et al. map the positive real values of this ratio onto the interval between $`1`$ and $`+1`$ using the transformation

$$X(x)=\frac{4}{\pi }\text{arctan}\left(\sqrt{\frac{p_+(x)}{p_{}(x)}}\right)\mathrm{\hspace{0.33em}1}.$$
(11)
In order to suppress the quantum fluctuations one evaluates $`X(x)`$ only for $`x`$ near instanton or anti-instanton peaks, i.e. for $`x`$ near the maxima of the scalar density $`p_0(x)`$. We computed our results using cutoffs of 1%, 6.25% and 12.5% for the number of lattice points $`x`$ (ordered according to $`p_0(x)`$) and averaged over $`X(x)`$. On the gauge field ensemble the local chirality variable is expected to show a double peak structure only for instanton dominance.

In Fig. 6 we show the local chirality of the near-zero modes as a function of the imaginary part of their eigenvalues. The data were computed from 200 configurations on $`16^4`$ lattices at $`\beta _1=8.45`$. When using a cut of 1% on the number of supporting lattice points (top plot) we find the best signal while for cuts of 6.25% (middle plot) and 12.5% (bottom plot) the signal becomes weaker due to quantum fluctuations. Furthermore it is obvious for all three plots that the strongest signal for local chirality is found for the eigenvectors with eigenvalues closest to the origin. For larger values of $`|\text{Im }\lambda |`$ the signal of local chirality for the corresponding eigenvectors weakens and the mode becomes dominated by quantum fluctuations.

## 6 The near-zero modes at finite temperature

Let us now discuss what happens to instantons at the QCD phase transition. At the critical temperature $`T_c`$ the theory enters the deconfined phase and chiral symmetry is restored. The spectrum develops a gap (see e.g.  for a plot of typical spectra) and the density of eigenvalues is zero near the origin. According to the Banks-Casher formula (9) this gap in the spectral density amounts to a vanishing chiral condensate.

The topological charge does not vanish abruptly at $`T_c`$ and topological objects are observed also for $`T>T_c`$. Thus the liquid of weakly interacting instantons responsible for chiral symmetry breaking at $`T<T_c`$ must undergo some change in its structure as the temperature is increased above $`T_c`$. It may be expected that above $`T_c`$ instead of interacting weakly, instantons and anti-instantons form tightly bound molecules . This prediction from instanton models can be tested on the lattice (compare  for a similar study using the staggered Dirac operator).

In Fig. 7 we show the inverse participation ratio as a function of the imaginary part of the corresponding eigenvalue, as in Fig. 5, but now for finite temperature. In particular we show the results for configurations where the Polyakov loop $`L`$ is in its real branch, which is the domain of the Polyakov loop value for the full, unquenched theory. For the equivalent plot with complex Polyakov loop see . Due to the spectral gap in the deconfined phase, there are no near-zero modes and the graph has support for $`|\text{Im }\lambda |>0.15`$. One finds that now the most localized states are near the edge of the spectrum. Localization decreases quickly as $`|\text{Im }\lambda |`$ is increased and the mode becomes a bulk mode dominated by quantum fluctuations. The inverse participation ratio, however, provides only information on the localization of the modes and does not test for other properties such as local chirality. Thus it makes sense to repeat last section’s study of Horvath et al.’s local chirality variable for $`T>T_c`$.

In Fig. 8 we compare the local chirality of the near-zero modes for an ensemble at $`T<T_c`$ (top plot, $`6\times 20^3`$ lattices, 400 configurations at $`\beta _1=8.10`$) with an ensemble at $`T>T_c`$ (bottom plot, $`6\times 20^3`$ lattices, 147 configurations with real branch Polyakov loop at $`\beta _1=8.45`$). The plots were computed using only the eigenvectors with $`|\text{Im }\lambda |`$ less than some cut corresponding to roughly the 4 eigenvalues closest to the origin, respectively the edge of the spectrum. Zero modes were left out. For $`T<T_c`$ we find a very pronounced double peak structure, i.e. a strong signal for local chirality. At $`T>T_c`$ the double peak structure is gone and we do not observe local chirality in the eigenmodes. It is obvious that the weakly interacting instantons and anti-instantons undergo a drastic change as $`T`$ increases above $`T_c`$. The resulting states are still localized (see Fig. 7) but the local chirality present for weakly interacting topological excitations has vanished, probably due to the formation of tightly bound molecules.

Acknowledgements: We thank the Leibniz Rechenzentrum in Munich for computer time on the Hitachi SR8000 and their operating team for support and training. C. Gattringer and C.B. Lang thank the DOE’s Institute of Nuclear Theory at the University of Washington for its hospitality and the DOE for partial support during the completion of this work. We have profited considerably from discussion with Pierre van Baal, Tom DeGrand, Stefan Dürr, Peter Hasenfratz, Ivan Hip and Ferenc Niedermayer.