text:
# 1 Hard probes at LHC with CMS detector

## 1 Hard probes at LHC with CMS detector

High-$`p_T`$ QCD phenomena are perspective tools to study the properties of the hot matter created in ultrarelativistic nuclear interactions . Hard probes such as quarkonia, hard jets and heavy quarks are not thermalized but carry information about the earliest stages of a system subject to numerous final-state effects in heavy ion collisions at SPS and RHIC. At the LHC, a new regime is reached where hard and semi-hard QCD multi-particle production can already dominate over underlying soft events.

The following signals of jet quenching due to medium-induced parton energy loss  can be identified as being observable in heavy ion collisions at LHC .
1) The suppression of high-$`p_T`$ jet pairs as compared to what is expected from independent nucleon-nucleon interactions pattern . It includes modification of jet impact parameter dependence  and azimuthal anisotropy of jet quenching in non-central collisions .
2) The $`p_T`$-imbalance in $`\gamma +`$jet  and $`Z(\mu ^+\mu ^{})+`$jet  channels.
3) The suppression of high-$`p_T`$ hadron and photon (leading particle) yields due to a modification of the jet fragmentation function .
4) The modification of the high mass dimuon spectra from semileptonic $`B`$ and $`D`$ meson decays and secondary charmonium production due to heavy quark energy loss .

All above signals can be studied in CMS Heavy Ion Programme . The Compact Muon Solenoid (CMS) is a general purpose detector designed primary to search for the Higgs boson in proton-proton collisions at LHC . Accordingly, the detector is optimized for accurate measurements of the characteristics of high-energy leptons, photons and hadronic jets in a large acceptance, providing unique capabilities for “hard probes” for both $`pp`$ and $`AA`$ collisions. The central element of CMS is the magnet, a $`13`$ m long and $`6`$ m diameter solenoid with an internal radius $`3`$ m, which will provide a strong $`4T`$ uniform magnetic field. The $`4\pi `$ detector consists of a $`6`$ m long and $`1.3`$ m radius central tracker, electromagnetic (ECAL) and hadronic (HCAL) calorimeters inside the magnet and muon stations outside. The tracker and muon chambers covers the pseudorapidity $`|\eta |<2.4`$, while the ECAL and HCAL calorimeters reach $`|\eta |=3`$. A pair of quartz-fiber very forward (HF) calorimeters, located $`\pm 11`$ m from the interaction point, cover the region $`3<|\eta |<5`$ and complement the energy measurement. The tracker is composed of pixel layers and silicon strip counters. The CMS muon stations consist of drift tubes in the barrel region and cathode strip chambers in the end-cap regions. The electromagnetic calorimeter is made of almost $`83000`$ scintillating $`\mathrm{PbWO}_4`$ crystals, and the hadronic calorimeter consists of scintillator inserted between copper absorber plates. The main characteristics of calorimeters are presented in Table 1.

Table 2 presents the event rates for various channels including hard jets in one month of Pb$``$Pb beams assuming two weeks of data taking and luminosity $`L=5\times 10^{26}`$cm-2s-1. The production cross sections in minimum bias $`AA`$ collisions were obtained from those in $`pp`$ interactions at the same energy, $`\sqrt{s}=5.5`$ TeV, using simple parameterization $`\sigma _{AA}^h=A^2\sigma _{pp}^h`$. The $`pp`$ cross sections were evaluated using the PYTHIA\_$`6.1`$  with the CTEQ5L parton distribution functions. The estimated for the CMS acceptance statistics will be large enough to carefully study the dijet rate as a function of impact parameter as well as the $`\phi `$ and $`\eta `$ distributions of jets. The statistics for the $`\gamma `$+jet channel are satisfactory for studying the $`E_T`$-imbalance of the process, but the large background from jet$`+`$jet$`(\pi ^0)`$ is still under investigation. The statistics for $`Z(\mu ^+\mu ^{})`$+jet channel are rather low, but the background is less than $`10\%`$ in this case. Note that jet reconstruction in heavy ion collisions at CMS with efficiency and purity close to $`100\%`$ is possible starting $`E_T^{jet}100`$ GeV with energy resolution worser by factor $`2`$ as compared to $`pp`$ case .

## 2 Jets versus impact parameter and azimuthal angle

We have analyzed the impact parameter dependence of jet rates in Pb$``$Pb collisions at LHC (see Ref.  for model details). Figure 2 shows dijet rates in different impact parameter bins ($`E_T^{jet}>100`$ GeV, $`|\eta ^{jet}|<2.5`$): $`(i)`$ without energy loss, $`(ii)`$ with collisional loss only, $`(iii)`$ with collisional and radiative loss from BDMS model . The bin width $`2`$ fm in Fig. 2 corresponds to estimated accuracy of impact parameter determination using correlation between $`b`$ and energy flow deposited in very forward CMS HF calorimeter . The maximum and mean values of $`dN^{dijet}/db`$ distribution get shifted towards the larger $`b`$ because jet quenching is stronger in central collisions than in peripheral one’s,. Since the coherent LPM radiation induces a strong dependence of the radiative energy loss of a jet on the angular cone size , the result for jets with non-zero cone size $`\theta _0`$ is expected to be between $`(iii)`$ ($`\theta _00`$) and $`(ii)`$ cases. Thus the observation of a dramatic change in the $`b`$-dependence of jet rates as compared to what is expected from independent nucleon-nucleon interactions pattern, would indicate the existence of medium-induced parton rescattering.

Other interesting feature of jet production in semi-central $`AA`$ collisions can be azimuthal anisotropy of high-$`p_T`$ hadrons  or jets  due to parton energy loss in azimuthally non-symmetric volume of quark-gluon plasma. Figure 2 shows the distribution of jets over azimuthal angle $`\phi `$ with collisional and radiative loss (a) and collisional loss only (b) for $`b=0`$, $`6`$ and $`10`$ fm. The same conditions as for Fig. 2 are fulfilled. The plot is normalized on the distributions of jets over $`\phi `$ in Pb$``$Pb collisions without energy loss. The azimuthal anisotropy gets stronger as going from central to semi-central collisions, but the absolute suppression factor reduces with increasing $`b`$. The methodical advantage of azimuthal jet observables is that one needs to reconstruct only azimuthal position of jet without measuring total jet energy. The summarized in papers  methods for determination of the event plane are applicable for studying anisotropic flow in current heavy ion experiments at SPS and RHIC, and might be also used at LHC . It was also suggested recently in  the method for measurement of jet azimuthal anisotropy without reconstruction of the event plane. This technique is based on the calculation of correlations between the azimuthal position of jet axis and the angles of particles (not incorporated in the jet) using as weights the particle momenta or energy deposition in the calorimetric sectors.

## 3 Heavy quarks: $`B\overline{B}\mu ^+\mu ^{}`$ and $`BJ/\psi `$ modes

The open heavy flavour production is important due to opportunity to study the massive colour charge in quark-gluon environment, while jets probe the behaviour of massless partons in QCD-matter. The charm and bottom production cross section are much larger at LHC as compared with RHIC. Thus the systematical studies of heavy flavour (especially $`B`$) physics can be performed in CMS. It was found in the high invariant mass dimuon range $`10`$ GeV/c$`{}_{}{}^{2}M70`$ GeV/c2 the dominant contribution comes from $`B\overline{B}`$ fragmentation, and its extraction from uncorrelated and correlated background is possible . In-medium gluon radiation and collisional energy loss of heavy quarks can result in suppression and modification of high mass dilepton  and secondary $`BJ/\psi `$  spectra. In addition, finite quark mass effects can lead to the suppression of medium-induced radiation of heavy quarks as compared with massless partons, enhancing e.g. the $`B/\pi `$ ratio .

The important predicted feature of heavy quark production at LHC is that the contribution of gluon splittings in initial- or final-state shower evolution to heavy flavour yield can be important , while at present accelerator energies the most of heavy quarks being produced due to direct hard scatterings (“pair creation”). These are so-called “flavour excitation” (one heavy quark is produced in vertex of hard process and another quark is created from initial state parton shower) and “gluon splitting” (both heavy quarks are produced from final state parton shower). It was shown in  that such “showering” $`b\overline{b}`$ pairs are an important contribution to the dilepton continuum at high masses and, in fact, are of the same order as direct $`q\overline{q}Q\overline{Q}X`$ and $`ggQ\overline{Q}X`$ pair creation at the LHC for dimuons with masses greater than $`10`$ GeV/c2. Moreover, most of secondary $`J/\psi `$’s produced in $`B`$ decays originate from the $`ggg(g^{}Q\overline{Q})`$ diagram .

Figure 3 shows the dimuon invariant mass spectra at CMS acceptance for several models: without loss; with rescattering and collisional loss; with rescattering, collisional and radiative loss calculated as incoherent limit of Zakharov’s formula , and also from “dead-cone” model . The suppression factor $`2`$ -$`4`$ can be clearly observed over “nuclear shadowing” effect , $`15\%`$. Estimated suppression by a factor of $`1.5`$-$`2`$ for secondary $`J/\psi `$ is less than for $`B\overline{B}`$ decays and is comparable with the shadowing corrections ($`30\%`$).

## 4 Conclusion

To summarize, with its large acceptance, almost hermetic hadronic and electromagnetic calorimetry of fine granularity and good energy and spatial jet resolution, and excellent muon capabilities, the Compact Muon Solenoid (CMS) detector is a promising device to study jet quenching at LHC in various channels. Jet reconstruction in heavy ion collisions at CMS is possible starting $`E_T^{jet}=100`$ GeV, and estimated statistics will be large enough to carefully study the jet rate as a function of impact parameter, as well as azimuthal and rapidity distributions of jets. Impact parameter and event plane determination are possible using energy flow measurements. In-medium energy loss of heavy quarks also can be observable in $`B\overline{B}\mu ^+\mu ^{}`$ and $`BJ/\psi \mu ^+\mu ^{}`$ modes. However, in order to determine the baseline rates precisely, measurements in $`pp`$ or $`dd`$ collisions at the same or similar energies per nucleon as in the heavy ion runs are strongly desirable.

Acknowledgements. I would like to thank the organizers of the Symposium for the warm welcome and stimulating atmosphere. Discussions with M. Bedjidian, D. Denegri, Yu.L. Dokhitzer, D. Kharzeev, O.L. Kodolova, A.N. Nikitenko, S.V. Petrushanko, L.I. Sarycheva, S.V. Shmatov, A.M. Snigirev, V.V. Uzhinskii, I.N. Vardanian, S.A. Voloshin, R. Vogt, U. Wiedemann, P.I. Zarubin and G.M. Zinoviev are gratefully acknowledged.


format: text

text:
# 1 Introduction

## 1 Introduction

The problem of short-range nucleon-nucleon correlations (SRC) in nuclei is still an open field in nuclear physics . An indirect evidence for SRC was obtained by the observation of a depletion of valence shells in $`(e,e^{}p)`$ reactions . It is also well known that SRC modify the nucleon momentum distribution , i.e. the Fourier transform of the one-body density matrix $`\rho _1(\stackrel{}{r}_1,\stackrel{}{r}_1^{})`$ in the variable $`\stackrel{}{r}_1\stackrel{}{r}_1^{}`$, by introducing significant contributions at values of momentum beyond the Fermi momentum. Experimental information from inclusive $`(e,e^{})`$ and exclusive $`(e,e^{}p)`$ reactions established the existence of a high-momentum component in the momentum distribution.

Beyond the one-body density matrix and momentum distribution, rich information on the nuclear ground state and nucleon-nucleon correlations is contained in the two-body density matrix (2DM) $`\rho _2(\stackrel{}{r}_1,\stackrel{}{r}_2;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})`$ and its various Fourier transforms in momentum space. Experimentally, exclusive two-nucleon knockout reactions can provide information on the relative behaviour of nucleon pairs embedded in the nuclear medium and give more direct insight into SRC . A very important role is played by electron-nucleus scattering experiments with the availability of high-energy 100% duty cycle electron beams which made possible to carry out double coincidence $`(e,e^{}NN)`$ experiments. Past, present and near future experiments at NIKHEF, MAMI and Jefferson Lab use 3He, 12C and 16O as targets .

The analysis of $`(e,e^{}2N)`$ experiments is a complicated task . Among the most important processes of longitudinal character contributing to electromagnetically induced two-nucleon knockout is the scattering off a correlated nucleon pair followed by two-nucleon knockout during which each of the two nucleons interacts with the residual nucleus (final-state interactions - FSI). The minimal starting point in the analysis is the plane-wave impulse approximation (PWIA), in which the FSI are neglected and the two-nucleon knock-out cross section can be expressed in terms of the two-hole spectral function $`S(\stackrel{}{p}_1,\stackrel{}{p}_2;E)`$. The latter gives the joint probability of removing from the target two nucleons with momenta $`\stackrel{}{p}_1`$ and $`\stackrel{}{p}_2`$ leaving the residual nucleus with energy $`E`$ with respect to the target nucleus ground state. Under certain conditions it seems possible to extract information on the 2DM and its Fourier transforms and on nucleon-nucleon correlations. Taking the above into account, the development of a simple model of the ground-state 2DM in order to clarify the two-body characteristics of the nuclear ground state and the extent to which it is affected by the SRC appears useful.

