text:
# 1 Introduction

## 1 Introduction

The dilaton $`\varphi `$ and the axion $`\chi `$ are scalar fields predicted by superstring theory . Both of them arise in a natural way in the massless spectrum of 10 dimensional (10d) type IIB superstring theory and its lower dimensional compactifications. In the language of 4d gauge theory of field strength $`F_{\mu \nu }`$ and its dual $`\stackrel{~}{F}_{\mu \nu }`$, the $`\varphi `$ and $`\chi `$ have very special couplings. The dilaton $`\varphi `$ couples to gauge fields through a term $`exp(\varphi )F^2`$ and the axion $`\chi `$ couples to the topological term. The $`\varphi `$ and $`\chi `$ fields play a central role in superstring dualities , F-theory compactifications  and in the derivation of the exact results in 4d N=2 supersymmetric gauge theories .
Recently it was observed in  that a string inspired coupling of a dilaton $`\varphi `$ to 4d $`SU(N_c)`$ gauge fields $`A_\mu =T^aA_\mu ^a`$, with $`T^a`$ the $`(N_c^21)`$ $`SU(N_c)`$ generators, yields a phenomenologically interesting potential V(r) for the quark-quark interactions. Following , this potential is obtained as follows: First we start from the following model for the scalar field-gluon coupling

$$L(\varphi ,A)=\frac{1}{4G(\varphi )}F_{\mu \nu }^aF_a^{\mu \nu }\frac{1}{2}(_\mu \varphi )^2+W(\varphi )+J_\mu ^aA_a^\mu .$$
(1)
Then we choose $`G(\varphi )`$, the coupling of the scalar field $`\varphi `$ to the $`SU(N_c)`$ field strength $`F_{\mu \nu }`$, and the interacting lagrangian $`W(\varphi )`$ as:

$$G(\varphi )=const.+\frac{f^2}{\varphi ^2}$$

$$W(\varphi )=\frac{1}{2}m^2\varphi ^2$$
(2)
where the parameter $`f`$ is a scale characterizing the strength of the scalar-gluon coupling and m is the mass of the scalar field $`\varphi `$. Next we consider the equations of motion of the fields $`A_\mu `$ and $`\varphi `$ and solve them for static point like color source of current density $`J_a^\mu =\rho _a\eta ^{\mu 0}`$. After some straightforward algebra, we find that the Dick quark interaction potential $`V_D(r)`$ is given by:

$$V_D(r)=\frac{1}{r}f\sqrt{\frac{N_c}{2(N_c1)}}\mathrm{ln}[exp(2mr)1]$$
(3)
Eq.(3) is very remarkable since for large values of $`r`$ it leads to a confining potential $`V_D(r)2fm\sqrt{\frac{N_c}{2(N_c1)}}r`$. In this regard, we will show in this paper that for a general gluon-dilaton coupling $`G(\varphi )`$, the quark interactions potential $`V(r)`$ reads as:

$$V(r)=𝑑r\frac{G[\varphi (r)]}{r^2}.$$
(4)
Such form of the potential is very attractive. On one hand it extends the usual Coulomb formula $`V_c1/r`$ which is recovered from eq.(4) by taking $`G=1`$. Moreover for $`Gr^2`$, which by the way corresponds to a coupling $`G(\varphi )\varphi ^2`$, and $`W(\varphi )=\frac{m^2}{2}\varphi ^2`$, $`m0`$, eq.(4) yields a linearly increasing interquark potential $`Vr`$ having the good behaviour to describe the $`SU(N_c)`$ quarks confinement . On the other hand eq.(4) may be also used to describe other non perturbative effects associated with higher dimension quark and gluon vacuum condensates. Following , see also , one may extract interesting phenomenological informations on the dilaton-gluon coupling $`G[\varphi ]`$ by comparing eq.(4) to the Bian-Huang-Shen’s potential $`V_{BHS}(r)`$ namely:

$$V_{BHS}(r)\frac{1}{r}\underset{n0}{}C_nr^n$$
(5)
where $`C_n`$’s are related to the quark and gluon vacuum condensates. In fact one can do better if one can put the coupling $`G(\varphi )`$ in the form $`G[\varphi (r)]`$. In this case one can predict the type of vacuum condensates of the $`SU(N_c)`$ gauge theory which contributes to the quark-quark interaction potential . Thus, although, the derivation of the formula (4) for the interquark potential from eq.(1) is by itself an important result, there remain however other steps to overcome before one can exploit eq.(4). As mentioned above, a crucial step is to determine what type of couplings $`G(\varphi )`$ which can be put in the form $`G[\varphi (r)]`$. In other words for what couplings $`G(\varphi )`$, one can solve the equation of motion of the scalar field $`\varphi `$. This is a technical problem without solving it one cannot get V(r). An other step which remains to clarify is to show how the effective model (1) may be got from a more fundamental theory. If this is possible, one may for instance justify the mass scale $`f`$ introduced by hand in eqs.(2,3). One might also get some informations on the axion field couplings and more generally on the moduli of 10d superstrings compactified on six dimensional compact manifolds and especially type IIB on Calabi-Yau Threefolds . In trying to explore eq.(4), we have observed some remarkable facts among which we quote the three following :
1) The functional $`G[\varphi (r)]`$, and then the potential V(r) of eq.(4) may be obtained from the following one dimensional lagrangian $`L_D`$

$$L_D=\frac{1}{2}(y^{})^2+r^2W(y/r)+\frac{\alpha }{2r^2}G(y/r)$$
(6)
where $`y=r\varphi `$, $`y^{}=(\frac{dy}{dr})`$and $`\alpha =\frac{g^2}{16\pi ^2}\frac{N_c1}{2N}`$ and where $`g`$ is the gluon coupling constant. In particular $`L_D`$ reads, for $`W(\varphi )`$ and $`G(\varphi )`$ like in eq(2),as:

$$2L_D=(y^{})^2+m^2y^2+\frac{\mu ^2}{y^2}$$
(7)
where $`\mu =\alpha f^2`$.
2) Eq.(7) has a striking resemblance with the following harmonic superspace lagrangian $`L_{EH}`$ used in  in the derivation of the 4d Eguchi-Hanson metric

$$2L_{EH}=(D^{++}\omega )^2+m_{}^{++}{}_{}{}^{2}\omega ^2+\frac{\mu _{}^{++}{}_{}{}^{2}}{\omega ^2}.$$
(8)
In this equation, $`\omega `$ is an analytic harmonic superspace (hs) superfield taken to be dimensionless, $`D^{++}`$ is the hs covariant derivative and $`m^{++}`$ and $`\mu ^{++}`$ are coupling constants. More details on hs tools will be described in section 3. Much more precision can be found in . For the moment note only the formal analogy between y, dy/dr, m and $`\mu `$ of eq.(7) with $`\omega `$, $`(D^{++}\omega )`$, $`\mu ^{++}`$ and $`m^{++}`$ respectively. Both of models eqs.(7) and (8), involve hermitian fields with a self interacting potential proportional to the inverse of the square of the scalar field variable.
3) The Dick potential (3) is viable only for non zero mass dilaton field exactly as in 4d N=2 supersymmetric theories where the scalar potential is proportional to the mass eigenvalues of the central charges of the 4d N=2 superalgebra . Recall by the way that in 4d N=2 supersymmetric QFT, mass terms are generated by central charges. We shall see in sections 3 and 4 that this formal analogy between the Dick model (1) and 4d N=2 QFT’s is much deeper since it allows us to derive a new model containing eq.(1) and where the symmetries behind the solvability of Dick equations as well as the couplings of both the dilaton and axion fields are manifest.

The aim of this paper is to generalize the Dick model (1) by exploiting the formal analogy with 4d N=2 supersymmetric theories formulated in hs  and using known 4d N=2 exact results. In addition to the derivation of new model exhibiting a U(1) gauge invariance, we give an interpretation of the mass scale $`f`$, introduced by hand in eq.(2), as a Kahler moduli of a blown up SU(2) singularity of Calabi-Yau threefold of type II superstring compactifications. The appearance of the local U(1) symmetry in the analysis of eqs.(1-3) has a quite interesting consequence as it offers a possibility to incorporate in the game the axion field $`\chi `$ couplings. Recall that in Dick model as formulated in , the role of the topological field $`\chi `$ is ignored. We shall show in section 4 how this field can be incorporated by going in a general gauge other than $`\varphi =\varphi ^{}`$.

The presentation of this paper is as follows: In section 2, we formulate the Dick problem as a one dimensional field theory. In section 3 we give general solutions including those of . In section 4, we review briefly the building of the Eguchi Hanson hyperkahler metric in harmonic superspace. In section 5 we use the formal analogy between the Eguchi Hanson model and our one dimensional field theoretical formulation of the Dick problem to determine the dilaton couplings, the axion ones and interprete the mass scale $`f`$ as a kind of Fayet-Iliopoulos coupling. Our conclusion is given in section 6.

## 2 The Dick model as a one dimensional field theory.

Following ref., the analysis of the Coulomb problem of the theory (1) is based on considering a point like static color source which in its rest frame is described by a current $`J_a^\mu =g\delta (r)C_a\eta _0^\mu `$ where $`C_a`$ is the expectation value of the $`SU(N_c)`$ generator for a normalized spinor in the color space. These $`C_a`$’s satisfy the algebraic identity

$$\underset{a=1}{\overset{N^21}{}}C_a^2=\frac{(N_c1)}{2N_c}$$
(9)
The next step is to use the residual SO(3) space symmetry, which remains after setting $`J_a^\mu =\rho _a\eta _0^\mu `$, to rewrite the equations of motion

$$[D_\mu ,G^1(\varphi )F^{\mu \nu }]=J^\nu $$

$$_\mu ^\mu \varphi =\frac{W}{\varphi }\frac{1}{4}F_{\mu \nu }^aF_a^{\mu \nu }\frac{G^1(\varphi )}{\varphi }$$
(10)
into a simple form. Indeed setting $`F_a^{0i}=\frac{gC_a}{4\pi }_iV`$; $`\alpha =\frac{g^2}{16\pi }\frac{(N_c1)}{2N_c}`$ one finds after some easy algebra:

$$\frac{dV}{dr}=r^2G[\varphi ](a)$$

(11)
$$\mathrm{\Delta }\varphi =\frac{W}{\varphi }+\frac{\alpha }{r^4}\frac{G(\varphi )}{\varphi }(b)$$

Note that eqs.(11) have four unknown field quantities; the field $`\varphi `$, the interacting color potential V(r), the dilaton-gluon coupling $`G(\varphi )`$ and the $`\varphi `$ potential $`W`$. To solve eqs.(11) one has to fix two of them. For example choosing $`2W=m\varphi ^2`$and $`G(\varphi )`$ as in equation (3), one finds :

