text:
# A Relativistic Mean Field Model for Entrainment in General Relativistic Superfluid Neutron Stars

## I Introduction

A new generation of gravitational wave detectors (LIGO, VIRGO, etc) are now working to detect gravitational waves from compact objects, such as black holes and neutron stars. With this detection we expect to have a unique probe of the physics that dictates their behavior. This is ushering in a new era where strong-field relativistic effects will play an increasingly important role. Only through their inclusion can we hope to accurately decipher what gravitational wave data will have to tell us. With that in mind, we present here a fully relativistic model of the so-called entrainment effect (to be described in some detail below) that is a necessary feature of the dynamics of superfluid neutron stars.

For the densities appropriate to neutron stars there are attractive components of the strong force that should lead, via BCS-like mechanisms, to nucleon superfluidity and superconductivity. Indeed, calculations of supra-nulcear gap energies consistently lead to the conclusion that superfluid neutrons should form in the inner crust of a mature neutron star, with superfluid neutrons and superconducting protons in the core. Even more exotic possibilities have been suggested, such as pion condensates, superfluid hyperons, and superconducting quark matter. Perhaps most important is the well-established glitch phenomenon in pulsars the best description of which is based on superfluidity and quantized vortices. Superfluidity should affect gravitational waves from neutron stars by modifying the rotational equilibria and the modes of oscillations that these objects support AC01c ; AC01b ; GC02 .

The success of superfluidity in describing the glitch phenomena is due in part to the fact that the superfluid neutrons of the inner crust represent a component that can move freely (for certain timescales) from the rest of the star. Explaining the glitch phenomena then becomes a question of how to transfer angular momentum between the various “rotationally decoupled” components. For the modes of oscillation, it is by now well established that a similar “decoupling,” this time between the superfluid neutrons of the inner crust and core and a conglomerate of the remaining charged constituents (e.g. crust nuclei, core superconducting protons, and crust and core electrons), leads to a mode spectrum for superfluid neutron stars that is quite different from that of their ordinary fluid counterparts (see GC02 , and references therein, for a complete review).

Several recent studies AC01a ; PCA02 ; ACL02 ; ACP02a ; ACP02b  have established that the entrainment effect is an important element in modelling the rotational equilibria and modes of oscillation of superfluid neutron stars. Sauls S89  describes the entrainment effect as a result of the quasiparticle nature of the excitation spectrum of the superfluid and superconducting nucleons. That is, the bare neutrons (or protons) are accompanied by a polarization cloud containing both neutrons and protons. Since both types of nucleon contribute to the cloud the momentum of the neutrons is modified so that it is a linear combination of both the neutron and proton particle number density currents, and similarly for the proton momentum. Thus when one species of nucleon acquires momentum, both types of nucleons will begin to flow.

In the core of a neutron star, the Fermi energies of nucleons (as well as some of the leptons) can become comparable to their mass-energies, because the Fermi energies are a function of the local particle number densities, and these can be quite high. This implies that any Newtonian model for entrainment must become less reliable as one probes deeper into the core of a neutron star, and thus a relativistic formulation is required. In fact, we will see that the Newtonian parameterized model of Prix et al. PCA02  does deviate most from the relativistic model in the core. There are two purposes for which a relativistic formulation is necessary. At the microscopic level, the nucleons will (locally) have average speeds that are comparable to the speed of light. As well, at a mesoscopic level, the fluid elements, which contain a large number of nucleons, could have average speeds that are also comparable to the speed of light. The formalism that we develop here will be relativistic in both respects. One should note, though, that in realistic astrophysical scenarios (e.g. when an isolated neutron star undergoes linearized oscillations, or a pulsar exhibits a glitch) the fluid element average speeds are typically only a few percent of that of light.

To date, studies of superfluid dynamics in neutron stars have relied on models of entrainment that are obtained in the Newtonian regime. For instance, a few of the most recent studies LM00 ; ACL02  have employed a parameterized model for entrainment that is inspired by the Newtonian, Fermi-liquid calculations of Borumand et al. BJK96 . An alternative formulation PCA02 —motivated by mathematical simplicity that allows for analytic solutions for slowly rotating Newtonian superfluid neutron stars—for parameterizing entrainment has been recently put forward. Here we take a different approach, and this is to use a $`\sigma \omega `$ relativistic mean field model, of the type that is described in detail by Glendenning G97 . Although a relativistic Fermi-liquid formalism exists CB76 , we prefer to use the mean field model because it is sufficiently simple that semi-analytical formulas result, and a clear connection between the coupling parameters at the microscopic level can be made to the macroscopic properties (such as mass and radius) of the star. An immediate consequence is the ability to compare the relativistic entrainment model with the two parameterized models. We will see that the model used by Prix et al. PCA02 , although limited, is a better fit than the other formulation.

The next section begins with a review of the $`\sigma \omega `$ model. That is followed by an application of the mean field approximation to obtain an equation of state that includes entrainment. In Sec. 3, we briefly review the general relativistic superfluid formalism and how it is used to describe slowly rotating configurations. We then use the mean field results to produce explicit models. After some concluding remarks, an appendix is given that contains some of the technical details and results. Throughout we will use “MTW” MTW  conventions, a consequence of which is that several equations will have minus sign differences with, for instance, those of G97 .

## II Relativistic Mean Field Theory of Coupled Fluids

To create a seamless conceptual basis for general relativistic calculations of dynamic processes in neutron stars, we need a covariant formalism that describes the strongly interacting coupled neutron and proton fluids. It should be sufficiently simple that it provides physical insight, yet accurate enough that it can serve as the basis for realistic numerical calculations. For static stars, this role is played by the $`\sigma \omega `$ effective mean-field theory G97 . Our task in this paper is to generalize this theory to dynamic stars. In particular, we are interested in situations where there is relative motion of the two fluids, since the entrainment of one by the other turns out to play a large role in the dynamics.

The Lagrangian density for the baryons and the mesons that the baryons exchange is as in the static case. It is

$$L=L_b+L_\sigma +L_\omega +L_{int},$$
(1)
with

$$L_b=\overline{\psi }(i\gamma _\mu ^\mu m)\psi $$
(2)
as the baryon Lagrangian. Here $`\psi `$ is an 8-component spinor with the proton components as the top 4 and the neutron components as the bottom 4. The $`\gamma _\mu `$ are the corresponding $`8\times 8`$ block diagonal Dirac matrices. The Lagrangian for the $`\sigma `$ mesons is

$$L_\sigma =\frac{1}{2}_\mu \sigma ^\mu \sigma \frac{1}{2}m_\sigma ^2\sigma ^2.$$
(3)
The Lagrangian for the $`\omega `$ mesons is

$$L_\omega =\frac{1}{4}\omega _{\mu \nu }\omega ^{\mu \nu }\frac{1}{2}m_\omega ^2\omega _\mu \omega ^\mu $$
(4)
where $`\omega _{\mu \nu }=_\mu \omega _\nu _\nu \omega _\mu .`$ The interaction Lagrangian density is

$$L_{int}=g_\sigma \sigma \overline{\psi }\psi g_\omega \omega _\mu \overline{\psi }\gamma ^\mu \psi .$$
(5)
The Euler-Lagrange equations are

$`\left(\mathrm{}+m_\sigma ^2\right)\sigma `$ $`=`$ $`g_\sigma \overline{\psi }\psi ,`$ (6)
$`\left(\mathrm{}+m_\omega ^2\right)\omega _\mu +_\mu ^\nu \omega _\nu `$ $`=`$ $`g_\omega \overline{\psi }\gamma _\mu \psi ,`$ (8)
$`(i\gamma _\mu ^\mu m)\psi `$ $`=`$ $`g_\omega \gamma _\mu \omega ^\mu \psi g_\sigma \sigma \psi .`$ (10)
Finally, the stress-energy tensor takes the form

$$T^{\mu \nu }=T_b^{\mu \nu }+T_\sigma ^{\mu \nu }+T_\omega ^{\mu \nu }+T_{int}^{\mu \nu }$$
(11)
containing contributions from the baryons $`(b)`$, the mesons $`(\sigma ,\omega )`$, and the interaction. Individually, these are

$`T_b^{\mu \nu }`$ $`=`$ $`i\overline{\psi }(\gamma ^\mu ^\nu \eta ^{\mu \nu }\gamma ^\alpha _\alpha )\psi m\eta ^{\mu \nu }\overline{\psi }\psi ,`$ (12)
$`T_\sigma ^{\mu \nu }`$ $`=`$ $`^\mu \sigma ^\nu \sigma {\displaystyle \frac{1}{2}}\eta ^{\mu \nu }m_\sigma ^2\sigma ^2{\displaystyle \frac{1}{2}}\eta ^{\mu \nu }^\alpha \sigma _\alpha \sigma ,`$ (14)
$`T_\omega ^{\mu \nu }`$ $`=`$ $`\left(^\mu \omega ^\alpha ^\alpha \omega ^\mu \right)^\nu \omega _\alpha {\displaystyle \frac{1}{2}}\eta ^{\mu \nu }m_\omega ^2\omega ^\alpha \omega _\alpha {\displaystyle \frac{1}{4}}\eta ^{\mu \nu }m_\omega ^2\omega ^{\alpha \beta }\omega _{\alpha \beta },`$ (16)
$`T_{\iota nt}^{\mu \nu }`$ $`=`$ $`\eta ^{\mu \nu }g_\sigma \sigma \overline{\psi }\psi \eta ^{\mu \nu }g_\omega \omega _\alpha \overline{\psi }\gamma ^\alpha \psi .`$ (18)
We now solve these equations in the mean field approximation, eventually in a frame in which the neutrons have zero spatial momentum while the protons have on average a wavevector $`K_\mu =(K_0,0,0,K_z).`$ In this approximation we ignore all gradients of the averaged sigma and omega fields, and the neutrons and protons are taken to be in plane-wave states. The problem simplifies considerably and we find for the $`\sigma `$ and $`\omega _\mu `$ fields and the stress-energy tensor $`T_\nu ^\mu `$ that

$`m_{}`$ $`=`$ $`mc_\sigma ^2\overline{\psi }\psi ,`$ (19)
$`g_\omega \omega _\mu `$ $`=`$ $`c_\omega ^2\overline{\psi }\gamma _\mu \psi ,`$ (21)
$`T_\nu ^\mu `$ $`=`$ $`{\displaystyle \frac{1}{2}}\left(c_\omega ^2g_\omega \omega ^\alpha g_\omega \omega _\alpha +c_\sigma ^2\left[mm_{}\right]^2\right)\delta _\nu ^\mu i\overline{\mathrm{\Psi }}\gamma ^\mu _\nu \mathrm{\Psi },`$ (23)
where, for later convenience, we have introduced the notation $`c_\sigma ^2=(g_\sigma /m_\sigma )^2`$ and $`c_\omega ^2=(g_\omega /m_\omega )^2`$ and the Dirac effective mass $`m_{}`$, i.e.

$$g_\sigma \sigma =mm_{}.$$
(24)
Restricting to the zero-momentum frame of the neutrons leads to a set of algebraic equations for the $`\omega _\mu `$ field:

$`g_\omega \omega _0`$ $`=`$ $`c_\omega ^2\overline{\psi }\gamma _0\psi ,`$ (25)
$`g_\omega \omega _z`$ $`=`$ $`c_\omega ^2\overline{\psi }\gamma _z\psi .`$ (27)
The final equation is not needed in the case where both neutrons and protons have zero average momentum, since $`\omega _z`$ then vanishes by isotropy. In this case, the neutrons and protons have a common rest frame and $`\overline{\psi }\gamma ^0\psi =\psi ^{}\psi =n+p`$ where $`n`$ and $`p`$ are the baryon number densities of the neutrons and protons, respectively. The addition of the spatial velocity component complicates the solution of the problem considerably, in part because there is no longer a common rest frame for all the baryons. Each expectation value on the RHS of these equations involves an integration over the Fermi spheres of the particles, whose radii can be shown (c.f. the next section) to be $`k_n=\left(3\pi ^2n^0\right)^{1/3}`$ and $`k_p=\left(3\pi ^2p^0\right)^{1/3}`$, where $`n^0`$ ($`p^0`$) is the zero-component of the conserved neutron (proton) number density current $`n^\mu `$ ($`p^\mu `$). The proton Fermi surface is displaced by $`K_z\widehat{z}`$. We are interested in the case $`K_z<<k_n,k_p`$, but the expressions for general $`K_z`$ are not more complicated than the power series expansion.

Noting that

$$\overline{\psi }\left[\gamma ^\mu (i_\mu g_\omega \omega _\mu )m_{}\right]\psi =0$$
(28)
we find

$$\left(k^0+g_\omega \omega ^0\right)^2=\left(\stackrel{}{k}+g_\omega \omega ^z\widehat{z}\right)^2+m_{}^2,$$
(29)
where we have dropped expectation value brackets for the mean values of the fields. The energy $`\epsilon `$ of a baryon in a plane-wave state is given by