Here we focus on the two-body momentum distribution (2bMD) $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$, which is the Fourier transform of the 2DM $`\rho _2(\stackrel{}{r}_1,\stackrel{}{r}_2;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})`$ in the variables $`\stackrel{}{r}_1\stackrel{}{r}_1^{}`$ and $`\stackrel{}{r}_2\stackrel{}{r}_2^{}`$. Another Fourier transform of the half-diagonal two-body density matrix $`\rho _{2h}(\stackrel{}{r}_1,\stackrel{}{r}_2,\stackrel{}{r}_1^{})=\rho _2(\stackrel{}{r}_1,\stackrel{}{r}_2;\stackrel{}{r}_1^{},\stackrel{}{r}_2)`$ , namely the generalized momentum distribution $`\eta (\stackrel{}{p},\stackrel{}{Q})`$, has been studied previously in the case of infinite nuclear matter  and finite nuclei . The generalized momentum distribution plays an important role in the study of the FSI. The 2bMD is straighforwardly connected to the two-nucleon spectral function $`S(\stackrel{}{p}_1,\stackrel{}{p}_2;E)`$ via an integration with respect to the energy (see Eq.(10) below). The two-nucleon spectral function in the case of infinite nuclear matter has been studied using correlated basis function perturbation theory including central and tensor correlations . In the case of 16O it has been evaluated including short-range correlations through defect functions calculated via the Bethe-Goldstone equation or correlation functions derived from variational calculations and long-range correlations using dressed RPA .

In this work, first an analytical expression is derived for the 2bMD of $`Z=N`$, $`\mathrm{}`$-closed nuclei in the independent-particle shell model with harmonic oscillator wave-functions and applied to the nucleus 16O. Our method is an extension of that developed in Refs.  for the study of the nuclear form factor, the nuclear charge, matter and momentum distributions, the one-body density matrix and the generalized momentum distribution in closed shell nuclei. This calculation is expected to reproduce the main features of the 2bMD, stemming from the finite nuclear size and Fermi statistics, being reliable at certain kinematical domains. By introducing the quantity $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)/\eta \left(\stackrel{}{p}_1\right)\eta \left(\stackrel{}{p}_2\right)`$, where $`\eta \left(\stackrel{}{p}\right)`$ is the nucleon momentum distribution, we pinpoint the effects of the finite size and of the statistical correlations.

Next, the effect of the SRC on the 2bMD of 4He is studied by including Jastrow-type correlations in our calculation, via the lowest term of a cluster expansion. This so-called low order approximation (LOA)  has been exploited for the one-body density matrix and the two-body density matrix by Bohigas and Stringari  and Dal Rì, Stringari and Bohigas  and has been widely used to study single-particle nuclear properties, such as the (point-nucleon or charge) density , the elastic or charge form factor , the one-body density matrix in coordinate  or momentum space  and the nucleon momentum distribution , as well as two-body nuclear quantities, namely the two-body density matrix in coordinate space , the static structure function $`S\left(Q\right)`$ , the relative pair density distribution $`\rho _2\left(r_{12}\right)`$ , the center-of-mass pair density distribution $`\rho _2\left(R_{CM}\right)`$ , the two-body relative momentum distribution $`\eta _2\left(\stackrel{}{q}_{\mathrm{r}el}\right)`$ , the two-body center-of-mass momentum distribution $`\eta _2\left(\stackrel{}{P}_{\mathrm{C}M}\right)`$  and the combined two-body center-of-mass and relative momentum distribution $`\eta _2(P_{\mathrm{C}M},q_{\mathrm{r}el})`$ . A two-gaussian central correlation function, $`f\left(r\right)=1c_1\mathrm{exp}\left(r^2/\beta _1^2\right)+c_2\mathrm{exp}\left(r^2/\beta _2^2\right)`$, is used in this work. We find significant deviations from the independent-particle picture at high momenta, especially when the vectors $`\stackrel{}{p}_1`$ and $`\stackrel{}{p}_2`$ are antiparallel or parallel, as it is clearly shown by the values of the quantity $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)/\eta \left(\stackrel{}{p}_1\right)\eta \left(\stackrel{}{p}_2\right)`$. There is, however, within the LOA, significant sensitivity of the results in some kinematical regions, mainly on the correlation function used to describe the SRC. For certain correlation functions one may obtain negative values of $`\eta _2`$ in few regions of momenta.

The definition and properties of the two-body momentum distribution are presented in detail in Section 2. In Section 3 an exact expression is derived for the 2bMD of $`Z=N`$, $`\mathrm{}`$-closed nuclei in the independent-particle shell model with harmonic oscillator wavefunctions and applied to the case of 16O. In Section 4 the effect of SRC on the 2bMD of 4He is explored by including Jastrow-type correlations in our calculation via the LOA, and the results as well as the reliability of our approximation are discussed. Finally, in Section 5 we give a summary of the results and hints for possible future developments.

## 2 The Two-Body Momentum Distribution

The two-body momentum distribution $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$ gives the combined probability density of finding in the nucleus a nucleon with momentum $`\stackrel{}{p}_1`$ and another one with momentum $`\stackrel{}{p}_2`$. In a system of $`A`$ $`\left(A2\right)`$ identical particles, described by a unit-normalized state $`|\mathrm{\Psi }`$, it can be defined as the expectation value

