text:
# Stability of Color-Flavor Locked Strangelets

## I Introduction

The properties of quark matter have attracted a good deal of interest since it has been suggested that strange quark matter could be a ground state of strongly interacting matter.EW ; SQM  The most likely place for strange quark matter, perhaps, is the interior of compact stars even if it is not absolutely stable. The existence of absolutely stable strange quark matter is still an open question and it may be realized in the form of strangelets (small lumps of strange quark matter with roughly the same amount of up, down and strange quarks). Theoretically, there has been various investigations of the stability of (non)strange quark matter. Within the MIT bag model, Farhi and Jaffe FJ  found a reasonable window of the model parameters (i.e. the bag constant, current quark masses and strong coupling constant) for which strange quark matter is stable, while nonstrange quark matter is unstable as compared to a gas of $`{}_{}{}^{56}\text{Fe}`$. By including finite size effects (surface tension and/or curvature energy), similar analyses have been done for finite lumps of quark matter.SQM ; FJ ; BJ ; GJ ; JM  It is now well established that finite size effects increase the energy of finite quark lumps.

During the last decade, significant advances have been made in our understandings of the phase structure of hot and/or dense QCD. CSC  At the present time, it is widely accepted that a color-flavor locked (CFL) phase Alford98  is the ground state of cold, dense quark matter. 111As discussed in Refs. IB ; AR ; SRP ; NBO , it is of importance to take account of the constraint imposed by color and electric neutrality. For neutrino free, color and electrically neutral quark matter at zero temperature, it has been shown that the CFL phase optimizes the pairing energy and, then, is favored over the 2-flavor superconducting (2SC) phase in all (or almost all) the range of densities.AR ; SRP  Madsen CFLS  has studied the stability of CFL strangelets and found that CFL strangelets are significantly more stable than normal (unpaired) strangelets (see also Ref. Lugones1 ). However, that work ignored issues of density dependence of the bag constant and the CFL gap. What remains a question is the possible effects of the phase structure. In the previous paper KH , we have studied the chiral and 2SC phase of finite quark lumps. We found that finite size effects enhance the restoration of chiral symmetry. In this case, it is likely that strangelets lie in the color superconducting phase.

The purpose of this work is to study the behavior of the CFL gap in strangelets and to look at its effects on the stability of CFL strangelets. In order to describe the density dependent bag constant and the density dependent CFL gap we choose to use the three-flavor Nambu–Jona-Lasinio (NJL) model.NJ  (The NJL model is a simple tractable model to investigate the stability of quark matter/droplets \[see, for example, Ref. eos \].) We consider up, down and strange quarks to be massless, for simplicity. In order to take account of finite size effects, we apply what is called multiple reflection expansion (MRE).MRE  The MRE has been used for calculating thermodynamic quantities such as the energy per baryon number and the free energy of finite quark lumps. As far as the general structure is concerned, the results are in good agreement with those given in the MIT bag model.SQM ; JM  Using the NJL model with the MRE, we formulate the thermodynamic potential for spherical CFL strangelets. We then investigate the stability of CFL strangelets, including the dynamical effects as well as the finite size effects.

The organization of this paper is as follows. In the next section, we briefly review the thermodynamic potential of the CFL phase in bulk. Then, we formulate the thermodynamic potential for CFL strangelets. The gap equation and related thermodynamic quantities are also derived. In Sec. 3, we present numerical results. Finally, Sec. 4 is devoted to the conclusions.

## II Thermodynamic potential

First, we briefly survey the thermodynamic potential of the CFL phase in bulk. To describe the CFL phase in the NJL model we follow Buballa and Oertel.BO2  The $`\text{U}(3)_L\times \text{U}(3)_R`$ symmetric NJL Lagrangian is

$`=\overline{q}i\gamma ^\mu _\mu q+G_1{\displaystyle \underset{i=0}{\overset{8}{}}}\left[(\overline{q}\tau ^iq)^2+(\overline{q}i\gamma _5\tau ^iq)^2\right],`$ (1)
where $`q`$ denotes a quark field with three flavors ($`N_f=3`$) and three colors ($`N_c=3`$), the coupling constant $`G_1`$ has a dimension $`[G_1]=[\text{mass}]^2`$ and the Gell-Mann matrices $`\tau ^i`$ $`(i=1,\mathrm{},8)`$ with $`\tau ^0=\sqrt{2/3}𝑰_f`$ act in the flavor space.

In the case of the CFL, the quarks form the following $`(\overline{3}_c,\overline{3}_f)`$ condensate

$`\mathrm{\Delta }_{ij}^{\alpha \beta }=q_i^\alpha q_j^\beta C\gamma _5ϵ^{\alpha \beta X}ϵ_{ijX},`$ (2)
where $`i,j`$ denote flavor indices, $`\alpha ,\beta `$ denote color indices and repeated indices are summed. In principle, there must be the additional $`(6_c,6_f)`$ condensate. However, we neglect such condensate because it is much smaller than $`(\overline{3}_c,\overline{3}_f)`$ condensate Alford98 ; Shovkovy99  and its contribution to the thermodynamic potential is negligible. The interaction corresponding to the $`(\overline{3}_c,\overline{3}_f)`$ CFL paring can be written as

$`_{\text{CFL}}=G_2{\displaystyle \underset{i=\alpha =2,5,7}{}}\left(\overline{q}i\gamma _5\tau ^i\lambda ^\alpha C\overline{q}^T\right)\left(qCi\gamma _5\tau ^i\lambda ^\alpha q\right),`$ (3)
where $`C`$ is a charge conjugation matrix, defined by $`C^1\gamma _\mu C=\gamma _\mu ^T`$ and $`C^T=C`$, and $`\tau ^i`$ $`(i=2,5,7)`$ \[$`\lambda ^\alpha `$ $`(\alpha =2,5,7)`$\] denote the antisymmetric generators of $`\text{SU}(3)_f`$ \[$`\text{SU}(3)_c`$\]. The coupling constant $`G_2`$ can be obtained from Eq. (1) by making use of the Fierz transformation. However, we leave $`G_2`$ a free parameter. In this approximation, all nine (three colors times three flavors) quarks participate in the pairing and yield eight quasiparticles with gap $`\mathrm{\Delta }`$ and one with gap $`2\mathrm{\Delta }`$. (These quasiparticles correspond to an octet and a singlet of the unbroken $`\text{SU}(3)`$, respectively.) Then, in the mean-field approximation, the thermodynamic potential $`\mathrm{\Omega }=\mathrm{\Omega }(\mathrm{\Delta };\mu ,T)`$ at finite quark chemical potential $`\mu `$ and temperature $`T`$ is given by

$`\mathrm{\Omega }`$ $`=`$ $`{\displaystyle \frac{3\mathrm{\Delta }^2}{4G_2}}8T{\displaystyle \underset{\pm }{}}{\displaystyle \underset{n=\mathrm{}}{\overset{\mathrm{}}{}}}{\displaystyle \frac{d^3k}{(2\pi )^3}\mathrm{ln}(\omega _n^2+ϵ_\pm ^2)}`$ (4)
$`T{\displaystyle \underset{\pm }{}}{\displaystyle \underset{n=\mathrm{}}{\overset{\mathrm{}}{}}}{\displaystyle \frac{d^3k}{(2\pi )^3}\mathrm{ln}(\omega _n^2+\xi _\pm ^2)},`$
where $`ϵ_\pm =\sqrt{(k\pm \mu )^2+\mathrm{\Delta }^2}`$, $`\xi _\pm =\sqrt{(k\pm \mu )^2+4\mathrm{\Delta }^2}`$ and $`\omega _n=(2n+1)\pi T(n=0,\pm 1,\pm 2,\mathrm{})`$ denote the fermionic Matsubara frequencies. The summation over $`\omega _n`$ is straightforward. Henceforth, we restrict ourselves to $`T=0`$. The $`T0`$ limit of $`\mathrm{\Omega }`$ is

$`\mathrm{\Omega }={\displaystyle \frac{3\mathrm{\Delta }^2}{4G_2}}{\displaystyle \underset{\pm }{}}{\displaystyle _0^{\mathrm{\Lambda }_{\text{UV}}}}{\displaystyle \frac{k^2dk}{2\pi ^2}}\left(8ϵ_\pm +\xi _\pm \right).`$ (5)
Here, we have introduced the ultraviolet cutoff $`\mathrm{\Lambda }_{\text{UV}}`$ in the three-dimensional momentum space.

Now we incorporate finite size effects into the thermodynamic potential. To this end, we use the density of states derived from the MRE.MRE ; JM ; BJ  In the MRE framework, the density of states for a spherical system is given by $`k^2\rho _{\text{MRE}}/(2\pi ^2)`$, where

$`\rho _{\text{MRE}}`$ $`=`$ $`\rho _{\text{MRE}}(k,m,R)`$ (6)
$`=`$ $`1+{\displaystyle \frac{6\pi ^2}{kR}}f_S\left({\displaystyle \frac{k}{m}}\right)+{\displaystyle \frac{12\pi ^2}{(kR)^2}}f_C\left({\displaystyle \frac{k}{m}}\right)+\mathrm{},`$
with $`m`$ being the Dirac mass of quark and $`R`$ being the radius of the sphere. The second (third) term on the right-hand side represent the surface (curvature) contribution to the fermionic density of states. The ellipsis implies higher order terms in $`1/R`$, which are neglected throughout. In the case of massless fermions, the functions $`f_S`$ and $`f_C`$ have the following limits,

$`\underset{m0}{lim}f_S(k/m)=0,\underset{m0}{lim}f_C(k/m)=1/(24\pi ^2).`$ (7)
Then, we use the following MRE density of states,

$`\rho _{\text{MRE}}=1{\displaystyle \frac{1}{2(kR)^2}}.`$ (8)
One can see that the finite size effects reduce the density of states (the decreasing tendency is more pronounced at low momenta) and $`\rho _{\text{MRE}}`$ becomes negative at small $`kR`$. To avoid the unphysical negative density of states we shall introduce an infrared cutoff $`\mathrm{\Lambda }_{\text{IR}}=\sqrt{2}/(2R)`$ in the momentum space.

Using the density of states (8), we express the effective potential of the spherical CFL strangelets as follows,

$`\mathrm{\Omega }_{\text{MRE}}={\displaystyle \frac{3\mathrm{\Delta }^2}{4G_2}}{\displaystyle \underset{\pm }{}}{\displaystyle _{\text{MRE}}}\left(8ϵ_\pm +\xi _\pm \right),`$ (9)
where we have introduced the following notation,