$$\varphi =\varphi _D(r)=r^1[\frac{\alpha f}{m}(1exp(2mr)]^{1/2}$$
(12)
$$V_D(r)=\frac{1}{r}f\sqrt{\frac{N_c}{2(N_c1)}}ln[exp(2mr)1].$$
(13)
In general given $`G(\varphi )`$ and $`W(\varphi )`$, the color potential V(r) can be exactly determined up on solving one equation namely equation (11-b). For later use, let us introduce the new dimensionless field $`y=r\varphi `$ and take the spherical coordinate frame $`(r,\theta ,\phi )`$ to rewrite the lagrangian (1) as:

$$L=\frac{r^2}{2G(\varphi )}F_{0r}^aF_a^{0r}\frac{r^2}{2}_r\varphi ^r\varphi +r^2W(\varphi )+F_{0r}^a\rho _a.$$
(14)
In deriving eq.(14), we have used the stationarity of the color source, the SO(3) symmetry and the identity $`\mathrm{\Delta }(1/r)=\delta (r)`$. Putting this equality back into eq.(14) and using the change of variable $`y=r\varphi `$ together with the convention notations $`y^{}=_ry`$; $`^r\varphi =r^2_r\varphi `$ as well as eq.(9), one gets the lagrangian form eq.(6). Consequently the coupling $`G(\varphi )`$ of eq.(1) appears as a part of interacting potential of the one dimensional field theory eq.(6). From this point of view, the finding of the interquark potential V(r) is equivalent to solve the equation of motion

$$y^{\prime \prime }\frac{L_D}{y^{}}+y^{}\frac{L_D}{y}+_r^{exp}L_D=0.$$
(15)
## 3 Solving the Dick model

First of all observe that the lagrangian (6) including the Dick model (7) is a particular one dimensional field theory of lagrangian

$$L=\frac{1}{2}(y^{})^2U(y,r)$$
(16)
where $`U(y,r)`$ is a priori an arbitrary potential. Though simple, this theory is not easy to solve except in some special cases. A class of solvable models is given by potentials of the form :

$$U(y)=\lambda ^2y^{2(n+p)}+\gamma ^2y^{2(qn)}+\delta y^k$$
(17)
where $`n,p,q`$ and $`k`$ are numbers and $`\lambda ^2`$ , $`\gamma ^2`$ and $`\delta `$ are coupling constants scaling as $`(lenght)^2`$. The next thing to note is that eq.(17) has no explicit dependence in $`r`$ and consequently the following identity usually hold :

$$y^2=U+c$$
(18)
where c is a constant. Actually eq.(18) is just an integral of motion which may be solved under some assumptions. Indeed by making appropriate choices of the coupling $`\lambda `$ as well as the integral constant c , one may linearise y’ in eq.(18) as follows :

$$y^{}=U_1+U_2.$$
(19)
Once the linearisation in y’ is achieved and the terms $`U_1`$ and $`U_2`$ are identified, we can show that the solutions of eq.(18) are classified by the product $`U_1U_2`$ and the ratio $`U_1/U_2`$. In what follows we discuss briefly some interesting examples. For convenience let us rewrite eq.(18) as:

$$y^2=w_0+w_1+C_0$$
(20)
where $`w_0=m^2y^2`$ and $`w_1`$ is the interaction term which we take for the moment to be the Dick interaction that is $`w_1=c_1^2y^2`$, where $`c_1`$ is a coupling constant. Starting from eq.(20), it is not difficult to see that there are two possibilities to put it in the form (19):

### 3.1 First possibility: the Dick solution

This corresponds to take $`w_0=U_1^2`$ that is $`U_1=my`$ and $`U_2=c_1y^1`$. Putting back into eq.(19) one gets the Dick solution given by eqs.(12,13).

### 3.2 Second possibility: New solutions

In this case the mass term is related to the product $`U_1U_2`$ as :

$$U_1U_2=\pm \frac{1}{2}m^2y^2$$
(21)
Eq.(21) cannot however determine $`U_1`$ and $`U_2`$ independently as in general the following realizations are all of them candidates,

$$U_1=\lambda y^{n+p};U_2=\gamma y^{qn}$$
(22)
where the integers $`p`$ and $`q`$ are such that $`p+q=2`$ and where $`\lambda \gamma =\pm m^2`$. A remarkable example corresponds to take $`p+q=1`$. In this case we distinguish two solutions according to the sign of the product of $`\lambda \gamma `$ . For $`\lambda \gamma =+m^2`$, the solution is

$$y(r)=[\frac{1}{\lambda }tan(\frac{nmr}{\sqrt{2}}+const)]^{\frac{1}{n}}.$$
(23)
For $`\lambda \gamma =m^2`$ , we have:

$$y(r)=[\frac{1}{\lambda }tanh(\frac{nmr}{\sqrt{2}}+const)]^{\frac{1}{n}}.$$
(24)
The solutions (23) and (24) have quite interesting features inherited essentially from the features of tan and tanh functions. We remark that for n=0 the solution is:

$$y(r)=const.exp(\frac{\lambda +\gamma }{\sqrt{2}}r).$$
(25)
In the end of this section, it should be noted that one can go beyond the above mentioned solutions which are just special case of general models involving interactions classified by the following constraint equations

$$U_1.U_2y^k$$
(26)
where $`U_1`$ and $`U_2`$ are as in eq.(19) and $`k`$ is an integer. For $`k=0`$, one gets the Dick model and for $`k=2`$ one has solutions described in subsection 3-2 . For general values of $`k`$, one has to know moreover the ratio $`U_1/U_2`$ in order to work out solutions. For the example where

$$U_1=\lambda y$$

$$U_2=\gamma y^{k1};kinteger$$
(27)
one can check, after some straightforward algebra, that the solution of y is just a generalization of eq.(12) that is

$$y_k(r)=[r\varphi _D]^{\frac{2}{(2k)}}.$$
(28)
For $`k=0`$, one discovers the solution (12).

## 4 The Eguchi Hanson hs model

To start recall that Eguchi Hanson metric is a vacuum solution of the self dual euclidean four dimensional gravity. It is a Ricci flat hyperkahler metric having an $`SU(2)\times U(1)`$ isometry. There are different, but equivalent, ways of writing this metric. A remarkable way of expressing this metric is that using a local coordinate system exhibiting manifestly the $`SU(2)\times U(1)`$ symmetry. The element of length $`ds^2`$ reads as:

$$ds^2=g_{iajb}df_1^{ia}df_1^{jb}+k_{iajb}df_2^{ia}df_2^{jb}+h_{iajb}df_2^{ia}df_1^{jb}$$
(29)
where the metric factors are given by:

$$g_{iajb}=ϵ_{ab}ϵ_{ij}\frac{4f_{2ia}f_{2jb}}{f_1^{kc}f_{1kc}+f_2^{kc}f_{2kc}}(a)$$

$$k_{iajb}=ϵ_{ab}ϵ_{ij}\frac{4f_{1ia}f_{1jb}}{f_1^{kc}f_{1kc}+f_2^{kc}f_{2kc}}(b)$$
(30)
$$h_{iajb}=\frac{4f_{1ia}f_{2jb}}{f_1^{kc}f_{1kc}+f_2^{kc}f_{2kc}}c)$$

together with the SU(2) isovector constraint

$$ϵ_{ab}(f_1^{ia}f_2^{jb}+f_1^{ja}f_2^{ib})\lambda ^{ij}=0.$$
(31)
A tricky way to derive this metric is to use results of 4d N=2 supersymmetric non linear $`\sigma `$ models. In the harmonic superspace approach where 4d N=2 supersymmetry is manifest, the field theoretical model giving the family of Eguchi Hanson metrics reads in the superfield language as :

$$S[\omega ]=\frac{1}{2k^2}𝑑z^{(4)}𝑑u[(D^{++}\omega )^2m^{++2}\omega ^2\frac{\lambda ^{++2}}{\omega ^2}]$$
(32)
In this equation $`\omega =\omega (x_A,\theta ^+,\overline{\theta }^+,u)`$, is an analytic hs superfield taken to be dimensionless. $`D^{++}=(u^{+i}\frac{}{u^i}2\theta ^+\sigma ^m\overline{\theta }^+_m)`$ is the hs covariant derivative; $`dz^4`$ is the analytic superspace measure with U(1) Cartan charge (-4) and the couplings $`m^{++}`$ and $`\lambda ^{++}`$ are given by

$$m^{++}=u_i^+u_j^+m^{ij};\lambda ^{++}=u_i^+u_j^+\lambda ^{ij}$$
(33)
where $`u_i^+`$ and $`u_i^{}`$ are the harmonic variables parameterizing the $`SU(2)/U(1)S^2`$ sphere. We shall not use here after these hs tools, we are only interested in the formal analogy with the Dick problem. This is why we shall give here after only the necessary material. For more details on the HS method and the derivation of the Eguchi Hanson metric see . Note also that the Eguchi Hanson metric with $`SU(2)\times U(1)`$ isometry corresponds to $`m^{++}=0`$. Metrics with $`m^{++}0`$ have a $`U(1)\times U(1)`$ symmetry and fall in the family of multicenter metrics . Let us take $`m^{++}=0`$ and sketch the main steps in putting eq.(32) in the form (29-31). In fact there are two possible paths one may follow: First, a direct method which consists to start from the superfield equation of motion of the hermitian hs superfield $`\omega `$:

$$D^{++2}\omega =\frac{\lambda ^{++2}}{\omega ^3}$$
(34)
and use the $`\theta `$expansion of the superfield $`\omega `$, that is

$$\omega =\varphi +\theta ^{+2}M^{(2)}+\overline{\theta }^{+2}\overline{N}^{(2)}+\theta ^+^m\overline{\theta }^+B_m^{(2)}+\theta ^{+2}\overline{\theta }^{+2}P^{(4)};$$
(35)
where we have ignored fermions. Then fix N =2 supersymmetry partially on shell by eliminating the auxiliary fields $`P^{(4)}`$ and $`B_m^{(2)}`$. The relevant equations are those corresponding to the projection of eq.(34) along the $`\theta ^+=0`$ and $`\theta ^+\sigma ^m\overline{\theta }^+`$ directions, i.e.

$$^{++2}\varphi =\lambda ^{++2}/\varphi ^3$$

$$^{++}B_m^2=2(_m\frac{3}{2}\frac{\lambda ^{++2}}{\varphi ^4}B_m^{(2)})\varphi .$$
(36)
The next thing to do is to find the explicit dependence of $`\varphi `$ and $`B_m^{(2)}`$ in harmonic variables $`u_i^\pm `$ by solving eqs.(36). Then put the solution into eq.(32) once the integrations with respect to $`\theta ^+`$ and $`\overline{\theta }^+`$ are performed. In other words put the solutions $`\varphi =\varphi (u_i^\pm )`$ ,$`B_m^{(2)}=B_m^{(2)}(u_i^\pm )`$ into the following component field action:

$$S[\omega ]\frac{1}{k^2}𝑑x^4𝑑u[^{++}B_m^{(2)}^m\varphi +^mB_m^{(2)}^{++}\varphi ].$$
(37)
The last step is to integrate with respect to harmonic variables. Once this is done, we get the bosonic part of the 4d N=2 supersymmetric non linear $`\sigma `$ model from which one can read the Eguchi Hanson metric in the $`\omega `$ representation. The second method, which interest us here, is indirect but it has the merit of being based on hs superfield theory exhibiting manifestly the $`SU(2)\times U(1)`$ symmetry. The main steps of this approach are as follows:

1. Instead of working with a real superfield $`\omega `$, we take a complex superfield $`\omega `$ : $`\overline{\omega }\omega `$ .

2. Modify the action eq.(32) as:

$$S[\omega ]\frac{1}{2k^2}𝑑z^{(4)}𝑑u[|(D^{++}+iV^{++})\omega |^2+\lambda ^{++}V^{++}],$$
(38)
where $`V^{++}`$ is a U(1) gauge superfield. Eq.(38) is invariant under the following gauge transformations of parameter $`\lambda `$.

$$\omega ^{}=\mathrm{exp}(i\lambda )\omega ;V^{++}=V^{++}+D^{++}\lambda .$$
(39)
Note that $`V^{++}`$ has no kinetic term. It is an auxiliary superfield which can be eliminated through its equation of motion namely

$$2V^{++}=\frac{1}{\omega \overline{\omega }}[i(\overline{\omega }D^{++}\omega \omega D^{++}\overline{\omega })\lambda ^{++}].$$
(40)
For the special case where $`\omega `$ is real; $`\overline{\omega }=\omega `$ , eq.(40) reduces to

$$V^{++}=\lambda ^{++}/\omega ^2$$
(41)
and consequently the action (38) coincides with eq.(32). Note by the way that the term $`\lambda ^{++}V^{++}`$ is a Fayet-Iliopoulos (FI) coupling.

3. Rewrite eq.(38) in an equivalent form by using O(2) notations i.e. express the complex supefield $`\omega =\omega _1+i\omega _2`$ as an O(2) doublet $`(\omega _1,\omega _2)`$ and introducing two other auxiliary superfield $`F_1^{++}`$ and $`F_2^{++}`$

$$S[\omega _1,\omega _2,F_1^{++},F_2^{++}]=\frac{1}{2k^2}dz^{(4)}du[(F_1^{++})^2+2F_1^{++}D^{++}\omega _1+(12)$$
(42)
$$V^{++}(\omega _1F_2^{++}\omega _2F_1^{++}+\lambda ^{++})]$$

Eliminating $`V^{++}`$,$`F_1^{++}`$ and $`F_2^{++}`$ and choosing the gauge $`\omega _2=0`$ we reproduce the action $`S_{EH}`$ (32) with $`m^{++}=0`$. The second order action(42) is interesting since it has a manifest SU(2) invariance rotating $`\omega _i`$ and $`F_i^{++}`$ . To make this invariance more explicit we make the following change for both $`(\omega _1,F_1^{++})`$ and $`(\omega _2,F_2^{++})`$.

$$\omega =U_a^{}q^{+a};F^{++}=U_a^+q^{+a}$$

$$q^{+a}=ϵ^{ab}q_b^+;q_a^+=(q^+,\overline{q}^+);ϵ^{12}=1.$$
(43)
Thus for both $`\omega _1`$ and$`\omega _2`$, we have $`\omega =U_a^{}q_I^{+a}`$ with $`I=1,2`$ and so on. Putting back into eq.(42), we get the following action,

$$S=\frac{1}{2k^2}dz^{(4)}du[\overline{q}_1^+D^{++}q_1^++\overline{q}_2^+D^{++}q_2^+++V^{++}(\overline{q}_1^+q_2^++\overline{q}_2^+q_1^++\lambda ^{++})].$$
(44)
This action has the invariance under the following groups:

i) O(2) gauge group acting as

$$\delta q_I^+=ϵ_{IJ}\lambda q_J^+;\delta V^{++}=D^{++}\lambda ;\overline{\lambda }=\lambda .$$
(45)
ii) U(1) subgroup of the rigid SU(2) automorphism group of supersymmetry that leaves $`\lambda ^{++}`$ invariant.

iii) The SU(2) Pauli Cursey symmetry rotating $`q_I^+`$ and $`ϵ_{IJ}\overline{q}_J^+`$.
Now starting from the last form of the E.H. action (42) and solving the $`\theta ^+=0`$ and $`\theta ^+\sigma ^m\overline{\theta }^+`$components of the equations of motion:

$$D^{++}q_Iϵ_{IJ}V^{++}q_J^+=0$$

$$ϵ^{IJ}\overline{q}_I^+q_J^++\lambda ^{++}=0$$
(46)
in the Wess-Zumino gauge, one gets by following the same lines described for the direct method, the E.H. metric (29-31).

## 5 The Dick model revisited

In section 2 we have learnt that the Dick problem may be formulated as a one dimensional field theory of lagrangian $`L_D`$ given by eq.(16) namely:

$$2L_D=(\frac{dy}{dr})^2m^2y^2\frac{\mu ^2}{y^2}$$
(47)
In section 3 we have showed that hyperkhaler metrics of the Eguchi Hanson family can be derived from the following 4d N=2 supersymmetric model,

$$2L_{EH}=(D^{++}\omega )^2m_{}^{++}{}_{}{}^{2}\omega ^2\frac{\mu _{}^{++}{}_{}{}^{2}}{\omega ^2}$$
(48)
In section 4 we have seen that this lagrangian is equivalent to the following first order one, once the auxiliary U(1) gauge superfield $`V^{++}`$ is eliminated through its equation of motion:

$$2L_{EH}^{}=|D^{++}\omega |^2m_{}^{++}{}_{}{}^{2}\omega \overline{\omega }V^{++}(\overline{\omega }D^{++}\omega \omega D^{++}\overline{\omega }\mu ^{++})+V_{}^{++}{}_{}{}^{2}\omega \overline{\omega }$$
(49)
This form of $`L_{EH}`$ may be also transformed into two other forms as shown on eq.(42) and eq.(44). The difference between $`L_{EH}`$ and $`L_{EH}^{}`$ is that in eq.(48) $`\omega `$ is hermitian whereas in eq.(49) $`\omega `$ is complex. As we have seen, we can go from eq.(49) to eq.(48) either by constraining the superfield to be real, that is :

$$\omega _2=0,$$
(50)
Or equivalently by keeping $`\omega _20`$ and working in the Wess-Zumino gauge:

$$D^{++}V^{++}=0$$
(51)
Eq.(51) turns on to be helpful in the derivation of the Eguchi Hanson metric. Now using the formal analogy between eqs.(47) and (48), it is not difficult to see that the $`L_D`$ lagrangian may be also formulated in terms of the auxiliary fields $`F`$ and $`\overline{F}`$ as:

$$L_2=F\overline{F}+\overline{F}y+F\overline{y}+V[\xi +i(y\overline{F}\overline{y}F)]$$
(52)
where V is a one dimensional U(1) gauge field and $`\xi `$ is a 1d constant vector breaking explicitly invariance under space translations. Note that in this formulation, the two scalars $`y_1`$ and $`y_2`$ of the complex field $`y=\frac{1}{\sqrt{2}}(y_1+iy_2)`$ represent respectively the dilaton $`\varphi `$ and the axion $`\chi `$ in agreement with the requirement of F-theory and 10d type IIB superstring and 4d N supersymmetric gauge theory. Eliminating the auxiliary fields F and through their equations of motion namely

$$F=(+iV)y=y(a)$$

(53)
$$\overline{F}=(+iV)\overline{y}=\overline{y}(a)$$

one obtains the following first order lagrangian $`L_1`$

$$L_1=|y|^2+m^2y\overline{y}+\xi V.$$
(54)
Moreover eliminating the auxiliary U(1) gauge field V through its equation of motion

$$V=\frac{1}{2y\overline{y}}[\xi +i(\overline{y}y+y\overline{y})],$$
(55)
one gets the one dimensional field theory of the dilaton-axion system extending eq.(7) which may be recovered from eqs.(54,55) by going to the gauge fixing $`y_2=\chi =0`$. However to exhibit the effect of the axion field $`\chi `$, one has to keep $`\chi 0`$ and imposes a constraint on the gauge field $`V`$ that we write as follows:

$$C(V,V)=0$$
(56)
Using this constraint, the first order lagrangian is no longer invariant under the change $`ye^{i\psi }y`$ and $`VV+\psi `$, where $`\psi `$ is the U(1) gauge parameter, but as a counterpart one can work out a non trivial solution for the axion field $`\chi =\chi (r)`$ by solving the conjugate where the field V should be substituted in equations of motion (55) $`^2y=m^2y`$ and its complex by the value $`V_0(r)`$ verifying the constraint eq.(56) and satisfying the identity $`(\overline{y}yV_0)=0`$. In the end of this study we would like to note that the term $`\xi V`$ appearing in eq.(52) plays a similar role as the Fayet Iliopoulos term $`m^{++}V^{++}`$ of the 4d N=2 supersymmetric Eguchi Hanson model (32) . Thus the mass scale $`f`$ introduced by hand in the Dick model may be viewed, under some assumptions, as the scale of breaking of the U(1) symmetry rotating the dilaton and axion fields. Recall by the way that in general supersymmetric gauge theories with a U(1) gauge invariance, the FI term is generally used to break supersymmetry and/or gauge invariance. The FI couplings are Kahler moduli of Calabi-Yau threefolds on which 10d type II superstrings are compactified on, and their magnitude are of order of the Calabi-Yau compactification scale.

## 6 Conclusion

Inspired from the dilaton-gluon coupling in superstring theory, Dick built a field theoretical model having the remarkable property of leading to a confining quark-quark interaction potential. The model is mainly a 4d $`SU(N_c)`$ gauge theory coupled to a massive scalar field $`\varphi `$ of lagrangian (1) and a dilaton-gluon coupling $`G(\varphi )=1+f^2/\varphi ^2`$, where $`f`$ is a mass scale introduced by hand. The parameter $`f`$ may be compared with the mass scale of the $`\sigma `$ model of mesonic theory . The confining phase of Dick model is parameterized by non zero mass for the dilaton and non vanishing $`f`$. In trying to analyze the potential $`V_D(r)`$ we have observed that the Dick problem has a perfect formal analogy with the problem of building the Eguchi-Hanson metric in 4d N=2 supersymmetric harmonic superspace. This formal similarity appears at several levels. In the introduction we have quoted some of these striking analogies. For example, vanishing masses for both Dick and Eguchi-Hanson scalar fields lead to trivial potentials. An other example is that the mass scale $`f`$ which is introduced by hand and interpreted as a compactification scale by Dick plays a similar role as the 4d N=2 FI coupling appearing in the Eguchi Hanson model (32). Recall by the way that now it is well established that the FI couplings are of order of the compactification scale since they are just the Kahler moduli of the Calabi Yau threefold on which the 10d type II superstring are compactified on. To understand the thing behind the striking similarity between the Dick problem and the Eguchi Hanson one, we have reformulated the Dick problem as a one dimensional field theory. As a consequence we have found a general formula for the inter-quark potential V(r) which of course depend on the nature of the dilaton-gluon $`G[\varphi ]`$ as shown on eq.(4). The beauty of this formula is not only because it extends the Coulomb and Dick theory but also because it can be compared with known parameterizations of the confinement especially the contribution of the quark and gluon vacuum condensates. From this point of view, the Dick model as we have formulated it, may be viewed as a phenomenological theory modeling the non-perturbative contributions responsible of confinement. In this regard a more explicit analysis will be presented in .
Moreover, having at hand the 1d field theoretical formulation of the Dick model and the analogy with 4d N=2 Eguchi Hanson model, we have shown how the axion field may be incorporated in the game in agreement with the requirement of F-theory and 10d type IIB superstrings according to which the dilaton and the axion form a complex field. In our formulation the dilaton-axion model is represented by a 1d U(1) gauge theory of lagrangian (52). The $`U(1)SO(2)`$ symmetry rotates the dilaton and the axion fields and allows to interpret the Dick mass scale as a kind of FI coupling.
This research work has been supported by the program PARS number 372-98 CNR.

References

1. \] M. Green, J. Schwartz and E. Witten, Superstring Theory , Cambrige University Press (1987).
2. \] C. Vafa, Lectures on String and Dualities HUTP-79/A009, hep-th/9702201.
3. \] A. Sen, Unification of string dualities , hep-th/9609176.
4. \] C. Vafa, Evidence for F-theory , HUTP-96/A004, hep-th/9602022.
5. \] N. Seiberg and E. Witten , hep-th/9602022; Nucl. Phys. B426(1994)19 and Nucl. Phys. B431(1994) 484.
6. \] R. Dick, ”Confinment from a massive scalar in QCD”, LMU-TPW-98/06, to appear in Eur. Phys. Journal. C6 (1999). Proceedings of the Miniworkshop on ”Recent Developements in Quantum Field Theory”, UFR-PHE, Rabat (October 1998) to appear.
7. \] R. Dick, ”Implication of dilaton in gauge theory and cosmology”, hep-th/9701047, LMU-TPW-97/10. Phys. lett. B397 (1997)193; Phys. lett. B409 (1997)321.
8. \] T. Huang, Z.Huang, Phys. Rev. D39 (1989)1213.
9. \] J.G. Bian, T. Huang and Q.X.Shen, The quarkonium potential within short and intermediate range in QCD , BIHEP-TH-92-38 and references therein.
10. \] S. Katz, P. Mayr and C. Vafa, Mirror summetry and exact solution of 4d N=2 gauge theories I , HUTP-97/A025, OSU-M.97-5, IASSNS-HEP 97/65, hep-th/9706110.
11. \] A. Galperin, E. Ivanov, V. Ogievetsky and P.K. Townsend, Class. Quantum .Grav. 3 (1986) 625.
12. \] A. Galperin, E. Ivanov, S. Kalitzin, V. Ogievetsky and E. Sokatchev, Class. Quantum . Grav. 1 (1985) 469.
13. \] J. Wess and J. Bagger, Supersymmetry and supergravity , Princeton Press (1983). S. Ferrara and C. Savoy, Supergravity 81 (Cambridge U.P; Cambridge 1982).
14. \] E.A. Gibbons and P.J. Ruback, Commun. Math. Phys. 115 (1988) 267.
15. \] G.W. Gibbons, D. Olivier, P.J. Ruback and G. Valent, Nucl . Phys. B296(1988)679.
E. Sahraoui and E.H. Saidi, Class. Quantum. Grav. 16 (1999) 1.
16. \] K. Wehrberger, Phys. Rep. 5-6 (1993) 273.
17. \] M. Chabab, R. Markazi and E.H. Saidi, (in progress).
18. \]M.A. Shifman, A.I. Vainshtein and V.I. Zakharov, Nucl. Phys. B147 (1979) 385.


format: text

text:
# Quasiparticle Localization Transition in Dirty Superconductors

## I Introduction

The problem of Anderson localization in normal disordered electronic systems has been studied for years, and has continued to pose intriguing puzzles . Recently, the question of whether an analogous transition can occur in the behaviour of low-energy quasiparticle excitations about a superconducting ground state has come to light (see for e.g. Ref..) In fact, the stability of two very distinct superconducting phases which may be characterized by their transport properties - one with extended quasiparticle states at the Fermi energy, and the other with localized states- has been established. Here, it is the transition point between these two phases, and its critical properties, that forms the subject of our attention.

While one can afford to ask of the superconducting systems the same questions that have been addressed in normal systems, the former displays refreshingly new physics that shows both conceptual and qualitative differences from the latter. To begin with, quasiparticle excitations in superconducting systems do not conserve charge, and thus one cannot study the transition through charge dynamics. This absence of $`U(1)`$ symmetry in the Bogoliubov deGennes Hamiltonian, which we employ to describe the excitations in question, has marked consequences. As emphasized by Ref., it gives rise to a host of new universality classes ,with critical exponents whose values are significantly different from those of their normal partners. In particular, the density of states(DOS) exhibits astonishing features. In this paper, we explore these surprises in the context of superconducting systems that respect spin-rotational ($`SU(2)`$) invariance and time-reversal (T) symmetry.

