text:
# Momentum flux density, kinetic energy density and their fluctuations for one-dimensional confined gases of non-interacting fermions

## I Introduction

The study of ultracold gases of fermionic alkali atoms is being pursued experimentally  by using basically the same trapping and cooling techniques which have led in 1995 to the achievement of Bose-Einstein condensation in gases of bosonic isotopes . The main focus of the experimental effort is towards the realization of novel superfluids from pairing between fermions in different hyperfine states. However, attention has also been given to dilute Fermi gases in a fully spin-polarized state under magnetic confinement. Such a system presents special motives of fundamental interest from the fact that the Pauli principle suppresses the interactions in the s-wave scattering channel and thereby essentially isolates the kinetic energy functional which is invoked in the density functional theory of inhomogeneous fluids .

A quasi-one dimensional Fermi gas can be realized in the experiments  inside a very anisotropic magnetic trap, producing approximately harmonic confinement which is very tight in two directions. The particle and kinetic energy densities of the corresponding one-dimensional (1D) model of non-interacting fermions in a harmonic potential well have been evaluated exactly up to mesoscopic numbers of particles . The shell structure in the particle density which was noticed in earlier studies of 3D Fermi gases  is greatly enhanced in lower dimensionality (see also Ref. ). As a further motive of interest, the 1D model of noninteracting fermions has the same particle density profile and thermodynamic properties as the 1D model of point-like bosons with hard core interactions . The system of interacting bosons in 1D is attracting attention due to its rich phase diagram .

The particle and kinetic energy densities are the zero-th moment and a second moment of the one-body density matrix. In our earlier study  we have outlined a general method by which one may calculate the successive moments of the density matrix in a system of $`N`$ Fermi particles moving independently in a 1D potential. In the fully quantum regime alternative definitions can be given for the higher-order moments, depending on how the operators are symmetrized . In particular, two different second-order moments can be defined for an inhomogeneous gas under confinement, yielding two different physical quantities which are the kinetic energy density and the momentum flux density. The method here proposed is suitable for the higher-order moments in any of their forms. It exploits the definition of the Green’s function operator $`\widehat{G}(x)`$ in coordinate space and allows one to use the numerical methods already developed in solid state physics  to evaluate the density of single-particle states from the Green’s function in the energy domain.

In this paper we give a detailed account of the general method in Sec. II and some new applications to the case of harmonic confinement in the following Sections. In Sec. III we contrast the kinetic energy density with the momentum flux density in the presence of inhomogeneity. In Sec. IV we assess the quantal (zero-temperature) fluctuations in the kinetic energy and momentum flux densities by means of a calculation of fourth moments of the one-body density matrix. Properties of the kinetic energy functional are discussed in Sec. V in the context of density functional theory. An exact local-density method exists for the system of present interest, as shown in a different approach by Lawes and March , and is explicitely exhibited here for various values of the number of Fermi particles. Finally, in Sec. VI we summarize our main conclusions.

## II The method

We consider a system of $`N`$ non-interacting Fermi particles in one spatial dimension subjected to an external confinement given by a real potential $`V_{ext}(x)`$. We may take the eigenfunctions $`\psi _i(x)`$ as real . The one-body density matrix at zero temperature is

$$\rho (x,x_1)=\underset{i=1}{\overset{N}{}}\psi _i(x)\psi _i(x_1).$$
(1)
Its zero-th moment is the equilibrium density profile $`n(x)`$, and its odd moments are zero due to the reality of the wavefunctions. For the even moments of higher order several definitions have been considered in the literature (for a review see e.g. Ziff et al. ), among which we shall focus on the following:

$$P_n(x)=(i)^n\left[\frac{^n}{x_1^n}\rho (x,x_1)\right]_{x_1=x}$$
(2)
and

$$S_n(x)=(i)^n\left[\frac{^n}{r^n}\rho (x+r/2,xr/2)\right]_{r=0}$$
(3)
(we have set $`\mathrm{}=1`$ and $`m=1`$). For $`n2`$ these definitions give rise to different equilibrium properties in an inhomogeneous system, although their semiclassical limits and their integrals coincide. In the specific case of $`n=2`$ Eqs. (2) and (3) give twice the kinetic energy density and the momentum flux density, respectively: these will be discussed in Sec. III.

The method developed in Ref.  is here extended to evaluate both $`P_n(x)`$ and $`S_n(x)`$. For this purpose we write Eqs. (2) and (3) in terms of the Green’s function $`\widehat{G}(x)=(x\widehat{x}+i\epsilon )^1`$ in coordinate space as

$$P_n(x)=\frac{1}{\pi }\underset{\epsilon 0^+}{lim}\mathrm{Im}\underset{i=1}{\overset{N}{}}\psi _i|\widehat{G}(x)\widehat{p}^n|\psi _i$$
(4)
and

$$S_n(x)=\frac{1}{2^n\pi }\underset{\sigma =0}{\overset{n}{}}\left(\begin{array}{c}n\\ \sigma \end{array}\right)\underset{\epsilon 0^+}{lim}\mathrm{Im}\underset{i=1}{\overset{N}{}}\psi _i|\widehat{p}^\sigma \widehat{G}(x)\widehat{p}^{n\sigma }|\psi _i.$$
(5)
Here, $`|\psi _i`$ are the eigenstates of the Hamiltonian in the coordinate representation and $`\widehat{p}`$ is the momentum operator. Equation (4) is easily proved by making use of the expression of $`\widehat{G}(x)`$ in the coordinate representation:

$`P_n(x)`$ $`=`$ $`{\displaystyle \frac{(i)^n}{\pi }}\underset{\epsilon 0^+}{lim}\mathrm{Im}{\displaystyle \underset{i=1}{\overset{N}{}}}{\displaystyle dx_i\psi _i(x_i)\frac{1}{xx_i+i\epsilon }\frac{^n}{x_i^n}\psi _i(x_i)}`$ (6)
$`=`$ $`(i)^n{\displaystyle \underset{i=1}{\overset{N}{}}}{\displaystyle dx_i\psi _i(x_i)\delta (xx_i)\frac{^n}{x_i^n}\psi _i(x_i)},`$ (7)
where $`x_i`$ are the eigenvalues of the position operator $`\widehat{x}`$. Equation (5) can be proved in a similar fashion.

The momentum distribution $`n(p)`$ can be analogously written in terms of the trace of the Green’s function $`\widehat{G}(p)=(p\widehat{p}+i\epsilon )^1`$ in momentum space, taken on the eigenvectors $`|\varphi _i`$ of the Hamiltonian in the momentum representation:

$$n(p)=\underset{i=1}{\overset{N}{}}dp_i|\varphi _i(p_i)|^2\delta (pp_i)$$
(8)
or

$$n(p)=\frac{1}{\pi }\underset{\epsilon 0^+}{lim}\mathrm{Im}\underset{i=1}{\overset{N}{}}\varphi _i|\widehat{G}(p)|\varphi _i,$$
(9)
where $`p_i`$ are the eigenvalues of the momentum operator.

To evaluate the Green’s functions in the specific case of harmonic confinement, we make use of the representation for the position and the momentum operators in the basis of the eigenstates of the harmonic oscillator, i.e. $`\widehat{x}=(a+a^{})/\sqrt{2}`$ and $`\widehat{p}=i(a^{}a)/\sqrt{2}`$ with $`a|\psi _i=\sqrt{i1}|\psi _{i1}`$ and $`a^{}|\psi _i=\sqrt{i}|\psi _{i+1}`$ (we have set the harmonic oscillator frequency $`\omega =1`$). Explicitly, the representation of $`\widehat{x}`$ in matrix form is given by

$$\widehat{x}=\left(\begin{array}{ccccccc}0& 1/\sqrt{2}& & & & & \\ 1/\sqrt{2}& 0& 1& & & & \\ & 1& 0& \sqrt{3/2}& & & \\ & & \sqrt{3/2}& 0& \sqrt{2}& & \\ & & & \sqrt{2}& 0& \sqrt{5/2}& \\ & & & & \mathrm{}& \mathrm{}& \mathrm{}\end{array}\right)$$
(10)
and similarly for $`\widehat{p}`$. The various moments of the one-body density matrix and the momentum distribution are then evaluated as the trace of the operators $`\widehat{p}^\sigma \widehat{G}(x)\widehat{p}^{n\sigma }`$ and of $`\widehat{G}(p)`$ on the first $`N`$ states (TrN). The trace is calculated by the general method given below. The quantities $`n(x)`$ and $`n(p)`$ could also be calculated from a suitably modified Kirkman-Pendry relation , as is shown in Appendix A.

### A Expression for the evaluation of the trace

The partial trace $`\mathrm{Tr}_N`$ of a generic matrix $`Q`$ is related to the determinant of the inverse matrix $`Q^1`$ by

$$\mathrm{Tr}_NQ=\left[\mathrm{ln}det(Q^1+\lambda 𝕀_N)\right]/\lambda |_{\lambda =0}.$$
(11)
This relation is demonstrated as follows:

$`\left[\mathrm{ln}det(Q^1+\lambda 𝕀_N)\right]/\lambda |_{\lambda =0}=\left[\mathrm{Tr}\mathrm{ln}(Q^1+\lambda 𝕀_N)\right]/\lambda |_{\lambda =0}`$ (12)
$`=\mathrm{Tr}\left[(Q^1+\lambda 𝕀_N)^1𝕀_N\right]_{\lambda =0}=\mathrm{Tr}(Q𝕀_N)=\mathrm{Tr}_NQ.`$ (13)
Here we take $`Q=\widehat{p}^\sigma \widehat{G}(x)\widehat{p}^{n\sigma }`$. In the case of harmonic confinement this matrix is the product of a $`(2n+1)`$-diagonal matrix and of the inverse of a tridiagonal matrix. Exploiting Eq. (11) we can write

$$P_n(x)=\frac{1}{\pi }\underset{\epsilon 0^+}{lim}\mathrm{Im}\frac{}{\lambda }\left[\mathrm{ln}det(x+i\epsilon \widehat{K}^{n,0})\right]_{\lambda =0},$$
(14)
and