$$\epsilon (\stackrel{}{k})=E(\stackrel{}{k})g_\omega \omega ^0=\sqrt{\left(\stackrel{}{k}+g_\omega \omega ^z\widehat{z}\right)^2+m_{}^2}g_\omega \omega ^0.$$
(30)
Thus we see that $`\omega _0`$ contributes a constant shift, $`\omega _z`$ gives a preferred frame for the momenta, and $`\sigma `$ renormalizes the mass to the Dirac mass.

As an example of how the expectation values are evaluated, we give the scalar density (letting $`K_z=K`$, for ease of notation):

$`\overline{\psi }\psi `$ $`=`$ $`{\displaystyle \frac{1}{(2\pi )^3}}{\displaystyle _{occ}}d^3k{\displaystyle \frac{E}{m}}`$ (31)
$`=`$ $`{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}|<k_n}}d^3k{\displaystyle \frac{m_{}}{\sqrt{(\stackrel{}{k}+g_\omega \omega _z\widehat{z})^2+m_{}^2}}}+{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}K\widehat{z}|<k_p}}d^3k{\displaystyle \frac{m_{}}{\sqrt{(\stackrel{}{k}+g_\omega \omega _z\widehat{z})^2+m_{}^2}}}`$ (33)
$`=`$ $`{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}|<k_n}}d^3k{\displaystyle \frac{m_{}}{\sqrt{(\stackrel{}{k}+g_\omega \omega _z\widehat{z})^2+m_{}^2}}}+{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}|<k_p}}d^3k{\displaystyle \frac{m_{}}{\sqrt{(\stackrel{}{k}+g_\omega \omega _z\widehat{z}+K\widehat{z})^2+m_{}^2}}},`$ (35)
and the average four-velocity components of the baryons:

$`\overline{\psi }\gamma ^0\psi `$ $`=`$ $`{\displaystyle \frac{1}{(2\pi )^3}}{\displaystyle _{occ}}d^3k{\displaystyle \frac{E}{k_0}}`$ (36)
$`=`$ $`{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}|<k_n}}d^3k+{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}K\widehat{z}|<k_p}}d^3k`$ (38)
$`=`$ $`{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}|<k_n}}d^3k+{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}|<k_p}}d^3k,`$ (40)
$`\overline{\psi }\gamma ^z\psi `$ $`=`$ $`{\displaystyle \frac{1}{(2\pi )^3}}{\displaystyle _{occ}}d^3k{\displaystyle \frac{E}{k_z}}`$ (42)
$`=`$ $`{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}|<k_n}}d^3k{\displaystyle \frac{k^z+g_\omega \omega ^z}{\sqrt{(\stackrel{}{k}+g_\omega \omega ^z\widehat{z})^2+m_{}^2}}}+{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}K\widehat{z}|<k_p}}d^3k{\displaystyle \frac{k^z+g_\omega \omega ^z}{\sqrt{(\stackrel{}{k}+g_\omega \omega ^z\widehat{z})^2+m_{}^2}}}`$ (44)
$`=`$ $`{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}|<k_n}}d^3k{\displaystyle \frac{k^z+g_\omega \omega ^z}{\sqrt{(\stackrel{}{k}+g_\omega \omega ^z\widehat{z})^2+m_{})^2}}}+{\displaystyle \frac{2}{\left(2\pi \right)^3}}{\displaystyle _{|\stackrel{}{k}|<k_p}}d^3k{\displaystyle \frac{k^z+g_\omega \omega ^z+K}{\sqrt{(\stackrel{}{k}+g_\omega \omega ^z\widehat{z}+K\widehat{z})^2+m_{}^2}}}.`$ (46)
Thus we have reduced the problem to a set of nonlinear equations for the $`m_{}`$, $`\omega _0`$, and $`\omega _z`$ fields that must be solved numerically. This can be done for any set of the input parameters $`k_n`$, $`k_p`$, and $`K`$. The interaction and mass parameters for the effective fields have been determined from nuclear physics, and they are discussed further below. Once this is done, we still need expressions for the stress-energy tensor, which is the input for the Einstein equations.

Again, specializing to the zero-momentum frame of the neutrons, the only non-zero stress-energy tensor components are

$`T_0^0`$ $`=`$ $`{\displaystyle \frac{1}{2}}c_\omega ^2\left(\overline{\psi }\gamma ^0\psi ^2\overline{\psi }\gamma ^z\psi ^2\right){\displaystyle \frac{1}{2}}c_\sigma ^2\left(m^2m_{}^2\right)\overline{\psi }\gamma ^ik_i\psi ,`$ (47)
$`T_z^0`$ $`=`$ $`\overline{\psi }\gamma ^0k_z\psi ,`$ (49)
$`T_x^x`$ $`=`$ $`T_y^y={\displaystyle \frac{1}{2}}c_\omega ^2\left(\overline{\psi }\gamma ^0\psi ^2\overline{\psi }\gamma ^z\psi ^2\right){\displaystyle \frac{1}{2}}c_\sigma ^2\left(mm_{}\right)^2+\overline{\psi }\gamma ^xk_x\psi ,`$ (51)
$`T_z^z`$ $`=`$ $`{\displaystyle \frac{1}{2}}c_\omega ^2\left(\overline{\psi }\gamma ^0\psi ^2\overline{\psi }\gamma ^z\psi ^2\right){\displaystyle \frac{1}{2}}c_\sigma ^2\left(mm_{}\right)^2+\overline{\psi }\gamma ^zk_z\psi .`$ (53)
Some of the expressions have been simplified using the equations of motion.

Each component of $`T_\nu ^\mu `$ again involves an integration over the Fermi surfaces, but now in terms of completely known parameters. For example, to determine $`T_z^z`$, we need

$`\overline{\psi }\gamma ^zk_z\psi `$ $`=`$ $`{\displaystyle \frac{2}{(2\pi )^3}}{\displaystyle _{|\stackrel{}{k}|<k_n}}d^3kk^z{\displaystyle \frac{E}{k^z}}+{\displaystyle \frac{2}{(2\pi )^3}}{\displaystyle _{|\stackrel{}{k}K\widehat{z}|<k_p}}d^3kk^z{\displaystyle \frac{E}{k^z}}`$ (54)
$`=`$ $`{\displaystyle \frac{2}{(2\pi )^3}}{\displaystyle _{|\stackrel{}{k}|<k_n}}d^3kk_z\left[k_z+g_\omega \omega _z\right]\left[\left(\stackrel{}{k}+g_\omega \omega ^z\widehat{z}\right)^2+m_{}^2\right]^{1/2}+`$ (58)
$`{\displaystyle \frac{2}{(2\pi )^3}}{\displaystyle _{|\stackrel{}{k}|<k_p}}d^3k\left[k_z+K\right]\left[k_z+g_\omega \omega _z+K\right]\left[\left(\stackrel{}{k}+g_\omega \omega ^z\widehat{z}+K\widehat{z}\right)^2+m_{}^2\right]^{1/2},`$
and for $`T_z^0`$

$$\overline{\psi }\gamma ^0k_z\psi =\frac{2}{(2\pi )^3}_{|\stackrel{}{k}|<k_n}d^3kk_z+\frac{2}{(2\pi )^3}_{|\stackrel{}{k}|<k_p}d^3k\left(k_z+K\right)=\frac{k_p^3}{3\pi ^2}K.$$
(59)
The main result of this section is thus a well-defined prescription for producing the functions $`T_\nu ^\mu (k_n,k_p,K_z)`$. In the next section we take this prescription and produce from it the so-called master function, including entrainment, that is used in the general relativistic superfluid field equations.

## III General Relativistic Superfluid Formalism

The formalism to be used here, and motivation for it, has been described in great detail elsewhere C89 ; CL94 ; CL95 ; CL98a ; CL98b ; LSC98 ; CLL99 ; P00 ; AC01c ; GC02 , and so we will review only the highlights. The central quantity of the superfluid formalism is the master function $`\mathrm{\Lambda }`$. It depends on the three scalars $`n^2=n_\mu n^\mu `$, $`p^2=p_\mu p^\mu `$ and $`x^2=p_\mu n^\mu `$ that can be formed from the conserved neutron ($`n^\mu `$) and proton ($`p^\mu `$) number density currents. Furthermore, the master function is such that $`\mathrm{\Lambda }(n^2,p^2,x^2)`$ corresponds to the total thermodynamic energy density if the neutrons and protons flow together (as measured in the comoving frame). Once the master function is provided the stress-energy tensor is given by

$$T_\nu ^\mu =\mathrm{\Psi }\delta _\nu ^\mu +n^\mu \mu _\nu +p^\mu \chi _\nu ,$$
(60)
where

$$\mathrm{\Psi }=\mathrm{\Lambda }n^\rho \mu _\rho p^\rho \chi _\rho $$
(61)
is the generalized pressure, and

$`\mu _\nu `$ $`=`$ $`n_\nu +𝒜p_\nu ,`$ (62)
$`\chi _\nu `$ $`=`$ $`𝒜n_\nu +𝒞p_\nu ,`$ (63)
are the chemical potential covectors. We also have

$$𝒜=\frac{\mathrm{\Lambda }}{x^2},=2\frac{\mathrm{\Lambda }}{n^2},𝒞=2\frac{\mathrm{\Lambda }}{p^2}.$$
(64)
The momentum covectors $`\mu _\nu `$ and $`\chi _\nu `$ are dynamically, and thermodynamically, conjugate to $`n^\nu `$ and $`p^\nu `$ and their magnitudes are the chemical potentials of the neutrons and the protons, respectively. The two covectors also make manifest the entrainment effect; that is, we see that the momentum of one constituent ($`\mu _\nu `$, say) carries along some of the mass current of the other constituent ($`\mu _\nu `$ is a linear combination of $`n_\nu `$ and $`p_\nu `$). We can also see that there is no entrainment unless the master function depends on $`x^2`$.

The field equations for this system take the form of two conservation equations for the neutrons and protons, i.e.

$$_\mu n^\mu =0,_\mu p^\mu =0,$$
(65)
which is a reasonable approximation given that the weak interaction timescale is much longer than the dynamical timescale of neutron stars for small amplitude deviations from equilibrium E88 , and two Euler equations, i.e.

$$n^\mu _{[\mu }\mu _{\nu ]}=0,p^\mu _{[\mu }\chi _{\nu ]}=0,$$
(66)
where the square braces means antisymmetrization of the enclosed indices.

### III.1 Extracting the Master Function from the Mean Field Results

The two scales that enter this problem are the microscopic, on the scale of the nucleons, and the mesoscopic where one speaks in terms of the two interpenetrating superfluids. The fundamental “particles” at the fluid level are the fluid elements which contain, say, an Avogadro’s number worth of nucleons. The connection between the micro- and meso-scopic levels is via the averaged stress-energy components calculated earlier. Consider a fluid element deep in the core of the neutron star and orient the local coordinate frame in such a way that the z-axis of the frame is in the same direction as the proton momentum with respect to the neutrons. As shown just below, a unique combination of the averaged stress-energy components determined via the mean field theory will yield the master function. As this quantity is a scalar, the functional relationship we obtain between $`\mathrm{\Lambda }`$ and the two particle number densities and the relative velocity of the protons with respect to the neutrons can then be applied anywhere in the star.

The key idea is to use the (local) relationship

$$T_\nu ^\mu =\mathrm{\Psi }\delta _\nu ^\mu +n^\mu \mu _\nu +p^\mu \chi _\nu $$
(67)
to obtain $`\mathrm{\Lambda }`$. In the perfect fluid case, the identification is made immediate by the fact that there is only one four-velocity $`u^\mu `$ for the system, and hence a preferred rest-frame for the particles. The local energy density of the fluid is thus uniquely obtained from $`\mathrm{\Lambda }=T_\nu ^\mu u^\nu u_\mu `$. In the superfluid case there are two four-velocities, and thus no preferred rest-frame. Fortunately, we can still obtain $`\mathrm{\Lambda }`$ in a unique, and covariant, way, by using the trace $`TT_\mu ^\mu `$ and the three scalars that can be formed from contracting $`T_\nu ^\mu `$ with $`n^\mu `$ and $`p^\mu `$, i.e. $`T_\nu ^\mu n_\mu n^\nu `$, $`T_\nu ^\mu n_\mu p^\nu `$, and $`T_\nu ^\mu p_\mu p^\nu `$. We thus find that $`\mathrm{\Lambda }`$ is given by

$$\mathrm{\Lambda }=\frac{1}{2}T+\frac{3}{2}\left(x^4n^2p^2\right)^1\left(n^2p^2\left[\frac{1}{n^2}n^\mu n^\nu +\frac{1}{p^2}p^\mu p^\nu \right]x^2\left[n^\mu p^\nu +p^\mu n^\nu \right]\right)T_{\mu \nu },$$
(68)
and the generalized pressure is

$$\mathrm{\Psi }=\frac{1}{3}\left(T\mathrm{\Lambda }\right).$$
(69)
In like manner we find that

