text:
# Loading a vapor cell magneto-optic trap using light-induced atom desorption

## I Introduction

The evaporative cooling techniques used to achieve Bose-Einstein condensation in atomic gases  rely on loading large numbers of atoms into magnetic traps with long trap lifetimes. The approach originally taken by Anderson et al.  was to load Rb atoms into a vapor cell  magneto-optic trap (MOT)  and to subsequently transfer atoms into a magnetic trap located in the same cell. Large numbers of atoms and long lifetimes were achieved by optimizing the Rb partial pressure and by working with long MOT loading times.

We have found a simple way to improve such a setup by modulating the vapor pressure such that it is high for initial trap loading and then low again in order to achieve long lifetimes in a magnetic trap. The technique requires the use of a white light source (WLS) with radiation incident upon the inner walls of the vapor cell. When such a light source is turned on, Rb atoms that coat the inner walls of the stainless steel vacuum chamber are quickly desorbed and the Rb vapor pressure suddenly increases. The vapor pressure soon returns to equilibrium after the WLS is turned off. This enables loading large numbers of atoms into the MOT in a relatively short amount of time, while preserving the low pressures required for long magnetic trap lifetimes. The WLS method that we describe here is a possible alternative to the double-chamber techniques  and Zeeman slowing techniques  currently used to capture atoms before evaporatively cooling in a magnetic trap. In our experiments, the WLS method is used in the manner described here for the achievement of BEC in a vapor cell, where the WLS frees us from environmentally induced variations in vapor pressure; for example, regardless of chamber temperature, we can load large numbers of atoms into our MOT and achieve BEC .

Light-induced atom desorption (LIAD) has previously been used to obtain optically thick Na and Rb vapors in cells made of glass, pyrex, and sapphire . In most of these experiments, the inner walls of the vapor cells were coated with paraffin or silane in order to enhance the LIAD efficiency by reducing the alkali atom adsorption energy . In our work, an optically thick vapor was not required. Since we did not need desorption rates characteristic of coated cells, we could desorb atoms directly from stainless steel.

## II Background

We first review the basic mechanisms involved in the operation of a vapor cell MOT, lucidly described in Ref. , in order to understand the gains available with the LIAD method. In a vapor cell MOT, atoms are loaded into the MOT at a rate $`R`$. This rate depends upon the size and intensity of the laser cooling and trapping beams and the Rb partial pressure. Atoms with velocities below a critical velocity will be captured by the trap. Atoms are also lost from the trap due to collisions, limiting the number that can be loaded into the MOT. The rate equation for the number, $`N`$, of trapped atoms is given by

$$\frac{dN}{dt}=RN\left(\frac{1}{\tau _\text{b}}+\frac{1}{\tau _{\text{Rb}}}\right)\beta n^2𝑑V,$$
(1)
where $`1/\tau _\text{b}`$ is the trap loss rate due to collisions with background gas atoms and $`1/\tau _{\text{Rb}}`$ is the loss rate determined by collisions with untrapped Rb atoms. The trap density, $`n`$ in the volume integral, contributes to density-dependent losses within the trap with a loss coefficient of $`\beta `$. The loss rate $`1/\tau _\text{b}`$ is proportional to the pressure of the background gas, and like $`R`$, $`1/\tau _{\text{Rb}}`$ is proportional to the Rb partial pressure. In the absence of density-dependent collisional losses ($`\beta =0`$) , and with a Rb partial pressure that is much higher than the background pressure ($`1/\tau _{\text{Rb}}1/\tau _\text{b}`$), the rate equation becomes

$$\frac{dN}{dt}=R\frac{N}{\tau _{\text{Rb}}}.$$
(2)
The limiting number, $`N_{\text{lim}}`$, that can be loaded into the MOT is obtained when the increase in number due to loading balances the loss due to collisions. At this point, $`N=N_{\text{lim}}`$ and $`dN/dt=0`$, yielding

$$N_{\text{lim}}=R\tau _{\text{Rb}},$$
(3)
independent of the Rb partial pressure .

Frequently, the background-gas collisions can not be neglected, and the total number reached will be less than $`N_{\text{lim}}`$. The maximum number $`N_{\text{max}}`$ that can be captured for a given Rb partial pressure will then be given by

$$N_{\text{max}}=N_{\text{lim}}\frac{\tau _{\text{MOT}}}{\tau _{\text{Rb}}}=R\tau _{\text{MOT}},$$
(4)
where

$$\frac{1}{\tau _{\text{MOT}}}=\frac{1}{\tau _\text{b}}+\frac{1}{\tau _{\text{Rb}}}.$$
(5)
If the trap starts filling at time $`t=0`$, the number of atoms in the MOT at any point in time is given by

$$N(t)=N_{\text{max}}\left[1\mathrm{exp}\left(\frac{t}{\tau _{\text{MOT}}}\right)\right].$$
(6)
Because of the appearance of $`\tau _{\text{MOT}}`$ as the time constant in the exponential, we define $`\tau _{\text{MOT}}`$ as the ”MOT loading time.”

The lifetime of a *magnetic trap* in the same chamber also depends upon the collision rate of trapped atoms with background atoms. Thus the magnetic trap lifetime $`\tau `$ is proportional to $`\tau _{\text{MOT}}`$. We express this proportionality as $`\tau =\tau _{\text{MOT}}/\alpha `$, where for our traps, $`\alpha 4`$. For evaporative cooling experiments, where large numbers of atoms and long magnetic trap lifetimes are both necessary, the product of total number $`N_{\text{max}}`$ and magnetic trap lifetime is the critical parameter to maximize . Because of the relationship between $`\tau `$ and $`\tau _{\text{MOT}}`$, we can alternatively view this requirement as maximizing the product of $`N_{\text{max}}`$ and $`\tau _{\text{MOT}}`$. We must therefore find the optimum Rb partial pressure for a given background pressure. Multiplying Eq. 4 by $`\tau _{\text{MOT}}`$ leads to maximization of $`\tau _{\text{MOT}}^2/\tau _\text{b}`$ ($`N_{\text{lim}}`$ is independent of vapor pressure). Under optimal conditions, with constant Rb partial pressure, $`N_{\text{max}}\tau _{\text{MOT}}`$ is maximized for $`\tau _\text{b}=\tau _{\text{Rb}}=\tau _{\text{MOT}}/2`$ and hence $`N_{\text{max}}=N_{\text{lim}}/2`$ and $`N_{\text{max}}\tau _{\text{MOT}}=N_{\text{lim}}\tau _\text{b}/4`$.

However, we can further improve the number-lifetime product (which from now on we will generally designate as $`N\tau `$) by permitting a modulation of the Rb vapor pressure. If the Rb partial pressure is *temporarily* increased until the trap contains the maximum possible number of atoms ($`N=N_{\text{max}}=N_{\text{lim}}`$), at which point the Rb vapor is suddenly reduced to a negligible level ($`\tau _{\text{MOT}}\tau _\text{b}`$), an increase of a factor of 4 in $`N\tau `$ will be realized. Furthermore, the time needed to load the MOT is significantly shortened when $`\tau _{\text{Rb}}\tau _\text{b}`$ during the loading interval, increasing the repetition rate of the experiment.

The goal of our experiment was to realize gains in $`N\tau `$ by modulating the Rb vapor pressure in the described manner with the white light source, thus improving conditions for evaporative cooling and obtaining BEC.

## III Experimental setup and measurement techniques

Our stainless steel vacuum chamber consisted of a vapor cell atom trapping chamber with indium-sealed windows, a liquid nitrogen-filled cold finger which protruded into the chamber, and a Rb cold finger at $`0^{}`$C, as shown in Fig. 1. The vacuum in the chamber was maintained by a Ti-sublimation pump and an ion pump. We maintained a Rb vapor in the chamber by slightly opening a valve between the chamber and a Rb cold finger. This replenished Rb that was pumped out of the chamber.

The MOT was constructed using a forced dark SPOT technique : a 4 mm opaque spot was placed in the center of the path of the repumping laser light, and was imaged onto the region in the chamber where the trap was formed. Another laser beam filled the hole in the repumping beam, and was used to optically pump trapped atoms into a dark state. This technique reduced the trap loss rate due to light-assisted, density-dependent collisions between trapped atoms. The Rb trapping light was tuned 13 MHz below the $`5S_{1/2},F=25P_{3/2},F^{}=3`$ transition, and was provided by six 23 mW/cm2, 1.2 cm diameter laser beams. The 2.7 mW/cm2 repumping laser beam was tuned 15 MHz below the $`F=1F^{}=2`$ transition, and the 9 mW/cm2 forced optical pumping light was tuned to the $`F=2F^{}=2`$ transition. The number of atoms in the trap was measured by detecting light scattered by the trapped atoms. This was done by turning off the $`F=2F^{}=2`$ light for $`50`$ ms and filling the hole in the repumping beam with a separate bypass repumping beam such that the trapped atoms were made bright by scattering light from the trapping beams. A fraction of the light scattered by the trapped atoms was collected and focused onto a calibrated photomultiplier tube. Loading rates ($`R`$) and MOT loading time constants ($`\tau _{\text{MOT}}`$) were measured by detecting the number of atoms at sequential points in time as the trap filled.

The white light used to enhance trap loading was provided by a fiber optic illuminator, consisting of a halogen bulb with variable power and a fiber bundle which pointed the light into the vapor cell. The coupling of the light from the bulb into the fiber gave a maximum intensity onto the inner vapor cell wall of $`10`$ W/cm2. The WLS was switched on and off electronically.

To measure $`\tau _{\text{MOT}}`$, we measured the number of atoms loaded into the trap as a function of time both with and without the WLS. The loading curves were exponential in time, as expected for a trap without light-assisted losses. Typical filling curves are shown in Fig. 2(a) for various WLS intensities. In the figure, the curve representing the fastest filling rate, with a WLS intensity of $`10`$ W/cm2, shows a loading time constant of $`\tau _{\text{MOT}}=67`$ s and a maximum number of $`1.3\times 10^8`$ atoms as determined by the exponential fit. Without the WLS, the loading time constant was $`\tau _{\text{MOT}}=538`$ s and the maximum number was $`2\times 10^7`$ atoms. Values of number loaded and loading time constants for the curves shown in Fig. 2(a) are given in the second and third columns of Table I
.

A key factor to consider in optimizing $`N\tau `$ using the WLS scheme is the time for the vapor pressure to return to lower equilibrium values once the WLS is switched off. We define this time as the vapor pressure recovery time. A liquid nitrogen cold finger in the main chamber was used to decrease the Rb vapor pressure and shorten the recovery time after the WLS was switched off. In our cell, the cold finger had little effect on the background gas pressure, but shortened the recovery time by a factor of $`2`$. Furthermore, our experimental timing sequence consisted of a MOT loading phase with the WLS switched on, followed by a MOT holding phase, during which the atoms were held in a MOT with the WLS switched off. This enabled us to keep a large number of atoms trapped while waiting for the vapor pressure to decrease before extinguishing the trapping light.