$`{\displaystyle _{\text{MRE}}}={\displaystyle _{\mathrm{\Lambda }_{\text{IR}}}^{\mathrm{\Lambda }_{\text{UV}}}}{\displaystyle \frac{k^2dk}{2\pi ^2}}\rho _{\text{MRE}}.`$ (10)
Here, we emphasize that the Fermi momentum $`k_F(=\mu )`$ in Eq. (9) is common to all nine quarks. In the case of the three flavors of massless quarks, a strangelet with the common Fermi momentum is automatically color and electrically neutral.

A strangelet must be in a color singlet state. As noted in Refs. AR ; ABMW , color neutrality is nothing but a prerequisite condition for color singletness. However, it has been shown that, as long as a quark lump of color superconductor is color neutral, the projection over a color singlet state makes a negligible contribution to the thermodynamic potential.ABMW  Thus, color neutrality is a good approximation to color singletness. A small strangelet $`(A10^7)`$ does not need to satisfy electric neutrality. This is because the electron Compton wavelength is larger than such a strangelet and, then, electrons mainly stay outside of the quark phase.SQM  Therefore the phases of small strangelets could be different from that of quark matter in bulk.

We focus on CFL strangelets embedded in a vacuum and derive a set of coupled equations with the MRE. For computation of a finite system, we choose a fixed radius $`R`$. First, we derive the gap equation, which is the extremum condition of $`\mathrm{\Omega }_{\text{MRE}}`$ with respect to $`\mathrm{\Delta }`$:

$`{\displaystyle \frac{\mathrm{\Omega }_{\text{MRE}}}{\mathrm{\Delta }}}={\displaystyle \frac{3\mathrm{\Delta }}{2G_2}}4\mathrm{\Delta }{\displaystyle \underset{\pm }{}}{\displaystyle _{\text{MRE}}}\left({\displaystyle \frac{2}{ϵ_\pm }}+{\displaystyle \frac{1}{\xi _\pm }}\right)=0.`$ (11)
In addtion, we need to take account of the pressure balence relation. The relative pressure inside strangelet $`P_{\text{MRE}}`$ is given by

$`P_{\text{MRE}}=p_{\text{MRE}}p_{\text{vac}},`$ (12)
where $`p_{\text{MRE}}=\mathrm{\Omega }_{\text{MRE}}`$ and $`p_{\text{vac}}`$ has been introduced to measure the pressure relative to outside of strangelet. Note that $`p_{\text{vac}}`$ is nothing but the pressure of chirally broken vacuum. In the mean-field approximation, its value can be determined by using the Lagrangian (1) as follows:

$`p_{\text{vac}}=2N_fN_c{\displaystyle _0^{\mathrm{\Lambda }_{\text{UV}}}}{\displaystyle \frac{d^3k}{(2\pi )^3}}\sqrt{k^2+M^2}{\displaystyle \frac{3M^2}{8G_1}},`$ (13)
Here, $`M`$ is the dynamically generated quark mass, which is a nontrivial solution to the following equation,

$`{\displaystyle \frac{3M}{4G_1}}2N_fN_c{\displaystyle _0^{\mathrm{\Lambda }_{\text{UV}}}}{\displaystyle \frac{d^3k}{(2\pi )^3}}{\displaystyle \frac{M}{\sqrt{k^2+M^2}}}=0.`$ (14)
We solve Eqs. (11) and (12) self-consistently and, then, compute the baryon number of the strangelet $`A`$ by making use of the following relation,

$`A=Vn_B,`$ (15)
where $`V=4\pi R^3/3`$ is the volume of the spherical strangelet and

$`n_B={\displaystyle \frac{1}{3}}{\displaystyle \frac{\mathrm{\Omega }_{\text{MRE}}}{\mu }}={\displaystyle \frac{1}{3}}{\displaystyle \underset{\pm }{}}{\displaystyle _{\text{MRE}}}\left[{\displaystyle \frac{8(k\pm \mu )}{ϵ_\pm }}{\displaystyle \frac{k\pm \mu }{\xi _\pm }}\right],`$ (16)
is the baryon number density of the strangelet.

## III Numerical results

In numerical calculation, we use the following set of parameters: $`\mathrm{\Lambda }_{\text{UV}}=0.6`$ GeV, $`G_1\mathrm{\Lambda }^2=2.311`$.BO2  At $`T=\mu =0`$, these parameters yield the constituent mass $`M=0.35`$ GeV and the bag constant $`B^{1/4}=0.182`$ GeV. The quark-quark coupling constant $`G_2`$ has been set to be $`G_2/G_1=3/4`$ in the literature.BO2 ; SRP ; GNA  This relation between $`G_1`$ and $`G_2`$ is obtained from a four-fermion interaction with a quantum number of a one-gluon exchange,

$`_{\text{OGE}}=g{\displaystyle \underset{\alpha =1}{\overset{8}{}}}(\overline{q}\gamma ^\mu \lambda ^\alpha q)^2,`$ (17)
by making use of the Firez transformation. Moreover, this is consistent with that determined by fitting the nucleon mass within a Fadeev approach.BT  However, we vary $`G_2`$ in the range of $`G_2/G_1=(1/21)`$ to see its effect on the stability of strangelets.

We begin with the discussion of the solution to the gap equation (11). Figure 1 shows the CFL gap for finite volumes as a function of the quark chemical potential. We find that the finite size effects reduce the gap. However, the decrease in the gap is not pronounced even at small radii; therefore, the $`R\mathrm{}`$ limit is a good approximation as to the size of the gap. 222In this respect the result is similar to that in Ref. ABMW  where finite size effects on the 2SC gap have been studied by confining the quarks to a cubic box. This is easily understood as follows. Color superconductivity is brought about by the quarks near the Fermi surface. Recall that at a fixed radius $`R`$ the decrease in the density of states is not pronounced at large momenta. Then, if the Fermi momentum is sufficiently large, the finite size effects do not greatly affect the gap.

In order to look at CFL strangelets, we solve the gap equation (11), the pressure balance relation (12) and the baryon number relation (15), self-consistently. By numerically solving these equations, we obtain the baryon number dependence of the gap, the quark chemical potential and the radius. After that the energy per baryon number of the pressure-balanced strangelet,

$`{\displaystyle \frac{E}{A}}|_{p=0}={\displaystyle \frac{}{n_B}}|_{p=0}=3\mu ,`$ (18)
can be evaluated. In Eq. (18), $`=\mathrm{\Omega }_{\text{MRE}}+\mu _an_aϵ_{\text{vac}}`$ denotes the energy density of the strangelet, where $`a(u,d,s)`$, $`n_u=n_d=n_s(=n_B)`$ are the quark number densities, and $`ϵ_{\text{vac}}(=p_{\text{vac}})`$ is the energy density of the chirally broken vacuum.

Figure 2 shows the CFL gap as a function of the baryon number for the cases of $`G_2/G_1=1/2,3/4,1`$. Each curve remains approximately constant as long as $`A`$ is not too small. The change of the gap at small baryon numbers $`(A100)`$ is a consequence of the fact that the actual quark chemical potential increases when $`A`$ decreases. We note that the results for $`A10`$ should not be considered robust ones, because in this regime the quark chemical potentials are close to the ultraviolet cutoff $`\mathrm{\Lambda }_{\text{UV}}`$.

In Fig. 3, the energy per baryon number $`E/A`$ for CFL strangelets is plotted as a function of the baryon number. For comparison, we also present the result for normal strangelets. As is evident from the figure, CFL strangelets are considerably more stable than normal strangelets. Needless to say, this decrease in $`E/A`$ is due to the pairing energy contribution to the thermodynamic potential $`\mathrm{\Omega }_{\text{MRE}}`$. Note that $`E/A`$ depends on the size of the gap; hence, the quark-quark coupling constant $`G_2`$. As $`G_2`$ is increased, the gap increases and then the contribution to $`\mathrm{\Omega }_{\text{MRE}}`$ from the pairing energy also increases. Thus, it is reasonable that the growth of $`G_2`$ tends to lower $`E/A`$ of CFL strangelets. The gap can be large enough to have a great effect on $`E/A`$. For instance, we observe that at $`G_2=G_1`$ CFL strangelets of $`A160`$ are absolutely stable. Within our model parameters ($`\mathrm{\Lambda }_{\text{UV}}`$ and $`G_1`$), $`G_2/G_1`$ needs to exceed $`0.88`$ for the existence of the absolutely stable CFL strangelets. Note also that each curve shows the typical behavior of $`E/A`$ of finite quark lumps. $`E/A(=3\mu )`$ increases with decreasing $`A`$. (As mentioned earlier, this is a consequence of the increase of the quark chemical potential by the finite size effects.) On the other hand, as $`A`$ grows, $`E/A`$ approaches the value without the finite size effects.

## IV Conclusions

In this work, we have explored the color-flavor locked phase in strangelets. We used the three-flavor Nambu–Jona-Lasinio model. To take account of finite size effects we applied the multiple reflection expansion. Using the NJL model with the MRE, we formulated the thermodynamic potential for CFL strangelets and derived related thermodynamic quantities. We clarified the behavior of the CFL gap and its effects on the stability of strangelets, although our model is nothing but a toy model, neglecting the finite strange quark mass. The CFL gap is almost independent of $`A`$ as long as $`A`$ is not too small. Due to the contribution from the formation of the CFL condensate, CFL strangelets are more stable than strangelets without the CFL. Further, if $`G_2`$ is larger than the critical value, we have a chance to find absolutely stable CFL strangerets. These results complement some of the conclusions made by Madsen.CFLS  Of course we cannot completely accept the existence of absolutely stable CFL strangelets, because we do not have a good knowledge of $`G_2`$. However, a Fadeev approach to baryons BT ; BHIT  would be a guide to the determination of $`G_2`$. Moreover, we had neglected many aspects of dense QCD. Most important is probably the finite strange quark mass $`m_s`$. When $`m_s`$ is nonzero, we have to take account of unlocking phase transition.Alford99  As pointed out in Ref. NBO , the main reason for the difficulty to find the stable 2SC matter arises from the constraint of electric neutrality. In contrast, small strangelets, which can be regarded as free from this constraint, could be in the 2SC phase.Kiriyama  It should be mentioned, however, that the increase of the chemical potential caused by the finite size effects might prefer the CFL phase to the 2SC phase. It would be of great importance to study competition with other phases (hadronic phase, 2SC phase, and so on), taking account of nonzero $`m_s`$ as well as color neutrality.IB ; AR ; SRP ; NBO  A more careful analysis is left to future studies.