$$S_n(x)=\frac{1}{2^n\pi }\underset{\sigma =0}{\overset{n}{}}\left(\begin{array}{c}n\\ \sigma \end{array}\right)\underset{\epsilon 0^+}{lim}\mathrm{Im}\frac{}{\lambda }\left[\mathrm{ln}det(x+i\epsilon \widehat{K}^{n,\sigma })\right]_{\lambda =0}$$
(15)
where $`\widehat{K}^{n,\sigma }=\widehat{x}\widehat{p}^{n\sigma }\lambda 𝕀_N\widehat{p}^\sigma `$. We have thus reduced the problem to the evaluation of the determinant of the matrix $`(x+i\epsilon \widehat{K}^{n,\sigma })`$, which is tridiagonal in the tail and whose first $`N`$ rows have only a few non-vanishing elements for low values of $`n`$.

If we make a partition in blocks of dimension $`n\times n`$ for the first part of the matrix, the whole operator can be written as

$$\widehat{K}^{n,\sigma }=\underset{M\mathrm{}}{lim}\widehat{K}_M^{n,\sigma }=\underset{M\mathrm{}}{lim}\left(\begin{array}{cccccc}𝒜_1& _{1,2}& & & & \\ _{2,1}& 𝒜_2& _{2,3}& & & \\ & _{3,2}& \mathrm{}& _{3,4}& & \\ & & \mathrm{}& \mathrm{}& \mathrm{}& \\ & & & \mathrm{}& \mathrm{}& _{M1,M}\hfill \\ & & & & _{M,M1}& 𝒜_M\hfill \end{array}\right).$$
(16)
Here, we have introduced the truncated matrix $`\widehat{K}_M^{n,\sigma }`$, which will be used below in actual numerical calculations.

### B The renormalization procedure

The determinant of the matrix $`(x+i\epsilon \widehat{K}^{n,\sigma })`$ with a tridiagonal representation such as that in Eq. (16) can be factorized into a product of determinants of matrices having the same dimension as that of the blocks of the partition .

As a first step, it is easy to show that if we partition the matrix $`\widehat{K}_M^{n,\sigma }`$ into four blocks, thereby defining $`\stackrel{~}{}_{i,j}`$ an $`\widehat{K}_i^{n,\sigma }`$ from

$$\widehat{K}_M^{n,\sigma }=\left(\begin{array}{cc}\widehat{K}_{M1}^{n,\sigma }& \stackrel{~}{}_{M1,M}\\ \stackrel{~}{}_{M,M1}& 𝒜_M\end{array}\right),$$
(17)
we can write

$$det(x+i\epsilon \widehat{K}_M^{n,\sigma })=det(x+i\epsilon \widehat{K}_{M1}^{n,\sigma })det(x+i\epsilon 𝒜_M\stackrel{~}{}_{M,M1}(x+i\epsilon \widehat{K}_{M1}^{n,\sigma })^1\stackrel{~}{}_{M1,M}).$$
(18)
Applying recursively this procedure and taking the limit $`M\mathrm{}`$ we obtain

$$det(x+i\epsilon \widehat{K}^{n,\sigma })=det(x+i\epsilon 𝒜_1)\underset{M\mathrm{}}{lim}\underset{j=2}{\overset{M}{}}det(x+i\epsilon 𝒜_j\stackrel{~}{}_{j,j1}(x+i\epsilon \widehat{K}_{j1}^{n,\sigma })^1\stackrel{~}{}_{j1,j}).$$
(19)
It is now important to notice that, owing to the particular form of the matrices $`\stackrel{~}{}_{j1,j}`$ and $`\stackrel{~}{}_{j,j1}`$, it is not necessary to explicitly invert the matrices $`(x+i\epsilon \widehat{K}_{j1}^{n,\sigma })`$. Rather, we may use a renormalization procedure for the operator $`\widehat{K}_{j1}^{n,\sigma }`$ to further simplify the expression (19) and to calculate the inverse of matrices having dimension at most equal to $`n\times n`$. Renormalization is based on the idea that if one only needs to describe a given subspace of the whole Hilbert space of the system, one can define an effective operator acting in the subspace and taking into account the contribution from the rest of the system. Renormalization allows us to write

$$det(x+i\epsilon \widehat{K}^{n,\sigma })=\underset{j=1}{\overset{\mathrm{}}{}}det(x+i\epsilon \stackrel{~}{𝒜}_j)$$
(20)
where $`\stackrel{~}{𝒜}_1=𝒜_1`$ and

$$\stackrel{~}{𝒜}_j=𝒜_j+_{j,j1}(x+i\epsilon \stackrel{~}{𝒜}_{j1})^1_{j1,j}$$
(21)
for $`j>1`$. The $`_{j,j+1}`$ have been introduced in Eq. (16). In the specific case of $`n=2`$ and $`\sigma =0`$ Eq. (20) yields back the expression obtained in Ref. .

Operatively, we have studied the convergence of the function $`det(x+i\epsilon \widehat{K}_M^{n,\sigma })`$ on increasing $`M`$. We have performed our calculations up to $`n=4`$ and $`N=1000`$ fermions with $`M10^7`$ and $`\epsilon =0.001`$.

## III Kinetic energy and momentum flux densities

The second moments $`P_2(x)`$ and $`S_2(x)`$ of the one-body density matrix from Eqs. (2) and (3) have different physical meanings and show different behaviours. The function $`P_2(x)=_i^N\psi _i|\delta (xx_i)\widehat{p}^2|\psi _i`$ is twice the kinetic energy density $`T(x)`$  and for harmonic confinement has already been evaluated in Ref.  together with the equilibrium density profile. This quantity is of main interest in the context of Density Functional Theory, as will be discussed in Sec. V. On the other hand, $`S_2(x)=(1/2)_i^N\psi _i|\delta (xx_i)\widehat{p}^2+\widehat{p}\delta (xx_i)\widehat{p}|\psi _i`$ is the momentum flux density $`\mathrm{\Pi }(x)`$ which enters the equations of generalized hydrodynamics . In 3D fluids this quantity becomes a tensor, known as the kinetic stress tensor.

The equation of motion of the particle density $`n(x,t)`$ for non-interacting fermions in 1D is

$$\frac{^2}{t^2}n(x,t)=\frac{^2}{x^2}\mathrm{\Pi }(x,t)+\frac{}{x}\left[n(x,t)\frac{}{x}V_{ext}(x,t)\right].$$
(22)
At equilibrium in a static potential we obtain from Eq. (22) an exact relation between the momentum flux density and the density profile,

$$\frac{d}{dx}\mathrm{\Pi }(x)=n(x)\frac{d}{dx}V_{ext}(x)+\mathrm{const}.$$
(23)
The integration constant can be fixed by imposing suitable boundary conditions for $`x\pm \mathrm{}`$.

Another exact relation can be derived by exploiting Eqs. (1)-(3) together with the definition of $`n(x)`$ as the zero-th moment of the density matrix. This reads:

$$\mathrm{\Pi }(x)=2T(x)+\frac{1}{4}\frac{d^2}{dx^2}n(x).$$
(24)
Equation (24) exhibits the relation which exists between the momentum flux density and the kinetic pressure $`P(x)=2T(x)`$ in the inhomogeneuos system. Evidently, in the homogeneous limit the momentum flux density and the kinetic pressure coincide – and therefore they also coincide within a local-density description for the confined gas. They instead differ within a quantum-mechanical description of the the inhomogeneous system, when one keeps exactly into account the role of confinement. This difference is readily illustrated in the fully quantum case of a single fermion in an harmonic oscillator potential, where we have

$$\mathrm{\Pi }(x)=T(x)+V(x)$$
(25)
with $`V(x)=n(x)V_{ext}(x)`$ being the potential energy density associated with the confinement.

For the case of harmonic confinement we have employed the method described in Sec. II to evaluate the exact profiles for both $`T(x)`$ and $`\mathrm{\Pi }(x)`$ at various numbers of particles. Figure 1 shows these results for $`N`$= 4, 12 and 24 fermions. While $`T(x)`$ shows $`N`$ oscillations and negative tails , $`\mathrm{\Pi }(x)`$ is everywhere positive and its shell structure is visible only in its first derivative (see inset in Fig. 1).

## IV Quantum fluctuations of the kinetic energy and momentum flux densities

From the fourth moments $`S_4(x)`$ and $`P_4(x)`$ of the density matrix we can evaluate the local fluctuations of the momentum flux and kinetic energy densities. For the former the fluctuations are most conveniently introduced through the Wigner distribution function,

$$f_W(p,x)=𝑑re^{ipr}\rho (x+r/2,xr/2).$$
(26)
In this representation the relevant moments are given by $`n(x)=_pf_W(p,x)`$ and by $`\mathrm{\Pi }(x)=_pp^2f_W(p,x)n(x)p^2`$. The mean square fluctuation of the momentum flux density is therefore given by

$$\mathrm{\Delta }\mathrm{\Pi }(x)=\underset{p}{}(p^2p^2)^2f_W(p,x)=S_4(x)\mathrm{\Pi }^2(x)/n(x),$$
(27)
where $`S_4(x)`$ is as defined in Eq. (3). Equation (27) in the coordinate representation for the density matrix can be rewritten as

$$\mathrm{\Delta }\mathrm{\Pi }(x)=\left[\frac{^2}{r^2}\frac{\mathrm{\Pi }(x)}{n(x)}\right]^2\rho (x+r/2,xr/2)|_{r=0}.$$
(28)
The fluctuations of the kinetic energy density are obtained in a similar fashion as

$`\mathrm{\Delta }T(x)`$ $`=`$ $`\left[{\displaystyle \frac{1}{2}}{\displaystyle \frac{^2}{x_1^2}}{\displaystyle \frac{T(x_1)}{n(x_1)}}\right]^2\rho (x,x_1)|_{x_1=x}`$ (29)
$`=`$ $`{\displaystyle \frac{1}{4}}P_4(x){\displaystyle \frac{T^2(x)}{n(x)}}+{\displaystyle \frac{1}{2}}{\displaystyle \frac{d}{dx}}\left[n(x){\displaystyle \frac{d}{dx}}\left({\displaystyle \frac{T(x)}{n(x)}}\right)\right],`$ (30)
with $`P_4(x)`$ as defined in Eq. (2).

For harmonic confinement in the local density approximation (LDA) we have $`n_{LDA}(x)=(2Nx^2)^{1/2}/\pi `$, $`\mathrm{\Pi }_{LDA}(x)=2T_{LDA}(x)=(2Nx^2)^{3/2}/3\pi `$ and $`S_4^{LDA}(x)=P_4^{LDA}(x)=(2Nx^2)^{5/2}/5\pi `$. It is straightforward to show from these expressions that the integrated mean square fluctuations, divided by the square of the corresponding integrated quantity, scale as $`1/N`$. The same scaling is found in the exact numerical results at large $`N`$.

