text:
# Conductivity of a 2DEG in Si/SiGe heterostructure near metal- insulator transition: role of the short and long range scattering potential

## I Introduction

According to the well established one-parameter scaling theory of conductivity  any two dimensional electron system can exist only in insulating state and no MIT should be possible. Thus, it came quite as a surprise when an apparent metallic state and an MIT were observed first in Si MOSFETs  and, later, in several other 2D systems, such as p-Si/SiGe, n- and p-GaAs/AlGaAs etc . The fact that the MITs have always been observed in systems with a large parameter $`r_s`$ (the ratio of the interaction energy to the kinetic energy) is an indication that electron-electron interaction plays an important role in the MIT. While it still remains unclear whether the observed MIT is a genuine transition and whether a zero temperature metallic state can really exist at certain levels of interaction and disorder, its discovery has stimulated an extensive experimental and theoretical work in the area. In particular, a considerable progress has been made  in our understanding of the temperature dependent corrections to the conductivity from electron-electron interaction not only in the diffusive regime $`k_BT\tau /\mathrm{}1`$, as previously, but also in the ballistic and the intermediate regimes $`k_BT\tau /\mathrm{}1`$, typical for systems showing a MIT. The theoretical results  have been reported to give a good description of the conductance temperature dependence and magneto- resistance in parallel magnetic field in several types of 2D systems on the metallic side of the apparent MIT, . It should, however, be noted, that the theory  relies on an approximation where the scattering potential is treated as point-like scatterers. While appropriate for Si MOSFETs, this approximation is definitely not true for high mobility heterostructures where the impurities are separated from the 2DEG by an undoped spacer and the scattering potential is predominantly long range. Therefore, a certain caution should be used when applying the results of Ref.  for the description of the experimental data in high mobility heterostructures.

Till recently, in the context of the new approach  no explicit discussion has been made of the interaction contribution to the magnetoresistance in a transverse magnetic field. The point has been recently addressed in Ref.  where the interaction contribution to the magnetoresistance in a classically strong transverse field is calculated for a smooth disorder and for arbitrary $`k_BT\tau /\mathrm{}`$. Soon afterwards the results of Ref.  have been tested in an experiment on magnetoresistance in a GaAs/AlGaAs heterostructure with a 2DEG and a good agreement with theory  has been found .

In the present work we attempt an experimental investigation of these questions in a Si/SiGe heterostructure with a 2DEG. Despite the abundance of information on the MIT in Si MOSFETS and different other types of heterostructures, so far there has been no observation of an MIT in a n-Si/SiGe. This may in part be explained by the difficulty in fabricating a gated n-Si/SiGe structure with a stable and controllable behavior. At the same time an experimental study of these questions in a n-Si/SiGe is expected to be of some interest. Such a study would provide an opportunity to compare the electron transport properties of a Si MOSFET and a Si/SiGe heterostructure in the vicinity of the MIT and in the “metallic” regime. This comparison may be instructive as these two silicon structures have very similar electron energy spectrum, but differ in the structure of the disorder potential, which is predominantly short range in Si MOSFETs but has an important long range component in a Si/SiGe heterostructure. In this work we report the first observation of the MIT in a Si/SiGe heterostructure with a 2D electron gas and perform an analysis of the magnetoresistance data and the conductance temperature dependence on the metallic side of the MIT using the results of Ref.  and Ref. . The theory  produces a good fit for the linear part of the conductance themperature dependence. However the value of the interaction constant $`F_0^\sigma `$ serving as the only fitting parameter turns out to be about two times smaller than that reported for Si MOSFETs with similar values of $`r_s`$, . Also, in contrast to  we find a rather big divergence between our magnetoresistance data and the behavior predicted in Ref. . We attribute all these results to the possibility that in the intermediate and ballistic regimes the case of a mixed mechanism of scattering, i.e. the situation where both the short and the long range scattering are equally important for electron transport, cannot be properly described by a model that ignores either one of the other component of disorder.

## II Experimental procedures

The n-type modulation doped Si/Si1-xGex samples used in this study were grown by MBE . On a high resistivity p-type Si substrate a 1.4$`\mu `$m graded Si/Si1-xGex layer, in which the Ge content increases linearly from 5 to 35%, was grown, followed by 0.5$`\mu `$m of constant composition Si0.7Ge0.3 alloy. These layers produced a strain-relaxed but low-dislocation-density surface onto which the tensilely strained Si channel, 110 $`\AA `$ thick, was deposited, followed by an undoped Si/Si1-xGex spacer layer, a doped Si/Si1-xGex supply layer and finally undoped Si/Si1-xGex and Si cap layers. The thickness and As doping levels in the supply layer were varied slightly between the samples. To measure the magnetoresistivity of the samples, the layers were photolithographically processed into 50 $`\mu `$m wide Hall bars with the distance between the four voltage probes on each side 100, 250 and 100 $`\mu `$m. The mesa pattern for the Hall bars was etched by reactive ion etching to define the channel geometry. The ohmic contacts to the Hall bars were prepared by ion doping of phosphorus.

The magnetotransport properties of the structures were measured at 0.4-4.2 K in magnetic fields up to 15 T using a superconducting solenoid and He3 and VTI cryostats. The data was acquired using standard a.c. lock-in techniques with a frequency around 10 Hz. Driving currents through the Hall bar of 0.1 $`\mu `$A were chosen to avoid the heating effects. All our samples behaved in a more or less similar way. In what follows we discuss the results obtained on a sample with the electron density $`n_s=(3.536.23)\times 10^{11}`$cm-2 depending on the prior illumination and the highest mobility $`\mu 6\times 10^4`$ cm2/Vs.

## III Results and Analysis

As has been mentioned in the introduction the nature of the so called MIT remains a topic of debate. Since all the experiments are necessarily conducted at finite temperatures it is impossible to say whether the metallic state observed in these experiments would persist down to zero temperature or whether it will ultimately evolve into an insulating state. Despite a large number of articles on the MIT in various types of 2D semiconductor systems, so far there has been no report of a MIT in a Si/SiGe heterostructure with a 2DEG. The MIT is typically observed by measuring the temperature dependence of a sample resistance for different carrier densities which are usually controlled by voltage applied to the electrostatic top gate. Reportedly, it has been found difficult to fabricate a n-Si/SiGe heterostructure with a gate that does not leak. The leakage from the gate to the 2DEG results in a sample instability that makes it unsuitable for experiment. In the present work, after several unsuccessful attempts to fabricate an operational gated n-Si/SiGe structure, we decided to perform our study using structures without gate. Among the several samples that we have tried there were some that, upon cooling down to low temperatures in the dark were found in high resistance low electron density states. The samples were also found to be very sensitive to illumination. Exposing them to a specially selected dose of red LED radiation permitted us to fine-tune the sample resistance with a good precision. The sample states of different resistivity thus obtained were characterized by a high stability yielding reproducible dependences of resistance versus temperature and magnetic field over a considerable, up to several hours, period of time. Shown in Fig. 1 is a series of resistance versus temperature traces for different electron densities obtained using the procedure described above. Similar to the previously reported MITs in other 2D semiconductor systems  we observe a transition from the insulator to metallic-like behavior as the electron density is varied from $`3.53\times 10^{11}`$cm-2 to $`4.87\times 10^{11}`$cm-2 . To our knowledge this is the first observation of this kind in a Si/SiGe heterostructure with a 2DEG. In Fig. 1,b,c,d some of the curves are given in more detail. With our method of varying the electron density and a rather limited temperature range used it is difficult to determine with a high precision the ”separatrix”, i.e. the curve corresponding to the transition from the metallic to insulator regime. Nevertheless we consider the curve shown in Fig. 1c to be close to that transition point. There are two things about this curve that are worth mentioning. Firstly, the resistivity of the sample corresponding to the separatrix is about $`0.3h/e^2`$ which is about ten times lower than that reported for Si MOSFETs . Secondly, the peculiar non-monotonic temperature dependence corresponding to our separatrix is also quite unusual. At present we have no explanation for these observations.

Before the first observation of the MIT in Si MOSFET  it has been widely accepted that the electrons confined to two dimensions have an insulating ground state at zero magnetic field and zero temperature. On the other hand the application of a strong perpendicular magnetic field changes the situation, resulting in a transition form the insulating state to a QH liquid state. This magnetic-field-driven transition from the insulator to the QH liquid state has attracted much attention  and some interesting hypothesis have been proposed for its explanation . Later, it was suggested  that there possibly exists a link between the MIT transition observed at zero magnetic field and the finite magnetic field transitions studied previously. In the present work we have made a similar observation in our samples. Fig. 2 shows the magnetoresistance traces taken at different temperatures corresponding to some of the zero field sample states presented in Fig. 1. For any given magnetic field it is the sign of the resistance temperature dependence that is the indicator of whether the state of the sample at this particular field is insulating or metallic. Let us start with a sample in an insulating state at zero field. In Fig. 2a we see that for $`n_s=3.53\times 10^{11}`$cm-2 the zero field insulating state persists up to $`B_c^L2.8`$T where the first insulator-to-QH liquid transition occurs. This QH liquid state corresponding to filling factor $`\nu =4`$ is then maintained until the QH liquid-to- insulator transition turns the sample again into insulator at $`B_c^H4.8`$T . According to Fig. 1, increasing the electron density from $`n_s=3.53\times 10^{11}`$cm-2 brings the sample closer to the zero field MIT. It is interesting to note, that in magnetic field this results in a corresponding gradual decrease of $`B_c^L`$ until, finally, at the carrier density $`n_s=4.05\times 10^{11}`$cm-2 appropriate for the zero field MIT, $`B_c^L`$ becomes exactly zero, Fig. 2b. On further increasing the carrier density, the zero field metallic state is maintained at finite field and the first magnetic-field-driven transition to occur is thus a QH liquid-to-insulator transition at higher fields, Fig. 2c. As the authors of Ref. , we find the behavior described above to be an indication that both the finite- and zero-field transitions are likely to share a common physical origin.