It should be also noted that the MRE contains several problems concerning its reliability. First, the density of state $`\rho _{\text{MRE}}`$ should have the terms proportional to $`1/R^3`$, $`1/R^4`$, and so on. These terms would be dominant at small radii, i.e., at small baryon numbers. Furthermore, $`\rho _{\text{MRE}}`$ in Eq. (6) causes the unphysical negative density of states at small $`kR`$. Although our qualitative results for relatively large baryon numbers $`(A100)`$ would not change, it is desirable to take more rigorous way of including finite size effects. In particular, as noted in Ref. JM , a detailed study of strangelet decay modes will have to rely on shell model calculations.

In this work, we restricted ourselves to zero temperature. However, for applications such as heavy ion collision experiments and the cosmological quark-hadron phase transition, it is interesting to investigate strangelets at finite temperature. The present approach is straightforwardly generalized to nonzero temperature. Studies along this line are now in progress.


format: text

text:
# Geometrical picture of dynamical facilitation

## 1 Introduction

The view of the glass transition (see  for general reviews) as a dynamical phenomenon emerges naturally from the real-space picture of supercooled liquids developed in . This work demonstrates that many of the properties of supercooled liquids can be simply explained by the existence and statistical properties of spatially correlated dynamics, or dynamical heterogeneity . Dynamical heterogeneity has been observed both in deeply supercooled liquids , in the mildly supercooled regime in simulations , and in experiments on colloidal suspensions .

Heterogeneous dynamics, in combination with structural homogeneity, as observed in glass formers, appears in systems possessing both localized excitations and facilitated dynamics . The concept of dynamic facilitation stipulates that the dynamics of a local region is suppressed if is bordered by immobile regions, and thus reflects the idea that particles in a jammed liquid cannot rearrange unless adjacent to mobile particles. Dynamic facilitation is realized in KCMs, which are simple models of non interacting particles or excitations in a lattice whose dynamics is subject to kinetic constraints (see  for a review). Examples of KCMs are spin facilitated models, like the Fredrickson-Andersen (FA) model  or the East model , and constrained lattice gases, such as the Kob-Andersen model .

In this paper we study the effect that dynamic facilitation has on the geometry of configuration space, and how in turn this affects the dynamics . We show, by generalizing KCMs to continuous degrees of freedom, that kinetic constraints impose a metric structure on configuration space which does not affect the static distribution of states, and therefore the thermodynamics, but which has a dramatic effect on the paths between states at low energies. We show that dynamical slow-down at low temperatures occurs because the system is forced to evolve along configuration space paths of least metric distance, or geodesics, which are distinct from the paths that relax the energy. This competition between available trajectories and energy relaxation gives rise to dynamical barriers, or bottlenecks, which are unrelated to the statics of the system. This view of the dynamics is distinct from that based on the idea that a static rugged energy landscape is responsible for the properties of glass formers (see e.g. ).

This paper is organised as follows. In section 2 we discuss the geometric interpretation of dynamic facilitation. In sections 3 and 4 we illustrate our ideas with simple examples. In section 5 we state our conclusions.

## 2 Geometrical interpretation of constrained dynamics

The real-space picture of supercooled liquids of Refs.  is based on two observations. First, at low temperature very few particles in the liquid are mobile, and these mobility excitations are localized in space. Second, regions of the liquid cannot become mobile unless their neighbouring regions are mobile. These observations are implemented as follows. A supercooled fluid in $`d`$ spatial dimensions is coarse-grained into cells of linear size of the order of the static correlation length as given by the pair correlation function. Cells are classified by a scalar mobility field, $`x`$, identified by coarse-graining the system on a microscopic time scale. Mobile regions carry a free energy cost, and when mobility is low interactions between cells are not important. Adopting a thermal language, we expect static equilibrium to be determined by the non-interacting Hamiltonian,

$$H[x]=\underset{\mu =1}{\overset{N}{}}f_\mu (x^\mu ),$$
(1)
where $`\mu =1,\mathrm{},N`$ labels the lattice site, $`x^\mu `$ is the mobility field at site $`\mu `$, and $`f_\mu `$ is an analytic function of $`x^\mu `$ chosen so that the ground state of the model is such that $`|x^\mu |`$ is small. At low mobility, the distinction between single and multiple occupancy is irrelevant, as is the distinction between discrete or continuous degrees of freedom. In what follows we assume $`x^\mu `$ to be real. The dynamics of the mobility field is given by a master equation,

$$_tP(x,t)=\underset{\mu }{}𝒞_\mu \left(x\right)\widehat{_\mu }P(x,t),$$
(2)
where $`P(x,t)`$ is the probability that the system has configuration $`x\{x^\mu |\mu =1,2,\mathrm{},N\}`$ at time $`t`$. The local operators $`\widehat{_\mu }`$ are the same as those for unconstrained local dynamics. Their action on $`P`$ can be written

$$\widehat{_\mu }P\left(x\right)=\underset{x^\mu }{}w(x^\mu x^\mu )P(x^\mu ,\{x^{\nu \mu }\})\underset{x^\mu }{}w(x^\mu x^\mu )P(x^\mu ,\{x^{\nu \mu }\}),$$
(3)
where $`w(x^\mu x^\mu )`$ is the probability of going from configuration $`(x^\mu ,\{x^{\nu \mu }\})`$ to configuration $`(x^\mu ,\{x^{\nu \mu }\})`$ in unit time. We ensure a unique equilibrium configuration exists by requiring (2) to obey detailed balance with respect to (1), i.e.

$$\frac{w(xx^{})}{w(x^{}x)}=e^{\beta (H[x^{}]H[x])}.$$
(4)
The kinetic constraint, $`𝒞_\mu \left(x\right)`$, is designed to suppress the dynamics of cell $`\mu `$ when surrounded by immobile regions. It must also be such that Eq. (2) satisfies detailed balance, which is achieved if $`𝒞_\mu \left(x\right)`$ does not depend on $`x^\mu `$ itself. To reflect the local nature of dynamic facilitation we allow $`𝒞_\mu `$ to depend only on the $`\{x^\nu \}`$ of nearest neighbours of $`\mu `$, and require that $`𝒞_\mu `$ is small when local mobility is scarce. Equations (1)-(4) define a generic KCM.

We now show that Eq. (2) has a simple geometrical interpretation. Using standard methods  we pass from the master equation (2) to the Fokker-Planck equation

$$_tP(x,t)=\underset{\mu }{}_\mu ^{(\mathrm{FP})}P(x,t),$$
(5)
where

$$_\mu ^{(\mathrm{FP})}=\frac{}{x^\mu }\left\{𝒞_\mu \left(x\right)f_\mu ^{}(x^\mu )+T\frac{}{x^\mu }𝒞_\mu \left(x\right)\right\}.$$
(6)
The prime on $`f_\mu `$ denotes differentiation with respect to $`x^\mu `$. Equation (5) has a simple physical interpretation: it describes single-particle driven diffusion on an $`N`$-dimensional curved space $`\{x^\mu |\mu =1,2,\mathrm{},N\}`$ . The form of the driving term in (6) comes from the Hamiltonian (1) and the kinetic constraint. But the curvature of this space is due solely to the kinetic constraint: from the form of a Fokker-Planck equation on a curved manifold  we identify the inverse metric tensor of this space as $`g^{\mu \nu }(x)=\delta ^{\mu \nu }𝒞_\mu \left(x\right)`$. Thus the metric reads

$$g_{\mu \nu }(x)=\delta _{\mu \nu }𝒞_\mu ^1\left(x\right).$$
(7)
No summation is implied. The effect of dynamic facilitation is to endow the configuration space $`x`$ with a nontrivial metric. Since dynamic facilitation is responsible for the interesting dynamical properties of these models, we expect that these dynamical properties can be extracted from the metric tensor.