In order to evaluate vapor pressure recovery times, we measured the dependence of loading rates on time just after the WLS was switched off. For the data shown in Fig. 2(b), the WLS was left on until the Rb partial pressure reached a saturated level. The WLS light was then turned off, and the number of atoms loaded into a MOT in 5 seconds was repeatedly measured. After each measurement, the MOT light was kept off for 5 s, and then the MOT started filling again for the subsequent 5 s filling rate measurement. This set of measurements indicated the speed at which $`\tau _{\text{MOT}}`$ and the Rb vapor pressure could recover after the WLS was turned off, and demonstrated that the recovery time was roughly independent of the WLS intensity. The fastest loading rate shown with the WLS on was $`8\times 10^5`$ atoms/s, and with the light off was $`2.7\times 10^4`$ atoms/s. Each loading rate vs. time curve in Fig. 2(b) was fit with a sum of two decaying exponential curves. The time constants for the loading rate to return to lower equilibrium values were $`8`$ s for the fast recovery time ($`\tau _{rec,1}`$), and between 113 s and 167 s for the slower recovery time ($`\tau _{rec,2}`$). Table I contains a list of recovery times.

To help evaluate the vapor cell performance, the values for $`N_{\text{lim}}`$ and $`\tau _\text{b}`$ were estimated by measuring $`N`$ and $`\tau _{\text{MOT}}`$ for various Rb partial pressures. Experimentally, we varied the Rb partial pressure by adjusting the intensity of the WLS. We estimated $`\tau _\text{b}700`$ s and $`N_{\text{lim}}=1.9\times 10^8`$ atoms for our operating parameters by linear extrapolation with our data.

## IV Model

We now describe a detailed model for determining the numbers and lifetimes of traps loaded with the WLS to demonstrate the possibility of increasing $`N\tau `$ under realistic experimental conditions. Specifically, this model includes the effects of finite vapor pressure recovery times and finite loading times. Without the use of the WLS, and with long loading times, $`N\tau `$ in a magnetic trap can obtain a maximum optimal value of

$$(N\tau )_{\text{opt}}N_{\text{lim}}\tau _\text{b}/4\alpha $$
(7)
with $`\tau _{\text{Rb}}=\tau _\text{b}`$. We will compare the performance of a WLS-loaded MOT to $`(N\tau )_{\text{opt}}`$ to demonstrate the effectiveness of a WLS-loaded MOT.

For a trap loaded with the WLS, calculating $`N\tau `$ is more complicated. We divide the experimental cycle into three time periods. During the first period, the MOT is loaded, and the WLS remains on for the duration of this period. We call this period the *MOT loading phase*, which has a duration of time $`t_1`$. The cycle then enters the *MOT holding phase*, which has a duration of time $`t_2`$. In the holding phase, the WLS remains off, allowing the vapor pressure to recover while continuing to hold a large fraction of the trapped atoms in the MOT. In the third period of the cycle, the MOT beams are also turned off and the atoms are trapped in a magnetic trap. This period begins at time $`t_\text{T}=t_1+t_2`$.

Variables for the number of atoms in the trap can be defined at the boundaries of the time periods. At the beginning of the loading phase, $`N=0`$. By the end of the loading phase at time $`t_1`$, $`N_1`$ atoms are in the MOT. The cycle then enters the holding phase, during which some atoms are lost from the trap due to collisions with other trapped atoms at a rate that is faster than the decreasing loading rate into the trap. We define $`N_{\text{WLS}}`$ to be the number of trapped atoms remaining at the end of this period. The “WLS” subscript emphasizes that this number was obtained using the WLS. The cycle then enters the magnetic trap phase, and $`N_{\text{WLS}}`$ atoms are loaded into the magnetic trap. Because of the continually decreasing vapor pressure (from having used the WLS and then turning it off), the number of atoms in the magnetic trap decays faster than exponentially. Since we desire to maximize the number-lifetime product for the magnetic trap, we define an *effective lifetime* $`\tau _{\text{WLS}}`$ as the time at which the number of atoms in the magnetic trap has reached $`(1/e)N_{\text{WLS}}`$. The entire cycle as described is illustrated in Fig. 3.

Our intent in this analysis is to compare $`N_{\text{WLS}}\tau _{\text{WLS}}`$ with both $`N\tau `$ for unmodulated Rb pressures at varying loading times and with $`(N\tau )_{\text{opt}}`$, as defined in Eq. 7. First, we calculate the number of atoms $`N_1`$ in the MOT at $`t_1`$. At the beginning of the MOT loading phase, the WLS is turned on, and the loading time constant associated with the Rb partial pressure quickly drops to a value of $`\tau _1`$. We thus obtain

$$N_1=\frac{N_{\text{lim}}}{1+\tau _1/\tau _\text{b}}\left[1\mathrm{exp}\left(t_1\left[1/\tau _1+1/\tau _\text{b}\right]\right)\right]$$
(8)
from the use of Eqs. 6 and 4.

The trapped atoms then enter the holding phase. The WLS is turned off, and the number of atoms in the MOT is governed by the rate equation $`dN/dt=R(t)N/\tau _2(t)N/\tau _\text{b}`$, where $`\tau _2`$ is the loading time constant associated with the decaying Rb vapor pressure. The time dependence of $`R`$ and $`\tau _2`$ is made explicit, since these values depend upon the decreasing Rb vapor pressure. The loading rate $`R(t)`$ and the loss rate $`1/\tau _2(t)`$ are assumed to decay exponentially (with a time constant of the vapor pressure recovery time) to the steady-state values $`R(t)0`$ and $`1/\tau _20`$ (negligible Rb vapor pressures) as the vapor pressure recovers. The rate equation is numerically integrated to determine the number of atoms, $`N_{\text{WLS}}`$, left in the MOT at time $`t_\text{T}`$, the point at which the MOT is turned off.

We finally must determine the effective lifetime $`\tau _{\text{WLS}}`$ of the magnetic trapping phase of the cycle by numerically solving the rate equation $`dN/dt=\alpha (N/\tau _3(t)N/\tau _\text{b})`$. Here, $`1/\tau _3(t)`$ has an initial value of $`1/\tau _2(t_2)`$ and decays exponentially in time to 0 as the vapor pressure continues to recover. Finally, we can write the number-lifetime product of the WLS-loaded MOT, designated by $`(N\tau )_{\text{WLS}}`$, as $`(N\tau )_{\text{WLS}}=N_{\text{WLS}}\tau _{\text{WLS}}`$.

## V Results

We numerically investigated the performance of the MOT loaded with the WLS by comparing $`(N\tau )_{\text{WLS}}`$ with $`N\tau `$ for unmodulated pressures (Fig. 4). Figure 4(a) shows the number-lifetime product due to trapping atoms in a MOT for a time $`t_\text{T}=\tau _\text{b}`$ as a function of the fraction of the loading cycle that the WLS is used. We assume that $`N_{\text{max}}=N_{\text{lim}}/2`$ for unmodulated partial pressures, and an arbitrarily chosen value of $`N_{\text{max}}=N_{\text{lim}}/(1+\tau _1/\tau _\text{b})=(0.9)N_{\text{lim}}`$ (see Eq. 8), or equivalently $`\tau _1=(0.1)\tau _\text{b}`$, for the modulated partial pressures. Here, the chosen value of $`N_{\text{max}}`$ can not be set to $`N_{\text{lim}}`$ due to limitations in the numerical calculations. Fig. 4(b) shows the same conditions as Fig. 4(a), but here we have plotted the ratio of $`(N\tau )_{\text{WLS}}`$ to $`N\tau `$ with unmodulated pressures after a total MOT trapping time of $`t_\text{T}`$.

As suggested by Fig. 4, the optimum time to leave on the WLS is determined by the maximum point on a given curve. In the calculations, the gain in $`N\tau `$ after using the WLS is less than the maximum possible value of 4 due to the need to allow the vapor pressure to recover before loading the atoms into a magnetic trap. The highest values that can be achieved for $`N\tau `$ with and without the WLS are plotted against total loading time $`t_\text{T}`$ in Fig. 5 for the same conditions as in Fig. 4. The gain in using the WLS is again less than the ideal maximum factor of 4 for long loading times. However, for short loading times, $`N\tau `$ for unmodulated pressures is much lower than $`(N\tau )_{\text{opt}}`$ as shown by the gray curve in Fig. 5(a). Modulated vapor pressures can give substantial benefits in this regime, as shown by the larger $`N\tau `$ ratios in Fig. 5(b).

As a concrete example of reading the plots given here, we assume that we have a system that has a vapor pressure recovery time of $`0.035\tau _\text{b}`$. Thus we are interested in the uppermost curves in Figs. 4(a-d) and 5(a,b). If we load the vapor cell MOT without modulating the Rb partial pressure, we can achieve a value of $`N\tau 1`$ (in units of $`(N\tau )_{\text{opt}}N_{\text{lim}}\tau _\text{b}/4\alpha `$) after loading the trap for a total time of $`t_\text{T}=2\tau _\text{b}`$, as shown in the lower (gray) curve of Fig. 5(a). However, if we modulate the Rb pressure with the WLS, we can triple the value of $`N\tau `$ for the same total MOT trapping time. To determine the proper time to remove the WLS light, Fig. 4(c) indicates beginning the MOT holding phase $`0.12\tau _\text{b}`$ before loading the magnetic trap (thus $`t_\text{1}=1.88\tau _\text{b}`$) for optimum trap loading.

Alternatively, we can shorten the loading time $`t_\text{T}`$ to $`0.5\tau _\text{b}`$, as demonstrated in Fig. 5(a), and maintain the same gain in $`N\tau `$. In doing so, we would not only gain a factor of 3 in $`N\tau `$, but we would also increase the repetition rate of the experiment by as much as a factor of 4. If instead we load the experiment for a fixed amount of time in either case (with or without the WLS), we should look at Fig. 5(b) to compare the $`N\tau `$ products. For a total MOT trapping time of $`t_\text{T}=0.4\tau _\text{b}`$, we would achieve over a five-fold gain in $`N\tau `$ by modulating the vapor pressure with the WLS.

Experimentally, we were able to obtain an $`N\tau `$ product of $`(N\tau )_{\text{WLS}}=5.0\times 10^9`$ atoms$``$s using the WLS technique, achieved with a MOT loading phase of duration $`t_1=100`$ s and a MOT holding phase of duration $`t_2`$ = 50 s. This value of $`(N\tau )_{\text{WLS}}`$ is 2.2 times larger than $`(N\tau )_{\text{opt}}=2.3\times 10^9`$ atoms$``$s, $`95\%`$ of the maximum value of $`(N\tau )_{\text{opt}}`$ for optimized Rb partial pressure, reached with a loading time of $`3\tau _\text{b}=466`$ s. Without the WLS, the Rb partial pressure was optimized when $`\tau _{\text{MOT}}=\tau _{\text{Rb}}/2=\tau _\text{b}/2`$ and $`N=N_{\text{lim}}/2`$. This number of $`N\tau `$ is inferred from measurements of $`N_{\text{lim}}`$ and $`\tau _\text{b}`$ mentioned previously.

Note that in addition to the gain in $`N\tau `$, the time to reach the above value of $`(N\tau )_{\text{WLS}}`$ is 3.1 times faster than the time to reach the above value of $`(N\tau )_{\text{opt}}`$ (without the WLS), tripling the repetition rate of experiments. The WLS experimental technique would be even more beneficial by shortening the recovery time of the vapor pressure. This might be accomplished by keeping a larger fraction of the inner surface of the vapor cell at cryogenic temperatures or through optimization of the surface adsorption chemistry.

## VI Conclusions and summary

The use of LIAD to enhance loading of vapor cell MOTs may be applicable to other atomic species. Lithium vapor cells, for instance, are difficult to work with due to the high temperatures needed to create a substantial Li vapor. Yet if LIAD were to work well with Li adsorbed on stainless steel, or between co-adsorbed Li atoms on a surface, a Li vapor cell MOT would be practical. Although the effect has not yet been quantitatively explored as it has been for Rb, we observed a LIAD induced increase in the loading rate into a Cs MOT in a Cs vapor cell with aluminum walls. In general, when first using the LIAD technique, the WLS intensity should be raised incrementally to monitor the loading time constant. When the WLS loading time constant becomes too short ($`<10`$ s) the vapor pressure can potentially become high enough that atoms may re-adsorb onto cold chamber windows and may possibly form small clusters of atoms.

