text:
# Scaling properties of granular materials

## Abstract

Given an assembly of viscoelastic spheres with certain material properties, we raise the question how the macroscopic properties of the assembly will change if all lengths of the system, i.e. radii, container size etc., are scaled by a constant. The result leads to a method to scale down experiments to lab-size.

Assume the dynamics of a certain granular system $`S`$ is known. Will the dynamics change if we rescale all sizes by a constant factor $`\alpha `$, i.e., $`R_i^{}\alpha R_i`$, but leaving the material properties unchanged? prime  If scaling affects the system properties, how do we have to modify the material properties to assure that the system $`S`$ and the scaled system $`S^{}`$ behave identically?

It is frequently desired to investigate large scale phenomena in granular systems experimentally, such as geophysical effects or industrial applications. To this end one has to rescale all lengths of the system to meet the restrictions of the laboratory size, i.e. big boulders in the original system become centimeter sized particles in the experiment. Of course, one wishes that the effects which occur in the large system occur equivalently in the scaled system too. With the assumption of viscoelastic particle deformation we will show that naive scaling will modify the properties of a granular system such that the original system and the scaled system might reveal quite different dynamic properties. To guarantee equivalent dynamical properties of the original and the scaled systems we have to modify the material properties in accordance with the scaling factor and we have to redefine the unit of time.

As another consequence of the scaling properties we claim that for numerical simulations of granular material it is not sufficient to provide relative data such as, e.g. to describe the container size in units of the particle diameter. We will show an example where the dynamics of a granular system changes significantly with system size, although all relative sizes are kept constant.

In a simple approximation, a granular system may be described as an assembly of spheres of radii $`R_i`$, $`i=1,\mathrm{},N`$. If two particles of radii $`R_i`$ and $`R_j`$ at positions $`\stackrel{}{r}_i`$ and $`\stackrel{}{r}_j`$ touch, i.e., if $`\xi _{ij}R_i+R_j\left|\stackrel{}{r}_i\stackrel{}{r}_j\right|>0`$, they feel an interaction force

$$\stackrel{}{F}_{ij}=F_{ij}^n\stackrel{}{n}_{ij}+F_{ij}^t\stackrel{}{t}_{ij},$$
(1)
with the unit vector in normal direction $`\stackrel{}{n}_{ij}\left(\stackrel{}{r}_j\stackrel{}{r}_i\right)/\left|\stackrel{}{r}_j\stackrel{}{r}_i\right|`$ and the respective unit vector in tangential direction $`\stackrel{}{t}_{ij}`$. Eventually, external forces as, e.g., gravity, may also act on the particles.

The normal force $`F^n`$ can be subdivided into elastic and dissipative parts $`F^n=F_{\mathrm{el}}^n+F_{\mathrm{dis}}^n`$ indices . The elastic force for colliding spheres is given by Hertz’s law Hertz

$$F_{\mathrm{el}}^n=\frac{2Y}{3\left(1\nu ^2\right)}\sqrt{R^{\mathrm{eff}}}\xi ^{3/2}\kappa \xi ^{3/2},$$
(2)
with $`R^{\mathrm{eff}}=R_iR_j/\left(R_i+R_j\right)`$ and $`Y`$, $`\nu `$ being the Young modulus and the Poisson ratio. Equation (2) also defines the prefactor $`\kappa `$ which we will need below.

The formulation of the dissipative part of the force $`F_{\mathrm{dis}}^n`$ depends on the mechanism of damping. Here we will focus on viscoelastic damping which is the most simple assumption for dissipatively colliding bodies brittle . It implies that the elastic part of the stress tensor is a linear function of the deformation tensor and the dissipative part of the stress tensor is a linear function of the deformation rate tensor. It is valid if the characteristic velocity (the impact rate $`g`$) is much smaller than the speed of sound $`c`$ in the material and the viscous relaxation time $`\tau _{\mathrm{vis}}`$ is much smaller than the duration of the collisions $`\tau _c`$ BSHP . The range of the viscoelastic model is, hence, limited from both sides: the collisions should not be too fast to assure $`gc`$, $`\tau _{\mathrm{vis}}\tau _c`$, and not too slow to avoid influences of surface effects as adhesion. For viscoelastically colliding spheres the dissipative part of the normal force reads Kuwabara ; BSHP ; Morgato

$$F_{\mathrm{dis}}^n=A\frac{d\xi }{dt}\frac{d}{d\xi }F_{\mathrm{el}}^n=\frac{3}{2}A\kappa \sqrt{\xi }\frac{d\xi }{dt}.$$
(3)
The dissipative material constant $`A`$ is a function of the viscous constants $`\eta _{1/2}`$, the Young modulus $`Y`$ and the Poisson ratio $`\nu `$ (for details see BSHP ):

$$A=\frac{1}{3}\frac{\left(3\eta _2\eta _1\right)^2}{3\eta _2+2\eta _1}\frac{\left(1\nu \right)\left(12\nu \right)}{Y\nu ^2}.$$
(4)
Combining the forces (2) and (3) one obtains the equation of motion

$$\frac{d^2\xi }{dt^2}+\frac{\kappa }{m^{\mathrm{eff}}}\left(\xi ^{3/2}+\frac{3}{2}A\sqrt{\xi }\frac{d\xi }{dt}\right)=0,$$
(5)
with $`m^{\mathrm{eff}}=m_im_j/(m_i+m_j)`$ and with the initial conditions $`\xi |_{t=0}=0`$ and $`d\xi /dt|_{t=0}=g`$.

In dimensionless variables, $`\widehat{\xi }\xi /\xi _0`$, $`\tau t/\tau _0`$, Eq. (5) reads DimAnalys

$$\frac{d^2\widehat{\xi }}{d\tau ^2}+\frac{5}{4}\widehat{\xi }^{3/2}+\frac{3}{2}\left(\frac{5}{4}\right)^{3/5}A\left(\frac{\kappa }{m^{\mathrm{eff}}}\right)^{2/5}g^{1/5}\sqrt{\widehat{\xi }}\frac{d\widehat{\xi }}{d\tau }=0,$$
(6)
with $`\widehat{\xi }|_{\tau =0}=0`$ and $`d\widehat{\xi }/d\tau |_{\tau =0}=1`$. The characteristic length $`\xi _0`$ is the maximal compression for the equivalent undamped (elastic) problem which can be found by equating the kinetic energy of the impact $`m^{\mathrm{eff}}g^2/2`$ and the elastic energy at the instant of maximal compression $`2\kappa \xi _0^{5/2}/5`$. As characteristic time $`\tau _0`$ we define the time in which the distance between the particles changes by the characteristic length just before the collision starts colltime

$$\xi _0\left(\frac{5}{4}\frac{m^{\mathrm{eff}}}{\kappa }\right)^{2/5}g^{4/5},\tau _0\xi _0/g.$$
(7)
The only term in Eq. (6) which depends explicitly on the system size and on material properties is the prefactor in front of the third term. If the scaling procedure affects the value of this term it will change the dynamics of the system. To ensure identical behavior of the scaled system, however, besides the identity of this prefactor, further requirements have to be met which will be discussed below.

Expanding our abbreviations we obtain

$`A\left({\displaystyle \frac{\kappa }{m^{\mathrm{eff}}}}\right)^{2/5}g^{1/5}`$
$`=\left(2\pi \right)^{2/5}AF(R_i,R_j)g^{1/5}\left({\displaystyle \frac{Y}{\rho \left(1\nu ^2\right)}}\right)^{2/5},`$ (8)
with $`\rho `$ being the material density. The function $`F(R_i,R_j)`$ collects all terms containing $`R_i`$ and $`R_j`$:

$$F(R_i,R_j)\frac{\left[R_iR_j/\left(R_i+R_j\right)\right]^{1/5}}{\left[R_i^3R_j^3/\left(R_i^3+R_j^3\right)\right]^{2/5}}$$
(9)
Scaling the radii by $`\alpha `$, the function $`F`$ scales $`F(R_i^{},R_j^{})=F(\alpha R_i,\alpha R_j)=\alpha ^1F(R_i,R_j)`$. Obviously, simple scaling of the system in general affects the prefactor Eq. (8) already via the scaling properties of $`F(R_i,R_j)`$, hence, in general the original system and the system where the lengths have been scaled by a factor $`\alpha `$ differ in their dynamic properties. More explicitly, one can show that naively scaling the system by a factor $`\alpha <1`$ will lead to a comparatively more damped dynamics.

To provide equivalent dynamical properties of the systems, therefore, we have to modify the material properties in a way to assure that the equations of motion of both systems are equivalent, which in turn assures that the prefactors (8) of the original system and the scaled system are identical.

One of the few things which cannot be modified in experiments with reasonable effort is the constant of gravity $`G`$. That implies that going from $`S`$ to $`S^{}`$ not only $`G`$ but all other accelerations must remain unaffected too

$$\left(\frac{d^{\mathrm{\hspace{0.17em}2}}x}{dt^2}\right)^{}=\frac{d^{\mathrm{\hspace{0.17em}2}}(\alpha x)}{d\left(t^{}\right)^2}=\frac{d^{\mathrm{\hspace{0.17em}2}}x}{dt^2},$$
(10)
yielding $`t^{}=\sqrt{\alpha }t`$. Hence, scaling all lengths $`x^{}=\alpha x`$ implies that times scale as $`t^{}=\sqrt{\alpha }t`$ if we require that the gravity constant stays unaffected. Thus, the clock in the scaled system $`S^{}`$ must run by a factor $`\sqrt{\alpha }`$ faster (or slower) than the clock in the original system. In other words, if in the original system we observe a phenomenon at time $`t=100`$ s, we will find the same effect in the scaled system at time $`t^{}=\sqrt{\alpha }\mathrm{\hspace{0.17em}100}`$ s. Scaling of time is a direct consequence of scaling the lengths if the constant of gravity has the same value in both systems.

In the scaled system the equation of motion of a particle contact reads

$$\frac{d^2\xi ^{}}{dt^2}+\frac{\kappa ^{}}{\left(m^{\mathrm{eff}}\right)^{}}\left(\xi ^{}\right)^{3/2}+\frac{3}{2}A^{^{}}\frac{\kappa ^{}}{\left(m^{\mathrm{eff}}\right)^{}}\sqrt{\xi ^{}}\frac{d\xi ^{}}{dt^{}}=0.$$
(11)
If we apply our scaling relations which were introduced above, i.e.

$$\xi ^{}=\alpha \xi ,\frac{d\xi ^{}}{dt^{}}=\sqrt{\alpha }\frac{d\xi }{dt},\frac{d^2\xi ^{}}{dt^2}=\frac{d^2\xi }{dt^2}$$
(12)
we obtain

$$\frac{d^2\xi }{dt^2}+\alpha ^{3/2}\frac{\kappa ^{}}{\left(m^{\mathrm{eff}}\right)^{}}\xi ^{3/2}+\frac{3}{2}\alpha A^{}\frac{\kappa ^{}}{\left(m^{\mathrm{eff}}\right)^{}}\sqrt{\xi }\frac{d\xi }{dt}=0.$$
(13)
Comparing with Eq. (5) we find the conditions to assure identity of the equations of motion:

$$\frac{\kappa ^{}}{\left(m^{\mathrm{eff}}\right)^{}}=\alpha ^{3/2}\frac{\kappa }{m^{\mathrm{eff}}},A^{}=\sqrt{\alpha }A.$$
(14)
Using the definitions of $`\kappa `$ (Eq. (2)) and $`m^{\mathrm{eff}}`$ yields finally

$$\left(\frac{Y}{\rho \left(1\nu ^2\right)}\right)^{}=\alpha \frac{Y}{\rho \left(1\nu ^2\right)},A^{}=\sqrt{\alpha }A.$$
(15)
If we choose material constants which obey Eqs. (15) we will obtain the original equation of motion after scaling the system back to its original size, i.e., both systems are equivalent.

When we incorporate the tangential force $`F^t`$ into the analysis, of course, we have to require that this force scales in the same way as any other force, namely

$$\frac{\left(F^t\right)^{}}{\left(m^{\mathrm{eff}}\right)^{}}=\frac{F^t}{m^{\mathrm{eff}}},$$
(16)
given that accelerations are invariant under scaling. This requirement has to be met by appropriately scaling the material constants, particularly the friction constant, resulting in an additional scaling equation. Its form depends on the underlying friction model, i.e. on the functional dependence of the tangential force on the geometry and the material properties as well as on the compression and the relative velocity. For instance, if we assume the most simple tangential force law