$`𝒜`$ $`=`$ $`\left(n_\mu p^\nu T_\nu ^\mu +x^2\mathrm{\Lambda }\right)/\left(x^4n^2p^2\right),`$ (70)
$``$ $`=`$ $`\left(p_\mu p^\nu T_\nu ^\mu +p^2\mathrm{\Lambda }\right)/\left(x^4n^2p^2\right),`$ (72)
$`𝒞`$ $`=`$ $`\left(n_\mu n^\nu T_\nu ^\mu +n^2\mathrm{\Lambda }\right)/\left(x^4n^2p^2\right).`$ (74)
One other necessary component of uniting the mean field theory with the superfluid formalism is to relate (locally) $`n^\mu `$ and $`p^\mu `$ to the mean particle flux of the neutrons and protons; i.e.

$`n^\mu `$ $``$ $`nu_n^\mu =\overline{\psi }_n\gamma ^\mu \psi _n,`$ (75)
$`p^\mu `$ $``$ $`pu_p^\mu =\overline{\psi }_p\gamma ^\mu \psi _p,`$ (77)
where $`\psi _n`$ and $`\psi _p`$ are the neutron and proton, respectively, components of the Dirac spinor $`\psi `$. Recall again that we have arranged that the average neutron and proton particle fluxes are in the z-direction. Thus, the unit vectors have only two components:

$`u_n^\mu `$ $`=`$ $`\{u_n^0,0,0,u_n^3\},u_n^0=\sqrt{1+\left(u_n^3\right)^2},`$ (78)
$`u_p^\mu `$ $`=`$ $`\{u_p^0,0,0,u_p^3\},u_p^0=\sqrt{1+\left(u_p^3\right)^2}.`$ (80)
It thus follows that

$$\mathrm{\Lambda }=T_0^0+T_z^zT_x^x.$$
(81)
We will use cylindrical coordinates and define $`\varphi _z=g_\omega \omega _z`$ so that

$`m_{}`$ $`=`$ $`m{\displaystyle \frac{c_\sigma ^2}{2\pi ^2}}m_{}({\displaystyle _{k_n}^{k_n}}dk_z[k_n^2+\varphi _z^2+m_{}^2+2\varphi _zk_z]^{1/2}+{\displaystyle _{k_p}^{k_p}}dk_z[k_p^2+(\varphi _z+K)^2+m_{}^2+2(\varphi _z+K)k_z]^{1/2}`$ (84)
$`{\displaystyle _{k_n}^{k_n}}dk_z[(k_z+\varphi _z)^2+m_{}^2]^{1/2}{\displaystyle _{k_p}^{k_p}}dk_z[(k_z+\varphi _z+K)^2+m_{}^2]^{1/2}),`$
$`n^0`$ $`=`$ $`{\displaystyle \frac{1}{3\pi ^2}}k_n^3,p^0={\displaystyle \frac{1}{3\pi ^2}}k_p^3,`$ (86)
$`n^z`$ $`=`$ $`{\displaystyle \frac{1}{2\pi ^2}}{\displaystyle _{k_n}^{k_n}}𝑑k_z\left(k_z+\varphi _z\right)\left(\left[k_n^2+m_{}^2+\varphi _z^2+2\varphi _zk_z\right]^{1/2}\left[\left(k_z+\varphi _z\right)^2+m_{}^2\right]^{1/2}\right),`$ (88)
$`p^z`$ $`=`$ $`{\displaystyle \frac{1}{2\pi ^2}}{\displaystyle _{k_p}^{k_p}}dk_z(k_z+\varphi _z+K)([k_n^2+m_{}^2+(\varphi _z+K)^2+2(\varphi _z+K)k_z]^{1/2}`$ (92)
$`[(k_z+\varphi _z+K)^2+m_{}^2]^{1/2}).`$
It is to be remarked that even though only the protons are given spatial momentum $`K`$ the neutron four-velocity nevertheless has a non-zero spatial component. This, in fact, is a signature of the entrainment effect, which is a momentum induced in one of the fluids will cause part of the other fluid to flow.

Of primary importance to the fluid equations are the $`𝒜`$, $``$, and $`𝒞`$ coefficients. We could, in principle, use Eq. (77) to express $`(n^2,p^2,x^2)`$ in terms of $`(k_n,k_p,K)`$, but practically speaking this is not possible. Fortunately we see from Eq. (74) that we can construct these coefficients algebraically from the mean field values of the stress-energy tensor components. Thus, when it comes to the numerical work, we use the set $`(k_n,k_p,K)`$ as the independent variables. Note that because the master function is a scalar, it must be invariant if $`KK`$, and is thus an even function of $`K`$.

### III.2 Equilibrium models

The equilibrium configurations are spherically symmetric and static, so the metric can be written in the Schwarzschild form

$$\mathrm{d}s^2=e^{\nu (r)}\mathrm{d}t^2+e^{\lambda (r)}\mathrm{d}r^2+r^2\left(\mathrm{d}\theta ^2+\mathrm{sin}^2\theta \mathrm{d}\varphi ^2\right).$$
(93)
The two metric coefficients are determined from two Einstein equations, which are written as

$$\lambda ^{}=\frac{1e^\lambda }{r}8\pi re^\lambda \mathrm{\Lambda }|_0,\nu ^{}=\frac{1e^\lambda }{r}+8\pi re^\lambda \mathrm{\Psi }|_0.$$
(94)
The equations that determine the radial profiles of $`n(r)`$ and $`p(r)`$ have been derived by Comer et al. CLL99  and they are

$$𝒜_0^0|_0p^{}+_0^0|_0n^{}+\frac{1}{2}\mu |_0\nu ^{}=0,𝒞_0^0|_0p^{}+𝒜_0^0|_0n^{}+\frac{1}{2}\chi |_0\nu ^{}=0,$$
(95)
where

$`𝒜_0^0`$ $`=`$ $`𝒜+2{\displaystyle \frac{}{p^2}}np+2{\displaystyle \frac{𝒜}{n^2}}n^2+2{\displaystyle \frac{𝒜}{p^2}}p^2+{\displaystyle \frac{𝒜}{x^2}}pn,`$ (96)
$`_0^0`$ $`=`$ $`+2{\displaystyle \frac{}{n^2}}n^2+4{\displaystyle \frac{𝒜}{n^2}}np+{\displaystyle \frac{𝒜}{x^2}}p^2,`$ (98)
$`𝒞_0^0`$ $`=`$ $`𝒞+2{\displaystyle \frac{𝒞}{p^2}}p^2+4{\displaystyle \frac{𝒜}{p^2}}np+{\displaystyle \frac{A}{x^2}}n^2,`$ (100)
and $`\mathrm{\Lambda }|_0`$, $`\mathrm{\Psi }|_0`$, $`\mu |_0`$, and $`\chi |_0`$ are given in the appendix. A zero subscript means that after the partial derivatives are taken, then one takes the limit $`K0`$.

Of course, since the variables that are more suited to the mean field theory are the two Fermi wavenumbers, $`k_n`$ and $`k_p`$, we replace everywhere $`n=k_n^3/3\pi ^2`$ and $`p=k_p^3/3\pi ^2`$ and solve for the wavenumbers instead. We have also found a more convenient way of determining the Dirac effective mass $`m_{}|_0(k_n,k_p)`$ and that is we have turned the transcendental algebraic relation in Eq. (153) of the appendix into a differential equation via

$$m_{}^{}|_0=\frac{m_{}}{k_n}|_0k_n^{}+\frac{m_{}}{k_p}|_0k_p^{},$$
(101)
where $`k_n^{}`$ and $`k_p^{}`$ are obtained from Eq. (95).

The “boundary” conditions that must be imposed include a set at the center and another at the surface of the star. Demanding a non-singular behavior at the center of the star imposes that $`\lambda (0)=0`$, and consequently that $`\lambda ^{}(0)`$ and $`\nu ^{}(0)`$ must also vanish. This and Eq. (95) imply that $`k_n^{}(0)`$ and $`k_p^{}(0)`$ have to vanish as well. A smooth joining of the interior spacetime to a Schwarzschild vacuum exterior at the surface of the star, i.e. $`r=R`$, implies that the total mass $`M`$ of the system is given by

$$M=4\pi _0^Rdrr^2\mathrm{\Lambda }|_0(r)$$
(102)
and that $`\mathrm{\Psi }|_0(R)=0`$.

### III.3 The low velocity limit for fluid elements

Two immediate applications of the formalism developed here are to model slowly rotating configurations AC01c  and linearized perturbations, or quasinormal modes CLL99 ; ACL02 . In both cases the fluid element velocities are small in the sense that they are typically only a few percent of the speed of light. The net effect is that these applications require only the first few terms from an expansion of the master function in terms of the entrainment parameter (i.e. $`x^2`$ in the canonical formulation, and $`K^2`$ in what we present here). Such an expansion has been described in AC01a ; ACL02 , and thus only the highlights will be reproduced here. It should be noted, however, that if one wanted to model rapidly rotating superfluid neutron stars PNC02 , say, then the expansion to be described below will be inappropriate.

For a region within the fluid that is small enough that the gravitational field does not change appreciably across it one can show that

$$x^2=np\left(\frac{1\stackrel{}{v}_n\stackrel{}{v}_p/c^2}{\sqrt{1(v_n/c)^2}\sqrt{1(v_p/c)^2}}\right).$$
(103)
If it is the case that the individual three-velocities $`\stackrel{}{v}_{n,p}`$ are small with respect to the speed of light, i.e.

$$\frac{v_{n,p}}{c}<<1,$$
(104)
then it will be true that $`x^2np`$ to leading order in the ratios $`u_n/c`$ and $`v_p/c`$. Thus, an appropriate expansion of the equation of state is

$$\mathrm{\Lambda }(n^2,p^2,x^2)=\underset{i=0}{\overset{\mathrm{}}{}}\lambda _i(n^2,p^2)\left(x^2np\right)^i,$$
(105)
since $`x^2np`$ is small with respect to $`np`$. In this case the $`𝒜`$, $`𝒜_0^0`$, etc. coefficients that appear in the field equations can be written as

$`𝒜`$ $`=`$ $`{\displaystyle \underset{i=1}{\overset{\mathrm{}}{}}}i\lambda _i(n^2,p^2)\left(x^2np\right)^{i1},`$ (106)
$``$ $`=`$ $`{\displaystyle \frac{1}{n}}{\displaystyle \frac{\lambda _0}{n}}{\displaystyle \frac{p}{n}}𝒜{\displaystyle \frac{1}{n}}{\displaystyle \underset{i=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\lambda _i}{n}}\left(x^2np\right)^i,`$ (108)
$`𝒞`$ $`=`$ $`{\displaystyle \frac{1}{p}}{\displaystyle \frac{\lambda _0}{p}}{\displaystyle \frac{n}{p}}𝒜{\displaystyle \frac{1}{p}}{\displaystyle \underset{i=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{\lambda _i}{p}}\left(x^2np\right)^i,`$ (110)
$`𝒜_0^0`$ $`=`$ $`{\displaystyle \frac{^2\lambda _0}{pn}}{\displaystyle \underset{i=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{^2\lambda _i}{pn}}\left(x^2np\right)^i,`$ (112)
$`_0^0`$ $`=`$ $`{\displaystyle \frac{^2\lambda _0}{n^2}}{\displaystyle \underset{i=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{^2\lambda _i}{n^2}}\left(x^2np\right)^i,`$ (114)
$`𝒞_0^0`$ $`=`$ $`{\displaystyle \frac{^2\lambda _0}{p^2}}{\displaystyle \underset{i=1}{\overset{\mathrm{}}{}}}{\displaystyle \frac{^2\lambda _i}{p^2}}\left(x^2np\right)^i.`$ (116)
For quasinormal mode and slow-rotation calculations, each of the coefficients are evaluated on the background, so that $`x^2=np`$, and thus only the first few $`\lambda _i`$ are needed. In fact, one needs to retain only $`\lambda _0`$ and $`\lambda _1`$, where the latter contains the information concerning the entrainment effect. Some details are given in the appendix, and the final results are