In summary, we have demonstrated that the technique of non-thermal light induced atom desorption can be used to effectively increase the number of atoms that can be loaded into a vapor cell MOT. This technique benefits atom trapping experiments where large numbers of atoms and long trap lifetimes are crucial.


format: text

text:
# Universal saturation of electron dephasing in three-dimensional disordered metals

## Abstract

We have systematically investigated the low-temperature electron dephasing times $`\tau _\varphi `$ in more than 40 three-dimensional polycrystalline impure metals with distinct material characteristics. In all cases, a saturation of the dephasing time is observed below about a (few) degree(s) Kelvin, depending on samples. The value of the saturated dephasing time $`\tau _0`$ \[$`\tau _\varphi (T0\mathrm{K})`$\] falls basically in the range 0.005 to 0.5 ns for all samples. Particularly, we find that $`\tau _0`$ scales with the electron diffusion constant $`D`$ as $`\tau _0D^\alpha `$, with $`\alpha `$ close to or slightly larger than 1, for over two decades of $`D`$ from about 0.1 to 10 cm2/s. Our observation suggests that the saturation behavior of $`\tau _\varphi `$ is universal and intrinsic in three-dimensional polycrystalline impure metals. A complete theoretical explanation is not yet available.

The issue concerning the saturation of conduction electron dephasing time $`\tau _\varphi (T)`$ in disordered conductors and mesoscopic systems has recently attracted renewed theoretical  and experimental  attention. In particular, it is of great interest in whether or not the saturation of $`\tau _\varphi `$ might be universal in all material (e.g., polycrystalline metal, amorphous metal, and MBE-grown semiconductor) systems and in all (zero, one, two, and three) dimensions. It is also of importance to ask whether or not there might be a common, intrinsic source that causes the (almost) universally observed saturation. The interest in this issue of the existence or not of a finite conduction electron dephasing time at absolute zero temperature is closely connected with the concern about the validity of the Fermi-liquid behavior, the possibility of the occurrence of a quantum phase transition, and also the feasibility of quantum computing . There are also works suggesting that this might explain the long-standing persistent current problem in metals . In practice, the value of $`\tau _\varphi `$ can be reliably determined using quantum-interference studies  such as weak-localization effects and universal conductance fluctuations. It can also be determined using studies of the shape change of the electron distribution function upon applying a bias voltage across the conductor .

In this work, we have systematically investigated the electron dephasing times $`\tau _\varphi `$ in numerous three-dimensional polycrystalline impure metals. Our samples were made of various materials using various fabrication techniques at different laboratories. Also, the samples were measured at different laboratories at very different times. Our samples considered in this work include dc sputtered thick Au50Pd50 films, dc and/or RF sputtered thick Ag40Pd60 films, dc sputtered thick Sb films, thermal-evaporation deposited thick AuxAl films ($`1.8x2.2`$), thermal-evaporation deposited thick Sc85Ag15 films, and arc-melted V100-xAlx alloys ($`20x24`$). The thick film samples are typically on the order 4000 $`\mathrm{\AA }`$ $`\times `$ 0.3 mm $`\times `$ 17 mm, while the arc-melted samples are typically on the order 0.1 $`\times `$ 0.1 $`\times `$ 10 mm3. Since more than 40 samples are studied, it would not be practical to list the material parameters for all the samples. However, we notice that all of the samples studied are made from very high-purity starting materials obtained from reputable suppliers such as Alfa Aesar, Cerac, and Goodfellow. Since the major aim of this work is to study experimentally whether there might exist a universal behavior of the saturated (or, the zero-temperature) dephasing time $`\tau _0`$ \[$`\tau _\varphi (T0\mathrm{K})`$\], the use of many kinds of samples with very distinct characteristics should therefore serve well this purpose. Any behavior of $`\tau _0`$ common for all these materials, if observed, must be a manifestation of the very general nature of the zero-temperature dephasing time. In addition, we notice two more practices regarding our experimental method. (a) Our measurements of $`\tau _\varphi `$ had been performed over a long time period of over three years (1997-2000). (b) Our measurements of $`\tau _\varphi `$ were carried out at two different laboratories located at two different sites 80 km apart. Nevertheless, the same 3He fridge and the same electronic measuring systems (except with different grounding) were used for all measurements. Surprisingly, regardless of all the above-mentioned very different preparation and measurement conditions, it is noteworthy that our $`\tau _0`$ measured in all samples varies with the electron diffusion constant $`D`$ with a simple power law as $`\tau _0D^\alpha `$, with $`\alpha `$ close to or slightly larger than 1 (see below). Such universal behavior could not be just “accidental” and certainly deserves serious theoretical and experimental attention. We emphasize that we are concerned with three-dimensional polycrystalline metals. It is our opinion that the saturation behavior (e.g., the functional form of $`\tau _0`$ on disorder) is universal for a given dimensionality and a given kind of sample (e.g., polycrystalline or well-textured semiconductor) structure, while it might not be universal over different dimensionalities and different sample structures.

We have measured the low-field magnetoresistances at liquid-helium temperatures and compared with three-dimensional weak-localization predictions to extract the values of $`\tau _\varphi (T)`$ for our samples. The details of the data analysis procedure had been discussed previously . Here we merely stress that, for every sample studied in this work, the three-dimensional weak-localization predictions can well describe our experimental measurements. Thus, $`\tau _\varphi `$ can be very reliably extracted. Figure 3 shows the total electron dephasing rate $`\tau _\varphi ^1`$ as a function of temperature for three representative thick Ag40Pd60 films with resistivities as indicated in the caption to Fig. 3. The symbols are experimental data and the solid curves (emphasizing particularly the low measurement temperatures) are the least-squares fits to the form

$$\tau _\varphi ^1(T)=\tau _0^1+\tau _\mathrm{i}^1(T),$$
(1)
where we assume that the total dephasing time comprises two terms, with $`\tau _0`$ being a temperature-independent (i.e., the saturated) dephasing time which dominates at the lowest measurement temperatures, and $`\tau _\mathrm{i}`$ being the relevant inelastic scattering time which is usually dominant at a (few) degree(s) Kelvin or so and higher, strongly depending on samples. One sees that, for every sample, $`\tau _\varphi ^1`$ first decreases with decreasing temperature in the high-temperature part of our measurements. This decrease in the total dephasing rate can readily be ascribed to the weakening of the inelastic scattering process with decreasing temperature. That is, $`\tau _\varphi ^1`$ is dominated by $`\tau _\mathrm{i}^1`$ and it decreases with weakening $`\tau _\mathrm{i}^1`$ above 3 K or so for these particular samples. At lower temperatures, the inelastic scattering becomes weaker and weaker so that $`\tau _\varphi ^1`$ becomes more and more dominated by the temperature-independent $`\tau _0^1`$. Inspection of Fig. 3 reveals that Eq. (1) can well describe the experimental data. In fact, in some cases (e.g., Fig. 3) when the saturation of $`\tau _\varphi `$ already sets in at a few degrees Kelvin high, $`\tau _0`$ can already be determined using essentially the measured value at the lowest temperatures. This “advantage” of (relatively) high-temperature saturation largely reduces the uncertainties in the least-squares fits, making the inferred values of the adjusting parameters $`\tau _0`$ and $`\tau _\mathrm{i}`$ very reliable.

Before the recent renewed interest in $`\tau _0`$, a saturation of the electron dephasing time as that shown in Fig. 3 has been observed in a good number of experiments, including reduced-dimensional metals  and semiconductors . One of the most adopted early explanations proposed for the observed saturation invokes hot-electron effects. However, this kind of explanation can be ruled out for our measurements in Fig. 3. In the hot-electron picture, it is argued that if the Joule heating caused by measurement currents were sufficiently large, the conduction electrons would be removed out of thermal equilibrium with the phonon bath. Then, the electron temperature would be higher than the phonon (lattice) temperature, resulting in a “roll-over” of $`\tau _\varphi `$ at the lowest measurement temperatures. Since we are concerned with three-dimensional systems in this work, Joule heating is negligible in our case. To confirm this assertion experimentally, we plot in Fig. 3 the variations of resistivity with temperature for those representative thick Ag40Pd60 films just considered in Fig. 3. Comparison examination of Figs. 3 and 3 indicates that the resistivities of the samples keep increasing with reducing temperature all the way down to our lowest measurement temperatures (400 mK), though a saturation in $`\tau _\varphi `$ is already seen at a much higher temperature ($``$ 3 K in these particular samples). Indeed, the resistivity rises with decreasing temperature given in Fig. 3 can essentially be described by disorder enhanced electron-electron interaction effects , suggesting that the temperature of the electron gas is the same to the temperature of the lattice. Therefore, no appreciable Joule heating was produced in our measurements. Consequently, the saturation of $`\tau _\varphi `$ shown in Fig. 3 can not be interpreted in terms of hot-electron effects. This is true for all samples studied.

The central result of this work is plotted in Fig. 3. Figure 3 shows a plot of our measured saturated dephasing time $`\tau _0`$ as a function of the electron diffusion constant $`D`$ for more than 40 samples made from very different materials and prepared with very different fabrication techniques. Different symbols label different material systems as indicated in the caption to Fig. 3. In all cases, the samples are three-dimensional with regard to the electron phase coherence length $`\sqrt{D\tau _\varphi }`$. Inspection of Fig. 3 reveals that the saturated dephasing time $`\tau _0`$ for all samples falls essentially on the same region in the coordinates. In particular, it reveals that all that matter in determining the value of $`\tau _0`$ is the diffusion constant $`D`$, regardless of the distinct material characteristics (e.g., electronic structures) of the various samples. (It is understood that the evaluation of the diffusion constant $`D`$ for some metal alloys could probably be subject to an uncertainty of, e.g., as large as a factor $``$ 2, due to the complex electronic structures. However, our experimental availability of a very wide range of $`D`$ ensures that even with a factor of 2 in the uncertainty of $`D`$ would not change the main conclusion of this work .) Close inspection of Fig. 3 indicates that $`\tau _0`$ varies essentially with the diffusion constant $`D`$ as $`\tau _0D^\alpha `$, with $`\alpha `$ close to or slightly larger than 1, for about two decades of $`D`$. As $`D`$ increases from about 0.1 to 10 cm2/s (corresponding to the product $`k_Fl`$ ranging from of order unity to of order 10, where $`k_F`$ is the Fermi wave number, and $`l`$ is the electron elastic mean free path), $`\tau _0`$ decreases correspondingly from about 0.5 to 0.005 ns. In this figure, the straight solid line is drawn proportional to $`D^1`$ and is a guide to the eye. This observed $`\tau _0D^\alpha `$, with $`\alpha `$ close to or slightly larger than 1, is totally unexpected. On the contrary, it is often conjectured that $`\tau _0`$ should increase with reducing disorder, at least in one and two dimensions . Until now, it is not known exactly how differently the saturated dephasing time should behave in different dimensionalities and in different sample structures. Surprisingly, we notice that our observation of an essentially inverse linear dependence of $`\tau _0`$ on $`D`$ actually implies an essentially “constant” saturated dephasing length of $`\sqrt{D\tau _0}`$ 1000 $`\mathrm{\AA }`$ in all samples. For comparison, the experimental values of $`\tau _0`$ previously reported by Lin and Giordano  in thin, polycrystalline AuPd films are also indicated in Fig. 3 by the two vertical bars: the vertical bar at $`D`$ 2.2 cm2/s stands for the range of $`\tau _0`$ observed in their dc sputtered films, while the vertical bar at $`D`$ 23 cm2/s stands for that in their thermally evaporated films. The observation of Fig. 3 is for the first time such a systematic study has been done on so many samples with a wide range of diffusion constant.