Up to date there is still no theoretical description available for the behavior of a 2D system in the vicinity of the MIT, i.e. when $`\rho _{xx}h/e^2`$ , and all the discussions of the MIT are, therefore, necessarily qualitative in character. On the contrary there is now a number of comprehensive theoretical models for a situation when $`\rho _{xx}h/e^2`$. These models offer predictions concerning the corrections to the zero temperature Drude conductivity $`\sigma _D=e^2n\tau /m^{}`$ under various physical conditions. We now turn to the discussion of these quantum corrections in the limit $`\rho _{xx}h/e^2`$, i.e. when our samples show metallic behavior.

On the whole, there are two distinct types of corrections to conductivity. First, there is the weak localization correction calculated in the approximation of non-interacting electrons and describing the variation of conductivity due to the electron coherent scattering by the random potential :

$$\mathrm{\Delta }\sigma _{xx}^{wl}=\alpha p(e^2/\pi h)ln(k_BT\tau /\mathrm{}),$$
(1)
where it is assumed that the phase breaking time $`\tau _\varphi `$ varies as $`T^p`$ . The amplitude $`\alpha `$ is expected to be 1 for normal scattering (-0.5 for pure spin-orbit scattering and 0 for scattering on magnetic impurities). In Si samples only the Coulomb scattering is present. The above expression for weak localization correction to conductivity remains valid for all 2D systems with $`\rho _{xx}h/e^2`$ . The situation is more complicated when the interaction between electrons is taken into account. The first theory to consider the interaction corrections to conductivity was developed for the diffusive regime $`k_BT\tau /\mathrm{}<1`$ i.e. for a situation where the electron experiences multiple scattering by impurities during its life time . Of the two components of the interaction correction in this regime, the correction in the Cooper channel is negligible when $`T_cT`$ , where $`T_c`$ is the critical temperature for superconductor transition, and the remaining correction in the Diffusion channel has the form:

$$\mathrm{\Delta }\sigma _{xx}^{ee}=\frac{e^2}{2\pi ^2\mathrm{}}ln\left[\frac{k_BT\tau }{\mathrm{}}\right]\left[1+3\left[1\frac{ln(1+F_0^\sigma )}{F_0^\sigma }\right]\right]$$
(2)
where $`F_0^\sigma `$ is the interaction constant in the triplet channel. This constant is negative if the electron-electron interaction tends to align the electron spins and its absolute value increases with the interaction strength characterized by parameter $`r_s`$ . Note that the sign of this logarithmic in temperature correction depends on the absolute value of $`F_0^\sigma `$ . It is believed that in the limit of very strong interaction the constant $`F_0^\sigma `$ approaches -1 leading to the Stoner instability with all electron spins parallel in the absence of magnetic field.

Recently, a theory has been proposed , that along with the diffusive regime covers the intermediate and ballistic regimes $`k_BT\tau /\mathrm{}1`$. It is shown that the interaction corrections in all these regimes are a consequence of the same physical process, namely the coherent electron scattering on the modulated density of other electrons (Friedel oscillations) caused by an impurity with a short range potential. The quantum interaction correction in this theory can be expressed as:

$$\mathrm{\Delta }\sigma _{xx}^{ee}=\delta \sigma _C+15\delta \sigma _T$$
(3)
where

$$\delta \sigma _C=\frac{e^2}{\pi \mathrm{}}\frac{k_BT\tau }{\mathrm{}}\left[1\frac{3}{8}f(k_BT\tau /\mathrm{})\right]\frac{e^2}{2\pi ^2\mathrm{}}ln\left[\frac{1}{k_BT\tau /\mathrm{}}\right];$$

is the charge channel correction and

$$\delta \sigma _T=\frac{F_0^\sigma }{\left[1+F_0^\sigma \right]}\frac{e^2}{\pi \mathrm{}}\frac{k_BT\tau }{\mathrm{}}\left[1\frac{3}{8}t(k_BT\tau /\mathrm{};F_0^\sigma )\right]$$

$$\left[1\frac{1}{F_0^\sigma }ln(1+F_0^\sigma )\right]\frac{e^2}{2\pi ^2\mathrm{}}ln\left[\frac{1}{k_BT\tau /\mathrm{}}\right];$$

is the correction in the triplet channel. We also take into account the two fold valley degeneracy at (100) Si surface that changes the triplet term prefactor from 3 to 15. The explicit expressions for the functions $`f(k_BT\tau /\mathrm{})`$ and $`t(k_BT\tau /\mathrm{};F_0^\sigma )`$ are given in Ref. . In the diffusive limit Eqs.(3) reproduce the well known logarithmic temperature correction Eq.(2). In the intermediate and ballistic regimes Eqs.(3) result in a linear with temperature correction to conductivity with the slope and the sign determined by the absolute value of the interaction constant $`F_0^\sigma `$. For small $`F_0^\sigma `$ the triplet channel contribution is small and the temperature dependence, governed mostly by the correction in the charge channel, is insulating, whereas for higher $`F_0^\sigma `$ values the negative triplet term increases and the behavior becomes metallic. Also, since the triplet term prefactor increases with the valley degeneracy, systems with a higher valley degeneracy show metallic behavior for lower values of $`F_0^\sigma `$ and hence for lower interaction strength than systems with lower or no valley degeneracy.

Fig. 3 shows the typical magnetoresistance traces obtained at different temperatures in our samples after the electron concentration had been saturated with a maximum LED illumination. The corresponding electron mobility is $`\mu =61800`$cm2/Vs and the electron density $`n_s=6.23\times 10^{11}`$cm-2. In Fig. 4 we plot the zero field conductivity values form Fig. 3a versus temperature. In line with the earlier report on a high mobility n-Si/SiGe structute, one can see in Fig. 4 that the sample behavior is markedly metallic. The temperature dependence appears to be linear for $`T1.25`$K and saturates at lower temperatures. By extrapolating the linear part of the curve to $`T=0`$ we find the Drude conductivity and the corresponding momentum relaxation time $`\tau =6.8\times 10^{12}`$s. Hence, we get $`k_BT\tau /\mathrm{}=0.89\times T`$ which means that in the experimental temperature range $`0.42.7`$K the sample operates in the intermediate and ballistic regimes. The solid line in Fig. 4 is a theoretical fit given by the sum of the weak localisation term (Eq.(1)) and the interaction correction (Eqs.(3)) with the only fitting parameter being the interaction constant $`F_0^\sigma =0.155`$. The weal localisation term prefactor $`\alpha p=1.5`$ used in the fit is obtained from the analysis of the magnetoresistance data in low fields.

It is impossible to compare the calculated value of $`F_0^\sigma `$ with the experimental parameter $`r_s`$ which is $`6.7`$ for our sample since the functional relation between $`F_0^\sigma `$ and $`r_s`$ is unknown for $`r_s>1`$. It is, however, interesting to note that systematically higher values of $`F_0^\sigma `$ have been reported for various other types of 2D systems and, in particular, for Si MOSFETs with $`r_s`$ values close to ours. For example, $`F_0^\sigma 0.27`$ has been found in an n-channel Si MOSFET  and in a p-type GaAs/AlGaAs heterostructure  with $`r_s7`$. The difference in the $`F_0^\sigma `$ values obtained in our Si/SiGe heterostructure and in an n-channel Si MOSFET  is, in our opinion, the consequence of a difference in the structure of the scattering potential in these two silicon systems. The linear-in-$`T`$ term of Ref.  arises only due to the scattering on a short range impurity potential (in a purely smooth disorder this term is exponentially suppressed). In a Si MOSFET with a predominantly short range scattering potential the results of Ref.  are, thus, quite appropriate. In our samples the single particle scattering time $`\tau _q`$ obtained from the magnitude of the Shubnikov-de Haas oscillations is about 6 times smaller than the momentum relaxation time $`\tau `$ indicating that apparently both types of disorder are present. Therefore the prefactor in front of the linear term would be reduced in our case as compared to the corresponding expression in Ref.  and, hence, a lower value of $`F_0^\sigma `$ would be obtained from fitting the theory  to experiment which, under these conditions, is unlikely to reflect the actual interaction intensity in the sample. This effect is expected to be even more pronounced in such high mobility structures as n- and p- GaAs/AlGaAs and, in our opinion, a certain caution should be used when comparing experiment with theory  in these cases.

On the whole the theory  produces a good fit for the experimental curve at temperatures $`T1.25`$K. On the lower temperature side the dependence $`\sigma _{xx}(T)`$ tends to saturate, in contrast to the theoretical prediction. Similar departures from theory have been observed in all our samples and also by other authors in a Si 2D system different from ours . A non-zero valley splitting and a strong sample specific inter-valley scattering have been suggested as a possible explanation of the weakening of the $`\sigma _{xx}(T)`$ dependence. At present there is no theory available that takes inter-valley scattering into account.

Let us now analyze the transverse magnetoresistance in our samples, Fig. 3a. In the absence of any quantum corrections to the conductivity one has: $`\sigma _{xx}^0(B)=\sigma _D/\left[1+(\mu B)^2\right]`$, $`\sigma _{xy}^0(B)=\mu B\sigma _{xx}^0(B)`$, $`\rho _{xy}^0(B)=\mu B/\sigma _D`$ and $`\rho _{xx}^0(B)=1/\sigma _D=const`$, i.e. the resistance of the sample is independent of magnetic field. If, on the other hand, there are some small corrections to the conductivity tensor components $`\sigma _{xx}(B)=\sigma _{xx}^0(B)+\mathrm{\Delta }\sigma _{xx}`$ and $`\sigma _{xy}(B)=\sigma _{xy}^0(B)+\mathrm{\Delta }\sigma _{xy}`$, then, converting the conductivity tensor into the resistivity tensor, one gets:

$$\rho _{xx}(B)=\rho _D+\rho _D^2\left[\left((\mu B)^21\right)\mathrm{\Delta }\sigma _{xx}+2\mu B\mathrm{\Delta }\sigma _{xy}\right]$$
(4)
In general there will be contributions to the magnetoresistance both from the weak localization and from the interaction. Weak localization gives rise to magnetic-field-dependent corrections to $`\sigma _{xx}^0(B)`$ and $`\sigma _{xy}^0(B)`$ but these are suppressed and can be neglected in magnetic fields higher than $`B_{tr}\mathrm{}/2el_e^2`$ which is $`1mT`$ for our sample. The situation with the interaction corrections to the magnetoresistance is more complicated. In the diffusive limit $`k_BT\tau /\mathrm{}1`$ it has been established that the longitudinal conductivity correction $`\mathrm{\Delta }\sigma _{xx}^{ee}`$, given by Eq.(2), remains unchanged at low and classically strong magnetic fields, while, at the same time, the Hall conductivity is not affected by interaction: $`\mathrm{\Delta }\sigma _{xy}^{ee}0`$ . As can be seen from Eq.(4) for $`\omega _c\tau >1`$ this will lead to a parabolic negative magnetoresistance:

$$\rho _{xx}(B)=\rho _D+\rho _D^2(\mu B)^2\mathrm{\Delta }\sigma _{xx}^{ee}(T)$$
(5)
with $`\mathrm{\Delta }\sigma _{xx}^{ee}(T)`$ given by Eq.(2). Finally, at still higher fields $`B_ZkT/g\mu _B`$ a positive magnetoresistance develops caused by the Zeeman effect on the interaction correction $`\mathrm{\Delta }\sigma _{xx}^{ee}`$, .

In contrast to the diffusive regime, the behavior of magnetoresistance is less well understood in the intermediate and ballistic regimes $`k_BT\tau /\mathrm{}1`$. The theory of interaction corrections at arbitrary relation between temperature and elastic mean-free time  makes certain predictions concerning the corrections to the Hall coefficient  and magnetoresistance in a parallel field  but leaves magnetoresistance in a perpendicular magnetic field out of consideration. Thus it remains unclear whether the corrections to conductivity obtained in this theory for zero field would persist in finite magnetic field and give rise to a magnetoresistance in a manner similar to the diffusive regime, Eq.(5), or whether they will be modified or suppressed by magnetic field. As has been mentioned above, the results of the theory  are obtained in the approximation where the impurities are treated as point-like scatterers. This is justified at zero magnetic field where a long-range potential reduces the backscattering and suppresses the interaction corrections in the ballistic regime. The situation, however, changes in a strong magnetic field which increases the probability of an electron to return back and, thus, restores the interaction correction even in the case of a purely smooth disorder. Recently, a theory  has been proposed that considers the interaction in the ballistic and intermediate regime in classically strong fields $`\omega _c\tau >1`$. The theory  takes into account only the long-range component of the scattering potential and shows that in this case the interaction will also lead to a parabolic negative magnetoresistance described by Eq.(5) with a distinct expression for $`\mathrm{\Delta }\sigma _{xx}^{ee}(T)`$. However, as will be shown below, while justified in the case of a predominantly long range scattering potential, the exclusion of the short range component of disorder from consideration may be inappropriate in a more general case of a mixed mechanism of scattering. Finally, it might be well to point out that unlike the diffusive regime the conductivity correction $`\mathrm{\Delta }\sigma _{xx}^{ee}(T)`$ extracted from magnetoresistance in the ballistic regime  has nothing to do with the correction at zero field (the latter being exponentially suppressed in the case of a smooth disorder).

Up to date there have been several experimental studies of the interaction-related transverse magnetoresistance in the intermediate magnetic field range where the weak localization is already suppressed while the Zeeman effect on the interaction correction has not yet set in. First to mention are the experiments  where the observed temperature dependent negative magnetoresistance was compared to the theoretical prediction for the diffusive regime, Eq.(5). In one of these works  the samples used were n-Si/SiGe heterostructures very similar to ours. However with the exception of Ref.  these experiments have been performed on high mobility structures in the ballistic regime where the validity of Eq.(5) for the description of the magnetoresistance is questionable. Recently an experimental work has been reported  whose aim was to study the magnetoresistance caused by electron-electron interactions in the intermediate regime in a n-GaAs/AlGaAs heterostructure with long range fluctuation potential and to compare the results with the theoretical prediction . A good agreement between the theory and experiment has been reported.

Similar to  in our case we have a relatively high mobility 2D electron gas in the intermediate and ballistic regimes. The important role of the long range scattering potential is evidenced by the fact that the transport time is significantly longer than the quantum single-particle time. The sharp change of the magnetoresistance observed in Fig. 3a in low magnetic field is caused by the WL effect suppressed at higher fields. Then follows a relatively flat region, which, according to , is supposed to correspond to the suppression of backscattering by the long-range potential in low fields. This effect is more pronounced for higher temperatures. At still higher field the backscattering is restored by the increasing magnetic field and a parabolic magnetoresistance sets in. In our sample the condition $`\omega _c\tau =1`$ is satisfied at $`B=0.16`$T and so the parabolic magnetoresistance in Fig. 3a is observed in classically strong fields. On the other hand, except for $`T=0.4`$K, the magnetoresistance traces in Fig. 3a are measured at fields well below $`B_Z`$\[T\]$`=0.74\times T`$\[K\] at which the Zeeman effect on the interaction correction should become important. That means that the parabolic magnetoresistance in Fig. 3a is observed in the magnetic field region where the results of Ref.  are expected to be applicable. In Fig. 5a we plot the resistivity from Fig. 3a as a function of $`B^2`$. Shown in Fig. 5b are $`\mathrm{\Delta }\sigma _{xx}^{ee}(T)`$ obtained using Eq.(5) from the slope of the straight part of the magnetoresistance traces in Fig. 5a. The theoretical curves plotted in Fig. 5b are

$$\mathrm{\Delta }\sigma _{xx}^{ee}=\frac{2}{\pi }G_F(k_BT\tau /\mathrm{})$$
(6)
and

$$\mathrm{\Delta }\sigma _{xx}^{ee}=\frac{2}{\pi }\left[G_F(k_BT\tau /\mathrm{})G_H(k_BT\tau /\mathrm{};F_0^\sigma )\right]$$

where $`G_F(k_BT\tau /\mathrm{})`$ is the exchange contribution to the interaction correction to conductivity in classically strong magnetic fields and in the case of a long range scattering potential. $`G_H(k_BT\tau /\mathrm{};F_0^\sigma )`$ is the triplet Hartree term. The explicit expressions for these functions are given in Ref. . For plotting in Fig. 5b we take these functions in the appropriate for our experimental situation limit $`\kappa k_F`$, where $`k_F`$ is the Fermi wave vector and $`\kappa =4\pi e^2\nu `$ is the inverse screening length. The Hartree term is plotted for $`F_0^\sigma =0.155`$, the interaction constant obtained from the fit of the zero-field resistance temperature dependence to theory . As can be seen there is a rather big divergence between the experimental points and the theoretical curves in Fig. 5b. The exchange correction independent of the interaction constant lies closer to the experiment. Taking into account the Hartree contribution with a sign opposite to the exchange term increases the discrepancy, which becomes even greater if it is multiplied by a factor $`5`$ corresponding to the valley degeneracy, (compare Eqs.(3)). In our opinion the explanation why in contrast to Ref.  we don’t find a good agreement with Ref.  in our experiment may lie with the fact that theory  was developed with only the long-range component of the disorder potential taken into account. The structures studied in Ref.  were of a higher mobility than ours and with a more important role of the long-range scattering potential compared to that of the point-like scatterers. Probably in that case the approximation made in Ref.  is adequate, whereas in our samples the short-range scattering can not be completely neglected and an extension of theory  to the case of a mixed disorder is needed to explain our results.

## IV Conclusions

To conclude, we have observed for the first time a MIT in a Si/SiGe heterostructure with a 2DEG. We have also observed the evolution of the magnetic-field-driven transition from an insulator to a QH liquid state as our system is step by step driven through the MIT by gradually increasing the electron density. Similar to Ref.  we find that the magnetic field corresponding to the insulator-to-QH liquid state transition becomes zero at the critical electron density appropriate for the MIT in zero field. On the metallic side of the MIT we find that the theory  of the interaction corrections to conductivity for an arbitrary relation between temperature and the momentum relaxation time provides a good description of the linear part of the conductivity temperature dependence at zero field but with the interaction constant $`F_0^\sigma `$ about two times smaller than in Si MOSFETs with a similar value of $`r_s`$. The recently developed theory  of the interaction contribution to magnetoresistance at classically strong magnetic fields and for arbitrary $`k_BT\tau /\mathrm{}`$ does not seem to be appropriate for the description of the parabolic negative magnetoresistance observed in our system. The theory  was built in the assumption that there is only a long-range scattering potential. The importance of a short-range disorder in our samples may be the reason for the observed discrepancy with the theory. Both these theories  neglect either the long or the short range component of disorder and for that reason do not seem to be quite adequate for the description of the case of a mixed disorder.

We are grateful to I.V. Gornyi and A.D. Mirlin for helpful discussion. This work was supported by PICS-RFBR (N 1577.), RFBR (N 02-02-16516), NATO Linkage (N CLG.978991), INTAS (N01-0014), programs ”Physics and Technology of Nanostructures” of the Russian ministry of Industry and Science and ”Low dimensional quantum structures” of RAS.


format: text

text:
# Metastable states, anomalous distributions and correlations in the HMF model

## 1 Introduction

There is by now a huge number of papers devoted to understanding the physical foundation and the possible applications of the nonextensive statistical mechanics formalism proposed by Constantino Tsallis in 1988 . In particular, although the topic is still a matter of debate , the last few years have seen an important turning point for what concerns dissipative systems both low-dimensional (for example unimodal maps)  and high-dimensional ones (for example turbulence, plasma and nuclear physics, econophysics) . Progress has been made also for conservative systems both for maps  and for high-dimensional hamiltonian systems with long-range interactions , but open problems still remain and need a more detailed investigation. With this scenario in mind, in this paper we present and discuss new numerical simulations relative to the Hamiltonian Mean Field (HMF) model, a system of $`N`$ fully coupled classical spins which has been intensively studied in the last years. It is nowadays well known that, for some values of the initial energy, below the critical point, the system does not immediately relax to the Boltzmann-Gibbs equilibrium, but remains trapped in anomalous metastable states and exhibits non-Gaussian velocity distributions that can be reproduced by the probability density functions (pdfs) predicted by the Tsallis nonextensive thermodynamics . Such states have been named Quasi-Stationary States (QSS) because they become stable in the thermodynamic limit if the $`N\mathrm{}`$ limit is taken before the $`t\mathrm{}`$ limit. In fact in such a case the force between spins tends to zero and the largest Lyapunov exponent vanish .