$`𝒜|_0`$ $`=`$ $`c_\omega ^2+{\displaystyle \frac{c_\omega ^2}{5\mu ^2|_0}}\left(2k_p^2{\displaystyle \frac{\sqrt{k_n^2+m_{}^2|_0}}{\sqrt{k_p^2+m_{}^2|_0}}}+{\displaystyle \frac{c_\omega ^2}{3\pi ^2}}\left[{\displaystyle \frac{k_n^2k_p^3}{\sqrt{k_n^2+m_{}^2|_0}}}+{\displaystyle \frac{k_p^2k_n^3}{\sqrt{k_p^2+m_{}^2|_0}}}\right]\right)+{\displaystyle \frac{3\pi ^2k_p^2}{5\mu ^2|_0k_n^3}}{\displaystyle \frac{k_n^2+m_{}^2|_0}{\sqrt{k_p^2+m_{}^2|_0}}},`$ (117)
$`|_0`$ $`=`$ $`{\displaystyle \frac{3\pi ^2\mu |_0}{k_n^3}}c_\omega ^2{\displaystyle \frac{k_p^3}{k_n^3}}{\displaystyle \frac{c_\omega ^2k_p^3}{5\mu ^2|_0k_n^3}}\left(2k_p^2{\displaystyle \frac{\sqrt{k_n^2+m_{}^2|_0}}{\sqrt{k_p^2+m_{}^2|_0}}}+{\displaystyle \frac{c_\omega ^2}{3\pi ^2}}\left[{\displaystyle \frac{k_n^2k_p^3}{\sqrt{k_n^2+m_{}^2|_0}}}+{\displaystyle \frac{k_p^2k_n^3}{\sqrt{k_p^2+m_{}^2|_0}}}\right]\right)`$ (121)
$`{\displaystyle \frac{3\pi ^2k_p^5}{5\mu ^2|_0k_n^6}}{\displaystyle \frac{k_n^2+m_{}^2|_0}{\sqrt{k_p^2+m_{}^2|_0}}},`$
$`𝒞|_0`$ $`=`$ $`{\displaystyle \frac{3\pi ^2\chi |_0}{k_p^3}}c_\omega ^2{\displaystyle \frac{k_n^3}{k_p^3}}{\displaystyle \frac{c_\omega ^2k_n^3}{5\mu ^2|_0k_p^3}}\left(2k_p^2{\displaystyle \frac{\sqrt{k_n^2+m_{}^2|_0}}{\sqrt{k_p^2+m_{}^2|_0}}}+{\displaystyle \frac{c_\omega ^2}{3\pi ^2}}\left[{\displaystyle \frac{k_n^2k_p^3}{\sqrt{k_n^2+m_{}^2|_0}}}+{\displaystyle \frac{k_p^2k_n^3}{\sqrt{k_p^2+m_{}^2|_0}}}\right]\right)`$ (125)
$`{\displaystyle \frac{3\pi ^2}{5\mu ^2|_0k_p}}{\displaystyle \frac{k_n^2+m_{}^2|_0}{\sqrt{k_p^2+m_{}^2|_0}}}+{\displaystyle \frac{3\pi ^2}{k_p^3}}\sqrt{k_p^2+m_e^2},`$
$`𝒜_0^0|_0`$ $`=`$ $`{\displaystyle \frac{\pi ^4}{k_n^2k_p^2}}{\displaystyle \frac{^2\mathrm{\Lambda }}{k_pk_n}}|_0=c_\omega ^2+{\displaystyle \frac{\pi ^2}{k_p^2}}{\displaystyle \frac{m_{}|_0\frac{m_{}}{k_p}|_0}{\sqrt{k_n^2+m_{}^2|_0}}},`$ (127)
$`_0^0|_0`$ $`=`$ $`{\displaystyle \frac{\pi ^4}{k_n^5}}\left(2{\displaystyle \frac{\mathrm{\Lambda }}{k_n}}|_0k_n{\displaystyle \frac{^2\mathrm{\Lambda }}{k_n^2}}|_0\right)=c_\omega ^2+{\displaystyle \frac{\pi ^2}{k_n^2}}{\displaystyle \frac{k_n+m_{}|_0\frac{m_{}}{k_n}|_0}{\sqrt{k_n^2+m_{}^2|_0}}},`$ (129)
$`𝒞_0^0|_0`$ $`=`$ $`{\displaystyle \frac{\pi ^4}{k_p^5}}\left(2{\displaystyle \frac{\mathrm{\Lambda }}{k_p}}|_0k_p{\displaystyle \frac{^2\mathrm{\Lambda }}{k_p^2}}|_0\right)=c_\omega ^2+{\displaystyle \frac{\pi ^2}{k_p^2}}{\displaystyle \frac{k_p+m_{}|_0\frac{m_{}}{k_p}|_0}{\sqrt{k_p^2+m_{}^2|_0}}}+{\displaystyle \frac{\pi ^2}{k_p}}{\displaystyle \frac{1}{\sqrt{k_p^2+m_e^2}}},`$ (131)
where $`m_{}/k_n|_0`$ and $`m_{}/k_p|_0`$ can be found in the appendix and $`m_e=m/1836`$. For reasons to be discussed below, we have included contributions due to a normal fluid of highly degenerate electrons.

### III.4 Equilibrium configurations

We now use our model to construct static and spherically symmetric configurations. A priori there are two input parameters, which are the neutron and proton Fermi wavenumbers at the center of the star. However, we can reduce this to just the neutron wavenumber by imposing the condition of chemical equilibrium between the nucleons. In order to have a chemical equilibrium that is believed to be representative of neutron stars (i.e. proton fractions $`x_p=p/(n+p)0.1`$), we have added to the master function a term (see, for instance, ST83 ; ACP02b ) that accounts for a highly degenerate gas of relativistic leptons (in our case, just electrons). Fig. 1 gives the mass $`M`$ as a function of the central neutron number density $`n(0)`$. We see the typical behaviour of general relativistic neutron stars, and that is a maximum value for the mass. Beyond this maximum, the stars will be in unstable equilibria. As canonical models of superfluid neutron stars, we have chosen configurations that are near to the maximum mass, but on the stable branch of the curves (cf. Table 1).

In several earlier studies LM00 ; PCA02 ; ACL02 , a parameterized model for entrainment has been used that is based on the Newtonian calculations of Borumand et al. BJK96 . The parameter, $`\epsilon _{mom}`$, in this model can be shown to be related to our coefficient $`𝒜|_0`$ via

$$\epsilon _{mom}=\frac{𝒜|_0p}{m𝒜|_0\left(n+p\right)}.$$
(132)
In the previous studies the “physical” range has been taken to be $`0.04<\epsilon _{mom}<0.2`$. But already with this formula we find that $`\epsilon _{mom}`$ has a singularity whenever $`m=𝒜|_0(n+p)`$. And in Fig. 2 we see that indeed $`\epsilon _{mom}`$ has a discontinuity in both curves (which are for Models I and II of Table 1). This is perhaps the most significant short term result of this work, and hence taking $`\epsilon _{mom}`$ to be a constant must be considered a non-viable option.

An alternative parameterization is that of Prix et al. PCA02 . They point out that there is some ambiguity in what is meant by the nucleon effective masses (i.e. the Landau, as opposed to the Dirac, effective masses G97 ), which can be traced to whether one chooses to define these masses with respect to the zero-momentum or zero-velocity frame of the nucleons. Their parameterization makes use of the zero-velocity frame. In this case one finds

$$\epsilon _{vel}=\frac{𝒜|_0n}{m}.$$
(133)
Immediately we see that this formulation does not have an obvious singularity anywhere, and this is reflected in Fig. 3 which contains plots of the radial profile of $`\epsilon _{vel}`$ for Models I and II of Table 1. The “physical” range has been taken to be $`0.4<\epsilon _{vel}<0.7`$. We see from the figure that this is the range in the outer portion of the star, but in the core the entrainment stretches well out of this range.

Also used in Prix et al. PCA02  is the so-called symmetry energy parameter

$$\sigma =\frac{\mu /p}{\chi /n}$$
(134)
that can be related to terms PAL88  in the equation of state that tend to force an equal number of protons and neutrons (as in most nuclei). For the relativistic mean field model used here, it is not difficult to show that $`\sigma =1`$ (which is consistent with the range of values used by PCA02 ).

## IV Slow Rotation Configurations: The Frame-Dragging

The key distinguishing feature of slowly rotating superfluid neutron stars is that the neutrons can rotate at rates different from that of the protons. The slow-rotation approximation is valid when the angular velocities are small enough that the fractional changes in pressure, energy density, and gravitational field due to the rotation are all relatively small. This translates into the inequalities (cf. H67 ; AC01c )

$$\mathrm{\Omega }_n^2\mathrm{or}\mathrm{\Omega }_p^2\mathrm{or}\mathrm{\Omega }_n\mathrm{\Omega }_p<<\left(\frac{c}{R}\right)^2\frac{GM}{Rc^2},$$
(135)
where the speed of light $`c`$ and Newton’s constant $`G`$ have been restored, and $`R`$ and $`M`$ are the radius and mass, respectively, of the non-rotating configuration. Since $`GM/c^2R<1`$, we also see that

$$\mathrm{\Omega }_nR<<c\mathrm{and}\mathrm{\Omega }_pR<<c,$$
(136)
and thus the slow-rotation approximation ought to be useful for most astrophysical neutron stars. In fact a comparison AC01c ; PCA02  of the above conditions to empirical estimates for the Kepler frequency (i.e. the rotation rate at which mass-shedding sets in at the equator) that can been obtained from calculations using realistic supranuclear equations of state reveals that even the fastest observed pulsars can be classified as slowly rotating.

The only quantities that contain terms linear in the angular velocities are the metric coefficient $`\omega (r)`$, that represents the dragging of inertial frames, and the fluid four-velocities. All other effects due to rotation enter at the second-order in the angular velocities. It is useful to define

$$\stackrel{~}{L}_n=\omega \mathrm{\Omega }_n,\stackrel{~}{L}_p=\omega \mathrm{\Omega }_p.$$
(137)
Up to an overall minus sign, these represent rotation frequencies as perceived by local zero-angular momentum observers. The Einstein equation that determines the frame-dragging has been shown to be AC01c

$$\frac{1}{r^4}\left(r^4e^{(\lambda +\nu )/2}\stackrel{~}{L}_p^{}\right)^{}16\pi e^{(\lambda \nu )/2}\left(\mathrm{\Psi }|_0\mathrm{\Lambda }|_0\right)\stackrel{~}{L}_p=16\pi e^{(\lambda \nu )/2}\mu |_0n\left(\mathrm{\Omega }_n\mathrm{\Omega }_p\right).$$
(138)
It is of the same form as that obtained by Hartle H67  except for the non-zero source term on the right-hand-side.

Exterior to the star, there is vacuum, and so the solution for the frame-dragging is the same as that considered by Hartle H67 , i.e.

$$\omega (r)=\frac{2J}{r^3}.$$
(139)
Assuming that the frame-dragging is continuous at the surface of the star, then

$$J=\frac{8\pi }{3}_0^Rdrr^4e^{(\lambda \nu )/2}\left[\mu |_0n\stackrel{~}{L}_n+\chi |_0p\stackrel{~}{L}_p\right],$$
(140)
where $`J`$ is the total angular momentum. Andersson and Comer AC01c  have furthermore shown that the neutron total angular momentum is

$$J_n=\frac{8\pi }{3}_0^Rdrr^4e^{(\lambda \nu )/2}\left[\mu |_0n\stackrel{~}{L}_n+𝒜|_0np\left(\mathrm{\Omega }_n\mathrm{\Omega }_p\right)\right]$$
(141)
and

$$J_p=\frac{8\pi }{3}_0^Rdrr^4e^{(\lambda \nu )/2}\left[\chi |_0p\stackrel{~}{L}_p+𝒜|_0np\left(\mathrm{\Omega }_p\mathrm{\Omega }_n\right)\right]$$
(142)
for the proton total angular momentum, from which it follows that

$$J=J_n+J_p.$$
(143)
Any solution of Eq. (138) for the frame-dragging is to be such that the interior matches smoothly onto the known vacuum solution in Eq. (139). This means that we must have, for instance,

$$\stackrel{~}{L}_p(R)=\mathrm{\Omega }_p+\frac{2J}{R^3}.$$
(144)
We can easily see that $`\stackrel{~}{L}_p`$ and its derivative are smooth provided that we have

$$\stackrel{~}{L}_p(R)=\mathrm{\Omega }_p\frac{R}{3}\frac{d\stackrel{~}{L}_p}{dr}|_{r=R}.$$
(145)
Having obtained a value for $`\stackrel{~}{L}_p(0)`$ that satisfies Eq. (145), an acceptable solution is in hand, and we can thus determine the angular momentum of the configuration from Eq. (144).

In Fig. 4 we have plots of the radial profile of the frame-dragging for Model I for a range of values of the ratio $`\mathrm{\Omega }_n/\mathrm{\Omega }_p`$. For the values considered in the left panel we see that the frame-dragging is much like that of an ordinary one-fluid star, and is consistent with solutions obtained by Andersson and Comer AC01c . For the negative ratios, we see that the frame-dragging is negative but increases monotonically towards zero. This is the behaviour we should expect, since the bulk of the matter is simply rotating the opposite way. There is some asymmetry between the negative and positive ratios, but that is due to the small number of protons that rotate oppositely to the neutrons when the ratio is negative. In the right panel, we examine the solutions near to a ratio of zero. The frame-dragging is no longer monotonic and actually becomes negative inside the star. An explanation of this can be understood as follows: in the interior the protons carry most of the angular momentum and thus have the largest impact on the frame-dragging, but further away from the center, the much larger mass contained in the neutrons begins to dominate AC01c .

Fig. 5 considers the same range for the ratio of the angular speeds, by showing how the total angular momentum $`J`$, and the neutron and proton angular momenta, $`J_n`$ and $`J_p`$ respectively, vary as $`\mathrm{\Omega }_n/\mathrm{\Omega }_p`$ is changed. As one might expect, when the ratio becomes greater than one, then the angular momentum in the neutrons is significantly greater than that of the protons. Likewise, as the ratio becomes smaller we find that the protons can dominate. Something of a surprise is the extreme right-hand-side of the curves, where $`J_p`$ actually becomes negative, and yet the angular speeds of the neutrons and protons have the same sign. We can explain this behaviour as a purely general relativistic effect that is intimately connected with the frame-dragging. With respect to infinity, particles that are rotating around at the same rate as the local inertial frames are found to have zero angular momentum. Thus, those particles that would be lagging behind the frames, even though their angular trajectories would be in the same direction as the frames, will nevertheless have negative angular momentum. Finally, one other feature is the configuration where $`J_n=J_p`$. In this case, the angular speeds are not equal, nor are the total neutron and proton particle numbers equal, and yet the angular momenta of both fluids are the same.