$$\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)=\mathrm{\Psi }\left|\underset{\stackrel{}{s},\stackrel{}{s^{}}}{}a_{\stackrel{}{p}_2,\stackrel{}{s^{}}}^{}a_{\stackrel{}{p}_1,\stackrel{}{s}}^{}a_{\stackrel{}{p}_1,\stackrel{}{s}}a_{\stackrel{}{p}_2,\stackrel{}{s^{}}}\right|\mathrm{\Psi }.$$
(1)
The 2bMD is the diagonal element of the 2DM in momentum space $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2;\stackrel{}{p}_1^{},\stackrel{}{p}_2^{})`$,

$$\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)=\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2;\stackrel{}{p}_1,\stackrel{}{p}_2)$$
(2)
and may also be defined as the Fourier transform of the 2DM $`\rho _2(\stackrel{}{r}_1,\stackrel{}{r}_2;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})`$ in the variables $`\stackrel{}{r}_1\stackrel{}{r}_1^{}`$ and $`\stackrel{}{r}_2\stackrel{}{r}_2^{}`$:

$$\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)=\rho _2(\stackrel{}{r}_1,\stackrel{}{r}_2;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})\mathrm{e}^{i\stackrel{}{p}_1\left(\stackrel{}{r}_1\stackrel{}{r}_1^{}\right)}\mathrm{e}^{i\stackrel{}{p}_2\left(\stackrel{}{r}_2\stackrel{}{r}_2^{}\right)}\mathrm{d}^3r_1\mathrm{d}^3r_2\mathrm{d}^3r_1^{}\mathrm{d}^3r_2^{}.$$
(3)
(By taking $`\mathrm{}=1`$, the 2bMD has the dimension (length)6.) The normalization adopted is such that

$$\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)\mathrm{d}^3p_1\mathrm{d}^3p_2=A\left(A1\right),$$
(4)
and the following sequential relation holds:

$$\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)\mathrm{d}^3p_2=\left(A1\right)\eta \left(\stackrel{}{p}_1\right),$$
(5)
where $`\eta \left(\stackrel{}{p}\right)`$ is the one-body momentum distribution, normalized to the number of particles

$$\eta \left(\stackrel{}{p}\right)\mathrm{d}^3p=A.$$
(6)
Equation (5) connects the 2bMD to the generalized momentum distribution  ($`\eta \left(\stackrel{}{p},\stackrel{}{Q}=0\right)=\left(A1\right)\eta \left(\stackrel{}{p}\right)`$).

The 2bMD $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$ is related to the two-body center-of-mass momentum distribution $`\eta _2\left(\stackrel{}{P}_{\mathrm{C}M}\right)`$ ($`\stackrel{}{P}_{\mathrm{C}M}=\stackrel{}{p}_1+\stackrel{}{p}_2`$ is the two-nucleon center-of-mass momentum), and to the two-body relative momentum distribution $`\eta _2\left(\stackrel{}{q}_{\mathrm{r}el}\right)`$ ($`\stackrel{}{q}_{\mathrm{r}el}=\frac{1}{2}\left(\stackrel{}{p}_1\stackrel{}{p}_2\right)`$ is the relative two-nucleon momentum)

$$\eta _2\left(\stackrel{}{P}_{\mathrm{C}M}\right)=\eta _2(\stackrel{}{P}_{\mathrm{C}M}/2+\stackrel{}{q}_{\mathrm{r}el},\stackrel{}{P}_{\mathrm{C}M}/2\stackrel{}{q}_{\mathrm{r}el})\mathrm{d}^3q_{\mathrm{r}el},$$
(7)
$$\eta _2\left(\stackrel{}{q}_{\mathrm{r}el}\right)=\eta _2(\stackrel{}{P}_{\mathrm{C}M}/2+\stackrel{}{q}_{\mathrm{r}el},\stackrel{}{P}_{\mathrm{C}M}/2\stackrel{}{q}_{\mathrm{r}el})\mathrm{d}^3P_{\mathrm{C}M}.$$
(8)
It is also related to the combined two-body center-of-mass and relative momentum distribution $`\eta _2(P_{\mathrm{C}M},q_{\mathrm{r}el})`$

$$\eta _2(P_{\mathrm{C}M},q_{\mathrm{r}el})=\eta _2(\stackrel{}{P}_{\mathrm{C}M}/2+\stackrel{}{q}_{\mathrm{r}el},\stackrel{}{P}_{\mathrm{C}M}/2\stackrel{}{q}_{\mathrm{r}el})d\mathrm{\Omega }_{P_{\mathrm{C}M}}d\mathrm{\Omega }_{q_{\mathrm{r}el}}.$$
(9)
Its relationship with the two-nucleon spectral function $`S(\stackrel{}{p}_1,\stackrel{}{p}_2;E)`$ is

$$\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)=dES(\stackrel{}{p}_1,\stackrel{}{p}_2;E).$$
(10)
In a totally uncorrelated system the 2bMD is factorized as

$$\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)=\left(1\frac{1}{A}\right)\eta \left(\stackrel{}{p}_1\right)\eta \left(\stackrel{}{p}_2\right).$$
(11)
The constant $`11/A`$ appears due to the normalization condition (4). For the ground state of a non-interacting Fermi system (a system of $`A`$ identical, non-interacting fermions) where only statistical correlations are present, starting from the 2DM , the 2bMD can be shown to have the form

$$\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)=\eta \left(\stackrel{}{p}_1\right)\eta \left(\stackrel{}{p}_2\right)\frac{1}{\nu }\left|\eta _1(\stackrel{}{p}_1,\stackrel{}{p}_2)\right|^2,$$
(12)
where $`\eta _1(\stackrel{}{p}_1,\stackrel{}{p}_2)`$ is the one-body density matrix in momentum space and $`\nu `$ is the level degeneracy. Note that in the special case of a non-interacting Fermi system with $`A=\nu `$ Eq. (12) becomes identical with Eq. (11).

As an example we consider the ground state of an infinitely extended ideal Fermi gas (volume $`\mathrm{\Omega }`$ and number of particles $`A`$ going to infinity while the particle density stays constant) with Fermi wave number $`k_\mathrm{F}`$ and degeneracy $`\nu `$. Equation (12) yields for the 2bMD per particle-pair $`\stackrel{~}{\eta }_2^\mathrm{F}(\stackrel{}{p}_1,\stackrel{}{p}_2)=\eta _2^\mathrm{F}(\stackrel{}{p}_1,\stackrel{}{p}_2)/\left[A\left(A1\right)\right]`$

$$\stackrel{~}{\eta }_2^\mathrm{F}(\stackrel{}{p}_1,\stackrel{}{p}_2)=\left(\frac{4}{3}\pi k_\mathrm{F}^3\right)^2\theta \left(k_\mathrm{F}p_1\right)\theta \left(k_\mathrm{F}p_2\right)\left[1\delta _{\stackrel{}{p}_1\stackrel{}{p}_2}/\nu \right].$$
(13)
The momentum distribution per particle, $`\stackrel{~}{\eta }^\mathrm{F}\left(p\right)=\eta ^\mathrm{F}\left(p\right)/A=\eta ^\mathrm{F}\left(\stackrel{}{p}\right)/A`$, is given by

$$\stackrel{~}{\eta }^\mathrm{F}\left(p\right)=\left(\frac{4}{3}\pi k_\mathrm{F}^3\right)^1\theta \left(k_\mathrm{F}p\right).$$
(14)
Let us define the quantity $`\xi `$ as a measure of finite-size effects, of statistical and short-range correlations of dynamical origin

$$\xi (\stackrel{}{p}_1,\stackrel{}{p}_2)\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)/\eta \left(\stackrel{}{p}_1\right)\eta \left(\stackrel{}{p}_2\right).$$
(15)
For an infinitely extended ideal Fermi gas $`\xi `$ is defined for $`p_1,p_2k_F`$ and equals

$$\xi ^\mathrm{F}(\stackrel{}{p}_1,\stackrel{}{p}_2)=\{\begin{array}{cc}1,\hfill & \stackrel{}{p}_1\stackrel{}{p}_2\hfill \\ 11/\nu ,\hfill & \stackrel{}{p}_1=\stackrel{}{p}_2\hfill \end{array}.$$
(16)
From Eqs. (11) and (16) we realize that, in the case of the infinitely extended ideal Fermi system, if $`\stackrel{}{p}_1\stackrel{}{p}_2`$, $`\xi `$ equals 1. This holds even in the presence of dynamical SRC, if long-range order does not exist . In the case of the infinitely extended ideal Fermi gas, $`\xi =11/\nu `$ if $`\stackrel{}{p}_1=\stackrel{}{p}_2`$. As for the finite, non-interacting Fermi system, we realize from Eq. (12) that, if $`\stackrel{}{p}_1=\stackrel{}{p}_2`$, $`\xi =11/\nu `$. Deviations of $`\xi `$ from this value show the effect of dynamical correlations. For $`\stackrel{}{p}_1\stackrel{}{p}_2`$ deviations of $`\xi `$ from $`11/A`$ is a measure of statistical and (or) dynamical correlations in a system of finite size.

## 3 Method of Calculation for $`\mathrm{}`$-Closed Nuclei \- Application to 16O

### 3.1 Method of Calculation

Let us now consider a system of $`A`$ identical non-interacting fermions in its ground state. The fermions occupy the lowest single-particle energy eigenstates $`|n_j`$ ($`j=1,2,\mathrm{},A/\nu `$) and the 2bMD is given by Eq. (12). In the case of the nucleus, $`A`$ is the mass number and $`\nu `$ is the degeneracy due to the nucleon spin and isospin.We consider $`Z=N`$, $`\mathrm{}`$-closed nuclei. In order to obtain closed analytical expressions for $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$, we have assumed that the nucleons move in an isotropic harmonic oscillator potential and that the spin-orbit coupling is negligible, the center-of-mass and finite nucleon size corrections are small and the Coulomb interaction (relevant for protons) is small. Then Eq. (12) gives:

$$\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)=\frac{b^6}{\pi ^3}\mathrm{e}^{p_1^2b^2}\mathrm{e}^{p_2^2b^2}\underset{\mu _1=0}{\overset{2N_{\mathrm{max}}}{}}\left(p_1b\right)^{\mu _1}\underset{\mu _2=0}{\overset{2N_{\mathrm{max}}}{}}\left(p_2b\right)^{\mu _2}𝒢_{\mu _1\mu _2}\left(\mathrm{cos}\gamma \right),$$
(17)
where the coefficients $`𝒢_{\mu _1\mu _2}`$ are given by

$$𝒢_{\mu _1\mu _2}\left(\mathrm{cos}\gamma \right)=4f_{\mu _1/2}f_{\mu _2/2}\underset{\lambda _1=0}{\overset{\mu _1}{}}\underset{\lambda _2=0}{\overset{\mu _2}{}}𝒦_{\lambda _1\lambda _2}\left(\mathrm{cos}\gamma \right)𝒦_{\left(\mu _1\lambda _1\right)\left(\mu _2\lambda _2\right)}\left(\mathrm{cos}\gamma \right)$$
(18)
and the coefficients $`f_{\mu /2}`$ and $`𝒦_{\lambda _1\lambda _2}\left(\mathrm{cos}\gamma \right)`$ are given by

$$f_{\mu /2}=\underset{n\mathrm{},\text{occ.}}{}f_n\mathrm{}^{\mu /2},$$
(19)
$$𝒦_{\lambda _1\lambda _2}\left(\mathrm{cos}\gamma \right)=\underset{n\mathrm{},\text{occ.}}{}K_n\mathrm{}^{\lambda _1\lambda _2}P_{\mathrm{}}\left(\mathrm{cos}\gamma \right).$$
(20)
In the above equations $`\gamma `$ is the angle between the vectors $`\stackrel{}{p}_1`$ and $`\stackrel{}{p}_2`$, b is the harmonic oscillator parameter, $`P_{\mathrm{}}\left(x\right)`$ is a Legendre polynomial and $`N_{\mathrm{max}}=\left(2n+\mathrm{}\right)_{\mathrm{max}}`$ is the number of energy quanta of the highest occupied $`n\mathrm{}`$-level. The coefficients $`f_n\mathrm{}^{\mu /2}`$, $`K_n\mathrm{}^{\lambda _1\lambda _2}`$ are rational numbers that enter the corresponding expressions of $`\eta \left(\stackrel{}{p}\right)`$ and $`\eta _1(\stackrel{}{p}_1,\stackrel{}{p}_1^{})`$ respectively. Analytical expressions for $`f_n\mathrm{}^{\mu /2}`$ and $`K_n\mathrm{}^{\lambda _1\lambda _2}`$ can be found in Ref. . They are different from zero only if the indices $`\mu `$ and $`\lambda _1+\lambda _2`$ are even.

The corresponding expression for the spherically symmetric nucleon momentum distribution $`\eta \left(\stackrel{}{p}\right)`$ is

$$\eta \left(\stackrel{}{p}\right)=\eta \left(p\right)=\frac{b^3}{\pi ^{3/2}}\mathrm{e}^{p^2b^2}\underset{\lambda =0}{\overset{N_{\mathrm{max}}}{}}\left(pb\right)^{2\lambda }2f_\lambda .$$
(21)
It is verified that the 2bMD $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$ given by Eq. (17) satisfies the property of Eq. (5).

### 3.2 Application to 16O

Using the expressions (17) and (21), we obtain results for the quantities $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$ and $`\xi (\stackrel{}{p}_1,\stackrel{}{p}_2)`$, in the case of the nucleus 16O. The harmonic oscillator parameter $`b=1.7825`$fm has been determined so as to reproduce the experimental value of the charge root mean square radius of 16O, $`r^{1/2}=2.737`$fm . In Fig. 1 the 2bMD of 16O in the harmonic oscillator model (continuous lines) is plotted for $`\stackrel{}{p}_1`$ parallel to $`\stackrel{}{p}_2`$ and for $`p_1=`$ 0, 1 and 1.5 fm-1, as a function of the variable p2 ($`\stackrel{}{p}_2=\mathrm{p}_2\widehat{p}_1`$). In the same graphs the 2bMD of an infinite ideal Fermi gas is shown (dashed lines), calculated via Eq. (13) with $`\nu =4`$ and with Fermi wave number $`k_\mathrm{F}`$ equal to the one corresponding to 16O, namely $`k_\mathrm{F}=1.1`$fm-1. This latter value has been calculated using local density approximation and coincides with the one used in Ref. . In order to make the comparison meaningful, the 2bMD of the ideal Fermi gas was normalized to the same number of particle pairs as those found in 16O, namely to $`A\left(A1\right)=240`$.

It is evident from Eq. (13) that the 2bMD of the ideal Fermi gas exhibits discontinuities at $`p_{1,2}=k_\mathrm{F}`$ and at $`\stackrel{}{p}_1=\stackrel{}{p}_2`$. This is observed in the first two panels of Fig. 1 ($`p_1=`$ 0, 1 fm-1) at p$`{}_{2}{}^{}=p_1`$ and at p$`{}_{2}{}^{}=\pm k_\mathrm{F}`$. In the case of the finite system 16O in the harmonic oscillator model a rather similar behaviour is observed, but the 2bMD changes smoothly with p2. For $`p_1=1.5`$fm-1 (which is higher than $`k_\mathrm{F}=1.1`$fm-1) the 2bMD of the ideal Fermi gas has vanished and that of 16O has dropped considerably.

In Fig. 2 the quantity $`\xi (\stackrel{}{p}_1,\stackrel{}{p}_2)`$ is plotted in the case of 16O as a function of $`\mathrm{cos}\gamma `$ and for $`p_1=p_2=`$ 0, 1 and 4 fm-1. The values of $`\xi `$ lie between 0.75($`=11/\nu `$) and 1. $`\xi `$ is equal to 0.75 when $`\stackrel{}{p}_1=\stackrel{}{p}_2`$ (as shown for example in the first panel, and in the other two panels when $`\mathrm{cos}\gamma =1`$). We realize that for high values of $`p_1`$ and $`p_2`$, $`\xi `$ tends to 1 and 0.75 for $`\mathrm{cos}\gamma =0`$ and $`\mathrm{cos}\gamma =1`$ respectively. The deviation of the values of $`\xi `$ from 1 or 0.75 is a finite-size effect in the presence of statistical correlations. In the case when $`p_1p_2`$ (not shown) the quantity $`\xi `$ as a function of $`\mathrm{cos}\gamma `$ remains higher than 0.75 and below or equal to 1, following a similar behaviour as in Fig. 2, namely $`\xi \left(1\right)>\xi \left(1\right)`$ and $`\xi \left(x_0\right)=1`$ for some value $`\mathrm{cos}\gamma =x_0<0`$.

## 4 The Effect of Short-Range Correlations

Up to now we have ignored the effect of dynamical, short-range correlations on the 2bMD. It is expected that such a calculation is reliable only in restricted kinematical domains, namely for $`p_1`$, $`p_2k_\mathrm{F}`$. In this Section we make an attempt to investigate the effect of the SRC by including Jastrow-type correlations in our study.

### 4.1 Method of Calculation

The approach adopted here is based on the Jastrow formalism with a state-independent central correlation function, which introduces the short range correlations. Employing the low-order approximation (LOA) of Ref.  for the 2DM and performing the spin-isospin summation ($`\nu =4`$) in Eq. (14) of Ref. , we obtain for the correlated 2DM $`\rho _2^{\mathrm{corr}}(\stackrel{}{r}_1,\stackrel{}{r}_2;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})`$ the following expression:

$`\rho _2^{\mathrm{corr}}`$ $`(\stackrel{}{r}_1,\stackrel{}{r}_2;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})=\left[1+g(r_{12},r_{1^{}2^{}})\right]\rho _2(\stackrel{}{r}_1,\stackrel{}{r}_2;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})`$ (22)
$`+{\displaystyle }[g(r_{13},r_{1^{}3})+g(r_{23},r_{2^{}3})][\rho _1(\stackrel{}{r}_1,\stackrel{}{r}_1^{})\rho _{2h}(\stackrel{}{r}_2,\stackrel{}{r}_3,\stackrel{}{r}_2^{})`$
$`\nu ^1\rho _1(\stackrel{}{r}_1,\stackrel{}{r}_2^{})\rho _{2h}(\stackrel{}{r}_2,\stackrel{}{r}_3,\stackrel{}{r}_1^{})\nu ^1\rho _1(\stackrel{}{r}_1,\stackrel{}{r}_3)\rho _2(\stackrel{}{r}_2,\stackrel{}{r}_3;\stackrel{}{r}_2^{},\stackrel{}{r}_1^{})]\mathrm{d}^3r_3`$
$`\nu ^1{\displaystyle }{\displaystyle }g(r_{34},r_{34})\{\rho _{2h}(\stackrel{}{r}_2,\stackrel{}{r}_4,\stackrel{}{r}_3)\rho _2(\stackrel{}{r}_1,\stackrel{}{r}_3;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})`$
$`+\rho _1(\stackrel{}{r}_1,\stackrel{}{r}_3)[\rho _1(\stackrel{}{r}_2,\stackrel{}{r}_2^{})\rho _{2h}(\stackrel{}{r}_3,\stackrel{}{r}_4,\stackrel{}{r}_1^{})\nu ^1\rho _1(\stackrel{}{r}_2,\stackrel{}{r}_1^{})\rho _{2h}(\stackrel{}{r}_3,\stackrel{}{r}_4,\stackrel{}{r}_2^{})`$
$`\nu ^1\rho _1(\stackrel{}{r}_2,\stackrel{}{r}_4)\rho _2(\stackrel{}{r}_3,\stackrel{}{r}_4;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})]\}\mathrm{d}^3r_3\mathrm{d}^3r_4,`$
where $`r_{ij}=\left|\stackrel{}{r}_i\stackrel{}{r}_j\right|`$, $`g(r,r^{})f\left(r\right)f\left(r^{}\right)1`$ and $`f\left(r\right)`$ is the correlation function, which has to obey the conditions $`f\left(0\right)<1`$ and for $`r\mathrm{}`$, $`f\left(r\right)1`$.

In Eq. (22) the one-body density matrix $`\rho _1(\stackrel{}{r},\stackrel{}{r^{}})`$, the 2DM $`\rho _2(\stackrel{}{r}_1,\stackrel{}{r}_2;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})`$ and the half-diagonal 2DM $`\rho _{2h}(\stackrel{}{r}_1,\stackrel{}{r}_2,\stackrel{}{r}_1^{})`$ are calculated in the harmonic-oscillator model. The correlated 2bMD, $`\eta _2^{\mathrm{corr}}(\stackrel{}{p}_1,\stackrel{}{p}_2)`$, is then calculated by Fourier transforming $`\rho _2^{\mathrm{corr}}(\stackrel{}{r}_1,\stackrel{}{r}_2;\stackrel{}{r}_1^{},\stackrel{}{r}_2^{})`$ according to Eq. (3). The correlated momentum distribution $`\eta ^{\mathrm{corr}}\left(\stackrel{}{p}\right)`$ is calculated likewise in the LOA by performing the spin-isospin summation in expression (13) of Ref.  for the one-body density matrix and by Fourier-transforming with respect to $`\stackrel{}{r}_1\stackrel{}{r}_1^{}`$.

The LOA preserves the normalization of the density matrices. As a consequence, the one- and two-body momentum distributions calculated in the LOA obey the normalization conditions (6) and (4), as well as the sequential relation (5). However, in the LOA some probability distributions (i.e. positive-definite quantities) may obtain negative values within regions of their domain. This problem is probably related to the $`A`$representability of the LOA two-body density matrices. For example, it was encountered in calculations of the relative distribution function $`\rho _2\left(r_{12}\right)`$  and the combined two-body center-of-mass and relative momentum distribution $`\eta _2(P_{\mathrm{C}M},q_{\mathrm{r}el})`$ , which implies that it may show up in the 2bMD too, since the latter is related to $`\eta _2(P_{\mathrm{C}M},q_{\mathrm{r}el})`$ via Eq. (9). In the case of the $`\rho _2\left(r_{12}\right)`$ the problem was subsequently solved by an appropriate modification of the correlation function . The above remarks are highly relevant to the choice of the correlation function used in this work, as we will discuss in due course.

### 4.2 Application to 4He

The 2bMD and the momentum distribution of 4He in the harmonic oscillator model, which we will hereafter denote by $`\eta _2^0(\stackrel{}{p}_1,\stackrel{}{p}_2;b)`$ and $`\eta ^0(p;b)`$ respectively, are given by (see Eqs. (17), (21))

$$\eta _2^0(\stackrel{}{p}_1,\stackrel{}{p}_2;b)=\frac{12b^6}{\pi ^3}\mathrm{e}^{b^2\left[p_1^2+p_2^2\right]},$$
(23)
$$\eta ^0(p;b)=\frac{4b^3}{\pi ^{3/2}}\mathrm{e}^{p^2b^2}.$$
(24)
yielding $`\xi (\stackrel{}{p}_1,\stackrel{}{p}_2)=0.75`$. The case falls into the category of non-interacting Fermi systems with $`A=\nu `$, hence Eq. (11) holds (see remark under Eq. (12)). The harmonic oscillator parameter $`b=1.382`$fm reproduces the experimental value of the charge root mean square radius of 4He, $`r^2^{1/2}=1.67`$fm  in this model (corrections due to center-of-mass motion and finite nucleon size are taken into account). We include SRC in the LOA using a two-gaussian (2G) correlation function, $`f\left(r\right)=1c_1\mathrm{exp}\left(r^2/\beta _1^2\right)+c_2\mathrm{exp}\left(r^2/\beta _2^2\right)`$. The expression of the 2bMD is

$`\eta _2^{\mathrm{corr}}(\stackrel{}{p}_1,\stackrel{}{p}_2)=\eta _2^0(\stackrel{}{p}_1,\stackrel{}{p}_2;b)+\delta \eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2;b,c_1,y_1)+\delta \eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2;b,c_2,y_2)`$ (25)
$``$ $`24c_1c_2{\displaystyle \frac{b^6}{\pi ^3}}\{{\displaystyle \frac{1}{\left[\left(1+4y_1\right)\left(1+4y_2\right)\right]^{3/2}}}\mathrm{e}^{b^2\left[\frac{1}{2}\left(\frac{1+2y_1}{1+4y_1}+\frac{1+2y_2}{1+4y_2}\right)\left(p_1^2+p_2^2\right)+\left(\frac{2y_1}{1+4y_1}+\frac{2y_2}{1+4y_2}\right)\stackrel{}{p}_1\stackrel{}{p}_2\right]}`$
$`+`$ $`{\displaystyle \frac{2}{\left(1+3y_1+3y_2+8y_1y_2\right)^{3/2}}}[\mathrm{e}^{\frac{1+2y_1+2y_2}{1+3y_1+3y_2+8y_1y_2}p_1^2b^2}\mathrm{e}^{p_2^2b^2}+(p_1p_2)]`$
$``$ $`{\displaystyle \frac{5}{\left(1+2y_1+2y_2\right)^{3/2}}}\mathrm{e}^{b^2\left(p_1^2+p_2^2\right)}\}`$
where $`y_ib^2/\beta _i^2`$ and

$`\delta \eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2;b,c,y)`$ $``$ $`{\displaystyle \frac{12b^6}{\pi ^3}}\{[{\displaystyle \frac{10c}{\left(1+2y\right)^{3/2}}}{\displaystyle \frac{5c^2}{\left(1+4y\right)^{3/2}}}]\mathrm{e}^{b^2\left(p_1^2+p_2^2\right)}`$ (26)
$``$ $`{\displaystyle \frac{2c}{\left(1+4y\right)^{3/2}}}\mathrm{e}^{b^2\left[\frac{1+3y}{1+4y}\left(p_1^2+p_2^2\right)+\frac{2y}{1+4y}\stackrel{}{p}_1\stackrel{}{p}_2\right]}`$
$`+`$ $`{\displaystyle \frac{c^2}{\left(1+4y\right)^{3/2}}}\mathrm{e}^{b^2\left[\frac{1+2y}{1+4y}\left(p_1^2+p_2^2\right)+\frac{4y}{1+4y}\stackrel{}{p}_1\stackrel{}{p}_2\right]}`$
$`+`$ $`{\displaystyle \frac{2c^2}{\left[\left(1+2y\right)\left(1+4y\right)\right]^{3/2}}}[\mathrm{e}^{\frac{1}{1+2y}p_1^2b^2}\mathrm{e}^{p_2^2b^2}+(p_1p_2)]`$
$``$ $`{\displaystyle \frac{4c}{\left(1+3y\right)^{3/2}}}[\mathrm{e}^{\frac{1+2y}{1+3y}p_1^2b^2}\mathrm{e}^{p_2^2b^2}+(p_1p_2)]\}.`$
The corresponding expression for the momentum distribution $`\eta ^{\mathrm{c}orr}\left(p\right)`$ of 4He is

$`\eta ^{\mathrm{c}orr}\left(p\right)=\eta ^0(p;b)+\delta \eta (p;b,c_1,y_1)+\delta \eta (p;b,c_2,y_2)`$
$`6c_1c_2{\displaystyle \frac{b^3}{\pi ^{3/2}}}\left[{\displaystyle \frac{\mathrm{e}^{p^2b^2\frac{1+2y_1+2y_2}{1+3y_1+3y_2+8y_1y_2}}}{\left(1+3y_1+3y_2+8y_1y_2\right)^{3/2}}}{\displaystyle \frac{\mathrm{e}^{p^2b^2}}{\left(1+2y_1+2y_2\right)^{3/2}}}\right],`$ (27)
where

$$\delta \eta (p;b,c,y)\frac{12b^3}{\pi ^{3/2}}\left[\frac{c^2\mathrm{e}^{p^2b^2\frac{1}{1+2y}}}{[(1+2y)(1+4y)^{3/2}}\frac{2c\mathrm{e}^{p^2b^2\frac{1+2y}{1+3y}}}{\left(1+3y\right)^{3/2}}\right].$$
(28)
If we put $`c_2=0`$ in the above equations we end up with the corresponding expressions for the single gaussian (1G) correlation function, $`f\left(r\right)=1c\mathrm{exp}\left(r^2/\beta ^2\right)`$ which has been used in evaluating several quantities of 4He . We have introduced the 2G correlation function, that can overshoot unity at intermediate distances. This was found necessary in order to eliminate the negative values of $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$ that result at some regions of intermediate values of momenta whenever a 1G correlation function is used.

The parameterization of the new, 2G correlation function by means of a rigorous procedure such as minimization and fitting to some set of experimental data (e.g. charge form factor data) lies beyond the scope of this work, therefore it was obtained somehow heuristically. First, it has been chosen to maintain the small-$`r`$ behaviour of the 1G correlation function of Ref. , by setting $`c_1c_2=0.8`$ and $`\beta _1=0.8`$fm. Subsequently, the extra parameters have been determined requiring that the healing distance is equal to or somewhat lower than 1fm and the relative pair density distribution $`\rho _2\left(r_{12}\right)`$ is approximately zero at $`r_{12}=0`$. The value of the harmonic oscillator parameter $`b`$ was again chosen so as to reproduce the experimental charge root mean square radius of 4He. We ended up with $`c_2=1.7`$ (and therefore $`c_1=0.9`$), $`\beta _2=1.1`$fm and $`b=1.4`$fm. In the left panel of Fig. 3 the resulting 2G correlation function, that overshoots unity, is plotted along with the 1G one of Ref. . The corresponding values of a wound parameter $`\kappa =\left(1f\left(r_{12}\right)\right)^2\rho _2^0\left(r_{12}\right)\mathrm{d}^3r_{12}`$ (where $`\rho _2^0\left(r_{12}\right)`$ is the relative pair density distribution function calculated in the harmonic oscillator model and normalized to unity) for the 2G and 1G correlations functions are equal to 0.008 and 0.018 respectively. In the center panel of the same figure the momentum distribution per particle $`\eta \left(p\right)/A`$ is plotted for both correlation functions, 2G and 1G, along with the result of a variational calculation . In the right panel of Fig. 3 the proton-proton relative pair density distribution (normalized to unity) $`\rho _{pp}\left(r_{12}\right)`$ is shown for both correlation functions 2G, 1G. The bars cover the results of several variational calculations using a variety of two- and three-nucleon interaction models . The values of $`\eta \left(p\right)`$ and $`\rho _{pp}\left(r_{12}\right)`$ calculated with the two-gaussian correlation function compare quite well with the ones of other methods.

Some numerical results for $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$ are presented in Figs. 4 and 5. In Fig. 4, the 2bMD of 4He is plotted in logarithmic scale, for $`\stackrel{}{p}_2`$ parallel to $`\stackrel{}{p}_1`$, $`\stackrel{}{p}_2=`$p$`{}_{2}{}^{}\widehat{p}_{1}^{}`$, as a function of p2, for $`p_1=`$ 0, 1, 2, 3 fm-1. The correlated 2bMD, as given by Eq. (25), is plotted with continuous lines and the harmonic-oscillator result, Eq. (23), with dashed lines. It seems that, within the present approximation, the correlated 2bMD shows deviations from the harmonic-oscillator picture mainly at high values of $`p_1`$ and/or $`p_2`$, developing high-momentum tails. In most cases, the deviations are larger when $`\stackrel{}{p}_1`$, $`\stackrel{}{p}_2`$ are antiparallel rather than parallel.

In Fig. 5 the quantity $`\xi (\stackrel{}{p}_1,\stackrel{}{p}_2)`$ for 4He is plotted in logarithmic scale, as a function of $`\mathrm{cos}\gamma `$ for $`p_1=p_2=`$ 2 and 4 fm-1. The continuous line includes SRC and the dashed line corresponds to the harmonic oscillator model and it equals 0.75. Deviations from this value shows the effect of SRC. The correlated $`\xi `$ in the present approximation grows larger than unity at backward angles, and more so for large values of $`p_1`$, $`p_2`$, as indicated by the $`p_1=p_2=4`$fm-1 graph. In this latter graph, deviations from 0.75 are large even at forward angles. We should note that for low values of $`p_1`$ and $`p_2`$ (not shown), the correlated $`\xi `$ takes values between about 0.7 and 1.

We have also compared our results for the 2bMD with those obtained by using the single-gaussian correlation function of Ref. . Although the overall behaviour of the 2bMD is not dramatically affected by the choice of $`f\left(r\right)`$, the two sets of results differ largely at certain kinematical regions, in particular at intermediate values of momenta where the single-gaussian correlation function yields negative results. (If the 2DM is approximated by only the first terms of the expression (22), i.e. $`\rho _2^{\mathrm{c}orr}=f\left(r_{12}\right)f\left(r_{1^{}2^{}}\right)\rho _2`$, these regions of negative values show up as exact zeros at certain points .) As anticipated in Sec. 4.1, a similar, unphysical situation has been encountered in the study of other two-body quantities, when the LOA is used . It has been argued  that such regions indicate where the SRC start to show up dramatically, and are most sensitive to the approximation employed to account for the SRC. The results in these regions are also sensitive to the particular correlation function used. Let us note here that even in the case of the 2G correlation function a somewhat lower value of the oscillator parameter (eg. $`b=1.35`$fm) retrieves the negative values of the 2bMD, indicating that the uncorrelated single-particle wavefunctions also play some role.

## 5 Summary and Conclusions

In summary, a compact analytical expression has been derived for the two-body momentum distribution, $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$, of $`Z=N`$, $`\mathrm{}`$-closed nuclei within the context of the independent-particle shell model, using harmonic-oscillator wavefunctions. An application to the 16O nucleus revealed interesting features arising from the Fermi statistics and the finite nuclear size. Next, the effect of short-range correlations was investigated by including Jastrow-type correlations in the calculation of $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$ of the 4He nucleus, using the low-order approximation of Ref.  and a two-gaussian correlation function. Significant deviations from the independent-particle picture were found for large values of $`p_1`$ and $`p_2`$ and for angles between the vectors $`\stackrel{}{p}_1`$ and $`\stackrel{}{p}_2`$ close to $`\gamma =`$180 and 0. The quality of our approximation was also explored. It was found that, at certain regions of momenta, the results are sensitive to the correlation function used.

Using the same approximation, one could evaluate the 2bMD and investigate the effects of short-range correlations in other interesting nuclei, for example in 16O and 40Ca. Furthermore, a general expression within the harmonic oscillator model applicable to every $`\mathrm{}`$closed nucleus could be derived for the two-body momentum distribution $`\eta _2(\stackrel{}{p}_1,\stackrel{}{p}_2)`$, the two-body center-of-mass momentum distribution $`\eta _2\left(\stackrel{}{P}_{\mathrm{C}M}\right)`$ and the two-body relative momentum distribution $`\eta _2\left(\stackrel{}{q}_{\mathrm{r}el}\right)`$, as well as for the two-hole spectral function $`S(\stackrel{}{p}_1,\stackrel{}{p}_2;E)`$. The above calculations are relevant to the study of the two-nucleon knock-out reactions.

Acknowledgements

Partial financial support from the University of Athens under Grant 70/4/3309 is acknowledged. P.P. acknowledges a scholarship from the Greek State Scholarship Foundation (IKY).


text:
# 1 Plot of 𝑤_𝐼-𝑤_𝑄 as a function of 𝛽² for 𝑛=2. The long dashed curve corresponds to 𝑡/𝜏=1, the short dashed curve corresponds to 𝑡/𝜏=5, and the solid curve corresponds to 𝑡/𝜏=10.

## Acknowledgments

We would like to thank J. Cline, T. Dent, G. Gabadadze, and R. Rattazzi for helpful discussions.


text:
# Two pairing parameters in superconducting grains

## I Introduction

The question of how superconductivity is affected by small sample size, raised already by Anderson in 1959, has experienced a recent revival of interest, which started with the experimental work of Ralph, Black and Tinkham. In these experiments it was shown that small superconducting grains, with size much smaller than the coherence length, but level spacing $`d`$ smaller than the bulk gap $`\mathrm{\Delta }`$, have a gap of order $`2\mathrm{\Delta }`$ in their tunneling spectrum, and grains in the regime $`d\mathrm{\Delta }`$ (to be called “ultrasmall” regime) do not. Since the pairing parameter is the basic quantity in bulk superconductivity, several efforts have been made to define pairing parameters which are relevant in the regime of small grains, and reduce to $`\mathrm{\Delta }`$ in the bulk limit. The need for such new definitions arises for two apparent reasons: (i) for ultrasmall grains, in which $`d>\mathrm{\Delta }`$, the quantity $`\mathrm{\Delta }`$ has no direct physical meaning; (ii) in both ultrasmall grains, and small grains for which $`d\mathrm{\Delta }`$, if the grains are isolated the appropriate ensemble is the canonical one, in which the usual definition of $`\mathrm{\Delta }`$ that is used in the grand canonical ensemble \[see $`\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}`$ in Eq. (3)\] is trivially zero. In this paper we show that there is a third, fundamental reason for the inadequacy of $`\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}`$ to describe pairing correlations in small superconducting grains. Unlike the situation in bulk superconductors, in which all the various superconducting properties can be characterized by one energy scale, $`\mathrm{\Delta }`$, in general there exist two distinct energy scales that characterize pairing correlations, and the difference between them becomes important particularly for small grains. The first such energy scale, which we denote by $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$, characterizes single-particle properties, such as excitation energies and even-odd effects. The second, denoted as $`\mathrm{\Delta }_{\mathrm{cu}}`$, characterizes cumulative properties, such as the condensation energy \[see Eqs. (12), (20)\], to which pairing correlations of all the levels up to $`\omega _\mathrm{D}`$ contribute. The reason why these two scales are in general distinct is that levels far from the Fermi energy, namely those with energy $`|\epsilon E_\mathrm{F}|`$ between $`\mathrm{\Delta }`$ and the Debye frequency $`\omega _\mathrm{D}`$ (to be called “not condensed” or “far” levels), make a more significant contribution to cumulative properties than to single-particle ones. (For a discussion on the role of the far levels in superconducting grains and persistent currents in normal metal rings see Refs. ). The contribution of the far levels to physical properties turns out to be proportional to the level spacing; in the large grain limit where the level spacing becomes exceedingly small, this contribution can thus be neglected. In this limit, the single-particle and cumulative properties are therefore both determined by the correlations of only the “condensed levels”, i.e. those within $`\mathrm{\Delta }`$ of $`E_\mathrm{F}`$, so that the two scales $`\mathrm{\Delta }_{\mathrm{cu}}`$ and $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ become identical.

Previous attempts to define, in terms of pair correlation functions, a pairing parameter adequate to describe small grains in the canonical ensemble resulted in parameters characterizing cumulative properties of the grain, but not single particle ones. We shall discuss a particular definition for such a cumulative parameter, denoted $`\mathrm{\Delta }_{\mathrm{cu}}`$, which is purposefully chosen such that $`\mathrm{\Delta }_{\mathrm{cu}}`$ reduces to $`\mathrm{\Delta }`$ in the bulk limit. For small grains, we show the correspondence of this definition to cumulative properties such as the condensation energy. We then define a single-particle parameter $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ and show the correspondence of this definition to single-particle properties of the grain. In the bulk limit both $`\mathrm{\Delta }_{\mathrm{cu}}`$ and $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ reduce to $`\mathrm{\Delta }`$ as expected, but at a different size scale (see Fig. 1). For the single-particle properties, we find $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}\mathrm{\Delta }`$ in the ultrasmall regime ($`d>\mathrm{\Delta }`$), and $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}\mathrm{\Delta }`$ for larger grains. In contrast, for the cumulative properties we find $`\mathrm{\Delta }_{\mathrm{cu}}\mathrm{\Delta }`$ not only in the ultrasmall regime, but also in an “intermediate regime” in which $`\mathrm{\Delta }>d>\mathrm{\Delta }^2/\omega _\mathrm{D}`$. Regarding the relation between cumulative and single-particle parameters, we find $`\mathrm{\Delta }_{\mathrm{cu}}\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ in both the ultrasmall and the intermediate regimes. Note that both the range of the intermediate regime, and the value of the pairing parameters and the various superconducting properties obtained below, depend explicitly on $`\omega _\mathrm{D}`$. This $`\omega _\mathrm{D}`$-dependence, as well as the very existence of the intermediate regime, are direct manifestations of the contribution of the pairing correlations of the far levels. Finally, we discuss some possible strategies for measuring the single-particle and cumulative quantities in small grains. Throughout the paper, except when discussing the even-odd effect, we assume for simplicity an even number of electrons in the grain.