In this work we study two classes of out-of-equilibrium initial conditions that lead to QSS, and we discuss how the different initial conditions may affect the dynamics of the relaxation towards equilibrium. In particular, we concentrate our attention on the existence of correlations in the orientation angles and velocity of the spins, on the presence of anomalous velocity distributions, and on the clusters formations. The main result of this paper is that Tsallis nonextensive thermodynamics formalism seems to provide a general framework, to interpret the statistical properties of the dynamically created QSS, for Hamiltonian systems with long-range interactions , only when long-range correlations and fractal structures in phase space are present. Recently some critical comments  have been raised on a previous study of two of us . In this paper we provide new calculations and detailed discussions which reply in part to this criticism and clarify also some points advanced in another work . A more detailed reply to the points raised in ref.  will be presented elsewhere .

The paper is organized as follows. In section II the details of the model are summarized and the initial conditions used are discussed. Numerical results for QSS are discussed in section III, while velocity pdfs and dynamical correlations are studied in a quantitative way in section IV. Conclusions are drawn in section V.

## 2 Model and initial conditions

The HMF model, consists of N planar classical inertial spins interacting through an infinite-range potential . The Hamiltonian can be written as:

$$H=K+V=\underset{i=1}{\overset{N}{}}\frac{p_{i}^{}{}_{}{}^{2}}{2}+\frac{1}{2N}\underset{i,j=1}{\overset{N}{}}[1cos(\theta _i\theta _j)].$$
(1)
In the latter $`\theta _i`$ is the orientation angle of the $`ith`$ spin and $`p_i`$ the conjugate variable representing the rotational velocity, considering the mass of the rotators equal to one. The summation in V is extended to all couples of spins and is not restricted to first neighbors. In order to make H formally extensive, i.e. $`VN`$ when $`N\mathrm{}`$, the coupling constant in the potential is divided by N. However, this often adopted prescription cures only a part of the problem since the energy remains non-additive, that is the system cannot be simply divided in two independent sub-systems. The magnetization M, i.e. the modulus of $`𝐌=\frac{1}{N}_{i=1}^N𝐦_i`$, whith $`𝐦_i=[cos(\theta _i),sin(\theta _i)]`$, is the order parameter of the model. The latter shows a a second-order phase transition from a low-energy clustered phase with magnetization $`M1`$, to a high-energy one, where the spins are homogeneously oriented on the unit circle and $`M0`$. The model can be solved exactly in the canonical ensemble and the dependence of the energy density $`U=E/N`$ on the temperature $`T`$, is given by

$$U=\frac{T}{2}+\frac{1}{2}\left(1M^2\right).$$
(2)
The critical point is at energy density $`U_c=\frac{3}{4}`$ corresponding to a critical temperature $`T_c=\frac{1}{2}`$ . The dynamics of HMF model can be investigated by starting the system with out-of-equilibrium initial conditions and integrating numerically the equations of motion . In a special region of energy values ($`\frac{1}{2}<U<U_c`$) the results of the simulations show, for a transient regime which depends on the system size, a disagreement with the canonical ensemble. In this region the dynamics is characterized by Lévy walks and anomalous diffusion, while in correspondence the system shows a negative specific heat . Ensemble inequivalence and negative specific heat have also been found in self-gravitating systems , nuclei and atomic clusters , though in the present model such anomalies emerge as dynamical metastable features  and vanish, if $`N`$ is finite, for $`t\mathrm{}`$. In this paper we focus on a particular energy value belonging to the anomalous region, namely $`U=0.69`$, and we study the time evolution of temperature, $`\mu `$-space structures, dynamical clusters in angle, velocity correlations and velocity distributions considering two out-of-equilibrium initial conditions in particular:

* the first one, used already in previous papers, consists of all orientation angles $`\theta _i=0`$ and velocity $`p_i`$ uniformly distributed accordingly to the total energy available. In this case the potential energy is initially zero and we have the maximum of kinetic energy, while the total magnetization M is equal to 1;
* in the second case considered, both angles and velocities are uniformly distributed. This implies that the total magnetization M is initially zero, so that the potential energy is $`\frac{1}{2}`$, see eq.(2). Therefore in this case only a smaller fraction of total energy is left for the kinetic part.

For both initial conditions the total momentum is initially put to zero. In the following, for brevity, we will indicate these two initial conditions respectively by $`M1`$ and $`M0`$, thus referring to the initial value of the magnetization.

## 3 Metastable equilibrium and Quasi-Stationary States

We discuss in this section the numerical results obtained starting from M1 and M0 initial conditions. Different sizes of the system are studied. The accuracy of the calculations is such that one gets an energy conservation $`\mathrm{\Delta }E/E=10^5`$. As in previous papers the fourth order Yoshida simplectic algorithm  has been used for the integration. At variance with the previous publications the averages where taken only over the ensemble and not also along the time evolution. This point will be discussed in detail later.

In fig.1 we plot, for system sizes N=500 and N=1000, the time evolution of twice the kinetic energy per particle, a quantity which is the temperature in stationarity situations. Different regimes exist according to the initial conditions considered. More specifically, when starting with $`M1`$ initial condition (IC) four regimes are observable:

1. a fast decay from $`M=1`$ to $`M0`$ in a time which is about 80-100. No sensible dependence on the size is observed;
2. a plateau region with $`M(N)<M_{eq}`$ and $`T(N)<T_{eq}`$ whose duration depends also on the size of the system as N ;
3. a power-law relaxation towards the true equilibrium value $`T_{eq}=0.476`$ ;
4. an equilibrium regime with $`M_{eq}=0.307`$.

In refs. the kinetic energy per particle was actually computed considering not only an ensemble average, but also a time average. This was done to improve statistics, though in principle time averages are justified only when correlations do not exist. In fig.2 the old procedure and the new one are compared. One observes that the two methods give the same plateau if a large number of events are considered in the ensemble average. The only difference is in the decay from the initial condition to the QSS regime. One can easily check that also the scaling laws found in ref.  continue to be valid. In particular the lifetime of the plateau diverges with $`N`$ and the value of the plateau temperature $`T_{QSS}(N)`$ converges towards the asymptotic temperature $`T_{\mathrm{}}=0.38`$ as $`N\mathrm{}`$ following the same laws of refs.. At variance with the criticism advanced in , it is evident that one can really talk of an approximate plateau for the average kinetic energy per particle, see fig.2. The log timescale used in fig.1 and in previous papers  is useful to have an overall view of the temporal behavior and does not change the macroscopic quasistationary features. This will be clearer when analizing also the M0 IC, for which the system is initially put in the state at $`M0`$ and $`T0.38`$ and there it remains for quite a long time which also in this case depends on the size studied. In fact when starting with $`M0`$ IC only three regimes are observable:

1. a plateau region with $`M0`$ and $`T0.38`$ whose length depends on the size as $`N`$. This has been checked numerically but it is not shown here for lack of space;
2. a power-law relaxation towards the true equilibrium value, but faster than the previous $`M1`$ case;
3. an equilibrium regime with $`M_{eq}=0.307`$.

In general, for both IC, the values obtained for the magnetization and the temperature in the plateau regime for $`N\mathrm{}`$ are consistent with the caloric curve (2) extended below the critical point. Therefore this average procedure is a-posteriori justified and consistent with a metastable macroscopic thermodynamic description, contrary to what claimed in ref..

Although both initial conditions drive the system towards the same out-of-equilibrium macroscopic values for the magnetization and the temperature in the $`N\mathrm{}`$ limit, for the QSS regime, the corresponding microscopic dynamics is very much different in the two cases. This fact can be seen when one plots the time evolution of the $`\mu `$-space, as done in fig.3 for a typical event and $`N=10000`$. In red we represent the rotators having at time $`t=0`$ the highest velocities. In this way it is possible to visualize the kind of mixing occurring for the initial conditions considered. It is immediately evident that there are correlations and structures in the plateau regime only for M1 IC, while for M0 IC the distributions remain homogeneous with a very slow mixing of the initially fast rotators. This different behavior can be observed also in figure 4, where we plot the pdf in angle at various times. A typical realization for N=10000 is plotted. One observes clusters formation in the plateau regime only for M1 while these are completely absent for M0 IC. However, at equilibrium, the pdfs are equivalent and the plots, for both M1 and M0, show only one rotating big cluster corresponding to $`M=0.307`$ which is the equilibrium magnetization (according to the canonical prediction eq.(2) with $`U=0.69`$ and $`T=0.476`$).

These numerical results reinforce our previous statement that although macroscopically similar the two metastable plateaux have a very different microscopic dynamics. Further evidence is provided if one studies the entropy and free energy as a function of time.

This is done in figs.5 and 6 respectively. To obtain the entropy we adopted a lattice over the $`\mu `$-space made up of $`n=10000`$ cells and we used the following formula:

$$S(t)=<\underset{k=1}{\overset{n}{}}f_k(t)ln(f_k(t))>,$$
(3)
where $`f_k`$ is the relative frequency of points for the $`kth`$ cell. The brackets $`<\mathrm{}>`$ denote ensemble averages. The curves plotted in the figures refer to an ensemble of $`1000`$ events for $`N=1000`$. This particular size was used in this case and in other calculations when a good average over the ensemble was essential in order to have a robust value, computer time being too long for greater sizes and such big ensemble averages. On the other hand no qualitative changes are anyhow observed for greater sizes. As intuitively expected from the temporal evolution of the $`\mu `$-space, the entropy grows with time towards equilibrium only for M1 IC due to the high degree of order of this initial configuration, while it remains fixed at the maximum value for the other, much more uniform, initial condition, $`M0`$. Correspondingly, we compute the free energy by extrapolating to a dynamical situation the equilibrium thermodynamic formula