$$F^t=\mu F^n,$$
(17)
with $`\mu `$ being the friction coefficient we can conclude that the friction coefficient has to be invariant with respect to scaling

$$\mu ^{}=\mu .$$
(18)
In literature, one can find a variety of different models which might be more realistic for the description of the dynamics of a granular material (for an overview see, e.g., wolf ). Applying the same procedure for these laws will result in different scaling relations and, hence, in different conditions of the type Eq. (18). A thorough discussion of the scaling of more realistic tangential friction laws will be published elsewhere.

To demonstrate the derived scaling laws, we present results of a 2D Molecular Dynamics simulations of stationary flow of 1000 particles down an inclined plane. In the “original system” the length of the plane is 10 m (periodic b.c.), its slope is 15 degrees, and the particle diameters have been randomly chosen with an average of 10 cm. The material parameters are given in Table 1. Particles are subjected to normal viscoelastic forces as well as tangential forces (of the simple type discussed above) which, although not very realistic, is enough for the purpose of showing the scaling properties. Figure 1 shows a snapshot of the flow.

The solid curve in Fig. 2 represents the density profile perpendicular to the inclined plane of the original system and of the correctly scaled systems I and II, taking into account the necessary modification of the material parameters according to Eqs. (15,18). The modified material parameters are summarized in Table 1. It is not surprising that all three lines collapse perfectly, since in scaled variables the equations of motion of the original system and the correctly scaled system are identical. Therefore, particle trajectories and, hence, the density plot must be identical.

Figure 2 also shows two curves for the stationary density profile of the systems scaled by $`\alpha =1/16`$ (system I) and $`\alpha =16`$ (system II) where only the lengths have been scaled but erroneously the material parameters have been kept constant (incorrect scaling). These curves deviate from the profile of the original systems since the original and the scaled systems obey different equations of motion. As anticipated, the scaling by $`\alpha =16`$ only applied to lengths leads to a more elastic dynamics, i.e., comparatively more dilute flow (dot-dashed curve in Fig. 2). This effect is amplified when larger scaling factors are applied ($`\alpha =49`$), where in twice the simulation time, no steady state is found while the density profile keeps relaxing towards more dilute states. The instantaneous density profile is shown as a thin-solid curve in Fig. 2, as an indication that incorrect scaling can ultimately lead to completely deviated results.

The simulations show that the simple length scaling of a granular system by a constant factor $`\alpha `$ changes the dynamical properties significantly, if the material parameters are kept constant. In order to obtain identical results for a scaled system one has also to modify the material constants and redefine the unit of time. These necessary scaling relations are summarized in Table 2.

The knowledge about these scaling relations offers the possibility to scale down real world systems, e.g. geophysical or industrial granular systems, to sizes where laboratory experiments can be performed. If one scales down such a granular system one has to replace the original material by a material which meets the scaling requirements discussed in the text.

We want to give an example: Assume in the original system one deals with steel spheres ($`Y=20.610^{10}`$ Nm-2, $`\nu =0.29`$ and $`\rho =7,700`$kg m-3) of average radius $`\overline{R}=10`$ cm and system size of $`L=10`$ m. The property whose scaling behavior is known is $`Y/(\rho (1\nu ^2))=2.9210^7`$ m2s-2. One wishes to know (to measure) a certain value at time $`t=100`$ s. In the lab we perform the experiment with an equivalent system of size $`L^{}=1`$ m, i.e. we scale the system by the factor $`\alpha =0.1`$, including all radii. From the scaling relations we see that we have to find a material whose scaled property is $`Y^{}/(\rho ^{}(1\nu ^2))0.310^7`$ m2s-2. From tables tables  we find that we can use plexiglass ($`Y=0.3210^{10}`$ Nm-2, $`\nu =0.35`$ and $`\rho =1,200`$ kg m-3) in order to obtain this value. Therefore, we have to perform the experiment with plexiglass spheres and have to measure the value of interest at time $`t^{}=31.6`$ s.

One can imagine that not for all scaling factors $`\alpha `$ one will find a proper material, however, nowadays it is possible to manufacture materials which can meet demanding requirements, such as high softness along with a custom-designed damping constant.

The scaling properties have also consequences for Molecular Dynamics simulations of granular systems: namely, it is not sufficient in simulations to provide relative parameters such as the container size in units of the particle radius. As demonstrated, the result of a simulation (and, of course, also of a real world experiment) depends on absolute values.

###### Acknowledgements.

This work has been supported by Deutsche Forschungsgemeinschaft via grants PO 472/6-1, PO 472/7-1.


format: text

text:
# The zero point field and the reduction of the wave packet in semiclassical electrodynamics.

## 1 Introduction.

Many physicists think that the zero point field is introduced by quantum electrodynamics (QED); historically, its knowledge is anterior to QED, and the computation of its mean amplitude by Planck  and Nernst  allowed the identification of the energy of the optical modes with the energy of quantum harmonic oscillators, which is the starting point of QED.

The interpretation of the zero point field is, however, fundamentally different in the classical theory and in QED: in the classical theory, the electromagnetic field is, in all points, a physical quantity and any solution of Maxwell’s equations represents precisely an optical mode; in QED, the electromagnetic field is a mathematical tool substituted to the wave function of the photon, and the postulate “reduction of the wave packet”changes the mode corresponding to a photon, thus the electromagnetic field. This fundamental difference is generally not clearly set, even in the best books (Milonni ).

For a practical use, a confusion between the quantum and classical zero point fields is not important because as long as the reduction of the wave packet is not used, the equations are exactly the same; this similitude was well studied by Marshall and Santos  and previous authors; our aim is seeking after an extension of their work to other optical problems, finding, in particular, a classical equivalent to the reduction of the wave packet, so that it becomes possible to use either interpretation in all problems.

The choice between the classical and QED interpretations depends on the spirit of the physicist. Using the interpretations of both the classical (or semi-classical) and quantum electromagnetic theories increases the reliability of the explanations, each theory having its specific advantages and traps; QED introduces more compact computations, more results use it; classical electrodynamics requires more complicated computations, but it seems next to the experiments; as QED is more used than the classical theory, its traps appear more in the literature :

Confusing the two interpretations of the electromagnetic field, thus using wrongly the postulate “reduction of the wave packet”, the best specialists of the quantum theory tried to discourage Townes from discovering the maser . To avoid such errors Lamb and al.  propose to reject the photon; in despite of their affirmation, is seems a de facto rejection of the principles of quantum electrodynamics.

For some physicists QED is useless in high field experiments: “The subtle interplay between real and imaginary parts of the complex linear and nonlinear susceptibilities follows quite naturally from the semiclassical treatment. How can this information be obtained from a theory in which the fields are quantified ? \[…\] The semiclassical theory which is used in this monograph will describe all situations correctly in a much simpler fashion”(Bloembergen ).

To describe correctly the semiclassical theory, we qualify the electromagnetic field as follows :

-“genuine”for the total electromagnetic field, with the genuine electric field noted $`\widehat{E}`$;

“zero point”for the field which remains after the largest, coherent, physically possible absorption of the genuine field; notation: $`E_0`$;

“conventional”for the genuine minus the zero point; as it is the commonly considered field, the conventional electric field will be simply written $`E`$;

“stochastic”for a stochastic zero point field.

Section 2 reminds the origin and properties of the classical fields, in particular the zero point fields.

Section 3 gives the correct relation between the signal of a photoreceiver and the electromagnetic field.

Section 4 applies the result of section 3 to the computation of the fourth order interferences.

Section 5 proposes a classical equivalent of the postulate “reduction of the wave packet”.

## 2 Recall of the properties of the classical electromagnetic field.

As Maxwell’s equations are linear in the vacuum, a linear combination of modes is a solution; with given limiting conditions, complete sets of modes allow to develop any mode. The energy of a mode is finite if its time and space extensions are limited; supposing that the mode is alone in the space, its energy is: $`W=w(\nu )d\nu `$ where $`w(\nu )`$ is the energy $`\frac{1}{2}(ϵ_0E^2+\mu _0H^2)dv`$ relative to the fields whose frequencies are between $`\nu `$ and $`\nu +\mathrm{d}\nu `$, at any instant. The usual normalisation is such that $`\frac{w(\nu )\mathrm{d}\nu }{\nu }=h`$. Two modes are orthogonal if the energy of their sum is the sum of their energies.

A punctual source $`S_O`$ of electromagnetic field (i. e. the most general multipole) in a point $`O`$ may be developed linearly using the derivatives of the three-dimensional Dirac’s distribution $`\delta _O`$:

$$S_O=\underset{p,q,r0}{}f_{t,O,p,q,r}_x^p_y^q_z^r\delta _O$$
(1)
where $`f_{t,O,p,q,r}`$ is a distribution relative to the variable $`t`$, and where $`f_{t,O,0,0,0}`$ is a constant.

Here the $`f_{t,O,p,q,r}`$ are supposed regular functions of time $`f_{O,p,q,r}(t)`$, and the number of non-zero $`f_{O,p,q,r}(t)`$ is supposed finite. For instance, the $`Oz`$ oriented dipole of a molecule in $`O`$ is proportional to the distribution $`f_{O,0,0,1}(t)_z\delta _O`$, with, for instance $`f_{O,0,0,1}(t)=\mathrm{sin}\omega t`$.

The fields radiated by these multipoles at the point source $`O`$ are mathematically singular in $`O`$, it is an approximation of the physical problem discussed later. These fields are often, improperly qualified “strictly spherical”although they are not invariant by a rotation around any axis $`Ou`$; the reason is that the fields on two spheres of centre $`O`$ have the same dependence on the Euler angles while it is not the case for a beam focussed in $`O`$ which diffracts.

Consider a particular source $`S_O`$, and the source $`S_O^{}`$ obtained replacing $`f_{O,p,q,r}(t)`$ by $`f_{O,p,q,r}(t)`$; by a time inversion, the source $`S_O^{}`$ becomes a source $`S_O^{}`$ which absorbs a “strictly spherical”field; summing the first and the last problems, $`S_O`$ and $`S_O^{}`$ cancel, it remains a “strictly spherical”field $`F_O`$ which converges to $`O`$ then diverges. Thus, the light emitted by a multipole will be considered as produced by the evolution of a “strictly spherical”converging beam, without the multipole. This field $`F_O`$ is regular everywhere, except in $`O`$.

The field $`F_O`$ corresponding to $`S_O`$ is absorbed if it is completely cancelled by an opposite field; this opposite field cannot be a linear combination of a finite number of fields $`F_A,F_B,F_C\mathrm{}`$ corresponding to multipoles placed in points $`A,B,C\mathrm{}`$ different of $`O`$, because these fields are regular in $`O`$. Thus the absorption of the electromagnetic field radiated by a point source requires an infinity of point absorbers, the residual field constitutes the zero point field. The building of the zero point field shows that it is an ordinary field.

An objection to the previous mathematical description may be set:

The dipoles introduced in spectroscopy are not punctual; but their dimensions are supposed small in comparison with their distances, so that the conclusion remains approximately valid: the absorption of the field radiated by a small source requires a lot of small absorbers, a long time during which it remains a part of the field radiated by the dipole.

Far from sources, the zero point field $`E_0`$ is generated by a large number of sources supposed incoherent, so that it is stochastic, it is characterised by its mean amplitude. If an atom emits a photon (that is, in the classical theory, a wave corresponding to an energy $`h\nu `$), the radiated field, large near the atom, is not immediately compensated by propagation, diffusion or absorption; thus the fluctuations of the field are shaped by next or coherent sources, so that “stochastic”is not a sure property of the zero point field.

A macroscopic consequence of the structuring of the zero point field by matter is observed in the Casimir effect : long wavelengths are rejected from inside two parallel plates, so that a lower pressure of radiation attracts the plates.

Following Einstein, the stimulated amplification of a mode by a source depends only on the amplitude of the true, genuine incident field $`\widehat{E}`$; the conventional, usual field $`E=\widehat{E}E_0`$ is a purely mathematical object, it has no physical, individual existence.