We have used the Green’s function method described in Sec. II to evaluate the exact profiles for both $`\mathrm{\Delta }\mathrm{\Pi }(x)`$ and $`\mathrm{\Delta }T(x)`$ in harmonic confinement, as shown in Fig. 2 for $`N`$ =4, 12 and 24 fermions. The positions of the $`N1`$ maxima in $`\mathrm{\Delta }T(x)`$ are located in correspondence to the minima of the kinetic energy density in Fig. 1 (a), and again negative tails are present at the boundaries. In the case of $`N=2`$ we have tested our numerical results against the analytic expression given for $`\mathrm{\Delta }T(x)`$ by low-order Hermite polynomials, the result being shown in Fig. 3. Evidently, we have obtained in this case a very accurate numerical evaluation of up to the fourth moment of the density matrix.

## V Properties of the kinetic energy functional

We have seen that the second moment $`P_2(x)`$ of the one-body density matrix yields the kinetic energy density $`T(x)=P_2(x)/2`$. The total kinetic energy is obtained by integration over the spatial coordinate and in the spirit of Density-Functional Theory (DFT)  can be viewed as a functional of the particle density $`n(x)`$:

$$E_{kin}[n(x)]=𝑑xT(x)=𝑑x\widehat{T}[n(x)].$$
(31)
Here we have introduced an operator $`\widehat{T}`$ which acts on the density profile to yield

$$\widehat{T}[n(x)]=T(x).$$
(32)
A main point of interest for DFT is to assess the properties of the unknown operator $`\widehat{T}`$.

Knowing exactly the function $`T(x)`$ the specific case of 1D harmonic confinement, we may try to construct a function $`\stackrel{~}{T}(n)=T(x(n))`$ in a local density approach, with $`x(n)`$ being given by inversion of the density profile. The inversion has to be performed locally, since the exact density profile is not globally invertible owing to its shell structure. The resulting $`\stackrel{~}{T}(n)`$, as reported in Fig 4, is a multi-valued function around the values of the particle density corresponding to each local maximum in the profile. While this shows that a local-density approach cannot characterize the operator $`\widehat{T}`$, further progress can be made by studying the derivative of $`\stackrel{~}{T}(n)`$ with respect to the density, which is the local chemical potential $`\stackrel{~}{\mu }_{loc}(n)=d\stackrel{~}{T}(n)/dn`$.

We show here that for the specific case of a 1D harmonic confinement the inverse $`\stackrel{~}{n}(\mu _{loc})`$ of the function $`\stackrel{~}{\mu }_{loc}(n)`$ is a single-valued function. This is most simply seen by inverting the Euler equation

$$\frac{\delta E_{kin}[n(x)]}{\delta n(x)}=\mu V_{ext}(x)\mu _{loc}(x)$$
(33)
in the domain $`x0`$ and by substituting its value into the density profile $`n(x)`$ to obtain

$$\stackrel{~}{n}(\mu _{loc})=n(x(\mu _{loc})).$$
(34)
Here, $`x(\mu _{loc})`$ is the inverse of the function $`\mu _{loc}(x)`$. Because of the symmetry $`xx`$ for the harmonic potential, the region $`x<0`$ is not needed to obtain $`\stackrel{~}{n}(\mu _{loc})`$.

An analytic calculation for $`\stackrel{~}{n}(\mu _{loc})`$ has been reported by Lawes and March  for $`N=1`$ and 2. The numerical method described in Sect. II allows us to easily compute the density profile for numbers of particles up to $`N=1000`$. By combining these results for $`n(x)`$ with the analytic expression for $`x(\mu _{loc})`$, we have obtained the function $`\stackrel{~}{n}(\mu _{loc})`$ which is shown in Fig. 5 for various values of $`N`$. At variance from the LDA prediction, the exact density profile is non-vanishing at negative values of $`\mu _{loc}`$ and shows a number of oscillations (equal to $`N/2`$ for even $`N`$ or to $`(N1)/2`$ for odd $`N`$), with decreasing amplitude as $`N`$ increases.

The property of single-valuedness of $`\stackrel{~}{n}(\mu _{loc})`$ is a characteristic of 1D harmonic confinement and can be understood to be a consequence of the local nature of the relation $`\stackrel{~}{n}(\mu _{loc})`$ vs. $`\mu _{loc}`$, which was already pointed out by Lawes and March . Anharmonicity would be sufficient to make it invalid.

## VI Summary and conclusions

In summary, we have presented a general method for evaluating the $`n`$-th moments of the one-body density matrix (defined either through the derivatives with respect to the relative coordinate or through the derivatives with respect to the second position variable) for a system of confined non-interacting 1D fermions. We have applied this method to the case of harmonic confinement to calculate up to fourth moments.

Regarding second moments, the momentum flux density $`\mathrm{\Pi }(x)`$ has been computed and compared with the kinetic energy density $`T(x)`$ at various numbers of particles. While $`T(x)`$ shows a prominent shell structure and negative tails, $`\mathrm{\Pi }(x)`$ is everywhere positive and has a less marked shell structure, which can be evidenced in its first spatial derivative. The relation between these two different physical quantities has been elucidated by giving an exact relation between them and the equilibrium density profile of the fermion cloud. The momentum flux density determines the dynamical properties of the cloud, whereas the kinetic energy density is relevant in the context of Density Functional Theory.

We have also calculated the fourth moments of the density matrix in order to display the quantal fluctuations of the kinetic energy and momentum flux densities. Finally, we have evaluated the local relationship which exists between density profile and local chemical potential for fermion clouds in 1D harmonic confinement. This relation has allowed us to display an exact property of the DFT kinetic energy functional for a 1D Fermi gas up to mesoscopic values of the number of particles.

###### Acknowledgements.

This work was partially supported by MURST through PRIN2000. One of us (A.M.) would like to thank Dr. F. Illuminati for useful discussions and for giving reference to the work of Ziff et al. .

## A Density profile from the Kirkman-Pendry relation

We derive in this Appendix an expression alternative to Eq. (14) for the evaluation of the density profile in a system of $`N`$ fermions under 1D harmonic confinement. The trace of the Green’s function $`\widehat{G}(x)`$ over its first N states in Eq. (4) is obtained through an extension of the Kirkman-Pendry relation , already used for the density of states in a quasi-1D system .

We define the Green’s function

$$\widehat{𝒢}(\delta ,x)=\frac{1}{x+\delta \widehat{\xi }(x)+i\epsilon }$$
(A1)
where $`\delta `$ is an auxiliary continuous variable and $`\widehat{\xi }(x)`$ is an effective position operator of dimension $`N\times N`$. This is defined by setting $`[\widehat{\xi }(x)]_{i,j}=[\widehat{x}]_{i,j}`$ if $`(i,j)(N,N)`$ and $`[\widehat{\xi }(x)]_{N,N}=\stackrel{~}{x}_{N,N}(x)`$, with

$$\stackrel{~}{x}_{N,N}(x)=\frac{N/2}{x+i\epsilon {\displaystyle \frac{(N+1)/2}{x+i\epsilon \mathrm{}}}}.$$
(A2)
The single term $`\stackrel{~}{x}_{N,N}(x)`$ contains the contribution of all the states which are not occupied by the fermions. Its asymptotic value is $`\stackrel{~}{x}_{N,N}(x)=i\sqrt{2Nx^2}`$ for $`N\mathrm{}`$.

An expression for $`\mathrm{Tr}_N\widehat{G}(x)`$ is obtained from the Green’s function element $`𝒢_{1,N}(\delta ,x)`$ between the first and the last occupied state by using the expression

$`{\displaystyle \frac{}{\delta }}\mathrm{ln}𝒢_{1,N}(\delta ,x)`$ $`=`$ $`{\displaystyle \frac{}{\delta }}\mathrm{ln}{\displaystyle \frac{_{i=1}^{N1}[\widehat{\xi }(x)]_{i,i+1}}{det[x+\delta \widehat{\xi }(x)+i\epsilon ]}}={\displaystyle \frac{}{\delta }}\mathrm{ln}{\displaystyle \frac{_{i=1}^{N1}[\widehat{\xi }(x)]_{i,i+1}}{_{i=1}^N[x+\delta \xi _i(x)+i\epsilon ]}}`$ (A3)
$`=`$ $`{\displaystyle \frac{}{\delta }}\mathrm{ln}{\displaystyle \underset{i=1}{\overset{N}{}}}{\displaystyle \frac{1}{x+\delta \xi _i(x)+i\epsilon }}={\displaystyle \underset{i=1}{\overset{N}{}}}{\displaystyle \frac{1}{x+\delta \xi _i(x)+i\epsilon }}`$ (A4)
where $`\xi _i(x)`$ are the eigenvalues of the operator $`\widehat{\xi }(x)`$. Therefore, the particle density is given by

$$n(x)=\frac{1}{\pi }\underset{\epsilon 0^+}{lim}\mathrm{ImTr}_N\widehat{G}(x)=\frac{1}{\pi }\underset{\epsilon 0^+}{lim}\mathrm{Im}\left[\frac{}{\delta }\mathrm{ln}𝒢_{1,N}(\delta ,x)\right]_{\delta =0}.$$
(A5)
The matrix element $`𝒢_{1,N}(\delta ,x)`$ can be calculated by a further renormalization  of $`\widehat{\xi }(x)`$, without explicitly inverting the matrix $`(x+\delta \widehat{\xi }(x)+i\epsilon )`$. As in the method presented in Sec. II, convergence is achieved by setting $`M=10^7`$ and $`\epsilon =0.001`$.

The same procedure can be used to evaluate the momentum distribution as the trace of an $`N\times N`$ operator $`\widehat{\pi }(p)`$, defined analogously to $`\widehat{\xi }(x)`$. Due to the symmetry between the operators $`\widehat{x}`$ and $`\widehat{p}`$ in the harmonic oscillator Hamiltonian, we obtain essentially the same final expression as for of the particle density $`n(x)`$. Thus, the momentum distribution for $`N`$ fermions under 1D harmonic confinement shows $`N`$ oscillations.