$$F(t)=UT(t)S(t).$$
(4)
As shown in fig.6 the free energy remains constant for $`M0`$ in the QSS regime. This implies that for this IC the system is consistently in a metastable regime, which is not however the true one that is reached only at greater times, when the free energy of the system has a real global minimum. On the other hand, for the $`M1`$ IC, the figure shows a different temporal evolution before equilibration. After a short rapid transient decay, an almost flat region, at a different free energy value, always emerges in the QSS regime, implying also in this case the presence of a metastable situation. The use of formula (4) is thus a-posteriori justified since the temperature is almost constant in the plateau regime. It is also very useful in clarifying the scenario that emerges which is consistent with a quasi-equilibrium situation. This results reinforce the validity of our procedure to extract average thermodynamic quantities from microscopic dynamics and can be considered a further answer to the criticism raised in ref.. On the other hand, the fact that there is no plateau for the temporal evolution of the entropy for $`M1`$ IC reflects the slow microscopic dynamics observed in the evolution of the $`\mu `$-space structures. This is perfectly in line with the quantitative analysis performed by means of the correlation dimension in ref.. Thus one can summarize this situation by saying that only at a macrosopic level the quasi-stationarity is valid, while microscopically the system never stops its slow temporal evolution. The evolution really stops only for $`N=\mathrm{}`$, since in this case the magnetization in the QSS is exactly zero and so is the force which is proportional to $`M`$. In this sense the inversion of the time limit with the size limit to infinity stabilize the QSS regime. However, when considering finite $`N`$ systems, the fluctuations due to finite-size noise, of the order $`1/\sqrt{N}`$, induce a force which is not exactly zero, and when they grow with time, the system is slowly forced towards equilibration. As shown in ref. the accuracy of our calculations is much smaller than this finite-size noise.

For the M0 IC no entropy increase is observed. But also in this case some microscopic dynamics exists. This is proved by the mixing of fast particles in fig.3. However no holes or structures are present in this case and this is the reason why the entropy temporal evolution is almost flat.

It is important to stress the fact that, for very long times, both the two initial conditions lead the system towards the equilibrium expected values.

We would like to notice at this point that the slope of the entropy growth in fig.5 does not give the Kolmogorov-Sinai entropy as demonstrated numerically for maps . In fact in our case we are considering the reduced $`\mu `$-space and not the complete phase-space which has to be used. Unfortunately in our case the phase space has too many dimensions and this fact prevents the direct investigation of the generalized q-entropy proposed by Tsallis along the same lines of ref. for maps.

## 4  Anomalous velocity pdfs, velocity correlations and aging

In this section we study in a quantitative way the eventual presence of velocity correlations and memory effects in the HMF dynamics. In Fig.7 we report the velocity pdfs for the plateau regime (at $`t=1000`$) and at equilibrium $`(t=1000000)`$ and for the case N=10000. For both IC, the initial velocity pdfs quickly acquire and maintain during the entire duration of the metastable state a non-Gaussian shape. But while for M1 IC they can be fitted with a truncated generalized Tsallis pdf , this is not possible for M0 IC. In fact in this case the pdf shows a flat almost rectangular shape, which clearly reveals the fact that structures and correlations, present in the other case, are here missing. Finally, for very long times the system reaches always the equilibrium distribution and the agreement with the gaussian is almost perfect.

A quantitative way to estimate the velocity correlations is the calculation of the two-time autocorrelation function, defined as

$$C(t,0)=\frac{<𝐏(t)𝐏(0)><𝐏(t)><𝐏(0)>}{\sigma _p(t)\sigma _p(0)},$$
(5)
where $`𝐏=(p_1,p_2,\mathrm{}p_N)`$ is the N-component velocity vector and the brackets $`<\mathrm{}>`$ indicate the average over the ensemble, while $`\sigma _p(t)`$ and $`\sigma _p(0)`$ are the standard deviations at time t and at time zero. In fig.8 we plot the velocity autocorrelation functions for M1 and M0 IC in the QSS regime, obtained for the case N=1000 and an average over 500 events. It is immediately evident that for M0 the autocorrelation decays faster than for M1 IC. This result clarifies the controversial numerical results obtained in ref. by using a slightly different, but equivalent, definition for the autocorrelation function. In fact we do find an almost exponential behavior only for the M0 IC, while in ref.  the author claimed to find an exponential in both cases. In general we can reproduce the autocorrelation decay by means of the q-exponential function

$$e_q(x)=A\left[1(1q)\frac{x}{\tau }\right]^{\frac{1}{(1q)}}$$
(6)
proposed by Tsallis in his thermodynamic nonextensive q-formalism , where $`A`$ is the saturation value and $`\tau `$ the characteristic time. One gets the usual exponential for $`q=1`$. In this way we can quantitatively discriminate between the two different initial conditions. In fact we get for M1 the value $`q=1.55`$, revealing a power law scaling, while for M0 we get a value of q much smaller, i.e. $`q=1.12`$. This means that for M0 the decay tends to be exponential (an exact exponential is obtained when $`q=1`$). This almost exponential decay is also evident in the inset of the figure, where we show the same data in linear-log scale. An exponential fit is also reported for comparison. Such numerical results may be considered a quantitative confirmation of the different microscopic dynamics of the system according to the particular initial condition, in spite of the macroscopic presence of QSS for both cases.

It is interesting to investigate also if the relaxation process observed in fig.8 depends on the history of the system, i.e the possible existence of strong long-term memory effects, the so-called aging phenomenon. This is typical of glassy systems but it has also been found in a great variety of complex systems . Recently it has been observed also in the HMF model by Montemurro, Tamarit and Anteneodo . In the following we perform a similar analysis to study this effect for both M1 and M0 IC, but considering only the N-component velocity vector defined previously. In their case also angles were taken into account in the state vector entering the correlation function and only M1 IC were considered. The autocorrelation function to be considered is that one defined before but calculated with respect to a set of four different waiting times $`t_w`$ (that plays the role of age of the system)

$$C(t+t_w,t_w)=\frac{<𝐏(t+t_w)𝐏(t_w)><𝐏(t+t_w)><𝐏(t_w)>}{\sigma _p(t+t_w)\sigma _p(t_w)}.$$
(7)
In Fig.9 (a) we show the autocorrelation decay for M1 IC, while in Fig.10 (a)that one for M0 IC. If aging were not present, all the curves should overlap into one, because the function $`C(t+t_w,t_w)`$ would depend only on the difference between the two times $`t+t_w`$ and $`t_w`$, i.e. $`C(t+t_w,t_w)=C(t)`$, and the decay would be exponential. Actually this is not what happens as shown in the figures.

For M1 IC aging properties surely exist, i.e. it can be observed a characteristic two times dependence on $`t`$ and $`t_w`$ and the data can be scaled onto one curve, panel(b), considering as discussed in ref.:

$$𝐂(t+t_w,t_w)=f(\frac{t}{t_{w}^{}{}_{}{}^{\beta }}),$$
(8)
with $`\beta =\frac{1}{4}`$. In such a way we obtain, for the tail a power law decay, i.e. $`f(t/t_w^\beta )(t/t_w^\beta )^\lambda `$ with scaling parameter $`\lambda =1.54`$ that is here characteristic of a slow relaxation dynamics and long-term memory effects. A q-exponential curve can reproduce very well this behavior with $`q=\frac{\lambda +1}{\lambda }=1.65`$, $`A=0.7`$ and $`\tau =60`$. Notice that the scaling indices are different from those found in ref.. A sistematic study of these scaling exponents would be desirable, but is beyond the scope of the present work and we postpone it to future investigations. As discussed in ref. the origin of aging in the HMF model is not completely clear. The most likely scenario is that of a weak ergodicity-breaking. As originally proposed by Bouchaud for glassy systems in ref., the latter occurs when the phase-space is a-priori not broken into mutually inaccessible regions, as in the true ergodicity-breaking case, but the system can remain trapped for very long times in some regions of the complex potential landscape. In HMF the weak ergodicity-breaking can be related to the complex dynamics generated by the vanishing of the largest Lyapunov exponent and by a sort of dynamical frustration due to the existence, in the QSS regime, of different small clusters (see fig.4). These clusters compete in trapping more and more particles until only one of them remain when the system reach the standard equilibrium.

On the other hand, a different behavior is obtained for M0 IC, as shown in fig. 10. The same scaling law used for the M1 IC case does not seem to apply, see panel (b), and the decay is now more rapid. In this case it is not so clear as before if aging does really exist, but again a quantitative difference emerges for the two dynamics. A more detailed study in this respect is left for the future.

Before closing, it is important to notice that the aging observed is not restricted to the QSS region, but extends also to the following relaxation regime towards equilibration.

## 5 Conclusions

In this paper we have studied the microscopic dynamics towards relaxation in the HMF model, a Hamiltonian system of fully coupled inertial spins. Two different out-of-equilibrium initial conditions which lead to metastable long-living quasi-stationary states (QSS) have been studied. Our numerical results indicate that the dynamics of the relaxation process is extremely rich and more complex than what previously thought. We have found that for a certain class of initial conditions, i.e. M1 IC, one can have metastable QSS with non Gaussian velocity pdfs which can be reproduced with truncated q-exponential curves and for which a finite size scaling seems to apply . We have shown that in such a regime there are also strong velocity correlations in time which decay as a q-exponential. Moreover, confirming similar studies , in correspondence of this metastable states and in the following relaxation regime, aging in the velocity correlations has been found. Such a strong long-term memory effect can be attributed to a glassy-like weak ergodicity-breaking due to a sort of dynamical frustration which disappears at equilibrium. These new numerical facts, added to the anomalous diffusion , the structures in angle and the weak mixing produced by a vanishing Lyapunov exponent , previously observed for this out-of-equilibrium regime, although do not constitute a rigorous proof, certainly provide with a strong indication that the nonextensive formalism proposed by Tsallis could probably be applied also for Hamiltonian many-body systems, expecially when the infinite size limit is performed before infinite time limit. A rigorous link between the entropic index $`q`$ and the dynamical properties of nonextensive Hamiltonian many-body systems is still missing and must be found in order to confirm definetely Tsallis formalism. Such metastable QSS provide an ideal benchmark for testing the nonextensive theory, a project that we plan to continue in future works. We have also presented numerical evidence that when adopting a different class of out-of-equilibrium initial conditions, i.e. M0 IC, one can have also metastable QSS, but with no structures in the $`\mu `$-space and a fast (almost exponential) decay of velocity correlations. In this case nonextensive statistics certainly does not apply. Open questions, concerning aging for example, for this M0 IC case still remain and will be studied in the future.