The effect of the metric on configuration space trajectories is made clearer by writing a formal solution to equation (5:

$$P(x_b,t_b|x_a,t_a)=_a^b\sqrt{|g|}𝒟xe^{𝒮[x]},$$
(8)
where $`P(x_b,t_b|x_a,t_a)`$ is the conditional probability that a system starting at point $`x_a`$ in configuration space at time $`t_a`$ is found at point $`x_b`$ at time $`t_b`$. We have introduced the determinant of the metric, $`|g|det\left(g_{\mu \nu }\right)=_\mu 𝒞_\mu ^1`$. The path integral (8) is taken over all configurations of the system, weighted by the dynamic action $`𝒮[x]`$. Up to constant terms this action reads

$$𝒮[x]=\frac{1}{T}\left(H[x_b]H[x_a]\right)+\frac{1}{T}_{t_a}^{t_b}𝑑t\left(\dot{x}^\mu g_{\mu \nu }\dot{x}^\nu +\frac{H}{x^\mu }g^{\mu \nu }\frac{H}{x^\nu }\right),$$
(9)
where dots denote differentiation with respect to time, and the summation convention for once-repeated upper and lower indices is implied.

The action (9) weights the configuration space trajectories of our system. The system’s classical trajectory, i.e. the path taken by the dynamics when $`T0`$, is one of extremum $`𝒮`$. At finite $`T`$, stochastic fluctuations cause perturbations about this path. Equation (9) shows clearly the way the various ingredients of a dynamically facilitated model influence its trajectory through configuration space. The first term depends only on the Hamiltonian (1) evaluated at the endpoints. It gives rise to Boltzmann weights in the path integral (8) and plays no role in choosing the trajectory of the system. The path is chosen by a competition between the two terms in the integral of Eq. (9). The first term, $`𝒮_1T^1𝑑t\dot{x}^\mu g_{\mu \nu }\dot{x}^\nu `$, is purely geometrical. It is the action for free motion in a curved background. It depends on the metric $`g_{\mu \nu }`$, and therefore on $`𝒞^1`$. The second term, $`𝒮_2T^1𝑑t_\mu Hg^{\mu \nu }_\nu H`$, depends on the Hamiltonian (1), but is proportional to the inverse metric $`g^{\mu \nu }`$, and so to $`𝒞`$. At low temperatures, we expect $`𝒮_2`$ to be important only far from equilibrium, when mobility is plentiful—or equivalently, when the magnitudes of the coordinates $`|x^\mu |`$ are large and $`𝒞`$ is large—and the system’s relaxation will take place by way of gradient descent on the energy surface. When mobility has diminished significantly the constraints $`𝒞`$ become small. We then expect $`𝒮_1`$ to dominate. Hence the classical trajectory of the system in the region of weak dynamic facilitation will be governed by extremum $`𝒮_1`$. This corresponds to the equation of a geodesic :

$$\ddot{x}^\mu +\mathrm{\Gamma }_{\alpha \beta }^\mu \dot{x}^\alpha \dot{x}^\beta =0.$$
(10)
The Christoffel symbols $`\mathrm{\Gamma }_{\alpha \beta }^\mu `$ are defined as $`\mathrm{\Gamma }_{\alpha \beta }^\mu \frac{1}{2}g^{\mu \nu }(g_{\nu \alpha ,\beta }+g_{\beta \nu ,\alpha }g_{\alpha \beta ,\nu })`$, and the metric tensor raises and lowers indices. Equation (10) tells us how the kinetic constraint controls the system’s trajectory in configuration space in the dynamically constrained regime. When mobility is scarce, the constraint $`𝒞`$ is small, and the metric distance between points in configuration space, $`s_{ab}_a^b\sqrt{g_{\mu \nu }dx^\mu dx^\nu }`$, is large, by virtue of (7). By minimising this distance, and so following Eq. (10), we expect the system to follow a path that does not lead directly to a configuration of lowest energy. This is a geometrical interpretation of dynamical arrest.

## 3 A simple example

We now illustrate these ideas with a simple example. We define a model of two variables, $`\{x^\mu \}=\{x,y\}`$, with each variable constraining the other. Although such a model cannot exhibit true glassiness it is useful because we can visualise directly the two-dimensional configuration space it inhabits. We choose the Hamiltonian (1) to be quadratic,

$$H=\frac{J}{2}\left(x^2+y^2\right),$$
(11)
where $`J`$ sets the energy scale. We choose the constraint functions to be quadratic also, with each coordinate constraining the other:

$$𝒞_x=y^2,𝒞_y=x^2.$$
(12)
Thus, the Fokker-Planck equation (5) is determined. Equation (7) gives the metric of this system as

$$g_{\mu \nu }=\left(\begin{array}{cc}1/y^2& 0\\ 0& 1/x^2\end{array}\right).$$
(13)
The metric is singular, and so distances in this configuration space become very large, when $`x`$ or $`y`$ vanish, reflecting the tendency of the kinetic constraint to suppress the dynamics of the system in the regime of low mobility. The space described by this metric is negatively curved with curvature scalar $`R=4\left(x^2/y^2+y^2/x^2\right)`$ . The geodesics follow from (10) which in this case read

$$y\ddot{x}+(y/x)^3\dot{y}^22\dot{x}\dot{y}=0,x\ddot{y}(x/y)^3\dot{x}^22\dot{x}\dot{y}=0.$$
(14)
To illustrate the effect of the kinetic constraint on the trajectories chosen by the system we first cast Equation (5) into a form more convenient for numerical simulation. In the usual way  we recognize that (5) is equivalent to the coupled Langevin equations

$$\left(\begin{array}{c}\dot{x}(t)\\ \dot{y}(t)\end{array}\right)=\left(\begin{array}{c}Jy^2(t)x(t)\\ Jx^2(t)y(t)\end{array}\right)+\left(\begin{array}{c}|y(t)|\eta _1(t)\\ |x(t)|\eta _2(t)\end{array}\right),$$
(15)
where the Gaussian white noise variables $`\eta _i(t)`$ have zero mean and variance $`\eta _i(t)\eta _j(t^{})=2T\delta (tt^{})\delta _{ij}`$. In this representation the kinetic constraint manifests itself as a state dependent rate: the driving term is linear in $`𝒞`$ while the noise term goes as $`\sqrt{𝒞}`$. The unconstrained version of the same model, i.e. one for which $`(𝒞_x,𝒞_y)=(1,1)`$, would have Langevin equation $`\dot{x}^\mu (t)=Jx^\mu +\eta ^\mu (t)`$.

We now show explicitly that low temperature trajectories with fixed initial and final configurations are approximated by the geodesic equation (14). We have used Transition Path Sampling (TPS)  to ensure that we obtained a set of trajectories which both obeyed the endpoint conditions are were well-sampled from the dynamical action (9). In order to generate out-of-equilibrium trajectories for the relaxation from high energy configurations we used the Crooks-Chandler  algorithm; alternatives such as the local algorithm of Ref.  are also suitable.

Figure 1 presents the comparison of the dynamical evolution of our toy KCM with the geometry of its configuration space. TPS was used on the system defined by Eq. (15) to generate trajectories between the configurations indicated by a star. In Fig. 1a we show both a characteristic trajectory (thin red path) and the average over the set of all trajectories obtained via TPS (thick blue path). We also show geodesics (black curves) joining the final configuration to the initial configuration (indicated by a star), and to other configurations along the gradient descent path (indicated by black squares). The solutions to (14) are degenerate: there is one geodesic going left from the initial point and another going right. The distinction depends on the initial conditions on the path. In the TPS simulations we force the initial displacement to be to the left to obtain a set of trajectories which follows the left hand geodesics. These are the paths shown in Fig. 1a and 1b. In Fig. 1b we also show the average over a small set of right hand TPS trajectories and the corresponding pairs of geodesics.

In Fig. 1 we see the features described in the previous section. Initially the system is in the high energy region and it relaxes approximately by gradient descent. As it reaches the region of weak dynamic facilitation at low energies, the system relaxes by following approximately the geodesics of configuration space. Shown for comparison in Fig. 1a are the trajectories followed by an unconstrained version of the same model (dashed green path). In this case relaxation is by gradient descent. Given that the rate for motion in the $`x`$ (resp. $`y`$) direction vanishes when $`y`$ (resp. $`x`$) vanishes, trajectories have to cross the coordinate axis orthogonally, which gives trajectories (and geodesics) a “rectangular” shape. Equipotential lines are concentric circles, and the conflict between allowed dynamical paths and energy relaxation is evident.

## 4 An even simpler example

We can gain more insight into the nature of dynamic facilitation by considering an even simpler, one-variable model. This can be thought of as modeling facilitated dynamics for an order parameter, rather than for microscopic variables. Equations (5) and (6) tell us that an $`N`$-variable, dynamically constrained system looks like a single particle diffusing in $`N`$ space dimensions. Thus a one-variable dynamically constrained model looks like a single particle diffusing in one space dimension. We can exploit the fact that one dimensional spaces are always flat  to ‘remove’ the kinetic constraint. We will see below that this can be done at the expense of introducing an effective free energy barrier.

Consider a model with one degree of freedom, $`x`$, with Hamiltonian (1) $`H=Jx^2/2`$. We choose the dynamical constraint to be $`𝒞(x)`$. This differs from the higher-dimensional cases in that the variable constrains itself. Consequently, the transformation from Fokker-Planck to Langevin form is ambiguous , and depends on whether one choses pre-point or mid-point discretisation of time-dependent quantities (Itô or Stratonovich calculus). We choose the former, so that our model is causal. The Langevin description of our one-variable model is then

$$\dot{x}(t)=Jx𝒞(x)+T𝒞^{}(x)+\sqrt{𝒞(x)}\eta (t),$$
(16)
where the prime denotes differentiation with respect to $`x`$, and the noise term $`\eta (t)`$ has zero mean and variance $`\eta (t)\eta (t^{})=2T\delta (tt^{})`$. Since in this case the constraint depends on the variable it constrains, the second term on the right hand side of (16) is needed to ensure an equilibrium density proportional to $`\mathrm{exp}(Jx^2/2T)`$ . In the higher-dimensional cases we considered above there was no self-facilitation and this term was absent. The corresponding unconstrained model has $`𝒞(x)=1`$, and its Langevin equation is $`\dot{x}(t)=_x(Jx^2/2)+\eta (t)`$.

The appearance of the multiplicative noise term $`\sqrt{𝒞(x)}\eta (t)`$ in Eq. (16) is due to the metric-augmented diffusion term in (6), and is therefore a consequence of the kinetic constraint. The statement that a one-dimensional space is flat is equivalent to the statement that the noise may be made additive by a local change of variables. If we make the change of variables

$$y(t)^{x(t)}\frac{dx^{}}{\sqrt{𝒞(x^{})}},$$
(17)
we can use Itô’s formula  to rewrite (16) as

$$\dot{y}(t)=Jx\sqrt{𝒞(x)}+\frac{1}{2}T\sqrt{𝒞(x)}_x\mathrm{ln}𝒞(x)+\stackrel{~}{\eta }(t),$$
(18)
where $`x`$ should be eliminated in favour of $`y`$. The noise $`\stackrel{~}{\eta }(t)`$ has the same mean and variance as $`\eta (t)`$. Provided we can write the deterministic terms in (18) as the gradient of an effective potential we can regard Equation (18) as unconstrained, in the sense that the prefactor of the noise does not contain the dependent variable.

Consider the specific case of a quadratic constraint, $`𝒞(x)=x^2`$. Equation (16) becomes

$$\dot{x}=Jx^3+2Tx+x\eta ,$$
(19)
where we have assumed $`x>0`$. The change of variables (17) now amounts to a Cole-Hopf transformation, $`y\mathrm{ln}x`$. Hence the regime where dynamic facilitation is relevant corresponds to $`y0`$. Equation (18) now reads

$$\dot{y}=_y\left(\frac{J}{2}e^{2y}Ty\right)+\stackrel{~}{\eta }.$$
(20)
Equation (20) looks like an unconstrained model with an effective entropic barrier $`\mathrm{\Delta }F=Ty=T\mathrm{ln}x`$, which prevents the particle from accessing the dynamically constrained regime of $`y0x0`$ . Figure 2 illustrates this.

## 5 Conclusions

In this paper we have presented a geometric interpretation of dynamic facilitation. By generalizing kinetically constrained models to continuous degrees of freedom we have shown that in these systems the kinetic constraints can be seen as endowing configuration space with a non-trivial metric structure. This metric structure in turn determines the dynamical trajectories. In the regime where the kinetic constraint is dominant, i.e. at low temperatures or energies, the geometry forces the dynamics to proceed through paths close to geodesics, rather than close to paths which relax the energy. This conflict between available paths and energy relaxation gives rise to dynamical bottlenecks or barriers, despite the fact that there are no static barriers in these systems.

Our view of dynamical arrest in terms of the geometry of configuration space is distinct from that based on the idea of motion in an effective rugged energy (or free energy) landscape (see  and references therein). In the landscape perspective it is often assumed that relevant features of the energy surface, like local minima and saddle points, determine the behaviour of the system, and that knowledge of the statistical properties of these features is sufficient to explain both dynamic and thermodynamic properties. In our approach the effective energy surface displays no particular features, the thermodynamics is uninteresting, and there are no special configurations, such as local minima or saddles, which play any major roles either in the statics or in the dynamics. All the interesting structure is in the paths between configurations, and therefore in the dynamics. In this picture it is the effective metric structure of configuration space that distinguishes a glassy system from a normal one, and not the ruggedness of its energy surface.

We are grateful to D. Chandler for discussions. We acknowledge financial support from EPSRC Grants No. GR/R83712/01 and GR/S54074/01, the Glasstone Fund, and Linacre College Oxford.


text:
# Atomistic simulations of the incipient ferroelectric KTaO3

## Abstract

A parameterized effective Hamiltonian approach is used to investigate KTaO3. We find that the experimentally observed anomalous dielectric response of this incipient ferroelectric is well reproduced by this approach, once quantum effects are accounted for. Quantum fluctuations suppress the paraelectric–to–ferroelectric phase transition; it is unnecessary to introduce defects to explain the dielectric behavior. The resulting quantum-induced local structure exhibits off-center atomic displacements that display longitudinal, needle-like correlations extending a few lattice constants

Numerous experimental and theoretical studies have been carried out on the perovskite KTaO3 over the last forty years (see, e.g, Refs Kvyatkovskii ; salce ; Wemple ; singhKTaO3 ; Kevin ; Tinte ; ref7 ; ref8 ; raman ; ref6 ; Ang ; comes  and references therein), making this material one of the most-studied “incipient ferrolectrics.” The main reason for this interest is that the dielectric constant of KTaO3 increases continuously with decreasing temperature down to $``$10K, but then saturates to a plateau at a large value ($``$ 4000) at lower temperatures while remaining paraelectric and cubic all the way down to zero Kelvin salce ; Wemple . These anomalous low-temperature features are usually thought to be caused by the suppression of a paraelectric–to–ferroelectric phase transition by zero-point quantum fluctuations salce ; Wemple  (hence the name “incipient ferroelectric” or “quantum paraelectric” used to describe KTaO3 and other materials, such as SrTiO3, exhibiting similar unusual dielectric and structural properties). Surprisingly, this generally-accepted picture is apparently not supported by various first-principles calculations, using density-functional theory (DFT) either in its local-density approximation (LDA)Kohn  or generalized-gradient approximation (GGA)pw91 ; pbe  form, since these simulations all predict that KTaO3 should be paraelectric at $`T`$=0 even when neglecting zero-point motion singhKTaO3 ; Kevin ; Tinte . This raises the possibility that LDA and GGA are not accurate enough to adequately reproduce the qualitative properties of incipient ferroelectrics. An alternate explanation for this discrepancy between first-principles calculations and experiments is that the simulations assume a perfect material while real samples may contain defects such as oxygen vacancies and Fe+3 ions salce ; Wemple ; ref7 ; ref8  that might lead to the observed anomalous properties of KTaO3. In fact, the interpretations of various experiments raman ; ref6  still remain controversial as to whether they are attributable to extrinsic effects (i.e., defects-induced) or intrinsic off-center atomic displacements. Furthermore, while previous studies invoke the existence of ferroelectric microregions inside the macroscopically-paraelectric KTaO3 system to explain some of its properties raman ; Ang , there has never been any direct determination of the size and shape of these proposed polar regions, to the best of our knowledge. For instance, the pioneering work of Ref. raman  made several assumptions in their analysis of low-temperature Raman spectra – such as isotropy of these microregions – to extract a characteristic size $``$ 16 Å for these polar regions.

In this Letter, we use large-scale atomistic simulations to shed light on the aforementioned long-standing problems. We report calculations on KTaO3 using a parameterized effective Hamiltonian approach. Our main findings are that (i) LDA and GGA are indeed not accurate enough to reproduce the observed anomalous properties of KTaO3, even qualitatively; (ii) these properties can be understood without the need of introducing defects, if quantum fluctuations are present to suppress the paraelectric–to–ferroelectric transition; (iii) the low-temperature local structure of KTaO3 is characterized by off-center atomic displacements that are longitudinally correlated, in a needle-like (and thus anisotropic) way, with a correlation length spanning a few 5-atom unit cells.

We use the effective Hamiltonian ($`H_{\mathrm{eff}}`$) approach developed in Ref. ZhongDavid  to investigate finite-temperature properties of KTaO3. Within this approach, the total energy $`E_{\mathrm{tot}}`$ is a function of three types of local degrees of freedom: (1) the ui (B-site centered) local soft-mode amplitude in each $`i`$ 5-atom cell, describing the local polarization in each cell; (2) the $`v_i`$ (A-site centered) inhomogeneous strain variables; and (3) the homogeneous strain tensor. $`E_{\mathrm{tot}}`$ contains 18 parameters and 5 different contributions: a local-mode self energy, a long-range dipole-dipole interaction, a short-range interaction between local modes, an elastic energy, and an interaction between the local modes and strains ZhongDavid . This effective Hamiltonian approach has been successfully used to model, understand, and design ferroelectric perovskites (see Refs. ZhongDavid ; DavidReview ; LaurentReview ; jorge ; kevin2  and references therein). $`E_{\mathrm{tot}}`$ is used in two different kinds of Monte-Carlo (MC) simulations: classical Monte Carlo (CMC) Metropolis , which does not take into account zero-point phonon vibrations, and path-integral quantum Monte Carlo (PI-QMC) jorge ; zhong ; ceperley  which includes purely quantum-mechanical zero-point motion. Consequently, comparing the results of these two different Monte-Carlo techniques allows a precise determination of quantum effects on macroscopic and microscopic properties of perovskites. 12$`\times `$12$`\times `$12 KTaO3 supercells (corresponding to 8,640 atoms) are used in all Monte-Carlo simulations. We typically perform 30,000 MC sweeps to thermalize the system and 70,000 more to compute averages, except at low temperatures in PI-QMC where more statistics is needed. For example, we use 180,000 MC sweeps for thermalization and 240,000 sweeps at 3K to accurately predict the dielectric response. (Note that we are not aware of any previous work reporting the dielectric response computed using PI-QMC)

In PI-QMC, each 5-atom cell interacts with its images at neighboring imaginary times through a spring-like potential (mimicking the zero-point phonon vibrations), while all the 5-atom cells interact with each other at the same imaginary time through the internal potential associated with $`E_{\mathrm{tot}}`$. The product $`TP`$, where $`T`$ is the simulated temperature and $`P`$ is the number of imaginary time slices (Trotter number), controls the accuracy of the PI-QMC calculation. In all our simulations we use $`TP`$=600, which we find leads to sufficiently converged results. Outputs of the PI-QMC simulations thus contain local-modes ui(t), where $`i`$ indexes the 5-atom unit cells of the studied supercell while the imaginary time $`t`$ ranges between 1 and $`P`$. Note that CMC simulations can be thought of as corresponding to $`P=1`$, so that they do not yield imaginary-time-dependent outputs.

Figure 1(a) shows the $`\chi _{33}`$ dielectric susceptibility – where the index $`3`$ refers to the  pseudo-cubic direction – as predicted by the $`H_{\mathrm{eff}}`$ approach, with all its parameters being derived from LDA calculations on small supercells of KTaO3 at its experimental lattice constant. (Technical details of these LDA calculations are similar to those of Ref. Kevin ). It can be clearly seen that CMC calculations yield a $`\chi _{33}`$ that is continuously increasing as the temperature is decreasing down to nearly zero Kelvin. Turning on quantum effects leads to the appearance of a plateau below $``$100K with a value of $``$100 for the dielectric constant. These CMC and PI-QMC simulations both predict a cubic paraelectric ground state. A plateau for the dielectric response has indeed been experimentally observed in KTaO3 salce ; Wemple , but reaching a much higher dielectric constant ($``$ 4,000) and over a much narrower temperature range (i.e., below 10K) than in Fig. 1(a).

In view of these two discrepancies, we have experimented with making minor adjustments in the LDA-fitted parameters in the hope of obtaining better agreement with experimental data. We have found that this can be done by adjusting just one of the 18 parameters, namely, the parameter denoted $`\kappa _2`$ in Ref.ZhongDavid , which describes the harmonic part of the local-mode self-energy. (In our model, reducing $`\kappa _2`$ favors ferroelectricity with respect to paraelectricity since it leads to a decrease of the zone-center transverse optical frequency by weakening short-range repulsions). Figure 1b shows that decreasing this single $`\kappa _2`$ parameter by $``$18 % from its LDA value of 0.0866 a.u. (atomic units) leads to reasonable agreement between our PI-QMC simulations and measurements, not only for the value of the dielectric constant plateau, but also at temperatures above 10K.

Furthermore, this modified $`\kappa _2`$ also results in a dramatic difference between the two kinds of Monte-Carlo calculations. CMC simulations yield a ferroelectric rhombohedral ground-state. The corresponding Curie temperature is around 30K, as evidenced by the peak in dielectric response displayed in Fig. 1(b). On the other hand, PI-QMC predicts a paraelectric ground state. In other words, quantum effects suppress the paraelectric–to–ferroelectric phase transition, which is consistent with the accepted picture salce ; Wemple . Figures 1(a-b) thus (i) reveal that extrinsic defects (such as impurities or vacancies), which have been proposed to be responsible for the anomalous properties of KTaO3 salce ; Wemple ; ref7 ; ref8 , are not needed to reproduce the experimental behavior of this material; and (ii) strongly suggest that, unlike in strongly ferroelectric perovskitesDavidReview ; LaurentReview , the LDA is not accurate enough for simulating KTaO3.

As for GGA, Tinte et al. Tinte  report zone-center optical frequencies in cubic KTaO3 that are all positive and very close to the LDA values. Consequently, according to Fig 1(a-b), we can conclude that a GGA effective Hamiltonian would not provide a significant improvement over our LDA one, and will also fail in reproducing experimental results. This may make KTaO3 a useful test-case for the development of new functionals within DFT or other ab-initio methods.

We now analyze the microscopic local structure of KTaO3 at low temperature. Figure 2 depicts the magnitude of the local modes ui inside each $`i`$ 5-atom cell versus the angle that these local modes make with the pseudo-cubic  direction, as obtained from a $`T`$=3K snapshot among the thermally equilibrated Monte-Carlo configurations using $`E_{\mathrm{tot}}`$ with the modified $`\kappa _2`$. (The magnitude of the local mode is directly proportional to the magnitude of the local polarization, e.g., $`\left|𝐮\right|=`$ 0.006 and 0.026 a.u. correspond to a local polarization $``$ 0.0583 and 0.253 $`C/m^2`$, respectively). Figure 2(a) displays the CMC results, while Fig. 2(b) corresponds to PI-QMC footnoteJorge . Comparing Figs. 2(a) and 2(b) reveals how quantum effects affect the microscopic structure of KTaO3: the local polarizations go from all lying close to the  direction (corresponding to an angle $``$ $`54^{}`$) and having a relatively large magnitude, to being heavily-scattered in direction and having a much smaller but non-zero magnitude. The fact that KTaO3 is predicted to exhibit non-zero local dipoles, even when quantum fluctuations are accounted for, is consistent with the first-order lines observed to appear in Raman spectra which are forbidden in the ideal cubic perovskite structure raman ; ref6 . Furthermore, an inspection of Fig. 2(b) does not reveal any obvious polar microregions. For instance, our quantum-statistical results do not show the local-mode distributions breaking up into clusters centered along $`111`$ directions (i.e., angles of $``$ $`54`$ and/or $`125^{}`$) as would be expected for such polar microregions.

To gain further insight into the local structure of KTaO3, we decided to compute an additional set of coefficients defined as

$$\theta _\mu (𝐫)=\frac{3}{N}\underset{i=1}{\overset{N}{}}\frac{u_{i,\mu }u_{i+𝐫,\mu }}{\left|𝐮_𝐢\right|\left|𝐮_{i+𝐫}\right|}.$$
(1)
Here $`\mu `$ denotes the $`x`$, $`y`$, or $`z`$ Cartesian axis chosen along the ,  or  cubic directions, respectively. The index $`i`$ runs over all the $`N`$ B-sites; $`u_{i,\mu }`$ and $`u_{i+𝐫,\mu }`$ are the $`\mu `$ components of the local modes in cell $`i`$, and in the cell centered at a distance $`𝐫`$ from cell $`i`$, respectively. The case in which the local dipoles all have the same (non-zero) magnitude and are all aligned along a given $`111`$ direction yields a value of 1 for $`\theta _\mu (𝐫)`$, for any $`𝐫`$ and for any $`\mu `$. This case corresponds to a ferroelectric rhombohedral state having identical local and average structures. On the other hand, the other limiting case — for which neighbors at a distance $`𝐫`$ do not exhibit any correlation between the $`\mu `$-components of their local modes — is associated with a zero value for $`\theta _\mu (𝐫)`$.

Figure 3 depicts $`\theta _x(𝐫)`$ (i.e., $`\mu =x`$) for $`𝐫`$ lying in the $`x`$-$`y`$ plane. The results correspond to one snapshot of a thermally equilibrated Monte Carlo configuration at $`T`$=3K, using the $`H_{\mathrm{eff}}`$ with the modified $`\kappa _2`$. Panels (a) and (b) correspond to CMC and PI-QMC simulations respectively. One can see that CMC technique leads to a $`\theta _x(𝐫)`$ close to unity for any $`𝐫`$, and thus generates a macroscopically- and microscopically-ferroelectric rhombohedral structure, as consistent with Fig. 2(a).

On the other hand, PI-QMC simulations give a more complex behavior for $`\theta _x(𝐫)`$ at low temperature. One can see that the $`x`$ components of the local modes are longitudinally correlated in a needle-like fashion: $`\theta _x(𝐫)`$ adopts large values only when $`𝐫`$ is along the  direction. These values decreasing as the magnitude of $`𝐫`$ increases. (The same result is obtained for all symmetry related cases, e.g., for $`\theta _x(𝐫)`$ in the $`x`$-$`z`$ plane, etc.) Figure 2(b) further reveals that $`\theta _x(𝐫)0.5`$ for neighbors at a distance of $`\pm 2a`$ (where $`a4`$Å is the cubic lattice constant) along the $`x`$ axis. This is in good agreement with the characteristic size of $`16`$Å extracted from low-temperature Raman spectra of KTaO3 raman . On the other hand, our simulations go against the hypothesis of isotropic correlation made in ref. raman . The longitudinal needle-like correlations depicted in Fig. 3(b) have also been predicted to occur in classical ferroelectrics just above the paraelectric–to–ferroelectric transition temperature DavidFE1998 . In fact, they are pretransitional effects that are probably common to most ferroelectric perovskites, the peculiarity of quantum paraelectric KTaO3 being that the phase transition does not actually occur. Finally, note that these needle-like correlations are consistent with the peculiar diffuse X-ray scattering observed in Ref. comes .

In summary, we have performed large-scale atomistic simulations to investigate the (defect-free) incipient ferroelectric KTaO3 system using a parameterized effective-Hamiltonian approach. The effect of quantum-mechanical zero-point motion is investigated by comparing the results of classical and path-integral Monte Carlo simulations. We find that the fitting of all the $`H_{\mathrm{eff}}`$ parameters within LDA yields a theoretical dielectric constant that is in poor quantitative agreement with experiment, strongly suggesting that LDA is inadequate for this material. Results in the literature also indicate that GGA will not improve the LDA result. On the other hand, a small modification of a single parameter in $`H_{\mathrm{eff}}`$ from its LDA value is enough to obtain reasonable agreement between theory and experiment for the dielectric constant over a wide temperature range. This modified $`H_{\mathrm{eff}}`$ leads to the predictions that (i) KTaO3 is ferroelectric classically, but becomes paraelectric once zero-point phonon vibrations are included, and (2) the quantum-induced local structure of KTaO3 is characterized by non-zero local dipoles that have longitudinal, needle-like correlations with a correlation length spanning a few unit cells.

We thank D. Ceperley, B. Dkhil, M. Itoh, J.M. Kiat, W. Kleemann, J. Kohanoff, I. Kornev, S.A. Prosandeev, G. Samara, J. Shumway, and E. Venturini for useful discussions. This work is supported by Office of Naval Research Grants N00014-01-1-0365 (CPD), N00014-01-1-0600 and N00014-97-1-0048, National Science Foundation Grant DMR-9983678 and Department of Energy Contract DE-AC04-94AL85000 at Sandia National Laboratories.


text:
# Low-Metallicity Gas Clouds in a Galaxy Proto-Cluster at Redshift 2.38

## 1 Introduction

There has been considerable recent work on investigating the inter-relationship between gas (traced by QSO absorption lines) and galaxies at high redshifts (eg. Pascarelle et al., 2001; Francis et al., 2001; Williger et al., 2002; Adelberger et al., 2003; Zhan, Fang & Burstein, 2003). A complex picture is emerging, with gas tracing galaxies on larger scales but some classes of absorption-line system avoiding galaxies on small scales.

QSO absorption-line measurements such as these could potentially shed light on the origin of the highly metal enriched intra-cluster medium of galaxy clusters today. Many models have been proposed (eg. Renzini, 1997; Chiosi, 2000; Martinelli, Matteucci and Colafrancesco, 2000; Aguirre et al., 2001; De Lucia, Kauffmann & White, 2003) to produce this gas, but all have problems, and there is currently no consensus.

In this paper, we present high resolution spectra of a QSO that lies behind a probable galaxy proto-cluster at redshift 2.38 (Fig 1). QSO LBQS 2139$``$4434 lies at redshift 3.23. Its sight line passes only 22″ (160 projected proper kpc, for $`H_0=70\mathrm{k}\mathrm{m}\mathrm{s}^1\mathrm{Mpc}^1`$, $`\mathrm{\Omega }_{\mathrm{matter}}=0.3`$ and $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$) from the extremely red object B1. Francis et al. (2001) demonstrated that B1 is actually a pair of large elliptical galaxies at redshift 2.38..

B1 is embedded within a luminous, extended Ly$`\alpha `$ emitting nebula or “blob”. Such Ly$`\alpha `$ nebulae are thought to be associated with massive dark halos, and may be the ancestors of cooling flows today (Steidel et al., 2000; Francis et al., 2001). It lies within a probable galaxy proto-cluster (Francis and Hewett, 1993; Francis et al., 1996; Francis, Woodgate and Danks, 1997; Francis, Wilson and Woodgate, 2001), which is itself embedded in a 80 Mpc scale “Great Wall” (Palunas et al., 2004; Francis et al., 2004).

This sight-line thus allows us to probe the gas in an unusual part of the early universe: the halo of an cluster elliptical galaxy.

## 2 Observations and Reduction

Our existing spectra of LBQS 2139$``$4434 (Francis and Hewett, 1993; Francis, Wilson and Woodgate, 2001) were of too low a resolution and narrow a wavelength coverage to seriously constrain the physical properties of the gas. We therefore re-observed it using the University College London Echelle Spectrograph (UCLES) on the Anglo-Australian Telescope (AAT). Observations were carried out on the nights of 2001 August 20 — 23 and the total usable integration time was 30,600 sec. D’Odorico, Petitjean & Cristiani (2002) independently obtained high resolution spectra of this QSO, using the UVES spectrograph. They kindly allowed us to use their spectrum in this analysis.

The data were reduced using standard procedures, and set to the vacuum heliocentric frame. Absorption-lines properties were measured interactively using the XVOIGT program (Mar and Bailey, 1995). We fit our UCLES spectrum, the UVES spectrum, and a lower resolution AAT spectrum which gave better coverage at UV wavelengths (Francis, Wilson and Woodgate, 2001), using for each line whichever spectrum gave the most reliable data. Where a given line was seen in multiple spectra, the derived parameters were always consistent within the error bars. Parameters we derived from the UVES spectra agree very well with those derived by D’Odorico, Petitjean & Cristiani (2002).

## 3 Results

QSO LBQS 2139$``$4434 was previously known to show high equivalent width Ly$`\alpha `$ absorption at the redshift of the elliptical galaxies making up B1 (and of the cluster). We confirm this: the absorption centroid redshift matches the C IV emission-line centroid of B1 to within $`200\mathrm{k}\mathrm{m}\mathrm{s}^1`$, the uncertainty being dominated by the width of the emission lines.

For the first time, we unequivocally detect metal-line absorption from this absorption-line system (Fig 2). No high ionization absorption is seen: the tentative detection of C IV by Francis, Wilson and Woodgate (2001) turns out, at higher resolution, to be caused by blends of unrelated lines at other redshifts.

The metal-line absorption is clearly resolved into two narrow components, separated in velocity by $`59\mathrm{k}\mathrm{m}\mathrm{s}^1`$: a lower redshift component (component 1) at $`z=2.3797`$ and a higher redshift component (component 2) at $`z=2.3804`$. Both components are narrow and spectrally unresolved: we place upper limits on the Doppler $`b`$ parameter of $`<6.0\mathrm{km}\mathrm{s}^1`$ in all transitions. For the Fe II transitions of component 1 we can place a tighter upper limit of $`b<5\mathrm{k}\mathrm{m}\mathrm{s}^1`$. There is no evidence for velocity substructure within either of these components.

If we assume that the neutral hydrogen is located in one or both of the metal-line clouds, this allows us to break the degeneracy (noted by Francis, Wilson and Woodgate (2001)) between velocity dispersion and column density in our fits to the Ly$`\alpha `$ emission profile, and measure the neutral hydrogen column density. The combined neutral atomic hydrogen column density is $`10^{19.7\pm 0.1}\mathrm{cm}^2`$. We are unable to constrain how it is divided up between the two components.

There is also a strong degeneracy between the assumed velocity dispersion $`b`$ and the inferred metal column densities. For component 1, we can place a lower limit of $`b>2\mathrm{k}\mathrm{m}\mathrm{s}^1`$ from the O I line: if $`b`$ were lower, then to get a good fit we would need such a high column density that we would start to see damping wings, which are not observed. We cannot place an observational lower limit on $`b`$ in component 2, but physically the sound speed of even a cold neutral medium should be at least $`1\mathrm{k}\mathrm{m}\mathrm{s}^1`$, so we take that as our lower limit. In Table 1, we use these limits on $`b`$ to bracket the possible column densities.

We list upper limits on non-detected lines by assuming the same limiting values of $`b`$, and seeing interactively how strong such an absorber could be before the residuals go clearly into emission.

## 4 The Physical State of the Absorbing Gas

The strength of the low ionization species, such as O I and Si II, compared with our upper limits on high ionization species such as C IV and Si IV clearly indicate that the absorbing gas is cool, and at least partially neutral.

There are two possibilities for the temperature of these gas clouds: we could have a warm neutral/partially ionised phase with $`T10^4`$K or a cold neutral phase with $`10^2<T<10^3`$ K (eg. Sutherland & Dopita, 1993).

We modeled the first possibility (warm neutral/partially ionised phase) using Gary Ferland’s CLOUDY code (Ferland, 1996). For plausible choices of the UV ionizing background (Haardt & Madau, 1996; Scott et al., 2000), we can reproduce the absorption-line columns measured with a cloud of density $`0.1<N_H<10\mathrm{c}\mathrm{m}^3`$, size 2 — 200 pc and mass 0.001 – 10 $`M_{\mathrm{}}`$. If the density is towards the lower end of this range, the hydrogen could be up to $`80`$% ionised.

The second possibility is almost impossible to model, due to the complex molecular cooling processes, but cold neutral gas in our own galaxy does show ionization states and column densities comparable to those we observe.

Which of these models is correct? Our stringent upper limits on the velocity dispersions of the two components provide some evidence in favor of the cold model. The sound speed in the warm model is $`10\mathrm{k}\mathrm{m}\mathrm{s}^1`$, while in the cold model it would be only 1 — 3 $`\mathrm{km}\mathrm{s}^1`$. One would expect the internal turbulent velocity dispersions to be comparable to the sound speed. The rapid movement of these clouds through B1’s dark matter halo should be ample to excite and sustain such turbulent motions. Our upper limit of $`B<6\mathrm{k}\mathrm{m}\mathrm{s}^1`$ is thus most consistent with a cool model.

On the other hand, if the velocity dispersion really is low ($`b2\mathrm{k}\mathrm{m}\mathrm{s}^1`$), we are forced to invoke extremely large O I column densities to fit the observed absorption line. This would imply an extraordinary overabundance of oxygen compared to silicon: \[O/Si\] $`>1.5`$ for $`b2`$. Larger values of $`b`$ bring the inferred column density of O I down to more reasonable values. A cold neutral phase with supersonic turbulence would be an interesting possibility.

If the gas were in the cold phase, we might expect dust to be present. We can place upper limits on the quantity of dust by noting, however, that the QSO has a UV continuum slope indistinguishable from other QSOs in the Large Bright QSO Sample (Francis et al., 1991). This allows us to place a rough upper limit on the extinction in these clouds of $`E(BV)<0.1`$.

We also note (Fig 2) that the O I line touches zero flux. Thus the absorbing cloud must be large enough to completely obscure the background QSO continuum emitting region.

In conclusion, we are unable to discriminate between the two physical models. We consider the oxygen overabundance argument somewhat stronger than the velocity dispersion argument, but neither is conclusive. Higher resolution spectra will be needed to resolve this puzzle.

## 5 Conclusions

Our line of sight though the halo of B1, and probably though the central regions of a galaxy proto-cluster, is thus intersecting two small gas clouds. Their properties are very similar to those of high velocity clouds (HVCs) in the halos of our own galaxy (eg. Gibson et al., 2001; Tripp et al., 2003). Their metallicities must be low: $`3.0<`$\[Fe/H\] $`<1.1`$, \[C/H\] $`<1.65`$ and $`2.9<`$ \[Si/H\] $`<1.65`$ (Table 1, and allowing for neutral hydrogen fractions as low as 20%). The warmer the gas, the lower our inferred metallicity. If there is substantial molecular hydrogen, the metallicities would be lower still.

Given the continued controversy over the nature even of the HVCs in our own galactic halo, it would be foolhardy to make dogmatic statements about the origins of these distant clouds. Many of the same models that have been proposed for our local HVCs (eg. tidally disrupted dwarf galaxies, in-falling primordial gas, gas trapped in dark matter sub-halos, debris from galactic fountains) could apply here.

The metallicities we measure for these clouds are lower than the typical intra-cluster medium (ICM) metallicities of galaxy clusters today. This implies either that this enrichment of this cluster had not yet taken place, or that the cold gas clouds we see are in some way isolated from the enriched gas. It is possible, for example, that these clouds are embedded in a very hot, ionized, X-ray emitting, metal-rich ICM. We would not detect such a hot phase in absorption. We do not detect B1 in X-rays (Williger et al. 2004, in preparation), but our current limits do not place tight constraints on the existence of such an ICM.

These observations may also shed light on the enigmatic Ly$`\alpha `$ blobs. While our sight-line does not pass through the blob, it does probe the same dark matter halo. It is quite possible that the Ly$`\alpha `$ emission region is made up of small clouds similar to those we are probing in absorption. The emission could be driven by photoionization and photo-evaporation of these clouds by some concealed AGN, or could be produced by fast shocks created at cloud-cloud collisions, as discussed by Francis et al. (2001).

We would like to thank D’Odorico, Petitjean and Cristiani for allowing us to use their UVES spectrum of QSO 2139$``$4434, and Ralph Sutherland for many illuminating discussions. We would also like to thank all the people who helped us get back to work at Mt Stromlo Observatory so soon after the January 2003 bushfires.


text:
# Resistance of superconducting nanowires connected to normal metal leads

Current address: \]Department of Physics and Astronomy, University of Basel, Klingelbergstrasse 82, CH-4056 Basel, Switzerland