The mapping between hard-core bosons and non-interacting fermions in 1D  allows us to evaluate the particle density profile for mesoscopic systems of strongly correlated bosons by the same method. However, the mapping is restricted to real-space properties and does not apply to the momentum distribution.


format: text

text:
# Summary: ZeV Air Showers Where do we Stand?

## Introduction

The discovery of events with energies above $`10^{20}`$ eV dates back to the 1960’s, to the early days of air shower detection experiments EeVents . In fact the cosmic ray spectrum has been observed as a continuum at all energies since the beginning of the 20th century, apparently only limited by the acceptance of the existing detectors. The prediction of the Greisen-Zatsepin-Kuz’min (GZK) cutoff dates back to the same decade of the 1960’s GZK . The cutoff should appear at energies just above $`410^{19}`$ eV because of proton attenuation in the Cosmic Microwave Background (CMB) . Heavy nuclei are also attenuated in both the infrared and cosmic microwave background radiation fields, at roughly the same energies through photo-disintegration and pair production. The mean free path for these processes is of order a few Mpc and even allowance of successive interactions with small fractional energy loss rises the attenuation length to roughly 50 Mpc. The attenuation distances for photons in the same background fields are even smaller.

If these events are coming from extragalactic distances, as suggested by the close to isotropic distribution in the relatively well known galactic fields, they should show the GZK cutoff just below the $`10^{20}`$ eV range. If on the contrary they are not extragalactic, we are facing with an unknown source of the highest energies particles ever discovered, which is quite close on cosmological scales, either challenging dimensional analysis of acceleration processes or opening up the way to new physics. A large number of hypotheses have been put forward including new astrophysical objects, new particles, new interactions or the violation of well established principles. The nature of most of these hypotheses, which span a number of research fields which are traditionally very far apart, is clearly a sign that we are debating a remarkable problem. Nonetheless this problem has resisted the efforts of theorists and experimentalists chasing for an acceptable solution for over 40 years.

In this very successful conference we have heard about recent developments in the field both from the theoretical and experimental sides. One of the main ideas that sprouts from it is that a new generation of large aperture experiments has just started. It will drastically speed up the remarkably low pace in building up statistics dictated by fluxes below one particle per square kilometer per century. These detectors will also help to determine the primary composition.

In this article I attempt to summarize the material presented in the conference. This is by no means an easy task and inevitably I will present the field in a subjective way and I will make omissions for which I apologize in advance. I first briefly comment on the mystery of the ultra high energy cosmic rays, to then refer to some of the alternative solutions which were discussed that I have divided in three groups, those models that require acceleration, those models that require fragmentation and decay of massive particles and those that avoid the GZK effect. I then discuss the importance of composition, one of the issues that was most addressed in the conference. Finally I refer to the experiments that were discussed and conclude.

## The Post GZK Flux Puzzle

The high energy tail of the measured cosmic ray flux, those particles arriving well above the energy of the Greisen-Zatsepin-Kuz’min cutoff, presents a complex challenge that is still unresolved as was pointed out in many of the talks biermann ; blasi ; cronin ; ellis ; gaisser ; harari ; selvon ; linsley ; masperi ; olinto ; sigl ; swain ; weiler . If we conservatively assume that the high energy end of the cosmic ray spectrum is due to well established particles such as protons or nuclei (that constitute the low energy end of the cosmic ray spectrum) or even photons, we can be fairly confident that their interactions with magnetic fields and background radiation fields are well understood. If these particles are coming from distances exceeding a few tens of Mpc the observed flux should have an imprint of the interactions with the background radiation fields, that is it should display the GZK cutoff at about $`410^{19}`$ eV.

The actual measurements are quite limited to the primary energy spectrum, the arrival directions and some information related to the nature of the arriving particles (mass composition) based on shower development. As more data has accumulated over the years from different experiments and efforts have been made to analyze the combined data alan  it has become rather clear that: a) The data exhibit a hardening of the spectral index at an energy of $`810^{18}`$ eV 111There is however data from the Fly’s Eye experiment, using the fluorescence technique, that suggests the change of slope occurs at about half this energy. and b) That there is no evidence of the GZK cutoff, with data reaching $`310^{20}`$eV. 222The highest energy event seen with the fluorescence technique provides a calorimetric measurement of shower energy, reinforcing the results obtained with particle detector arrays. The change of spectral index is very suggestive of a different component of the spectrum above $`10^{19}`$ eV, the Ultra High Energy Cosmic Rays (UHECR) which provide the central topic of this conference. Throughout this article I will also refer to such a component as the post GZK particle flux.

According to the observations these particles are unlikely to be coming from distances exceeding about 50 Mpc. This distance scale was discussed as an observational ”horizon” for protons (or nuclei) because of their interactions with the background fields linsley . So far only two mechanisms have been suggested by which particles can attain the highest energies, namely: a) Direct acceleration of charged particles and b) Fragmentation of the decay products of other particles.

If the particles are accelerated, the exceptional energies achieved are very demanding on possible sources that could accelerate them. The constraint comes from basic dimensional arguments. The acceleration of a particle of charge $`Ze`$ to an energy $`E`$ requires a minimum value for the product of the size of the accelerator $`L`$ and its magnetic field $`B`$, namely:

$$E<ZeBcL$$
(1)
where $`c`$ is the speed of light. This is traditionally illustrated by the Hillas plot hillas  and several versions of it have been discussed at the conference cronin ; olinto ; sigl ; selvon . Few of the astrophysical objects and structures known satisfy the constraint. Moreover those that do must be very efficient in reaching the maximum energy and many of them are either too large or too distant compared to the 50 Mpc horizon. There are however possible acceleration sites that cannot be excluded yet.

The difficulties in acceleration scenarios has opened the way to alternative proposals in which the post GZK particles either avoid the interactions with the CMB or they are postulated to be a secondary flux produced locally from the decay of other particles.

On pretty general grounds all models for the origin of the post GZK particle flux are very dependent on at least three assumptions in a very interrelated fashion, namely:

* Hypothesis about the possible production sites for these particles which determine both the source distributions and the distance travelled before reaching us.
* Assumptions about the nature of the particles themselves which determine their interactions with the magnetic fields and background radiation fields in their path to us.
* Models for magnetic field distributions in the Galaxy, galactic halo, clusters of galaxies and intergalactic space that condition the distribution of arrival directions.

Part of the puzzle lies in the fact that it is very difficult to extract unconditional conclusions with the limited available data. Each model has to be tested against observations. The interpretation of the information on the post GZK particle arrival directions is completely dependent on these three hypotheses. Moreover the knowledge of the magnetic fields are pretty limited outside the light-matter distribution of our galaxy. As an example it was pointed out that it is possible that the magnetic field distribution provides strong flux magnification and depletions in preferred directions, that could explain the apparent absence of the GZK cutoff and the observed multiplicities in the arrival directions harari .

In a way the difficulties associated with this intricate interdependence of assumptions makes these particles so attractive to different fields of research including particle physics at extreme energies, astrophysical objects and the study of extragalactic magnetic fields.

## Models

The different types of alternatives for the origin of the UHECR were reviewed by Angela Olinto olinto  who classified them in two main groups: Those that push these conventional acceleration ideas to extremes in order to accommodate the data or Zevatrons and those that invoke new physics. Particular attention was paid in her review to the differences in the arrival direction distributions, the spectral features and primary composition that can be expected from different models. It was stressed that most models imply strong requirements on the magnetic fields to fit the observations olinto . A number of proposals were critically discussed by several speakers. I will select a few.

### Standard Fermi Acceleration Models

Stochastic particle acceleration in the interaction of particles with astrophysical shocks is the conventional astrophysical answer to the question about the origin of the cosmic rays below the GZK cutoff. The mechanism can be extended to the Zevatron models to also explain the origin of the post GZK flux. These models are also called ”bottom-up” scenarios in contrast to a solutions which avoid acceleration by assuming particles are created with high energy already, the so called ”top-down” models.

A first general approach to these models was made by Thomas Gaisser gaisser  on a energy balance argument. The power density needed to be injected in cosmic rays to produce the higher energy region of the cosmic ray spectrum can be computed for a given type of object with a known distribution and assuming a nominal spectral index of $`\gamma =2`$ characteristic of Fermi acceleration, mimicking the well known argument supporting Fermi acceleration in Super Nova Remnants for the origin of the bulk of cosmic rays. It turns out that for galaxies, clusters of galaxies, Active Galactic Nuclei (AGN) and Gamma Ray Bursts (GRB) the UHECR power requirement is a reasonably small fraction of the power density emitted by each of these classes of objects. The argument can be reconsidered for strong relativistic shocks such as those expected to be found in AGN’s and GRB’s, which typically result in steeper spectral indices $`\gamma 2.22.3`$. In these shocks the particle can be accelerated provided it achieves a minimum energy by some other means. The resulting power balance is very dependent on both the spectral index and the injection energy. Part of the extra power that would be required because of the steeper spectrum is compensated by an increased injection energy. Simple energetic considerations open up these possibilities and do not allow much discrimination between them gaisser .

Fermi Acceleration in Gamma Ray Bursts has been suggested as a possible UHECR source waxman . Although the energetics may be adequate, the implied assumption that the bursts have a comoving density which is independent of redshift is somewhat contrived stecker . We should probably have to wait for a better understanding of these interesting phenomena before this possibility can be critically revised. In this respect we heard of an interesting proposal for GRB observations which could shed some light into the general GRB problem oscar .

Powerful radiogalaxies are certain candidates for post GZK acceleration. It was stressed that these objects provide the largest shock waves known and that the standard radioastronomical observations already demand highly energetic particles to explain the energy transport along the jets biermann . It has been claimed that a single source could provide the solution of the UHECR problem ahn . It is possible to devise reasonable magnetic field models in which backtracking of post GZK particles under the assumption of mostly proton primaries and one helium nucleus, leads to M87, a nearby radiogalaxy in the Virgo cluster ahn ; biermann . These hypotheses will be further tested once the statistics builds up. Very interesting results for magnetic field flux magnifications and spectral deformations in general and for this particular model were presented harari . These local magnifications could significantly affect power requirement estimates and the arguments based on them.