## V Conclusions

We have developed a formalism that uses relativistic mean field theory for supra-nuclear density matter that can be applied to general relativistic superfluid neutron stars. In this formalism we have also allowed for the entrainment effect between the various superfluid species. We have also shown how to use our formalism in the relativistic superfluid field equations that have been recently developed for modelling slowly rotating equilbirium configurations AC01c , as well as the linearized oscillations CLL99 ; ACL02 . Our results indicate that parameterized models of entrainment that use a zero-momentum rest-frame definition for nucleon effective masses have a singularity, whereas those based on zero-velocity rest-frames do not.

Our results should find a wide range of applications, not the least of which is to understand better the role of entrainment in the superfluid modes of oscillation (e.g. the avoided crossings described by AC01a ) and subsequent imprints AC01b ; GC02  that may be left in neutron star gravitational waves (emitted, for instance, during glitches). Applications planned for the near future will include numerical studies of rapidly rotating superfluid neutron stars (using an adaptation of the very accurate LORENE code PNC02 ) and continued research on the newly discovered two-stream instability ACP02a ; ACP02b , which has been proposed as a trigger mechanism for glitches in pulsars. For the rapid rotation calculations one must necessarily imploy the full formalism discussed here in the sense that $`K`$ will no longer be kept small, since the LORENE code is specifically designed to accurately handle relative velocities of the neutrons with respect to the protons that approach the speed of light. And for the two-stream instability entrainment provides one of the main couplings between the two fluids.

###### Acknowledgements.

We like to thank N. Andersson and R. Prix for useful input at various stages of this work. GLC gratefully acknowledges partial support from NSF Gravitational Theory, Grant No. PHYS-0140138, a Saint Louis University SLU2000 Faculty Research Leave award, and EPSRC grant GR/R52169/01 in the UK. RJ gratefully acknowledges partial support from the NSF Materials Theory Program, Grant No. DMR-0081039.

## Appendix: Limiting Forms

The slow-rotation approximation is such that only terms up to and including $`𝒪(\mathrm{\Omega }_{n,p}^2)`$ are required. This translates into keeping only those terms in the mean field theory up to and including $`𝒪(K^2)`$. This is because those quantities like $`m_{}`$ and $`\mathrm{\Lambda }`$ are scalars, and can only depend on terms that are even in $`K`$. Likewise, those quantities that are like vectors, e.g. $`\omega _z`$, can only depend on terms that are odd in $`K`$. Because $`m_{}`$ and $`\omega _z`$ are known only implicitly, we determine their expansion coefficients by assuming they take the form

$`\varphi _z`$ $`=`$ $`{\displaystyle \frac{\varphi _z}{K}}|_0K,`$ (146)
$`m_{}`$ $`=`$ $`m_{}|_0+{\displaystyle \frac{m_{}}{K^2}}|_0K^2,`$ (148)
where

$`m_{}|_0`$ $`=`$ $`m_{}(k_n,k_p,0)`$ (149)
$`=`$ $`mm_{}|_0{\displaystyle \frac{c_\sigma ^2}{2\pi ^2}}(k_n\sqrt{k_n^2+m_{}^2|_0}+k_p\sqrt{k_p^2+m_{}^2|_0}+{\displaystyle \frac{1}{2}}m_{}^2|_0\mathrm{ln}\left[{\displaystyle \frac{k_n+\sqrt{k_n^2+m_{}^2|_0}}{k_n+\sqrt{k_n^2+m_{}^2|_0}}}\right]+`$ (153)
$`{\displaystyle \frac{1}{2}}m_{}^2|_0\mathrm{ln}\left[{\displaystyle \frac{k_p+\sqrt{k_p^2+m_{}^2|_0}}{k_p+\sqrt{k_p^2+m_{}^2|_0}}}\right]).`$
By inserting Eq. (148) into Eq. (92), and expanding and keeping terms to the appropriate orders, we find

$`{\displaystyle \frac{\varphi _z}{K}}|_0`$ $`=`$ $`{\displaystyle \frac{c_\omega ^2}{3\pi ^2}}{\displaystyle \frac{k_p^3}{\sqrt{k_p^2+m_{}^2|_0}}}\left(1+{\displaystyle \frac{c_\omega ^2}{3\pi ^2}}\left[{\displaystyle \frac{k_n^3}{\sqrt{k_n^2+m_{}^2|_0}}}+{\displaystyle \frac{k_p^3}{\sqrt{k_p^2+m_{}^2|_0}}}\right]\right)^1,`$ (154)
$`{\displaystyle \frac{m_{}}{k_n}}|_0`$ $`=`$ $`{\displaystyle \frac{c_\sigma ^2}{\pi ^2}}{\displaystyle \frac{m_{}|_0k_n^2}{\sqrt{k_n^2+m_{}^2|_0}}}\left({\displaystyle \frac{3m2m_{}|_0}{m_{}|_0}}{\displaystyle \frac{c_\sigma ^2}{\pi ^2}}\left[{\displaystyle \frac{k_n^3}{\sqrt{k_n^2+m_{}^2|_0}}}+{\displaystyle \frac{k_p^3}{\sqrt{k_p^2+m_{}^2|_0}}}\right]\right)^1,`$ (156)
$`{\displaystyle \frac{m_{}}{k_p}}|_0`$ $`=`$ $`{\displaystyle \frac{c_\sigma ^2}{\pi ^2}}{\displaystyle \frac{m_{}|_0k_p^2}{\sqrt{k_p^2+m_{}^2|_0}}}\left({\displaystyle \frac{3m2m_{}|_0}{m_{}|_0}}{\displaystyle \frac{c_\sigma ^2}{\pi ^2}}\left[{\displaystyle \frac{k_n^3}{\sqrt{k_n^2+m_{}^2|_0}}}+{\displaystyle \frac{k_p^3}{\sqrt{k_p^2+m_{}^2|_0}}}\right]\right)^1,`$ (158)
$`n^z`$ $`=`$ $`{\displaystyle \frac{1}{3\pi ^2}}{\displaystyle \frac{k_n^3}{\sqrt{k_n^2+m_{}^2|_0}}}{\displaystyle \frac{\varphi _z}{K}}|_0K,`$ (160)
$`p^z`$ $`=`$ $`{\displaystyle \frac{1}{3\pi ^2}}{\displaystyle \frac{k_p^3}{\sqrt{k_p^2+m_{}^2|_0}}}\left({\displaystyle \frac{\varphi _z}{K}}|_0+1\right)K.`$ (162)
We note that the coefficient $`m_{}/K^2|_0`$ cancels everywhere, which is why it is not written here. Also, we find

$`\mathrm{\Lambda }|_0`$ $`=`$ $`{\displaystyle \frac{1}{4\pi ^2}}\left(k_n^3\sqrt{k_n^2+m_{}^2|_0}+k_p^3\sqrt{k_p^2+m_{}^2|_0}\right){\displaystyle \frac{1}{2}}m_\omega ^2\omega _0^2{\displaystyle \frac{1}{4}}c_\sigma ^2\left(2mm_{}|_0\right)\left(mm_{}|_0\right)`$ (165)
$`{\displaystyle \frac{1}{8\pi ^2}}\left(m_ek_p\left[2k_p+m_e\right]\sqrt{k_p^2+m_e^2}m_e^4\mathrm{ln}\left[{\displaystyle \frac{k_p+\sqrt{k_p^2+m_e^2}}{m_e}}\right]\right),`$
$`\mu |_0`$ $`=`$ $`g_\omega \omega _0+\sqrt{k_n^2+m_{}^2|_0},`$ (167)
$`\chi |_0`$ $`=`$ $`g_\omega \omega _0+\sqrt{k_p^2+m_{}^2|_0},`$ (169)
$`\mathrm{\Psi }|_0`$ $`=`$ $`\mathrm{\Lambda }|_0+{\displaystyle \frac{1}{3\pi ^2}}\left(\mu |_0k_n^3+\left[\chi |_0+\sqrt{k_p^2+m_e^2}\right]k_p^3\right).`$ (171)
The condition of chemical equilibrium is that $`\mu |_0=\chi |_0+\sqrt{k_p^2+m_e^2}`$.


format: text

text:
# The spinor bundle of Riemannian products

## Abstract.

In this note we compare the spinor bundle of a Riemannian manifold $`(M=M_1\times \mathrm{}\times M_N,g)`$ with the spinor bundles of the Riemannian factors $`(M_i,g_i)`$. We show that - without any holonomy conditions - the spinor bundle of $`(M,g)`$ for a special class of metrics is isomorphic to a bundle obtained by tensoring the spinor bundles of $`(M_i,g_i)`$ in an appropriate way. For $`N=2`$ and an one dimensional factor this construction was developed in \[Baum 1989a\]. Although the fact for general factors is frequently used in (at least physics) literature, a proof was missing.

I would like to thank Shahram Biglari, Mario Listing, Marc Nardmann and Hans-Bert Rademacher for helpful comments. Special thanks go to Helga Baum, who pointed out some difficulties arising in the pseudo-Riemannian case.

This text is part of the authors PhD. thesis \[Klinker 2003\]

We consider a Riemannian manifold $`(M=M_1\times \mathrm{}\times M_N,g)`$, which is a product of Riemannian spin manifolds $`(M_i,g_i)`$ and denote the projections on the respective factors by $`p_i`$. Furthermore the dimension of $`M_i`$ is $`D_i`$ such that the dimension of $`M`$ is given by $`D=_{i=1}^ND_i`$.

The tangent bundle of $`M`$ is decomposed as

(1)
$$T_{(x_0,\mathrm{},x_N)}M=p_1^{}T_{x_1}M_1\mathrm{}p_N^{}T_{x_N}M_N.$$

We omit the projections and write $`TM=_{i=1}^NTM_i`$.

The metric $`g`$ of $`M`$ need not be the product metric of the metrics $`g_i`$ on $`M_i`$, but is assumed to be of the form

(2)
$$\begin{array}{c}g_{ab}(x)|_{TM_i}=A_{i}^{}{}_{a}{}^{c}(x)g_{i}^{}{}_{cd}{}^{}(x_i)A_{i}^{}{}_{b}{}^{d}(x),\hfill \\ \hfill \text{for }D_1+\mathrm{}+D_{i1}+1a,bD_1+\mathrm{}+D_i,1iN\end{array}$$

In particular, for those metrics the splitting (1) is orthogonal, i.e. the frame bundle of $`M`$ can be reduced to a $`SO(D_1)\times \mathrm{}\times SO(D_N)`$-principal bundle, and this is isomorphic to the product of the frame bundles over $`M_i`$. The explicit form of the isomorphism is

(3)
$$\begin{array}{c}P_{SO}(M_1)\times P_{SO}(M_N)(E_1(x_1),\mathrm{},E_N(x_N))\hfill \\ \hfill (A_1(x)E_1(x_1),\mathrm{},A_N(x)E_N(x_N))P_{SO}(M)\end{array}$$

It is clear that such a manifold need not have a splitting of the holonomy group into subgroups of $`SO(D_i)`$, which would lead to an, at least local, Riemannian product (cf. \[Joyce 2000, sect. 3.2\]). Examples for such spaces are the Eguchi Hanson space – where we have $`M\times S^3`$ with metric $`\alpha (r)dr^2+h_r(S^3)`$ – or warped products of metrics.

We consider spinor bundles $`S_i`$ over $`M_i`$ and we are going to construct a bundle $`S`$ over $`M`$ from these spinor bundles. We will discuss what conditions are necessary for the bundle $`S`$ to be the spinor bundle over $`M`$ and how we have to modify the given Clifford multiplication on $`S_i`$ to extend it to $`S`$.

As is well known, the Clifford algebra of the sum of two vector spaces is the ($`_2`$-graded) tensor product of the Clifford algebras of the two summands (This is denoted by $`\mathrm{C}\mathrm{}(VW)=\mathrm{C}\mathrm{}(V)\widehat{}\mathrm{C}\mathrm{}(W)`$ in \[Lawson and Michelsohn 1989\]). Of course this is not restricted to two factors, but can be iterated. The case of two factors, i.e. $`N=2`$, will always be emphasized.

We consider the pullback bundles $`p_i^{}S_iM`$ of $`S_i`$ over $`M`$ and once more we omit the projections in our further notation. From $`S_i`$ we construct a bundle $`W`$ on $`M`$

$`W`$ $`=(S_1^+S_1^{})\mathrm{}(S_N^+S_N^{})`$
(4) $`={\displaystyle \underset{ϵ_1+\mathrm{}+ϵ_N=0}{}}{\displaystyle \underset{i=1}{\overset{N}{}}}S_i^{ϵ_i}{\displaystyle \underset{ϵ_1+\mathrm{}+ϵ_N=1}{}}{\displaystyle \underset{i=1}{\overset{N}{}}}S_i^{ϵ_i}.`$
For $`N=2`$ this is