Remark that the field radiated by a dipole is the same for an emission, a refraction or an absorption; in the last case the dipolar field cancels a part of an external field; the equality of Einstein’s $`B`$ coefficients for stimulated emission and absorption (demonstrated by thermodynamics) appears natural. The absorption of light is considered now as a decrease of the energy in a mode down to the zero point energy, the emission is an amplification of the energy in a mode; thus all systems are connected at least by the zero point field, even in classical physics isolated systems do not exist.

The oscillating dipole is a particular system of moving charges; any moving electron radiates a field, but if it belongs to a stationary system the interference of the radiated field with the zero point field does not change the energy of this last field, in the average. Sommerfeld’s electron does not fall on the proton, but the fluctuations of its interaction with the zero point field produce the Lamb shift .

Remark that the existence of the zero point field depends on the emission of an electromagnetic field by charged particles; if the charged particles are unable to emit a very high frequency wave, there is no zero point for high frequencies, no UV divergence.

Planck’s constant $`h`$ connects the density of zero point energy in the universe to microphysics; is $`h`$ a cosmological or a microphysical constant?

## 3 Absorption and detection.

The signal of a photoelectric cell cannot be a function of the conventional field $`E`$ which has no physical existence; for a long time, the sensitivity of the detectors of light was bad so that the zero point field could be neglected and the genuine field $`\widehat{E}`$ could be replaced by the conventional field

A light receiver is excited by an attenuation of the energy of a mode from a value higher than $`h\nu /2`$ to nearly $`h\nu /2`$. Generally this is possible if the mode was amplified by a source, but in the dark, cold, good photocells generate a noise signal which seems produced by the particularly large fluctuations of the zero point field.

The net available energy on a receiver is proportional to the difference between the input and output energies in the exciting mode; set $`E_0`$ the total zero point field received by the cell, and $`E_0^{}`$ the small fraction of this field amplified by the source, with an amplification coefficient $`1+\gamma `$; the available energy on the source is

$$\widehat{E}^2E_0^2=(\gamma E_0^{}+E_0)^2E_0^2=2\gamma E_0^{}E_0+(\gamma E_0^{})^2.$$
(2)
Usually $`E_0`$ is much smaller than $`\gamma E_0^{}`$, $`2\gamma E_0^{}E_0`$ is neglected, the usual rule is got; on the contrary, supposing that $`\gamma `$ is small, $`(\gamma E_0)^2`$ is neglected; for long observations, the time- average of the total zero point amplitude is zero, so that no signal is observed; for a short observation the phases are constant so that the detected signal is proportional to $`\gamma E_0^{}`$ that is to the amplitude , not the intensity of the conventional field, with, however, the additional factor $`E_0`$ which has an arbitrary phase. The phase factors of $`E_0`$ and $`E_0^{}`$ fluctuate almost independently: without a sophisticated detection nothing appears.

The classical emission or absorption is modelled by an excitation of a mono- or polyatomic molecule by an electromagnetic field up to a barrier between the two involved relative minimums of potential $`u_1`$ and $`u_2`$; if the initial and final states are stationary $`|u_2u_1|=h\nu `$ . But our macroscopic experiments usually use nearly plane modes while the molecules emit light generally through dipoles or quadrupoles. Quantum mechanics transforms the geometry of the waves using the “reduction of the wave packet”. Section 5 gives the classical reduction of the wave packet.

## 4 Fourth order interferences.

A sophisticated detection is performed in the fourth order interference experiments with photon counting: two elementary measurements are done while the phases of $`E_0`$ and $`E_0^{}`$ are constant (see, for instance, ). The result of these experiments is easily got qualitatively using the conventional rules , but the contrast of the computed fringes is lower than shown by the experiments. In the simplest experiment  two small photoelectric cells are put in the interference fringes produced by two point sources; the interferences are not visible because they depend on the fast changing difference of phase $`\varphi `$ of the modes of the zero point field amplified by the sources. The sources are weak; the signal is the correlation of the counts of the cells, so that the phases of $`E_0`$ and $`E_0^{}`$ are eliminated.

Distinguishing the photoelectric cells by an index j equal to 1 or 2, set $`\delta _j`$ the difference of paths for the light received by the cells. The amplitude of the conventional field received by a cell is proportional to $`\mathrm{cos}(\pi \delta _j/\lambda +\varphi /2)`$, so that, assuming the linearity in the conventional field, the probability of a simultaneous detection is proportional to

$$E_0^2E_{0}^{}{}_{}{}^{2}\mathrm{cos}(\frac{\pi \delta _1}{\lambda }+\frac{\varphi }{2})\mathrm{cos}(\frac{\pi \delta _2}{\lambda }+\frac{\varphi }{2})\mathrm{cos}(\frac{\pi \delta _1}{\lambda }+\frac{\varphi }{2})\mathrm{cos}(\frac{\pi \delta _2}{\lambda }+\frac{\varphi }{2}).$$
(3)
The mean value of this probability got by an integration over $`\varphi `$ is zero for $`\delta _1\delta _2=\lambda /2`$, so that the visibility has the right value 1. Assuming the usual response of the cells proportional to the square of the conventional field, the visibility would have the wrong value 1/2 . Equation 3 may be found directly computing the interference of the zero point modes amplified by the sources.

## 5 Classical reduction of the wave packet.

The reduction of the wave packet breaks the symmetry of the waves, transforming, in particular, a local wave into a beam, for instance a spherical dipolar wave into a plane wave.

The polarisation of a transparent matter by a light beam may be observed by a variation of the energy levels, or detecting Kerr effect,…Thus the beginning of a pulse of light must transfer energy to matter, and this energy is recovered in the tail of the pulse (except for a small incoherent Rayleigh scattering). In the tail, the field is amplified, although there is no transition, no inversion of population, the polarisation mixing only slightly the initial state of the molecules with other states. This power of amplification applies not only to the exciting mode, so that many modes, usually initially at the zero point, are amplified, later reabsorbed in the medium : there is a dynamical equilibrium between the exciting field, the other modes and the polarisation of the molecules. The modes which are excited are dipolar or quadrupolar, they radiate far only the small incoherent Rayleigh field : they may be qualified “local”. On the contrary, the interactions with the light pulse are strong because they are coherent.

As the local modes are amplified, the strongest and longest fluctuations of their field may be able to excite molecules up to a barrier, such that an absorbing transition occurs; the mean energy of the local field, then of the molecules, is decreased, the amplification of the tail of the pulse is decreased, the medium has absorbed the light. In a laser, a similar process explains the coherent amplification by incoherently pumped molecules.

The emission of a field during the absorption of a quantum of energy is not instantaneous; during this emission and a short time after it, the probability for a strong and long fluctuation of the field is lowered, so that a sub-poissonian photon statistic appears; neglecting the space-time structuring of the field leads to the poissonian statistic .

## 6 Low level “Impulsive Stimulated Raman Scattering”

Consider a nonlinear property which may be developed as a series of the electric field of an electromagnetic field. It must be a function of the genuine field. The development of this function using the conventional field is linear for low values of this field: there is a linearisation of the properties at low values of the conventional field.

If the property cannot be written as a series development, a specific treatment is necessary; consider, for instance the adaptation to the natural incoherent light of the “Impulsive stimulated Raman Scattering”(ISRS); ISRS is known since 1968  and now commonly used . It is not a simple Raman scattering, but a parametric effect, combination of two space-coherent Raman scatterings, so that the state of the interacting molecules is not changed. The hot222The temperature of a spectral line is deduced from Planck’s laws exciting beam and its scattered beams interfere into a single frequency, redshifted, beam; the cold beam is blueshifted.

The interference of two collinear light beams into a single frequency beam is often performed involuntarily, for instance using a Michelson interferometer in with a moving mirror produces a Doppler frequency shift. The mixture of the frequencies requires an observation shorter than the beats of the two beams. The resulting frequency is intermediate, in proportion of the conventional amplitudes of the initial frequencies. In ISRS, the scattered amplitude which is permanently mixed with the incident amplitude remains extremely low and the Raman frequency is much lower than the exciting frequency, so that the relative frequency shift $`\mathrm{\Delta }\nu /\nu `$ is nearly proportional to the conventional scattered amplitude and to the Raman frequency333Computing this interference process with the genuine amplitudes is more difficult because it is necessary to take into account the zero point field of the three involved modes. .

The conventional scattered amplitude is proportional to the square of the incident genuine amplitude, square which, at a low level, is proportional to the conventional amplitude. More simply, while at high intensities the scattered amplitude is stimulated, at low intensity it is spontaneous, but it remains coherent in the absence of collisions because the behaviour of all molecules on a wave surface is the same (just as in refraction). Thus, the relative frequency shift is nearly constant, only slightly modified by the dispersion of the tensor of polarisability.

ISRS is obtained using ultrashort light pulses, that is “pulses shorter than all relevant time constants”, usually femtosecond laser pulses. The relevant time constants in a gas may be adapted to incoherent light :

i) to avoid that the collisions destroy the coherence of the excitation of the molecules during the pulses, the pressure must be very low;

ii) to obtain an interference of the scattered and incident lights into a single frequency light, the period which corresponds to the virtual Raman transition must be larger than the length of the impulsions. The molecules must have transitions in the radiofrequencies, generally hyperfine transitions.

The Universe, provides good experimental conditions for a confusion of this interaction with a Doppler effect: the paths are long and the pressures often low, a lot of observed polyatomic molecules or atoms perturbed by Zeeman effect have hyperfine structures. The absorption spectra of the molecules which are destroyed at their first collision, H$`{}_{}{}^{+}{}_{2}{}^{}`$ for instance, cannot be seen because the redshift simultaneous with the absorption of their spectra widens, thus weakens, their lines.

## 7 Conclusion.

The electromagnetic fields, in particular the zero point field, often, improperly qualified “stochastic”, obey the same equations in quantum and classical theories. Both, very different interpretations, are useful, giving the same final results, having their specific advantages : Quantum electrodynamics provides ready to use properties or postulates, but a common improper use of some of its concepts, the photon for instance, leads to wrong conclusions ; classical electrodynamics is more intuitive, but it requires often more complicated demonstrations.

An isolated system cannot exist in classical optics because the electromagnetic fields expand everywhere, although they are amplified or attenuated by matter.

The teachers should point the approximation made neglecting the zero point field and replace unnecessarily approximate rules, for instance the first Planck’s law, by the rigorous rules.


text:
# Transport in random quantum dot superlattices

## I Introduction

Latest advances in nanotechnology make it possible to growth quantum dot (QD) superlattices . In view of the analogy between atoms and QDs, it is expected that strongly confined levels overlap when QDs are closely packed. Although this analogy cannot be complete since carriers in QDs are influenced by phonons, defects and interface states, the resulting collective levels will depend on the particular arrangement of QDs. In this sense, QD arrays grown by molecular beam epitaxy can be completely random , partially regimented  or may be regularly stacked (high regimentation) . Electronic states in highly regimented QDs are adequately described within envelope-function approximation with a three dimensional Kronig-Penney model and the occurrence of miniband have been established . However, the lack of periodicity observed in random QDs superlattices demands different approaches.

In the present work we introduce a two-dimensional effective-mass model to study the effects of scattering by disorder on electron transmission through random quantum-dots superlattices. To this end, the Ben Daniel-Duke equation is discretized, boundary conditions are discussed and scattering solutions are found by means of the transfer-matrix method for any arbitrary array of QDs. The model is worked out in a two-dimensional space for computational limitations, although it will be clear that generalization to three dimensions is rather straightforward. Finally, we present the numerical results for the transmission coefficient through random arrays of coupled QDs and the main conclusions of the work.

## II Model

We consider the Ben Daniel-Duke equation with constant effective mass $`m^{}`$ at the $`\mathrm{\Gamma }`$ valley for the electron envelope function $`\chi (y,z)`$. In order to find numerically the single electron states, the whole sample is divided into three different regions, namely left (I) and right (III) contacts and the random QDs (II), where scattering by disorder takes place. Figure 1 shows a schematic view of the three spatial regions. We then consider a mesh with lattice spacings $`a_y`$ and $`a_z`$ in the $`y`$ and $`z`$ directions, respectively. Defining $`t_y\mathrm{}^2/(2m^{}a_y^2)`$ and $`t_z\mathrm{}^2/(2m^{}a_z^2)`$, we obtain the following discrete effective-mass equation