## II Pairing correlations in the grand-canonical ensemble

We consider the reduced BCS Hamiltonian

$$\widehat{H}=\underset{j,\sigma =\pm }{}ϵ_jc_{j\sigma }^{}c_{j\sigma }\lambda d\underset{i,j}{\overset{}{}}c_{i+}^{}c_i^{}c_jc_{j+}.$$
(1)
where the second sum (the pairing interaction) is restricted to energies within $`\omega _\mathrm{D}`$ of $`E_\mathrm{F}`$, and $`+()`$ denote spin up (down). The Hamiltonian (1) is the usual BCS Hamiltonian used when discussing superconducting grains  and its validity is discussed in, e.g., Refs. . For a grain with a given, finite number of electrons this Hamiltonian has an exact solution, obtained by Richardson, and independently by Gaudin. In the macroscopic limit, where the canonical and grand-canonical ensembles produce the same results, Richardson’s solution reduces to the BCS solution. As a result of the pairing interaction, the ground state of a superconductor is different from the noninteracting Fermi state. It is a coherent superposition of various paired many-body noninteracting states, defined as Slater determinant of real one electron wavefunctions. The coherence means that the amplitudes for all these states in the superposition are real, up to an overall, global phase factor. This is true both for the BCS wave function in the grand canonical formalism