$`W=`$ $`(S_1^+S_1^{})(S_2^+S_2^{})`$
$`=`$ $`\left(S_1^+S_2^+S_1^{}S_2^{}\right)\left(S_1^+S_2^{}S_1^{}S_2^+\right)`$
The bundles $`S_i^{ϵ_i}`$ are subbundles of $`S_i`$ with the $`_2`$-degree defined by the label111We identify $`_2=\{0=+,1=\}`$. $`ϵ_i=_2`$. This induces on $`W`$ the natural $`_2`$-grading given by the decomposition (4). An element $`\mathrm{\Xi }W`$ is called totally homogenous, if $`\mathrm{\Xi }^ϵ=\xi _1^{ϵ_1}\mathrm{}\xi _N^{ϵ_n}`$, i.e. it has only contributions from one of the summands in the decomposition (4). We will use the multi index notation $`ϵ=(ϵ_1,\mathrm{},ϵ_N)`$.

The subbundles of $`S_i`$ are chosen in such a way that the Clifford multiplication acts via $`TM_iS_i^\pm S_i^{}`$, i.e.

(5)
$$X_i\xi _i=X_i\left(\genfrac{}{}{0pt}{}{\xi _i^+}{\xi _i^{}}\right)=\left(\genfrac{}{}{0pt}{}{X_i\xi _i^{}}{X_i\xi _i^+}\right)=\left(\begin{array}{cc}0& X_i\\ X_i& 0\end{array}\right)\xi _i.$$

Clifford multiplication is an odd operation on $`S_i^+S_i^{}`$, i.e. $`|X_i\varphi _i|=|\varphi _i|+1`$.

We will consider the continuation of this to the bundle $`W`$. Therefore we introduce the linear map $`\delta _k:WW`$, which is defined by its action on totally homogenous elements $`\mathrm{\Xi }^ϵ=\xi _1^{ϵ_1}\mathrm{}\xi _N^{ϵ_N}`$

(6)
$$\delta _k(\mathrm{\Xi }^ϵ)=()^{ϵ_1+\mathrm{}+ϵ_{k1}}\mathrm{\Xi }^ϵ$$

$`\delta _k`$ can be seen as the sign, we get by “putting an odd operator acting on $`S_k`$ at the right place in the tensor product”. We define $`\delta _1=\mathrm{id}`$.

For $`X=X_1+\mathrm{}+X_NTM`$ and totally homogenous $`\mathrm{\Xi }^ϵ=\xi _1^{ϵ_1}\mathrm{}\xi _N^{ϵ_N}W`$ we define

(7)
$$\begin{array}{c}X\mathrm{\Xi }^ϵ:=_{i=1}^N(\mathrm{𝟏}\mathrm{}A_iX_i\mathrm{}\mathrm{𝟏})\mathrm{\Xi }^ϵ\hfill \\ \hfill =_{i=1}^N\delta _i(\mathrm{\Xi }^ϵ)\xi _1^{ϵ_1}\mathrm{}A_iX_i\xi _i^{ϵ_i}\mathrm{}\xi _N^{ϵ_N}.\end{array}$$

In particular for $`N=2`$

$$\begin{array}{c}(X_1+X_2)\mathrm{\Xi }^{(ϵ_1,ϵ_2)}=(\mathrm{𝟏}A_1X_1+A_2X_2\mathrm{𝟏})\mathrm{\Xi }^{(ϵ_1,ϵ_2)}\hfill \\ \hfill =A_1X_1\xi _1^{ϵ_1}\xi _2^{ϵ_2}+()^{ϵ_1}\xi _1^{ϵ_1}A_2X_2\xi _2^{ϵ_2}\end{array}$$

###### Proposition 1.

Let $`(M,g)`$, $`(M_i,g_i)`$, $`S_i`$ and $`W`$ be as before. The Clifford relation

(8)
$$(XY+YX)\mathrm{\Xi }=2g(X,Y)\mathrm{\Xi }$$

holds for all $`X,YTM`$ and $`\mathrm{\Xi }W`$.

###### Proof.

Because of linearity it is sufficient to prove the statement for $`X_iTM_i`$ and $`Y_jTM_j`$ and totally homogenous $`\mathrm{\Xi }^ϵ`$. We have to distinguish the cases $`i=j`$ and $`i<j`$. We recall the property $`\delta _k(X_i\mathrm{\Xi }^ϵ)=\{\begin{array}{cc}\delta _k(\mathrm{\Xi }^ϵ)\hfill & i<k\hfill \\ \delta _k(\mathrm{\Xi }^ϵ)\hfill & ik\hfill \end{array}`$.

With these informations we get

$`(X_i`$ $`Y_j+Y_jX_i)\mathrm{\Xi }^ϵ=`$
$`=`$ $`X_i\delta _j(\mathrm{\Xi }^ϵ)\xi _1^{ϵ_1}\mathrm{}A_jY_j\xi _j^{ϵ_j}\mathrm{}\xi _N^{ϵ_N}`$
$`+Y_j\delta _i(\mathrm{\Xi }^ϵ)\xi _1^{ϵ_1}\mathrm{}A_iX_i\xi _i^{ϵ_i}\mathrm{}\xi _N^{ϵ_N}`$
$`=`$ $`\{\begin{array}{cc}\delta _i(A_jY_j\mathrm{\Xi }^ϵ)\delta _j(\mathrm{\Xi }^ϵ)\xi _1^{ϵ_1}\mathrm{}A_iX_i\xi _i^{ϵ_i}\mathrm{}A_jY_j\xi _j^{ϵ_j}\mathrm{}\xi _N^{ϵ_N}\hfill & \\ +\delta _j(A_iX_i\mathrm{\Xi }^ϵ)\delta _i(\mathrm{\Xi }^ϵ)\xi _1^{ϵ_1}\mathrm{}A_iX_i\xi _i^{ϵ_i}\mathrm{}A_jY_j\xi _j^{ϵ_j}\mathrm{}\xi _N^{ϵ_N}\hfill & i<j\hfill \\ \delta _i(\mathrm{\Xi }^ϵ)\delta _i(\mathrm{\Xi }^ϵ)\xi _1^{ϵ_1}\mathrm{}(A_iX_i)(A_iY_i)\xi _i^{ϵ_i}\mathrm{}\xi _N^{ϵ_N}\hfill & \\ +\delta _i(\mathrm{\Xi }^ϵ)\delta _i(\mathrm{\Xi }^ϵ)\xi _1^{ϵ_1}\mathrm{}(A_iY_i)(A_iX_i)\xi _i^{ϵ_i}\mathrm{}\xi _N^{ϵ_N}\hfill & i=j\hfill \end{array}`$
$`=`$ $`\{\begin{array}{cc}0\hfill & i<j\hfill \\ \xi _1^{ϵ_1}\mathrm{}(A_iX_iA_iY_i+A_iY_iA_iX_i)\xi _i^{ϵ_i}\mathrm{}\xi _N^{ϵ_N}=g(X_i,Y_i)\mathrm{\Xi }^ϵ\hfill & i=j\hfill \end{array}`$
$`=`$ $`2g_i(A_iX_i,A_iY_j)\mathrm{\Xi }^ϵ=2g(X_i,Y_i)\mathrm{\Xi }^ϵ`$
Up to now we did not specify the subbundles of $`S_i^{ϵ_i}S_i`$. We have to distinguish two different situations:

* If $`dimM_i=2n_i`$ is even, the bundle $`S_i`$ (of rank $`2^{n_i}`$) itself admits a natural $`_2`$-grading induced by the volume element. And we take exactly this one.
* If $`dimM_i=2n_i+1`$ is odd, the spinor bundle $`S_i`$ (of rank $`2^{n_i}`$) does not admit a natural $`_2`$-grading. In this case we double the bundle and define $`S_i^+:=S_i`$ and $`S_i^{}:=\mathrm{\Pi }_iS_i`$.

In this definition of the subbundles, $`\mathrm{\Pi }_i:S_i^+S_i^{}S_i^{}S_i^+`$ denotes the parity operator with $`(\mathrm{\Pi }_i(S_i^+S_i^{}))^\pm =S_i^{}`$ and $`\mathrm{\Pi }_i^2=\mathrm{id}`$. Explicitly we have

(9)
$$\mathrm{\Pi }_i\left(\genfrac{}{}{0pt}{}{\xi _i^+}{\xi _i^{}}\right)=\left(\genfrac{}{}{0pt}{}{\xi _i^{}}{\xi _i^+}\right).$$

The parity operator is naturally extended to $`W`$ by its action on totally homogenous elements

(10)
$$\mathrm{\Pi }(\mathrm{\Xi })=\frac{1}{\sqrt{N}}\underset{i=1}{\overset{N}{}}()^{ϵ_1+\mathrm{}+ϵ_{i1}}\xi _1^{ϵ_1}\mathrm{}\mathrm{\Pi }_i\xi _i^{ϵ_i}\mathrm{}\xi _N^{ϵ_N}$$

where we use that the parity operator is formally odd. The proof for $`\mathrm{\Pi }^2=\mathrm{id}`$ is similar to that for the Clifford relation.

Let $`N_e`$ and $`N_o`$ be the number of even and odd dimensional manifolds in the product $`M=M_1\times \mathrm{}\times M_N`$, respectively (i.e. $`N_e+N_o=N`$). The dimension of the product manifold is $`dimM=2n+N_o`$ with $`n:=_{i=1}^Nn_i`$ and the rank of the bundle $`W`$ is $`2^{n+N_o}`$.

The spinor bundle of $`M`$ should have the rank $`2^{n+[\frac{N_o}{2}]}`$. We will construct a subbundle $`S`$ of $`W`$ of this rank.

The trick is to diagonalize some of the bundles which come from the odd dimensional manifolds. The diagonalization is denoted by $`\mathrm{\Delta }`$ and the bundle $`S`$ is constructed as follows:

Choose $`N_o[\frac{N_o}{2}]`$ of the odd dimensional bundles222This does not depend on the choice, because the resulting (non graded) bundles are isomorphic. and consider the subbundle $`SW`$ given by

(11)
$$S:=\underset{i=1}{\overset{N_e}{}}S_i\underset{j=N_e+1}{\overset{NN_o+[\frac{N_o}{2}]}{}}(S_j\mathrm{\Pi }_jS_j)\underset{k=NN_o+[\frac{N_o}{2}]+1}{\overset{N}{}}\mathrm{\Delta }(S_k\mathrm{\Pi }_kS_k)$$

which has the rank $`2^K`$ with

$`K`$ $`={\displaystyle \underset{i=1}{\overset{N_e}{}}}n_i+{\displaystyle \underset{j=N_e+1}{\overset{NN_o+[\frac{N_o}{2}]}{}}}n_j+(NN_o+[{\displaystyle \frac{N_o}{2}}]N_e)+{\displaystyle \underset{k=NN_o+[\frac{N_o}{2}]+1}{\overset{N}{}}}n_k`$
$`={\displaystyle \underset{i=1}{\overset{N}{}}}n_i+[{\displaystyle \frac{N_o}{2}}]=n+[{\displaystyle \frac{N_o}{2}}]`$
We specialize this to the case $`N=2`$:

* $`dimM_1=2n`$, $`dimM_2=2m`$: $`S`$ over $`M=M_1\times M_2`$ is given by

(12)
$$S=S_1S_2$$
* $`dimM_1=2n`$, $`dimM_2=2m+1`$: The subbundle $`S`$ over $`M`$ is given by

(13)
$$S=S_1\mathrm{\Delta }(S_2\mathrm{\Pi }_2S_2)S_1S_2.$$
* $`dimM_1=2n+1`$, $`dimM_2=2m+1`$: In this case $`S`$ is defined by

$$\begin{array}{c}S=(S_1\mathrm{\Pi }_1S_1)\mathrm{\Delta }(S_2\mathrm{\Pi }_2S_2)(S_1S_1)S_2\hfill \\ \hfill S_1S_2S_1S_2(S_2S_2)S_1,\end{array}$$

where we have to emphasize that the equivalences may not respect the $`_2`$-grading.

That this bundle is indeed a spinor bundle will be clear from the following remark. Although we have established the Clifford multiplication, the bundle is not a priori a spinor bundle. For that it should be constructed as an associated vector bundle to the $`Spin`$-principal bundle $`P_{Spin(D)}(M)`$ or of one reduction of this. If such a reduction does not exist in an appropriate way, then $`W`$ is not of this kind. The reason is that we do not have an action of $`Spin(D)`$ on the standard fibre of $`W`$.

The maximal subgroup of $`Spin(D)`$ which is compatible with the structure of the standard fibre is

(14)
$$\begin{array}{c}S(Pin(D_1)\times \mathrm{}\times Pin(D_N)):=Pin(D_1)\times \mathrm{}\times Pin(D_N)\mathrm{C}\mathrm{}^+(D)\hfill \\ \hfill _{i=1}^N\mathrm{C}\mathrm{}(D_i).\end{array}$$