One of the most popular (early) explanations for the saturation of $`\tau _\varphi `$ invokes spin-spin scattering due to the presence of a very minor amount of magnetic impurities in the sample. This explanation was refuted in many experiments . In the present work, this explanation can also be ruled out simply because our samples were made from very different high-purity sources obtained from various reputable suppliers, and measured at very different times. It is hard to conceive that magnetic spin-spin scattering due to “accidental” contamination could have caused the “systematic” variation of $`\tau _0`$ with $`D`$ as observed in Fig. 3. Moreover, since we are concerned with bulk samples, any spin-flip scattering that might result from surface effects, such as interfaces, substrates, and paramagnetic surface oxidation, can largely be minimized (while surface effects could be significantly more important in reduced-dimensional systems). Indeed, as mentioned, we do not see a logarithmic Kondo type of dependence in resistivity (Fig. 3) in our temperature range, strongly ruling out magnetic impurities playing any role in our samples. It is also worthwhile pointing out in passing that, although it has been widely taken for granted for over years that magnetic impurities will cause dephasing and produce a temperature independent $`\tau _\varphi `$, Webb and coworkers  have recently ruled out magnetic impurities as the cause of the observed saturation.

The origin for the observed saturated dephasing time given in Fig. 3 is not exactly clear. To our knowledge, Altshuler, Aronov, and Khmelnitsky  have considered the dephasing of electron wave amplitudes by high-frequency electromagnetic noises. They found that an electromagnetic noise can be, on the one hand, already large enough to cause dephasing while, on the other hand, still too small to cause significant Joule heating of the conduction electrons. Unfortunately, in the most effective frequency range in causing dephasing, their theory predicts a form $`\tau _0D^{1/3}`$ which is in disagreement with our experimental observation. Besides, recent experiments  has explicitly demonstrated that direct dephasing due to radiation could not be the cause of the observed saturation. Therefore, we have to turn to other explanations, such as, among others, zero-point fluctuations of the electrons , dynamical two level systems , and coherent charge transfer between crystallites . Unfortunately, in terms of the zero-point fluctuations theory, we observe, in almost all samples at 500 mK, a $`\tau _0`$ already exceeds 3 to 5 orders in magnitude above the upper bound set by the theory. \[The saturation time in three dimensions is given by the zero-point fluctuations theory as $`\tau _{\mathrm{GZ}}^{(3\mathrm{D})}=(22\pi ^2\mathrm{}/e^2)(D^2/\rho _0v_F^3)`$, where $`\rho _0`$ is the residual impurity, and $`v_F`$ is the Fermi velocity. See, e.g., Refs. .\] Even worse, the zero-point fluctuations theory predicts increasing $`\tau _0`$ with increasing $`D`$, opposite to our observation in polycrystalline impure metals. Although an increase of $`\tau _0`$ with increasing mobility has been observed in, e.g., semiconductor quantum wires  and quantum dots , it is conceived  that the saturation behavior in polycrystalline impure metals could be different from that in well-textured semiconductor quantum wires and dots. Indeed, it is our opinion that the saturation behavior should be universal among polycrystalline impure metals, while different functional forms of $`\tau _0`$ on disorder might be expected for different sample structures. If this were the case, then the current concept for zero-temperature dephasing would need major redirection.

In conclusion, we have observed a saturated dephasing time $`\tau _0`$ on the order of 0.005 to 0.5 ns in more than 40 three-dimensional polycrystalline impure metals with distinct characteristics. Taken many metals and alloys together, our result indicates that $`\tau _0`$ varies essentially with the electron diffusion constant $`D`$ as $`\tau _0D^\alpha `$, with $`\alpha `$ close to or slightly larger than 1, for about two decades of $`D`$. It is for the first time such a systematic study has been done on so many samples with a wide range of diffusion constant. A complete theoretical explanation for this observation is not yet available. On the other hand, our result suggests that the origin for the observed saturation at low temperatures might not be universal over all dimensionalities and all kinds of sample structures. If this were the case, then the current concept for zero-temperature dephasing would need major redirection.

The authors are grateful to T. J. Li, S. B. Chiou, S. F. Chang, Y. L. Zhong, and A. K. Meikap for help in the experiment. This work was supported by Taiwan National Science Council through Grant No. NSC 89-2112-M-009-033.


text:
# 1 Introduction

## 1 Introduction

In order to have a precise function in the biological “machinery”, it is important for proteins to have an unique conformation at physiological temperatures. This is termed the native state. Anfinsen  proved in his famous experiment with ribonuclease the important fact that the folding of the polypeptide chain is thermodynamically determined.

One simplified view of protein folding is that the protein is supposed to follow a specific folding pathway of conformations in a descending landscape of Gibbs free energy . This is a picture of a folding protein that is forced to follow a specific “path” of successive conformational steps of increasing structural order. We will use this pathway assumption in this paper.

A protein in physiological environments (pH, ionic strength etc.), and temperatures is packed in a very compact way. It is then termed folded. An increase of the temperature will eventually denaturate the protein, i.e. it unfolds. Other ways to unfold the protein are for instance to change the pressure, denaturant concentration or the pH. The fact that proteins also unfold at low temperatures, termed as cold unfolding, makes the system very unusual. A major difficulty in experiments of cold unfolding is that the temperature is around and below the freezing point of water. In a frozen aqueous solution, one cannot observe any conformational transitions .

A general feature of small globular proteins is that they thermodynamically seems to unfold in an “all-or-none” manner. This means that they unfold cooperatively without noticeable intermediates , with a deviation from a two-state system not exceeding 5 %. The deviation from a single macroscopic system can be explained by presence of unstable intermediates . It is worth noting that all these experiments have been done only for the warm unfolding. The occurrence of intermediate states in larger proteins  is not a contradiction to the two-state behavior in the experiments in Refs. , because the latter only considers small globular proteins.

The van’t Hoff enthalpy relation (for heat of reaction)

$$\mathrm{\Delta }H=\alpha k_BT_c^2\frac{\mathrm{\Delta }C}{Q},$$
(1)
is a powerful way to quantify the sharpness of a smoothed out first order phase transition. As shown in Fig. 1, $`T_c`$ is the transition temperature (at the middle of the peak), $`Q`$, which is the same as $`\mathrm{\Delta }H`$ (no pressure), is the released energy (latent heat), and $`\mathrm{\Delta }C`$ is the peak height of the transition. $`\alpha `$ is a dimensionless proportionality factor. For a given $`\mathrm{\Delta }H`$ and $`Q`$, then the value of $`\alpha `$ is inversely proportional to $`\mathrm{\Delta }C`$. In this respect a smaller $`\alpha `$ corresponds to a sharper transition.

In this article we will explain the physical basis of a protein model, that reformulates the water interactions proposed in earlier models by Hansen et al. and Bakk et al.. We will compare thermodynamical quantities, as the heat capacity, to experiments. The protein is also investigated in a temperature region below accessible experimental data.

## 2 Modeling the protein

### 2.1 The polypeptide chain

The polypeptide chain is modeled as in earlier articles by Hansen et al. and Bakk et al., where the protein is supposed to follow a pathway as described in Section 1. The protein is equipped with $`N`$ contact points, which we here call nodes. A node in contact means that the protein has a “correct” conformation on the folding pathway, and we say that the node is folded. Due to the fact that a protein is a complex 3-dimensional system, a folded node likely has non-local contacts with respect to the amino acid sequence in the polypeptide chain.

Each node is assigned only two energies $`ϵ_0`$ or 0 . By using the binary variables $`\varphi _i\{0,1\}`$, the energy associated to each individual node is written
$`E_i=ϵ_0\varphi _1\varphi _2\mathrm{}\varphi _{i1}\varphi _i.`$ The value $`\varphi _i=0`$ means an unfolded node, while $`\varphi _i=1`$ is equivalent to a folded node. The product term meets the assumption about a folding pathway, because $`E_i=ϵ_0`$ only if all nodes before $`i`$ is folded, in addition to a folded node $`i`$ itself. For a system of $`N`$ nodes the Hamiltonian of the polypeptide chain (vacuum energy of the protein) is

$$_c=ϵ_0(\varphi _1+\varphi _1\varphi _2+\mathrm{}+\varphi _1\varphi _2\mathrm{}\varphi _N).$$
(2)
Folding of node $`i`$ is only proper in the case of an unique conformation before this folding step. According to Eq. 2, every attempt to fold $`i`$ while $`j`$ (describing one or several nodes $`<i`$) is unfolded will not gain energy in the system, because then the protein is supposed to be in an energetical unfavorable conformation.

The unfolded protein has more degrees of freedom relative to the folded, because the unfolded polypeptide backbone will have rotational freedom. We incorporate this by assigning each unfolded node $`f`$ degrees of freedom. The parameter $`f`$ is interpreted as the relative increase in the degrees of freedom for an unfolded node compared to a folded node.

### 2.2 Water interactions

Interactions between water and protein surface is important. Proteins are during the evolution “designed” to interact with water, simply because they are exposed to water in vivo. Makhatadze and Privalov  states that in sum hydration effects destabilize the native state, and decreasing temperature implies increasing destabilizing action. The water that access the hydrophobic protein interior during unfolding is supposed to obtain an “ice-like” structure around the apolar surfaces. Hence, this structured water will both decrease the entropy and the energy compared to “free” water , and thus impacts the thermodynamics of the system.

Hansen et al. proposed a simple model where each water molecule interacts with a node, and not with other water molecules, by a “ladder” of $`g`$ equidistant energies accessible