The issue of a delocalization-localization transition within the superconducting state has been a recurrent theme (see for e.g. Ref). But it is only in recent years that analyses of both microscopic and continuum models have paid focused attention from a variety of different avenues, and delved into the prospect of making this transition realizable in physical systems. To repeat the example offered of the dirty d-wave superconductor, consider a system of impure superconducting sheets with d-wave pairing coupled to one another. At the nodal points, one has low-energy quasiparticle excitations, and in fact, due to disorder, one even expects a finite DOS at the Fermi energy. For low interplane coupling strength or high impurity concentration, one would expect these states to be localized, and upon increasing the coupling or lowering the disorder, one could conceive of accessing a critical point beyond which these states become extended. Here, we study such a transition on more generic grounds by exploiting the field-theoretic set-up offered by Ref., and through numerical analyses, both of which serve to bring out the novel features of superconducting systems quite dramatically.

We review the framework used to describe quasiparticle excitations in the absence of interactions. We expand on some details of the two phases, the ’thermal metal’ and the ’thermal insulator’. Then, in order to keep our paper self-contained, we elaborate on the field-theoretic and numerical methods. We proceed to discuss the field-theoretic and numerical predictions for the localization length exponent $`\nu `$ associated with the thermal metal-thermal insulator transition, and show that they both confirm the existence of a new universality class with $`\nu <\nu _n`$, where $`\nu _n`$ is the corresponding exponent for normal systems. We then study the unique properties of the DOS at criticality. Finally, we mention the characteristics of transitions in superconducting systems besides those that respect $`SU(2)`$ and T symmetry, and possibilities for experiment.

## II Models and Symmetries

### A The Superconducting Hamiltonian

Quasiparticle excitations about the superconducting ground state are well described within the framework of the Bogoliubov deGennes(BdG) Hamiltonian. For the spin-singlet paired superconductor, which we focus on here, the BdG Hamiltonian has the form

$`_0=H_1+H_2,`$ (1)
$`H_1`$ (2)
$`=`$ $`{\displaystyle d^dxc_\sigma ^{}(x)\left(\frac{\left(i\mathrm{}\stackrel{}{}\frac{e}{c}\stackrel{}{A}(x)\right)^2}{2m}E_F+V(x)\right)c_\sigma (x)},`$ (3)
$`H_2`$ (4)
$`=`$ $`{\displaystyle d^dxd^dx^{}(c_{}^{}(x)\mathrm{\Delta }(x,x^{})c_{}^{}(x^{})+c_{}(x)\mathrm{\Delta }^{}(x,x^{})c_{}(x^{}))},`$ (5)
where $`c^{}`$ and $`c`$ are electron creation and annhilation operators respectively, $`m`$ the mass , $`E_F`$ the Fermi Energy , $`V(x)`$ a random potential describing impurities in the system, and $`A(x)`$ a vector potential describing any external magnetic field. The lattice version of $`_0`$, which is more tractable for numerics, has the form

$$_{0L}=\underset{ij}{}[t_{ij}\underset{\alpha }{}c_{i\alpha }^{}c_{j\alpha }+(\mathrm{\Delta }_{ij}c_i^{}c_j^{}+h.c)].$$
(6)
From Hermiticity, one requires the condition $`t_{ij}=t_{ji}^{}`$, while spin rotation invariance requires $`\mathrm{\Delta }_{ij}=\mathrm{\Delta }_{ji}`$. The BdG Hamiltonian contains anamolous terms that reflect the fact that the excitations do not conserve charge, and hence break the associated $`U(1)`$ symmetry. As a result, the BdG Hamiltonian lives in an ’extended particle-hole space’ with twice the degrees of freedom of normal electronic systems. More importantly, the absence of $`U(1)`$ charge conservation endows the BdG Hamiltonian with symmetry properties which are completely different from normal systems. It is this difference that plays the key role in giving rise to properties that are unique to superconducting systems.

The symmetries of the BdG Hamiltonian have recently been studied in the context of mesoscopic systems and random matrix theory. Within the class of models described in Ref, the Hamiltonians that we study in this paper have SU(2) and T symmetry; they describe systems that are singlet-paired and have no spin-orbit scattering, thus having spin-rotational invariance, and in addition, they have time-reversal invariance.

Since the quasiparticles conserve both spin and energy, one can explicitly write Eq.6 in terms of conserved quantities by defining a new set of fermionic $`d`$operators:

$$d_i=c_i,d_i=c_i^{}.$$
(7)
Thus one makes a particle-hole transformation on the spin-down operators, but leaves the spin-up operators unchanged. The Hamiltonian of Eq.6 now takes the form

$$_L=\underset{ij}{}d_i^{}\left(\begin{array}{cc}t_{ij}& \mathrm{\Delta }_{ij}\\ \mathrm{\Delta }_{ij}^{}& t_{ij}^{}\end{array}\right)d_j\underset{ij}{}d_i^{}H_{ij}d_j,$$
(8)
where time reversal invariance requires $`t_{ij}`$ and $`\mathrm{\Delta }_{ij}`$ be real. Spin rotational invariance now requires

$$\sigma ^yH_{ij}\sigma ^y=H_{ij}^{},$$
(9)
where $`\sigma ^y`$ is the standard Pauli matrix. Spin conservation along the z-direction is evident from the fact that the physical spin is related to the number operator for the ’d-particles’ by

$$S_i^z=\frac{\mathrm{}}{2}\left(d_i^{}d_i1\right),$$
(10)
and that $`_L`$ of Eq.8 conserves particle number.

The Hamiltonian $`_L`$, in principle, may be diagonalized by solving the eigenvalue equation

$$\underset{i}{}H_{ij}\left[\begin{array}{c}u(j)\\ v(j)\end{array}\right]=E\left[\begin{array}{c}u(i)\\ v(i)\end{array}\right],$$
(11)
where $`u`$ and $`v`$ describe the wave-function amplitudes at each site. Given Eq.11, one can construct the following state

$$i\sigma ^y\left[\begin{array}{c}u(j)^{}\\ v(j)^{}\end{array}\right]=\left[\begin{array}{c}v^{}(i)\\ u^{}(i)\end{array}\right],$$
(12)
with eigenvalue $`E`$; $`SU(2)`$ invariance requires a symmetric dispersion about the Fermi energy, and eigenvalues of the BdG Hamiltonian come in pairs $`(E,E)`$.

We now draw attention to the physical situation described by the wavefunctions of Eq.11, and the relation they bear to the phase transition of interest. To begin with, in gapless superconductors, which we focus on here, one finds states at and about the Fermi energy. These states have a profound impact on thermodynamic and transport properties, as discussed previously. Moreover, the transport properties depend crucially on the nature of the eigenstates of Eq.11 at $`E=0`$, which in turn is determined by the spatial configuration of the random potential $`V(x)`$ and the gap-function $`\mathrm{\Delta }(x)`$. In particular, depending on whether these eigenstates are extended or localized, one can conceive of two very different superconducting phases - the ’thermal metal’ which is capable of transporting energy, and the ’thermal insulator’ which cannot conduct energy over large length scales. In fact, it has been shown that both phases are stable in three dimensions. Thus, one can characterize the two phases by the thermal conductivity $`\kappa `$. In systems with SU(2) symmetry, since the same quasiparticles that carry energy also carry spin, the spin conductivity $`\sigma _s`$, may also be used to describe the two phases, and is related to $`\kappa `$ via an analog of the Weidemann-Franz law in the thermal metal:

$$\frac{\kappa }{T\sigma _s}=const.$$
(13)
Note however, that as the quasiparticles do not carry well defined charge, the Weidemann-Franz law breaks down with regards to electrical conductivity.

The two distinct superconducting phases at hand are in analogy with, but quite different from, the metallic and insulating phases in normal systems. We study the critical properties of the transition between the two phases using the formalism and physical set-up described below.

### B A Field-theoretic Formulation

As with the case of Anderson localization, the problem of quasiparticle localization within the superconducting state can be described within a field-theoretic framework. To briefly review its key features, in previous work, starting with the Hamiltonian of Eq.6, T.Senthil.et.al derived a useful field theoretic action from which one can extract a rich variety of properties of superconducting systems. Coupling the fermionic degrees of freedom in Eq.6 to an infinitesimal Zeeman field $`\eta `$ (which acts as a chemical potential for the d-quasiparticles), and employing the replica method to calculate disorder averages, they obtained an effective action. Fluctuations about the saddle-point of this action were captured near two spatial dimensions by a non-linear sigma model($`NL\sigma M`$) treatment, yielding for the final form of the action,

$$S_{NL\sigma M}=d^dx[\frac{1}{2t}Tr(\stackrel{}{U}.\stackrel{}{U^{}})\eta Tr(U+U^{})].$$
(14)
Here, ’t’ is a dimensionless coupling constant that has the physical interpretation of inverse spin conductance, i.e. $`\frac{1}{t}=\frac{\pi }{2}\sigma _s`$. $`U(x)`$ is a matrix with symplectic $`Sp(2n)`$ group structure, where ’n’ is the number of replicas.

The action $`S_{NL\sigma M}`$ given above is quite different from the analogous action describing normal systems. The field theory , referred to as the ’the principle chiral $`Sp(2n)`$ model’, has its first term invariant under the global ’rotation’ $`UA^{}UB`$, where $`A^{},BSp(2n)`$ , thus possessesing $`Sp(2n)\times Sp(2n)`$ symmetry. The second term in the action reduces the symmetry to an $`Sp(2n)`$ symmetry as it only allows invariance under $`UA^{}UA`$. A knowledge of this symmetry structure proves to be very useful in deriving critical properties.

As done for electrical conductivity in normal systems, we employ a scaling theory for the inverse spin conductivity $`t`$, and analyze the critical point seperating the thermal metal and the thermal insulator. Specifically, we extract the localization length exponent and the unusual singular behaviour of the DOS at the Fermi energy.

### C Hamiltonian for Numerics

The predictive power of the effective field theory of the previous section lies in the fact that some results derived from it are universal to all Hamiltonians satisfying the appropriate symmetries. In practice, we find that in obtaining the localization length exponent, an analog of the tight-binding Anderson model commonly used for numerics shows crisp data with much less noise than other models that we have studied. We focus on the same model for the density of states since the localisation length numerics enables us to identify the critical point quite accurately. With reference to Eq.6, the couplings take the form