$$|\mathrm{BCS}=\underset{j}{}(u_j+v_jb_j^{})|\mathrm{Vac}$$
(2)
where $`b_jc_jc_{j+}`$, and for the exact wave function given by Richardson’s solution. (In the grand canonical formalism, the above coherence relates to noninteracting many-body states with a given number of pairs. One can add a constant phase between states of different number of pairs, which is referred to as the superconducting phase). A particular characteristic of the structure of the ground state of a superconductor is that the occupation probability of levels above the Fermi energy is *nonzero*, and that of levels below the Fermi energy is smaller than unity. A measure of the pairing correlations, which exploits both the non Fermi like occupation probability and the phase coherence mentioned above, is given by the pair amplitude, which in the g.c. ensemble is given by

$$\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}\lambda d\underset{j}{}c_jc_{j+}.$$
(3)
For any many-body BCS-like eigenstate, characterized by the set $`\{f_{j\sigma }\}`$ of the occupation probabilities of the BCS quasiparticles, one obtains

$$\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}(\{f_{j\sigma }\})=\lambda d\underset{j}{}u_jv_j^{}(1f_{j+}f_j).$$
(4)
Specifically, for the ground state given in Eq. (2), for which $`f_{j\sigma }=0`$ for all $`j`$, we obtain

$$\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}(\mathrm{g}.\mathrm{s}.)=\lambda d\underset{j}{}u_jv_j,$$
(5)
which gives in the bulk limit $`\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}(\mathrm{g}.\mathrm{s}.)=\omega _\mathrm{D}/\mathrm{sinh}(1/\lambda )\mathrm{\Delta }`$. The mean occupation of level $`j`$ is given by $`v_j^2`$, and $`u_j^2+v_j^2=1`$. Note that the nature of $`\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}`$ is cumulative, being a sum of the contributions of all levels. As a result of the phase coherence mentioned above, all the contributions to the sum in Eq. (5) come with the same phase.

## III Cumulative pairing parameter $`\mathrm{\Delta }_{\mathrm{cu}}`$

We now turn to the definition of the cumulative canonical pairing parameter. A natural extension of the concept of long range order in the bulk, suggests the following definition for a canonical pairing parameter:

$$|\mathrm{\Delta }_{\mathrm{cu}}|^2(\lambda d)^2𝑑r_1𝑑r_2F(r_1,r_2)$$
(6)
where $`F(r_1,r_2)`$ is a function characterizing pairing correlations, which is given by

$`F(r_1,r_2)\psi _+^{}(r_1)\psi _{}^{}(r_1)\psi _{}(r_2)\psi _+(r_2)`$ (7)
$`\psi _+^{}(r_1)\psi _+(r_2)\psi _{}^{}(r_1)\psi _{}(r_2).`$ (8)
We will show that this definition for the cumulative canonical pairing parameter is adequate, as it has the following properties: (i) it is meaningful for a canonical ensemble; (ii) in the bulk limit it is equivalent to the grand canonical definition (3) for any given BCS eigenstate; (iii) it is related to the condensation energy and to the mean occupation of the noninteracting levels not only in the bulk limit, but also in the opposite limit of ultrasmall grains \[see Eqs. (1213) below\].

Expanding each of the $`\psi `$ operators in the basis of the noninteracting single-particle eigenstates, we obtain

$$|\mathrm{\Delta }_{\mathrm{cu}}|^2=(\lambda d)^2\underset{ij}{}(b_i^{}b_jc_{i+}^{}c_{j+}c_i^{}c_j).$$
(9)
Since the terms in each of the brackets are number conserving, they are meaningful in the canonical ensemble, which is used below for the evaluation of $`|\mathrm{\Delta }_{\mathrm{cu}}|^2`$ for small grains. However, in the macroscopic limit, one can use the grand canonical ensemble to evaluate $`|\mathrm{\Delta }_{\mathrm{cu}}|^2`$ within the BCS approximation. By using the Bogoliubov transformation and the BCS wavefunctions, we find that for any many-body BCS-like eigenstate

$`|\mathrm{\Delta }_{\mathrm{cu}}|^2=`$ $`(\lambda d)^2{\displaystyle \underset{ij}{}}u_iv_i^{}u_j^{}v_j(1f_{i+}f_i)(1f_{j+}f_j)`$ (10)
$`=|\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}|^2,`$ (11)
where the last equality is a result of Eq. (4). The canonical pairing parameter is equal, in the bulk limit, to the grand canonical pairing parameter for any many-body state, and therefore both definitions are equivalent in this limit. Our definition for $`\mathrm{\Delta }_{\mathrm{cu}}`$ differs slightly from that given in Ref. . The difference lies in the last term in Eq. (9), which results in the exact equivalence to $`\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}`$ in Eq. (11).