This contains the subgroup $`Spin(D_1)\times \mathrm{}\times Spin(D_N)`$, which will be of interest soon.

In the next proposition we show that we are indeed able to write $`W`$ as an associated bundle, if we demand the reduction of the structure group of $`M`$ to $`SO(D_1)\times \mathrm{}\times SO(D_N)`$ compatible with the natural splitting (1). This is a weaker condition than (2), which was at least necessary to get the Clifford multiplication.

###### Proposition 2.

Let $`M`$ be a Riemannian spin manifold. If the structure group of $`M`$ can be reduced to $`SO(D_1)\times \mathrm{}\times SO(D_N)`$, we have a reduction of the spin principal bundle to $`Spin(D_1)\times \mathrm{}\times Spin(D_N)`$.

###### Proof.

We use the following notations: $`G:=Spin(D)`$, $`\stackrel{~}{G}:=Spin(D_1)\times \mathrm{}\times Spin(D_N)`$, $`H:=SO(D)`$ and $`\stackrel{~}{H}:=SO(D_1)\times \mathrm{}\times SO(D_N)`$.

Let $`P_H`$ be the principal bundle of orthonormal frames of $`(M,g)`$ and $`P_{\stackrel{~}{H}}`$ be the reduction to $`\stackrel{~}{H}`$, which we denote by $`ı`$. Furthermore let $`P_G`$ the $`Spin`$-principal bundle over $`M`$ with the two fold covering $`\lambda :P_GP_H`$, compatible with the right action of $`G`$ and $`H`$ and the cover $`\lambda :GH`$ for which we take the same symbol.

We collect this by writing $`P_G\stackrel{\lambda }{}P_H\stackrel{ı}{}P_{\stackrel{~}{H}}`$. We use the pullback construction for fibre bundles cf. (\[Kolář et al. 1993\]) to complete this edge to a commutative diagram and we denote the pullback by $`(P_G\times P_{\stackrel{~}{H}})\mathrm{}P_H`$.

(15)
Its total space is given by

(16)
$$\left\{(p_G,p_{\stackrel{~}{H}})P_G\times P_{\stackrel{~}{H}}\right|\lambda (p_G)=ı(p_{\stackrel{~}{H}})\}$$

The bundle $`(P_G\times P_{\stackrel{~}{H}})\mathrm{}P_H`$ is a principal bundle with the right action of its standard fibre $`G\times \stackrel{~}{H}\mathrm{}H=\{(g,\stackrel{~}{h})|\lambda (g)=ı(\stackrel{~}{h})=\stackrel{~}{h}\}\lambda ^1(\stackrel{~}{H})\stackrel{~}{G}`$ on the total space defined in the obvious way. We take the action on the Cartesian product $`(p_G,p_{\stackrel{~}{H}})a=(p_Ga,p_{\stackrel{~}{H}}\lambda (a))`$ which is compatible with the quotient:

(17)
$$\lambda (p_Ga)=\lambda (p_G)\lambda a=ı(p_{\stackrel{~}{H}})\lambda a=ı(p_{\stackrel{~}{H}}\lambda (a))$$

This shows that the $`Spin`$-principal bundle $`P_G`$ is reduced to a $`Spin(D_1)\times \mathrm{}\times Spin(D_N)`$-principal bundle. The reduction is just the left vertical arrow in the diagram. ∎

From this we get

###### Theorem 3.

Let $`(M=M_1\times \mathrm{}\times M_N,g)`$ be a Riemannian manifold, given as a (not necessarily Riemannian) product of the simply connected Riemannian spin manifolds $`(M_i,g_i)`$. The metric $`g`$ is connected to the metrics $`g_i`$ via (2).

Then $`M`$ is spin and the spinor bundle is isomorphic to the subbundle $`S`$ of $`W`$ constructed in this section.

###### Proof.

That $`M`$ is spin, follows immediately from the fact that the second Stiefel-Whitney class behaves additively under products of manifolds. The ON-frame bundles w.r.t. $`g`$ and $`g_1+\mathrm{}+g_N`$ are isomorphic cf. (3). On $`M`$ there exists exactly one spin structure, because $`M`$ is simply connected. So the spin principal bundles obtained by $`g`$ and $`g_1+\mathrm{}+g_N`$ are isomorphic. From the last proposition we get that the $`Spin`$-principal bundle $`P_G`$ over $`(M,g)`$ is reducible to $`Spin(D_1)\times \mathrm{}\times Spin(D_N)`$.

With the notations from the previous proposition and the construction for $`SW`$ we have established the following isomorphism

(18)
$$S\left((P_G\times P_{\stackrel{~}{H}})\mathrm{}P_H\right)\times _{\stackrel{~}{G}}\widehat{S}$$

where $`\widehat{S}`$ denotes the standard fibre of $`S`$. ∎

###### Remark 4.

An important example of this construction is the case $`N=2`$ with one of the factors being one dimensional and the metric is a warped product. This has been discussed in detail in \[Baum 1989a\] and \[Baum 1989b\].

From the construction it is clear that the spin connection obtained from the Levi-Civita connection of $`g`$ need not be compatible with the tensor structure of $`S`$. Explicit formulas for the connection in the case of an one dimensional factor and warped products can also be found in \[Baum 1989b\]. By claiming that the holonomy of $`M`$ is contained in $`SO(D_1)\times \mathrm{}\times SO(D_N)`$ – which is the same as to say that the projections $`p_i`$ are parallel – we make sure this further compatibility. As we mentioned above, this further assumption forces $`M`$ to be a local Riemannian product (in the case of $`M`$ being simply connected and complete this decomposition is global).

This constructions yields for $`A=id`$ the following

###### Proposition 5.

Let $`M`$ be the Riemannian product of $`(M_i,g_i)`$. Then the spinor bundle of $`M`$ with respect to the induced spin structure is given by $`S`$ from the construction above.

We add some remarks, which explain our further notation and draw the attention to the pseudo-Riemannian case.

###### Remark 6.

1. Proposition 5 is also true in the case of metrics which are not of Euclidean signature.
2. For the special case $`N=D`$, i.e. $`D_i=1`$, our constructions ends up with the Clifford action on $`^2\mathrm{}^2`$ as given in \[Baum et al. 1991\].
3. For another special case $`N=2`$ and $`D_2=1`$ this construction leads to the discussion in \[Baum 1989a\].
4. For $`N=2`$ in the cases $`_1`$ and $`_2`$, our construction yields the – at least in physics literature – frequently used decomposition of the $`\gamma `$-matrices $`\mathrm{\Gamma }^A`$, for $`1AD`$ in tensor products of $`\gamma `$-matrices of the respective factors $`\gamma ^a`$, $`\gamma ^\alpha `$ for $`1aD2n=D_2`$, $`1\alpha 2n=D_1`$. The $`_2`$-grading is ensured by using the volume element $`\widehat{\gamma }:=\gamma ^1\mathrm{}\gamma ^{2n}`$ which anticommutes with all $`\gamma ^a`$ and define

(19)
$$\mathrm{\Gamma }^a=\gamma ^a\widehat{\gamma },\mathrm{\Gamma }^\alpha =\mathrm{𝟏}\gamma ^\alpha ,$$

compare e.g. \[Duff et al. 1986\].


text:
# A “Little Big Bang” Scenario of Multifragmentation

## 1 Introduction

Recent theoretical studies of the morphology of simple fluids have suggested the presence of self-bound clusters of particles . The (time averaged) mass distributions of these clusters depend essentially on the energy and density of the system and resemble those found in nuclear multifragmentation ($`U`$-shapes, power law, exponentials). Furthermore, these distributions (event averaged) remain invariant if the system is allowed to expand freely.

The above observations suggest the following scenario for multifragmentation: In a first step, the nucleus is excited (and possibly compressed). In this ”hot” and dense phase, at least part of it reaches thermal equilibrium. Clusters, defined as self-bound ensembles of particles, are present in this medium. Immediately after, the system starts to expand, out of equilibrium, as *an ensemble of interacting clusters*. Asymptotically, these *clusters* cease to interact with each other, becoming the observed *fragments*. The sudden expansion reveals the cluster distribution of the primordial hot and dense phase of the system (*a Little Big Bang*) and “freezes” it. The long range Coulomb repulsion between clusters helps this freezing process. The smallness of the system also helps, because most of the clusters are close to the surface and escape freely into the vacuum. The clusters have ramified shapes in the dense medium and become, at asymptotic times, spherical fragments.

The ”standard” mechanism  of multifragmentation, adopted in the so-called ”Statistical Equilibrium Models” (SEM), is different: Once the system is excited, it undergoes an homogeneous and quasi-static expansion. At some sufficiently low density (freeze-out), it ”recondensates” as *an ensemble of spherical non-overlapping and non-interacting fragments*. At this stage, the models calculate the fragment mass distributions as those of a gas of non-interacting fragments with internal nuclear structure, at equilibrium. These models have been extremely successful in describing the observed fragment mass distributions . It is however important to notice that the main hypothesis (equilibrium at low freeze-out density, spherical fragments) required in these models to make the calculations feasible, are not substantiated by other, more microscopic, approaches of the multifragmentation phenomenon . In fact, despite the success in describing the mass distributions, this ”standard” scenario leads to some recurrent problems. Among them, the most important is the prediction of the kinetic energies of the fragments that are too low compared with experimental observations . This problem is ”corrected” by adding an *ad hoc* radial flow to the fragments. However the origin of this flow is unclear and seems inconsistent with the hypothesis of equilibrium and the low freeze-out density. Another problem concerns the internal temperature of the fragments, which in the experiments seems to be much lower than the temperature of the gas of fragments . We will see that the proposed high density scenario with pre-existing clusters helps to solve these problems in a very natural way.

This approach is substantiated by Classical Molecular Dynamics (CMD) calculations  of an ensemble of particles either confined in a container or freely expanding. Our goal is to interpret the (well known) results obtained with these methods at the light of our present knowledge on clusters in hot and dense fluids. We are thus not proposing a well-finished model of nuclear multifragmentation which could be directly compared with experimental data, but rather a generic scenario.

This paper is organized as follows. In section 2, we present results on the mapping between thermodynamics and clustering in dense and confined medium, in large and in nuclear-size systems. In section 3 we discuss the free expansion of these latter systems. Some final remarks will be found in Section 4.

## 2 Clusters in Confined Systems

In order to investigate the mapping between thermodynamics and clustering, we first have to define what a physical cluster is, over the whole range of densities and temperatures. A definition solely based on a proximity criterium (interaction radius) in configuration space is obviously not adequate at high densities and/or temperatures.

At low density, a cluster is naturally defined as a self bound ensemble of particles. By continuity we shall extend this definition to a dense medium. Several approaches can be used to identify self-bound clusters among the particles of the system. For instance, they can be defined as a partition of particles which minimizes the interaction energy between clusters . Clusters can also be built, bond by bond, using Hill’s criterium : Two particles are bound if their relative kinetic energy is less than their potential energy. In other words, a cluster is a set of particles close to each other in phase space and *not only* in configuration space. It must be emphasized that, with this model, knowing the positions and velocities of the particles, no adjusting parameter is necessary to recognize the clusters whatever the density or the temperature. In previous works , we checked that the cluster size distribution does not depend significantly on the specific definition of self-bound clusters mentioned above.

The mapping between thermodynamics and clustering can be studied by performing classical microcanonical molecular dynamic simulations of a simple fluid. Particles confined in a cubic container interact through a Lennard-Jones potential ($`V(r)=4ϵ((\sigma /r)^{12}(\sigma /r)^6)`$), where $`ϵ`$ and $`\sigma `$ define the energy and length scales respectively. The natural time-unit is defined by $`\tau _0=\sqrt{m\sigma ^2/(48ϵ)}`$, where $`m`$ is the mass of the particles. When Coulomb interaction is included, it is done with the prescription of reference . Self-bound clusters are recognized by means of Hill’s criterium. The details of the calculations can be found in references .

In order to establish an accurate mapping between thermal and geometrical properties of this fluid, we first present results for the case of a large system ($`N12000`$ particles) with periodic boundary conditions  (and, of course without Coulomb interaction):

* When crossing the condensation curve at sub-critical densities, a macroscopic cluster appears in the liquid-gas coexistence region444This correspondence between thermodynamics and clustering is not observed when clusters are simply defined by a constant interaction radius in the configurational space.. Although defined in a complete different theoretical framework, our clusters mark the liquid-gas coexistence line as Fisher’s clusters do 555Notice however that Fisher’s clusters are non-interacting and not microscopically defined..
* A percolation line starts at the thermodynamical critical point (within the uncertainties inherent to the critical slowing down) and goes through the supercritical region of the phase diagram. Thus, a power law behavior of the cluster size distribution does not necessarily imply that the system is in its thermodynamical critical state. Also notice that this percolation line does not extend below the critical temperature.
* It can be shown that the Hill prescription we are using to define self-bound clusters in a realistic fluid is equivalent  to the prescription introduced by Coniglio and Klein  to define physical clusters in the lattice-gas model. This gives a theoretical framework to understand the existence of this critical line. Along this line (usually called Kertész line), the critical exponents are those of the random percolation, except at the thermodynamical critical point, where the exponents are those of the correlated percolation .
* The energy of the system is nearly constant along the critical percolation line.