## Abstract

We study experimentally the low temperature resistance of superconducting nanowires connected to normal metal reservoirs. We find that a substantial fraction of the nanowires is resistive, down to the lowest temperature measured, indicative of an intrinsic boundary resistance due to the Andreev-conversion of normal current to supercurrent. The results are successfully analyzed in terms of the kinetic equations for diffusive superconductors.

Superconducting nanowires are resistive down to very low temperatures due to dynamical changes of the macroscopic phase (phase-slip). Thermally activated phase-slip (TAP) becomes more likely for reduced cross-sectional dimensions because the free-energy barrier scales with the area of the wire. Upon approaching $`T0`$, phase slip due to thermal activation disappears and resistivity persists only by macroscopic quantum tunneling through the free-energy barrier. These processes have recently been studied in suspended carbon nanotubes coated with a thin layer of a superconducting molybdenum-germanium (MoGe-)alloy,MTNature2000 ; MTPRL2001  and receive increased theoretical attention.Sachdev2004 ; BuchlerPRL2004  In a separate experimentkociak01  ropes of single-walled carbon nanotubes show signs of superconductivity, which should be strongly influenced by phase-slips as well.

A second potential cause of low temperature resistance in superconducting wires is the penetration of a static electric field at normal-metal–superconductor interfaces. It reflects the conversion of current carried by normal electrons into one carried by Cooper-pairs via Andreev-reflection. It has been studied extensively close to the critical temperature, where it is related to quasiparticle charge imbalance.clarke72 ; hsiangclarke1980  Although hardly experimentally studied, at very low temperatures a similar resistive contribution is expected to be present, reflecting a length of the order of the coherence length. Since the coherence length is a sizable portion of the resistive length of the nanowires it may contribute significantly to the measured two-point resistance. Here we report experimental results on the resistance of narrow superconducting wires connected to normal metal leads(for short NSN). We find a strong contribution to the resistance due to the conversion processes, which is analyzed and understood in terms of the non-equilibrium theory for dirty superconductors.