The shape of the spectrum at the cutoff is a signature of the source distribution and an indirect handle on the distance to the production sites. If the post GZK particles are of extragalactic origin but they are produced at higher rates in a region close to us relative to the 50 Mpc proton horizon, the GZK cutoff effect is mitigated berezinskii ; bahcall . A recent simulation of this effect was discussed at the conference blasi  assuming the post GZK flux is correlated with the dark matter distribution. The spectral index assumed for the injection spectrum becomes important and the observed flux above $`10^{20}`$eV would require a hard spectral index for injection with $`\gamma 2`$. A distinct signature for this scenario is the anisotropy of arrival directions. Certainly if local density enhancements are responsible for the post GZK particles the arrival directions should map these density enhancements. Large magnetic fields would be needed to explain observed distribution of arrival directions but this is consistent with what we know about magnetic fields biermann . Some claims for an excess in the direction of the supergalactic plane have been made stanev , but the level of confidence is still low because the statistics is poor.

A recent model for acceleration of iron like nuclei in the galaxy was addressed olintons . These nuclei are stripped off the neutron star surface with high magnetic fields that allow acceleration to the observed high energies. As the magnetic flux of a rotating neutron star reaches the light cylinder it can be converted into kinetic energy of the particles in a relativistic wind. If the magnetic field is high and the neutron star has a rapid rotation the acceleration can reach GZK energies. This alternative demands large galactic magnetic fields to reproduce the small enhancement in the direction of the galactic plane. The enhancement should increase as the energy of the observed particles rises. A particularly distinctive characteristic is the spectral index at production which is $`\gamma 1`$ olinto . Such a hard spectral index puts most of the energy in the high energy end of the spectrum and new experiments extending the energies to the ZeV region should be able to measure the spectral index.

### Alternative Solutions

Most of the solutions avoiding Fermi acceleration in nearby objects are motivated by physics beyond the standard model and would correspond to new physics in Angela Olintos’s review olinto . As regards the observations from the expected fluxes two very different categories emerge. In one class of alternatives the post GZK flux involves standard model particles from the decay and fragmentation of other particles. A second class avoids the cutoff mechanism either with new particles or with standard model particles that have unexpected behaviours.

#### Fragmentation Origin

Most of the alternative models share the feature of producing the bulk of the post GZK particle flux by fragmentation into pions. These models span very different scenarios, for instance the pions can be due to quarks which are in turn the decay products of a more massive particle such as a light electroweak $`Z`$ boson weiler , a Wimpzilla wimpzilla  which is a non thermal long lived particle, or a much heavier $`X`$ particle hillshramm  from a possible extension of the Standard Model into a Grand Unification Theory (GUT). These decaying particles are conjectured to be produced in a variety of mechanisms including, primordial origin, couplings to gravity, annihilations of topological defects and local interactions of distant ultra high energy neutrinos, some of which were addressed by several speakers olinto ; sigl ; weiler ; ellis ; masperi . The produced massive particles usually decay into quarks which fragment into hadrons, mostly pions that in turn decay into photons and neutrinos. There are models however in which the quarks can be emitted directly, for instance in primordial black hole evaporation macgibbon . As a final result the fragmentation and decay chain eventually ends up as photons, neutrinos and protons and no heavy nuclei can be expected. Such different scenarios ultimately predict roughly similar relative rates for the different particles because they rely on the same production mechanism.

Annihilation of topological defects were historically one of the first mechanisms postulated as a source of massive particles that give rise to post GZK particles as fragmentation products hillshramm . Under the heading of topological defect a large variety of objects can be included such as monopoles, cosmic strings bhatacharjee , vortons which are classically stable loops of superconducting cosmic strings masperi ; masperi2 , and combinations of these such as necklaces that combine strings and monopoles vilenkin , or monopole-antimonopole pairs connected by a string blanco  just to name a few. A number of reviews have dealt with these objects sigl ; siglbhat . The possibilities are already heavily constrained by various arguments and observations, such as the expected abundances, their effects on Big Bang Nucleosynthesis, and the expected secondary fluxes of photons and neutrinos sigl ; siglbhat .

A particularly important feature of these alternative sources concerns their cosmological distributions. The possibility that the post GZK particles are locally produced in our galactic halo following the fate of cold matter is recently gaining more attention blasi ; berezinskycluster ; tanco . If the clustering is on scales smaller than the 50 Mpc scale the GZK cutoff effect on the spectrum will be mitigated. In many of the proposed models involving fragmentation the sources are expected to be clustered. An interesting signature of the clustering sources is the expected dipole asymmetry associated to our position not being central in the halo and clearly related to the halo size.

##### Clustered sources:

One possibility is that the post GZK particles are due to decays of heavy metastable particles that cluster in our galactic halo following the fate of cold matter and eventually decay. They are mostly motivated by and different aspects of string and M-theories. Some of these were addressed in the conference for instance pentons motivated by M-Theory compactifications and cryptons that are supermassive bound states from string hidden sectors ellis .

Another model predicting a locally enhanced cosmic ray density is motivated by the recent findings about neutrino oscillations, and gives qualitatively similar predictions. The model was discussed in several talks and involves the production of $`Z`$-bosons in neutrino-antineutrino annihilations as discussed by several speakers weiler ; sigl . The model requires a very high energy neutrino flux that could originate far beyond the proton horizon of 50 Mpc and a target which is provided by massive relic neutrinos that could cluster in our galactic halo. No assumptions are made on the origin of the neutrino beam but some suggestions have been made postGZKnu . The maximum neutrino energy has to be pushed even further firstly because the observed particles are secondary products of the neutrino interactions and secondly because the neutrinos themselves can also be expected to come from decays of other particles. This model is subject to important restrictions because of neutrino flux bounds mainly from the absence of observations of horizontal or upward going air showers BlancoPRD .

All these models have become less likely in the view of very recent studies of composition for energies above $`10^{19}`$eV, that put the first limits on photon abundance zas ; zasbound . The restriction is strongest for models in which the sources are clustered in the vicinity of our galaxy so that the photon flux does not travel long enough distances to become relatively suppressed with respect to the proton flux.

#### Solutions without Fragmentation

A number of solutions to the Ultra High Energy puzzle has been suggested invoking particle behaviours that avoid the GZK cutoff. Some are familiar standard model particles that develop unexpected behaviours at large energies. They include hadronic-like cross sections for neutrinos that can mimic proton interactions in the atmosphere bordes ; ralston . Alternatively it has also been suggested that Lorentz invariance could be broken mestres ; glashow . Other solutions invoke new stable particles usually motivated by supersymmetry, also called uhecrons which have higher threshold energies in their interaction with the background photon fields but but do interact with matter farrar ; ellis .

Particular attention was paid to the recent suggestion that Lorentz invariance could be broken. At these extreme energies space-time can have a non trivial structure because of quantum fluctuations due to the recoil of the vacuum. High energy particles see distortions of the metric, a spacetime foam in which their speed becomes reduced. The Lorentz invariance break implied prevents particle production and particularly the photoproduction process that is mostly responsible for the GZK cutoff ellis . Most interesting where the first experimental constraints on these ideas searching for correlations of delays and distances in the radiation received form GRB’s in different energy bins. The observations imply that the possible quantum-gravity scale $`M`$ has a lower bound $`M>10^{15}`$GeV. Lorentz symmetry violations were also suggested from a CP-violating kinematic structure ahluwalia .

Lastly there are solutions involving other particles such as heavy or light monopolia kephart ; weiler ; swain  or dust grains. In the monopole models the monopoles themselves would induce the observed air showers, which have however a distribution of arrival directions which strongly disfavours such a hypothesis escobar . The dust grain alternative is ruled out by the shower development curve observed for the highest energy event kolb  but there is no knowledge of the behaviour of monopole induced showers weiler .

## Composition

Composition was discussed by several speakers as a fundamental tool to distinguish between different models which is crucial for the correct interpretation of the distribution of arrival directions which is mass dependent.

Particle production in models relying on decays of more massive particles is governed by the fragmentation processes. If the bulk of the post GZK particles are due to fragmented hadrons from heavier particle decays, the relative fluxes of different particle species can be extrapolated from lower energy fragmentation processes that are well known from $`e^+e^{}`$ collisions in accelerators. In fragmentation processes pions are typically produced at a rate which is of order ten times higher than protons. Photons arise from neutral pion decays while neutrinos are mostly produced in the decays of charged pions. The fragmentation and decay chain process is thus expected to produce a significantly larger fraction of photons and neutrinos than protons (neutrons decay and end up as protons) and no heavy nuclei. 333Incidentally it was pointed out during this conference that these extrapolations are not free of uncertainty and that they should be carefully reconsidered for instance in connection to supersymmetric extensions of the standard model sigl ; blasi ; ellis .

One must note that measurements of these fluxes have already proven to be very useful in constraining production models sigl ; gaisser . These neutral particles if observed will provide non deflected information on the fluxes, thus directly reflecting the source distribution anisotropy independently on assumptions about the magnetic fields which would be a most valuable piece of information.

As a result one can roughly expect flux ratios of order ten for photons and neutrinos with respect to protons. The photon ratio is further modulated by interactions with the background radiation fields. If the model produces more particles near the Galaxy, the photon fraction is of order 10 ellis ; sigl ; olinto ; selvon  while this value can be modified depending on distance to the source and fragmentation details to drop to values about 1-3 when the average distance travelled by the particles is in the scale of tens of Mpc kachelriess ; sarkar ; sigl . On the other hand the neutrino flux will not be attenuated and the neutrino to proton flux ratio can be expected to be of order ten or higher if there is proton attenuation.

One must remark that the predictions for the secondary fluxes are rather sensitive to proton, electron and photon absorption both at the production site and during transport. If photons are attenuated, pair production and synchrotron losses dump the photon energy density into the MeV region of the gamma ray spectrum, linking it to the unknown magnetic fields protstan ; sigl .

Alternatively stochastic acceleration mechanisms require charged particles. If the accelerated particles are heavy nuclei higher energies can be achieved. The possibility that these ultraenergetic particles are heavy nuclei is attractive on a double basis because these nuclei are more easily isotropized in a given magnetic field and because acceleration models are less constrained by Eq. 1. Nuclei heavier than iron such as gold have been also proposed as a plausible solution swain . Although the mean free path for interactions with the background radiation fields is small, the emerging depleted nuclei can undergo successive interactions and reach distances in the 100 Mpc range keeping a substantial fraction of the original energy stecker . Heavy nuclei thus relax the constrains on both source distance and intensity of the magnetic fields that observations impose.