$$t_{ij}=\{\begin{array}{cc}\hfill \frac{1}{\sqrt{2}},& ij,n.n.\hfill \\ \hfill V_{it},& i=j\hfill \\ \hfill 0,& otherwise\hfill \end{array}$$
(15)
$$\mathrm{\Delta }_{ij}=\{\begin{array}{cc}\hfill \frac{1}{\sqrt{2}},& ij,n.n.\hfill \\ \hfill V_{i\mathrm{\Delta }},& i=j\hfill \\ \hfill 0,& otherwise\hfill \end{array}$$
(16)
where n.n denotes nearest neighbours, and $`V_{it}`$ and $`V_{i\mathrm{\Delta }}`$ are on-site random variables chosen from a uniform probability distribution ranging from $`W`$ to $`+W`$. We work with a three-dimensional cubic lattice described by Eq.6, and having real couplings with the specific form

$$_L=\frac{1}{\sqrt{2}}[t\sigma ^z+\mathrm{\Delta }\sigma ^x].$$
(17)
Here, the $`\sigma `$’s denote Pauli matrices, and $`t`$ and $`\mathrm{\Delta }`$ are matrices with off-diagonal terms set to unity, and diagonal terms taking on values $`\sqrt{2}V_{it}`$ and $`\sqrt{2}V_{i\mathrm{\Delta }}`$.

### D The Transfer Matrix

While we obtain the DOS from the Hamiltonian of Eq.17 by the straightforward process of exact diagonalization, we extract the localization length exponent by means of a transfer matrix formulation that caters to Eq.17 and tremendously reduces the dimensions of the matrices involved in numerical work.

To describe its principle, consider a quasi 1-dimensional strip in d-dimensions with cross-sectional area $`L^{d1}`$, and in-going and out-going states at either end of this long strip as shown in Fig.1. One might formally obtain the scattering matrix S for the in-going and out-going states using the definition

$$S=\underset{T\mathrm{}}{lim}\mathrm{exp}\frac{iHT}{\mathrm{}}.$$
(18)
The scattering matrix, in the specific basis of in-going and out-going states can be written in terms of reflection and transmission matrices ’r’ and ’t’ respectively:

$$S=\left(\begin{array}{cc}r& t^{}\\ t& r^{}\end{array}\right).$$
(19)
Given this form of the scattering matrix, one can easily derive the transfer matrix T which has the property

$$\stackrel{}{\psi _L}=T\stackrel{}{\psi _R}.$$
(20)
The symmetries of the BdG Hamiltonian imply that the transfer matrix too has very specific symmetry properties which distinguish it from those of normal systems. Of late, these symmetries have been explored on group theoretic grounds.

The transfer matrix T can be constructed by multiplying a set of transfer matrices, each appropriate for a slice of the strip shown in Fig.1:

$$\stackrel{}{\psi }_{i+1,R}=T_i\stackrel{}{\psi }_{i,R},$$
(21)
$$T=\underset{i=1}{\overset{N}{}}T_i.$$
(22)
The form of the transfer matrix that we use for numerical calculation does not make the symmetry of the BdG Hamiltonian manifest, but it is tailored specifically for a tight-binding Hamiltonian such as the one described in Eq.6. We begin by writing the Schrodinger Eq.11 as a difference equation

$$A_i\stackrel{}{D}_i+B_{i,i+1}\stackrel{}{D}_{i+1}+B_{i,i1}\stackrel{}{D}_{i1}=E\stackrel{}{D}_i,$$
(23)
where $`\stackrel{}{D}_i`$ denotes the wave-function amplitudes on each slice in the spinful d-quasiparticle eigenbasis, and $`A_i`$ and $`B_{i,i+1}`$ are $`2L^2\times 2L^2`$ matrices of the form

$$A_i=\left[\begin{array}{cc}t_{ii}& \mathrm{\Delta }_{ii}\\ \mathrm{\Delta }_{ii}& t_{ii}\end{array}\right],$$
(24)
$$B_{i,i+1}=\left[\begin{array}{cc}t_{i,i+1}& \mathrm{\Delta }_{i,i+1}\\ \mathrm{\Delta }_{i,i+1}& t_{i,i+1}\end{array}\right],$$
(25)
where the $`t_{ij}`$ and $`\mathrm{\Delta }_{ij}`$’s are now matrices coupling slices ’i’ and ’j’ in the manner described by Eq.23. To obtain the transfer matrix, we rewrite Eq.23 as

$$\stackrel{}{D}_{i+1}=B_{i,i+1}^1(EA_i)\stackrel{}{D}_iB_{i,i+1}^1B_{i1,i}\stackrel{}{D}_{i1}.$$
(26)
For our model, the inter-slice coupling has the simple form

$$B_{i,i+1}=\frac{1}{\sqrt{2}}[I\sigma ^x+I\sigma ^z],$$
(27)
which satisfies the special property

$$B_{i,i+1}=B_{i,i+1}^1.$$
(28)
Using Eq.26, we are now in a position to define a transfer matrix as follows:

$$\left(\begin{array}{c}\stackrel{}{D}_{i+1}\\ \stackrel{}{D}_i\end{array}\right)=T_i\left(\begin{array}{c}\stackrel{}{D}_i\\ \stackrel{}{D}_{i1}\end{array}\right),$$
(29)
where $`T_i`$ has the relatively simple form

$$T_i=\left(\begin{array}{cc}B_{i,i+1}(EA_i)& I\\ I& 0\end{array}\right),$$
(30)
and the multiplicative property of Eq.22. With this transfer matrix at hand, which very closely resembles standard ones used for the Anderson model, we can extract the localization length for different values of energy E and disorder W in a manner completely analogous to the numerical treatment of the Anderson model. As we are interested in the behaviour of states at the Fermi energy, we set the energy E, to zero.

The procedure for extracting the localization length is quite standard, and has been elaborated on in great depth in many works. But to briefly outline the method, one begins with an orthonormal basis of vectors $`\widehat{o_i}(0)`$ in the space of the transfer matrix, representing the right-most states in Fig.1. One then assumes that for a given disorder strength W, and width L, there exists a set of eigenmodes $`\widehat{w_i}(L,W)`$ that describes typical eigenmodes for a quasi 1-dimensional system of N slices, where ’N’ is large enough to represent the average behaviour of the random disorder, and that these modes decay or grow in magnitude as $`\mathrm{exp}(\pm \gamma _iN)`$ upon multiplication with the corresponding transfer matrix $`T=_{i=1}^NT_i`$. One further assumes that when the set of basis vectors $`\widehat{o_i}(0)`$ , each of which may be represented as a linear combination of the $`\widehat{w_i}(L,W)`$,is multiplied by T, the resulting vectors $`\stackrel{}{v_i}(1)`$ (where the dependence on L and W is implicit), are each composed of the appropriate linear combination of modes $`\widehat{w_i}(L,W)`$ now weighted by the corresponding growth factors $`\mathrm{exp}(\pm \gamma _iN)`$. The length ’N’ is numerically restricted by the exponentially growing magnitude of the vectors $`\stackrel{}{v}`$.

Since all vectors have now grown fastest along the fastest growing mode, say $`\widehat{w_1}(L,W)`$, the magnitude and direction of any one of these modes, say $`\stackrel{}{v_1}(1)`$, are $`\mathrm{exp}\gamma _1N`$ and along $`\widehat{w_1}(L,W)`$ respectively. Projecting out the component $`\stackrel{}{v_1}(1)`$ from the next vector $`\stackrel{}{v_2}(1)`$ gives a resultant vector $`\stackrel{}{o_2}(1)`$ whose magnitude is roughly $`\mathrm{exp}\gamma _2N`$. Thus, by such a ’Gram-Schmidt’ orthogonalization procedure for the whole set of vectors $`\stackrel{}{v_i}(1)`$, one obtains a set of orthogonal vectors $`\stackrel{}{o_i}(1)`$ with associated Lyapunov exponents

$$\gamma _i(1)=\frac{\mathrm{ln}|\stackrel{}{o_i}(1)|}{N},$$
(31)
which give the characteristic inverse localization lengths associated with each mode. To reduce computational effort, we consider just the positive Lyapunov exponents corresponding to exponentially growing states. As we are restricted in our length size N, to obtain a fair estimate of the typical $`\gamma _i`$’s, we repeat the procedure of transfer matrix multiplication, now using as our initial basis vectors the normalized set $`\widehat{o_i}(1)`$ which more or less point along the ideal basis $`\widehat{w_i}(L,W)`$. An average value $`\frac{1}{M}_{j=1}^M\gamma _i(j)`$ obtained from M such iterations provides the desired estimate of the ideal $`\gamma _i`$’s. We associate the characteristic localization length $`\lambda (L,W)`$ with the slowest decaying mode, and thus with the inverse of the smallest positive Lyapunov exponent $`\gamma _{min}(L,W)`$.

In the quasi 1-dimensional case, all modes are exponentially decaying or growing (corresponding to in-going or out-going states respectively) since even the slightest disorder is enough to localize states. But in the 3-dimensional limit, where the 2-dimensional cross-sectional area becomes large, we know that the modes ought to experience a transition from extended to localized behaviour as a function of disorder. To determine the critical disorder strength $`W_c`$ for this transition, the nature of the modes in 3-dimensions and the localization length for an infinite size system of given disorder strength, one can use a finite size scaling analysis of the quasi 1-dimensional system. The scaling function that we will use to do so is the dimensionless parameter

$$\mathrm{\Lambda }(L,W)=\frac{\lambda (L,W)}{L}.$$
(32)
We now turn to the critical properties of the phase transition between the thermal insulator and thermal metal.

## III Critical Properties

As in disordered electronic systems, we have seen in Sec.II B that the thermal metal, the thermal insulator ,and the critical point seprating them may be characterized by their transport properties. The replica field theory of Eq.14, with its dimensionless coupling t, provides a powerful means of studying this transition. An analysis of the action in Eq.14 shows that in $`2+ϵ`$ dimensions, where $`ϵ=1`$ for our system, an unstable fixed point $`t_c`$ describes the critical point between the thermal metal and thermal insulator. One can study the scaling behaviour of t with system size L(see Fig.2) by deriving a form for the scaling function $`\beta (t)=\frac{dt}{d\mathrm{ln}L}`$, which we present explicitly in the next section. Near $`t_c`$, for $`t>t_c`$, the coupling $`t`$ grows larger with L, and thus exhibits a stable thermal insulator, while for $`t<t_c`$, a smaller and smaller value of $`t`$ with increasing length scale signals a thermal metal.

### A Localization Length Exponent

The $`\beta `$ function for the $`Sp(2n)`$ sigma model of the action of Eq.14 can be found in Ref., and it is given to cubic order(’two’ loop) in coupling t by

$$\beta (t;Sp(2n))=ϵt+(2n+1)+\frac{1}{2}(2n+1)^2t^3+𝒪(t^4),$$
(33)
where $`ϵ=d2`$, and n denotes the number of replicas. In the limit $`n0`$, at the critical point where the $`\beta `$ function vanishes, its derivative gives the inverse localization length exponent:

$$\frac{1}{\nu }=ϵ+\frac{ϵ^2}{2}+𝒪(ϵ^3).$$
(34)
In contrast, normal systems with time reversal symmetry and spin-rotational invariance may be identified with the $`Sp(4n)/Sp(2n)\times Sp(2n)`$ model of Ref., and its associated $`\beta _n`$ function has the form

$$\beta _n=ϵt+(4n+1)t^2+(8n^2+2n)t^3+𝒪(t^4).$$
(35)
The localization length exponent derived from Eq.35 has the value

$$\frac{1}{\nu _n}=ϵ+𝒪(ϵ^3),$$
(36)
quite different from the value of $`\nu `$ in Eq.34. The numerical evidence to follow supports the field theoretic result that the value of the localization length exponent $`\nu `$ for the superconducting system with T and $`SU(2)`$ is considerably lower than the analogous exponent $`\nu _n`$ for normal systems.

#### 1 Numerical Treatment

The standard numerical technique that we employ for extracting the localization length exponent $`\nu `$, shows that in three dimensions it takes on the value $`1.15\pm 0.15`$.

The finite-size scaling technique can be summarized as follows: scaling arguments require that the only relevant length scale in the system be the localization length $`\xi (W)`$ of the infinite sized system, and thus we have

$$\frac{\lambda (L,W)}{L}=\mathrm{\Lambda }(L,W)=h(\frac{\xi (W)}{L}),$$
(37)
where h is a scaling function yet to be determined. Close to the critical point $`W_c`$, we have the localization length $`\xi `$ of the infinite system behaving as

$$\xi |WW_c|^\nu ,$$
(38)
which means that the argument ’$`x`$’ of $`h(x)`$ in Eq.37 blows up at the critical point. However, $`\mathrm{\Lambda }`$ is well-behaved and finite; this is only possible if we have

$$\underset{x\mathrm{}}{lim}h(x)=const.,$$
(39)
where the constant refers to independence with respect to L for large L. Thus, the critical value $`\mathrm{\Lambda }_c`$ is common to all sufficiently large system sizes. Using Eq.38, we rewrite Eq.37 as

$$\mathrm{ln}\mathrm{\Lambda }(L,W)=f[L^{\frac{1}{\nu }}(WW_c)].$$
(40)
Linearizing the function ’f’ about the critical fixed point $`(W_c,\mathrm{ln}\mathrm{\Lambda }_c)`$ yields

$$\mathrm{ln}\mathrm{\Lambda }(L,W)=\mathrm{ln}\mathrm{\Lambda }_c+A(WW_c)L^{1/\nu }.$$
(41)
To procure the value of $`\nu `$, we use an iterative procedure which is equivalent to the widely used procedure of Ref. of performing a least square fit to obtain actual values of $`\xi (W)`$. We rewrite Eq.41 as

$`\mathrm{ln}\mathrm{\Lambda }`$ $`=`$ $`AL^{1/\nu }W(\mathrm{ln}\mathrm{\Lambda }_c+AL^{1/\nu }W_c)`$
$`=`$ $`a(L)Wb(L,\mathrm{\Lambda }_c,W_c),`$
assume an initial value for the critical point $`(W_c,\mathrm{ln}\mathrm{\Lambda }_c)`$ from Fig.4, obtain the functions $`a(L)`$ and $`b(L,\mathrm{\Lambda }_c)`$ by curve-fitting, extract $`A`$ and $`\nu `$, determine the value of the critical point thus obtained, and repeat the procedure till convergence is achieved (which happens rather quickly).If scaling is valid, we can collapse our data onto curves of $`\mathrm{ln}\mathrm{\Lambda }(W,L)`$ vs $`L^{1/\nu }(WW_c)`$.

In our transfer matrix calculations, we choose systems whose cross-sectional areas have linear dimensions of $`L=4,6,8,10`$, and our transfer matrices have dimensions $`4L^2\times 4L^2`$ with the given values of $`L`$. We choose the number of transfer matrices to be multiplied together by ensuring that the basis vectors do not grow upto a magnitude greater than $`10^5`$ upon being multiplied by the set of transfer matrices. We utilize a total of 2000 slices in the quasi 1-dimensional system for each value of L and W.

Fig.3 shows the plots for $`\mathrm{\Lambda }(L,W)`$ as a function of disorder for $`L=4,6,8,10`$. For fixed disorder W, an increasing $`\mathrm{\Lambda }(L,W)`$ with increasing system size $`L`$ indicates the extended regime, while a decreasing $`\mathrm{\Lambda }(L,W)`$ shows that the system is in the localized regime. In comparison to the Anderson model for normal systems, our simulations require a much smaller number of transfer matrices for relatively noise-free data, and we believe that this really is a consequence of the relatively low critical disorder strength.

Fig.4 shows the data for the iterative procedure which gives the value for the localization length exponent

$$\nu =1.15\pm 0.15.$$
(42)
Finally, Fig.5 indeed demonstrates clean data collapse close to the critical point.

#### 2 Summary of Results

Both field theory and numerics concur with the fact that the localization length exponent in the superconducting systems is significantly lower than that of their normal partners, clearly indicating a new universality class. An $`ϵ`$ expansion in $`d=2+ϵ`$ dimensions of the action in Eq.14 shows that the superconducting system has a localization length exponent $`\nu `$ of $`(ϵ+\frac{ϵ^2}{2}+𝒪(ϵ^3))^1`$ in contrast to a $`\nu _n`$ of $`(ϵ+𝒪(ϵ^3))^1`$ for normal systems. The field theory would thus predict $`\nu =\frac{2}{3}`$ versus $`\nu _n=1`$ in 3-dimensions. In comparison, one obtains the numerical estimate $`\nu =1.15\pm 0.15`$ versus $`\nu _n=1.54\pm 0.08`$ in 3-dimensions.

We must remark that the system sizes and computing power utilized in our numerical studies were relatively low compared to the current cutting edge procedures. As a lot of work has gone into refining techniques with regard to normal systems(for e.g. Ref), it is well worth employing them to study this novel phase transition and analogous ones in superconducting systems with other symmetries.

### B Density of States

The quasiparticle DOS in dirty superconducting systems exhibits some of the most stunning differences from normal systems. In normal systems, both in the Anderson metal and the Anderson insulator, i.e., in the absence of interactions, the DOS remains a smooth continuous function across the Fermi energy. In contrast, in gapless superconductors that respect $`SU(2)`$, well within the thermal metal, quantum interference effects cause a singularity at the Fermi energy which manifests itself as a $`\sqrt{E}`$ cusp in 3-dimensional systems. Deep in the thermal insulator, the density of states exhibits a power-law that vanishes at the Fermi energy with the form $`\rho |E|^\alpha `$, where $`\alpha =1`$ for systems possessing time-reversal invariance. About the critical point, the DOS once again shows power-law singularities, the details of which we discuss below. The curious form of the DOS (shown in Fig.6) has profound impact on thermodynamic properties, and in particular, manifests itself in quantities such as specific heat and spin susceptibility.

#### 1 Discussion of Critical Behaviour

The field theoretic action of Eq.14 not only offers concrete predictions for the DOS, if one were to use Wegner’s analogy with magnetic systems, it provides an intuitive picture for the behaviour about criticality. To elaborate, the quasiparticle DOS at the Fermi energy, which also gives a measure of the magnetization, acts as the ’order parameter’ of the field theory. It is given by

$$\rho =\underset{n0}{lim}\frac{\rho _0}{4n}<Tr(U^{}+U)>,$$
(43)
where $`\rho _0`$ is the bare DOS, and n the number of replicas. The field $`\eta `$, which has units of energy E, couples to the DOS in the action of Eq.14, and might be equated with the magnetic field in the magnetic analog.

With a little indulgence, one can go further with parallels between the field theory and the magnetic systems, as first suggested for normal systems:

| F.T. of $`S_{NL\sigma M}`$ | Magnetic Systems |
| --- | --- |
| Distance from criticality | Reduced Temperature |
| $`\mathrm{\Delta }=\frac{WW_C}{W_C}`$ | t |
| DOS, $`\rho `$ | Magnetization, m |
| Energy, E | Magnetic Field, h |
| $`\stackrel{~}{\chi }=\frac{d\rho }{dE}`$ | Magnetic Susceptibility,$`\chi `$ |

In normal systems, the analogy is clouded by the fact that the DOS is a continuous function of energy and disorder respectively. One can reconcile with this if the DOS obeys a power-law form with exponent zero, and in fact, one can show this to be the case on field theoretic grounds. But in superconducting systems, as we shall see, the analogy goes through in quite a striking manner, with a whole slew of nontrivial critical exponents:

| $`\rho (\mathrm{\Delta },E=0)|\mathrm{\Delta }|^\beta `$ |
| --- |
| $`\rho (\mathrm{\Delta }=0,E)|E|^{\frac{1}{\delta }}`$ |
| $`\stackrel{~}{\chi }|\mathrm{\Delta }|^\gamma `$ |
| $`\xi |\mathrm{\Delta }|^\nu `$ |

where we have used the notation of the tables above,and $`\xi `$ is the localization length and $`\nu `$ the associated exponent described in the previous section.

In order to derive expressions for $`\beta `$, $`\delta `$ and $`\nu `$, we start with the ’free energy density’ f, obtained from the action of Eq.14 :

$$f=\frac{1}{L^d}\underset{n0}{lim}\frac{ln𝒵_n}{n},$$
(44)
$$𝒵_n=\overline{Z^n}=𝑑\stackrel{~}{U}e^{S_{NL\sigma M}},$$
(45)
where $`d\stackrel{~}{U}`$,the integral volume element, takes into account the symplectic group structure of the matrices in the action $`S_{NL\sigma M}`$, n denotes the number of replicas, $`Z`$ denotes the partition function of single system, and the overbar above $`Z^n`$ refers to the average over disorder. Near criticality, $`f_s`$, the singular part of the free energy density, is expected to scale as follows:

$$f_s(\mathrm{\Delta },E)=\xi ^d\stackrel{~}{f}(\xi ^yE),$$
(46)
where ’y’ describes the scaling form of E, and Eq.38 gives the behaviour of the correlation length $`\xi |\mathrm{\Delta }|^\nu `$. Differentiating the free energy with respect to E results in the following form for the DOS:

$$\rho (\mathrm{\Delta },E)=\xi ^{d+y}F^\pm (\xi ^yE),$$
(47)
where $`F^+`$ corresponds to behaviour for $`\mathrm{\Delta }>0`$, and $`F^{}`$ for $`\mathrm{\Delta }<0`$. To obtain $`\beta `$, we set $`E=0`$, and compare the form of the resulting order parameter $`\rho (\mathrm{\Delta },E=0)`$ in the above table, yielding

$$\beta =\nu (dy).$$
(48)
Here, we require $`F^{}(0)`$ to be finite, and $`F^+(0)=0`$.

To extract $`\delta `$, we impose the physical constraint that $`\rho `$ be well-behaved and finite at criticality. This requires that $`\rho (\mathrm{\Delta }0,E)`$ be independant of the diverging correlation length, and thus yields

$$\frac{1}{\delta }=\frac{d}{y}1.$$
(49)
Finally, taking a derivative of $`\rho `$ in Eq.47 with respect to E gives us the following expression for $`\gamma `$:

$$\gamma =\nu (2yd).$$
(50)
#### 2 Results

To obtain estimates of critical exponents from field theoretic results in $`2+ϵ`$ dimensions, we use the value $`(ϵ+\frac{ϵ^2}{2})^1+𝒪(ϵ^3)`$ obtained for $`\nu `$ in the previous section, and the value $`\frac{1}{\delta }=\frac{ϵ}{4}+𝒪(ϵ^3)`$ from Ref.. Eq.48-50 then enable us to determine the critical exponents $`\beta `$ and $`\gamma `$ via the relationship $`\frac{4d}{ϵ+4}`$. Specifically, in the case of 3-dimensions, substituting the value $`ϵ=1`$, we obtain the rough estimates $`y=12/5`$, $`\delta =4`$, $`\nu =2/3`$, $`\beta =2/5`$ and $`\gamma =6/5`$.

Shifting our focus to numerical results, the method of exact diagonalization reveals that the superconducting system at hand does indeed show the novel singular behaviour in the DOS at the Fermi energy. In the data shown below, we have once more modelled the superconducting Hamiltonian after Eq.8 using periodic boundary conditions. Systems of linear dimension ’L’ have required matrices of dimension $`2L^3X2L^3`$, and we have explored system sizes with linear dimensions $`L=4,6,8`$.

Fig.7 shows the progression of the behaviour of the DOS with increasing disorder strength. As seen in the last panel of Fig.7, the DOS shows a power-law behaviour of the form $`\rho |E|`$ about the Fermi energy , $`E_F`$, consistent with expectations for the thermal insulator.

Fig.8 shows a zoom of the DOS about $`E_F`$ for disorder close to the critical strength $`W_C`$, for which we have an estimate from the localization length study of the previous section. One can easily discern that the DOS plummets down quite markedly, and does indeed exhibit singular power-law behaviour.

A plot of the DOS at the Fermi energy(Fig.9) shows that even relatively small system sizes provide numerical confirmation of the fact that $`\rho (E=0)`$ acts as the order parameter for the field theory of Eq.14; the DOS at $`E_F`$ is finite for low disorder, and it slowly drops to zero beyond a critical disorder strength. As discussed in the previous section, one would in fact expect the DOS for an infinite sized system to behave as $`\rho (\mathrm{\Delta },E=0)|\mathrm{\Delta }|^\beta `$, where $`\mathrm{\Delta }`$ is the distance from criticality within the thermal metal. Scaling arguments for extracting $`\beta `$ require

$$\rho _L(\mathrm{\Delta },E=0)=|\mathrm{\Delta }|^\beta Y(L\mathrm{\Delta }^\nu ),\mathrm{\Delta }<0,$$
(51)
where $`\rho _L`$ is the DOS associated with a system of linear dimension L, and $`Y`$ is a scaling function. One can rewrite the above equation in a form more conducive to numerics as follows:

$$\rho _L(\mathrm{\Delta },E=0)=L^{\frac{\beta }{\nu }}\stackrel{~}{Y}(\mathrm{\Delta }L^{\frac{1}{\nu }}),$$
(52)
where $`\stackrel{~}{Y}`$ is yet another scaling function with limiting behaviour $`\stackrel{~}{Y}(x\mathrm{})=|x|^\beta `$, reproducing the required dependence of $`\rho (\mathrm{\Delta },E=0)`$ on $`\mathrm{\Delta }`$ for infinite system size.

Fig.10 exhibits the plots of $`\rho _L`$ as a function of disorder for different system sizes L, and we make use of this data to procure the value of $`\beta `$ in Eq.52. To extract $`\beta `$, we perform a fit taking $`W_c`$, $`\nu `$ and $`\beta `$ as variable parameters. Exploiting the universal nature of the function $`\stackrel{~}{Y}`$ in Eq.52, we find the appropriate values of $`\rho `$ obtained by linear interpolation for a given set of system sizes and fixed argument in $`\stackrel{~}{Y}`$, and plot these on a log-log scale versus system size; the slope for a linear fit of such a set of points determines $`\beta `$. The actual value of $`\beta `$ is obtained by performing the above procedure for different values of the argument of $`\stackrel{~}{Y}`$ and taking the average of the $`\beta `$’s thus obtained.

The set of values $`\beta =0.15`$, $`W_c=4.67`$ and $`\nu =1.25`$ result in the data collapse shown in Fig.11. In comparison, as remarked at the beginning of this section, the field theoretic result predicts that $`\beta =0.4`$. Once more, as in the case of the localization length exponent, we comment on numerical accuracy; other simulations using exact diagonalization, for instance, those catering to specific physical situations, have used larger sytems sizes and number of realizations which would be well-worth employed here. However, the above numerics conveys quite clearly that the DOS at criticality exhibits a power law suppression about the Fermi energy, and that $`\rho (E=0,W)`$ acts as an ’order parameter’ with a non-trivial exponent $`\beta `$ in surprising contrast to normal systems which have $`\rho (E,W_c)`$ smooth about $`E_F`$, and a vanishing exponent $`\beta `$.

## IV Other Systems

We have studied the thermal metal-thermal insulator transition in superconducting systems with $`SU(2)`$ and $`T`$, and discussed the dirty d-wave superconductor as a possible physical realization. Superconducting systems with other symmetries too promise such a transition.

The properties and phases of superconducting systems with spin-rotational invariance, but no time reversal invariance are rather similar to our case which preserves both symmetries. The thermal metal and the thermal insulator are both stable only in 3-dimensions, and the transition cannot occur in lower dimensions where quasi-particle excitations are always found to be localized at large enough length scales. Due to the absence of time-reversal symmetry, the Hamiltonian of Eq.6 now has imaginary couplings. As described in detail in previous work, in this symmetry category, the pinned vortex state of a type II s-wave superconductor appears to be a fine candidate for exhibiting the thermal metal-thermal insulator transition. Low energy quasiparticles bound to the core of the vortices can tunnel from one vortex to another, and as the magnetic field is increased, the density and tunneling strength also increase. It is conjectured that there could exist a critical magnetic field $`H_{c4}`$, within the vortex phase at which the low energy states can permeate through the medium to form extended states. Field theoretic methods have analyzed the properties of these systems, and in particular, have shown that in parallel to the case with T, the density of states at the critical point has the power-law behaviour $`\rho |E|^{ϵ/2}`$, where $`E`$ is the energy and $`ϵ=d2`$.

The presence of spin-orbit scattering or of triplet pairing introduces new ingredients. It breaks spin-rotational symmetry, and as in the case of normal systems, field theoretic arguments predict the presence of a delocalization-localization transition not only in 3-dimensions but also in 2-dimensions.In Ref., T.Senthil et.al. cast their Hamiltonian in terms of Majorana fermions, and their formulation, among other things, is highly conducive to numerical work. They discuss possibilities for experiment, and mention heavy fermion systems where spin-orbit scattering is prominent, as a physical realization.

Finally, the tantalizing prospect of a Hall effect in superconducting systems has been explored in systems with other symmetries as well,i.e., those with SU(2). Ref shows that a superconductor with unconventional $`d_{x^2y^2}+id_{xy}`$ pairing symmetry is capable of exhibiting a phase with non-vanishing spin and thermal Hall conductances. Indeed, as in normal systems, sophisticated methods such as employing the network model and supersymmetric spin chains, have shed light on these systems.

## V Experiment

The transport properties of normal systems have been probed in great depth, and now we see that superconducting systems could potentially offer an equally rich range of experiments in the thermal metal-thermal insulator transition. In the previous section, we have mentioned a variety of experimental candidates for study, such as the dirty d-wave, the type II s-wave, heavy fermion and other superconducting systems. These systems must share the feature of gapless superconductivity; one requires states at the Fermi energy since it is these states that determine transport properties, and distinguish the thermal insulator from the thermal metal. Associated with each system, a tunable parameter such as disorder or magnetic field ought to be able to access the phases. As was previously discussed in detail, the type II s-wave superconductor in the vortex phase offers promise as a likely candidate for observing the novel transition since in principle, one need only tune the magnetic field, and generally the vortex phase exists over a large range of field.

As seen earlier, in contrast to normal systems, the density of states shows singular behaviour about the Fermi energy for both phases and at the critical point. We saw that for systems with T and SU(2), it obeys a power-law behaviour of the form $`\rho (E)|E|^\alpha `$, where $`\alpha =\frac{1}{2}`$ well within the thermal metal, $`\alpha =1`$ deep in the thermal insulator, and field theory predicts $`\alpha =\frac{ϵ}{4}`$ at the transition, with $`ϵ=1`$ for 3-dimensional systems. This singular behaviour ought to be reflected in thermodynamic quantities such as specific heat and spin susceptibility, and in tunneling experiments. In particular, the temperature dependence of the thermodynamic quantities would have a form $`CT^{1+\alpha }`$ for the specific heat, and $`\chi T^\alpha `$ for the spin susceptibility. Recent experiments of cuprate superconductors have observed a suppression of the specific heat close to the Fermi energy. However, measuring the differing behaviours to detemine the phase might prove tricky, especially since we have neglected various effects such as interactions which could come into play.

For systems that preserve spin-rotational invariance, quasiparticle excitations about the superconducting ground state not only conserve energy, but also spin; the spin-conductance $`\sigma _s`$ could be employed to determine whether the system inhabits the thermal metal or thermal insulator. In the thermal metallic phase, a magnetic field gradient would cause the spinful quasiparticle to diffuse across the sample, while in the thermal insulator, they would be unable to conduct spin. A variety of refined spin-injection techniques have been developed in semi-conductors to measure spin dynamics (see for e.g. Ref.), but by no means would it be a simple task to cater these experiments to superconductors.

We believe that by far, thermal conductivity measurements would offer most promise in probing the thermal metal-thermal insulator transition. For all superconducting systems with their differing symmetries, the thermal conductivity $`\kappa `$, distinguishes the two phases in that the ratio $`\kappa /T`$ tends to a finite constant in the thermal metal, and to zero as $`T0`$ in the thermal insulator. Along the lines of previous experiments, it would be extremely interesting to observe the transition by applying a small thermal gradient across a superconducting sample and measuring $`\kappa /T`$ as a function of a tunable parameter.

We thank T.Senthil and I.A. Gruzberg for many an illuminating conversation, and D. Whysong for indispensible advice on numerical work. This research was supported by NSF Grants DMR-97-04005, DMR95-28578 and PHY94-07194.


text:
# References

Neutron matter with a model interaction

M. Ya. Amusiaa,b, V. R. Shaginyanc,d 111E–mail: vrshag@thd.pnpi.spb.ru

aThe Racah Institute of Physics, Hebrew University, Jerusalem 91904, Israel;

bA.F. Ioffe Physical-Technical Institute, 194021 St.Petersburg, Russia;

cPetersburg Institute of Nuclear Physics, 188350 Gatchina, Russia;

d Department of Physics, University of Washington, Seattle, WA 98195, USA

## Abstract

An infinite system of neutrons interacting by a model pair potential is considered. We investigate a case when this potential is sufficiently strong attractive, so that its scattering length $`a`$ tends to infinity, $`a\mathrm{}`$. It appeared, that if the structure of the potential is simple enough, including no finite parameters, reliable evidences can be presented that such a system is completely unstable at any finite density. The incompressibility as a function of the density is negative, reaching zero value when the density tends to zero. If the potential contains a sufficiently strong repulsive core then the system possesses an equilibrium density. The main features of a theory describing such systems are considered.

PACS. 24.10.Cn Many-body theory - 21.10.Dr Binding energies

There exists a well known problem in many body physics related to the description of the ground state properties of an infinite system composed of interacting fermions. In general, this description is based usually on tedious numerical calculations, particularly when the interaction is rather strong. The well known exceptions from this situation, when it is possible to calculate the ground state properties analytically, are the Random Phase Approximation (RPA) for a high density electron gas  and the low density approximation for dilute gases . In both cases the kinetic energy $`T_k`$ is much bigger then the interaction energy $`E_{int}`$ of the system. This allows to apply some kind of a perturbation theory. In the case of homogeneous electron liquid it turns out that the analytical RPA-like description is also possible not only at very high but medium densities when $`T_kE_{int}`$ . Similar extension of the range of validity is impossible in the case of fermion systems at low densities $`\rho `$: there the gas approximation is not applicable if $`T_kE_{int}`$. In the cases when the pair interaction is attractive and sufficiently strong, the system can have a quasi equilibrium or equilibrium states in which $`T_kE_{int}`$. On the other hand, these states are preceded by special points with density $`\rho `$ values at which the incompressibility $`K(\rho )`$ of the system tends to zero. Thus, if it would be possible to predict the existence of such points then in principle it would become possible to conclude that the system has at least a quasi equilibrium state.

In this Short note we address the above mentioned problem and consider the ground state properties of the infinitely extended multi-fermion system. We demonstrate that it can be done analytically provided that the pair interaction between fermions is characterized only by the scattering length $`a\mathrm{}`$. One can say in this case that the scattering length is the dominant parameter of the problem under consideration. Such an investigation is of great importance since it can be applied to fermion systems interacting via potentials with not only infinite, but also sufficiently large $`a`$. For instance, the scattering length $`a`$ of neutron-neutron interaction is about $`20`$ fm, thus greatly surpassing the radius of the interaction $`r_0`$. On the other hand, it is possible now to prepare artificially a system composed of Fermi atoms interacting by an artificially constructed potential with almost any desirable scattering length, similarly to that how it is done for the trapped Bose gases, see e.g. . An experimental study, performed on such Fermi-system would be of great importance presenting new information on the behavior of dilute gases and the gas-liquid phase transition.

Let us consider the interaction of two isolated particles. We assume that this interaction is of finite radius $`r_0`$, which is small, so that $`p_Fr_01`$ ($`p_F=(3\pi ^2\rho )^{1/3}`$ is the Fermi momentum), but its strength is such that the scattering length is negative and infinitely big, $`a\mathrm{}`$. We assume also, that the density $`\rho `$ of the system under consideration is homogeneous. As it will be demonstrated below, in such a case the system is located in the vicinity of a phase transition, transforming it into a strongly correlated one. Therefore, the problem of calculating its ground state properties has to be treated for the most part qualitatively.

Let us start considering general properties of a Fermi system with some attractive two-particle bare interaction $`V(r)`$, which is sufficiently weak to create a two-particle bound state. Assume, that the scattering length $`a`$ corresponding to this potential is negative and finite. The ground state energy density $`E(\rho )`$ can be approximated by Skyrme-like expression ,

$$E(\rho )=\frac{3p_F^2}{10M}\rho +t_0\rho ^2+t_3\rho ^{7/3},$$
(1)
Here $`M`$ is the particle mass, and $`\rho `$ is their density in the system. The first term of eq. (1) is the kinetic energy $`T_k`$, while the second and the third are related to the interaction energy $`E_{int}`$ determined by the potential $`V(r)`$. The second term which is proportional to $`t_0`$ gives a proper description in the gas limit. The third term provides the behavior of $`E(\rho )`$ at higher densities, including that of the equilibrium density. Such structure of $`E(\rho )`$ appears if the interaction is sufficiently attractive so that $`t_0<0`$, and $`t_3>0`$. Note, that eq. (1) presents at least a qualitative description of the system under consideration giving a rather simple and reasonable presentation of the function $`E(\rho )`$. A more precise picture of the energy dependence upon density can be obtained using more sophisticated functionals for the ground state energy .

If the potential $`V`$ is of short range and purely attractive, then in the Hartree-Fock approximation the ground state energy $`E_{HF}`$ is given by the following expression

$$E_{HF}(\rho )=\frac{3p_F^2}{10M}\rho +t_{HF}\rho ^2,$$
(2)
where the parameter $`t_0=t_{HF}`$ , being negative, is entirely determined by the potential $`V(r)`$. For instance, in the case of a short range $`\delta `$-type interaction one has $`t_{HF}=v_0/4`$, with $`v_0`$ being the corresponding strength of the potential. Eq. (2) shows that at small densities $`E_{HF}>0`$ due to the kinetic energy term, but at sufficiently high densities $`\rho \mathrm{}`$ the Hartree-Fock energy becomes dominating, leading to the collapse of the system, with $`E_{HF}\mathrm{}`$. Keeping in mind that the Hartree-Fock approximation gives the upper limit to the binding energy $`E_{HF}E`$, one can conclude that the system does not have, in this case, an equilibrium density $`\rho _e`$ and energy $`E_e`$ since $`E_e\mathrm{}`$ when $`\rho \mathrm{}`$ . Note, that for a given and finite total number of particles $`N`$, the HF energy is not going to infinity and the system collapses into a small volume with the radius $`r_0`$, with the density $`\rho N/r_0^3`$. It is evident that the function $`E(\rho )`$ is positive at small densities, if the parameter $`t_0`$ is finite. Therefore, it must have at least one maximum at the density $`\rho _m`$ before it becomes negative, on the way to $`E\mathrm{}`$. If the potential $`V(r)`$ includes a repulsive core at sufficiently short distances, then $`t_3>0`$ . As a result, the system has an equilibrium density and energy, $`\rho _e`$ and $`E_e`$, respectively, determined by the repulsive core strength and its radius $`r_cr_0`$.

Now let us apply eq. (1) to demonstrate the most important features of the system under consideration:

a) when $`\rho 0`$ the third term on r.h.s. in eq. (1) can be omitted. The kinetic energy is relatively very big, $`T_kE_{int}`$, and $`t_0a`$, with $`a<0`$ being the scattering length. In that case we have a dilute Fermi gas with positive pressure $`P`$ and incompressibility $`K`$, the latter being determined by the equation, see e.g. ,

$$K(\rho )=\rho ^2\frac{dE^2(\rho )}{d\rho ^2}.$$
(3)
b) on the way to higher densities, which can be achieved by applying an external pressure, the system reaches the density $`\rho _{c1}<\rho _m`$ at which the incompressibility is equal to zero, $`K(\rho _{c1})=0`$. Remembering that at the maximum the second derivative is negative, one can conclude, as it is seen from eq. (3), that $`K(\rho _m)<0`$. In the range $`\rho _{c2}\rho \rho _{c1}`$ the incompressibility is negative, $`K<0`$, and as a result the system becomes totally unstable. In fact, in this density range all calculations of the ground state energy are meaningless since such a system cannot exist and thus be observed experimentally ;