$$\omega _i=\{\begin{array}{cc}\epsilon _w+(g1)\delta \hfill & \\ \mathrm{}\hfill & \\ \epsilon _w+2\delta \hfill & \\ \epsilon _w+\delta \hfill & \\ \epsilon _w,\hfill & \end{array}$$
(3)
which we will also apply in the model considered in this text. The interpretation of $`\omega _i`$ is the energy difference between a “structured” water molecule, associated to the unfolded parts of the protein, and a “free” water molecule in the bulk.

The observable states in a small globular protein is either the unfolded ($`\varphi _1=0`$), with water bounded to the surface that uncovers during unfolding of the protein, or the folded state ($`\varphi _1\mathrm{}\varphi _N=1`$) with no water in the protein interior. No intermediate states are detected for small globular proteins , hence one cannot know for sure how the water enters the protein interior during the unfolding. Hansen et al. and Bakk et al. have earlier only considered that the amount of water interactions increase proportional to the number of unfolded nodes, and with that the contact energy of the chain. In this paper we study, as a more speculative assumption, the case when a macroscopic contribution of water enters the protein surface when the last node is unfolded.

We note that Eq. 3 is the quantized energy levels of a magnetic dipole in an external field. In the continuum limit where $`g\mathrm{}`$ (with $`g\delta `$ finite), a classical magnetic dipole in an external field is obtained, and this again is analogous to an electrical dipole in an external electrical field. The dipolar water molecules are exposed to an electrical field from the permanent, and induced charges on the protein surface, thus Eq. 3 is a representation of that.

By using the same notation as in Eq. 2, we propose the Hamiltonian that corresponds to the water-protein interactions

$$_w=(1\varphi _1\varphi _2\mathrm{}\varphi _N)(\omega _1+\omega _2+\mathrm{}+\omega _M),$$
(4)
where $`M`$ is the number of water molecules. The folded protein is a highly ordered and dense packed structure where no water can access the interior. Due to Eq. 2, unfolding of the last node $`(\varphi _N=0)`$ implies a less dense pacing of the protein, and the cavities are now supposed to be big enough to let water access the interior of the protein. The next step, unfolding of node $`N1`$, implies likely an even lesser dense packing, and allows more water in the protein interior. We assume in this text that the water entering upon unfolding of node $`N1`$, will not interact with the protein surface, because it is regarded as a second layer of “free”water. Cohn and Edsall  states that roughly a monolayer of water is bounded to the protein, implying that the protein is only interacting with the first monolayer, thus the second, and third etc. water layers, successively entering the protein during unfolding, are regarded as “free water”. Hence, according to the latter possible (but somewhat speculative) explanation of how the water access the apolar interior of the protein, unfolding of nodes $`i<N`$ does not contribute energetical to the water Hamiltonian ($`_w`$) and thus not to the thermodynamics.

### 2.3 The statistical framework

The system Hamiltonian ($``$) describing both chain specific energy ($`_c`$) and water interactions ($`_w`$) is

$$\begin{array}{cc}\hfill =_c+_w=& ϵ_0(\varphi _1+\varphi _1\varphi _2\mathrm{}+\varphi _1\varphi _2\mathrm{}\varphi _N)\hfill \\ & +(1\varphi _1\varphi _2\mathrm{}\varphi _N)(\omega _1+\omega _2+\mathrm{}+\omega _M)\hfill \end{array}$$
(5)
Let now $`Z_i`$ be term number $`i`$ in the partition function which corresponds to folding of all nodes $`i`$ (pathway assumption), thus

$$Z_i=f^{Ni}e^{iϵ_0\beta }\left(e^{\epsilon _w\beta }\frac{1e^{g\delta \beta }}{1e^{\delta \beta }}\right)^M(i<N).$$
(6)
$`\beta 1/T`$ is a rescaled inverse absolute temperature in which the Boltzmann constant is absorbed. $`Z_0`$ is the term where all nodes are zero, i.e. a complete unfolded protein, while $`Z_N`$ corresponds to a folded protein. The factor $`f^{Ni}`$ in Eq. 6 is the degrees of freedom in the polypeptide chain that is available in the $`Ni`$ unfolded nodes. $`e^{iϵ_0\beta }`$ is the Boltzmann factor from $`i`$ contact energies $`ϵ_0`$ in the polypeptide chain. The last term in brackets is simply the sum over all distinct levels in one water molecule raised to the power of the number of water molecules $`M`$ bounded to the unfolded parts of the protein. We assume that $`\delta \beta 1`$ (i.e.$`g\mathrm{}`$), which is equal to an infinite small level spacing in Eq. 3. A first order Taylor expansion of the denominator in Eq. 6 yields

$$Z_i=f^N\left(\frac{e^{\epsilon _w\beta }}{\delta \beta }\right)^Me^{i(ϵ_0\beta \mathrm{ln}f)}(i<N),$$
(7)
assuming $`1e^{g\delta \beta }1`$ when $`g\mathrm{}`$. The last term in the partition function ($`Z_N`$) corresponds to a complete folded protein, where there are $`g^M`$ degrees of freedom from $`M`$ unbounded water molecules and $`N`$ contact energies $`ϵ_0`$, hence

$$Z_N=g^Me^{Nϵ_0\beta }.$$
(8)
By summing up the $`Z_i`$ terms in Eqs. 7 and 8, we obtain the partition function

$$Z=\underset{i=0}{\overset{N}{}}Z_i=f^Ng^M\left[\left(\frac{ae^{\mu \beta }}{\beta }\right)^M\frac{1r^N}{1r}+r^N\right],$$
(9)
where $`re^{\beta \mathrm{ln}f}`$, $`aϵ_0/(g\delta )`$ and $`\mu \epsilon _w/ϵ_0`$. The inverse temperature is here rescaled by $`ϵ_0\beta \beta `$. The parameter $`\mu `$ measures the strength of the water interactions relative the chain contact energy, thus $`\mu `$ is interpreted as an effective chemical potential. Changing $`\mu `$ means adding denaturants, changing pH or salt concentration etc.

The order parameter $`n`$ in this system measures the degree of folding, i.e. the mean number of folded nodes divided by $`N`$

$$n=\frac{_{i=0}^NiZ_i}{Z}=\frac{1}{N}\frac{\left(\frac{ae^{\mu \beta }}{\beta }\right)^M\frac{(N1)r^{N+1}Nr^N+r}{(1r)^2}+Nr^N}{\left(\frac{ae^{\mu \beta }}{\beta }\right)^M\frac{1r^N}{1r}+r^N}.$$
(10)
$`n=0`$ corresponds to an unfolded protein, while $`n=1`$ is interpreted as a completely folded protein.

## 3 Calculations and discussion

The heat capacity is $`C=\beta ^2^2(\mathrm{ln}Z)/\beta ^2`$. Fig. 2 shows a typical plot of the heat capacity $`C(T)`$ with three peaks (numbered 1, 2 and 3 from left). These characteristic peaks corresponds to three critical transition temperatures: $`T_1`$, $`T_2`$ and $`T_3`$, measuring the temperatures at the respective peak maxima. The corresponding order parameter $`n(T)`$ in Fig. 3, calculated from Eq. 10, shows that the protein is essential folded for $`T<T_1`$ and $`T_2<T<T_3`$, while the protein is nearly unfolded in the temperature intervals $`T_1<T<T_2`$ and $`T>T_3`$. From this picture it is reasonable to state that the physiological temperature interval is between peak 2 and 3. Accordingly, with reference to this temperature region, we call peak 1 for cold folding and peak 2 and 3 respectively for cold and warm unfolding. Peak 2 and 3 are both observed in experiments  and are also seen in the model of Hansen et al. and Bakk et al.. The model considered in this paper has, in addition to the cold and warm unfolding, the peculiarity of cold folding.

Experiments on cold unfolding are very difficult because most proteins unfolds below the freezing point of water. Chen and Chellman  and Privalov et al. have all done experiments where the cold unfolding temperature is elevated by denaturants, but denaturants make the interpretation of the data very difficult. However, Privalov  did experiments in super cooled water, which is easier to interpret. Unfortunately he was not able to detect the sharpness of the cold unfolding, and not at all the heat capacity below the cold unfolding. This means that our model may predict a cold folding transition at a temperature below the cold unfolding transition.

For temperatures below the cold folding ($`T0`$) analysis of the model gives $`n0.99`$, i.e. only the last node is unfolded, and corresponds to the global energy minimum. From Fig. 2 it is seen that $`T_11.45`$. In Eq. 9 the critical
$`re^{1/T\mathrm{ln}f}=1`$, caused by the contact energies of the polypeptide chain, implies $`T_1=1/\mathrm{ln}21.44`$. Hence, this is nothing but a transition initiated of the chain entropy. An increase of the temperature from $`T_1`$ takes the protein through a nearly unfolded state, whereupon the protein folds again at $`T_21.8`$. The temperature is now so high that the energy of water (Eq. 3) is small (thermal exited) compared to the chain contact energy $`ϵ_0`$, thus the protein prefers to fold again. Further increase of the temperature causes warm unfolding at $`T_33.0`$, because then the entropy of the chain again dominates the Gibbs free energy. It is interesting to note that the entropy of the chain causes two transitions, the cold folding and warm unfolding.

We now turn our interest to the sharpness of the transitions, i.e. the value of the parameter $`\alpha `$ in the van’t Hoff enthalpy relation (Eq. 1). For $`M=200`$ is $`\alpha 4`$ both for the cold and warm unfolding. This means that the protein is thermodynamically regarded as a two-state system that folds in an “all-or-none” manner. Privalov  has measured $`\alpha =4.0`$ for the warm unfolding. As far as we know there are no experiments on the sharpness of the cold unfolding, but Privalov  indicates a sharpness for the cold unfolding as well, thus according to our model. The cold folding transition has $`\alpha 12`$. This value is typical for a transition where one has small energy differences between the folded/unfolded states and the intermediate states. Remember that the “folded” state for $`T<T_1`$ is actually the first node unfolded, thus the unfolding will essential depend on the polypeptide chain with the Hamiltonian in Eq. 2, which can be shown corresponds to $`\alpha =12`$.

Finally we note the consequence of a decreased $`\mu `$ is an increasing separation between the cold and warm unfolding as seen in Fig. 4. This makes sense because a smaller $`\mu `$ is equivalent to a relatively smaller $`\epsilon _w`$ compared to $`ϵ_0`$ (see Eq. 9), i.e. it is less favorable for the protein to be bounded to water. The consequence is that the protein prefers to be folded in a larger temperature interval, in where the water is expelled to the bulk. However, the transition temperature $`T_1`$ is not changed because this transition is given by the value $`T_1=1/\mathrm{ln}f`$. It is also seen that a smaller $`\mu `$ is qualitatively equivalent to a smaller $`a`$.

An increase in $`M`$ is the same as a decrease in $`N`$, because then the water becomes more important relative to the chain, and will again allow a broader separation between $`T_2`$ and $`T_3`$. The broader separation is also seen for a decreasing $`f`$, because this is equivalent to a larger $`M`$.

Further increase of $`\mu `$ will eventually merge peak 1 and 2. It is interesting to note that $`\alpha =4`$ for the merged peaks, because then the transition is energetically dominated of the $`M`$ water molecules which caused the transition at $`T_2`$ in Fig. 2.

In sum the qualitative change from Fig. 3 to Fig. 4, by a decreasing $`\mu `$, is also obtained by an decrease of $`a`$, $`f`$ or $`N`$ or an increase of $`M`$.

## 4 Conclusion

We have in this paper studied a protein model with water interactions. The model is based on earlier models by Hansen et al. and Bakk et al.. In contrast to these similar models, where the water amount was supposed to increase linearly to the degree of unfolding of the polypeptide chain, we have, with a more speculative assumption, studied the situation where a macroscopic amount of water access the protein interior during unfolding of the last node is the only contribution to the water-protein Hamiltonian.

With reference to physiological temperatures we find that the protein exhibits cold and warm unfolding transitions, which is an experimental fact . These transitions are associated by a sharpness indicating, from at thermodynamical point of view, a two-state system, which is also experimental established . Decreasing the temperature further below the cold unfolding region the protein folds again (cold folding). This folding, caused by the chain entropy, has a less sharp transition, which corresponds to a transition where the intermediate folding energies does not differ significant from the folded/unfolded energies. In sum the model exhibits three unfolding/folding transitions.

It is interesting to note that both the cold folding and the warm unfolding is due to the polypeptide chain entropy.

Acknowledgements

A. B. thanks the Norwegian Research Council for financial support. A. H. thanks NORDITA and Niels Bohr Institute for warm hospitality and support. We thank J. S. Høye for enthusiastic and enlightening discussions.

Figure captions

Fig. 1. Schematic illustration of the heat capacity around an unfolding transition showing the parameters in the van’t Hoff enthalpy relation (Eq. 1). $`T_c`$ is the transition temperature, $`Q`$ (area of the peak) is the released energy (latent heat) and $`\mathrm{\Delta }C`$ is the peak height of the transition.

Fig. 2. Heat capacity $`C(T)`$ with the parameters $`a=0.077`$, $`\mu =3.3`$, $`f=2`$, and $`N`$=$`M`$=$`\mathrm{\hspace{0.17em}200}`$, showing three peaks. With reference to the temperature region between peak 2 and 3 (physiological temperatures) we call the transitions: 1) cold folding, 2) cold unfolding and 3) warm unfolding.