We chose to study samples (Fig.1) made of superconducting (S) aluminium (Al) because of its long coherence length. Our main interest is in the two-point resistance of the S-wire. Hence, we have chosen to work with thick and wide normal (N) contacts with a negligible contribution to the normal state resistance. To minimize interface resistances due to electronic mismatch of both materials, we have chosen to work for N with bilayers of aluminium covered with thick normal metal (Cu). In such a geometry the superconducting aluminium wire is directly connected to normal aluminium.

The S-wire of 100 nm thick Al is made by evaporating 99.999% purity Al at a rate of $`1`$ nm/s in a vacuum of $`1\times 10^8`$ mbar during evaporation. Films made in the same way, have a residual resistance ratio, $`RRR=R_{300K}/R_{4.2K}`$ of $`7.5`$ indicative of the level of impurity scattering. Taking the phonon resistivity of $`\rho _{ph}(Al)=2.7\mu \mathrm{\Omega }`$cm, the impurity resistivity is $`\rho _0=0.4\mu \mathrm{\Omega }`$cm. Using $`\sigma _0=N(0)e^2D`$ and renormalized free-electron parameters: $`N(0)=2.2\times 10^{47}`$ J-1m-3 and $`v_\mathrm{F}=1.3\times 10^6`$ m/s (Ref. ashcroftmermin, ), we find the elastic mean free path $`\mathrm{}`$ of 100 nm, presumably limited by the thickness. The superconducting transition temperature of the 100 nm film is 1.26 K, the usually enhanced value for aluminium thin films.