c) at some point $`\rho =\rho _{c2}>\rho _m`$ the contribution due to the repulsive core becomes sufficiently strong to prevent the further collapse of the system. The incompressibility attains $`K=0`$ at $`\rho _{c2}<\rho _e`$, being positive at the higher densities. Finally, the system becomes stable at $`\rho >\rho _{c2}`$, reaching equilibrium density at $`\rho _e`$ with equilibrium energy equal to $`E_e`$. It is obvious that $`K(\rho _e)>0`$ being proportional to the second derivative at the minimum, see eq. (3). It should be kept in mind that in this density domain, $`\rho \rho _{c2}`$, the function $`E(\rho )`$ is determined by the repulsive part of the potential which makes $`t_3>0`$. As it was mentioned above, without this component of $`V(r)`$ the system’s energy would infinitely increase, $`E_{HF}\mathrm{}`$, with density growth, $`\rho \mathrm{}`$, thus inevitably collapsing.

One could expect in principle the existence of metastable states at $`\rho >\rho _{c1}`$ if the potential $`V(r),`$ even being pure attractive, would have a complicated structure. It can be said that there could exist parameters of $`V(r)`$, which are able to open the possibility for the metastable states to be formed. On the other hand, a system of fermions interacting via a short-range, finite scattering-length, $`\delta `$-type potential $`V_\delta `$ , can be stable only in the dilute gas regime. While at the densities $`\rho \rho _{c1}`$ the incompressibility $`K`$ becomes negative, the system collapses. Indeed, the potential $`V_\delta `$ has no structure to ensure any metastable states at the densities $`\rho \rho _{c1}`$. As a result, one can write down a dimensionless expression for the ground state energy as a function of the only variable $`z=p_Fa`$ ,