We now turn to the opposite limit of ultrasmall grains. By examining Richardson’s exact solution, it was shown that pairing correlations, however small \[i.e. even for $`\lambda 1/\mathrm{ln}(\omega _\mathrm{D}/d)`$\] manifest themselves both in the form of the ground state wave function and in a finite condensation energy $`E_{\mathrm{cond}}`$. Here we make the connection between these two effects and the canonical pairing parameter. In particular, we shall show that to leading order in $`\lambda `$

$$E_{\mathrm{cond}}=\frac{|\mathrm{\Delta }_{\mathrm{cu}}|^2}{2\lambda d},$$
(12)
where $`|\mathrm{\Delta }_{\mathrm{cu}}|^2`$ is evaluated for the exact ground state (this result differs from the known bulk result $`E_{\mathrm{cond}}=\mathrm{\Delta }^2/(2d)`$ by the occurrence of $`\lambda `$ in the denominator, which is discussed below). Also, $`\mathrm{\Delta }_{\mathrm{cu}}`$ is related to the sum of the pairing correlations in all the levels, as are reflected in their occupation probabilities, since we shall show that

$$|\mathrm{\Delta }_{\mathrm{cu}}|^2=\frac{2\mathrm{ln}2}{\lambda }\mathrm{\Delta }_{\mathrm{occ}}^2,$$
(13)
where

$$\mathrm{\Delta }_{\mathrm{occ}}\lambda d\underset{j}{}\overline{u}_j\overline{v}_j,$$
(14)
with $`\overline{v}_j^2b_j^{}b_j`$ and $`\overline{u}_j^2=1\overline{v}_j^2`$. The quantity $`\mathrm{\Delta }_{\mathrm{occ}}`$, which has been defined in analogy to $`\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}`$, reflects the non-Fermi-like mean occupation probability in the many-body ground state of the noninteracting single-particle levels.

These results are found as follows: To obtain the value of $`|\mathrm{\Delta }_{\mathrm{cu}}|^2`$ in the ground state we use Richardson’s equations and expressions for the wavefunctions , to write the ground state to leading order in $`\lambda `$, in terms of the amplitudes of the various noninteracting many-body states appearing in it:

$`\varphi _{\mathrm{g}.\mathrm{s}.}(1,\mathrm{},N)`$ = $`1`$ (15)
$`\varphi _{\mathrm{g}.\mathrm{s}.}(1,\mathrm{},N;j,k)`$ = $`{\displaystyle \frac{\lambda d}{2(ϵ_kϵ_j)}}.`$ (16)
The first line refers to the amplitude that all the levels below the Fermi energy are occupied by pairs (noninteracting ground state of the system). The second line is the amplitude for the noninteracting many-body state which is the same as the noninteracting ground state, except for one pair excitation from level $`j`$ to level $`k`$. The amplitudes of all the other many-body noninteracting states are zero to first order in $`\lambda `$.

From Eq. (15), to zeroth order in $`\lambda `$, the two terms in Eq. (9) cancel each other. The second term in Eq. (9) has no contribution to first order in $`\lambda `$. The contribution to the first term which is first order in $`\lambda `$ comes from the fact that there is finite amplitude for levels above the Fermi energy to be occupied, and is given, using Eq. (15), by $`\lambda d/2(ϵ_jϵ_i)`$ for each $`j>,i<`$ (by $`<(>)`$ we refer to levels below (above) the Fermi energy). Since the sum in Eq. (9) is unrestricted, we get a factor of two, and to leading order in $`\lambda `$

$$|\mathrm{\Delta }_{\mathrm{cu}}|^2=(\lambda d)^3\underset{i<,j>}{}\frac{1}{ϵ_jϵ_i}.$$
(17)
The condensation energy was calculated in Ref. . In order to compare it with $`\mathrm{\Delta }_{\mathrm{cu}}`$ of Eq. (17) we present here its value for a general spectrum, which to leading order in $`\lambda `$ is given by

$$E_{\mathrm{cond}}=(\lambda d)^2\underset{i<,j>}{}\frac{1}{2ϵ_j2ϵ_i}.$$
(18)
This result is obtained directly from Richardson’s equations, and leads to Eq. (12).

Evaluating Eq. (17) for equally-spaced spectrum, we obtain

$$|\mathrm{\Delta }_{\mathrm{cu}}|^2=2\mathrm{ln}2\lambda ^3d\omega _\mathrm{D}.$$
(19)
The large magnitude (linear in $`\omega _\mathrm{D}`$) of this result is due to the fact that all the amplitudes in (15) have the same phase (which is a consequence of the coherence discussed in Sec. II), so that all the terms in Eq. (17) are added with the same sign.

In Ref.  we have shown that $`E_{\mathrm{cond}}E_{\mathrm{cond}}^{\mathrm{BCS}}+E_{\mathrm{cond}}^{\mathrm{pert}}`$, where $`E_{\mathrm{cond}}^{\mathrm{BCS}}=\mathrm{\Delta }^2/2d`$ is the contribution of the condensed levels, and $`E_{\mathrm{cond}}^{\mathrm{pert}}\lambda ^2\omega _\mathrm{D}`$ is the contribution of the far (not-condensed) levels, which can be calculated perturbatively. Similarly, $`|\mathrm{\Delta }_{\mathrm{cu}}|^2`$ is a sum of the contributions of all levels below $`\omega _\mathrm{D}`$, too, and is related to the condensation energy in both the BCS and the perturbative regimes. Therefore, it is natural to hypothesize that $`|\mathrm{\Delta }_{\mathrm{cu}}|^2`$ (similarly to $`E_{\mathrm{cond}}`$) can to a good approximation be written as the sum of a bulk contribution from the condensed levels and a perturbative contribution from the far levels, i.e. $`|\mathrm{\Delta }_{\mathrm{cu}}|^22\mathrm{ln}2\lambda ^3d\omega _\mathrm{D}+\mathrm{\Delta }^2`$. This would imply that $`\mathrm{\Delta }_{\mathrm{cu}}\mathrm{\Delta }`$ in both the ultrasmall and intermediate regimes (see Fig 1). Note, that the expression for the condensation energy in Eq. (12) is different from the bulk expression of the condensation energy in terms of $`\mathrm{\Delta }`$, namely $`E_{\mathrm{cond}}^{\mathrm{bulk}}=\mathrm{\Delta }^2/2d`$, by a factor of $`\lambda `$ in the denominator. Since we have shown that in the bulk limit $`\mathrm{\Delta }_{\mathrm{cu}}=\mathrm{\Delta }_{\mathrm{g}.\mathrm{c}.}`$, we conclude that $`\mathrm{\Delta }_{\mathrm{cu}}`$ and the condensation energy have a different functional dependence on $`\lambda `$. The above hypothesis results in

$$E_{\mathrm{cond}}=\frac{|\mathrm{\Delta }_{\mathrm{cu}}|^2}{2f(\lambda )d},$$
(20)
where $`f(\lambda )`$ is a monotonic function of $`\lambda `$, which equals $`\lambda `$ for $`\lambda 1/\mathrm{ln}N`$ and unity for $`\lambda 1/\mathrm{ln}N`$ (bulk limit).

We now turn to derive the relation (13) between $`\mathrm{\Delta }_{\mathrm{cu}}`$ and the mean occupation probabilities. In the exact ground state, $`\overline{v}_j^2`$ is given by

$$\underset{\{f_1\mathrm{}f_Nj\}}{}|\varphi (f_1,\mathrm{},j,\mathrm{}f_N)|^2,$$
(21)
where $`\{f_1\mathrm{}f_Nj\}`$ is any configuration of $`N`$ levels out of the $`2N`$ noninteracting single-particle levels, that includes the level $`j`$. To leading order in $`\lambda `$, using Eq. (15), we obtain for levels $`j`$ above $`E_\mathrm{F}`$

$$\overline{v}_j^2=(\lambda d)^2\underset{i<}{}\frac{1}{4(ϵ_jϵ_i)^2}.$$
(22)
Using the approximation of constant level spacing, we find that

$$\overline{v}_j^2=\frac{\lambda ^2d}{4ϵ_j},$$
(23)
where $`ϵ_j`$ is measured from the Fermi energy. This important result shows that the mean occupation $`\overline{v}^2(ϵ)`$ is proportional to $`ϵ^1`$, unlike the usual BCS result, where for $`ϵ\mathrm{\Delta }`$ the mean occupation is proportional to $`ϵ^2`$. Since the occupation probability for a single level is, by (23), proportional to the level spacing, this term can be neglected in the bulk limit. However, in finite-size grains, we find that in both the ultrasmall and intermediate regimes (i.e. $`\mathrm{\Delta }<\lambda \sqrt{d\omega _\mathrm{D}}`$), the occupation probability for energies of the order of $`\omega _\mathrm{D}`$ is in fact larger than that given by the BCS approximation.

As a result of Eq. (22), to first order in $`\lambda `$ we find that $`\overline{v}_j=\lambda d\sqrt{_{i<}1/4(ϵ_jϵ_i)^2}`$ and $`\overline{u}_j=1`$. For equally spaced levels, one therefore obtains from Eq. (23) that

$$\mathrm{\Delta }_{\mathrm{occ}}=\lambda d\underset{j}{}\overline{u}_j\overline{v}_j=\lambda ^2\sqrt{d\omega _\mathrm{D}},$$
(24)
which yields the relation between $`\mathrm{\Delta }_{\mathrm{occ}}`$ and $`\mathrm{\Delta }_{\mathrm{cu}}`$ given in Eq. (13). To summarize: the condensation energy, pairing parameter, and mean level occupation all receive significant contributions from the weak pairing correlation of the far (“not-condensed”) levels up to the Debye frequency. As a result, their magnitude is much larger than that given by the BCS approximation, not only in the ultrasmall regime, but also in the intermediate regime where $`\mathrm{\Delta }^2/\omega _\mathrm{D}<d<\mathrm{\Delta }`$.

## IV Single-particle pairing parameter $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$

Unlike the cumulative properties considered above, other superconducting properties, such as the excitation energy and the even-odd effect (quantified by the Matveev–Larkin (ML) parameter, see Ref. ), are related to the blocking of, say, only one or two levels to pairing correlations. Therefore, they do not depend on the correlations between all possible pairs of levels in the grain, but only on the correlations between these selected levels and all the other levels. As a result, *their values in small grains are much smaller than one would get by using the bulk analogy with $`\mathrm{\Delta }_{\mathrm{cu}}`$ as the pairing parameter*. We therefore define

$$\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}^i\lambda d\underset{j}{}(b_i^{}b_j+b_ib_j^{})$$
(25)
where the sum is over the noninteracting single-particle levels and $`i`$ is a selected such level. Since one frequently deals with the lowest energy levels, we define

$$\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}^{\overline{i}}$$
(26)
for $`\overline{i}`$ being the level closest to the Fermi energy (for our considerations the cases that $`\overline{i}`$ is below or above $`E_\mathrm{F}`$ are equivalent, and we take $`\overline{i}`$ to be below $`E_\mathrm{F}`$).

In the bulk limit, using the BCS approximation, we obtain for the ground state

$$\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}^i=2\lambda du_iv_i\underset{j}{}u_jv_j,$$
(27)
and specifically

$$\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}=\lambda d\underset{j}{}u_jv_j=\mathrm{\Delta }.$$
(28)
We now turn to the ultrasmall regime, and evaluate $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ in the ground state to second order in $`\lambda `$. Using Eq. (15) we obtain

$$\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}^i=\lambda d+\underset{j>}{}\frac{(\lambda d)^2}{ϵ_jϵ_i},$$
(29)
and for the equally spaced spectrum

$$\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}=\lambda d+\lambda ^2d\mathrm{ln}\frac{\omega _\mathrm{D}}{d}.$$
(30)
Let us now consider the excitation energy, say $`E_{\mathrm{exc}}`$, of the first excited state of an (even) superconducting grain. This state can be described by having the two levels nearest to $`E_\mathrm{F}`$ singly occupied with probability unity, thus breaking one pair, and the other $`N1`$ pairs occupying the remaining $`2N2`$ levels according to Richardson’s exact solution . The excitation energy has two different contributions. (i) The kinetic energy cost $`d`$ of occupying a level of higher energy with probability unity, and (ii) a pairing energy cost, which is given by $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$. Note that the latter has (iia) a diagonal part, given by $`\lambda d`$, which is related to the excess energy of two electrons occupying the same spatial noninteracting eigenstate, and (iib) an off-diagonal part due to the blocking of the two singly occupied levels to pairing correlations. Adding these contributions gives

$$E_{\mathrm{exc}}=d+\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}d+\lambda d+\lambda ^2d\mathrm{ln}\frac{\omega _\mathrm{D}}{d},$$
(31)
where the last expression is valid for a grain with equally spaced spectrum. Similarly, the ML parameter is given by