$$t_z(\chi _{n+1,m}+\chi _{n1,m})+t_y(\chi _{n,m+1}+\chi _{n,m1})+(U_{n,m}2t_z2t_y)\chi _{n,m}=E\chi _{n,m}.$$
(1)
The potential term $`U_{n,m}`$ in Eq. (1) is given by the conduction-band edge energy at the point $`(na_y,ma_z)`$. Therefore, disorder enters the equation through this diagonal term. Contacts are characterized by flat band conditions, $`U(n,m)=0`$, in the absence of applied electric field. Clearly, the effects of the applied field can be easily taken into account within the present approach by adding a linearly varying potential of the form $`U_\mathrm{F}(n,m)=eFma_z`$, $`F`$ being the applied field.

## III Transfer matrix formalism

In order to solve the tight binding-like equation (1) we use the transfer matrix method based on the solutions calculated for each slide of the system along the $`z`$ direction. For the sake of simplicity in the calculation, we define $`\stackrel{}{\varphi }_n`$ ($`n=0,1,\mathrm{},N+1`$) as a vector whose components are $`\varphi _n^m\chi _{n,m}`$ (m=1,…,M). Here $`M`$ and $`N+1`$ are the number of mesh divisions in the $`y`$ and $`z`$ directions, respectively. Thus, Eq. (1) can be cast in a more compact form

$$\left(\begin{array}{c}\stackrel{}{\varphi }_{n1}\\ \stackrel{}{\varphi }_n\end{array}\right)=\left(\begin{array}{cc}t_z^1(E_n)& \\ & 𝒪\end{array}\right)\left(\begin{array}{c}\stackrel{}{\varphi }_n\\ \stackrel{}{\varphi }_{n+1}\end{array}\right),$$
(2)
where $``$ and $`𝒪`$ are the $`M\times M`$ identity and null matrices respectively. The matrix $`_n`$ splits into the form $`_n=_n+_n`$. The diagonal elements of the tridiagonal matrix $`_n`$ are $`(_n)_{mm}=U_{n,m}2t_z2t_y`$ while nonvanishing off-diagonal elements equal $`t_y`$. The matrix $`_n`$ depends on the boundary conditions to be specified later.

We can obtain the expression of the envelope function amplitudes in the left contact as a function of the amplitudes in the right one

$$\left(\begin{array}{c}\stackrel{}{\varphi }_0\\ \stackrel{}{\varphi }_1\end{array}\right)=𝒯^{(N)}\left(\begin{array}{c}\stackrel{}{\varphi }_N\\ \stackrel{}{\varphi }_{N+1}\end{array}\right),𝒯^{(N)}\underset{n=1}{\overset{N}{}}\left(\begin{array}{cc}t_z^1(E_n)& \\ & 𝒪\end{array}\right),$$
(3)
where $`𝒯^{(N)}`$ is the transfer matrix for the heterostructure.

## IV Scattering solutions

The envelope functions within the contacts will be determined by the boundary conditions. These boundary conditions are open in the $`z`$ direction, and periodic on each slide, namely in the $`y`$ direction. Consequently, all elements of $`_n`$ vanish except $`(_n)_{1M}=(_n)_{M1}=t_y`$. The former conditions imply plane wave solutions in the $`z`$ axis, and the latter yield an energy discretization on $`y`$. As a consequence, this discretization results in a number of transverse channels equal to the number of points in the transverse mesh direction. Considering $`U_{n,m}`$ to be constant at the contacts, that is, considering perfect leads so no voltage drop occurs within regions I and III, and setting $`\chi _{n,1}=\chi _{n,M+1}`$, a particular solution of Eq. (1) is given by