$$\alpha E(z)=z^5(1+\beta (z)),$$
(4)
with $`\alpha =10\pi ^2Ma^5`$. In the low density limit, $`|ap_F|1`$ and when the interaction has the radius $`r_0`$, eq. (4) reads ,

$$\alpha E(z)=z^5\left[1+\frac{10}{9\pi }z+\frac{4}{21\pi ^2}(112\mathrm{ln}2)z^2+\left(\frac{r_0}{a}\right)^3z^3\gamma (\frac{r_0}{a},z)+\mathrm{}\right].$$
(5)
Here the function $`\gamma (y,z)`$ is of the order of one, $`\gamma (y,z)1`$. It is seen from eq. (5) that as soon as the scattering length becomes big enough, $`|a|r_0,`$ one can omit the contribution coming from the function $`\gamma `$ and neglect all the term proportional to $`(r_0/a)^3`$. Then eq. (5) reduces to eq. (4). Thus, in the case $`|a|\mathrm{}`$ we can use eq. (4) to determine the ground state energy $`E`$. Eq. (4) is valid up to the density $`\rho _{c1}`$ which is a singular point of the function $`\beta (z)`$, since beyond this point $`K<0`$, and the system is completely unstable. On the other hand, there is no physical reasons to have another irregular point in the region $`0\rho \rho _{c1}`$. We note, that eq. (5) is as well valid up to its own density $`\rho _{c1}^{^{}}\rho _{c1}`$, provided $`|a|r_0`$. Using eq. (3) for the incompressibility, one can calculate the position of the point $`z_{c1}`$ where $`K=0`$. Denoting the corresponding $`z`$ as $`z_{c1}=c_0`$, where $`c_0`$ is a dimensionless number, one is led to the conclusion that $`\rho _{c1}|a|^3`$ provided $`a`$ is sufficiently large to be the only dominating parameter. The system has only one stable region at small densities $`\rho \rho _{c1}`$ which decreases and even vanishes as soon as $`a\mathrm{}`$. One could expect that $`|c_0|\mathrm{}`$ as soon as $`a`$ becomes the dominant parameter so that the above given expression for $`\beta (z)`$ is valid in the whole domain $`|z|\mathrm{}`$. On the other hand, there exists another singular point $`z_{c2}`$ in the function $`E(\rho )`$ and the position of this point which corresponds to $`\rho _{c2}`$ depends mainly on the parameters such as $`r_0,`$ core radius $`r_c`$ of $`V(r)`$ rather than on the scattering length $`a`$. As to the function $`\beta (z),`$ by definition it does not contain any information about $`\rho _{c2}`$. Therefore, in order to take into account the existence of $`\rho _{c2}`$, one has to recognize that $`c_0`$ is finite, and the densities $`\rho _{c1}`$ and $`\rho _{c2}`$ are different. As a result, the function $`\beta `$ is determined in fact only in the region $`|z||z_{c1}|`$.