Fig. 3. The corresponding order parameter $`n(T)`$ to Fig. 2, describing the degree of folding. The chosen parameters are as in Fig. 2. $`n=0`$ corresponds to an unfolded protein, while $`n=1`$ is interpreted as a completely folded protein.

Fig. 4. Heat capacity where the effective chemical potential is $`\mu =3.2`$, slightly decreased from the value in Fig. 2 ($`\mu =3.3`$). All other parameters are chosen as in Fig. 2. The qualitative picture is a broader separation between $`T_2`$ and $`T_3`$ compared to Fig. 2.

A. Bakk, A. Hansen and K. Sneppen
A protein model exhibiting three folding transitions

A. Bakk, A. Hansen and K. Sneppen
A protein model exhibiting three folding transitions

A. Bakk, A. Hansen and K. Sneppen
A protein model exhibiting three folding transitions

A. Bakk, A. Hansen and K. Sneppen
A protein model exhibiting three folding transitions


text:
# 1 (a) 𝑃(𝑠) for the model () compared with the the semi-Poisson distribution, (b) the same as in (a) for short distances. In addition the theoretical result () is drawn as a dashed line and the GOE result as a dashed dotted curve.

Comment on ”Models of Intermediate Spectral Statistics”

T. Gorin1, M. Müller2 and P. Seba3,4

1 Centro de Ciencias Físicas, UNAM, Campus Morelos, Cuernavaca, Mexico, 62251, A.P. 48-3

2 Facultad de Ciencias, UAEM, C.P. 62210, Cuernavaca, Morelos, Mexico

3 Institute of Physics, Czech Academy of Science, Cukrovarnicka 10, Prague, Czech Republic

4 Pedagogical University, Department of Physics, Hradec Kralove, Czech Republic

PACS: 03.65, 05.45

## Abstract

In  it was proposed, that the nearest-neighbor distribution $`P(s)`$ of the spectrum of the Bohr-Mottelson model is similar to the semi-Poisson distribution. We show however, that $`P(s)`$ of this model differs considerably in many aspects from semi-Poisson. In addition we give an asymptotic formula for $`P(s)`$ as $`s0`$, which gives $`P^{}(0)=\pi \sqrt{3}/2`$ for the slope at $`s=0`$. This is different not only from the GOE case but also from the semi-Poisson prediction that leads to $`P^{}(0)=4`$

The motivation for this comment stems from the impression, that the article “models of intermediate spectral statistics” can easily be misinterpreted in two ways: One might be led to believe that (i) the semi-Poisson distribution is universal, and (ii) the universality class of “intermediate statistics” is as well defined and established as for example the Poisson ensemble or the GOE. In this comment, we will argue, that both statements are wrong.

The purpose of  is to present models which could constitute a “third” universality class of systems which show so called “intermediate statistics”, previously introduced by Shklovskii in . The Poissonian and the Gaussian ensemble (for definiteness, consider orthogonal ensembles only) are considered as the first two universality classes in this list.

As in the Poissonian and in the GOE case, where the respective members have common and unique statistical properties, one would expect the same to hold for the models with intermediate statistics. In  the authors concentrate on the distribution of nearest neighbor spacings. In the Poissonian case it is given by $`P(s)=\mathrm{exp}(s)`$, in the GOE case it is close to the well known Wigner surmise $`P(s)(\pi /2)\mathrm{exp}(\pi /4s^2)`$, whereas in the case of the “intermediate statistics” the candidate proposed in  is the semi-Poisson distribution $`P(s)=4s\mathrm{exp}(2s)`$.

In what follows we will show, that the level spacing distribution in the case of the Bohr-Mottelson model, that represents one of the candidate systems for the intermediate statistics mentioned in  is in fact very different from the proposed semi-Poisson distribution. This discrepancy can be found in fact already on a figure published in  but without discussing the problem. However,  gives an overview over the statistical properties of different variants of the Bohr-Mottelson model.

In  the following matrix model (originally introduced by Bohr and Mottelson) is presented as a possible candidate showing statistical properties similar to the semi-Poisson distribution:

$$H_{mn}=e_n\delta _{mn}+t_mt_n.$$
(1)
$`H`$ is a $`N\times N`$-matrix, $`e_n`$ are mutually independent random variables uniformly distributed over a finite interval, and the $`t_n`$ are chosen with equal absolute value squared $`t_n^2=r`$. The authors of  sketch a procedure for calculating analytically the 2-point correlation function. For small distances it should agree with $`P(s)`$, so that one can derive the slope $`P^{}(0)=\pi \sqrt{3}/2`$ of the spacing distribution at $`s=0`$. This slope is different from the GOE case, where it equals to $`\pi ^2/6`$ as well as from the slope of the semi Poisson distribution that equals to $`4`$. The difference to the GOE case is remarked, but a similar remark on the difference to the semi-Poisson is avoided. An unprejudiced reader might believe, that the correlation properties of the matrix model are similar to semi-Poisson, even though that the spectral statistics of this model differ remarkably.

In Figure 1 we present the numerical result for $`P(s)`$ obtained for an ensemble of 1000 matrices of dimension 750. For the statistical analysis we used only one third of the states in the center of the spectral region. The numbers $`e_n`$ are uniformly distributed over an interval $`[1,1]`$ and the elements of the vector $`\stackrel{}{t}`$ are chosen as $`t_i=\sqrt{\alpha /(\pi \rho N)}`$, where $`\rho `$ is the level density in the center of the spectrum, $`N`$ is the dimension of the matrix and $`\alpha =10`$ is the coupling constant. (We checked, that larger coupling does not change the numerical results). Figure 1.a demonstrates the qualitative differences in the behavior of $`P(s)`$ between the random matrix model and the semi-Poisson distribution. The slope at $`s=0`$ is smaller, the maximum of $`P(s)`$ is slightly shifted to the right, and for values $`s>1`$ it shows significant deviations well above of the statistical error. Hence the spacing distribution of the Bohr-Mottelson model is not close to semi-Poison.

Figure 1.b. shows a magnification of the interval $`0s1/2`$ using the same data as in figure 1.a. In this figure we additionally plotted the asymptotic result (1) for the present model as a dashed line. The basic idea of the approximation is the following: In order to get a short distance between two neighbored levels in the spectrum of $`H`$, three eigenvalues of $`H_0`$ have to come close together. Then the levels which are farther away can be neglected. Therefore we can restrict the sum

$$K(E)=\underset{i=1}{\overset{N}{}}\frac{t_i^2}{Ee_i}$$
(2)
whose roots define the eigenvalues of $`H`$, to those terms with the three consecutive eigenvalues. Resolving for the two roots, calculating their distance, and averaging over the levels $`e_i`$ leads to the following formula

$$P(s)=\frac{9s}{4}\underset{0}{\overset{\pi /2}{}}𝑑\varphi \frac{\mathrm{exp}\left[\frac{3s}{2}(\mathrm{cos}\varphi +\mathrm{sin}\varphi )/\sqrt{1+\mathrm{sin}(2\varphi )/2}\right]}{1+\mathrm{sin}(2\varphi )/2}.$$
(3)
The dashed curve in Figure 1.b. is obtained from a numerical integration of (3). At short distances, this approximation describes the numerical data much better than the semi-Poisson. A Taylor expansion of the integrand of (3) gives $`P^{}(0)=\pi \sqrt{3}/2`$ for the slope at $`s=0`$, the same result as found in .

A detailed numerical investigation of several statistical properties of the type of models can be found also in .