$$\mathrm{\Delta }_{ML}=\frac{1}{2}\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}.$$
(32)
Note that the pairing contribution to the energy cost of the first excitation equals $`2\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ in the bulk limit, whereas it equals $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ in an ultrasmall grain, and the same factor of two appears in the ML parameter. This reflects the fact that in the bulk limit each level is correlated with all the other levels, while in the ultrasmall grain, the dominant correlations are between levels on different sides of the Fermi energy.

## V Measurement strategies for $`\mathrm{\Delta }_{\mathrm{cu}}`$ and $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$

A possible measurement of the cumulative correlations of a superconducting grain was discussed in Ref. . It was shown that the condensation energy can be obtained from specific heat or spin magnetization measurements, where an explicit calculation was done for the latter.

In distinction with the case of bulk superconductivity, where the single-particle properties are easiest to measure through the lowest excitations of the system, in small grains the energy of the first excited state is given by Eq. (31), in which $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ is only a small correction to the level spacing. Matveev and Larkin suggested to measure their pairing parameter by the parity-induced alternation of Coulomb blockade peak spacings in small grains. Here we suggest the possibility that for an ensemble of small, weakly-coupled grains, $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ could be measured through their spin magnetization as a function of magnetic field at zero temperature. We assume that the coupling between the grains is weak enough such that the equilibrium properties of the system can be approximated by summing over the individual grains, but strong enough so tunneling between grains in a large enough ensemble (of say 50 grains) occurs within the time of the experiment.

For a single-grain, the spin magnetization shows steps at values that correspond to the Zeeman energy being equal the energy required to break a pair, of which the first is at $`E_{\mathrm{exc}}=d+\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ \[by Eq. (31)\], which is dominated by the large *single-grain* level spacing $`d`$. In order to avoid having to worry about the latter, we consider an ensemble of weakly connected superconducting grains, which would have an effective joint level spacing $`d_{\mathrm{ens}}`$ that is much smaller than $`d`$ if the ensemble is large, $`d_{\mathrm{ens}}d`$ (we neglect the charging energy. Note, that the charging energy was found to be much smaller than what naive estimates give. Since the tunneling between the grains is weak, the relevant ensemble is the canonical one). Then the ground state of the unconnected system is given by each of the grains having an even number of electrons, being in its ground state. The first excited state of the system is given by moving one electron between two grains, so that afterwards the first has a single-occupied level just below $`E_\mathrm{F}`$, and the second has a singly-occupied level just above $`E_\mathrm{F}`$. For a system of “normal” noninteracting grains the energy of such an excitation is of the order of the single-particle level spacing of the whole ensemble of grains, $`d_{\mathrm{ens}}`$. However, in a system of ultrasmall superconducting grains, since two singly occupied levels are created in two different grains, the energy of such an excitation would be given by $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ (which is $`d_{\mathrm{ens}}`$). Since the first step in the spin magnetization is at the energy of the first excited state of the system, in such a measurement we predict a gap of value $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$. (In a similar measurement in grains with $`d<\mathrm{\Delta }`$, the measured gap would be given by $`\mathrm{\Delta }`$). Thus, $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ can be interpreted as the superconducting gap as measured by single particle properties. In principle, performing such a measurement as a function of the grain size would monitor the change of $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$, from being much larger than $`\mathrm{\Delta }`$ for ultrasmall grains to equaling $`\mathrm{\Delta }`$ for large ones. Note by that reducing the grain size one increases this gap $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$.

Finally, we would like to note that the logarithmic dependence of the correlation (second) term of $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ in Eq. (30) is manifested in the interaction correction to the ensemble averaged magnetic response of small metallic grains, when considered within the BCS model.

## VI Summary

We have shown that various superconducting properties can be classified to two groups, one containing those properties which are single particle in nature, and the other containing those properties which are cumulative, a result of the summed contributions of many levels. Unlike the case in bulk superconductivity, where both these properties are characterized by one energy parameter $`\mathrm{\Delta }`$, in general, and in particular in small grains, two different energy parameters, $`\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}`$ for the former group and $`\mathrm{\Delta }_{\mathrm{cu}}`$ for the latter one, characterize superconducting properties.

In bulk superconductivity the relevant contribution to the various superconducting properties comes from the “condensed” levels, within $`\mathrm{\Delta }`$ of $`E_\mathrm{F}`$. However, in small grains the contribution of all levels up to $`\omega _\mathrm{D}`$ is significant, and this results in the existence of these two different parameters, as well as in the fact that $`\mathrm{\Delta }_{\mathrm{cu}}\mathrm{\Delta }_{\mathrm{s}.\mathrm{p}.}\mathrm{\Delta }`$ in the ultrasmall regime and $`\mathrm{\Delta }_{\mathrm{cu}}\mathrm{\Delta }`$ also in the intermediate regime. Experimental possibilities to measure the two pairing parameters were discussed.

## ACKNOWLEDGMENTS

M.S. is thankful for the support by the Lady Davis fund. This work was supported by a Center of Excellence of the Israel Science Foundation, Jerusalem and by the German Federal Ministry of Education and Research (BMBF) within the Framework of the German-Israeli Project Cooperation (DIP) and by the German-Israeli Foundation(GIF).


text:
# CPT-2002/P.4424 October 2002 Hadronic light-by-light scattering contribution to 𝑔_𝜇-2Talk presented at the 6th International Symposium on Radiative Corrections (RADCOR 2002) and the 6th Zeuthen Workshop on Elementary Particle Theory (Loops and Legs in Quantum Field Theory), Kloster Banz, Germany, 8-13 September, 2002.

## 1 Introduction

The muon $`g2`$ is an important quantity that provides a stringent test of the Standard Model and which is potentially sensitive to new physics. However, for this purpose one first needs to well understand the hadronic contributions, i.e. vacuum polarization effects and light-by-light scattering.

The present picture of hadronic light-by-light scattering is shown in Fig. 1 and the corresponding contributions to $`a_\mu `$ are listed in Table 1, taking into account the corrections made in the two full evaluations , after we had discovered the sign error in the pion-pole contribution .

There are three classes of contributions to the hadronic four-point function \[Fig. 1(a)\], which can be understood from an effective field theory (EFT) analysis of hadronic light-by-light scattering : (1) a charged pion loop \[Fig. 1(b)\], where the coupling to photons is dressed by some form factor ($`\rho `$-meson exchange, e.g. via vector meson dominance (VMD)), (2) the pseudoscalar pole diagrams \[Fig. 1(c)\] together with the exchange of heavier resonances ($`f_0,a_1,\mathrm{}`$) and, finally, (3) the irreducible part of the four-point function which was modeled in Refs.  by a constituent quark loop dressed again with VMD form factors \[Fig. 1(d)\]. The latter can be viewed as a local contribution $`\overline{\psi }\sigma ^{\mu \nu }\psi F_{\mu \nu }`$ to $`a_\mu `$. The two groups  used similar, but not identical models which explains the slightly different results for the dressed charged pion and the dressed constituent quark loop, although their sum seems to cancel to a large extent and the final result is essentially given by the pseudoscalar exchange diagrams. We take the difference of the results as indication of the error due to the model dependence.

Our approach to this problem consists of making an ansatz for the relevant Green’s functions in the framework of large-$`N_C`$ QCD. In this limit, an infinite set of narrow resonance states contributes in each channel. Then we perform a matching of the ansatz with chiral perturbation theory (ChPT) at low energies and with the operator product expansion (OPE) at high momenta in order to reduce the model dependence. In practice, it is sufficient to keep a few resonance states to reproduce the leading behavior in ChPT and the OPE. In this way we show in Section 2 that the pseudoscalar contribution now seems under control, due to our semi-analytical calculation , using a pion-photon-photon form factor $`_{\pi ^0\gamma ^{}\gamma ^{}}`$ which fulfills the relevant QCD short-distance constraints , in contrast to the form factors used in Refs. 111Furthermore, the calculations in Refs.  were based purely on numerical approaches. These findings are also corroborated by an EFT and large-$`N_C`$ analysis  which allows to calculate the leading and next-to-leading logarithms in $`a_\mu ^{\text{LbyL;had}}`$ (Sec. 3).

## 2 Pion-pole contribution

The contribution from the neutral pion intermediate state is given by the following two-loop integral (see Ref.  for all the details)

$`a_\mu ^{\text{LbyL;}\pi ^0}=e^6{\displaystyle \frac{d^4q_1}{(2\pi )^4}\frac{d^4q_2}{(2\pi )^4}\frac{1}{q_1^2q_2^2(q_1+q_2)^2}}`$
$`\times {\displaystyle \frac{1}{[(p+q_1)^2m^2][(pq_2)^2m^2]}}`$
$`\times [{\displaystyle \frac{_{\pi ^0\gamma ^{}\gamma ^{}}(q_1^2,(q_1+q_2)^2)_{\pi ^0\gamma ^{}\gamma ^{}}(q_2^2,0)}{q_2^2M_{\pi ^0}^2}}T_1`$
$`+{\displaystyle \frac{_{\pi ^0\gamma ^{}\gamma ^{}}(q_1^2,q_2^2)_{\pi ^0\gamma ^{}\gamma ^{}}((q_1+q_2)^2,0)}{(q_1+q_2)^2M_{\pi ^0}^2}}T_2],`$
(1)
that involves the convolution of two pion-photon-photon transition form factors, see Fig. 1(c). The $`T_i`$ are polynomials of up to sixth order in the momenta $`p,q_1,`$ and $`q_2`$, with $`p^2=m^2`$.

Since no data on the doubly off-shell form factor $`_{\pi ^0\gamma ^{}\gamma ^{}}(q_1^2,q_2^2)`$ is available, one has to resort to models. In order to proceed with the analytical evaluation of the two-loop integrals, we considered a certain class of form factors which includes the ones based on large-$`N_C`$ QCD that we studied in Ref. . For comparison, we have also used a vector meson dominance (VMD) and a constant form factor, derived from the Wess-Zumino-Witten (WZW) term.

In large-$`N_C`$ QCD, the pion-photon-photon form factor is described by a sum over an infinite set of narrow vector resonances, involving arbitrary couplings, although there are constraints at long and short distances. The normalization is given by the WZW term, $`_{\pi ^0\gamma ^{}\gamma ^{}}(0,0)=N_C/(12\pi ^2F_\pi )`$, whereas the OPE tells us that

$`\underset{\lambda \mathrm{}}{lim}_{\pi ^0\gamma ^{}\gamma ^{}}(\lambda ^2q^2,(p\lambda q)^2)`$
$`={\displaystyle \frac{2}{3}}{\displaystyle \frac{F_\pi }{q^2}}\left\{{\displaystyle \frac{1}{\lambda ^2}}+{\displaystyle \frac{1}{\lambda ^3}}{\displaystyle \frac{qp}{q^2}}+𝒪\left({\displaystyle \frac{1}{\lambda ^4}}\right)\right\}.`$ (2)
In the following, we consider the form factors that are obtained by truncation of the infinite sum in large-$`N_C`$ QCD to one (lowest meson dominance, LMD), and two (LMD+V), vector resonances per channel, respectively:

$$_{\pi ^0\gamma ^{}\gamma ^{}}^{\text{LMD}}(q_1^2,q_2^2)=\frac{F_\pi }{3}\frac{q_1^2+q_2^2c_V}{(q_1^2M_V^2)(q_2^2M_V^2)},$$
(3)
$`_{\pi ^0\gamma ^{}\gamma ^{}}^{\text{LMD+V}}(q_1^2,q_2^2)={\displaystyle \frac{F_\pi }{3}}\{(q_1^2q_2^2(q_1^2+q_2^2)+h_1(q_1^2+q_2^2)^2`$
$`+h_2q_1^2q_2^2+h_5(q_1^2+q_2^2)+h_7)/((q_1^2M_{V_1}^2)`$
$`\times (q_1^2M_{V_2}^2)(q_2^2M_{V_1}^2)(q_2^2M_{V_2}^2))\},`$ (4)
with the constants $`c_V=N_CM_V^4/(4\pi ^2F_\pi ^2)`$ and $`h_7=N_CM_{V_1}^4M_{V_2}^4/(4\pi ^2F_\pi ^2)`$. The parameters $`h_1,h_2,`$ and $`h_5`$ in the LMD+V form factor are not fixed by the normalization and the leading term in the OPE. We have determined these coefficients phenomenologically . In particular, $`_{\pi ^0\gamma ^{}\gamma ^{}}(Q^2,0)`$ with one photon on-shell behaves like $`1/Q^2`$ for large spacelike momenta, $`Q^2=q^2`$. Whereas the LMD form factor does not have such a behavior, it can be reproduced with the LMD+V ansatz, provided that $`h_1=0`$. A fit to the data yields moreover $`h_5=6.93\pm 0.26\text{GeV}^4`$. Analyzing the experimental data for the decay $`\pi ^0e^+e^{}`$ leads to the loose bound $`|h_2|<20\text{GeV}^2`$.