Now let us consider the calculation of the ground state energy $`E`$ of the system when the density approaches $`\rho \rho _{c1}|a|^3`$ from the low density side. As a rule, points $`\rho _{c1}`$ and $`\rho _{c2}`$ are missed in calculations because of the lack of the self consistency , which relates the linear response function of system with its incompressibility $`K`$,

$$\chi (q0,i\omega 0)=\left(\frac{d^2E}{d\rho ^2}\right)^1.$$
(6)
As we shall see below, these points can give important contributions to the ground state energy. To see it we express the energy of a system in the following form (see e.g. ),

$$E(\rho )=T_k+E_H\frac{1}{2}\left[\chi (q,i\omega ,g)+2\pi \rho \delta (\omega )\right]v(q)\frac{d𝐪d\omega dg}{g(2\pi )^4},$$
(7)
where $`\chi (q,i\omega ,g)`$ is the linear response function on the imaginary axis and $`v(q)=gV(q)`$, with $`V(q)`$ being the Fourier image of $`V(r)`$. The integration over $`\omega `$ goes from $`\mathrm{}`$ to $`+\mathrm{}`$, while the integration over the coupling constant $`g`$ runs from zero to the real value of the coupling constant, i.e. to $`g=1`$. At the point $`\rho =\rho _{c1}`$ the linear response function has a pole at the origin of coordinates $`q=0,\omega =0`$ due to eq. (6). At the densities $`\rho >\rho _{c1}`$ the function $`\chi (q,i\omega )`$ has poles at finite values of the momentum $`q`$ and frequencies $`i\omega `$. This prevents the integration over $`i\omega `$, making the integral in eq. (7) divergent. Thus, we conclude that it is the contribution of these poles that reflects the system’s instability in the density range $`\rho _{c1}\rho \rho _{c2}`$. Note, that violations of eq. (6) lead to serious errors in the calculation of the ground state energy. Eq. (7) can be rewritten, explicitly accounting for the effective interparticle interaction $`R(q,i\omega ,g)`$, (see e.g. ), in the following form

$$E(\rho )=T_k+E_H\frac{1}{2}\left[\frac{\chi _0(q,i\omega )}{1R(q,i\omega ,g)\chi _0(q,i\omega )}+2\pi \rho \delta (\omega )\right]v(q)\frac{d𝐪d\omega dg}{g(2\pi )^4}.$$
(8)
Here $`\chi _0`$ is the linear response function of noninteracting particles, while $`\chi `$ is given by the following equation

$$\chi (q,\omega )=\frac{\chi _0(q,\omega )}{1R(q,\omega )\chi _0(q,\omega )}.$$
(9)
It is seen from eqs. (6) and (9) that the denominator $`(1R\chi _0)`$ vanishes at $`\rho \rho _{c1}`$ while the radius of correlation tends to infinity . Thus, it is impossible to present the denominator as a power series in $`R\chi _0`$ approximating the expansion by the finite number of terms. This result is quite obvious since $`\rho _{c1}`$ is a singular point in the function $`E(\rho )`$ which makes it impossible to expand that function in the vicinity of this point. Therefore, one should try to satisfy eq. (6) in order to get proper results for the ground state calculations in the vicinity of the instability points. Such an approach was suggested in  and is based on the exact functional equation for the effective interaction $`R(q,\omega ,g)`$,

$$R(q,\omega ,g_0)=g_0v(q)\frac{1}{2}\frac{\delta ^2}{\delta \rho ^2(q,\omega )}\frac{\chi _0(k,iw)}{1R(k,iw,g)\chi _0(k,iw)}v(q)\frac{d𝐤dwdg}{g(2\pi )^4}.$$
(10)
As a result, the linear response function $`\chi `$ given by eq. (9) automatically satisfies eq. (6) . Our preliminary calculations , based on eq. (9) and applied to the case when the scattering length is sufficiently large but finite, confirm the result that $`\rho _{c1}|a|^3`$.

Let us suppose for a while that the bare potential is pure attractive. Then, the interval of the densities $`[0,\rho _{c1}]`$ within which the system is stable vanishes with the growth of $`|a|`$. As a result, in the limit $`a=\mathrm{}`$ the incompressibility becomes negative $`K0`$, making the considered system completely unstable at any density. Thus, the point at which $`a=\mathrm{}`$ is the only point of the system’s instability at all the densities. As soon as the scattering length deviates from its infinite value, that is $`+\mathrm{}>a>\mathrm{}`$ the system comes back to its stable state at list in the range of the density values $`\rho <\rho _{c1}|a|^3`$. It is of interest to understand whether it is possible to prove by e.g. numerical calculations, that $`\rho _{c1}|a|^3`$ when $`a\mathrm{}`$. From our point of view, at least at this moment, the answer is “no”. We are dealing with a system located in the vicinity of a phase transition, which transforms it into a strongly correlated one. As a result, it is hard to believe that the numerical calculations could be reliable. On the other hand, it is not really necessary to carry out numerical calculations if the problem allows a qualitative analysis. It was argued above, that there exists the only parameter to characterize the system which is the scattering length $`a`$. In fact the scattering length determines only the specific point $`\rho _{c1}`$ at which the incompressibility vanishes, separating the region of a dilute gas from the region of the system’s instability. As soon as $`a\mathrm{}`$ this last and the only parameters vanishes, driving the point $`\rho _{c1}|a|^3`$ of the curve $`E(\rho )`$ to the origin of coordinates. Thereafter, the system becomes unstable at all densities. And vice versa, as soon as the scattering length becomes finite the system is stable at list within the interval $`\rho \rho _{c1}|a|^3`$.

Note, that as it follows from our consideration, any Fermi system possesses an equilibrium density and energy if the bare particle-particle interaction contains a repulsive core and its attractive part is strong enough, so that $`a\mathrm{}`$. Indeed, at sufficiently small densities the ground state energy is negative (since the incompressibility $`K0)`$ and the system will collapse until the core stops the density growth. Therefore, the minimal value of the ground state energy must be negative when the repulsive core will enter the play to prevent the system from the further collapse. It is worth to remark, that superfluid correlations cannot stop the system squeezing, since their contribution to the ground state energy being negative increases in the absolute value with the growth of the density.

A liquid similar to the model one considered in this paper exists in Nature. This is liquid $`{}_{}{}^{3}He`$. If a helium dimer exists, its bound energy does not exceed $`10^4`$ meV while the ground state energy of helium liquid is about $`210^1`$ meV per atom . Because of this huge difference in binding energies, it is evident that there is no essential contribution coming from the binding energy of the dimer to the ground state energy of the liquid. In fact, the numerical calculations show that the pair potential is rather weak to produce the dimer $`He_2`$ . Thus, one can reliably consider an infinite homogeneous system of Helium atoms as consisting of particles interacting via pair potential, characterized by a very big but finite scattering length $`|a|r_0`$. Let us make also the following additional remark. It seems quite probable that the neutron-neutron scattering length ($`a20`$ fm) is sufficiently large to permit the neutron matter to have an equilibrium energy and density . Therefore, calculations of a neutron matter satisfying eq. (7) are quite desirable.

In summary, the homogeneous system of interacting fermions was considered. It was shown that when the scattering length $`a`$ is negative and sufficiently large the fermion matter becomes a strongly correlated system at the densities $`\rho |a|^3`$. Therefore, the consideration of such a system is connected to a number of problems which yet persist and have to be resolved. At the same time, the qualitative consideration presented above gives strong evidences that the point $`\rho _{c1}`$ at which the incompressibility vanishes is defined by $`\rho _{c1}|a|^3`$ provided the scattering length is the dominant parameter of the problem. Thus, a homogeneous system composed of fermions, interacting via a pure attractive potential, at $`a\mathrm{}`$ is completely unstable at all the densities, with the incompressibility as a function of the density being always negative. As soon as the density $`\rho `$ goes to zero the incompressibility goes to zero as well.

We thank G.F. Bertsch for attracting our attention to the discussed above many body problem which he has in fact suggested . One of us (VRS) is grateful to A. Bulgac for valuable discussions, and to the Department of Physics of the University of Washington where part of this work was done, for hospitality. This research was funded in part by INTAS under Grant No. INTAS-OPEN-97-603.