These results remain valid in small systems . We have studied the behaviour of small systems like atomic nuclei by performing similar calculations with $`N=189`$ particles confined in a cubic container with perfectly reflecting walls and interacting *via* a Lennard-Jones potential plus Coulomb. The “phase diagram” is represented in figure 1. Although the sharp effects found in the thermodynamical limit are smoothed by the finite size and surface effects, we reach basically the same conclusions.

* As is shown in figure 2, the shape of the cluster size distribution changes suddenly when crossing the condensation curve for a very small variation of temperature.
* The cluster size distribution presents a power law behavior along a percolation “line” (more precisely a band due to the finite size effects) which starts at the thermodynamical critical point (region) (see figure 1). Along this line, the slope of the cluster size distribution ($`\tau 2.5`$) is in accordance with percolation theory666Due to the range of our potential, in a small system, all particles interact with each other, so the mean field result is expected..

We note also that, in the supercritical region, a given cluster size distribution corresponds qualitatively to a given energy of the system. We stress that the converse is not true, notably below the coexistence line.

Studying the degrees of freedom associated to the internal motion of the particles inside the clusters, and the degrees of freedom associated to the center of mass motion of these clusters, we define two ”effective temperatures”, $`T^{}(s)`$ and $`T^{cm}(s)`$, as two third of the average kinetic energy of the corresponding ensemble 777However, if we were able to isolate a cluster from the others (this is what happens during the expansion), the internal effective temperature would become the real thermodynamic temperature of the isolated fragment.. Let us define the internal velocity $`v^{}`$ as the velocity of a particle calculated in the center of mass system of the cluster it belongs to. As is shown in Fig 3, the internal velocity distribution for cluster of size $`s=2`$ contrasts strongly with a Maxwell distribution. As it should be, by definition of the self-bound clusters, the distribution is truncated at $`v^{}=1`$. However, increasing the cluster size, the velocity distribution tends to a Maxwell distribution characterized by the asymptotic temperature $`T^{}(s)`$ for $`s>50`$ (see figure 4). Similar results  have been recently obtained by calculating the first order corrections (namely the interaction between a monomer and a self-bound cluster) to a perfect gas of clusters at low density. Whatever the size of the clusters and the thermodynamic state of the system, the internal effective temperature $`T^{}(s)`$ is always less than the (real) temperature of the fluid. On the other hand, the effective temperature associated to the ”gas of clusters”, $`T^{cm}`$ is greater than the (real) temperature of the fluid. In figure 4, we illustrate this inequality which was already mentioned . In brief, the system can be seen as a ”hot” gas of ”cold” physical clusters.

The present results show that various basic hypothesis of the SEM  are not fulfilled. These models consider an ensemble of spherical, non-interacting clusters, with an internal temperature equal to the one of the global system, confined in a volume at low densities ($`\rho 0.2`$). On the contrary, we find that:

1. The clusters exhibit fractal shapes in the confined system888As it should be, in large systems ($`N=12000`$), the fractal dimension ($`D_f=2.55`$) of the self-bound clusters along the percolation line has been checked to be in agreement with the one of random percolation . (see figure 5, left). We must point out that it is not incompatible with the cold spherical fragments we expect at the end of the expansion. Indeed, by following the evolution of a self-bound cluster during the expansion (cf section 3) we find that it becomes spherical as it should be (see figure 5, right).
2. Even at densities as low as $`\rho =0.10.2`$, the strong interaction (without Coulomb) between clusters still represents about 40 % of the total potential energy.
3. The internal temperature of the clusters is lower than the temperature of the system. Thus, this fragments will become in the course of the expansion (see next section) the observed fragments, without significant particle evaporation.

Usually three observables are used to extract the temperature from nuclear experimental data: i) the kinetic energies of the fragments, ii) the population of the excited states of fragments  and iii) the ratio of isotopic yields . All these methods assume that these observables will yield the thermodynamic temperature of the system, $`T`$. Our study contradicts this. As shown in figure 5, the internal and “kinetic energy” temperatures of the fragments are not a measure of $`T`$.

## 3 Fragments in Expanding Systems

Having studied the properties of clusters in confined systems, we now turn to the problem of their evolution once the system is allowed to expand freely. The questions we want to address are :

* How and when do the fragments appear ? Does the asymptotic fragment size distribution and kinetic energies relate to the initial configuration or to an intermediate freeze-out configuration ?
* Does the system follow a quasi-static path up to some, low density, freeze-out configuration where the final fragment distribution would be determined, as assumed in the standard  models ?
* How can the temperature and density of this initial configuration be inferred from the information measured at asymptotic times ?

In order to answer these questions, we have analyzed a number of expansions of the system using the CMD with Lennard-Jones plus Coulomb potentials . The system is characterized, as in section 2, by its energy (microcanonical ensemble), its density and the number of particles. It is enclosed in a cubic container with perfectly reflective walls. The calculations proceed in two steps :

* After a thermalisation period, the container is removed at a time defined as $`t=0`$.
* The system is then allowed to expand freely during a period of time sufficient to establish the final (asymptotic) fragment size distribution.

Let us first examine the time evolution of the components of the energy. Figure 6 displays, the evolution of the L-J potential energy, the Coulomb energy, the total kinetic energy and the total energy.

At $`t=0`$, the effect of the opening of the container is clearly visible by a rapid decrease of the potential energy. This corresponds to the fact that the mean distance between clusters increases rapidly. It can be noted that the asymptotic potential energy does not go to zero, indicating the existence of fragments of finite size. The kinetic energy shows a complementary behaviour : after a sharp decrease, the Coulomb acceleration is clearly observed.

We shall now study the evolution of the size of the largest fragment (Figure 7). The line gives the result of a single event whereas the full dots represent an average over 100 expansions. The striking feature of this figure is the observation that the average value of the largest *fragment* at asymptotic times is very close to the average value of the largest *cluster* identified in the container ($`t<0`$). The correspondence between the largest cluster in the container (time average) and the asymptotic largest fragment (event average) is systematically observed. This conclusion is valid for all the calculations that we have performed, independently of the initial density and energy.

Figure 8 shows confirmation of this conclusion by the examination of more detailed quantities. It presents the fragment size distributions as a function of time. To within statistical fluctuations, one observes that they are almost identical. Figure 9 shows the evolution in time of the multiplicity of fragments of size $`s=5`$. There again, the asymptotic value is very close to the initial one. In view of the rapid decrease of the fragment size distribution observed (figure 8), the stability of this quantity is impressive.

The correlation between the initial and final size distributions is, we believe, a consequence of the violent expansion phase (which follows the opening of the container), the small size of the system and the presence of the Coulomb repulsion between the clusters. Because of the small size of the system, the fragments are always close to the surface and will escape freely into the vacuum. For a very large system, the expansion phase would give rise, in the central region, to substantial thermal and chemical activity and hence modify the thermodynamic parameters of the system. In a small system, this modification appears almost unobservable. The influence of the expansion phase is enhanced by the presence of the Coulomb force. The Coulomb potential not only accelerates the expansion but also introduces, between the fragments, a Coulomb barrier, which inhibits particle exchanges.

The violence of the expansion phase is illustrated by the figure 10 which shows the *persistence* of fragments and the density of the system, as a function of time. The persistence is defined as the fraction of the number of particles present in a given fragment which will remain in this fragment at asymptotic times ($`t=400\tau _0`$). The density, a poorly defined quantity in such an inhomogeneous expansions, is deduced from the r.m.s. radius of an uniformly charged sphere of the same Coulomb energy. Figure 10 shows that this persistence evolves very rapidly and almost reaches its asymptotic value at densities of the order of $`0.8`$: An exchange of particles between fragments is clearly observed during this (short) period, but this limited chemical activity does not induce a modification of the average size of the fragments.

The picture that comes out from these calculations therefore corresponds to a fast, non equilibrium, expansion of a gas of clusters. This expansion acts as a ”developer” of the presence of these clusters, the Coulomb repulsion helping to ”fix” them. Thus, the (event) average distribution of fragments reflects the (time) average distribution of clusters in the dense and hot initial state.

We now study the kinetic energies of the fragments. Figure 12 shows the asymptotic kinetic energy of expansions *without* Coulomb interaction ($`E_{tot}=0.2`$). As can be seen, the kinetic energies obtained from different starting densities (0.3 and 1.3) are quantitatively different. This demonstrates once more that the expansion is a non-equilibrium phenomenon. If equilibrium, at a given energy, was maintained during the expansion, the system starting from $`\rho =1.3`$ would, at a point in time, reach the same point in $`\{T\rho \}`$ space as the one from which the system at $`\rho =0.3`$ starts and hence the two systems would be indistinguishable in all respect: Figure 12 demonstrates that it is not the case.

Figure 12 shows the kinetic energies in the case *with* Coulomb interaction. The shape of these distributions is characteristic of the presence of Coulomb forces and is close to what is observed in nuclear reactions at Fermi bombarding energies ($`50`$ MeV/A) for symmetric nuclei . The solid lines show a qualitative (analytical) prediction of what would be expected in the case of a uniformly charged medium.

The important feature that is revealed by this figure is the dependence of the kinetic energies on the initial density. This result appears to confirm the evolution suggested by the analysis of the fragment size distribution. Very rapidly after the removal of the container, the fragments cease to strongly interact with the medium and are accelerated by the (Coulomb) repulsive potential. Because this expansion takes place out of equilibrium, the Coulomb energy does not transform into thermal energy but in fragment kinetic energy. The higher the density of the initial configuration, the greater the final kinetic energies are. We can therefore consider that the asymptotic kinetic energies are a “measure” of the density at which the systems departs from equilibrium.

This scenario accounts, at least partially, for the radial flow often associated with multifragmentation processes . As stated above (section 1), the standard method of analysis, using the SEM, assumes a low density configuration and therefore naturally underpredicts the kinetic energies of the fragments when the system has evolved (out of equilibrium) from a high density state. To compensate for this, an extra component of “radial flow” is artificially introduced. In the scenario we suggest, such an extra flow is not necessary when the calculation is performed at the correct initial density. It is important to remark that our scenario does not necessarily imply large values for the kinetic energy. For example, for peripheral heavy ion reactions, no initial compression is expected. Even more, for proton-nucleus or $`\pi `$-nucleus collisions, initial densities lower than the normal equilibrium nuclear density are plausible because of the formation of holes inside a nucleus with normal volume .

Having “measured” the density, we now turn to the estimation of the temperature. The sensitivity of the kinetic energies to the density, provides a “method” to determine the thermodynamic parameters of the point at which the system departs from equilibrium. The analysis of the fragment size distribution will select which iso-$`n(s)`$ line is to be considered. The analysis of the kinetic energies will determine the density. The intersection of both will yield $`T`$ and $`\rho `$. Evidently, this method is physics (model) dependent, i.e. it is only if one is able to calculate the relationship $`E(T,\rho )`$ and the fragment size distribution at each $`(T,\rho )`$ point that the thermodynamical parameters can be extracted from the analysis of the data measured at asymptotic times.

## 4 Final remarks

We have explored a new scenario of multifragmentation, based on the observation that a dense and hot fluid, at equilibrium, can be viewed as a hot gas of cold clusters. These clusters are defined as self-bound ensembles of particles. When the system is allowed to expand freely, it proceeds, *out of equilibrium*, as an ensemble of interacting clusters. Once these clusters cease to interact with each other, they become the observable fragments.

We believe that this ”Little Big Bang” scenario is not only more realistic than the ”standard” low density freeze-out model but that it also allows to solve a number of recurrent problems:

* Equilibrium is assumed only at the beginning of the process, when the high density and temperature make this assumption more likely. No ”quasi-static” hypothesis of expansion up to a low freeze-out density is needed.
* Due to this compacity, the Coulomb repulsion between fragments generates larger kinetic energies. Clusters in the dense medium have ”ramified” (fractal in large systems) forms, allowing more compact configurations in the initial stage. The observed radial flows can thus be explained, at least partially, in a natural way.
* Isolated fragments are cold, because the precursor clusters are already cooler than the ensemble of the system.

The aim of this work was to present a new scenario of multifragmentation, not a well finished model producing results directly comparable with experimental data. One has to keep in mind that the present scenario is supported by classical molecular dynamics simulations of the expansion of a system of particles interacting through Lennard-Jones plus Coulomb potentials, initially confined in a container. Although we believe that these results are generic for describing the free expansion of any simple fluid (and that nuclear matter at high temperatures behaves as a simple fluid), it is necessary to confront our results with other calculations. Namely, with Quantum Molecular Dynamics  calculations to test the equilibrium hypothesis in the first stages of the reaction, and with Fermionic Molecular Dynamics  to check the importance of quantum effects, particularly at sub-critical temperatures.

We would like to thank F. Lavaud for his contribution to the expansion phase calculations. One of us (N.S.) wishes to acknowledge financial support of the European TMR Network-Fractals (Contract number: FMRXCT-980183).