To conclude, if one wants to insist on the introduction of a “third” universal ensemble, one should possibly use a criterium similar to the following (cited from : “…the main features are (i) the existence of level repulsion (as in random matrix ensembles), and (ii) slow (approximately exponential) fall-off …”. Even though this definition is quite “spongy”, it seems to be the only way to make sure, that the systems discussed fit into this class.

Acknowledgements M. Müller acknowledges financial support from the CONNACyT (No.32101-E).


text:
# Near-IR Imaging of LSB Galaxies: the High and Low HI Content Cases

## 1. Introduction: Background and Rationale

Low surface brightness galaxies (LSBGs) have been the subject of interest since Zwicky (1957) and Disney (1976) emphasized the fact that the central surface brightnesses of disk galaxies in the Hubble sequence Sa-Sb-Sc (which fall in a rather narrow range) might be the result of a selection effect, which is due to the difficulty in discovering galaxies of very low surface brightness. In practice, LSB means galaxies whose central surface brightness is fainter than 22.0 mag arcsec-2 in the $`B`$ band, and seem to have different properties than those of brighter galaxies, which might imply a different evolutionary track.

Previous studies (Longmore et al. 1982) has shown that LSBGs are much more gas rich than high surface brightness ones, and bluer than “normal” late type galaxies (de Blok, van der Hulst & Bothun 1995; Sprayberry et al. 1995; Dalcanton, Spergel & Summers 1997). These facts combined with their low metallicity content (McGaugh 1994; de Blok & McGaugh 1997), result in a low star formation rate in an unevolved system. However, because of the low surface brightness of the underlying population, only a small fraction of the total number of stars is needed to make the colors significantly blue. Therefore, the study of the M/L ratio in the near-IR for these galaxies would lead to a clearer picture. Also, more observational constraints are needed in order to explain the blue colors and low metallicity and to disentangle the role played by age and metallicity. Combined near-IR and optical observations are a perfect probe to clarify these issues.

Knezek & Wroten (1994) observed a small sample of LSBGs in $`J`$, $`H`$ and $`K`$, but the sample was very small and biased towards massive galaxies.

The fact that the integrated colors of some LSBGs are quite red (like Malin-type objects) suggests that they contain a significant old population of stars and were therefore not formed recently. The near-IR colors are good tracers of the mass, where additionally the evolution factor is minimized by the weakly evolving dominant population of red giant stars and low mass stars. Moreover, integrated color observations by McGaugh (1992) show that only the reddest LSBGs have bulges. The super-giant galaxies Malin 1 and Malin 2 appear to be redder than the majority of LSBGs with luminosities that are comparable to those of normal spirals. Our preliminary results in the near-IR seem to agree with this picture, in the sense that more massive LSBGs (high HI mass) are redder than less massive (small HI masses), low luminosity LSBGs.

With these caveats in mind, this systematic study will provide fundamental properties of LSBGs both in the optical and in the near-IR (for example, to derive the size of these galaxies in terms of the surface brightness in both the near-IR and in the optical). This offers the possibility of tracing accurately the mass distribution using near-IR imaging, and the luminosity distribution using optical colors, helping to trace better the evolutionary tracks of these galaxies. The spread in the observed integrated optical colors of LSBGs is large, showing that such objects have had a wide range of evolutionary histories. Near-IR colors prevent the melting of the light of all the star populations to infer the evolutionary (and metallicity) path of these galaxies. Optical and near-IR colors combined give a complete picture of the spatial distribution of the young and old stellar populations, as well as allow to quantify the degree of internal absorption by dust. Moreover, recent modeling of LSBGs predict that, from the assumption that blue LSBGs are currently undergoing a period of enhanced star formation, there should exist a population of red, non-bursting, quiescent LSBGs (Gerritsen & de Blok 1999; Bell et al. 2000; Beijersbergen, de Blok, & van der Hulst 1999). These galaxies should then also be metal-poor and gas-rich, and share many of the properties of the LSBGs observed in the optical. Near-IR imaging is a valuable window to test this last hypothesis and to study the stellar content of LSBGs, allowing to compare the properties of their old stellar population with the features observed for the young populations using optical data, particularly from blue bands.

## 2. The Sample

In order to investigate these issues, I have selected two subsamples of LSBGs from the catalogue of Impey et al. (1996). The two subsamples differ in the HI mass of their galaxies. The first subsample contains LSBGs with log\[MHI/M\] $``$ 9.0 (52 galaxies, sample A) and the second subsample galaxies with log\[MHI/M\] $``$ 9.5 (58 galaxies, sample B). In the final catalogue, only galaxies with $`\mu _0(B)23.5`$ mag arcsec-2 are included. This $`B`$ limiting central magnitude corresponds roughly to the equivalent limiting magnitude in $`K`$ for the IR Classic Camera at the Las Campanas 2.5m du Pont telescope for a stellar image in one hour integration time (5$`\sigma `$ and 1.5 arcsec aperture), given that $`<BK>2.53.0`$ for LSBGs (see for example Tully & Verheijen 1997). All the galaxies already have relevant information, such as radial velocities (i.e. distance), HI content, $`B_{total}`$ magnitudes, and the central blue surface magnitude $`\mu _0`$. Also some morphological information is available (Hubble types).

## 3. Status of the Project and Prospects

Currently, most of the near-IR observations are completed and 90% of them are reduced. All the sample will be imaged in $`B`$ and $`I`$ during 2000. Therefore, the final catalogue will include high S/N ($`1015`$) $`B`$, $`I`$, $`J`$ and $`K_s`$ imaging of 110 LSBGs from the catalogue of Impey et al. (1996), as described in §2. The current work with the data includes (1) structural characterization of the near-IR images (sizes, shapes and morphology), (2) computation of different kind of calibrated magnitudes for the whole sample, (3) characterization of the photometric profiles in terms of magnitudes and colors, and (4), comparison of derived photometric results with spectrophotometric models of galaxy evolution, considering the role of metallicity and dust, e.g. using PEGASE (Fioc & Rocca-Volmerange 1997) and GISSEL96 (Charlot, Worthey, & Bressan 1996). The approach will be similar to that applied by Galaz (2000) for a sample of E+A galaxies and by Bell et al. (2000). I emphasize that this is currently a unique large and homogeneous database of LSBGs with reduced near-IR imaging.

Figures 1 and 2 show $`J`$ and $`K_s`$ mosaiced images of 38 reduced LSBGs from sample A. Figures 3 and 4 mosaiced images of 47 reduced LSBGs from sample B.

The core of the reduction pipeline uses a modified version of DIMSUM, and has been applied successfully to all the data (see also Galaz 2000). In the future, we expect to observe fainter LSBGs in the near-IR and in the optical, using larger telescopes.

### Acknowledgments.

I acknowledge Carnegie Observatories for the huge (but necessary) telescope time allocated for this project. This work is supported under agreement between Fundación Andes and Carnegie Institution of Washington.

## References

Beijersbergen, M., de Blok, W., & van der Hulst, J. 1999, A&A, 351, 903

Bell, E., Barnaby, D., Bower, R. et al. 2000, MNRAS, 312, 470

de Blok, W., van der Hulst, J., & Bothun, D. 1995, MNRAS, 274, 235

Charlot, S., Worthey, G., & Bressan, A. 1996, ApJ, 457, 625 (GISSEL96)

Dalcanton, J., Spergel, D., & Summers, F. 1997, ApJ, 482, 659

Disney, M. 1976, Nature, 263, 573

Fioc, M., & Rocca-Volmerange, B. 1997, A&A, 326, 950

Galaz, G. 2000, AJ, 119, 2118

Gerritsen, P., & de Blok, W. 1999, A&A, 342, 655

Impey, C., Sprayberry D., Irwin M., & Bothun G. 1996, ApJS, 105, 209

Knezek, P., & Wroten, J. 1994, BAAS, 185, 77.13

Longmore, A., Hawarden, T., Goss, W., Mebold, U., & Webster, B. 1982, MNRAS, 200, 325

McGaugh, S., & de Blok, W. 1997, ApJ, 481, 689

McGaugh, S. 1994, ApJ, 426, 135

Sprayberry, D., Bernstein, G., Impey, C., & Bothun, G. 1995, ApJ, 438, 72

Tully, B., & Verheijen, M. 1997, ApJ, 484, 145

Zwicky, F. 1957, Morphological Astronomy (Berlin: Springer)


text:
# Neutron Stars and Quantum Billiards

## Abstract

Homogeneous neutron matter at subnuclear densities becomes unstable towards the formation of inhomogeneities. Depending on the average value of the neutron density one can observe the appearance of either bubbles, rods, tubes or plates embeded in a neutron gas. We estimate the quantum corrections to the ground state energy (which could be termed either shell correction or Casimir energy) of such phases of neutron matter. The calculations are performed by evaluating the contribution of the shortest periodic orbits in the Gutzwiller trace formula for the density of states. The magnitude of the quantum corrections to the ground state energy of neutron matter are of the same order as the energy differences between various phases.

PACS numbers: 21.10.Dr, 21.65.+f, 97.60.Jd, 05.45.Mt

It is somewhat surprising to find out that the problem we are going to present here, which could have naturally emerged as a typical study case in quantum chaos, has somehow eluded the attention of this community. Various billiard problems, in particular the Sinai billiard in 2D, are paradigmatic in quantum chaos studies. The structures believed to appear in the crust of neutron stars, namely bubbles, tubes, rods and plates, known also as the ”pasta phase” , are nothing else but one of several natural realizations of Sinai billiards. The energetics of these and similar phases at higher densities involving quarks has been studied for many years  (this reference list is not meant to be exhaustive).

All the studies mentioned above have estimated the ground state energy of neutron matter only in the liquid drop or Thomas–Fermi approximation. In this approximation applied to finite fermi systems – such as an atom, a nucleus, a metallic cluster or a quantum dot – one fails to put in evidence the fact that in the case of closed shells these systems show an enhanced stability. Apparently, an agreement has been reached in literature concerning the existence of the following chain of phase changes as the density is increasing: nuclei $``$ rods $``$ plates $``$ tubes $``$ bubbles $``$ uniform matter. There are a couple of studies, where the role of quantum corrections has been taken into account, one being the Hartree–Fock calculations of Refs. . In Ref.  the shell effects due to the bound nucleons in nuclei immersed in a neutron gas only (mainly protons) have been computed. There it was concluded that this type of quantum corrections to the ground state energy will not lead to any qualitative changes in the sequence of the nuclear shape transitions in the neutron star crust. However, no other type of quantum corrections was suggested. We have shown that the quantum corrections to the ground state energy of neutron matter arising from the unbounded motion of fermions outside the bubbles in particular, are more important quantitatively and thus can lead to qualitative changes in the physics of neutron matter at subnuclear densities . In particular, since the magnitude of the quantum corrections to the ground state energy is of the same order of magnitude as the energy differences between various “pasta phases”, there is no reason to expect that the sequence “nuclei $``$ rods $``$ plates $``$ tubes $``$ bubbles $``$ uniform matter” is occurring and instead strong disorder could dominate over various regular lattices.

Let us start our analysis with reviewing a somewhat simpler case, that of systems with only one bubble . About half a century ago it was predicted that if a large nucleus was ever to be observed, most likely it would have either the shape of a torus or of a sphere with an empty spherical cavity inside (this is what is typically referred to as a bubble nucleus) . Such geometries lead to a minimum of the sum of Coulomb and surface energies. Recently it was suggested that many fermion systems with a similar geometry can be created in the case of highly charged atomic clusters . In order to focus our attention on the quantum corrections, we shall discard here the role of long range Coulomb interaction, which is easy to account for, and we shall sidestep the issue of bubble stability as well and refer the interested reader to earlier literature. In the case of atomic clusters the Coulomb energy plays a secondary role anyway . For the sake of simplicity, we shall consider only one type of fermions with no electric charge.

It is known that in the case of a saturating $`N`$–fermion system the total energy has the general structure

$$E(N)=e_vN+e_sN^{2/3}+e_cN^{1/3}+E_{sc}(N).$$
(1)
The first three terms represent the smooth, liquid drop part of the total energy. $`E_{sc}(N)`$ is the pure quantum contribution, known as shell correction energy, the amplitude of which grows in magnitude approximately as $`N^{1/6}`$, see Ref. . The liquid drop part of the energy in the case of billiard systems is referred to as the Weyl energy .

One of the simplest questions one can raise is: How should one determine the position of a bubble or void? Considerations based on naive symmetry arguments have been used in the past, without any critical consideration, to simply position the bubble in the very center of the system. Recently we have shown that moving the bubble off–center can often lead to a greater stability of the system, due to shell correction energy effects alone . It is worth noting at this point that the shell correction energy bears a remarkable similarity to the Casimir energy in quantum field theory and in critical phenomena . The Casimir energy can be thought of as a specific measure of the magnitude and nature of the fluctuations induced in the energy spectrum by the presence of various “obstacles” and is computed with a formula very similar to the one used to estimate the shell correction energy . There is a large number of studies of various fermion systems, which shows that the character of the shell corrections is strongly correlated with the existence of regular or chaotic motion in the classical limit . Integrable systems have conspicuous large gaps in their spectra, while chaotic systems have a spectrum characterized by rather small fluctuations and essentially no large gaps. If a spherical bubble appears in a spherical system and if the bubble is positioned at the center, then for certain magic fermion numbers the shell correction energy $`E_{sc}(N)`$, and hence the total energy $`E(N)`$, has a very deep minimum. However, if the number of particles is not magic, in order to become more stable the system will in general tend to deform. Real deformations lead to an increased surface area and thus to an increased liquid drop energy. In the case of systems with bubbles/voids there is a somewhat unexpected extra soft deformation mode – bubble positioning. Merely shifting a bubble off–center deforms neither the bubble nor the external surface and therefore, the liquid drop part of the total energy of the system remains unchanged. Therefore, one would expect that classically this deformation mode has exactly zero energy. Quantum mechanics at first glance seems to lend some support to this idea. In recent years it was shown that in a 2–dimensional annular billiard, which is the 2–dimensional analog of spherical bubble nuclei, the motion becomes more chaotic as the bubble is moved further from the center . One might anticipate then that the importance of the shell corrections diminishes when the bubble is off–center and thus becomes easier to move a bubble, the further from the center it is. One might also naively assume that the most favorable arrangement is the one with the bubble at the very center, when the system has the highest symmetry. It came as somewhat of a surprise to find out that this is not the case most of the time and that the profile of the energy surface has a very unexpected structure .

We shall refer the reader to our previous work for more details concerning the physics of systems like bubble nuclei and revert now to the case of bubbles, rods, tubes and plates floating in an infinite neutron gas. All calculations we present here for infinite systems were performed at constant chemical potential. In order to better appreciate the nature of the problem we present here, let us consider the following situation. Let us imagine that two spherical identical bubbles have been formed in an otherwise homogeneous neutron matter. One can then ask the following apparently innocuous question: “What determines the most energetically favorable arrangement of the two bubbles?” According to a liquid drop model approach the energy of the system should be insensitive to the relative positioning of the two bubbles. A similar question was raised in condensed matter studies, concerning the interaction between two impurities in an electron gas. In the case of two “weak” and point–like impurities the dependence of the energy of the system as a function of the relative distance between the two impurities $`𝐚`$ is given by

$$E(𝐚)=\frac{1}{2}𝑑𝐫_1𝑑𝐫_2V_1(𝐫_1)\chi (𝐫_1𝐫_2𝐚)V_2(𝐫_2),$$
(2)
where $`\chi (𝐫_1𝐫_2𝐚)`$ is the Lindhard response function of a homogeneous Fermi gas and $`V_1(𝐫_1)`$ and $`V_2(𝐫_2)`$ are the potentials describing the interaction between impurities and the surrounding electron gas. At large distances, when $`k_Fa1`$, this expression leads to an interaction first derived by Ruderman and Kittel :

$$E(𝐚)\frac{\mathrm{cos}(2k_Fa)}{a^3},$$
(3)
where $`k_F`$ is the Fermi wave vector and $`m`$ is the fermion mass. This asymptotic behavior is valid only for point–like impurities, when $`k_FR1`$, where $`R`$ stands for the radius of the two impurities. This condition is typically violated for nuclei and bubbles immersed in a neutron gas, for which $`k_FR1`$. We have shown in Ref.  that in the case of large “impurities” ($`k_FR1`$) and at constant chemical potential the interaction energy at large separations $`a`$ becomes

$$E_{}\frac{\mathrm{}^2k_F^2}{2m}\left(\frac{R}{a}\right)^2\frac{2\mathrm{sin}(2k_Fa)}{\pi ^2}.$$
(4)
The most striking aspect of this interaction is the fact that it is inversely proportional to the square of the separation, instead of $`1/a^3`$ as in Eq. (3). The bubble–bubble interaction has a surprinsingly long range. Apart from Coulomb and gravitational interactions, there are no other longer ranged interactions known between two objects. One can advance a simple qualitative explanation of the difference in the power law exponent of the impurity–impurity interaction (3) and the bubble–bubble interaction (4) in terms of waves. A large “obstacle” reflects back more of the incident wave than a point object, which acts like a pure $`s`$–wave scatterer. We shall see below that perfectly planar surfaces lead to an interaction which decays with distance even slower, inversely proportional to the separation between surfaces. Naturally, when the incident and reflected waves interfere constructively, one does expect a more energetically favorable geometric configuration. Obviously, the mere fact that the interaction (4) oscillates shows that for each value of $`a=\pi (n+3/4)/k_F`$, where $`n`$ is an integer, the potential has a minimum and one thus expects that a bubble–bubble molecule with such a radius can be formed. Moreover, for each such size various vibrational and rotational states, corresponding to small fluctuations of the bubble–bubble distance and angular velocity are most likely to exist. Tunneling between various molecular sizes should also appear. However, the treatment of such molecular systems and tunneling should be deferred until the bubble inertia is also estimated.

The arguments of the cosine in Eq. (3) and of the sine in Eq. (4) are nothing else but the classical action in units of $`\mathrm{}`$ of the bouncing periodic orbit between the two impurities and this suggests the most likely method of attack for this problem, the semiclassical approximation based on the Gutzwiller trace formula for the density of states. As a matter of fact this is the way Eq. (4) has been derived . The formation of various inhomogeneities in an otherwise uniform Fermi gas can be characterized by several natural dimensionless parameters: $`k_Fa1`$, where as above $`a`$ is a characteristic separation distance between two such inhomogeneities; $`k_FR1`$, where $`R`$ is a characteristic size of such an inhomogeneity; and $`k_Fs=𝒪(1)`$, where $`s`$ is a typical “skin” thickness of such objects. The fact that the first two parameters, $`k_Fa`$ and $`k_FR`$, are both very large makes the adoption of the semiclassical approach natural. Since the third parameter, $`k_Fs`$, is never too large or too small, one might be tempted to discard the semiclassical treatment in the case of saturating fermion systems altogether as unreasonable. However, there is a large body of evidence pointing towards the fact that even though this last parameter in real systems is of order unity, the seemingly unreasonable approximation $`k_Fs1`$, which we shall adopt in this work, is surprisingly accurate . The corrections arising from considering $`k_Fs=𝒪(1)`$ lead mainly to an overall energy shift, which is largely independent of the separation among various objects embedded in a Fermi gas. On one hand, this type of shift can be accounted for in principle in a suitably implemented liquid drop model or Thomas–Fermi approximation. On the other hand, the semiclassical corrections to the ground state energy arising from the relative arrangement of various inhomogeneities have to be computed separately, as they have a different physical nature. We are thus lead to the natural assumption that a simple hard–wall potential model for various types of inhomogeneities appearing in a neutron fermi gas is a reasonable starting point to estimate quantum corrections to the ground state energy, see Refs.  and earlier references therein. One might expect that such simplifications will result in an overestimation of the magnitude of quantum corrections to the ground state energy, but that the qualitative effect would be reproduced.

One can compute the quantum corrections to the ground state energy analytically for plate–like nuclei or voids with the neutron gas filling the space between slabs. The quantum correction to the energy for such a system is given by:

$$\frac{E_{shell}}{L^3}=\frac{EE_{Weyl}}{L^3},$$
(5)
where the exact and the Weyl (smooth) energy  per unit volume are given by

$`{\displaystyle \frac{E}{L^3}}`$ $`=`$ $`{\displaystyle \frac{2}{L^3}}{\displaystyle \frac{\mathrm{}^2}{2ma^2}}{\displaystyle \frac{\pi ^3}{2}}\left({\displaystyle \frac{L}{a}}\right)^2\left[{\displaystyle \frac{1}{4}}\left({\displaystyle \frac{k_Fa}{\pi }}\right)^4N{\displaystyle \frac{N(N+1)(2N+1)(3N^2+3N1)}{120}}\right],`$ (6)
$`{\displaystyle \frac{E_{Weyl}}{L^3}}`$ $`=`$ $`{\displaystyle \frac{2}{L^3}}{\displaystyle \frac{\mathrm{}^2}{2ma^2}}{\displaystyle \frac{\pi ^3}{2}}\left({\displaystyle \frac{L}{a}}\right)^2\left[{\displaystyle \frac{1}{5}}\left({\displaystyle \frac{k_Fa}{\pi }}\right)^5{\displaystyle \frac{1}{8}}\left({\displaystyle \frac{k_Fa}{\pi }}\right)^4\right],`$ (7)
and where

$$N=\mathrm{Int}\left[\frac{k_Fa}{\pi }\right]$$
(8)
stands for the integer part of the argument in the square brackets, and $`a=L2R`$ is the distance between slabs and $`R`$ is half of the width of the slab. Here $`L^3`$ is the volume of an elementary (cubic) cell and the factor “2” in front stands for the two spin states. Using these formulas, one can show that the shell correction energy behaves as

$$\frac{E_{shell}}{L^3}=\frac{\mathrm{}^2k_F^3}{48aLm}G\left(\frac{k_Fa}{\pi }\right),$$
(9)
where $`G(x)`$ is an approximately periodic function of its argument. For $`x1`$, $`G(x+1)G(x)`$, with properties $`G(x=n)1`$ and $`1G(x)0.5`$. Note that at small separations $`a`$ the shell correction energy is attractive in character.

In the case of rods or tubes and bubbles we resorted to the Gutzwiller trace formula to estimate the corrections to the density of states. We were interested in what in the nuclear physics lingo would be called the “gross shell structure” and we discarded the fine details. As this approach implies a certain amount of spectral averaging one needs to account only for the existence of the shortest periodic orbits. As we are dealing with an infinite system there are no discrete levels and correspondingly no shells. However, the presence of “obstacles” induces fluctuations in both spatial matter distribution and density of states, which thus lead to the appearance of what we term shell correction or Casimir energy. For the case of spherical voids there are $`26`$ periodic orbits between the nearest neighbors of three different lengths $`2L_1=2(L2R)`$, $`2L_2=2(L\sqrt{2}R)`$ and $`2L_3=2(L\sqrt{3}R)`$ in the case of a simple cubic lattice. Thus the shell correction/Casimir energy and density are equal to:

$`{\displaystyle \frac{E_{shell}}{L^3}}`$ $`=`$ $`{\displaystyle \frac{1}{L^3}}{\displaystyle _0^\mu }\epsilon {\displaystyle \underset{i=1}{\overset{3}{}}}A_ig_{shell}(\epsilon ,L_i)d\epsilon ,`$ (10)
$`\rho _{out}`$ $`=`$ $`\rho _{Weyl}+{\displaystyle \frac{1}{L^3}}{\displaystyle _0^\mu }{\displaystyle \underset{i=1}{\overset{3}{}}}A_ig_{shell}(\epsilon ,L_i)d\epsilon .`$ (11)
The contribution due to one periodic orbit to the fluctuating part of the level density reads:

$$g_{shell}(\epsilon ,L_i)=\frac{mL_i}{2\pi \mathrm{}^2k}\underset{n=1}{\overset{\mathrm{}}{}}\frac{\mathrm{cos}(2nkL_i)}{\mathrm{sinh}^2(n\kappa _i)},$$
(12)
where

$`\epsilon `$ $`=`$ $`{\displaystyle \frac{\mathrm{}^2k^2}{2m}},`$ (13)
$`\kappa _i`$ $`=`$ $`\mathrm{ln}\left[1+{\displaystyle \frac{L_i}{R}}+\sqrt{{\displaystyle \frac{L_i}{R}}\left({\displaystyle \frac{L_i}{R}}+2\right)}\right].`$ (14)
and hence we get:

$`{\displaystyle \frac{E_{shell}}{L^3}}={\displaystyle \underset{i=1}{\overset{3}{}}}{\displaystyle \frac{\mathrm{}^2A_i}{16m\pi ^2L_i^2L^3}}{\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{[2nk_FL_i\mathrm{cos}(2nk_FL_i)+(2n^2k_F^2L_i^21)\mathrm{sin}(2nk_FL_i)]}{n^3\mathrm{sinh}^2(n\kappa _i)}},`$ (15)
$`\rho _{out}=\rho _{Weyl}+{\displaystyle \frac{1}{L^3}}{\displaystyle \frac{1}{4\pi }}{\displaystyle \underset{i=1}{\overset{3}{}}}A_i{\displaystyle \underset{n=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\mathrm{sin}(2nk_FL_i)}{n\mathrm{sinh}^2(n\kappa _i)}},`$ (16)
where $`A_1=6,A_2=12,A_3=8`$ respectively. Above, the summation over $`n`$ is over repetitions of a given orbit and such summation becomes superfluous when one performs a spectral averaging in order to extract the “gross shell structure”. This shell energy can be regarded as the interaction energy between bubbles immersed in neutron gas. The bubble–bubble interaction (4) was obtained from the above formula for the case of two bubbles only in the asymptotic limit. One can derive similar formulas for the case of rods . The fact that the interaction energy between various homogeneities is described in terms of periodic orbits, points to another nontrivial aspect: the interaction between three or more such objects cannot be reduced to simple pairwise interactions. Orbits which bounce among three or more objects would lead to three–body, four–body and so forth genuine interactions. Such orbits were not yet included in our analysis, even though they can be accounted for in a straight forward manner. Since periodic orbits of this kind are typically longer, one expects their contribution to be masked when one performs a spectral averaging in order to extract the “gross shell structure”.

Lack of space prevents us from presenting more detailed results and discussing many other aspects of the energetics of inhomogeneous fermi systems. Hopefully, the few results we have highlighted here will convince the reader of the richness of these systems and of the spectacular role played by geometry and the chaotic versus integrable character of the single–particle dynamics. Besides static properties, one should expect that the dynamics of such systems will be extremely rich and challenging to describe.