The sample is made in one evaporation run using shadow evaporation. A 300 nm PMMA/ 500 nm PMMA-MAA lift-off mask defines the width and the length of the S-wire. The width is approximately 250 nm and the length is varied from 1 $`\mu `$m to 4 $`\mu `$m. The wire and the aluminium of the normal reservoirs is deposited in one run, both 100 nm thick. The reservoirs are made normal by a 2$`^{\text{nd}}`$ evaporation of 470 nm copper(Cu) on top. The Cu is evaporated at a rate of $`2.5`$ nm/s at a pressure of $`1\times 10^7`$ mbar. The material-parameters are $`\rho _0=0.4\mu \mathrm{\Omega }`$cm, $`\mathrm{}=165`$ nm, $`D=66\times 10^3`$ m2/s using $`N(0)=1.5\times 10^{47}`$ J-1m-3, and $`v_\mathrm{F}=1.2\times 10^6`$ m/s. The non-superconducting properties of the final reservoirs are confirmed by measuring the resistance of a 100 nm Al/470 nm Cu bilayer down to the lowest temperature measured: 600 mK. No sign of superconductivity is observed. This is in agreement with the analysis of $`T_c`$ for a bilayer using the Usadel equationsaarts97 . Only for a transparency of the Al/Cu interface lower than 0.1 the $`T_c`$ would reach values above 600 mK. 0.1 is an unrealistically low transparency for an in vacuum prepared Al/Cu interface.