Note that the usual VMD form factor $`_{\pi ^0\gamma ^{}\gamma ^{}}^{\text{VMD}}(q_1^2,q_2^2)1/[(q_1^2M_V^2)(q_2^2M_V^2)]`$ does not correctly reproduce the OPE in Eq. (2).

For the form factors discussed above one can perform all angular integrations in the two-loop integrals analytically  using the method of Gegenbauer polynomials . The key observation is that all form factors can be written as follows

$$_{\pi ^0\gamma ^{}\gamma ^{}}(q_1^2,q_2^2)=\stackrel{~}{f}(q_1^2)\underset{M_{V_i}}{}\frac{\stackrel{~}{g}_{M_{V_i}}(q_1^2)}{q_2^2M_{V_i}^2}.$$
(5)
This allows to cancel all dependences on $`q_1q_2`$ in the numerators in $`_{\pi ^0\gamma ^{}\gamma ^{}}(q_1^2,(q_1+q_2)^2)`$ in the loop integrals in Eq. (1). Then one writes the propagators as (for Euclidean momenta $`K,L`$)

$$\frac{1}{(KL)^2+M^2}=\frac{Z_{KL}^M}{|K||L|}\underset{n=0}{\overset{\mathrm{}}{}}(Z_{KL}^M)^nC_n(\widehat{K}\widehat{L}),$$
(6)
with $`Z_{KL}^M=(K^2+L^2+M^2[(K^2+L^2+M^2)^24K^2L^2]^{1/2})/(2|K||L|)`$, and uses the orthogonality properties of the Gegenbauer polynomials

$`{\displaystyle 𝑑\mathrm{\Omega }(\widehat{K})C_n(\widehat{Q}_1\widehat{K})C_m(\widehat{K}\widehat{Q}_2)}`$
$`=`$ $`2\pi ^2{\displaystyle \frac{\delta _{nm}}{n+1}}C_n(\widehat{Q}_1\widehat{Q}_2),`$ (7)
where for instance $`\widehat{Q}_1\widehat{K}`$ denotes the cosine of the angle between the four-dimensional vectors $`Q_1`$ and $`K`$. After performing the angular integrations in this way, the pion-exchange contribution to $`a_\mu `$ can be written as a two-dimensional integral representation, where the integration runs over the moduli of the Euclidean momenta

$`a_\mu ^{\text{LbyL;}\pi ^0}`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}𝑑Q_1{\displaystyle _0^{\mathrm{}}}𝑑Q_2`$ (8)
$`\times {\displaystyle \underset{i}{}}w_i(Q_1,Q_2)f_i(Q_1,Q_2),`$
with universal \[for the above class of form factors\] weight functions $`w_i`$ (rational functions, square roots and logarithms) . The dependence on the form factors resides in $`f_i`$. In this way we could separate the generic features of the pion-pole contribution from the model dependence and thus better control the latter. This is not possible anymore in the final analytical result (as a series expansion) for $`a_\mu ^{\text{LbyL;}\pi ^0}`$ derived in Ref. . Note that the analytical result has not the same status here as for instance in QED. One has to keep in mind that there is an intrinsic uncertainty in the form factor of $`1030\%`$, furthermore the VMD form factor used in that reference has the wrong high-energy behavior.

The weight functions $`w_i`$ in the main contribution are positive and peaked around momenta of the order of $`0.5\text{GeV}`$. There is, however, a tail in one of these functions, which produces for the constant WZW form factor a divergence of the form $`\mathrm{ln}^2\mathrm{\Lambda }`$ for some UV-cutoff $`\mathrm{\Lambda }`$. Other weight functions have positive and negative contributions in the low-energy region, which lead to a strong cancellation in the corresponding integrals.

In Table 2 we present the numerical results for the different form factors.

All form factors lead to very similar results (apart from WZW). Judging from the shape of the weight functions described above, it seems more important to correctly reproduce the slope of the form factor at the origin and the available data at intermediate energies. On the other hand, the asymptotic behavior at large $`Q_i`$ seems not very relevant. The results for the LMD+V form factor are rather stable under the variation of the parameters, except for $`h_2`$. If all other parameters are kept fixed, our result changes in the range $`|h_2|<20\text{GeV}^2`$ by $`\pm 0.9\times 10^{10}`$ from the value for $`h_2=0`$.

Thus, with the LMD+V form factor, we get

$$a_\mu ^{\text{LbyL;}\pi ^0}=+5.8(1.0)\times 10^{10},$$
(9)
where the error includes the variation of the parameters and the intrinsic model dependence. A similar short-distance analysis in the framework of large-$`N_C`$ QCD and including quark mass corrections for the form factors for the $`\eta `$ and $`\eta ^{}`$ was beyond the scope of Ref. . We therefore used VMD form factors fitted to the available data for $`_{\pi ^0\gamma ^{}\gamma ^{}}(Q^2,0)`$ to obtain our final estimate

$`a_\mu ^{\text{LbyL;PS}}`$ $``$ $`a_\mu ^{\text{LbyL;}\pi ^0}+a_\mu ^{\text{LbyL;}\eta }|_{\text{VMD}}+a_\mu ^{\text{LbyL;}\eta ^{}}|_{\text{VMD}}`$ (10)
$`=`$ $`+8.3(1.2)\times 10^{10}.`$
An error of 15 % for the pseudoscalar pole contribution seems reasonable, since we impose many theoretical constraints from long and short distances on the form factors. Furthermore, we use experimental information whenever available.

## 3 EFT approach to $`a_\mu ^{\text{LbyL;had}}`$

In Ref.  we discussed an EFT approach to hadronic light-by-light scattering based on an effective Lagrangian that describes the physics of the Standard Model well below 1 GeV. It includes photons, light leptons, and the pseudoscalar mesons and obeys chiral symmetry and $`U(1)`$ gauge invariance.

The leading contribution to $`a_\mu ^{\text{LbyL;had}}`$, of order $`p^6`$, is given by a finite loop of charged pions with point-like electromagnetic vertices, see Fig. 1(b). Since this contribution involves a loop of hadrons, it is subleading in the large-$`N_C`$ expansion.

At order $`p^8`$ and at leading order in $`N_C`$, we encounter the divergent pion-pole contribution, diagrams (a) and (b) of Fig. 2, involving two WZW vertices. The diagram (c) is actually finite. The divergences of the triangular subgraphs in the diagrams (a) and (b) are removed by inserting the counterterm $`\chi `$ from the Lagrangian $`^{(6)}=(\alpha ^2/4\pi ^2F_0)\chi \overline{\psi }\gamma _\mu \gamma _5\psi ^\mu \pi ^0+\mathrm{}`$, see the one-loop diagrams (d) and (e). Finally, there is an overall divergence of the two-loop diagrams (a) and (b) that is removed by a local counterterm, diagram (f). Since the EFT involves such a local contribution, we will not be able to give a precise numerical prediction for $`a_\mu ^{\text{LbyL;had}}`$.

Nevertheless, it is interesting to consider the leading and next-to-leading logarithms that are in addition enhanced by a factor $`N_C`$ and which can be calculated using the renormalization group . The EFT and large-$`N_C`$ analysis tells us that

$`a_\mu ^{\text{LbyL;had}}=\left({\displaystyle \frac{\alpha }{\pi }}\right)^3\{f({\displaystyle \frac{M_{\pi ^\pm }}{m_\mu }},{\displaystyle \frac{M_{K^\pm }}{m_\mu }})`$
$`+N_C\left({\displaystyle \frac{m_\mu ^2}{16\pi ^2F_\pi ^2}}{\displaystyle \frac{N_C}{3}}\right)\left[\mathrm{ln}^2{\displaystyle \frac{\mu _0}{m_\mu }}+c_1\mathrm{ln}{\displaystyle \frac{\mu _0}{m_\mu }}+c_0\right]`$
$`+𝒪({\displaystyle \frac{m_\mu ^2}{\mu _0^2}}\times \text{log’s})+𝒪({\displaystyle \frac{m_\mu ^4}{\mu _0^4}}N_C\times \text{log’s})\},`$ (11)
where $`f(M_{\pi ^\pm }/m_\mu ,M_{K^\pm }/m_\mu )=0.038`$ represents the charged pion and kaon-loop that is formally of order one in the chiral and $`N_C`$ counting and $`\mu _0`$ denotes some hadronic scale, e.g. $`M_\rho `$. The coefficient $`𝒞`$ of the log-square term in the second line is universal and of order $`N_C`$, since $`F_\pi =𝒪(\sqrt{N_C})`$.

Unfortunately, although the logarithm is sizeable, in $`a_\mu ^{\text{LbyL;}\pi ^0}`$ there occurs a cancellation between the log-square and the log-term. If we fit our result for the VMD form factor for large $`M_\rho `$ to an expression as given in Eq. (11), we obtain

$`a_{\mu ;\text{VMD}}^{\text{LbyL};\pi ^0}`$ $``$ $`\left({\displaystyle \frac{\alpha }{\pi }}\right)^3𝒞\left[\mathrm{ln}^2{\displaystyle \frac{M_\rho }{m_\mu }}+c_1\mathrm{ln}{\displaystyle \frac{M_\rho }{m_\mu }}+c_0\right]`$ (12)
$`\stackrel{\text{Fit}}{=}`$ $`\left({\displaystyle \frac{\alpha }{\pi }}\right)^3𝒞\left[3.943.30+1.08\right]`$
$`=`$ $`\left[12.310.3+3.4\right]\times 10^{10}`$
$`=`$ $`5.4\times 10^{10},`$
which is confirmed by the analytical result of Ref.  (setting for simplicity $`M_{\pi ^0}=m_\mu `$): $`a_{\mu ;\text{VMD}}^{\text{LbyL};\pi ^0}=[128.0+1.7]\times 10^{10}=5.7\times 10^{10}`$. This cancellation is now also visible in the revised version of Ref. . In that paper the remaining parts of $`c_1`$ have been calculated: $`c_1=2\chi (\mu _0)/3+0.237=0.93_{0.83}^{+0.67}`$, with our conventions for $`\chi `$ and $`\chi (M_\rho )_{\mathrm{exp}}=1.75_{1.00}^{+1.25}`$.

Finally, the EFT analysis shows that the modeling of hadronic light-by-light scattering by a constituent quark loop is not consistent with QCD. The latter has a priori nothing to do with the full quark loop in QCD which is dual to the corresponding contribution in terms of hadronic degrees of freedom. Equation (11) tells us that at leading order in $`N_C`$ any model of QCD has to show the behavior $`a_\mu ^{\text{LbyL;had}}(\alpha /\pi )^3N_C[N_Cm_\mu ^2/(48\pi ^2F_\pi ^2)]\mathrm{ln}^2\mathrm{\Lambda }`$, with a universal coefficient, if one sends the cutoff $`\mathrm{\Lambda }`$ to infinity. From the analytical result given in Ref.  for the quark loop, one obtains the behavior $`a_\mu ^{\text{LbyL;CQM}}(\alpha /\pi )^3N_C(m_\mu ^2/M_Q^2)+\mathrm{}`$, for $`M_Qm_\mu `$, if we interpret the constituent quark mass $`M_Q`$ as a hadronic cutoff. Even though one may argue that $`N_C/48\pi ^2F_\pi ^2`$ can be replaced by $`1/M_Q^2`$, the log-square term is not correctly reproduced with this model. Therefore, the constituent quark model (CQM) cannot serve as a reliable description for the dominant contribution to $`a_\mu ^{\text{LbyL;had}}`$, in particular, its sign. Moreover, we note that the pion-pole contribution is infrared finite in the chiral limit222This can be shown by studying the low momentum behavior of the weight functions $`w_i`$ corresponding to the two-loop diagrams 2(a)–(c) and the one-loop diagrams 2(d)+(e) (given in Ref. ) for $`M_{\pi ^0}0`$, see also Ref. ., whereas the quark loop shows an infrared divergence $`\mathrm{ln}(M_Q/m_\mu )`$ for $`M_Q0`$ .

## 4 Conclusions

The pseudoscalar pole contribution $`a_\mu ^{\text{LbyL;PS}}`$ seems to be under control at the 15 % level. Moreover, the EFT and large-$`N_C`$ analysis provides a systematic approach to $`a_\mu ^{\text{LbyL;had}}`$ and yields the leading and next-to-leading logarithmic terms, enhanced by a factor $`N_C`$. It also shows that the modeling of hadronic light-by-light scattering by a constituent quark loop is not consistent with QCD. Since model calculations for the dressed charged pion and the dressed constituent quark loop yield slightly different results, our present estimate reads (by adding the errors linearly)

$$a_\mu ^{\text{LbyL;had}}=+\mathrm{\hspace{0.17em}8}(4)\times 10^{10}.$$
(13)
## Acknowledgments

I am grateful to my collaborators M. Knecht, M. Perrottet and E. de Rafael. I thank the organizers of RADCOR / Loops & Legs 2002 for providing a pleasant atmosphere and the Schweizerischer Nationalfonds for financial support.