If the bulk of the post GZK particles are charged, the study of the arrival directions will allow improvements of the current bounds on the poorly known magnetic fields outside the galaxy. These bounds will be stricter if composition can be determined. There are models with heavy nuclei involving both galactic and extragalactic sources for which magnetic effects should be most important. They will be tested against future data on both the energy spectrum and anisotropy measurements.

Acceleration models are also expected to produce secondary fluxes of photons and neutrinos. The mechanism is the interaction of the protons with matter or radiation and this can happen both at the source or during transport. Depending on details of the models themselves the secondary fluxes will have however different magnitudes gaisser ; sigl . Even if no interactions take place at the source itself the GZK cutoff mechanism should produce a minimum flux of secondaries in these models. In practically all cases however the relative fluxes of neutrinos and photons at energies above the GZK cutoff is expected to be less than that of protons, unless there is a high attenuation of the protons themselves at origin.

In models relying on fragmentation, photons are expected to present a significant fraction of photons in the post GZK flux. In this respect the measurements of photon composition will be crucial because they provide an indirect handle on the clustering process for the post GZK sources. Neutrino searches however will extend the observation horizon to distances well beyond the 100 Mpc range and will confirm clustering scenarios. Interesting new ideas were reported for establishing the photon composition of the post GZK flux which will be further discussed below. These ideas will certainly prove extremely valuable in selecting between acceleration and fragmentation mechanisms.

## EXPERIMENT: Exploring the Post GZK Flux

Altogether only seventeen events with energies exceeding $`10^{20}`$eV have been reported. Sixteen of them are from four array experiments, one from Volcano Ranch volcano , six from Haverah Park HP , one from Yakutsk array yakutsk  and eight from AGASA AGASA . The other event, which is the highest energy event observed, was detected with Fly’s Eye bird , a different detector concept that measures the fluorescence light from nitrogen as the shower develops in the atmosphere. The accumulation of events as a function of time is shown in Fig. 1 including an ultrahigh energy event detected in 1999 by AGASA lateagasa , and two events in a new analysis of the Haverah Park data for zenith angles above $`60^{}`$ zasrate  which have been discussed in this conference.

Although the statistics are now enough to convince the few remaining skeptics, the field clearly demands more data. On the one hand more and more precise data on the flux measurements themselves is needed to build up statistics for the spectral and anisotropy studies. This need is the drive for a number of projects that were presented in this conference. Some are already in construction but many presentations were concerned about experiments in planning. The importance of additional information in the search of the origin of the post GZK particles has already been stressed. Experiments and techniques that are sensitive to primary composition and to secondary fluxes of gamma rays and neutrinos played a central part in the conference.

### Existing Detectors

We heard reports on data obtained with three cosmic ray detectors: AGASA, Haverah Park and HiRes.

##### AGASA:

The status report of this ongoing experiment teshima ; lateagasa  was centered on the 8 events detected above $`10^{20}`$ eV and on the anisotropy of arrival directions. The search of coincidences in arrival directions within the angular resolution of the experiment for events with energies above $`410^{19}`$eV results in one 3-fold coincidence and four 2-fold coincidences with an estimated $`0.3\%`$ chance probability. The recently reported $`10\%`$ anisotropy observed in the region of the galactic center and anti-center for energies above $`10^{18}`$eV hayashida  was also addressed. Although this effect suggests a galactic origin, it is observed at a threshold energy well below the GZK cutoff and may not have much to do with the origin of the UHECR if as expected they indeed are a different component.

##### Haverah Park:

The interest of Horizontal Air Showers (HAS) was shown not to be limited to neutrino detection. We also heard a report on the first analysis of inclined HAS induced by UHECR zas . A recently developed model for muon density maps at ground level zasmodel  allows event reconstruction. As a first application, the analysis of the Haverah Park data has resulted in 7 (2) new events above $`410^{19}`$eV ($`10^{20}`$eV) zasrate .

##### HiRes:

HiRes hires  sokolsky  is the first large aperture detector that has recently started operation based on the nitrogen fluorescence technique. The detector has two wide aperture eyes and when operating in the stereo mode it is expected to measure the energy of cosmic rays between 0.1 and 200 EeV with a $`20\%`$ resolution, with improved sensitivity to the position of shower maximum and to the determination of the arrival directions. About five more events with energy above $`10^{20}`$ eV can be expected from the preliminary results of the first six months of data hires . The forthcoming data analysis of HiRes data is eagerly expected because there is only one calorimetric measurement of an event with energy above $`10^{20}`$ eV. It will also help to resolve small discrepancies between different experiments.

### Detectors in Planning

Most advanced plans are for ground detectors using well established technologies. The Auger observatory will be the first one and the Telescope Array and IceCube may be the following ones, but new projects involving detection from satellites are already quite advanced in planning and could start operating in the second half of the 2000’s decade.

#### Ground Experiments:

##### Auger:

The southern Auger observatory is now in construction in Pampa Amarilla, Mendoza, Argentina and will be the next large aperture detector to become operative cronin . The engineering array is expected to be finished in less than a year time. Detailed progress on tank monitoring and calibration salazar , tyvek reflectivity zepeda , photomultipliers fernandez  and on Fresnel mirrors for the fluorescence cordero  was reported. The analysis of HAS induced by cosmic rays zas  has shown that the commonly estimated acceptance for cosmic rays with zenith angles below $`45^{}`$ can be doubled for the highest energy end of the spectrum. As it is a hybrid detector designed to combine the fluorescence and the particle density sampling techniques, and it is the first one that will perform such hybrid measurements for showers with energy exceeding $`10^{19}`$eV, it will be most adequate for cross calibration of the two different techniques and will provide key information to the development of the field.

##### Telescope Array:

The Telescope Array (T.A.) is a project consisting on an array of eight stations each with 42 3-m diameter mirrors covering the solid angle from the horizon to $`58^{}`$ of zenith angle and having a similar aperture to the Auger observatories. It will have $`6\%`$ energy resolution and angular resolution of $`0.6^{}`$ teshima . Decisions on its funding are expected within a year from now.

##### IceCube:

The IceCube proposal gaisser ; icecube  is a well advanced project aiming to expand the proven AMANDA amanda  technology to a full scale detector instrumenting one km3 of deep clear ice. The expected rates of neutrinos depend both on the production mechanism and on the attenuation of cosmic rays during transport. The mechanism determines the ratio of neutrinos to protons and photons. Depending on the maximum energy of the neutrinos the flux normalization can span a wide range of values in plausible scenarios. The expected rates in IceCube in many of these scenarios are at the detectable rate gaisser .

#### Detectors from Space:

Observation of the fluorescence light induced by extensive air showers from a satellite is a promising alternative. The thin upper atmosphere allows distant detection and a large aperture can be achieved by placing the satellite at the planned orbit of the International Space Station (ISS) of about 400 km. Several initiatives were discussed in a number of presentations. These systems are sensitive both to UHECR and to UHE neutrinos by looking for deeply penetrating air showers. The Čerenkov light reflected from the clouds, earth or the sea can be used for distance calibration scarsi ; khrenov .

##### EUSO:

The Extreme Universe Space Observatory (EUSO) scarsi ; catalano ; cordero  is a multinational joint effort which has been approved by the European Space Agency for an accommodation study. It is planned to have a field of view of $`60^{}`$ spanning $`200\times 200`$ km2, Fresnel lens optics catalano ; cordero  and a finely segmented focal plane that will allow the registration of the fluorescence light from air showers. It could start operation in 2006 and is scheduled to operate for three years.

##### KLYPVE:

Another version of the same technique adapted to the Russian section of the ISS is provided by the KLYPVE project khrenov ; panasyuk ; garipov . It has a smaller field of view of $`15^{}`$ that will cover an area of $`100\times 100`$ km2 panasyuk , 10 m Fresnel mirrors optics and 5 mrad pixels khrenov .

This class of experiments is based on a well understood technique so they can be expected to be viable but it relies on significant technological developments. The technique has not yet been proven and in that respect the TUS (Track Ultraviolet from Space) project is a preparatory smaller Russian project using 1.5 m mirrors. It is expected to be hosted by the Resource-DK satellite to be launched in 2002-2003 panasyuk . Upgrades of these detectors were also discussed such as MULTIOWL and SUPEROWL as part of a Grand Observatory in Space takahashi . They intend to enhance the performance of EUSO or KLYPVE by tilting the system with respect to nadir, having a higher altitude orbit and using more mirrors to measure ZeV ($`E>10^{21}`$eV) cosmic rays and neutrinos.

### New Ideas for Composition

As was pointed out by many speakers composition is a fundamental clue in the search of the origin of the post GZK flux olinto ; weiler ; sigl . The unavoidable need of using indirect detection techniques makes it hard to accurately determine the nature of the particles themselves because the differences in atmospheric showers developed by different primaries are rather subtle. Most of them can be related to either the shower development or to the muon content in the showers which can be explored by different techniques. Intresting new ideas were discussed for the study of the nature of the post-GZK particles.

##### Neutrinos:

All the detectors for measuring the post GZK flux in construction or planning are sensitive to neutrinos by searching for deeply penetrating showers croninzas ; athar  and much use has been made of this technique to constrain high energy neutrino fluxes BlancoPRL . An interesting new idea for the detection of tau decays in the southern Auger observatory selvon . Tau neutrinos can be expected in flavour mixing scenarios which are motivated by atmospheric neutrino measurements. The taus are produced in charged current tau neutrino interactions with the Andes mountains. For neutrinos in the energy range $`10^{17}\mathrm{5\hspace{0.17em}10}^{18}`$ eV the tau decays can give a detectable signal in the Auger tanks in spite of the solid angle being small. It was pointed out that it may not be possible to search for similar effects with fluorescence satellite detectors because their energy threshold would prevent the tau from decaying in the detector field of view selvon .

##### Photons:

Two new ideas were discussed for the identification of a photon component in the post GZK particle flux:

The first one is due to a competition between the Landau-Pomeranchuck-Migdal (LPM) suppression of pair production and bremsstrahlung because of collective effects of the atmospheric nuclei and the interaction of photons with the geomagnetic field of the Earth. The LPM effect changes shower development in a dramatic fashion above an effective threshold. The shower develops much later becoming more elongated and a ground particle array would detect a much younger shower with a characteristic steep lateral profile. On the other hand the interaction of a photon with the Earth magnetic field can happen long before the start of shower development and has the effect of distributing the photon energy among lower energy secondary photons, electrons and positrons that add up to the primary photon energy. Provided the photon does not travel parallel to the magnetic field this interaction prevents the development of LPM showers and the atmospheric shower that follows is a scaled up version of the well established lower energy showers. There are however preferred directions along which the photon is unlikely to interact with the magnetic field of the Earth and the shower develops these strange LPM behaviour that have been carefully studied for the southern Auger observatory bertou . The alignment of the magnetic fields gives “holes” for LPM showers. This possibility was discussed as a means to establish the photon fraction in cosmic rays for both particle arrays teshima  and for the fluorescence technique scarsi .

A second idea was discussed in connection to the detection of horizontal air showers (HAS). Inclined showers induced by protons, nuclei or photons develop early in the atmosphere and only their muon component can be detected with particle arrays. The recent study of geomagnetic effects of the muon distributions at ground level zasmodel  allows the study of events above $`60^{}`$. The expected HAS rates have been shown to be very dependent on composition because of the relative content of muons in the showers induced by different primaries, particularly between photons and protons or nuclei. As a result the combined measurement of vertical and horizontal air shower rates provides an new handle on composition. These ideas have been applied to Haverah Park data which are inconsistent with a model in which photons are more than $`42\%`$ of the cosmic ray flux above $`10^{19}`$ eV, which provides a severe constrain on models which have locally enhanced cosmic ray production by fragmentation zasbound ; zas . Since the number of muons in a photon shower is expected to be smaller by a large factor with respect to a proton (or nucleus) of the same total energy the results are fairly robust. These ideas can also be used for setting bounds on the heavy ion fraction, but the results are more dependent on the muon production of the assumed interaction model.

### The Intermediate Energy Range

Some of the talks referred to experiments designed to measure lower energy particles. Low energy fluxes can be related to UHECR although the bulk of the radiation received in the low energy region is expected to be from a different origin from the UHECR gaisser .

##### AMS:

The Anti Matter Spectrometer (AMS) also discussed in the conference menhaca . It is an advanced detector to be flown in the ISS in 2003 for precision measurements of cosmic ray composition in the GeV energy range that has an important mexican involvement AMSmex . A successful prototype has already been flown.

##### Ground Based GRB watch:

Gamma Ray Bursts may be the source of UHECR. A review was presented of the unique Chacaltaya observatory located 5230 m above sea level at an atmospheric depth of only 540 g cm-2. It included a recent proposal to transport an extensive air shower array combining scintillator detectors and air Čerenkov collectors to be used in combination with a 100 m2 calorimeter oscar . The arrangement could allow measurements of inelasticity by identifying the low energy subshowers of quasielastic events with Čerenkov collectors on the surface and could have a sensitivity to gamma rays of 16 GeV from Gamma Ray Bursts. An interesting proposal was made for a gamma ray watch in the 1-100 GeV range combining different ground experiments such as Chacaltaya, MILAGRO and the Auger Observatory. These experiments would be operated in a low threshold energy mode using correlated excess in the single particle counts oscar .

##### Dark Matter Searches:

Searches for supersymmetric dark matter were reviewed by John Ellis ellis . These particles can be detected via their annihilation in the galactic halo producing antiparticles which could be detected by satellites such as AMS and other satellites searching for gamma rays. Alternatively this possibility can be explored by searching for excess neutrinos in the directions of the Sun or the Earth center produced by annihilations of the supersymmetric particles trapped in the corresponding gravitational potentials. This method is the most sensitive and many specific models can be ruled out with 10 km2 y of neutrino data. Lastly they can be searched for directly by looking for their elastic scattering in a detector and there is a recent unconfirmed claim of evidence for such scattering DAMA , but the final answer may come after LHC starts running in 2006.

## Epilogue: The Future of the Field

The frustratingly inconclusive statistics reached at present will cease to be a problem in the near future as it is illustrated in Fig. 2 where the number of events detected above $`10^{20}`$eV is plotted as a function of time. The statistics will be more than doubled in two years with forthcoming HiRes data. By the end of the year 2005 we should have two years of data from the Auger detector and more than ten times the statistics we now have. The rising statistics trend will be continued with even more ambitious projects such as EUSO, KLIPVE and/or the Telescope Array.

The new experiments will not only bring more accurate measurements of the post GZK spectrum, they will also bring about new handles on composition. All the experiments in planning or construction are also sensitive to high energy neutrinos by looking for horizontal air showers. In the EeV energy range the neutrino flux is expected to be related to the production mechanisms of the cosmic rays themselves. Such flux will be further explored by dedicated under-water or under-ice neutrino experiments in construction amanda  or development antares , quite possibly with the IceCube proposal to extend AMANDA to a km3 size detector and also by searching for horizontal air showers with all the projects that are sensitive to ZeV Air Showers.

The coming decade of the 2000’s will bring the eagerly expected data that will answer many of the open questions of the post GZK flux. All these experiments must be regarded as a whole effort towards the understanding of the highest energy particles in the Universe. Different techniques for the detection of Ultra High Energy particles will complement each other, resolve uncertainties which are still present and will open the way to many different and interesting aspects of high energy particle physics, astrophysics and cosmology.

Acknowledgements: I thank the organizers of this very fruitful conference in a fabulous setting under the spell of an impressive volcano. I also want to thank J.J. Blanco-Pillado, G. Parente, R.A. Vázquez and A.A. Watson for many comments and suggestions after reading this manuscript. This work was supported in part by CICYT (AEN99-0589-C02-02) and by Xunta de Galicia (PGIDT00PXI20615PR).


text:
# Electronic response of aligned multishell carbon nanotubes

\[

## Abstract

We report calculations of the effective electronic response of aligned multishell carbon nanotubes. A local graphite-like dielectric tensor is assigned to every point of the multishell tubules, and the effective transverse dielectric function of the composite is computed by solving Maxwell’s equations. Calculations of both real and imaginary parts of the effective dielectric function are presented, for various values of the filling fraction and the ratio of the internal and external radii of hollow tubules. Our full calculations indicate that the experimentally measured macroscopic dielectric function of carbon nanotube materials is the result of a strong electromagnetic coupling between the tubes, which cannot be accounted for with the use of simplified effective medium theories. The presence of surface plasmons is investigated, and both optical absorption cross sections and energy-loss spectra of aligned tubules are calculated.

\]

The so-called multishell carbon nanotubes, with diameters on the order of $`10\mathrm{nm}`$ and lengths reaching up to several microns, are needle-like tubes consisting of a finite number of concentric cylindrical two-dimensional layers of graphite sheet arranged around a central hollow with a constant separation between the layers. Since their discovery, these tubules have attracted much attention because of their size and their novel structural and electronic properties. An important step in the investigation of these new materials was achieved when composites were manufactured with nanotubes organized into well-aligned arrays, since alignment of nanotubes was expected to be important not only to test properties but also for potential applications such as atomic-scale field emitters. With the availability of aligned carbon nanotube films, optical measurements were carried out with polarized light, thereby evaluating the frequency-dependent effective dielectric function of the composite and showing that carbon nanotubes have an intrinsic and anisotropic metallic behaviour. Accurate calculations of the effective dielectric function and the optical response of densely packed carbon nanotubes, as obtained by solving Maxwell’s equations with the use of tensor-like dielectric functions, have been carried out only very recently.

In this paper, we extend the calculations reported in Ref. to the case of hollow-core nanotubes, and emphasize the effects of the inner radius on the effective dielectric function of the nanotube array. Calculations of the so-called energy-loss function, i.e., the imaginary part of the effective inverse dielectric function, are also presented, for various values of the filling fraction and the ratio of the internal and external radii of the tubules.

We take a periodic array of infinitely long multishell nanotubes of inner and outer radii $`r`$ and $`R`$, respectively, arranged in a square array with lattice constant $`a=2xR`$, as shown in Fig. 1. These tubules are assumed to be embedded in vacuum. In the energy range $`0`$ to $`10\mathrm{eV}`$ of interest in the interpretation of absorption cross sections and energy-loss spectra, the diameter of typical multishell carbon nanotubes ($`2R10\mathrm{nm}`$) is small in comparison to the wavelength of light. Hence, our results will be independent of the actual value of the outer radius of the cylinders. We also assume that this radius is large enough that a macroscopic dielectric function is ascribable to the coaxial cylinders. These cylinders occupy a fraction $`f=\pi (1\rho ^2)/(4x^2)`$ of space, with $`\rho =r/R`$ being the ratio between the inner and outer radii of the tubes. For simplicity, the magnetic permeabilities will be assumed to be equal to unity in all media.

Planar graphite is a highly anisotropic material, and the dielectric function is a tensor. This tensor may be diagonalized, by choosing cartesian coordinates with two of the axis lying in the basal plane and the third axis being the so-called $`c`$-axis. One defines the dielectric function $`\epsilon _{}(\omega )`$ perpendicular to the $`c`$-axis and the dielectric function $`\epsilon _{}(\omega )`$ for the electric field parallel to the $`c`$-axis. For carbon nanotubes, we assume full transferability of the dielectric tensor of planar graphite to the curved geometry of carbon tubules, as suggested by Lucas et al  for the case of multishell fullerenes. Hence, we simply assign a local graphite-like dielectric tensor to every point inside the nanotube and outside the inner core,

$$\widehat{\epsilon }(\omega )=\epsilon _{}(\omega )(\theta \theta +\mathrm{𝐳𝐳})+\epsilon _{}(\omega )\mathrm{𝐫𝐫},$$
(1)
$`\theta \theta `$, $`\mathrm{𝐳𝐳}`$, and $`\mathrm{𝐫𝐫}`$ being the unitary basis vectors of cylindrical coordinates, and take the principal dielectric functions $`\epsilon _{}(\omega )`$ and $`\epsilon _{}(\omega )`$ of graphite from Ref..

In the long-wavelength limit, a composite material may be treated as if it were homogeneous, with the use of an effective dielectric function $`\epsilon _{eff}`$. The optical absorption cross section of the composite is then directly given by $`\mathrm{Im}\epsilon _{eff}(\omega )`$. Also, for small values of the adimensional parameter $`qR`$ ($`qR<1`$), $`q`$ being the momentum transfer, the energy-loss spectra of a broad beam of swift electrons penetrating the composite is found to be well described by the $`q0`$ limit of the imaginary part of the effective longitudinal dielectric function, i.e., the so-called energy-loss function, $`\mathrm{Im}[\epsilon _{eff}^1(\omega )]`$.