The two-point resistance is probed with a current of 1 $`\mu `$A (linear regime) in a 3He system down to 600 mK. The voltage over the wire is measured with a lock-in amplifier at 133 Hz.

We find that the normal state resistance of nominally identical wires scatters substantially (10% or more), presumably due to grain-sizes compared to wire-widths. To allow a quantitative evaluation we have selected a set of wires with identical values of $`RRR=3.2\pm 0.1`$. Figure 2 shows the $`RT`$ measurements for this particular set of wires. Samples with different $`RRR`$ values show qualitatively identical behavior. All wires show a finite remaining resistance as most striking result.siddiqi00

Obviously, despite the differences in length, the resistances at low temperatures have identical values and follow the same trace. It indicates that the origin of this remaining resistance is likely due to the region in the S-wire next to the interface with the normal reservoir. The resistance at 600 mK is equal to a normal segment of the superconductor with a length of about 200 nm.

Figure 2 (Inset) also shows that the critical temperature of the wire decreases linearly with increasing the inverse square of the wire length.

The critical temperature of the wire should follow a straightforward Ginzburg-Landau analysis. For $`\xi _{\mathrm{GL}}(T)\pi L`$ the bridge should become superconducting. This leads to an $`1/L^2`$ onset-temperature according to:

$$\frac{T_c}{T_{c0}}=12.2\pi \frac{\mathrm{}D}{\mathrm{\Delta }(0)}\frac{1}{L^2}$$
(1)
For $`L\mathrm{}`$ the normal contact can no longer depress $`T_c`$ and we find the intrinsic critical temperature. For the studied wires it is found to be $`T_{c0}=1.26`$ K, identical to the independently determined values of the 100 nm Al film.

From the normal state resistance we infer an impurity resistance of $`\rho _0=1.1\mu \mathrm{\Omega }`$cm, in accordance with the $`RRR=3.2`$. It is significantly higher than the $`\rho _0`$ of the 100 nm Al film, leading to a lower diffusivity and elastic mean free path: $`D=160`$ cm2/s, and $`\mathrm{}=37`$ nm, respectively. The resulting coherence length $`\xi =\sqrt{\frac{\mathrm{}D}{2\pi k_\mathrm{B}T_c}}=124`$ nm.

Finally, we estimate a resistance contribution of 11 m$`\mathrm{\Omega }`$ due to the spreading resistance in the normal reservoirs at low temperatures, considerably less than the value we measure in Fig. 2.

Theoretically, since the studied nanowires show diffusive transport, the Usadel equations should apply to the system. diffusive  It is most convenient to calculate the normal current for a given applied voltage difference, assuming linear response. Schmid and Schönschmid75  have shown that within this limit the normal current can be described with a variation, $`\delta f(E,x)`$, in the electronic distribution function $`f(E,x)`$:

$$I=\frac{A\sigma }{e}_{\mathrm{}}^+\mathrm{}M(E,0)\frac{\delta f(E,0)}{x}𝑑E.$$
(2)
with $`x`$ the coordinate along the length of the wire. $`\delta f(E,x)`$ is determined from a Boltzmann-equation, which includes the conversion of electrons into Cooper-pairs, but ignores inelastic electron-phonon scattering (only relevant close to $`T_c`$):

$$\mathrm{}D\frac{}{x}\left(M\frac{\delta f}{x}\right)2\mathrm{\Delta }N_2\delta f=0.$$
(3)
The position-dependent spectral conductivity $`M(E,x)`$ consist of two parts $`N_1(E,x)=Re(G)`$ and $`N_2(E,x)=Re(F)`$, with $`G`$ the normal and $`F`$ the anomalous Greens function:

$$M(E,x)=\left(N_1(E,x)\right)^2+\left(N_2(E,x)\right)^2,$$
(4)
$`N_1(E,x)`$ is comparable to the standard BCS density of states. The applied voltage $`V`$, is taken into account via the boundary conditions for Eq. 3:

$$\delta f(E;0,L)=\frac{\pm eV/2}{4k_BT\mathrm{cosh}^2(E/2k_BT)}.$$
(5)
Hence, the normal metal leads are taken as equilibrium reservoirs. The strength of the pairing interaction $`\mathrm{\Delta }(x)`$ in Eq. 3 is determined by solving the Usadel equation in the Matsubara representation:

$`{\displaystyle \frac{1}{2}}\mathrm{}D{\displaystyle \frac{^2\theta }{x^2}}`$ $`=`$ $`\mathrm{\Delta }(x)\mathrm{cos}\theta +\omega _n\mathrm{sin}\theta `$
$`\omega _n`$ $`=`$ $`(2n+1)\pi k_BT,n=0,1,\mathrm{}`$ (6)
in conjunction with the self-consistency equation

$$\mathrm{\Delta }\mathrm{ln}\left(\frac{T_c}{T}\right)=2\pi k_BT\underset{n=0}{\overset{\mathrm{}}{}}\left(\frac{\mathrm{\Delta }}{\omega _n}\mathrm{sin}(\theta )\right).$$
(7)
The so-called proximity angle, $`\theta `$, parameterizes the normal Green’s function $`G=\mathrm{cos}\theta `$ and the anomalous Green’s function $`F=\mathrm{sin}\theta `$. The spectral functions are calculated by again solving the Usadel equation but now for $`\omega iE`$. The large normal metal reservoirs impose the boundary condition $`\theta (x=0,L)=0`$

Our main interest is the question how the current conversion process contributes to the resistance. First of all, the presence of decaying normal electron states suppresses the gap in the density of states.

In Fig. 3, the density-of-states $`N_1(E)`$ is given for several positions along the wire of $`L=2`$ $`\mu `$m, $`D=160`$ cm2/s, $`\mathrm{\Delta }_0=192`$ $`\mu `$eV, and $`T/T_c=0.4`$. Clearly, moving away form the normal leads the density of states resembles more and more the well-known BCS density of states. Note however, that a finite sub-gap value remains in the middle ($`x=1\mu `$m) even for very long wires. This is an intrinsic result for any NSN system and it is not due to current-flow, since this result is calculated in thermal equilibrium. (The back-action of the current-flow on the spectral properties can be neglected in the linear response regime).

In Fig. 4, we show the results of a calculation of the voltage as a function of position along the wire for two different temperatures: $`t=0.4`$, and $`t=0.9`$ with $`D=160`$ cm2/s, and $`\mathrm{\Delta }_0=192`$ $`\mu `$eV. At the temperature close to the transition temperature, the electric field penetrates the sample completely and the resistance is close to the normal state value. At low temperatures, the electric field still penetrates the superconductor over a finite length, leaving a middle piece with hardly any voltage drop. localvoltage  The penetration length is of the order of the coherence length. The inset shows the position dependent normal currents and supercurrents illustrating the current conversion processes.

In Fig. 5, a comparison is made between the calculated resistance as a function of temperature and the measurement for a $`L=2`$ $`\mu `$m wire. The calculation is done with parameters $`D=160`$ cm2/s, as determined from the impurity resistivity, and $`\mathrm{\Delta }(0)=1.764k_\mathrm{B}T_c=192`$ $`\mu `$eV with $`T_c=1.26`$ K determined from the length dependence. These parameters have been determined independently. Without any fitting parameter, the agreement between the model (dots) and the experiment (data points: open symbols) is encouragingly good. Only at the lower temperatures the observed resistance is slightly less than the theoretically predicted values.

Apparently, the model is overestimating the remaining resistance below $`T=0.8`$ K . Since there is little freedom left, we have very few options to remedy this discrepancy. The most likely option is that the rigid boundary conditions imposed at the interfaces should be relaxed. There is a finite possibility for superconducting correlations to extend into the normal metal reservoirs, which would mean that the boundary condition $`\theta (x=0,L)=0`$ is too rigid. Since the correlations will extend into a 3-dimensional volume we assume that using the boundary conditions $`_x\theta (x=0,L)=\theta (x=0,L)/a`$ *i.e.* a decay over a fixed characteristic length $`a`$, is a realistic assumption. It assumes a geometric dilution of the correlations. The result is shown in Fig. 5 by filled triangles. The best agreement between measurement and calculation model is obtained for $`a=75`$ nm. This value appears reasonable for a decay length since it is comparable to the dimensions (100 nm by 250 nm) of the wire which emits into the reservoirs. The shortest bridge shown in Fig. 2 is found to have a further reduction in resistance, which we attribute to a small misalignment as shown in Fig. 1. Note however that the model predicts that NSN devices will continue to be resistive down to $`T=0`$ K. For the rigid boundary conditions we find $`R/R_{4.2K}=0.3255`$ and for the soft boundary conditions we find $`R/R_{4.2K}=0.2609`$

An early indication of the low temperature boundary reported here is given by Harding *et al.*harding74 . They studied the resistance of thick sandwiches of Pb-Cu-Pb. By subtracting the contribution to the resistance of the Cu-layers they identified a remaining boundary resistance which depended on the mean free path in the Pb layers. In later work by Hsiang and Clarkehsiangclarke1980  such a resistance appeared to be unobservable, in contrast however to more recent work by Gu *et al.* PrattPRB2002  The geometry of our sample, with a negligible contribution to the resistance from the reservoirs, allows us to measure only the resistance due to the conversion processes in the superconductor. In response to the work of Harding *et al.* Krähenbühl and Watts-Tobin krahenbuhl79  derived an analytical expression for the effective length of the boundary resistance: $`x_0=\sqrt{\frac{1}{6}\pi \xi _0l}f(T)`$ with $`\mathrm{}`$ the mean free path for elastic scattering, $`\xi _0`$ the BCS coherence length and $`f`$ a function of temperature. For $`T=0`$ K, the function $`f`$ is of order 1. In contrast to our model Ref. krahenbuhl79,  allows for one-dimensional diffusion in N, rather than treating N as an equilibrium reservoir.

In conclusion, we have shown that the resistance of a superconducting nano-wire connected to normal leads has a finite DC resistance down to very low temperatures. The microscopic theory describes the data very well and provide a detailed image of the conversion from normal current into supercurrent, over a few coherence lengths. The results emphasize the important role of the length of the wires in relation to the nature of the contacts.Sachdev2004  It also explains results obtained with diffusion-cooled hot-electron bolometers in which normal leads are used to provide rapid out-diffusion of hot electrons from a superconducting wiresiddiqi00 .

We thank the Stichting voor Fundamenteel Onderzoek der Materie(FOM) for financial support. We thank R.S. Keizer and A.A. Golubov for stimulating and clarifying discussions.