$`\chi _{n,m}^l`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{𝒩_l}}}\mathrm{exp}\left(i{\displaystyle \frac{2\pi l}{M}}m\right)\mathrm{exp}\left(ik_la_zn\right)`$ (5)
$`+`$ $`{\displaystyle \underset{j=1}{\overset{M}{}}}\widehat{r}_{lj}{\displaystyle \frac{1}{\sqrt{𝒩_j}}}\mathrm{exp}\left(i{\displaystyle \frac{2\pi j}{M}}m\right)\mathrm{exp}\left(ik_ja_zn\right),(m,n)\mathrm{I}.`$ (6)
at the left contact while at the right contact we have

$$\chi _{n,m}^l=\underset{j=1}{\overset{M}{}}\widehat{t}_{lj}\frac{1}{\sqrt{𝒩_j}}\mathrm{exp}\left(i\frac{2\pi j}{M}m\right)\mathrm{exp}\left(ik_ja_zn\right),(m,n)\mathrm{III}.$$
(7)
where the normalization constant is choosen to ensure that all the propagating modes carry the same current

$$𝒩_j=\frac{1}{a_y^2}\mathrm{sin}^2\left(\frac{2\pi j}{M}\right)+\frac{1}{a_z^2}\mathrm{sin}^2\left(k_ja_z\right)$$

and

$$k_j=\frac{1}{a_z}\mathrm{cos}^1\left\{\frac{1}{2t_z}\left[E2t_y\left(\mathrm{cos}\frac{2\pi j}{M}1\right)\right]+1\right\}.$$

These expressions remain valid for an applied electric field $`F`$ provided the electronic momentum $`k_j`$ is substituted within region $`III`$ in the following way

$`k_j`$ $``$ $`q_j`$ (8)
$`q_j`$ $`=`$ $`{\displaystyle \frac{1}{a_z}}\mathrm{cos}^1\left\{{\displaystyle \frac{1}{2t_z}}\left[E+eV2t_y\left(\mathrm{cos}{\displaystyle \frac{2\pi j}{M}}1\right)\right]+1\right\}.`$ (9)
Here $`V=FL`$, where $`L`$ is the lenght of region $`II`$.

The matrices $`\widehat{r}`$ and $`\widehat{t}`$ appearing in equations (5) and (7) are the reflection and transmission matrices and they are responsible for the channel mixing due to scattering events. Thus, $`\widehat{r}_{ij}`$ represents the probability amplitude for an electron impinging in channel $`i`$ to be reflected into the chanel $`j`$. Note that the solution within region II is unknown. Actually, we are not interested in this solution since all the physics of the scattering problem is contained in the mixing matrices $`\widehat{t}`$ and $`\widehat{r}`$. In particular, we can compute the conductance. From the Landauer-Büttiker formalism, the zero temperature two-leads multichannel conductance can be calculated using de Fisher-Lee formula

$$G=\frac{2e^2}{h}\text{Tr}(\widehat{t}^{}\widehat{t}),$$
(10)
where Tr stands for the trace of the matrix.

## V Configurational and morphological disorder

In order to describe epitaxially random QD superlattices, we consider they are arranged on a nonregular lattice (configurational disorder). But, in addition, there exists another source of disorder due to spatial inhomogeneities that make the shape to be slightly different (morphological disorder). To mimic both types of disorder we assume an array of rectangular QDs, randomly displaced from the regular lattice sites, whose size also change randomly from dot to dot (see Fig. 1). To avoid the profusion of free parameters, we consider that the energy of the confinement potential $`\mathrm{\Delta }E_c`$ provided by the high bandgap semiconductor is the same for every QD. This is not a serious shortcoming since spatial inhomogeneities of the conduction-band offset and fluctuations of individual QDs shapes yield essentially the same results. In addition, our model could easily deal with nonconstant values of $`\mathrm{\Delta }E_c`$ if further improvements are required.

We will separate the effects of configurational and morphological disorder for the sake of clarity. Configurational disorder means that every QD shifts its position $`\delta 𝐫=(\delta y,\delta z)`$ while its size $`d_y\times d_z`$ remains unchanged (see Fig. 2). Here $`\delta y`$ and $`\delta z`$ are random uncorrelated variables with zero mean and distributed according to box probability functions of width $`W_y`$ and $`W_z`$, respectively. To simulate the change of the shape (morphological disorder) we consider that the QD is enlarged along the $`Z`$ axis an amount $`\delta \zeta `$ while its center stays on a regular lattice (see Fig. 2), where $`\delta \zeta `$ is uniformly distributed with zero mean and width $`W_\zeta `$.

## VI Results

As a working example, we have performed several numerical calculations in order to study the effect of both configurational and morphological disorder over the conductance of random QDs superlattices made of GaAs-InxGa1-xAs heterojunctions. We have considered regimented and disordered $`4\times 4`$ arrays of QDs to elucidate the effects of randomness. The conduction-band offset, $`\mathrm{\Delta }E_c`$, is taken to be $`70\%`$ of the difference of the gaps $`\mathrm{\Delta }E_g`$ in strained GaAs-InxGa1-xAs heterojunctions, where $`\mathrm{\Delta }E_g=1.45x`$eV. For definiteness we set $`x=0.35`$ and consequently $`\mathrm{\Delta }E_c=0.35`$eV. In addition, since we are mainly interested in the effects of the coupling between the QDs through the high bandgap semiconductor rather in the confined levels of individual QD, we have taken $`m^{}=0.067`$ in units of the bare electron mass, corresponding to the embedding semiconductor. Let us mention that the the model can be easily generalized to include a different effective mass inside the QDs. The size of the regimented QDs is $`d_y\times d_z`$ with $`d_y=8.0`$nm and $`d_z=1.6`$nm. The separation between centers in the regimented array is $`14.0`$ and $`6.8`$nm along the $`Y`$ and $`Z`$ axes, respectively. The number of mesh points along the two spatial directions are $`M=50`$ and $`N+1=39`$.

As typical results of our simulations, Fig. 3 shows the conductance versus Fermi energy, measured from the conduction-band edge in InxGa1-xAs, in the absence of applied field. Solid line corresponds to regimented QDs. The coupling between the QDs splits the energy levels and results in the formation of minibands . From Fig. 3 we observe the occurrence of two well-defined minibands below the barrier when disorder is absent, in agreement with previous results . Each band is characterized by four main conductance peaks and each peak is the convolution of four peaks that cannot be resolved except for the lower one in the higher miniband (see inset of Fig. 3).

Transport though the miniband changes as soon as some degree of randomness is considered in the model, as expected. Mean size and separation between QDs as well as fluctuations around those values strongly depend on the growth conditions (e.g. growth temperature) and subsequent thermal treatments . As an example, Fig. 3 shows the conductance for configurational ($`W_y=2.0`$nm and $`W_z=1.2`$nm) and morphological ($`W_\zeta =0.8`$nm) disorder, obtained by averaging over $`100`$ realizations of the disorder. As a main point, we notice that electronic states in the array of random QDs behaves like an amorphous material since the conductance strongly decreases while the allowed energies broadens due to the fluctuations of QD energy levels for each realization of the disorder.

Concerning the effects of a uniform applied electric field, we have computed the conductance for a given Fermi energy ($`E=0.27`$eV) as a function of the applied bias $`V`$. The potential was assumed to drop uniformly across region II. Figure 4 shows the conductance versus applied bias for regimented as well as random QD superlattices. Regimented QD superlattices present three well-defined negative $`dG/dV`$ regions due to resonant tunneling through the QDs. The observed peak-to-valley ratios become worse as soon as disorder is included in the calculation, especially in the case of configurational disorder.

## VII Conclusions

In summary, we have presented a method to study electron transport through random QD superlattices. The method is based on the transfer matrix formalism applied to the discrete Ben Daniel-Duke Hamiltonian for the electron envelope function. A careful analysis of the scattering solutions under appropriate boundary conditions (periodic and open along the lateral and longitudinal directions, respectively) allows us to obtain the two-channel conductance. For regimented QDs (regular array of QDs) the conductance shows clear signatures of the miniband structure, as previously predicted by a Kronig-Penney model for strongly coupled QDs . However, the novelty of the model lies in the fact that random QD superlattices can also be studied. To this end, two models of disorder (configurational and morphological) have been introduced. Disorder reduces the conductance due to Anderson localization of the envelope functions and broadens the allowed energy bands. The characteristic $`G`$$`V`$ presents several regions of negative $`dG/dV`$, although the peak-to-valley ratios strongly decrease due to disorder.

###### Acknowledgements.

Work in Madrid was supported by DGI-MCyT (Project MAT2000-0734) and CAM (Project 07N/0075/2001). P. Orellana would like to thank Milenio ICM P99-135-F and Cátedra Presidencial de Ciencias for financial support.


text:
# Anomalous and dimensional scaling in anisotropic turbulence

## Abstract

We present a numerical study of anisotropic statistical fluctuations in homogeneous turbulent flows. We give an argument to predict the dimensional scaling exponents, $`\zeta _d^j(p)=(p+j)/3`$, for the projections of $`p`$-th order structure function in the $`j`$-th sector of the rotational group. We show that measured exponents are anomalous, showing a clear deviation from the dimensional prediction. Dimensional scaling is subleading and it is recovered only after a random reshuffling of all velocity phases, in the stationary ensemble. This supports the idea that anomalous scaling is the result of a genuine inertial evolution, independent of large-scale behavior.

In recent years a huge amount of theoretical, numerical and experimental work has been done in order to study anisotropic turbulent fluctuations . Typical questions go from the theoretical point of calculating and measuring anomalous scaling exponents in anisotropic sectors , to the more applied problem of quantifying the rate of recovery of isotropy at scales small enough . Another important issue is the universality of anisotropic scaling exponents, i.e. whether they are an intrinsic characteristic of the Navier-Stokes non-linear evolution or they are fixed by a dimensional matching with the external anisotropic forcing.

Important steps forward in the analysis of anisotropic fluctuations have recently been done in Kraichnan models, i.e. passive scalars/vectors advected by isotropic, Gaussian and white-in-time velocity fields , with a large-scale anisotropic forcing . In those models, anomalous scaling arises as the result of a non-trivial null-space structure for the advecting operator (zero modes). Also, correlation functions in different sectors of the rotational group show different scaling properties. Scaling exponents are universal: they do not depend on the actual value of the forcing and boundary conditions, and they are fully characterized by the order of the anisotropy. Non-universal effects are felt only in coefficients multiplying the power laws.

Similar problems, like the existence of scaling laws in anisotropic sectors and, if any, the values of the corresponding scaling exponents are at the forefront of experimental, numerical and theoretical research for real turbulent flows. Only few indirect experimental investigations of scaling in different sectors  and direct decomposition in numerical simulations  have been attempted up to now.

The question is still open, evidences of a clear improving of scaling laws by isolating the isotropic sector have been reported, supporting the idea that the undecomposed correlations are strongly affected by the superposition of isotropic and anisotropic fluctuations . On a theoretical ground, only recently it has been highlighted the potentiality of $`SO(3)`$ decomposition to quantify different degrees of anisotropies for any correlation function . On the basis of this analysis, preliminary experimental evidences of the existence of a scaling law also in sectors with total angular momentum $`j=2`$ have been reported . The value of the exponent for the second order correlation function being close to the dimensional estimate $`\zeta _d^{j=2}(2)=4/3`$,  (where, from now on, subscript $`d`$ denotes the dimensional value).

Tipically, experimental investigations in real turbulent flows are flawed by the contemporary presence of anisotropies and strong non-homogeneities. The meaning of $`SO(3)`$ decomposition becomes opaque in presence of strong non-homogeneities and also the very existence of scaling laws cannot be given for granted .
To overcome such difficulties, some of us performed (see ) the numerical investigation of a “Random-Kolmogorov Flow” (RKF), a fully periodic Kolmogorov flow with random forcing phases, $`\delta `$-correlated in time.

In this Letter, we present a more extended analysis of the same data set, but focusing on new evidences that anisotropic scaling exponents are indeed universal and anomalous, i.e. they do not follow simple dimensional scaling. In order to do this, we also give a clear phenomenological background able to predict the dimensional scaling in anisotropic sectors. The Letter is organized as follows. First, we present a simple dimensional argument for all anisotropic sectors of structure functions of any order. With respect to this dimensional prediction, we show that anisotropic exponents are indeed anomalous. Moreover, we show that by performing a random re-shuffling of the velocity phases (in the stationary ensemble of anisotropic configurations), the leading anomalous scaling is filtered out and the sub-leading dimensional prediction is recovered. This is both a test of our dimensional prediction and a clean indication that forced velocity correlations are dominated by inertial terms; in this sense one may refer to them as the “equivalent” of the zero-modes responsible for anomalous scaling and universality in linear hydrodynamical problems . These findings leads to conclude that anisotropic fluctuations in turbulence are anomalous and universal.
We recall few details on the numerical simulations . The RKF is fully periodic; the large-scale anisotropic random forcing points in one direction, $`\widehat{𝒛}`$, has a spatial dependency only from the $`x`$ coordinate and it is different from zero at the two wavenumbers: $`𝒌_1=(1,0,0),𝒌_2=(2,0,0)`$. Namely, $`𝒇_i(𝒌_{\{1,2\}})=\delta _{i,3}f_{\{1,2\}}\mathrm{exp}(\theta _{\{1,2\}})`$, where $`f_{\{1,2\}}`$ are fixed amplitudes and $`\theta _{\{1,2\}}`$ are independent random phases, $`\delta `$-correlated in time. Random phases gives an homogeneous statistics, without destroying the high anisotropy introduced by the chosen forced wavenumber. We simulated the RKF at resolution $`256^3`$ and collected up to $`70`$ eddy turn over times.

Anisotropy is studied by means of $`SO(3)`$ decomposition of longitudinal structure functions :

$$S_p\left(𝑹\right)=\left[\left(𝒗\left(𝒙+𝑹\right)𝒗\left(𝒙\right)\right)\widehat{𝑹}\right]^p,$$
(1)
where we have kept only the dependency on $`𝑹`$ and neglected the small non-homogeneous fluctuations. We expect that the undecomposed structure functions are not the real “scaling” bricks of the theory. Theoretical and numerical analysis showed  that one must first decompose the structure functions onto irreducible representations of the rotational group and then study the scaling behavior of the projections. In practice, being the longitudinal structure functions scalar objects, their decomposition reduces to the projections on the spherical harmonics:

$$S_p(𝑹)=\underset{j=0}{\overset{\mathrm{}}{}}\underset{m=j}{\overset{j}{}}S_p^{jm}\left(\left|𝑹\right|\right)Y_{jm}(\widehat{𝑹}).$$
(2)
As usual, we use indexes $`(j,m)`$ to label, respectively, the total angular momentum and its projection along a reference axis, say $`\widehat{𝒛}`$. The whole physical information is hidden in the functions $`S_p^{jm}(R)`$. In particular, the main question we want to address here concerns their scaling properties:

$$S_p^{jm}(|𝑹|)A_{jm}|𝑹|^{\zeta ^j(p)}.$$
(3)
First, we need an estimate for the “dimensional” values of the exponents $`\zeta _d^j(p)`$ in all sectors. Our argument is based on the idea that large-scale energy pumping and/or boundary conditions are such as to enforce a large-scale, anisotropic, driving velocity field $`𝑼`$. Dimensional predictions for intermediate (small) scales anisotropic fluctuations may then be obtained by studying the influence of the large-scale $`𝑼`$ on the inertial range. Let us therefore evaluate the weight of anisotropic contributions as it comes out from a balance between inertial advection of the small scales and the “shear effect”, induced by the instantaneous large scale velocity configuration. Decomposing the velocity field in a small scale component, $`𝒖`$, and a large scale, strongly anisotropic component, $`𝑼`$, one finds the following equation for the time evolution of $`𝒖`$:

$$_tu_i+u_k_ku_i+U_k_ku_i+u_k_kU_i=_ip+\nu \mathrm{\Delta }u_i.$$
(4)
The major effect of the large-scale field is given by the instantaneous shear $`𝒮_{ik}=_kU_i`$ which acts as an anisotropic forcing term on small scales. A simple dimensional reasoning can be done as follows. Let us first consider the equation of motion for two point quantities $`u_l(𝒙^{})u_i(𝒙)`$ in the stationary regime; we may balance inertial terms and shear-induced terms as follows:

$$u_l(𝒙^{})u_k(𝒙)_ku_i(𝒙)𝒮_{ik}(𝒙)u_l(𝒙^{})u_k(𝒙),$$
(5)
which allows for a dimensional estimate of the anisotropic components of the LHS in terms of the RHS shear intensity and of the $`uu`$ isotropic part. Similarly for three point quantities we have (neglecting tensorial details): $`uuuu𝒮uuu`$ which can be easily generalized to any order velocity correlations. The shear term is a large-scale “slow” quantity and therefore, as far as scaling properties are concerned, we may safely estimate: $`𝒮_{ik}(𝒙)u_l(𝒙^{})u_k(𝒙)D_{ik}u_l(𝒙^{})u_k(𝒙)`$. Here the matrix $`D_{ik}`$ is associated to the combined probability to have a given shear and a given small scale velocity configuration. Clearly the $`D_{ik}`$ tensor brings angular momenta only up to $`j=2`$. One may therefore argue, by using simple composition of angular momenta, the following dimensional matching :

$$\stackrel{j}{S_p^j(R)}\stackrel{2j2}{R𝒮S_{p1}^{j2}(R)},$$
(6)
where on top of each term we have written the total angular momentum of that contribution. $`S_p^j(R)`$ is the projection on the $`j`$-th sector of the $`p`$-th order correlation function at scale $`R`$ (see equation 2), and $`𝒮`$ is the intensity of the shear term, $`D_{ik}`$, in the $`j=2`$ sector.

For instance, the leading behaviour of the $`j=2`$ anisotropic sector of the $`3`$-th order correlation function in the LHS of (5) is given by the coupling between the $`j=2`$ components of $`D_{ik}`$ and the $`j=0`$ sector of the $`2`$-th order velocity correlation in the RHS of (5): $`S_3^2(R)R𝒮S_2^0(R)R^{\zeta _d^{j=2}(3)}`$. By using the same argument and considering that now we know the scaling of $`j=0`$ and $`j=2`$ sectors of the third order correlation, we may estimate the scaling exponents of the fourth order correlation for $`j=2,4`$. From equation (6), we have the dimensional matching in the $`j=2`$ sector: $`S_4^2(R)R𝒮S_3^0(R)R^{\zeta _d^{j=2}(4)}`$ and in the $`j=4`$ sector: $`S_4^4(R)R𝒮S_3^2(R)R^2𝒮^2S_2^0(R)R^{\zeta _d^{j=4}(4)}`$. The procedure is easily extended to higher orders:

$`\zeta _d^{j=2}(p)`$ $`=`$ $`\zeta _d^{j=0}(p1)+1=(p+2)/3,p>2;`$ (7)
$`\zeta _d^{j=4}(p)`$ $`=`$ $`\zeta _d^{j=2}(p1)+1=(p+4)/3,p>3;`$ (8)
$`\zeta _d^{j=6}(p)`$ $`=`$ $`\zeta _d^{j=4}(p1)+1=(p+6)/3,p>4;`$ (9)
which can be summarized as

$$\zeta _d^j(p)=\frac{(p+j)}{3},$$

where intermittency effects in the isotropic sector have been neglected for simplicity. In this way, giving as input only the isotropic exponents, $`\zeta _d^{j=0}(p)`$, we are able to predict the scaling exponents up to $`j=2`$ for the third order structure functions, to $`j=4`$ for the fourth order, to $`j=6`$ for the fifth order and so on. We may do a little better by giving a prediction also for anisotropic fluctuations of second order correlation functions. This cannot be simply obtained by using the equations of motion, because the first one involving velocity correlations at different spatial locations, i.e. inertial range quantities, is that for $`_tu_i(x)u_j(x^{})`$, which fixes a constraint only for the third order correlation function (5). A way out is to ask the second order anisotropic fluctuations to be analytic in the shear intensity, $`𝒮`$, consistently with what one finds for higher order structure functions by the above dimensional estimate. With this assumption, we recover for $`j=2`$ Lumley prediction , $`\zeta _d^{j=2}(2)=4/3`$ by simply writing the first two terms dimensionally consistent with an expansion in the shear intensity: $`uu(ϵR)^{2/3}+𝒮R^{4/3}+\mathrm{}`$ where the first corresponds to the isotropic scaling, while the second captures anisotropies up to $`j=2`$ (higher $`j`$-sectors could be captured by adding other terms in the expansions). By using this argument, we may now remove the limit of validity of the dimensional prediction, (7-9), and extend it to all $`p`$ values.
We now come to our numerical results for the $`SO(3)`$ decomposition of longitudinal structure functions. In Figure 1, we present for the $`4`$-th order longitudinal structure function, an overview of all sectors $`(j,m)`$ which have a signal-to-noise ratio high enough to ensure stable results. Sectors with odd $`j`$s are absent due to the parity symmetry of our observable. We measured anisotropic fluctuations up to $`j=6`$. Scaling exponents can be measured in almost all sectors except for $`j=2`$ where an annoying oscillation in the sign of $`S_4^{2,m}(|R|)`$ prevents us from giving a quantitative statement.
We notice, as it is also summarized in Figure 2, that all measured exponents show a clear departure from the dimensional prediction. For example we measure in the $`j=4`$ sectors the values: $`\zeta ^4(2)=1.65(5)`$, $`\zeta ^4(4)=2.20(5)`$, $`\zeta ^4(6)=2.55(10)`$, and in the $`j=6`$ sector: $`\zeta ^6(2)=3.2(2)`$, $`\zeta ^6(4)=3.1(2)`$, $`\zeta ^6(6)=3.3(2)`$. This is a first clear sign that anisotropic scaling exponents are intermittent.
The importance of being anomalous does not stand on the exact values of the exponents, but on the connection between anomalous scaling and universality. Indeed, if correlation functions in the inertial range are not given by a dimensional matching with the large-scale shear, it means that they are fixed only by the inertial part of the Navier-Stokes evolution. In other words, they should enjoy strong universality properties with respect to changes of the large-scale physics, similarly to what happens to “zero-modes” responsible for anomalous scaling and universality in linear hydrodynamical problems .
Such a statement can even be tested in a different way.

We have taken the stationary configurations of the RKF and randomly re-shuffled all velocity phases : $`\widehat{u}_i(𝒌)P_{il}(𝒌)\widehat{u}_l(𝒌)exp(i\theta _l(𝒌))`$, where $`P_{il}(𝒌)`$ is the incompressibility projector and $`\theta _l(𝒌)=\theta _l(𝒌)`$. In this way we expect to filter out the dominant intermittent fluctuations coming from the inertial evolution, or at least those intermittent contributions connected to non-trivial phase organization. The rationale of the above statement comes from the observation that anomalous scaling, in linearly advected hydrodynamical models, is connected to the existence of statistically preserved structures with highly complex geometrical properties .

We imagine that once canceled the anomalous fluctuations, the subdominant fluctuations, due to the dimensional balancing with the forcing-shear terms, should show up. Still, it is worth to remark, the statistics of the velocity field stays non-gaussian.
In Figure 3 we show the results for the decomposition of $`4`$-th order structure functions (after phase randomization) in the $`j=4`$ anisotropic sector. As it can be seen, scaling properties change significantly going from the anomalous value (before randomization) to the dimensional predictions (after randomization).
This happens for all sectors and moments we have measured, as it is summarized in Figure 2, with the notable exception of the second order structure function where phase randomization has almost no effect. An interesting fact which can have two explanations. Phases randomization is not enough to completely filter out intermittency, especially for two points quantities which should be less sensible to phase correlation. Or, as noticed before, because second order correlation function is not constrained by any equation of motion, dimensional scaling may never exist for it even not as a sub-leading contribution. This is an important point which certainly deserves further numerical, experimental tests.

In conclusions we have presented a dimensional argument able to predict, by means of a matching between inertial correlations and shear-induced inertial terms, scaling exponents for all structure functions in any anisotropic sector. We have shown by a direct numerical simulation that anisotropic scaling exponents deviate from the previous dimensional prediction, showing anomalous values. When performing a random re-shuffling of all velocity phases, the dimensional scaling comes out as a sub-leading contribution. Everything points toward the conclusion that anisotropic fluctuations are anomalous and universal.

We thank I. Arad, M. Cencini and M. Vergassola for useful comments. This research was supported in part by the EU under the Grant No. HPRN-CT 2000-00162 “Non Ideal Turbulence” and by the INFM (Iniziativa di Calcolo Parallelo).


text:
# Positivity bounds on parton distributions at large 𝑁_𝑐

## 1 Introduction

The spectrum of the applications of the $`1/N_c`$ expansion  to the analysis of the dynamics of strong interactions is rather rich. However, the parton distributions of hadrons remained a sort of exception for a long time. Historically the large $`N_c`$ behavior of parton distributions was first analyzed in the context of the chiral quark soliton model  but the general structure of the $`N_c`$ counting in this model is the same as in QCD.

Among recent applications of the large $`N_c`$ expansion to parton distributions one should mention the phenomenological analysis of the antiquark polarized distribution , of the polarized gluon distribution  and the large $`N_c`$ relations between the quark distributions in the nucleon and in the $`\mathrm{\Delta }`$ resonance . Here we shall consider another interesting result: combining the large $`N_c`$ expansion with the probabilistic interpretation of parton distributions one can derive certain inequalities for parton distributions  which are *stronger* than the corresponding inequalities known in the case of finite $`N_c`$.

Our analysis will be based on the standard picture of the large $`N_c`$ baryons  as states that can be described by mean field (Hartree) equations of some effective theory. Although the exact action of this effective theory equivalent to the large $`N_c`$ QCD is not known certain conclusions can be made using only the assumption about the symmetry  of the solutions of the corresponding Hartree equations.

The large $`N_c`$ inequalities described here are not mere products of the large $`N_c`$ pure art per se. For example, they include the enhanced large $`N_c`$ version of Soffer inequality  for the transversity distribution that can be of interest for the phenomenological applications. Actually some of the results considered here were noticed earlier  in the context of the quark-soliton model . Since the derivation was based only on the spin-flavor symmetry of the large $`N_c`$ baryons the generalization to the case of the large $`N_c`$ QCD is rather straightforward.

## 2 Positivity bounds at finite $`N_c`$

Let us first recall the well-known positivity properties of parton distributions in the real world at finite number of colors. The unpolarized parton distribution $`q_f(x)`$ is defined as the probability to find a quark of type $`f`$ carrying fraction $`x`$ of the nucleon momentum in the infinite momentum frame. The positivity of this probability leads us to the simplest inequality $`q_f(x)0`$. If one considers the longitudinally polarized nucleon then $`q_f(x)`$ can be decomposed into two pieces coming from quarks with the longitudinal spin collinear with nucleon momentum ($`q_f^{}(x)`$) and with anticollinear spin ($`q_f^{}(x)`$)

$$q_f(x)=q_f^{}(x)+q_f^{}(x).$$
(1)
Since both functions $`q_f^{}(x)`$ and $`q_f^{}(x)`$ have a probabilistic interpretation they are positive

$$q_f^{}(x)0,q_f^{}(x)0.$$
(2)
The (longitudinally) polarized parton distribution is defined as

$$\mathrm{\Delta }_Lq_f(x)=q_f^{}(x)q_f^{}(x).$$
(3)
Now it follows from (1), (2), (3)

$$|\mathrm{\Delta }_Lq_f(x)|q_f.$$
(4)
Next, one can consider the nucleon and quarks with the transverse polarization and define by analogy with $`\mathrm{\Delta }_Lq_f`$ the transversity distribution $`\mathrm{\Delta }_Tq_f`$ which will obey inequality

$$|\mathrm{\Delta }_Tq_f(x)|q_f.$$
(5)
One more inequality can be obtained by choosing arbitrary directions for the polarizations of quarks and nucleons. The positivity of the corresponding probabilities leads to Soffer inequality

$$q_f+\mathrm{\Delta }_Lq_f2|\mathrm{\Delta }_Tq_f|.$$
(6)
## 3 Positivity bounds at large $`N_c`$

Thus in the real world with $`N_c=3`$ (and at any other finite $`N_c`$) we have the following set of the positivity bounds for the unpolarized ($`q_f`$), longitudinally polarized ($`\mathrm{\Delta }_Lq_f`$) and transversity ($`\mathrm{\Delta }_Tq_f`$) quark distributions

$$\begin{array}{c}|\mathrm{\Delta }_Lq_f|q_f,\\ |\mathrm{\Delta }_Tq_f|q_f,\\ 2|\mathrm{\Delta }_Tq_f|q_f+\mathrm{\Delta }_Lq_f.\end{array}(\mathrm{finite}N_c)$$
(7)
Below we shall derive the following analogs of these inequalities in the large $`N_c`$ limit:

$$\begin{array}{c}\left|\mathrm{\Delta }_Lq_f\right|\frac{1}{3}q_f\left[1+O(N_c^1)\right],\hfill \\ |\mathrm{\Delta }_Tq_f|\frac{1}{2}q_f\left[1+O(N_c^1)\right],\hfill \\ 2|\mathrm{\Delta }_Tq_f|\left(\frac{1}{3}q_f+\mathrm{\Delta }_Lq_f\right)\left[1+O(N_c^1)\right].\hfill \end{array}(N_c\mathrm{})$$
(8)
Note that these inequalities for the large $`N_c`$ case are *stronger* than the corresponding inequalities at finite $`N_c`$ (7) because of the factors $`1/3`$ and $`1/2`$ appearing in the right-hand sides. The origin of these factors can be traced to the spin-flavor symmetry of the large $`N_c`$ baryons; in particular, the factor of $`1/3`$ has nothing to do with the fact that $`1/N_c=1/3`$ in the real world.

This enhancement of the positivity bounds in the large $`N_c`$ limit is rather a striking result. Naively one could expect that in the leading order of the large $`N_c`$ limit some part of information about parton distributions would be lost but instead in this limit we have new stronger bounds. The reason of this phenomenon can be explained in simple terms without appealing to the technical derivation of (8). In the leading order of the large $`N_c`$ approximation the nucleon is degenerate with the $`\mathrm{\Delta }`$ resonance and other baryonic resonances with $`T=J=1/2,3/2,5/2,\mathrm{}`$ For all of these particles one can analyze the twist-two parton distributions and derive the corresponding positivity bounds. Moreover, since all these particles are degenerate at the level of the leading order of the $`1/N_c`$ expansion one can consider the formal superposition of states with different values of $`T=J`$ and derive the positivity bounds for the parton distribution function of these formal superpositions. Although the superpositions of states with different $`T=J`$ make no physical sense this formal trick allows us to derive the new enhanced positivity bounds (8) valid in the large $`N_c`$ limit.

## 4 Large $`N_c`$ counting for parton distributions

In order to derive the large $`N_c`$ inequalities (8) we should start from the analysis of the large $`N_c`$ behavior of the parton distributions. This analysis can be performed using standard large $`N_c`$ counting rules of QCD (see e.g. ). On the other hand, it is well known that the same large $`N_c`$ behavior of various quantities can be also extracted from the nonrelativistic quark model, Skyrme model, quark soliton model etc. Below we briefly describe the large $`N_c`$ behavior of parton distributions appealing mainly to the naive quark model and using sum rules for parton distributions.

Note the twist-two parton distributions $`q_f(x),\mathrm{\Delta }_Lq_f(x),\mathrm{\Delta }_Tq_f(x)`$ depend on Bjorken variable $`x`$ so that the analysis of the large $`N_c`$ behavior should involve a convention about the behavior of $`x`$ in this limit. From the physical point of view the nucleon consists of $`N_c`$ quarks (with an admixture of quark-antiquark pairs) so that it is natural to expect that at large $`N_c`$ limit the nucleon momentum in the infinite momentum frame is more or less uniformly distributed between $`O(N_c)`$ quarks and antiquarks. Then in the large $`N_c`$ limit it makes sense to concentrate on $`x=O(1/N_c)`$. This physical argument shows that that the large $`N_c`$ behavior of various parton distributions should have the form

$$q_f(x)=N_c^2\left[\varphi _f(N_cx)+O(N_c^1)\right],$$
(9)
$$\mathrm{\Delta }_{L,T}q_f(x)=N_c^2\left[\mathrm{\Delta }_{L,T}\varphi _f(N_cx)+O(N_c^1)\right].$$
(10)
Functions $`\varphi _f,\mathrm{\Delta }_L\varphi _f,\mathrm{\Delta }_T\varphi _f`$ are stable in the large $`N_c`$ limit so that the $`x`$ dependence of the form $`\varphi _f(N_cx)`$ implies that $`xN_c^1`$ at large $`x`$. The factor of $`N_c^2`$ ensures the consistency of the large $`N_c`$ behavior for various sum rules for parton distributions. For example for the baryon charge sum rule

$$\underset{f}{}\underset{0}{\overset{1}{}}𝑑x\left[q_f(x)\overline{q}_f(x)\right]=\underset{f}{}\underset{1}{\overset{1}{}}𝑑xq_f(x)=N_c$$
(11)
we find in the large $`N_c`$ limit

$$\underset{1}{\overset{1}{}}𝑑xq_f(x)\stackrel{N_c\mathrm{}}{}N_c^2\underset{1}{\overset{1}{}}𝑑x\varphi _f(N_cx),$$

$$\stackrel{N_c\mathrm{}}{}N_c\underset{N_c}{\overset{N_c}{}}𝑑y\varphi _f(y)\stackrel{N_c\mathrm{}}{}N_c\underset{\mathrm{}}{\overset{\mathrm{}}{}}𝑑y\varphi _f(y)=O(N_c)$$
(12)
which really agrees with the sum rule (11).

In the leading order of the large $`N_c`$ limit we have for the $`u`$ and $`d`$ distribution functions of the proton :

$$q_u=q_d\left[1+O(N_c^1)\right],\mathrm{\Delta }_Lq_u=\mathrm{\Delta }_Lq_d\left[1+O(N_c^1)\right],$$

$$\mathrm{\Delta }_Tq_u=\mathrm{\Delta }_Tq_d\left[1+O(N_c^1)\right].$$
(13)
These large $`N_c`$ relations have a simple interpretation in the nonrelativistic quark model: at large odd $`N_c`$ the proton consists of $`(N_c+1)/2`$ $`u`$-quarks and $`(N_c1)/2`$ $`d`$-quarks. In the leading order of large $`N_c`$ we have the same amount of $`u`$ and $`d`$ quarks which means that $`q_u=q_d\left[1+O(1/N_c)\right]`$. Next, since we want to have the spin of the proton $`1/2`$, at large $`N_c`$ there must be a cancelation between the $`O(N_c)`$ parts of the spin of $`u`$ and $`d`$ quarks, which is expressed by the last two equations (13).

## 5 Derivation of large $`N_c`$ inequalities

Below we sketch the main ideas of the derivation of inequalities (8). The complete description can be found in ref. .

Although we cannot compute parton distributions in the large $`N_c`$ QCD, still it is possible to write certain general relations using only the spin-flavor symmetry of the large $`N_c`$ baryons. The twist-two quark distributions $`q_u,\mathrm{\Delta }_Lq_u,\mathrm{\Delta }_Tq_u`$ in proton can be represented in the form

$`q_u(x)`$ $`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{Sp}\rho (x),`$
$`\mathrm{\Delta }_Lq_u(x)`$ $`=`$ $`{\displaystyle \frac{1}{6}}\mathrm{Sp}\left[\gamma ^5\tau ^3\rho (x)\right],`$ (14)
$`\mathrm{\Delta }_Tq_u(x)`$ $`=`$ $`{\displaystyle \frac{1}{12}}\mathrm{Sp}\left[\gamma ^5\left(\gamma ^1\tau ^1+\gamma ^2\tau ^2\right)\rho (x)\right].`$
Here $`\rho _{f^{}s^{},fs}(x)`$ is an $`8\times 8`$ matrix with $`SU(2)`$ flavor indices $`f^{},f`$ and with Dirac spin indices $`s^{},s`$. Matrix $`\rho (x)`$ depending on $`x`$ is determined by the dynamics of the large $`N_c`$ QCD and is not known. Nevertheless we are aware of the following general properties:

1) $`\rho `$ is a hermitean positive matrix (i.e. all eigenvalues are positive or zero)

$$\rho ^+=\rho ,$$
(15)
$$\rho 0,$$
(16)
2) $`\rho `$ lives in the subspace of the projector $`(1+\gamma ^0\gamma ^3)/2`$

$$\rho (\gamma ^0\gamma ^3)=(\gamma ^0\gamma ^3)\rho =\rho ,$$
(17)
3) $`\rho `$ commutes with $`i\gamma ^1\gamma ^2+\tau ^3i\gamma ^1\gamma ^21_{flavor}+1_{spin}\tau _3`$

$$[\rho ,(i\gamma ^1\gamma ^2+\tau ^3)]=0.$$
(18)
These properties allow us to write the following general representation for $`\rho (x)`$

$$\rho (x)=\frac{1+\gamma ^0\gamma ^3}{2}\left[c_1(x)\mathrm{𝟏}+c_2(x)\gamma ^5\tau ^3+c_3(x)\gamma ^5(\gamma ^1\tau ^1+\gamma ^2\tau ^2)\right]$$
(19)
where $`c_i(x)`$ are some scalar functions of $`x`$.

The positivity property (16) leads to the following constraints on coefficients $`c_i`$

$$c_1c_22|c_3|,$$
(20)
$$c_1+c_20.$$
(21)
Inserting the explicit representation for $`\rho `$ (19) into the general equations for parton distributions (14) we immediately express the parton distributions $`f_i(x)`$ in terms the coefficients $`c_i(x)`$. Then the constraints (20), (21) on these coefficients immediately lead to bounds (8) on parton distributions.

## 6 Phenomenological status

Concerning the practical applications of the large $`N_c`$ inequalities we have to keep in mind two competing factors:

1) large $`N_c`$ inequalities (8) are stronger than the standard inequalities (7) valid at finite $`N_c`$,

2) at finite $`N_c`$ inequalities (8) can be violated by large $`1/N_c`$ corrections.

In the case of $`\mathrm{\Delta }_Lq_f`$ the $`1/N_c`$ corrections are known to be rather large. Indeed, the Bjorken sum rule relates $`\mathrm{\Delta }_Lq_u\mathrm{\Delta }_Lq_d`$ with the axial constant $`g_A`$ and it is well known that $`g_A`$ is usually underestimated in various model calculations based on the leading order of the $`1/N_c`$ expansion, which is usually attributed to sizable $`1/N_c`$ corrections. Hence the $`1/N_c`$ corrections to $`\mathrm{\Delta }_Lq_f`$ should be also large. Therefore one should not wonder that our large $`N_c`$ inequality for $`\mathrm{\Delta }_Lq_u\mathrm{\Delta }_Lq_d`$

$$\frac{\left|\mathrm{\Delta }_Lq_u\mathrm{\Delta }_Lq_d\right|}{q_u+q_d}\frac{1}{3}\left[1+O(N_c^1)\right]$$
(22)
following from (8), (13) is in conflict with experimental data. For example, in the GRSV parametrization

$$\underset{x}{\mathrm{max}}\frac{\left|\mathrm{\Delta }_Lq_u\mathrm{\Delta }_Lq_d\right|}{q_u+q_d}0.6.$$
(23)
On the other hand, the large $`N_c`$ bounds (7) can be of use in the case of the transversity distribution and polarized antiquark distribution where the experimental knowledge is rather poor.

## 7 Concluding remarks

In the large $`N_c`$ limit one has the positivity bounds (7) for twist-two parton distributions that are stronger than the standard positivity bounds at finite $`N_c`$. However, the practical applications of the large $`N_c`$ inequalities are sensitive to the size of the $`1/N_c`$ corrections.

In the case of finite $`N_c`$ it is well known that positivity bounds (7) are compatible with the DGLAP evolution equations in the following sense. If inequalities (7) hold at some normalization point then the one-loop evolution to higher normalization points will not violate these inequalities. In ref.  it was shown that the large $`N_c`$ inequalities (8) have the same property: if one uses the evolution kernels only in the leading order of the large $`N_c`$ expansion then inequalities (8) are preserved by the evolution upwards.

Acknowledgements.

The results described here are obtained in collaboration with M.V. Polyakov. I am thankful to A. Belitsky, D.I. Diakonov, A.V. Efremov, L. Frankfurt, K. Goeke, N. Kivel, V.Yu. Petrov, A.V. Radyushkin, A. Shuvaev, M. Strikman and O. Teryaev for the discussions. It is a pleasure to thank R. Lebed for the organization of the workshop. The author appreciates the support of INT, DFG and BMBF.


text:
# Electron-phonon interaction in the solid form of the smallest fullerene C20

## Abstract

The electron-phonon coupling of a theoretically devised carbon phase made by assembling the smallest fullerenes C20 is calculated from first principles. The structure consists of C20 cages in an fcc lattice interlinked by two bridging carbon atoms in the interstitial tetrahedral sites (fcc-C22). The crystal is insulating but can be made metallic by doping with interstitial alkali atoms. In the compound NaC22 the calculated coupling constant $`\lambda /N(0)`$ is 0.28 eV, a value much larger than in C60, as expected from the larger curvature of C20. On the basis of the McMillan’s formula, the calculated $`\lambda `$=1.12 and a $`\mu ^{}`$ assumed in the range 0.3-0.1 a superconducting Tc in the range 15-55 K is predicted.

The discovery of high superconducting transition temperature in alkali-doped C60 ($`T_c`$=10-40 K) has attracted much attention on these materials as representatives of a new class of superconductors C60 ; gunnarsson . Superconductivity, albeit at lower temperature ($`T_c`$=7 K), has been recently observed in another member of this class, C70, which indicates the possibility that other fullerenes might also become superconducting batlogg\_C70 . In this respect, the smaller fullerenes are the most interesting since the electron-phonon coupling (and hopefully $`T_c`$) has been calculated to become larger with increasing curvature of the molecule tomanek ; cote2 ; lannoo ; colo . Therefore, the recent synthesis of the smallest fullerene, having the largest curvature, C20, starting from the stable C20H20 molecules prinzbach  represents a step further along the way towards high-$`T_c`$ superconductivity in carbon-based materials. Experimental evidence for the realization of the next step, the formation of a crystalline solid based on C20 clusters, has been recently provided as well Iqbal . Crystalline micro-domains of the new carbon phase have been detected by Trasmission Electron Microscopy (TEM) in carbon films deposited on a nickel substrate by the ultraviolet laser ablation of diamond. A face-centered-cubic (fcc) symmetry have been inferred from TEM data. The Raman spectrum shows a strong peak at 1618 cm-1 typical of the stretching of the ethylenic C=C double bond and the presence of C20 clusters as building blocks of the new crystal has been determined by mass spectroscopy of laser-desorbed fragments Iqbal . Based on this information we have theoretically conceived a C20-based solid whose calculated structural and vibrational properties are suitable to reproduce the experimental TEM and Raman data nostro . The structure comprises C20 cages in a fcc lattice interconnected by two bridging atoms per unit cell in the tetrahedral interstitial sites (Fig. 1). Hereafter it is referred to as fcc-C22, the unit cell containing 22 atoms.

In this letter, we report on the ab-initio study of the electronic properties of fcc-C22 doped with alkali metals (Na, Li) in the search for metallization and superconductivity. The calculated electron-phonon coupling constant $`\lambda `$ for the alkali-doped fcc-C22 is found to be substantially enhanced compared to C60 leading to the possibility of larger $`T_c`$ than in alkali-doped fullerite.

We use density functional theory in the local density approximation (LDA), norm-conserving pseudopotentials and plane-wave expansion of the Kohn-Sham orbitals up to a kinetic cutoff of 40 Ry notaLDA . Non-linear-core corrections have been added for the alkali ions. Phonon frequencies and electron-phonon coupling constants have been computed within the linear response theory RMP\_baroni ; Pwscf . This approach has been successfully applied to the calculation of phonon dispersion relations in insulator and metals RMP\_baroni ; degironc  and of the electron-phonon coupling constant $`\lambda `$ in simple and transition metals pavone . Brillouin Zone (BZ) integrations are performed using Hermite-Gaussian smearing of order one degironc  with a linewidth of 0.01 Ry and 18 special points in the irreducible wedge of the BZ (IBZ).

The calculation of the equation of state of fcc-C22 yields a cohesive energy of 0.63 eV/atom lower than that of diamond and a lattice constant of 8.61 $`\mathrm{\AA }`$. The density is 1.78 g/cm3, close to that of fullerite (1.68 g/cm3), but the bulk modulus is as large as 239 GPa due to the presence of covalent interlinks between adjacent C20 cages. The crystal belongs to the Fm3̄ space group and contains only three independent atoms in the unit cells, the other being obtained by symmetry operations. Their coordinates in units of the lattice constant are $`A`$ (0.25, 0.25, 0.25), $`B`$ (0.148, -0.148, 0.148), $`C`$ (0.078, -0.217, 0.0) nostro ; notasaito  (Fig. 1). Eight atoms out of twenty of the cage (atoms $`B`$) and the interstitial carbon atoms ($`A`$) are sp3 hydridized. The other twelve atoms of the cage ($`C`$) are sp2 hybridized forming C=C ethylenic bonds. There are three bond lengths in the crystal: the single sp3-like $`BC`$ bond in the pentagonal ring (1.531 $`\mathrm{\AA }`$), the C=C ethylenic-like bond (1.346 $`\mathrm{\AA }`$) long and the single sp3-like bond between a cage (equivalent to $`B`$) and an interstitial atom ($`A`$) (1.517 $`\mathrm{\AA }`$).

The electronic band structure of fcc-C22 is reported in Fig. 2. The system is an insulator with a slightly indirect $`XL`$ band gap of 2.47 eV. Projection of the Kohn-Sham states on the minimal atomic basis set shows that the lowest conduction bands and the highest valence bands are mainly formed by $`\pi ^{}`$ and $`\pi `$ states of the ethylenic-like bonds.

It is suggestive to wonder if doping might turn fcc-C22 metallic by electron transfer into the molecular-like lowest conduction bands, as occurs for instance in fullerite. The analogy with C60 is evident by considering that both C60 molecules in fullerite and C20 cages in fcc-C22 form a fcc lattice. In the A3C60 (A=K, Cs, Rb) compounds the three alkali metal ions per unit cell occupy the two tetrahedral sites and the octahedral site gunnarsson . In fcc-C22 the two tetrahedral sites are occupied by the interstitial carbon atoms, but the octahedral site is empty and might host a guest ion. In the search for metallization we have inserted a Na atom per unit cell in the octahedral site. Calculation of the equation of state of the compound NaC22 shows an increase in the lattice parameter ($`a`$=8.765 $`\mathrm{\AA }`$) with respect to the pure phase and a larger (2 $`\%`$) increase in the C=C bond length due to population of antibonding $`\pi ^{}`$ states. In fact, Na is ionized by donating its outermost electron to the host which becomes metallic. The Na-host interaction is weak and the band structure of NaC22 is very close to that of the pure system, but for the shift of the Fermi level inside the lowest conduction bands (Fig. 2b). Therefore the rigid band approximation holds to a large extent in the description of Na insertion in fcc-C22. Although the lowest conduction band is 3-fold degenerate at the $`\mathrm{\Gamma }`$ point (Fig. 2a), the Fm3̄ symmetry is preserved in NaC22 and no static Jahn-Teller distorsion is observed. This is due to the fact that the minima of the lowest conduction bands are located at the $`X`$ and $`L`$ points where the degeneracy is already lifted.

To study the superconducting properties we have computed the electron phonon coupling constant $`\lambda `$ as

$`{\displaystyle \frac{\lambda }{N(0)}}`$ $`=`$ $`{\displaystyle \underset{\alpha }{}}{\displaystyle \frac{d𝐪}{\mathrm{\Omega }_{\mathrm{BZ}}}\frac{\lambda _{\alpha 𝐪}}{N(0)}}`$
$`=`$ $`{\displaystyle \underset{\alpha }{}}{\displaystyle \frac{d𝐪}{\mathrm{\Omega }_{\mathrm{BZ}}}\frac{1}{\omega _{\alpha 𝐪}^2}\underset{n,m}{}\frac{d𝐤}{\mathrm{\Omega }_{BZ}}\frac{\delta (E_{𝐤n})\delta (E_{𝐤+𝐪m})}{N^2(0)}u_{𝐤+𝐪,m}𝐌^{\frac{1}{2}}ϵ_{\alpha 𝐪}V_{eff}^𝐪u_{𝐤,n}^2}`$
where $`N(0)`$ is the electronic density of states at the Fermi level (zero of energy) per spin and per cell, $`\lambda _{\alpha 𝐪}/N(0)`$ is the partial electron-phonon coupling due the $`\alpha `$-th phonon at the point $`𝐪`$ in the BZ (of volume $`\mathrm{\Omega }_{\mathrm{BZ}}`$), $`ϵ_{\alpha 𝐪}`$ is the phonon polarization vector, $`𝐌`$ is the atomic mass matrix, $`n`$ and $`m`$ are electronic band indices, $`u_{𝐤,n}`$ is the periodic part of the Kohn-Sham state with energy $`E_{𝐤n}`$, and $`V_{eff}^𝐪`$ is the derivative of the Kohn-Sham effective potential with respect to the atomic displacement caused by a phonon with wavevector $`𝐪`$. We replaced the $`\delta `$-functions in Eq. Electron-phonon interaction in the solid form of the smallest fullerene C20 with Hermite-Gaussian forms of 0.01 Ry linewidth. The BZ integration over $`𝐤`$ is performed on a grid of 18 special points in the IBZ. We computed $`\lambda _{\alpha 𝐪}/N(0)`$ for $`𝐪`$ along the high symmetry directions  and . The integration over $`𝐪`$ in Eq. Electron-phonon interaction in the solid form of the smallest fullerene C20 is performed over a uniform mesh in the IBZ by assigning to each point the value of $`\lambda _{\alpha 𝐪}/N(0)`$ computed for the closest point along a high symmetry directions nota\_converg . Finally, for NaC22 we obtain $`\lambda /N(0)`$=0.28 eV, $`N(0)`$=4.0 states/eV/spin/cell (as computed by the tetrahedron method with 213 points in the IBZ) and $`\lambda `$=1.12. A similar calculation for the charged fullerite C$`{}_{60}{}^{}{}_{}{}^{3}`$ at the $`\mathrm{\Gamma }`$ point only in the molecular-like approximation nota\_molecular  gives $`\lambda /N(0)`$=0.076 eV, in agreement with previous works gunnarsson  (or 0.064 eV by including the set of phonons considered in table IV of Gunnarsson’s review gunnarsson ). Therefore, although fcc-C22 is not strictly a fullerite, the C20 clusters being covalently interlinked via interstitial atoms, the solid retains some features of its molecular building blocks, giving rise to a substantial enhancement of $`\lambda /N(0)`$ with respect to C60 as expected from the larger curvature of C20 lannoo . The persistence of molecular features in the crystal is also confirmed by the calculation of $`\lambda /N(0)`$ in the molecular-like approximation nota\_molecular  for the pure fcc-C22 doped with an additional electron (equally shared by the three lowest $`t_u`$ states). This latter calculation gives $`\lambda /N(0)`$=0.31 eV, a value very close to the result for NaC22 at convergence in BZ integration (0.28 eV) nota\_litio . A tight-binding calculation of the electron-phonon interaction in the isolated C20 cluster gave a somehow lower value $`\lambda /N(0)=0.18`$ eV lannoo  which may be partially due to the semiempirical method used.

In order to identify the phonons which mostly contribute to the elctron-phonon coupling, we report in Table I the partial electron-phonon coupling constants $`\lambda _\alpha /N(0)`$ of the individual $`\alpha `$-th phononic band, obtained by averaging $`\lambda _{\alpha 𝐪}/N(0)`$ over the BZ. The frequency and symmetry of the $`\alpha `$-th phononic bands at the $`\mathrm{\Gamma }`$ point are also given in Table I for the pure fcc-C22 and NaC22 phases. Although the phonons with the largest matrix elements in Eq. Electron-phonon interaction in the solid form of the smallest fullerene C20 correspond to the stretching modes of the C=C ethylenic bonds (1359-1362 cm-1), their contribution to $`\lambda `$ is attenuated by the factor $`1/\omega ^2`$. As a consequence the largest $`\lambda _\alpha /N(0)`$ is due to the $`E_g`$ mode at 334 cm-1 (cfr. Table I); its displacement pattern (Fig. 3) produces the largest modulation of the dihedral angles which control the $`\pi `$ conjugation of the ethylenic bonds. This is at variance with C60 where the largest electron-phonon coupling constants are due to the phonons with the highest frequencies gunnarsson .

In system where the superconductivity is driven by the electron-phonon interaction an estimate of $`T_c`$ can be obtained by making use of the McMillan’s solution carbotte ; allen  of the Eliashberg equations:

$$T_c=\frac{\omega _{ln}}{1.2}\mathrm{exp}\left(\frac{1.04(1+\lambda )}{\lambda \mu ^{}(1+0.62\lambda )}\right)$$
(2)
where $`\omega _{ln}`$ is the weighted logarithmic average of the phonon frequencies carbotte ; allen  (470 cm-1 in NaC22) and $`\mu ^{}`$ is the Coulomb pseudopotential $`\mu =UN(0)`$ renormalized by retardation effects. $`U`$ is a typical screened Coulomb interaction among states at the Fermi level allen . Latest estimates give $`\mu ^{}`$ 0.3 for C60 although, due to the similarity between the electronic bandwidth (0.5 eV) and phonon energies (up to 0.2 eV), the applicability of Eq. 2 is questionable for fullerite gunnarsson  . We have not attempted to estimate numerically $`\mu ^{}`$ for NaC22, so our discussion now turns on a very speculative level. Although $`U`$ might be larger in NaC22 than in fullerite due to the smaller size of C20, $`N(0)`$ turns out to be lower (4.0 and 6.6-9.8 states/eV/spin/cell in NaC22 and K3C60 from band structure calculations gunnarsson , respectively) due to the larger dispersion of the conduction bands, three times larger than in C60. The estimate of the renormalization of $`\mu `$ leading to $`\mu ^{}`$ is difficult gunnarsson , however the renormalization is perhaps larger than in C60 due to the larger bandwidth. Thus, if we assume that the value of 0.3 estimated for fullerite is an upper bound for $`\mu ^{}`$ in NaC22 we obtain $`T_c`$=15K from Eq. 2 and $`\lambda `$=1.12. For the more optimistic case of $`\mu ^{}`$ in the range 0.1-0.2 typical of simple and transition metals carbotte ; allen , we obtain the remarkable result of $`T_c`$=55-33 K ($`T_c`$ changes linearly with $`\mu ^{}`$ in this range).

Since $`\lambda /N(0)`$ is mainly a molecular property, we have also attempted to increase $`\lambda `$ by increasing $`N(0)`$ with a lattice expansion. This has been achieved by substituting C by Si in the interstitial sites producing the compound NaSi2C20. This latter system has a substantially expanded lattice parameter ($`a`$=9.49 $`\mathrm{\AA }`$, $`N(0)`$=4.44 states/eV/spin/cell) but only a marginally enhanced $`\lambda `$=1.49 with respect to NaC22. In fact the increase in $`N(0)`$ is lower than expected from the lattice expansion due to the reduction of bands at the Fermi level from two in NaC22 (cfr. Fig. 2) to one in NaSi2C20. Moreover, a decrease in $`\omega _{ln}`$ due to the lattice expansion compensates the increase in $`\lambda `$ and finally yields a $`T_c`$ close (within 10 K) to that of NaC20. However, further increase in $`\lambda `$ is foreseeable for the compounds NaGe2C20 and NaSn2C20.

In summary, we have carried out ab-initio calculations on the electronic and superconductive properties of a new carbon phase fcc-C22 formed by assembling C20 clusters in a fcc lattice. The system is insulating but can be made metallic by doping with interstitial alkali metals. Although the structure is not a molecular crystal, the C20 cages being covalently interlinked via interstitial carbon atoms, the crystal retains some molecular features. In fact, as expected from the molecular properties of the isolated C20 cluster, the compound NaC22 exhibits a large electron-phonon coupling constant $`\lambda /N(0)=0.28`$ eV, strongly enhanced with respect to C60 due to the larger curvature of C20. Estimate of $`N(0)`$ from band structure calculation yields $`\lambda =`$1.12 which allows a rough estimate of $`T_c`$ from McMillan’s formula. For $`\mu ^{}`$ in the reasonable range of 0.3-0.1 we predict $`T_c`$=15-55 K for NaC22. Although a more precise estimate of $`T_c`$ requires the difficult calculation of $`\mu ^{}`$, the results presented here shows that NaC22 is a potentially interesting superconductor. We hope that these theoretical results could stimulate further experimental work on the synthesis and characterization of this new structure.

This work is partially supported by the INFM Parallel Computing Initiative, and by MURST, through COFIN99 and COFIN01 under contract 2001021133. We gratefully aknowledge Z. Iqbal for discussion, information and for providing us with his experimental data prior to publication.