In order to compute, with full inclusion of the electromagnetic interaction between the tubules, the effective dielectric function of our periodic system, we first solve Maxwell’s equations. The local dielectric function of our composite material depends on frequency. Thus, on-shell methods for solving Maxwell’s equations, which proceed at constant frequency, are particularly well suited for our calculations. Though the methodology of Ref. was originally reported to solve Maxwell’s equations with scalar dielectric functions, an extension of this formalism was developed in Ref., which enables us to consider tensor-like dielectric functions like the one of Eq. (1).

We first fix the frequency, thereby the local tensor-like dielectric function at any point of the composite material being specified, and approximate the continuous electromagnetic field by its values at a series of mesh points located on a simple square lattice with unit cell of side $`a`$, as shown in Fig. 1. The corresponding discretized Maxwell’s equations provide a relationship between the electromagnetic fields on either side of the unit cell, i.e., $`x`$ and $`x+a`$ (see Fig. 1). We exploit this to calculate, on applying Bloch’s theorem, the eigenvalues of the so-called transfer matrix, which will give the band structure of the system, i.e., the dispersion relations $`k(\omega )`$ of the Bloch-state wave vector of propagation versus the frequency. The effective transverse dielectric function is then obtained as follows,

$$\epsilon _{eff}(\omega )=\frac{k^2(\omega )c^2}{\omega ^2},$$
(2)
where $`c`$ represents the speed of light.

Fo an e.m. wave normally incident on the structure ($`k_y=k_z=0`$) there are two different values of $`\epsilon _{eff}(\omega )`$ corresponding to $`s`$ and $`p`$ polarizations. For e.m. waves polarized along the tubes ($`s`$ polarization), the electric field is parallel to the cylinders at every point, and is not modified by the presence of the interfaces. Hence, the effective dielectric function of the composite is simply the weighted average of the dielectric functions of the constituents,

$$\epsilon _{eff}(\omega )=f\epsilon _{}(\omega )+(1f).$$
(3)
Our numerical results for the effective dielectric function, as obtained from Eq. (2) for this polarization, accurately reproduce the exact results predicted by Eq. (3), which represents a good check for our scheme. According to Eq. (3), $`\mathrm{Im}\epsilon _{eff}(\omega )=f\mathrm{Im}\epsilon _{}(\omega )`$; therefore, the shape of the experimentally determined $`\mathrm{Im}\epsilon _{eff}(\omega )`$ must coincide for this polarization with $`\mathrm{Im}\epsilon _{}(\omega )`$. The measurements of $`\epsilon _{eff}(\omega )`$ reported in Ref. for light polarized along the tubes are plotted in Fig. 2 by solid circles, together with the effective dielectric function predicted by Eq. (3), which roughly reproduces with $`f0.5`$ the experimentally determined $`\mathrm{Im}\epsilon _{eff}(\omega )`$. We note, however, that neither the exact position of the maximum nor the details of the real part of the measured effective dielectric function can be reproduced with the use of Eq. (3). This may reflect existing differences between carbon structures in planar graphite and in nanotubes. Also, optical measurements may be affected by the presence of tubes with different orientations.

The effective energy-loss function, $`\mathrm{Im}[\epsilon _{eff}^1(\omega )]`$, corresponding to $`s`$ polarization is also plotted in Fig. 2. For this polarization, the effective inverse dielectric function has a single pole, which occurs at a reduced plasma frequency where $`\mathrm{Re}\epsilon _{}(1f)/f`$ and $`\mathrm{Im}\epsilon _{}0`$. For graphite, at low concentrations of tubules ($`f<0.2`$) the resonance condition is never satisfied, and in the limit as $`f0`$ the effective energy-loss function simply coincides with $`\mathrm{Im}\epsilon _{}(\omega )`$, thereby showing the characteristic peak at $`4.6\mathrm{eV}`$. This peak is associated with the maximum in the joint density of states of the $`\pi `$ valence and conduction bands in graphite.

Calculations of the $`p`$ component of the effective dielectric function of an array of plain ($`\rho =0`$) carbon nanotubes were reported in Ref. for various values of the ratio $`x`$ between the lattice constant and the outer diameter of the cylinders, showing that the trend with increasing the concentration of tubules is for the actual dipolar peak in $`\mathrm{Im}\epsilon _{eff}(\omega )`$ to be shifted from the isolated-cylinder dipole mode at $`6.5\mathrm{eV}`$ to lower energies. When the nanoparticles are brought into close contact, electromagnetic coupling between them converts the dipolar surface mode into a very localized one, trapped in the region between the nanostructures. At higher concentrations of tubules, when graphite forms a connected medium ($`x1`$), dipolar modes cannot be excited, and the optical absorption exhibits a single peak originated in the maximum of $`\mathrm{Im}\epsilon (\omega )`$ at $`4.6\mathrm{eV}`$. We have also calculated the energy-loss function of an array of plain carbon nanotubes, and have found that for $`x1`$ it shows a single peak at the bulk plasmon resonance at $`7\mathrm{eV}`$, where $`\mathrm{Re}\epsilon _{}(\omega )`$ and $`\mathrm{Im}\epsilon _{}(\omega )`$ are both small.

For hollow tubes ($`\rho 0`$), there are two distinct dipolar modes with either tangential or radial symmetry, similar to those present in the case of a thin planar film and a spherical shell. With the presence of anisotropy this two-mode structure is replaced by a more complicated spectral representation. Our full calculations of the $`p`$ effective dielectric function of a periodic array of hollow carbon nanotubes, as obtained for various ratios between inner and outer radii, are shown in Figs. 3 and 4 with $`x=2.0`$ and $`x=1.03`$, respectively. As in the case of plain cylinders ($`\rho =0`$), our full calculations nearly coincide for $`x=2.0`$ with a generalized Maxwell-Garnett (MG) effective dielectric function appropriate for anisotropic hollow tubes. However, in the close-packed regime ($`x=1.03`$) the strong electromagnetic coupling between the tubes yields non-negligible contributions from multipolar resonances. This multipolar coupling provokes a redshift and a blueshift of the MG dipolar resonances that are visible in the optical spectra and the energy loss, respectively. These are the low-energy dipolar mode with tangential symmetry, which absorbs light ($`\mathrm{Im}\epsilon _{eff}(\omega )`$ is maximum at the corresponding energy), and the high-energy dipolar mode with radial symmetry, which can be excited by moving charged particles ($`\mathrm{Im}[\epsilon _{eff}^1(\omega )]`$ is maximum at the corresponding energy).

The experimentally determined macroscopic dielectric function of close-packed carbon nanotubes, as reported in Ref. for $`p`$-polarized light, is exhibited in Fig. 4 by solid circles. A comparison between these measurements and the $`p`$ effective dielectric function of closed-packed ($`x=1.03`$) plain ($`\rho =0`$) carbon nanotubes was presented in Ref., showing a nice agreement between theory and experiment, though a shift of approximately 1 was needed to obtain agreement with the real part of $`\epsilon _{eff}(\omega )`$. If one chooses the ratio $`\rho `$ between the inner and outer radii of the tubes to be 0.6, which yields in the close-packed regime ($`x=1.03`$) a filling fraction $`f0.5`$, both real and imaginary parts of our calculated effective dielectric function show an excellent agreement with the experiment. This filling fraction ($`f0.5`$) is precisely what is required to also reproduce the dielectric function experimentally determined with the use of light polarized along the cylinders, as discussed above (see Fig. 2).

In conclusion, we have shown that the effective dielectric function of densely packed carbon nanotubes is the result of a strong electromagnetic coupling between the tubes. We have chosen the inner and outer diameters of hollow tubules to be $`6`$ and $`10\mathrm{nm}`$, respectively ($`\rho =0.6`$), and have obtained both $`s`$ and $`p`$ effective dielectric functions to be in the close-packed regime ($`a=10.3\mathrm{nm}`$, $`x=1.03`$) in excellent agreement with the experimental measurements. Small discrepancies have been observed in the case of light polarized along the tubes ($`s`$ polarization), which may reflect existing differences between the dielectric function of carbon structures in planar and curved graphite.

J.M.P. gratefully acknowledges partial support by the University of the Basque Country, the Basque Hezkuntza, Unibertsitate eta Ikerketa Saila, and the Spanish Ministerio de Educación y Cultura.


text:
# Erratum Is it possible to infer the equation of state of a mixture of hard discs from that of the one-component system?

((Molec. Phys., 1999, 96, 1185–1188))

The numerical values in the sixth and seventh columns of table 1 of the paper are not correct. Consequently, the comments made in the paper about the better performance of the equation of state (12) proposed in  over the equation of state (6) proposed in the paper are wrong. The corrected version of table 1 is reprinted here. In view of this table, it is quite apparent that the results obtained from equation (6) are practically indistinguishable from those obtained from equation (12). For the cases considered in the table, both equations of state differ by less than 0.05% when combined with equation (14) for the one-component system and by less than 0.02% when combined with equation (15) for the one-component system. These differences increase, but are kept relatively small, when mixtures more asymmetric are considered. For instance, in the case $`\alpha =0.1`$, $`x_1=0.01`$ and $`\eta =0.6`$, the values of $`p\sigma _1^2/k_BT`$ predicted by (6)+(14), (6)+(15), (12)+(14) and (12)+(15) are 183.87, 186.26, 185.40 and 186.53, respectively. Table 1 also shows that, in general, the use of equation (15) gives a better agreement with the simulation data typically for $`\eta <0.55`$, while the use of equation (14) is preferable for $`\eta >0.55`$.


text:
# Multi-wavelength observations of the red giant R Doradus with the MAPPIT interferometer

## References

Hofmann, K. H., Scholz, M., & Wood, P. R., 1998, A&A, 339, 846.

Jacob, A. P., Bedding, T. R., Robertson, J. G., & Scholz, M., 2000a, MNRAS, 312, 733.

Jacob, A. P., Bedding, T. R., Robertson, J. G., Barton, J. R., Haniff, C. A., Marson, R. G., & Scholz, M., 2000b, In Lena, P. J., & Quirrenbach, A., editors, Interferometry in Optical Astronomy, Proc. SPIE vol. 4406, page 723 (astro-ph/0005195).