## 6 Acknowlegments

We thank Celia Anteneodo and Constantino Tsallis for stimulating discussions. One of us (A.R.) would like to thank the Los Alamos National Laboratory for the financial support.


text:
# Experimental Vacuum Squeezing in Rubidium Vapor via Self-Rotation

## Introduction

Squeezed states of a radiation field are nonclassical states whose field-quadrature noise is phase dependent with a minimum below the shot noise level. First investigated theoretically Walls83  and generated experimentally Slusher85 ; Wu86  in the 1980’s, these states awoke initial interest as a way of overcoming the shot-noise precision restrictions in optical measurements Grangier87 ; Xiao87  and enhancing the capacity of communication channels YuenShapiro80 . Further interest to squeezed light was triggered by recent advances in the field of quantum information and communication. In this context they are viewed primarily as a means of generating quadrature entanglement Ou92 , a valuable resource that can be applied, for example, for quantum teleportation Furusawa98 , computation LLoyd99 , cryptography Ralph00  and manipulation of atomic quantum states Kuzmich00 . Protocols of error correction Braunstein982 ; Lloyd98  and entanglement purification Duan00  have been elaborated, establishing continuous variable entanglement as one of the most promising tools for quantum information technology.

A number of techniques for producing squeezed states of light have been developed, using for instance four wave mixing in atomic vapor Maeda87 , Kerr effect in cold atoms Lambrecht96  or optical fibers Shelby86 ; Silberhorn01 , optical parametric amplifiers Wu86 . All these methods require rather complicated setups that cannot be scaled easily, effectively preventing their practical application for the quantum information purposes. In the present paper we demonstrate a simple, scalable squeezed vacuum source which consists in essence of a continuous-wave (cw) diode laser and an atomic rubidium vapor cell. Easy and inexpensive generation of multiple Einstein-Podolsky-Rosen-entangled modes lies at hand, making them available to a variety of experiments and applications.

We make use of a phenomenon known as cross-phase modulation squeezing Haus95 : if linearly polarized light propagates through a nonlinear medium in which elliptical light would undergo polarization self-rotation, the orthogonal vacuum mode gets squeezed. Self-rotation and the associated squeezing effect are present in any Kerr medium but are especially strong in atomic vapor due to its large optical nonlinearity near resonance. To achieve squeezing it is sufficient to use a moderately intense cw diode laser as a pump. Pulsed lasers or cavities are not needed and the diode laser does not even have to be frequency stabilized.

## Theory

The theory of squeezing via polarization self-rotation for the single mode case has been given in Ref. Matsko02 . Here we present its brief quasiclassical overview. A detailed multimode description including the degrading effects of resonance fluorescence, phase mismatch and absorption will be published elsewhere.

Interaction with an elliptically polarized light field can cause an initially isotropic $`\chi ^{(3)}`$-nonlinear medium to become circularly birefringent: the two circular components of different intensities propagate with different phase velocities. Specifically in atomic vapor, the mechanisms leading to this birefringence include optical pumping, ac Stark shifts, as well as ac Stark shift induced quantum beats and Bennett structures Rochester01 . Propagating through a cell of length $`l`$, the polarization ellipse of light will rotate by an angle $`\phi `$ which is proportional to the ellipticity $`ϵ1`$:

$$\phi =glϵ,$$
(1)
with the self-rotation parameter $`g`$ dependent on the incident light intensity and frequency.

Consider a monochromatic, elliptically polarized light field of frequency $`\omega `$ propagating in the $`z`$-direction. The complex field amplitudes in the two linear polarization components are given by $`_x`$ and $`_y`$, with $`|_x||_y|`$. Assuming the “pump” amplitude $`_y`$ real, we decompose the “probe” $`x`$-field into two real quadrature components in and out of phase with the pump field: $`_x(z)=X_x(z)+iP_x(z)`$. Then the ellipticity $`ϵ`$ of the laser light is given by

$$ϵ\frac{P_x(z)}{_y(z)}.$$
(2)
The resulting self-rotation by $`\phi 1`$ causes the pump field to contribute its fraction $`\phi _y`$ into the $`x`$-polarization, namely into the $`X_x`$ quadrature which is in phase with $`_y`$:

$`X_x(l)`$ $`=`$ $`X_x(0)+\phi _y=X_x(0)+glP_x(0)`$ (3)
$`P_x(l)`$ $`=`$ $`P_x(0),`$ (4)
$`_y`$ and $`ϵ`$ remain effectively unchanged.

We now assume the $`x`$-polarized light field as classical noise described by a probability distribution in the quadrature plane (Fig. 1). To describe slow fluctuations, we regard the propagation of a light field described by $`X_x(z),P_x(z)`$ and $`_y`$ through the medium under stationary conditions. A random quadrature component $`P_x(0)`$ results in a small random ellipticity according to Eq. (2). The propagation through the cell will displace the point ($`X_x,P_x`$) in the phase space as described in Eqs. (3) and (4). The phase space probability distribution experiences a linear shear, resulting in reduction of the noise in some phase quadratures below the original level.

Although the above treatment is purely classical, it is completely identical in application to the quantum quadrature noise. Considering the Heisenberg evolution of the field-quadrature operators $`\widehat{X}_x(z)`$ and $`\widehat{P}_x(z)`$ we notice that the evolution described by Eqs. (3) and (4) can be associated with the interaction Hamiltonian $`\widehat{H}=\mathrm{}gc\widehat{P}_x^2/2`$, $`c`$ being the speed of light in the medium. Under a second order Hamiltonian, each point of the Wigner function — the quantum-mechanical analogue of the phase-space probability distribution — evolves according to the classical equations of motion Schleich . The incoming symmetric Gaussian Wigner function of the vacuum is therefore transformed into an elliptical Wigner function of a minimum uncertainty squeezed vacuum state.

In a homodyne measurement such a state results in a quadrature noise depending on the local oscillator phase $`\chi `$ (cf. Matsko02 ):

$`\mathrm{\Delta }E_x(\chi ,l)^2`$
$`={\displaystyle \frac{_0^2}{4}}\left((1\alpha l)(12gl\mathrm{sin}\chi \mathrm{cos}\chi +g^2l^2\mathrm{cos}^2\chi )+\alpha l\right).`$
Here $`\frac{_0^2}{4}`$ is the standard quantum noise limit (SQL) associated with the vacuum state, and we have taken into account absorption of a small fraction $`\alpha l1`$ of the squeezed light after passing the medium. This treatment of absorption is approximate because some of the light energy is actually lost in the atomic vapor along the propagation path. However, it permits to account for one essential feature of our experimental result: the minimum uncertainty property is lost when $`gl>0`$ and $`\alpha l>0`$.

## Experimental apparatus

The scheme of our experimental setup is shown in Fig. 2. A free running external cavity diode laser (Toptica DL-100) delivered 35 mW of optical power in a single mode at a wavelength of 780 nm, tunable across the rubidium $`D2`$ line. The beam was shaped with a pair of cylindrical lenses and a polarizing beam splitter (PBS) purified the linear polarization. The beam was focused by a lens ($`f=800`$mm) into a heated and magnetically shielded $`75`$mm rubidium cell which had optical quality anti-reflection coated windows and was filled with the nearly pure isotope 87Rb. With a focal beam diameter of 600 $`\mu `$m we operated in a high saturation regime.

To determine the polarization self-rotation parameter $`gl`$ a fine-adjustable $`\lambda /4`$ plate was placed in front of the cell. Its rotation angle in radians was equal to the generated ellipticity $`ϵ`$. The outputs of the second PBS were directly fed into two photodiodes. A $`\lambda /2`$ plate placed after the cell was adjusted so that the DC photodiode signals $`S_1`$ and $`S_2`$ were exactly equal when the laser was detuned far from the resonance. This “balanced polarimeter” was sensitive to the polarization ellipse rotation while being insensitive to changes in ellipticity. The rotation angle $`\phi 1`$ can be determined from the photodiode signal difference Rochester01 :

$$\phi =\frac{S_1S_s}{2(S_1+S_2)}$$
(6)
By recording $`S_1`$ and $`S_2`$ while scanning over the rubidium $`D2`$ line for many different initial small ellipticities $`ϵ1`$, the self-rotation parameter $`gl`$ for each laser detuning was determined by polynomial fitting. Fig. 3 shows $`gl`$ and the absorption coefficient $`\alpha l`$ for our working temperature of $`70^{}`$C. Note that a relatively high transmission on resonance is due to very high saturation of the atomic transition — in the linear small signal regime, the rubidium vapor is optically dense.

To measure the quantum field quadrature noise we employed the standard balanced homodyne detection technique. The $`y`$-polarized component of the laser beam emerging from the cell was used as the local oscillator. To this end it was first separated from the $`x`$-polarized probe mode by means of a polarizing beam splitter and its polarization was rotated by $`90^{}`$. It was then overlapped with the squeezed mode on a symmetric beam splitter. To achieve good mode matching, an auxiliary $`\lambda /2`$ plate placed directly after the cell was set to provide equal splitting into both polarization modes. The resulting Mach-Zehnder interferometer was carefully aligned for a maximum visibility of $`V=0.98\pm 0.01`$ corresponding to a mode matching efficiency of $`V^20.96`$. The auxiliary $`\lambda /2`$ plate was then removed.

A piezoelectric element was used to scan the local oscillator phase and thus the field quadrature whose fluctuations were measured. The intensities at the beam splitter outputs were measured with Si PIN photodiodes (Hamamatsu S3883 with a nominal quantum efficiency of 91%). The AC signals of the diodes were amplified with a low noise amplifier with a cutoff frequency of approximately 50 MHz and subtracted by a hybrid junction (H-9 from M/A-COM). The noise was measured with a spectrum analyzer operated in the zero span mode at a set of radio frequencies $`\mathrm{\Omega }_{RF}`$ between 3 and 30 MHz with a resolution bandwidth of 1 MHz.

## Results and discussion

The squeezing homodyne measurements were performed at a finite radio-frequency $`\mathrm{\Omega }_{RF}`$, as usual, because of the high optical and electronic noise near DC. Fig. 4 shows the phase dependent noise of the $`x`$-polarization mode which contains the squeezed vacuum. It is compared to the SQL which can be measured by blocking the squeezed vacuum path. The phase dependent noise minima fall below the SQL exhibiting a squeezing of $`(0.85\pm 0.05)`$dB. This corresponds to a squeezing of $`(1.23\pm 0.07)`$dB at the cell output when corrections for linear losses, inefficient photodiodes and electronic photodetector noise are made. The geometric mean noise of the squeezed mode, however, strongly exceeds the SQL so the observed ensemble is not a minimum uncertainty state. The excess noise cannot be explained by absorption alone and is presumably due to resonance fluorescence into the optical mode probed.

Fig. 6 shows the squeezing as a function of the laser detuning and the radio frequency $`\mathrm{\Omega }_{RF}`$. The variation of squeezing with the laser detuning follows approximately the self-rotation curve in Fig. 3. For small $`\mathrm{\Omega }_{RF}`$ the squeezing is reduced due to extra noise from resonance fluorescence which scatters light into our vacuum mode. This extra noise was measured by misaligning the mode matching until the phase dependence vanished and maximizing the spatial overlap between the pump and probe modes in the cell by maximizing the noise. Fig. 6 shows the extra noise normalized to the vacuum noise. Its spectral width approximates the natural linewidth of $`\mathrm{\Gamma }=6`$MHz resonancefluorescence , its dependence on the laser detuning is similar to the absorption line (Fig. 3). For $`\mathrm{\Omega }_{RF}\mathrm{\Gamma }`$ the excess noise becomes negligible. Maximum squeezing is achieved for $`\mathrm{\Omega }_{RF}5`$MHz to $`10`$MHz. The reason for the squeezing to degrade with higher $`\mathrm{\Omega }_{RF}`$ might be a growing phase mismatch between the pump and the probe expected in the neighborhood of an atomic resonance.

## Summary and outlook

We reported the generation of squeezed optical vacuum states in rubidium vapor with a simple setup. A squeezing of 0.85 dB was achieved which is among the best atomic vapor squeezing results accomplished. The squeezing of 6 to 8 dB predicted in Matsko02  was not reached due to the degrading effects of self-focusing, phase mismatch and resonance fluorescence. Further theoretical investigation of these effects will help us optimize the experimental parameters. The self-rotation parameter and the associated squeezing might increase by using a cell with a moderate amount of buffer gas Novikova02 . We are optimistic that a noticeable improvement of the squeezing is possible. By splitting the laser and directing it into several rubidium cells, numerous squeezed states, which are phase locked with respect to each other, can be generated easily.

## Acknowledgement

We acknowledge helpful discussions with K. P. Marzlin, J. H. Shapiro and I. Novikova, as well as financial support from the Deutsche Forschungsgemeinschaft and the Optik-Zentrum Konstanz.


text:
# The Vortex State in a Strongly Coupled Dilute Atomic Fermionic Superfluid

## Abstract

We show that in a dilute fermionic superfluid, when the Fermions interact with an infinite scattering length, a vortex state is characterized by a strong density depletion along the vortex core. This feature can make a direct visualization of vortices in Fermionic superfluids possible.

The existence of stable vortex states is one of the most spectacular manifestations of superfluidity in both Bose and Fermi systems. In Bose dilute atomic gases isolated quantized vortices  and subsequently large arrays of vortices  have been observed already. The quest for the observation of superfluidity in dilute Fermi gases started as soon as the first experimental evidence of an atomic dilute Fermi degenerate gas was published . In Bose systems methods for establishing unambigously both the superfluidity and the quantization of vortices exists , while most of the suggestions made so far for Fermi systems are rather indirect. The presence of vortices in BEC was confirmed by analyzing the density variations of an expanding BEC cloud . The Bose-Einstein condensation (BEC) creation was confirmed by studying the character of the expansion of the atomic cloud after the trap was removed . A recent experimental result in Fermi systems  suggests that a hydrodynamic expansion of a Fermi superfluid is a plausible scenario .

In Fermi systems significant density variations due to the presence of vortices are not expected . We have reasons to expect however, that under certain conditions the density variations induced by the presence of one or more vortices could be akin to those in the case of BEC, as a similar feature was recently established in the analysis of the spatial structure of a vortex in low density nuclear matter . At densities significantly smaller than nuclear saturation densities, the superfluid gap in homogeneous neutron matter can attain values rather large by normal standards, $`\mathrm{\Delta }0.25\epsilon _F`$, where $`\mathrm{\Delta }`$ is the value of the gap and $`\epsilon _F`$ is the Fermi energy. Even though the gap is still smaller than the Fermi energy, such values proved sufficient in the case of low density neutron matter to lead to major density depletions in the vortex core. As Pitaevskii and Stringari note , the observation of quantized vortices in a dilute Fermi gas would provide the ultimate proof that the system has undergone a transition to a superfluid state. A vortex is just about the only phenomenon in which a true stable superflow is created in a neutral system. Other phenomena would only somewhat indirectly be affected by the onset of superfluidity.

Recently, a new ab initio calculation of the properties of low density nuclear matter became available . A particular result of this analysis, which is valid in principle for any Fermion system, concerns the properties of a two fermion species interacting with an infinite scattering length. These authors have shown that the energy per particle of such a normal Fermi system is

$$_N=\alpha _N\frac{3}{5}\epsilon _F,$$
(1)
with $`\alpha _N0.54`$. This result was obtained for a system of fermions interacting with a short range attractive potential with a zero–energy bound state. As long as $`k_Fr_01`$, where $`r_0`$ is of the order of the radius of the potential, one expects on general grounds that the energy per particle of such a system is proportional to the energy per particle of a non–interacting system. The scattering length $`a`$ and the effective range $`r_0`$ parameterize the low energy behavior of the $`s`$–wave scattering phase of two particles $`k\mathrm{cot}\delta (k)1/a+r_0k^2/2,`$ where $`k`$ is the wave vector of the relative motion. A second result of this ab initio calculation concerns the energy per particle of the superfluid phase of such a system , namely that

$$_S=\alpha _S\frac{3}{5}\epsilon _F,$$
(2)
with $`\alpha _S0.44`$. A simple estimate of the corresponding value of the superfluid gap, using for the condensation energy the weak coupling BCS value $`3\mathrm{\Delta }^2/8\epsilon _F`$ leads to $`\mathrm{\Delta }0.4\epsilon _F`$. Such large values for the gap $`\mathrm{\Delta }`$ are not compatible with the BCS weak coupling limit, when $`k_F|a|1`$ and $`a<0`$. It is well known that in this limit the BCS approximation  leads to a too large value of the gap and that polarization corrections lead to a reduction of the gap , namely to

$$\mathrm{\Delta }=\left(\frac{2}{e}\right)^{7/3}\epsilon _F\mathrm{exp}\left(\frac{\pi }{2k_Fa}\right).$$

A recent analysis  of the vortex state in a dilute superfluid Fermi gas, using the simple BCS value for the gap (which exceeds by a factor of $`2.2`$ the true gap value) shows a relatively modest density depletion in the vortex core of about 10% at most.

The possibility that the value of the superfluid gap can attain large values was raised more than two decades ago in connection with the BCS $``$ BEC crossover . One can imagine that one can increase the strength of the two–particle interaction in such a manner that at some point a real two–bound state forms, and in that case $`a\mathrm{}`$. By continuing to increase the strength of the two–particle interaction, the scattering length becomes positive and starts decreasing. When $`ar_0`$ and $`a>0`$ the energy of the two–particle bound state is $`ϵ_2\mathrm{}^2/ma^2.`$ A dilute system of fermions, when $`nr_0^31`$, will thus undergo a transition from a weakly coupled BCS system, when $`a<0`$ and $`a=𝒪(r_0)`$, to a BEC system of tightly bound Fermion pairs, when $`a>0`$ and $`a=𝒪(r_0)`$ again. In the weakly coupled BCS limit the size of the Cooper pair is given by the so called coherence length $`\xi \frac{\mathrm{}^2k_F}{m\mathrm{\Delta }},`$ which is much larger than the interparticle separation $`\lambda _F=2\pi /k_F`$. In the opposite limit, when $`k_Fa1`$ and $`a>0`$, and when tightly bound pairs/dimers of size $`a`$ are formed, the dimers are widely separated from one another. Surprisingly, these dimers also repel each other with an estimated scattering length $`0.6\mathrm{}2a`$  and thus the BEC phase is also (meta)stable. The bulk of the theoretical analysis in the intermediate region where $`k_F|a|>1`$ was based on the BCS formalism  and thus is highly questionable. All these authors have considered only the simple ladder diagrams in the particle–particle channel. The inclusion of the additional ”bosonic” degrees of freedom , which represent nothing else but a two–atom bound state, falls into the same approximation scheme, which includes only ladder diagrams in the atom–atom channel. Even the simplest polarization corrections have not been included into this type of analysis so far. In particular, it is well known that in the low density region, where $`a<0`$ and $`k_F|a|1`$ the polarization corrections to the BCS equations lead to a noticeable reduction of the gap . Only a truly ab initio calculation could really describe the structure of a many Fermion system with $`k_F|a|1`$. In the limit $`a=\pm \mathrm{}`$, when the two–body bound state has exactly zero energy, and if $`k_Fr_01`$, one can expect that the energy per particle of the system is proportional to $`\epsilon _F=\mathrm{}^2k_F^2/2m`$, as it was recently confirmed by the variational calculations of Refs. .

As in any other Fermion systems , the knowledge of the energy per particle as a function of density in the homogeneous phase, Eq. (1), allows us to construct the normal part of the Local Density Approximation (LDA) for the energy density functional (EDF). The extension of the LDA to superfluid systems presented in Refs.  in conjunction with Eq. (2) lead us to a rather well defined EDF, appropriate for Fermion systems with infinite scattering length, namely

$`(𝒓)n(𝒓)={\displaystyle \frac{\mathrm{}^2}{m}}\left[{\displaystyle \frac{m}{2m^{}}}\tau (𝒓)+\beta n(𝒓)^{5/3}+\gamma {\displaystyle \frac{|\nu (𝒓)|^2}{n(𝒓)^{1/3}}}\right],`$ (3)
$`n(𝒓)={\displaystyle \underset{\alpha }{}}|v_\alpha (𝒓)|^2,\tau (𝒓)={\displaystyle \underset{\alpha }{}}|\mathbf{}v_\alpha (𝒓)|^2,`$ (4)
$`\nu (𝒓)={\displaystyle \underset{\alpha }{}}v_\alpha ^{}(𝒓)u_\alpha (𝒓).`$ (5)
Here $`n(𝒓)`$ and $`\nu (𝒓)`$ are the normal and anomalous densities and $`\tau (𝒓)`$ is the kinetic energy density, all of them expressed through the quasiparticle wave functions $`(u_\alpha (𝒓),v_\alpha (𝒓))`$. The summation over $`\alpha `$ should be interpreted as a sum or integral when appropriate, over all quasiparticle states with the Bogoliubov–de Gennes eigenvalues $`E_\alpha >0`$. The dimensionless couplings $`m/m^{}`$, $`\beta `$ and $`\gamma `$ have to be chosen so as to reproduce in the case of infinite homogeneous matter the expressions (12) for the energy per particle. Unfortunately the parameterization of the EDF is not unique and its form has to be restricted using some additional considerations. One possible choice corresponds to $`m/m^{}=\alpha _N`$, $`\beta =0`$ and $`\gamma 6.64`$. We shall refer to this case as the parameter set I. Another possible choice of parameters, referred to as the parameter set II, corresponds to $`m/m^{}=1`$, $`\beta =3(\alpha _N1)(3\pi ^2)^{2/3}/5`$ and $`\gamma 7.20`$. Since the anomalous density $`\nu (𝒓)`$ is a diverging quantity, a regularization procedure of the above expression (3) is required. This was performed according to the formalism described in great detail by us elsewhere . The renormalization procedure amounts to replacing $`\gamma `$ with a new well defined running coupling $`\gamma _{\mathrm{𝑒𝑓𝑓}}(𝒓)`$. Even though in this renormalization procedure there is an explicit energy cutoff appearing in the formalism, no observable is affected by its presence, when this energy cutoff is chosen appropriately. Up to the overall factor $`\mathrm{}^2/m`$, the energy density in Eq. (3) has the overall scaling $`(𝒓)[n(𝒓)]^{5/3}`$, which is expected from dimensional arguments for a system with an infinite scattering length and $`k_Fr_01`$. In principle, there are an infinite number of possible local energy densities satisfying this requirement. In particular, one could have considered for the superfluid part of the functional for example a contribution of the form $`|\nu (𝒓)|^{5/3}`$. Since a term $`|\nu (𝒓)|^2`$ works in both BCS and BEC limits and since the emerging pairing field $`\mathrm{\Delta }(𝒓)`$ has the simplest form and the renormalization procedure is also the simplest when using the EDF Eq. (3), we did not consider explicitly any other forms. Moreover, nobody ever suggested other forms for the contribution of the pairing correlations to the energy density in the past, as far as we are aware. We do not expect however any qualitative changes in our results resulting from considering other forms for the superfluid contribution to the energy density.

We looked for a self–consistent solution of the Bogoliubov–de Gennes equations corresponding to the EDF in Rel. (3), with a vortex in the pairing field along the $`Oz`$–axis and $`\mathrm{}/2`$ of net angular momentum per particle, namely $`\mathrm{\Delta }(𝒓)=\mathrm{\Delta }(\rho )\mathrm{exp}(i\varphi )`$, where $`𝒓=(\rho ,\varphi ,z)`$ are cylindrical coordinates. The mean–field depends only on the radial coordinate $`\rho `$. The quasiparticle wave functions have the structure

$$\left(\begin{array}{c}u_\alpha (𝒓)\hfill \\ v_\alpha (𝒓)\hfill \end{array}\right)=\left(\begin{array}{c}u_s(\rho )\mathrm{exp}\left[in\varphi +ikz\right]\hfill \\ v_s(\rho )\mathrm{exp}\left[i(n1)\varphi +ikz\right]\hfill \end{array}\right),$$
(6)
where $`s`$ labels the radial part of the quasiparticle states, $`n`$ is an integer, $`k`$ is the wave vector of the quasiparticle state along the vortex symmetry axis $`Oz`$ and $`\alpha =(s,k,n)`$. The presence of the vortex implies a net (super)flow of the Fermion superfluid around the vortex axis with the velocity profile

$$\mathrm{V}_s(\rho )\widehat{𝒆}_\mathit{\varphi }=\frac{\mathrm{}}{m\rho }\frac{1}{n(\rho )}\underset{\alpha }{}v_\alpha ^{}(𝒓)i\widehat{𝒆}_\mathit{\varphi }\frac{}{\varphi }v_\alpha (𝒓),$$
(7)
where $`\widehat{𝒆}_\mathit{\varphi }=(y,x,0)/\rho `$,

The most salient feature emerging from this analysis is the behavior of the density as a function of the radial coordinate $`\rho `$, the prominent density depletion on the vortex axis, see Fig. 1. The density at the vortex core is lowered to a value of $`0.3\rho _{\mathrm{}}`$, while the pairing field vanishes at the vortex axis as expected, see Fig. 2. By comparing the actual flow profile with that of an ideal vortex, see Fig. 3, one can see that only the fraction $`f_s(\rho )=\mathrm{V}_s(\rho )/V_v(\rho )`$ of the system is superfluid at a distance $`\rho `$ from the vortex axis. We have performed also calculations by imposing the artificial constraint that the pairing field in the homogeneous phase has a larger value than the one expected from the self–consistent solution, corresponding to the EDF Eq. (3). In that case the density depletion at the vortex core becames even larger. In hindsight this result could have been expected. Large values of the pairing field correspond to the formation of atom pairs/dimers of relatively small sizes. When these dimers are relatively strongly bound and when they are also widely separated from one another, they undergo a Bose–Einstein condensation. For a vortex state in a 100% BEC system the density at the vortex axis vanishes identically. Therefore, by increasing the strength of the two–particle interaction, the Fermion system simply approaches more and more an ideal BEC system, for which a density depletion of the vortex core is expected.

The value of the coupling constant, and correspondingly the value of the scattering length can nowadays be routinely controlled in experiments , following earlier theoretical suggestions . One can then expect that it would be feasible to create a BCS state as in the experiment of Ref. , stir the system so as to create one or several vortices in a manner similar to the experiments with dilute Bose gases  and change the ambient magnetic field so that the system is brought close to the Feshbach resonance. The first quantized vortex will appear only at some critical angular velocity, if the system is gradually spun up. The vortex quantization could be demonstrated for example as in Ref. . At zero temperature, the single vortex state is the lowest energy state with the total angular momentum $`N\mathrm{}/2`$ (here $`N`$ is the total number of particles) and its decay into lower energy states with total angular momenta differing by a few number of $`\mathrm{}/2`$ is strongly suppressed. The next vortex will appear at a larger critical angular velocity, and so forth. The density profiles of the vortices now develop significant density depletions at the cores, which could hopefully be subsequently visualized, if they survive the expansion when the trap is removed. As discussed in Ref. , a fluid with a power law equation of state (the present case) will have a simple scaling expansion, see also Refs. . As a consequence, the density depletion along the vortex core should survive the atomic cloud expansion upon trap release. Note that a system in which the scattering length is infinite is highly unusual, as, even when it expands and its density decreases, the relative role of the interaction does not change, unlike any other interacting system, and such a system remains strongly coupled at all densities.

By increasing further the strength of the interaction, two atoms will form a bound state. In the limit when the size of this bound state becomes significantly smaller then the mean interparticle separation the system becomes a BEC of atom pairs. These pairs repel each other with an estimated scattering length $`0.6\mathrm{}2a`$ . The BEC phase can in principle be described as well within the framework of the same formalism outline above. Such an approach would be unreasonable, since an accurate description of the relatively strongly bound two–atom system in terms of extended mean–field quasiparticle states would require an exceedingly large number of quasiparticle states. The ratio of the quantum states required for a reasonable accuracy to the number of particles per particle would be at least of order $`1/(na^3)1`$. The natural description in this limit is in terms of bosonic dimers, which weakly repel each other. The magnitude of the coherence length in the BEC phase is exceeding considerably the mean separation between pairs, $`\xi _b\sqrt{na}1/n^{1/3}`$. This coherence length $`\xi _b`$ determines the radius of the vortex core. Thus along with the two–atom scattering length $`a`$, the size of the vortex core could be controlled as well by means of the Feshbach resonance.

In summary, we have shown that the vortex in a dilute atomic superfluid Fermi gas interacting with an infinite scattering length develops a significant density depletion along its core, which could be visualized after expanding the atomic cloud.

We thank P.F. Bedaque, G.F. Bertsch, J. Carlson, E.K.U. Gross, N. Hayashi, V.R. Pandharipande, D.S. Petrov, B. Spivak, S. Stringari and A. Svidzinsky for discussions and/or correspondence. This work benefited from partial financial support under DOE contract DE–FG03–97ER41014.