text:
# Gravitational instability of slowly rotating isothermal spheres

## 1 Introduction

Recently, the statistical mechanics of self-gravitating systems has attracted considerable attention (Chavanis, Sommeria & Robert 1996, Chavanis & Sommeria 1998, de Vega, Sanchez & Combes 1998, Youngkins & Miller 2000, Follana & Laliena 2000, Semelin, Sanchez & de Vega 2001, Cerruti-Sola, Cipriani & Pettini 2001, Ispolatov & Cohen 2001, Chavanis 2002a, Taruya & Sakagami 2002, de Vega & Sanchez 2002, Huber & Pfenniger 2002…). This topic was introduced in the $`1960`$s by Antonov (1962) and Lynden-Bell & Wood (1968) and further developed by Hertel & Thirring (1971), Horwitz & Katz (1978), Katz (1978) and Padmanabhan (1989) among others (see a complete list of references in the review of Padmanabhan, 1990). These authors pointed out the particularity of self-gravitating systems to possess negative specific heats. They showed that this strange property is responsible for the inequivalence of statistical ensembles (microcanonical/canonical) and the occurence of giant phase transitions associated with gravitational collapse. For a long time, these topics were only discussed in the astrophysical literature and were considered as a curiosity (not to say a fallacy) from statistical mechanicians. The situation is changing lately as these properties are re-discovered for other physical systems with long-range interactions which can be studied in the laboratory (see, e.g., Gross 2001). For that reason, maybe, the statistical mechanics of self-gravitating systems comes back to fashion with new perspectives.

From our point of view, the statistical mechanics of self-gravitating systems is far from being completely understood and rests on simplifying idealizations. The first idealization is to enclose the system within a box so as to prevent evaporation. It is only under this condition (or by introducing more realistic truncated models) that a rigorous statistical mechanics of self-gravitating systems can be carried out. Second, for most astrophysical systems, the relaxation time by two-body encounters is much larger than the age of the universe so that a more subtle, collisionless, relaxation must be advocated to explain the structure of galaxies. This is the concept of violent relaxation formalized by Lynden-Bell in 1967. Then, it is implicitly assumed that the relaxation towards statistical equilibrium proceeds to completion, which is not necessarily the case in reality. Indeed, it is possible that the relaxation stops before the maximum entropy state is attained (see Lynden-Bell 1967, Tremaine, Hénon & Lynden-Bell 1987, Chavanis et al. 1996). This problem of incomplete relaxation must be approached with extensive numerical simulations. Of course, this program was started long ago (e.g., van Albada 1982) but only recently are $`N`$-body simulations carefully compared with the predictions of the statistical mechanics approach (Cerriti-Sola et al. 2001, Huber & Pfeninger 2002) with variable success.

The statistical equilibrium of a non-rotating classical gas enclosed within a box was first investigated by Antonov (1962). He worked in the microcanonical ensemble and found that thermodynamical equilibrium exists only above a critical energy $`E_{min}=0.335GM^2/R`$. Below that energy, the system is expected to collapse and overheat; this is the so-called “gravothermal catastrophe” (Lynden-Bell & Wood 1968). An isothermal collapse also occurs below a critical temperature $`T_{min}=\frac{GMm}{2.52R}`$ in the canonical ensemble. This thermodynamical instability is closely related to the dynamical Jeans instability (Semelin et al. 2001, Chavanis 2002a). In a recent series of papers, we considered some extensions of the Antonov problem to the case of self-gravitating fermions (Chavanis & Sommeria 1998, Chavanis 2002c), special and general relativity (Chavanis 2002b) and confined polytropes (Chavanis 2002d). These studies complete previous works on the subject. We also introduced a simple dynamical model of Brownian particles in gravitational interaction (Chavanis, Rosier & Sire 2002). By introducing by hands a friction and a noise, we force the system to increase entropy continuously, thereby avoiding the problem of incomplete relaxation. This model can be used to test precisely the ideas of equilibrium statistical mechanics (inequivalence of ensembles, phase transitions, gravitational instabilities…) and is sufficiently simple to allow for a thorough analytical investigation of the collapse regime when an equilibrium state does not exist.

In all these studies, the system is assumed to be non-rotating so that the conservation of angular momentum is trivially satisfied. The object of the present paper is to extend the statistical mechanical approach to the case of rotating self-gravitating systems. This problem has been considered previously by Lagoute & Longaretti (1996), Laliena (1999), Lynden-Bell (2000) and Fliegans & Gross (2002) with different types of models. Clearly, the most interesting situation is the case of rapidly rotating systems since a wide variety of structures can emerge as maximum entropy states (Votyakov et al. 2002). However, we shall restrict ourselves in the present paper to the case of slowly rotating systems for which the problem can be tackled analytically. We shall adapt to the case of isothermal spheres the classical procedure developed by Chandrasekhar (1933) for distorted polytropes, i.e. we shall expand the solutions of the Boltzmann-Poisson equation in terms of Legendre polynomials. A similar procedure was performed by Lagoute & Longaretti (1996) for rotating globular clusters subject to tidal forces and described by an extended Michie-King model. We believe that it is useful to consider the case of a gas enclosed within a box so as to make a clear connexion with the Antonov model when the rotation is set to zero. In particular, we shall derive the expression of the thermodynamical parameters for slowly rotating isothermal spheres and show how the classical spiral of Lynden-Bell & Wood (1968) in the $`ET`$ plane is modified by rotation. We shall show that rotation avances the onset of gravothermal catastrophe in the microcanonical ensemble and delays the isothermal collapse in the canonical ensemble. Using the turning point criterion of Katz (1978), we argue that the series of equilibrium becomes unstable at the point of minimum energy (in the microcanonical ensemble) or at the point of minimum temperature (in the canonical ensemble). Interestingly, these instabilities happen to be close, in each ensemble, to the point of maximum flattening. We establish the generalization of the singular isothermal solution to the case of a slowly rotating configuration. We also consider the case of slowly rotating self-gravitating fermions. This system exhibits phase transitions between “gaseous” states with an almost uniform distribution of matter and “condensed” states with a core-halo structure. By cooling below a critical temperature, an almost nonrotating gaseous medium can collapse into a rotating “fermion ball” containing a large fraction of mass and angular momentum.

## 2 Statistical mechanics of rotating self-gravitating systems

### 2.1 The mean-field approach

Consider a system of $`N`$ particles, each of mass $`m`$, interacting via Newtonian gravity. We allow the system to have a non vanishing angular momentum. Let $`f(𝐫,𝐯,t)`$ denote the distribution function of the system, i.e. $`f(𝐫,𝐯,t)d^3𝐫d^3𝐯`$ gives the mass of particles whose position and velocity are in the cell $`(𝐫,𝐯;𝐫+d^3𝐫,𝐯+d^3𝐯)`$ at time $`t`$. The integral of $`f`$ over the velocity determines the spatial density

$$\rho =fd^3𝐯.$$
(1)
The mass and angular momentum of the configuration are given by

$$M=\rho d^3𝐫,$$
(2)
$$𝐋=f𝐫\times 𝐯d^3𝐫d^3𝐯.$$
(3)
On the other hand, in the mean-field approximation, the energy can be expressed as

$$E=\frac{1}{2}fv^2d^3𝐫d^3𝐯+\frac{1}{2}\rho \mathrm{\Phi }d^3𝐫=K+W,$$
(4)
where $`K`$ is the kinetic energy and $`W`$ the potential energy. The gravitational potential $`\mathrm{\Phi }`$ is related to the star density by the Newton-Poisson equation

$$\mathrm{\Delta }\mathrm{\Phi }=4\pi G\rho .$$
(5)
The equilibrium configuration of the system is determined by maximizing the Boltzmann entropy

$$S=k\frac{f}{m}\mathrm{ln}\frac{f}{m}d^3𝐫d^3𝐯,$$
(6)
while conserving mass, angular momentum and energy. We stress that for systems interacting via a long-range potential, like gravity, the mean-field approximation is exact so that this procedure is entirely rigorous and provides a simple approach to the problem.

Introducing Lagrange multipliers $`\alpha `$, $`\beta `$ and $`\beta \mathrm{\Omega }\mathrm{\Omega }\mathrm{\Omega }`$ for each constraint, we find that the critical points of entropy are given by

$$f=A^{}e^{\beta (\frac{v^2}{2}+\mathrm{\Phi }𝛀(𝐫\times 𝐯))},$$
(7)
where

$$\beta =\frac{m}{kT}.$$
(8)
We can rewrite the foregoing expression for $`f`$ in the more suggestive form

$$f=A^{}e^{\frac{1}{2}\beta (𝐯𝛀\times 𝐫)^2}e^{\beta (\mathrm{\Phi }\frac{1}{2}(𝛀\times 𝐫)^2)}.$$
(9)
We note that the most probable distribution of a rotating self-gravitating system is a Maxwell-Boltzmann law with a dispersion the same at every point (isothermal distribution). Moreover, the most probable form of rotation is a rigid rotation $`𝛀`$. By integrating Eq. (9) over the velocity, we find that the density is given by the Boltzmann factor

$$\rho =Ae^{\beta (\mathrm{\Phi }\frac{1}{2}(𝛀\times 𝐫)^2)}.$$
(10)
The quantity in parenthesis is the “effective” potential in the rotating frame

$$\mathrm{\Phi }_{eff}=\mathrm{\Phi }\frac{1}{2}(𝛀\times 𝐫)^2,$$
(11)
accounting for inertial forces. The equilibrium configuration is obtained by solving the Boltzmann-Poisson system

$$\mathrm{\Delta }\mathrm{\Phi }=4\pi GAe^{\beta (\mathrm{\Phi }\frac{1}{2}(𝛀\times 𝐫)^2)}.$$
(12)
and relating the Lagrange multipliers $`A`$, $`\beta `$ and $`𝛀`$ to the constraints $`M`$, $`E`$ and $`𝐋`$. The configuration is thermodynamically stable if the second order variations of entropy are negative for any perturbation that satisfies the constraints on mass, energy and angular momentum. This corresponds to the microcanonical description. Alternatively, we could describe the system in the canonical ensemble assuming that $`\beta `$ and $`𝛀`$ are given, instead of $`E`$ and $`𝐋`$. In that case, the equilibrium configuration is obtained by maximizing the free energy $`J=S\beta E\beta 𝛀𝐋`$ at fixed $`M`$, $`\beta `$ and $`𝛀`$. Of course, the microcanonical and canonical ensembles yield the same critical points, i.e. the critical points of entropy at fixed mass, energy and angular momentum and the critical points of free energy at fixed mass, temperature and angular velocity coincide. Only the onset of instability, regarding the second order variations of $`S`$ or $`J`$, will differ from an ensemble to the other. In the non-rotating case, the thermodynamical stability analysis was performed by Padmanabhan (1989) in the microcanonical ensemble and by Chavanis (2002a) in the canonical ensemble, by solving an eigenvalue equation. The onset of instability can also be determined by the turning point criterion of Katz (1978) who extended the theory of Poincaré on linear series of equilibrium.

### 2.2 The rotating isothermal sphere

To determine the structure of rotating isothermal spheres, we first introduce the function $`\mathrm{\Psi }=\beta (\mathrm{\Phi }_{eff}\mathrm{\Phi }_0)`$, where $`\mathrm{\Phi }_0`$ is the gravitational potential at $`r=0`$. Then, the density field can be written

$$\rho =\rho _0e^\mathrm{\Psi },$$
(13)
where $`\rho _0`$ is the central density. Introducing the notations $`\xi =(4\pi G\beta \rho _0)^{1/2}r`$, $`v=\mathrm{\Omega }^2/2\pi G\rho _0`$ and using spherical coordinates $`r,\theta ,\varphi `$, we can rewrite the Boltzmann-Poisson equation (12) in the form

$$\frac{1}{\xi ^2}\frac{}{\xi }\left(\xi ^2\frac{\mathrm{\Psi }}{\xi }\right)+\frac{1}{\xi ^2}\frac{}{\mu }\left((1\mu ^2)\frac{\mathrm{\Psi }}{\mu }\right)=e^\mathrm{\Psi }v,$$
(14)
where $`\mu =\mathrm{cos}\theta `$ and we have neglected $`\varphi `$ by assuming that the configuration is symmetric with respect to the axis $`𝛀=\mathrm{\Omega }𝐳`$. The true gravitational potential $`\mathrm{\Phi }`$ is related to $`\mathrm{\Psi }`$ by the relation

$$\beta \mathrm{\Phi }=\mathrm{\Psi }+\frac{1}{4}v\xi ^2(1\mu ^2)+\beta \mathrm{\Phi }_0.$$
(15)
Eq. (14) is the fundamental equation of the problem. For $`v=0`$ (no rotation), we recover the Emden equation

$$\frac{1}{\xi ^2}\frac{d}{d\xi }\left(\xi ^2\frac{d\psi }{d\xi }\right)=e^\psi .$$
(16)
with $`\psi =\psi ^{}=0`$ at $`\xi =0`$ (Chandrasekhar 1942).

So far, we have made no approximation regarding the value of the angular velocity. We shall now consider the case of slowly rotating structures and let $`v0`$. Assuming the following form for our solution

$$\mathrm{\Psi }=\psi +v\mathrm{\Theta }+v^2\mathrm{\Theta }_2+\mathrm{}$$
(17)
and keeping terms only to first order in $`v`$, we find that $`\mathrm{\Theta }`$ satisfies the equation

$`{\displaystyle \frac{1}{\xi ^2}}{\displaystyle \frac{}{\xi }}\left(\xi ^2{\displaystyle \frac{\mathrm{\Theta }}{\xi }}\right)+{\displaystyle \frac{1}{\xi ^2}}{\displaystyle \frac{}{\mu }}\left((1\mu ^2){\displaystyle \frac{\mathrm{\Theta }}{\mu }}\right)=\mathrm{\Theta }e^\mathrm{\Psi }1.`$
(18)
Now, following a procedure that dates back to Milne (1923) and Chandrasekhar (1933), we shall assume for $`\mathrm{\Theta }`$ the following form

$$\mathrm{\Theta }=\varphi _0(\xi )+\underset{j=1}{\overset{+\mathrm{}}{}}A_j\varphi _j(\xi )P_j(\mu ),$$
(19)
where $`P_j(\mu )`$ are the Legendre polynomials satisfying the differential equation

$$\frac{}{\mu }\left((1\mu ^2)\frac{P_j}{\mu }\right)+j(j+1)P_j=0.$$
(20)
Substituting for $`\mathrm{\Theta }`$ from Eq. (19) in Eq. (18) and equating coefficients of $`P_j`$, we get

$$\frac{1}{\xi ^2}\frac{d}{d\xi }\left(\xi ^2\frac{d\varphi _0}{d\xi }\right)=e^\psi \varphi _01,$$
(21)
$`{\displaystyle \frac{1}{\xi ^2}}{\displaystyle \frac{d}{d\xi }}\left(\xi ^2{\displaystyle \frac{d\varphi _j}{d\xi }}\right)=\left({\displaystyle \frac{j(j+1)}{\xi ^2}}e^\psi \right)\varphi _j,`$
$`(j=1,2,\mathrm{})`$ (22)
with $`\varphi _j=\varphi _j^{}=0`$ at $`\xi =0`$.

So far, the $`A_j`$ are arbitrary. They will be determined by requiring that the gravitational potential and its radial derivative are continuous across the sphere at $`r=R`$. Now, outside the sphere the potential is given by the Laplace equation

$$\mathrm{\Delta }\mathrm{\Phi }_{ext}=0.$$
(23)
The physically acceptable solution of this equation tending to zero at infinity can be written as

$$\beta \mathrm{\Phi }_{ext}=\frac{B_0}{\xi }+v\underset{j=1}{\overset{+\mathrm{}}{}}\frac{B_j}{\xi ^{j+1}}P_j(\mu )$$
(24)
On, the other hand, according to Eq. (15), we have inside the spherical box

$$\beta \mathrm{\Phi }_{int}=\mathrm{\Psi }+\frac{1}{6}v\xi ^2(1P_2(\mu ))+\beta \mathrm{\Phi }_0,$$
(25)
where $`\mathrm{\Psi }(\xi ,\mu )`$ is given by Eqs. (17)(19) and use has been made of the identity $`1\mu ^2=\frac{2}{3}(1P_2(\mu ))`$. Let us denote by $`\alpha =(4\pi G\beta \rho _0)^{1/2}R`$ the value of $`\xi `$ at $`r=R`$. Comparing the inner and the external potentials at $`\xi =\alpha `$, and also their derivative, we find that $`A_j=B_j=0`$ if $`j2`$ and, for $`j=2`$,

$$\frac{B_2}{\alpha ^3}=A_2\varphi _2(\alpha )\frac{1}{6}\alpha ^2,$$
(26)
$$\frac{3B_2}{\alpha ^4}=A_2\varphi _2^{}(\alpha )\frac{1}{3}\alpha .$$
(27)
Solving for $`A_2`$, we get

$$A_2=\frac{5}{6}\frac{\alpha ^2}{3\varphi _2(\alpha )+\alpha \varphi _2^{}(\alpha )}.$$
(28)
Therefore, the solution of the Boltzmann-Poisson equation (14) to first order in $`v`$ is given by

$`\mathrm{\Psi }=\psi (\xi )+v\left\{\varphi _0(\xi )+{\displaystyle \frac{5}{6}}{\displaystyle \frac{\alpha ^2}{3\varphi _2(\alpha )+\alpha \varphi _2^{}(\alpha )}}\varphi _2(\xi )P_2(\mu )\right\},`$
(29)
with

$$\frac{1}{\xi ^2}\frac{d}{d\xi }\left(\xi ^2\frac{d\varphi _0}{d\xi }\right)=e^\psi \varphi _01,$$
(30)
$`{\displaystyle \frac{1}{\xi ^2}}{\displaystyle \frac{d}{d\xi }}\left(\xi ^2{\displaystyle \frac{d\varphi _2}{d\xi }}\right)=\left({\displaystyle \frac{6}{\xi ^2}}e^\psi \right)\varphi _2,`$ (31)
with $`\varphi _j=\varphi _j^{}=0`$ at $`\xi =0`$.

### 2.3 The slowly rotating singular isothermal sphere

For $`\xi +\mathrm{}`$, the solution of the Emden equation (16) behaves like (Chandrasekhar 1942)

$$e^\psi \frac{2}{\xi ^2}.$$
(32)
Substituting this asymptotic behaviour in Eq. (30) and making the change of variables $`\xi =e^t`$, we obtain

$$\frac{d^2\varphi _0}{dt^2}+\frac{d\varphi _0}{dt}+2\varphi _0=e^{2t}.$$
(33)
The solution of this equation is readily found. After returning to original variables, we get for large $`\xi `$’s:

$$\varphi _0=\frac{1}{8}\xi ^2+\frac{A}{\xi ^{1/2}}\mathrm{cos}\left(\frac{\sqrt{7}}{2}\mathrm{ln}\xi +\delta \right).$$
(34)
Following the same procedure for the function $`\varphi _2`$, we have

$$\frac{d^2\varphi _2}{dt^2}+\frac{d\varphi _2}{dt}4\varphi _2=0,$$
(35)
leading to

$$\varphi _2=\frac{A^{}}{\xi ^{1/2}}\mathrm{cosh}\left(\frac{\sqrt{17}}{2}\mathrm{ln}\xi +\delta ^{}\right).$$
(36)
Keeping only the leading term in Eqs. (34)(36) in the limit $`\xi +\mathrm{}`$, we find that

$$\varphi _0\frac{1}{8}\xi ^2;\varphi _2K\xi ^\gamma ,\gamma =\frac{\sqrt{17}}{2}\frac{1}{2}.$$
(37)
Substituting these results in Eq. (29) and returning to dimensional variables, we finally obtain

$`\rho ={\displaystyle \frac{1}{2\pi G\beta r^2}}+{\displaystyle \frac{\mathrm{\Omega }^2}{8\pi G}}\left\{1{\displaystyle \frac{20}{3(3+\gamma )}}\left({\displaystyle \frac{R}{r}}\right)^{2\gamma }P_2(\mu )\right\},`$ (38)
which is the expression of the slowly rotating singular isothermal sphere.

## 3 The thermodynamical parameters

### 3.1 The mass

In spherical coordinates, the mass is given by

$$M=2\pi _1^{+1}_0^R\rho r^2𝑑r𝑑\mu .$$
(39)
Introducing the dimensionless variables previously defined, we obtain, using Eq. (13),

$$\eta \frac{\beta GM}{R}=\frac{1}{2\alpha }_1^{+1}_0^\alpha e^\mathrm{\Psi }\xi ^2𝑑\xi 𝑑\mu .$$
(40)
Substituting for $`\mathrm{\Psi }`$ from Eq. (29) in Eq. (40) and recalling that $`_1^{+1}P_j(\mu )𝑑\mu =0`$, we obtain to first order in $`v`$:

$$\eta =\frac{1}{\alpha }_0^\alpha e^\psi (1v\varphi _0(\xi ))\xi ^2𝑑\xi .$$
(41)
Now,

$$_0^\alpha e^\psi \xi ^2𝑑\xi =_0^\alpha 𝑑\xi \frac{d}{d\xi }\left(\xi ^2\frac{d\psi }{d\xi }\right)=\alpha ^2\psi ^{}(\alpha ),$$
(42)
$`{\displaystyle _0^\alpha }e^\psi \varphi _0\xi ^2𝑑\xi ={\displaystyle _0^\alpha }𝑑\xi \left\{\xi ^2+{\displaystyle \frac{d}{d\xi }}\left(\xi ^2{\displaystyle \frac{d\varphi _0}{d\xi }}\right)\right\}`$
$`={\displaystyle \frac{\alpha ^3}{3}}\alpha ^2\varphi _0^{}(\alpha ).`$ (43)
Therefore, the normalized inverse temperature $`\eta `$ is expressed in terms of $`\alpha `$ and $`v`$ by the relation

$$\eta =\alpha \psi ^{}(\alpha )+v\left[\frac{\alpha ^2}{3}+\alpha \varphi _0^{}(\alpha )\right].$$
(44)
Expressing the central density in terms of $`\alpha `$, we find that

$$v=\frac{2\eta }{\alpha ^2}\omega ^2,\mathrm{with}\omega =\mathrm{\Omega }\left(\frac{R^3}{GM}\right)^{1/2}.$$
(45)
Therefore, the foregoing expression for $`\eta `$ can be rewritten to order $`\omega ^2`$:

$$\eta =\alpha \psi ^{}(\alpha )\left\{1+2\omega ^2\left[\frac{1}{3}+\frac{1}{\alpha }\varphi _0^{}(\alpha )\right]\right\}.$$
(46)
The validity of our perturbative approach requires that

$$ϵ2\omega ^2\left[\frac{1}{3}+\frac{1}{\alpha }\varphi _0^{}(\alpha )\right]1.$$
(47)
### 3.2 The angular momentum

Using the Maxwell-Boltzmann distribution (9), we can rewrite the angular momentum (3) in the form

$$𝐋=\rho 𝐫\times 𝐮d^3𝐫,\mathrm{with}𝐮=𝛀\times 𝐫.$$
(48)
After straightforward manipulations, we find that the angular momentum is related to the angular velocity by

$$L_i=I_{ij}\mathrm{\Omega }_j,$$
(49)
where $`I_{ij}=\rho (r^2\delta _{ij}r_ir_j)d^3𝐫`$ is the moment of inertia tensor. If the density depends only on $`r`$ and $`\theta `$, we have $`𝐋=I𝛀`$ with

$$I=2\pi _1^{+1}_0^R\rho r^4(1\mu ^2)𝑑r𝑑\mu .$$
(50)
For our purpose, it is sufficient to determine the angular momentum to first order in $`\mathrm{\Omega }`$. Therefore, we need just the expression of the moment of inertia $`I`$ for a non-rotating isothermal sphere. Introducing the dimensionless variables previously defined, and using the expression (46) for $`\eta `$, we get

$$\frac{I}{MR^2}=\frac{2}{3\alpha ^4\psi ^{}(\alpha )}_0^\alpha e^\psi \xi ^4𝑑\xi .$$
(51)
Then,

$`\lambda =(\alpha )\omega ,\mathrm{with}\lambda ={\displaystyle \frac{L}{\sqrt{GM^3R}}}.`$ (52)
### 3.3 The energy

Quite generally, the potential energy of a self-gravitating system can be written in the form (see, e.g., Binney & Tremaine 1987)

$$W=\rho 𝐫\mathrm{\Phi }d^3𝐫.$$
(53)
It can be readily verified that the condition $`f=f(\frac{w^2}{2}+\mathrm{\Phi }_{eff})`$, where $`𝐰=𝐯𝐮`$ and $`\mathrm{\Phi }_{eff}`$ is the effective potential (11), is equivalent to the condition of hydrostatic equilibrium in the rotating frame

$$p=\rho \mathrm{\Phi }\rho 𝛀\times (𝛀\times 𝐫),$$
(54)
where

$$p=\frac{1}{3}fw^2d^3𝐰,$$
(55)
is the local pressure. For an isothermal gas $`p=\frac{k}{m}\rho T`$. Inserting Eq. (54) in Eq. (53), we obtain after straightforward calculations

$$W+3pd^3𝐫+𝐋𝛀p𝐫𝑑𝐒=0,$$
(56)
where $`d𝐒`$ is a surface element normal to the spherical box. On the other hand, we can write the kinetic energy in the form

$$K=K_{rot}+K_{th},$$
(57)
where

$$K_{rot}=\frac{1}{2}\rho u^2d^3𝐫=\frac{1}{2}𝐋𝛀,$$
(58)
is the rotational energy and

$$K_{th}=\frac{1}{2}fw^2d^3𝐫d^3𝐰=\frac{3}{2}pd^3𝐫,$$
(59)
is the thermal energy. For an isothermal gas $`K_{th}=\frac{3}{2}NkT`$. Therefore, Eq. (56) becomes

$$W+2K=p𝐫𝑑𝐒,$$
(60)
which is the Virial theorem for a self-gravitating gas enclosed within a box. For the total energy $`E=K+W`$, we have

$$E=K+p𝐫𝑑𝐒.$$
(61)
This expression is valid even if the system is not axisymmetric. If, now, the density depends only on $`r`$ and $`\theta `$, and if the gas is isothermal, it is possible to combine the formulae derived in Sect. 2.2 to obtain the expression

$`\mathrm{\Lambda }{\displaystyle \frac{ER}{GM^2}}={\displaystyle \frac{3}{2\eta }}+{\displaystyle \frac{1}{2}}\lambda \omega {\displaystyle \frac{\alpha ^2}{\eta ^2}}e^{\psi (\alpha )}(1v\varphi _0(\alpha )),`$ (62)
where $`\eta `$ and $`v`$ can be expressed in terms of $`\alpha `$ and $`\omega `$ by Eqs. (46)(45). We also recall that our theory is valid to order $`\omega ^2`$. For $`\omega =0`$, Eqs. (46)(62) reduce to the equations of state obtained by Lynden-Bell & Wood (1968) for non-rotating isothermal spheres.

### 3.4 The entropy

Using Eqs. (9)(10) and (13), the distribution function of a rotating isothermal gas can be written

$`f=\left({\displaystyle \frac{\beta }{2\pi }}\right)^{3/2}e^{\beta \frac{w^2}{2}}\rho _0e^\mathrm{\Psi }.`$ (63)
Substituting this expression in Eq. (6) and expressing $`\rho _0`$ in terms of $`\alpha `$, we get

$`{\displaystyle \frac{mS}{k}}={\displaystyle \frac{1}{2}}M\mathrm{ln}\beta 2M\mathrm{ln}\alpha +{\displaystyle \rho \mathrm{\Psi }d^3𝐫}.`$ (64)
Throughout this paper, we shall not write the constant terms (depending on the fixed parameters $`M`$ and $`R`$) which enter in the expression of the entropy. Therefore, a term $`\beta K_{th}=\frac{3}{2}M`$ has been ignored in Eq. (64). Using the definition of $`\mathrm{\Psi }`$, the last integral can be rewritten

$`{\displaystyle \rho \mathrm{\Psi }d^3𝐫}=2\beta W\beta K_{rot}\beta M\mathrm{\Phi }_0.`$ (65)
Using $`W=EK_{rot}K_{th}`$, we obtain

$`{\displaystyle \frac{mS}{k}}={\displaystyle \frac{1}{2}}M\mathrm{ln}\beta {\displaystyle \frac{3}{2}}\beta 𝐋𝛀`$
$`2M\mathrm{ln}\alpha +2\beta E\beta M\mathrm{\Phi }_0.`$ (66)
Introducing the dimensionless parameters previously defined, we get

$`{\displaystyle \frac{S}{Nk}}={\displaystyle \frac{1}{2}}\mathrm{ln}\eta {\displaystyle \frac{3}{2}}\eta \lambda \omega 2\mathrm{ln}\alpha 2\eta \mathrm{\Lambda }\beta \mathrm{\Phi }_0.`$ (67)
We now need to determine the central potential $`\mathrm{\Phi }_0`$. The condition that $`\mathrm{\Phi }`$ and $`\mathrm{\Phi }/\xi `$ are continuous at $`\xi =\alpha `$ implies for $`j=0`$:

$`{\displaystyle \frac{B_0}{\alpha }}=\psi (\alpha )+v\varphi _0(\alpha )+{\displaystyle \frac{1}{6}}v\alpha ^2+\beta \mathrm{\Phi }_0,`$ (68)
$`{\displaystyle \frac{B_0}{\alpha ^2}}=\psi ^{}(\alpha )+v\varphi _0^{}(\alpha )+{\displaystyle \frac{1}{3}}v\alpha .`$ (69)
Comparing Eq. (69) with Eq. (44), we see that $`B_0/\alpha =\eta `$. Inserting this result in Eq. (68), we obtain

$`\beta \mathrm{\Phi }_0=\eta +\psi (\alpha )+v\left[{\displaystyle \frac{\alpha ^2}{6}}+\varphi _0(\alpha )\right].`$ (70)
With this new relation, the entropy (67) becomes

$`{\displaystyle \frac{S}{Nk}}={\displaystyle \frac{1}{2}}\mathrm{ln}\eta 2\mathrm{ln}\alpha +\psi (\alpha )+\eta 2\mathrm{\Lambda }\eta `$
$`{\displaystyle \frac{3}{2}}\eta \lambda \omega +v\left[{\displaystyle \frac{\alpha ^2}{6}}+\varphi _0(\alpha )\right].`$ (71)
### 3.5 The flattening function

The value of the potential at $`\xi `$ for a non-rotating configuration is $`\psi (\xi )`$. For a rotating configuration, the equation of the surface with the same value of the potential is given, to first order in $`v`$, by

$`\xi ^{}=\xi {\displaystyle \frac{v}{\psi ^{}(\xi )}}\left\{\varphi _0(\xi )+{\displaystyle \frac{5}{6}}{\displaystyle \frac{\alpha ^2}{3\varphi _2(\alpha )+\alpha \varphi _2^{}(\alpha )}}\varphi _2(\xi )P_2(\mu )\right\}.`$
(72)
If $`a`$ denotes the largest radius of the isodensity surface (at the equator $`\mu =0`$) and $`b`$ the smallest radius (at the pole $`\mu =1`$), and if we define the flattening by $`f=1b/a`$ (see, e.g., Lagoute & Longaretti 1996), we get

$`f(\xi )={\displaystyle \frac{5}{4}}{\displaystyle \frac{v}{\xi \psi ^{}(\xi )}}{\displaystyle \frac{\alpha ^2}{3\varphi _2(\alpha )+\alpha \varphi _2^{}(\alpha )}}\varphi _2(\xi ).`$ (73)
Considering the limit $`\xi +\mathrm{}`$ and returning to dimensional variables, we find that the flattening function of the rotating singular isothermal sphere is

$`f(r)={\displaystyle \frac{5}{4}}{\displaystyle \frac{\beta \mathrm{\Omega }^2R^2}{3+\gamma }}\left({\displaystyle \frac{r}{R}}\right)^\gamma ,`$ (74)
a result which can also be derived directly from Eq. (38).

Coming back to Eq. (73), we find that the flattening behaves with the distance (for a given value of $`\alpha `$) as

$`F(\xi )={\displaystyle \frac{\varphi _2(\xi )}{\xi \psi ^{}(\xi )}}.`$ (75)
On the other hand, the flattening at the edge of the configuration ($`\xi =\alpha `$) is given by

$`f(\alpha )={\displaystyle \frac{5}{2}}\omega ^2{\displaystyle \frac{\varphi _2(\alpha )}{3\varphi _2(\alpha )+\alpha \varphi _2^{}(\alpha )}},`$ (76)
where we have used Eqs. (45)(46) to eliminate the variable $`v`$ in profit of $`\omega `$. In the canonical ensemble (fixed $`\omega `$), this function depends on $`\alpha `$ as

$`F_{CE}(\alpha )={\displaystyle \frac{\varphi _2(\alpha )}{3\varphi _2(\alpha )+\alpha \varphi _2^{}(\alpha )}},`$ (77)
In the microcanonical ensemble (fixed $`\lambda `$), we must express the angular velocity in terms of the angular momentum using the relation (52). In that case, the flattening at the edge of the configuration depends on $`\alpha `$ as

$`F_{MCE}(\alpha )={\displaystyle \frac{1}{(\alpha )^2}}F_{CE}(\alpha ).`$ (78)
We shall come back to these results in the following section.

## 4 Equilibrium phase diagram

### 4.1 Microcanonical ensemble

The microcanonical ensemble corresponds to isolated systems characterized by their energy $`\mathrm{\Lambda }`$ and their angular momentum $`\lambda `$. In order to determine the equilibrium phase diagram $`\eta (\mathrm{\Lambda })`$ for different values of $`\lambda `$, we need to solve Eqs. (16)(30)(31) numerically. Expanding the functions $`\psi `$, $`\varphi _0`$ and $`\varphi _2`$ in Taylor series for $`\xi 0`$, we get

$$\psi =\frac{1}{6}\xi ^2\frac{1}{120}\xi ^4+\mathrm{},$$
(79)
$$\varphi _0=\frac{1}{6}\xi ^2+\frac{1}{120}\xi ^4+\mathrm{},$$
(80)
$$\varphi _2=\xi ^2+\frac{1}{14}\xi ^4+\mathrm{}$$
(81)
so that $`\psi ^{\prime \prime }(0)=1/3`$, $`\varphi _0^{\prime \prime }(0)=1/3`$ and $`\varphi _2^{\prime \prime }(0)=2`$. The integration can be continued numerically by a standard Runge-Kutta routine. In Figs. 1,2,3, we plot the curves $`(\alpha )`$, $`\mathrm{\Lambda }(\alpha )`$ and $`\eta (\alpha )`$ defined by Eqs. (51)(62)(46). These curves exhibit damped oscillations and tend to the values

$$_s=\frac{2}{9},\eta _s=2+\frac{27}{4}\lambda ^2,\mathrm{\Lambda }_s=\frac{1}{4}\frac{63}{32}\lambda ^2.$$
(82)
as $`\alpha +\mathrm{}`$. This asymptotic limit corresponds to the singular solution (38). The iso-density contours of the rotating singular isothermal sphere are represented in Fig. 4 (for $`\lambda =0.17`$).

In Fig. 5, we have represented the curve $`\eta (\mathrm{\Lambda })`$. It has a classical spiral behaviour as noted by a number of authors in the non-rotating case. There is no equilibrium state (i.e., no critical point of entropy) above the value $`\mathrm{\Lambda }_c(\lambda )`$. In that case, the system will collapse and overheat (gravothermal catastrophe). It is also at this point that the critical points of entropy become unstable (saddle points) in the series of equilibria (Katz 1978). We see that rotation tends to favour the instability, i.e., the gravothermal catastrophe occurs sooner than in the non-rotating case.

In Fig. 6, we have represented the angular velocity of the system as a function of energy for different values of the angular momentum. We observe that the curve has a spiral behavior similar to the $`\eta \mathrm{\Lambda }`$ diagram, but reversed. We might expect that the moment of inertia decreases as the system becomes more and more concentrated, resulting in an increase of angular velocity along the series of equilibria. This is true for moderate density contrasts (up to $`709`$), coinciding with the region of stability, but not for larger density contrasts. Indeed, although the central density tends to diverge, the mass contained in the core is low and does not dominate the moment of inertia. Therefore, the moment of inertia and the angular velocity have a non monotonous (in fact oscillatory) behaviour with $`\alpha `$ and saturate to finite values $`_s`$ and $`\omega _s=\lambda /_s`$ as $`\alpha +\mathrm{}`$.

In Fig. 7, we plot the flattening function $`F(\xi )`$ defined by Eq. (75). As expected, the flattening is is a monotonous function of the distance. For $`\xi 0`$, $`F3`$ and for $`\xi +\mathrm{}`$, $`F\xi ^\gamma `$. In Fig. 8, we plot the flattening at the edge of the configuration $`F_{MCE}(\alpha )`$ as a function of $`\alpha `$ in the microcanonical ensemble (see Eq. (78)). We observe that the curve displays damped oscillations towards the asymptotic value $`\frac{1}{3+\gamma }(9/2)^2`$. In particular the flattening (by unit of $`\lambda ^2`$) is maximum for $`\alpha =32.6\mathrm{}`$. Interestingly, this value lies precisely in the range of values at which the gravothermal catastrophe sets in (compare with Fig. 2).

### 4.2 Canonical ensemble

The canonical ensemble is characterized by the specification of the inverse temperature $`\eta `$ and the angular velocity $`\omega `$. In Fig. 9, we have represented the curve $`\mathrm{\Lambda }(\eta )`$ for different values of $`\omega `$. There is no equilibrium state (i.e., no critical point of free energy) above the value $`\eta _c(\omega )`$. In that case, the system will undergo an isothermal collapse. It is also at this point that the solutions become unstable (saddle points of free energy) in the series of equilibria (Katz 1978). We see that rotation tends to delay the instability, i.e., the isothermal collapse occurs later than in the non-rotating case.

In Fig. 10, we plot the flattening at the edge of the configuration $`F_{CE}(\alpha )`$ as a function of $`\alpha `$ in the canonical ensemble. The curve displays damped oscillations towards the value $`\frac{1}{3+\gamma }`$. The flattening (by unit of $`\omega ^2`$) is maximum for $`\alpha =5.4\mathrm{}`$. Interestingly, this value is close to the typical values at which the isothermal collapse sets in (compare with Fig. 11 giving the inverse temperature along the series of equilibria).

## 5 Rotating self-gravitating fermions

The previous results can be easily generalized to the case of self-gravitating fermions. This extension is relatively straightforward and we shall just give the main steps of the calculations. The thermodynamical parameters of a non-rotating Fermi gas at finite temperature have been calculated by Chavanis & Sommeria (1998) and we shall adopt a similar presentation. For a rotating configuration, the Fermi-Dirac distribution can be written

$$f=\frac{\eta _0}{1+\lambda e^{\beta \mathrm{\Phi }_{eff}}e^{\beta \frac{w^2}{2}}},$$
(83)
where $`\eta _0`$ is the maximum allowable value of the distribution function and $`\lambda >0`$ a strictly positive parameter (inverse of the fugacity) insuring that $`f\eta _0`$ (the other quantities have been defined previously). For quantum particles with spin $`s`$, $`\eta _0=(2s+1)m^4/(2\pi \mathrm{})^3`$. The distribution function (83) can be obtained by maximizing the Fermi-Dirac entropy

$$S=\left\{\frac{f}{\eta _0}\mathrm{ln}\frac{f}{\eta _0}+\left(1\frac{f}{\eta _0}\right)\mathrm{ln}\left(1\frac{f}{\eta _0}\right)\right\}d^3𝐫d^3𝐯,$$

at fixed mass, energy and angular momentum. This form of entropy also occurs in the context of violent relaxation for collisionless self-gravitating systems (Lynden-Bell 1967, Chavanis, Robert & Sommeria 1996, Chavanis & Sommeria 1998). Defining $`\mathrm{\Psi }=\beta (\mathrm{\Phi }_{eff}\mathrm{\Phi }_0)`$ and $`k=\lambda e^{\beta \mathrm{\Phi }_0}`$, the spatial density can be written

$$\rho =\rho _0\frac{I_{1/2}(ke^\mathrm{\Psi })}{I_{1/2}(k)},\mathrm{with}\rho _0=\frac{4\pi \sqrt{2}\eta _0}{\beta ^{3/2}}I_{1/2}(k),$$
(84)
where $`I_{1/2}`$ is the Fermi integral

$`I_n(t)={\displaystyle _0^+\mathrm{}}{\displaystyle \frac{x^n}{1+te^x}}𝑑x`$ (85)
of order $`n=1/2`$. Substituting the relation (84) in the Poisson equation (5) and introducing the dimensionless parameters $`\xi =(4\pi G\beta \rho _0)^{1/2}r`$ and $`v=\frac{\mathrm{\Omega }^2}{2\pi G\rho _0}`$, we obtain

$`{\displaystyle \frac{1}{\xi ^2}}{\displaystyle \frac{}{\xi }}\left(\xi ^2{\displaystyle \frac{\mathrm{\Psi }}{\xi }}\right)+{\displaystyle \frac{1}{\xi ^2}}{\displaystyle \frac{}{\mu }}\left((1\mu ^2){\displaystyle \frac{\mathrm{\Psi }}{\mu }}\right)={\displaystyle \frac{I_{1/2}(ke^\mathrm{\Psi })}{I_{1/2}(k)}}v.`$
(86)
We can also check that Eq. (15) keeps the same form. For $`v=0`$, Eq. (86) reduces to

$`{\displaystyle \frac{1}{\xi ^2}}{\displaystyle \frac{d}{d\xi }}\left(\xi ^2{\displaystyle \frac{d\psi }{d\xi }}\right)={\displaystyle \frac{I_{1/2}(ke^\psi )}{I_{1/2}(k)}},`$ (87)
with $`\psi =\psi ^{}=0`$ at $`\xi =0`$. For $`k+\mathrm{}`$ (classical limit), we can use the limiting form of the Fermi integral

$`I_n(t){\displaystyle \frac{1}{t}}\mathrm{\Gamma }(n+1),(t+\mathrm{})`$ (88)
and we recover the classical Emden equation (16). For $`k0`$ (completely degenerate limit), we have

$`I_n(t){\displaystyle \frac{(\mathrm{ln}t)^{n+1}}{n+1}},(t0)`$ (89)
and the differential equation (87) becomes equivalent to the Lane-Emden equation for a polytrope of index $`3/2`$.

We shall now consider the case of slowly rotating structures and let $`v0`$. Repeating the steps of Sec. 2.2, we find that Eq. (29) remains valid with the new functions $`\varphi _0`$ and $`\varphi _2`$ defined by

$`{\displaystyle \frac{1}{\xi ^2}}{\displaystyle \frac{d}{d\xi }}\left(\xi ^2{\displaystyle \frac{d\varphi _0}{d\xi }}\right)={\displaystyle \frac{1}{2}}{\displaystyle \frac{I_{1/2}(ke^\psi )}{I_{1/2}(k)}}\varphi _01,`$ (90)
$`{\displaystyle \frac{1}{\xi ^2}}{\displaystyle \frac{d}{d\xi }}\left(\xi ^2{\displaystyle \frac{d\varphi _2}{d\xi }}\right)=\left({\displaystyle \frac{6}{\xi ^2}}{\displaystyle \frac{1}{2}}{\displaystyle \frac{I_{1/2}(ke^\psi )}{I_{1/2}(k)}}\right)\varphi _2,`$ (91)
with $`\varphi _j=\varphi _j^{}=0`$ at $`\xi =0`$. In arriving at (90)-(91), we have used the identity

$`I_n^{}(t)={\displaystyle \frac{n}{t}}I_{n1}(t),(n>0),`$ (92)
for $`n=1/2`$, which can be easily established from Eq. (85).

We now determine the thermodynamical parameters of a slowly rotating self-gravitating Fermi gas. We can check that the inverse temperature is still given by Eq. (46) and that the relation (45) remains valid. On the other hand, eliminating the central density between Eq. (84) and $`\alpha =(4\pi G\beta \rho _0)^{1/2}R`$, we find that

$`\eta ={\displaystyle \frac{\mu ^2}{\alpha ^4}}I_{1/2}(k)^2,\mathrm{where}\mu =\eta _0\sqrt{512\pi ^4G^3MR^3},`$ (93)
is the “degeneracy parameter” (Chavanis & Sommeria 1998). Using Eq. (46), we obtain

$`\mu ^2I_{1/2}(k)^2=\alpha ^5\psi ^{}(\alpha )\left\{1+2\omega ^2\left[{\displaystyle \frac{1}{3}}+{\displaystyle \frac{1}{\alpha }}\varphi _0^{}(\alpha )\right]\right\}.`$ (94)
The relation between the angular momentum and the angular velocity is still given by Eq. (52) with the moment of inertia

$`={\displaystyle \frac{2}{3\alpha ^4\psi ^{}(\alpha )}}{\displaystyle _0^\alpha }{\displaystyle \frac{I_{1/2}(ke^\psi )}{I_{1/2}(k)}}\xi ^4𝑑\xi .`$ (95)
Using Eq. (61), the energy can be put in the form

$`\mathrm{\Lambda }={\displaystyle \frac{\alpha ^7}{\mu ^4}}{\displaystyle _0^\alpha }\left\{{\displaystyle \frac{I_{3/2}(ke^\psi )}{I_{1/2}(k)^5}}{\displaystyle \frac{3}{2}}v{\displaystyle \frac{I_{1/2}(ke^\psi )}{I_{1/2}(k)^5}}\varphi _0(\xi )\right\}\xi ^2𝑑\xi `$
$`+{\displaystyle \frac{1}{2}}\lambda \omega {\displaystyle \frac{2}{3}}{\displaystyle \frac{\alpha ^{10}}{\mu ^4}}\left\{{\displaystyle \frac{I_{3/2}(ke^{\psi (\alpha )})}{I_{1/2}(k)^5}}{\displaystyle \frac{3}{2}}v{\displaystyle \frac{I_{1/2}(ke^{\psi (\alpha )})}{I_{1/2}(k)^5}}\varphi _0(\alpha )\right\}.`$
(96)
Finally, for the entropy we obtain (see Appendix A)

$`{\displaystyle \frac{\eta _0S}{M}}=\mathrm{ln}k+\eta +\psi (\alpha ){\displaystyle \frac{7}{3}}\eta \mathrm{\Lambda }{\displaystyle \frac{4}{3}}\eta \lambda \omega +v\left[{\displaystyle \frac{\alpha ^6}{6}}+\varphi _0(\alpha )\right]`$
$`{\displaystyle \frac{2}{9}}{\displaystyle \frac{\alpha ^6}{\mu ^2}}\left\{{\displaystyle \frac{I_{3/2}(ke^{\psi (\alpha )})}{I_{1/2}(k)^3}}{\displaystyle \frac{3}{2}}v{\displaystyle \frac{I_{1/2}(ke^{\psi (\alpha )})}{I_{1/2}(k)^3}}\varphi _0(\alpha )\right\}.`$ (97)
For $`v=0`$, Eqs. (46) (94)(96) and (97) reduce to the equations of state obtained by Chavanis & Sommeria (1998) for the non-rotating Fermi gas. For $`k+\mathrm{}`$ (non degenerate limit), we recover the equations derived in Sec. 3.

The equilibrium phase diagram can be obtained in the following manner. For given $`k`$, $`\mu `$ and $`\omega `$, we can solve Eqs. (87) (90) and (91) until the value $`\xi =\alpha `$ for which the relation (94) is satisfied. Then, Eqs. (46) and (96) determine the temperature and the energy of the configuration. If the angular momentum is fixed instead of the angular velocity, we must use Eq. (52) with Eq. (95) to express $`\omega `$ in terms of $`\lambda `$. By varying the parameter $`k`$ (for a fixed value of the degeneracy parameter $`\mu `$), we can cover the whole diagram in parameter space. A complete description of this diagram has been given by Chavanis (2002c) in the non-rotating case.

In Fig. 12 we represent the equilibrium phase diagram of self-gravitating fermions for a degeneracy parameter $`\mu =10^5`$ and for different values of angular momentum. We observe that degeneracy has the effect of unwinding the spiral of Fig. 5. For sufficiently large values of the degeneracy parameter, there is still gravitational collapse at $`\mathrm{\Lambda }_c`$ accompanied by a rise of temperature, but this “gravothermal catastrophe” stops when the core of the system becomes degenerate. This leads to the formation of a “fermion ball” which contains a moderately large fraction of mass $`\alpha M`$ (at point $`D`$, we typically have $`\alpha 0.2`$). In the microcanonical ensemble, the decrease of potential energy in the core is compensated by an increase of temperature. Therefore, the mass $`(1\alpha )M`$ contained in the halo undergoes an expansion which, in our model, is arrested by the walls of the box. As a result, the density of the halo is almost uniform. Typical density profiles are given by Chavanis & Sommeria (1998) in the non rotating case. The expansion of the halo explains why the moment of inertia of the system increases despite the formation of a massive nucleus. Therefore, the angular velocity decreases during the collapse contrary to what might be expected. The angular velocity is represented as a function of energy in Fig. 13. It has a complicated behaviour which corresponds to the unwindement of the spiral of Fig. 6.

For smaller values of the degeneracy parameter, the gravitational phase transition is suppressed (Chavanis & Sommeria 1998, Chavanis 2002c) and the equilibrium phase diagram has the structure of Fig. 14. This diagram is similar to the one found by Fliegans & Gross (2002) in their two-dimensional model of rotating self-gravitating systems. Had they used a smaller value of their cut-off radius (which plays the role of the inverse of our degeneracy parameter), they would have probably obtained a diagram similar to the one of Fig. 12. However, when the cut-off radius (or degeneracy) is sufficiently large, the spiral unwinds and the $`\eta `$-$`\mathrm{\Lambda }`$ curve is univalued like in Fig. 14. For low energies, the equilibrium states have a core-halo structure with a partially degenerate nucleus and a dilute envelope. As energy decreases any further, the nucleus contains more and more mass and becomes smaller and smaller (a property of the $`RM^{1/3}`$ law of degenerate configurations). For rotating systems described in the microcanonical ensemble (fixed $`E`$ and $`L`$), the resulting decrease of moment of inertia is accompanied by an increase of angular velocity as shown in Fig. 15. In this diagram, the $`\omega \mathrm{\Lambda }`$ spiral of Fig. 6 is completely unwound and the angular velocity increases monotonically with $`\mathrm{\Lambda }`$.

In Fig. 16 we have focused our attention to what happens close to the minimum energy. For $`EE_{min}`$, $`T0`$ and the system has the same structure as a cold white dwarf star (see Appendix B). For our perturbative analysis to be valid, we have taken a very small angular velocity $`\lambda =0.01`$. We shall describe the diagram of Fig. 16 in the $`(L,T)`$ ensemble in which the angular momentum and the temperature are assumed given (this is a situation intermediate between microcanonical and canonical ensembles). For $`T>T_c`$, the system is in a gaseous phase with a smoothly decreasing density profile. For $`T<T_c`$, the system undergoes an isothermal collapse that only stops when gravity is balanced by the degeneracy pressure. The result of this phase transition is a “fermion ball” which contains almost all the mass, unlike in the microcanonical ensemble at the point of gravothermal catastrophe. For rotating systems, this isothermal collapse is accompanied by a discontinuous rise of angular velocity at the critical temperature $`T_c`$ (see Fig. 17). Therefore, even if the initial rotation of the system is negligible in the gaseous phase, after collapse the “fermion ball” can have appreciable rotation as suggested in Fig. 17. Its structure is then similar to a distorted polytrope of index $`n=3/2`$ as computed by Chandrasekhar (1933) in the limit of slow rotation.

## 6 Conclusion

In this paper, we have considered the effect of a small rotation on the thermodynamical stability of self-gravitating systems. We have worked in a finite box in order to make a clear connexion with the Antonov problem for non-rotating systems and render the statistical mechanics of these objects rigorous. Physically, this idealization means that our isothermal system is surrounded by a medium which fixes its size. We have found that for rotating systems, the well-known inequivalence of statistical ensembles for self-gravitating systems manifests itself in a striking manner: the instability is advanced in the microcanonical ensemble and delayed in the canonical one. In addition, we have found a connexion between the onset of instability and the configuration of maximum flattening in the series of equilibria. These results have be generalized to the case of self-gravitating fermions.

On the other hand, the case of rapidly rotating isothermal configurations is interesting because new, non trivial, structures can emerge as maximum entropy states (i.e., most probable states). The classification of such structures is complicated because many bifurcations can occur depending on the values of the control parameters $`(M,E,L)`$. In particular, we must be careful to select only entropy maxima, discarding the critical points of entropy which are only saddle points. This can be obtained either by solving the Boltzmann-Poisson equation (Votyakov et al. 2002) and checking the stability of the solutions or by using relaxation equations towards the maximum entropy state (Chavanis et al. 1996). These studies are important in order to obtain a classification of the most probable configurations of self-gravitating systems. The second approach (relaxation methods) will be explored in a future paper (Rosier & Chavanis, in preparation).

###### Acknowledgements.

I acknowledge interesting discussions with J. Katz. I have also benefited from stimulating interactions with the participants of the Les Houches School of Physics on “Dynamics and thermodynamics of systems with long range interactions” (February 2002).

## Appendix A The entropy of a rotating self-gravitating Fermi gas

The Fermi-Dirac entropy (5) can be written in the equivalent form

$`S={\displaystyle \left\{\frac{f}{\eta _0}\mathrm{ln}\left(\frac{f/\eta _0}{1f/\eta _0}\right)+\mathrm{ln}\left(1\frac{f}{\eta _0}\right)\right\}d^3𝐫d^3𝐯}.`$
(98)
Inserting the distribution function

$`f={\displaystyle \frac{\eta _0}{1+ke^\mathrm{\Psi }e^{\beta \frac{w^2}{2}}}},`$ (99)
in Eq. (98), we obtain

$`\eta _0S=M\mathrm{ln}k+{\displaystyle \rho \mathrm{\Psi }d^3𝐫}+\beta K_{th}`$
$`+\eta _0{\displaystyle \mathrm{ln}\left(1+\frac{1}{k}e^\mathrm{\Psi }e^{\beta \frac{w^2}{2}}\right)d^3𝐫d^3𝐯}.`$ (100)
The last integral can be reduced by an integration by parts to the expression $`(2/3)\beta K_{th}`$. Therefore,

$`\eta _0S=M\mathrm{ln}k+{\displaystyle \rho \mathrm{\Psi }d^3𝐫}+{\displaystyle \frac{5}{3}}\beta K_{th}.`$ (101)
Now, using the definition of $`\mathrm{\Psi }`$, we get

$`\eta _0S=M\mathrm{ln}k+2\beta W\beta K_{rot}M\beta \mathrm{\Phi }_0+{\displaystyle \frac{5}{3}}\beta K_{th}.`$ (102)
Since $`W=EK`$ and $`K_{th}=KK_{rot}`$, the foregoing expression can be rewritten

$`\eta _0S=M\mathrm{ln}k+2\beta E{\displaystyle \frac{1}{3}}\beta KM\beta \mathrm{\Phi }_0{\displaystyle \frac{8}{3}}\beta K_{rot},`$ (103)
or, according to Eqs. (58) and (61),

$`\eta _0S=M\mathrm{ln}k+{\displaystyle \frac{7}{3}}\beta E{\displaystyle \frac{1}{3}}\beta {\displaystyle p𝐫𝑑𝐒}M\beta \mathrm{\Phi }_0{\displaystyle \frac{4}{3}}\beta 𝐋𝛀.`$
(104)
The pressure at the surface of the sphere can be calculated with Eqs. (55) and (99). After simplification, we get

$`{\displaystyle \frac{\eta _0S}{M}}=\mathrm{ln}k{\displaystyle \frac{7}{3}}\eta \mathrm{\Lambda }\beta \mathrm{\Phi }_0{\displaystyle \frac{4}{3}}\eta \lambda \omega `$
$`{\displaystyle \frac{2}{9}}\eta {\displaystyle \frac{\alpha ^{10}}{\mu ^4}}\left\{{\displaystyle \frac{I_{3/2}(ke^{\psi (\alpha )})}{I_{1/2}(k)^5}}{\displaystyle \frac{3}{2}}v{\displaystyle \frac{I_{1/2}(ke^{\psi (\alpha )})}{I_{1/2}(k)^5}}\varphi _0(\alpha )\right\}.`$ (105)
Using Eq. (70), which remains valid for self-gravitating fermions, and Eq. (93) we finally arrive at the form (97).

## Appendix B The energy of a rotating polytrope

In this Appendix, we derive a simple analytic formula for the potential energy of a slowly rotating polytrope of index $`n`$. The index $`n=3/2`$ describes a completely degenerate Fermi gas at zero temperature, which is a particular limit of the model studied in Sec. 5.

For a self-gravitating system rotating with constant angular velocity $`\mathrm{\Omega }`$, the condition of hydrostatic equilibrium in the rotating frame can be written

$`p=\rho \mathrm{\Phi }_{eff},`$ (106)
where $`\mathrm{\Phi }_{eff}`$ is the effective potential defined in Eq. (11). Now, for a polytropic gas,

$`p=K\rho ^{1+\frac{1}{n}},`$ (107)
and the equation of hydrostatic equilibrium can be rewritten as

$`(n+1)\left({\displaystyle \frac{p}{\rho }}\right)=\mathrm{\Phi }_{eff}.`$ (108)
This equation integrates to give

$`(n+1)p=\rho \left(\mathrm{\Phi }_{eff}^b\mathrm{\Phi }+{\displaystyle \frac{1}{2}}|𝛀\times 𝐫|^2\right),`$ (109)
where $`\mathrm{\Phi }_{eff}^b`$ is the effective potential at the surface of the polytrope. On integrating Eq. (109) over the volume of the configuration, we obtain

$`(n+1){\displaystyle pd^3𝐫}=M\mathrm{\Phi }_{eff}2W+{\displaystyle \frac{1}{2}}𝐋𝛀.`$ (110)
Inserting this last relation in Eq. (56), and recalling that the pressure vanishes on the surface of a polytrope, we get

$`W={\displaystyle \frac{1}{5n}}\left[3M\mathrm{\Phi }_{eff}^b+\left(n+{\displaystyle \frac{5}{2}}\right)𝐋𝛀\right].`$ (111)
This general expression was derived by Chandrasekhar (1961) in a slightly different manner. We now consider the case of slowly rotating polytropes. In that case, an explicit expression for $`\mathrm{\Phi }_{eff}^b`$ can be deduced from the study of Chandrasekhar (1933) on distorted polytropes. The density profile can be written

$`\rho =\rho _0\mathrm{\Theta }^n,`$ (112)
with

$`\mathrm{\Theta }=\theta +v\left\{\psi _0(\xi )+{\displaystyle \underset{j=1}{\overset{+\mathrm{}}{}}}A_j\psi _j(\xi )P_j(\mu )\right\}.`$ (113)
The quantities have their usual meaning: $`\rho _0`$ is the central density, $`v=\mathrm{\Omega }^2/2\pi G\rho _0`$, $`\theta `$ is the solution of the Lane-Emden equation of index $`n`$ and $`P_j(\mu )`$ are Legendre polynomials. The functions $`\psi _j(\xi )`$ have been defined and tabulated by Chandrasekhar (1933). The inner and outer potentials are given by

$`\mathrm{\Phi }_{int}=\mathrm{\Phi }_{eff}^bR\left[\mathrm{\Theta }{\displaystyle \frac{1}{6}}v\xi ^2(1P_2(\mu ))\right],`$ (114)
$`\mathrm{\Phi }_{ext}=R\left[{\displaystyle \frac{C_0}{\xi }}+v{\displaystyle \underset{j=1}{\overset{+\mathrm{}}{}}}{\displaystyle \frac{C_j}{\xi ^{j+1}}}P_j(\mu )\right],`$ (115)
with $`R=(n+1)K\rho _0^{1/n}`$. The constant $`A_j`$, $`C_j`$ are determined by requiring the continuity of $`\mathrm{\Phi }`$ and its derivative on a sphere of radius $`\xi _1`$, the first zero of the Emden’s function with index $`n`$. For our purposes, we need to consider the equations obtained for $`j=0`$:

$`R{\displaystyle \frac{C_0}{\xi _1}}=Rv\psi _0(\xi _1){\displaystyle \frac{1}{6}}Rv\xi _1^2\mathrm{\Phi }_{eff}^b,`$ (116)
$`R{\displaystyle \frac{C_0}{\xi _1^2}}=R\theta _1^{}+Rv\psi _0^{}(\xi _1){\displaystyle \frac{1}{3}}Rv\xi _1.`$ (117)
In Eq. (116), we have used $`\theta (\xi _1)=0`$. Now, the relation between the mass and the central density is given by (Chandrasekhar 1933, Eq. (40)):

$`M=4\pi \left[{\displaystyle \frac{(n+1)K}{4\pi G}}\rho _0^{\frac{3n}{3n}}\right]^{3/2}\xi _1^2\theta _1^{}\left[1+v{\displaystyle \frac{\frac{1}{3}\xi _1\psi _0^{}(\xi _1)}{|\theta _1^{}|}}\right].`$
(118)
On the other hand, we define the characteristic radius $`R_\mathrm{\Omega }`$ of a rotating polytrope by the relation (see Chandrasekhar 1933, Eq. (9)):

$`R_\mathrm{\Omega }=\left[{\displaystyle \frac{(n+1)K}{4\pi G}}\rho _0^{\frac{1}{n}1}\right]^{1/2}\xi _1.`$ (119)
Inserting Eq. (117) in Eq. (118) and using Eq. (119), we obtain

$`R{\displaystyle \frac{C_0}{\xi _1}}={\displaystyle \frac{GM}{R_\mathrm{\Omega }}}.`$ (120)
On the other hand, we can check that

$`Rv={\displaystyle \frac{2\mathrm{\Omega }^2}{\xi _1^2}}R_\mathrm{\Omega }^2.`$ (121)
Substituting these results in Eq. (116), we find that the effective potential at the surface of a slowly rotating polytrope is given by

$`\mathrm{\Phi }_{eff}^b={\displaystyle \frac{GM}{R_\mathrm{\Omega }}}+2\mathrm{\Omega }^2R_\mathrm{\Omega }^2\left({\displaystyle \frac{\psi _0(\xi _1)}{\xi _1^2}}{\displaystyle \frac{1}{6}}\right).`$ (122)
We now need to relate the characteristic radius $`R_\mathrm{\Omega }`$ to the total mass $`M`$ of the system. Eliminating the central density $`\rho _0`$ between Eq. (118) and Eq. (119), and expressing $`v`$ in terms of $`R_\mathrm{\Omega }`$ with the aid of Eq. (121), we obtain to order $`\mathrm{\Omega }^2`$:

$`R_\mathrm{\Omega }=R_0\left[1+\mathrm{\Omega }^2{\displaystyle \frac{S_nK^{\frac{3n}{3n}}}{2\pi G^{\frac{2n+3}{3n}}M^{\frac{2n}{3n}}}}\right],`$ (123)
where $`R_0`$ is the radius of a non-rotating polytrope. It is related to its mass $`M`$ by the relation

$`GN_nM^{\frac{n1}{n}}R_0^{\frac{3n}{n}}=K.`$ (124)
The constant $`N_n`$ is defined and tabulated in Chandrasekhar (1942) for different values of the polytropic index $`n`$. The constant $`S_n`$ appearing in Eq. (123) is given by

$`S_n={\displaystyle \frac{n1}{3n}}\left({\displaystyle \frac{n+1}{4\pi }}\right)^{\frac{n}{1n}}N_n^{\frac{2n^2}{(1n)(3n)}}\xi _1^{\frac{2n}{1n}}{\displaystyle \frac{\frac{1}{3}\xi _1\psi _0^{}(\xi _1)}{|\theta _1^{}|}},`$
(125)
and typical values are listed in Table 1.

From the above results, we can obtain an explicit expression for the potential energy of a rotating polytrope. For an axisymmetrical system, we have the relation $`𝐋=I𝛀`$ where $`I`$ is the axial moment of inertia (50). To our order of approximation, we just need to determine the value of $`I`$ for a non-rotating polytrope. Using the relations $`\rho =\rho _0\theta ^n`$ and $`r=\xi R_0/\xi _1`$, expressing the central density as a function of $`R_0`$ by the relation (119) and using the mass-radius relation (124), we obtain

$`I=_nMR_0^2,`$ (126)
with

$`_n={\displaystyle \frac{8\pi }{3}}\left[{\displaystyle \frac{4\pi }{(n+1)N_n}}\right]^{\frac{n}{1n}}\xi _1^{\frac{53n}{n1}}{\displaystyle _0^{\xi _1}}\theta ^n\xi ^4𝑑\xi .`$ (127)
Therefore, according to Eqs. (111) (122) and (126), the potential energy of a slowly rotating polytrope of index $`n`$ is given by

$`W={\displaystyle \frac{1}{5n}}\left({\displaystyle \frac{3GM^2}{R_\mathrm{\Omega }}}+Q_nMR_0^2\mathrm{\Omega }^2\right),`$ (128)
with

$`Q_n=6{\displaystyle \frac{\psi _0(\xi _1)}{\xi _1^2}}1+\left(n+{\displaystyle \frac{5}{2}}\right)_n,`$ (129)
where $`R_0`$ and $`R_\mathrm{\Omega }`$ are given by Eqs. (124) and (123) respectively. For $`\mathrm{\Omega }=0`$, Eq. (128) reduces to the well-known Ritter’s formula (see Chandrasekhar 1942). Note also that the dimensionless rotation parameter can be written

$`v=Z_n{\displaystyle \frac{\mathrm{\Omega }^2R_0^3}{GM}},`$ (130)
with

$`Z_n={\displaystyle \frac{1}{2\pi }}\left({\displaystyle \frac{n+1}{4\pi }}\right)^{\frac{n}{1n}}N_n^{\frac{n}{1n}}\xi _1^{\frac{2n}{1n}}.`$ (131)
The polytropic indices $`n=1`$ and $`n=3`$ are special. For $`n=1`$, the problem can be solved analytically. The following results are well-known (Chandrasekhar 1942, 1933): $`\xi _1=\pi `$, $`\theta _1^{}=1/\pi `$, $`\psi _0(\xi _1)=1`$, $`\psi _0^{}(\xi _1)=1/\pi `$, $`N_1=2/\pi `$. According to Eq. (120), the radius $`R_\mathrm{\Omega }=R_0=(K/2\pi G)^{1/2}\xi _1`$ is independant of the central density and is uniquely determined by the value of $`K`$. The mass $`M`$ is arbitrary and independant of the radius. Eqs. (126) (130) remain valid with $`_1=2(\pi ^26)/3\pi ^3`$ and $`Z_1=2/\pi ^2`$. The other quantities can be obtained from these results (see Table 1). The case $`n=3`$ corresponds to the reverse situation. For $`\mathrm{\Omega }=0`$, the mass is given by $`M_0=(K/GN_3)^{3/2}`$ and the radius is arbitrary. For $`\mathrm{\Omega }0`$, we can still consider that the radius $`R_\mathrm{\Omega }=R_0`$ is a free parameter and determine the mass according to Eq. (118), which can be written

$`M=M_0\left[1+\mathrm{\Omega }^2S_3^{}{\displaystyle \frac{G^{1/2}R_\mathrm{\Omega }^3}{2K^{3/2}}}\right],`$ (132)
with

$`S_3^{}={\displaystyle \frac{\pi ^{1/2}}{\xi _1^3}}{\displaystyle \frac{\frac{1}{3}\xi _1\psi _0^{}(\xi _1)}{|\theta _1^{}|}}0.03309.`$ (133)
We now consider the case of a degenerate Fermi gas at zero temperature. As is well-known, this system is equivalent to a polytrope of index $`n=3/2`$. In addition, the constant $`K`$ is explicitly given by

$`K={\displaystyle \frac{1}{5}}\left({\displaystyle \frac{3}{4\pi \eta _0}}\right)^{2/3},`$ (134)
where $`\eta _0`$ is the distribution function of a completely degenerate Fermi gas. According to the Virial theorem (61), the total energy $`E=K+W`$ of this system is $`E=\frac{W}{2}`$. Therefore, using Eq. (128), we have

$`E={\displaystyle \frac{3GM^2}{7R_\mathrm{\Omega }}}+{\displaystyle \frac{1}{7}}Q_{3/2}MR_0^2\mathrm{\Omega }^2.`$ (135)
According to Eq. (124), the mass-radius relation of a non rotating “fermion ball” is

$`MR_0^3={\displaystyle \frac{\chi }{\eta _0^2G^3}},`$ (136)
with

$`\chi {\displaystyle \frac{1}{125}}\left({\displaystyle \frac{3}{4\pi }}\right)^2{\displaystyle \frac{1}{N_{3/2}^3}}\mathrm{5.9723\; 10}^3.`$ (137)
For a rotating fermion ball we have, according to Eq. (123),

$`R_\mathrm{\Omega }=R_0\left(1+{\displaystyle \frac{\kappa \mathrm{\Omega }^2}{2\pi G^4M^2\eta _0^2}}\right),`$ (138)
with

$`\kappa {\displaystyle \frac{1}{125}}\left({\displaystyle \frac{3}{4\pi }}\right)^2S_{3/2}\mathrm{3.9813\; 10}^3.`$ (139)
Then, the energy (135) can be expressed in terms of $`M`$, $`G`$, $`\eta _0`$ and $`\mathrm{\Omega }`$. In order to make the link with the variables introduced in Sec. 5, which are normalized by the box radius $`R`$, we note the relation

$`{\displaystyle \frac{R_0}{R}}=(512\pi ^4\chi )^{1/3}\mu ^{2/3}={\displaystyle \frac{6.6784}{\mu ^{2/3}}},`$ (140)
which directly results from Eq. (136) and the definition of the degeneracy parameter $`\mu `$.


format: text

text:
# Gamma-Ray Burst Afterglow emission with a decaying magnetic field

## 1 Introduction

It is widely accepted that GRBs afterglow emission is due to non-thermal radiation mechanisms produced in relativistic shocks, as the ejected matter from a highly energetic explosion expands into the surrounding medium. With current facilities we can observe at photon energies up to 10 keV and the modeling of the observed lightcurves shows in general a predominance of synchrotron emission over inverse Compton scattering on the synchrotron photons. (Panaitescu & Kumar 2001 (thereafter PK01)). Only in one case, GRB000926, has a Compton component been inferred (Harrison et al. 2001). This fact sets an upper limit on the density of the medium surrounding the explosion, whose properties can help us to unveil GRBs progenitors. In the standard framework the dynamics of the remnant expansion follows an adiabatic law. This behavior is expected when the electrons cannot cool on the expansion timescale or (even when they can) the bulk of the energy remains in relativistic protons and magnetic fields. The post shock magnetic field, necessary for the synchrotron emission, is mainly generated during the shock itself, because the shock compression of the pre-existing field alone would lead to a negligible magnetic energy per particle (e.g. Gruzinov 2001). The magnetic field is normally assumed to be uniform and its lengthscale to be of the order of the remnant scale. The electron energy distribution behind the expanding shock front is given by a power law, down to a minimum injection energy. Modeling eight afterglow lightcurves within this framework has resulted in a wide range of density from $`10^3`$ to $`27`$ cm-3, strongly indicating a low density environment (PK01). This is inconsistent with other independent observations placing GRBs birth in a star forming region: SN bumps in afterglow lightcurves (e.g. Bloom et al. 1999, Lazzati et al 2001, Bloom et al. 2002), metal lines (Antonelli et al 2000, Piro et al 2000, Reeves et al 2002), afterglow localizations close to star forming regions (Bloom, Kulkarni, Djorgovski 2002, Jaunsen et al 2002) and estimations of high density environment from the temporal evolution of external column density (Lazzati & Perna 2002). Another remarkable inference from this model is that the fraction of internal energy given to electrons ($`ϵ_e`$) is always much larger then the fraction given to magnetic field ($`ϵ_b`$). This implies that the magnetic field is relatively weak and the Compton emission dominates the overall cooling of the electrons.

Recently many publications have investigated the somewhat poorly understood generation of magnetic field in strong relativistic collisionless shocks (e.g. Kazimura et al 1998, Gruzinov and Waxman 1998, Medvedev & Loeb 1999, Gruzinov 2001). Kazimura et al (1998) found generation of a small-scale quasi-static magnetic field; similar results were found by Gruzinov (2001) at small simulation times but his longer run shows that the field quickly decays by Landau damping; as a result synchrotron emission could be insignificant. Medvedev & Loeb (1999) found more reassuring results. Relativistic two-stream magnetic instability can, they claim, generate stable randomly oriented strong magnetic field in the plane of the collisionless shock front, so synchrotron emission is possible. These contradictory results show our ignorance of the possible acceleration process, and of how (or if) a stable magnetic field is created there. In the standard model this ignorance is ”parameterized “ through $`ϵ_e`$ and $`ϵ_B`$: the simplistic assumption is made that a persistent magnetic field makes these parameters uniform in space and time. An observation that maybe suggests a time dependence of $`ϵ_e`$ and $`ϵ_B`$ is the independent determination of these parameters in GRB 970508 afterglow by Wijers and Galama (1999) at 12 days and Frail, Waxman and Kulkarni (2000) at $`1`$ year: while the former found $`ϵ_e=0.12`$ and $`ϵ_B=0.089`$ the latter inferred $`ϵ_eϵ_B=0.5`$.

These uncertainties, coupled with the apparent problems with the standard model, motivate us to explore the consequences for afterglow emission of a different scenario, where the propagation of the magnetic field to large scale is disfavored. In this contest we particularly focus on the possibility that external density could be larger, keeping the total energy to a reasonable value, without loosing a general agreement with data. In fact, at this stage, we deal with sparse set of measurements, almost never simultaneous at all wavelengths. These data can be modeled in different frameworks, (e.g. GRB 970508 (Frail et al 2000, Chevalier & Li 2000), GRB 980519 (Frail et al 2000b, Chevalier & Li 1999), GRB 000301C (Berger et al 2000, Li & Chevalier 2001), GRB 991208 (Galama et al 2000, Li & Chevalier 2001 and GRB 000418 (Berger et al 2001)). This analysis is meant to investigate how robust are the current estimation of shock and density parameters.

In this paper we present this more general theory for the GRB afterglow emission (§2); our conclusion (§3) is that a qualitative comparison of this model with data suggests the possibility that a shorter magnetic lengthscale can be involved, affecting the estimates of parameters, but only a quantitative modeling of data can give the conclusive answer.

## 2 Short magnetic lengthscale theory

The assumption of a post shock magnetic field persisting during all the expansion time has three main consequences for the synchrotron (S) emission and spectrum. First the emitting region linear dimension is $`Rc\times t_{exp}`$, where $`t_{exp}`$ is the expansion time in the lab frame; this enters in the calculation of the observed peak flux $`F_p`$. Then the observed cooling frequency $`\nu _c`$ is computed using the Lorenz factor, $`\gamma _c`$, of the electrons that cool radiatively on a timescale equal to the remnant age. Finally the lengthscale for synchrotron self-absorption has to be compared with the dimension of the fireball in order to compute the observed absorption frequency $`\nu _a`$. The only break frequency which is not affected by any assumption about the extension of the magnetic field, is the observed peak frequency $`\nu _m`$, corresponding to the minimum random Lorentz factor $`\gamma _m`$: it depends crucially only on the fractions of internal energy given to electrons, $`ϵ_e`$, and given to the magnetic field $`ϵ_b`$. Throughout the entire emission volume electrons cool also via inverse Compton (IC) on the synchrotron photons (Self Synchrotron Compton: SSC). This cooling process dominates the electrons cooling if the Compton parameter $`Y=L_{IC}/L_S`$ (IC and S luminosity ratio) is greater than one. This happens if $`ϵ_bϵ_e`$ and current estimates find in most objects $`ϵ_b/ϵ_e10^110^2`$, (Wijers & Galama 1999, Granot et al 1999, PK01).

If the magnetic field behind the shock persists for an average timescale of $`t_b<t_{exp}`$, its lengthscale is $`\lambda c\times t_b<R`$. In principle $`\lambda `$ could depend on time, however we do not have any consolidated understanding of this process, therefore for simplicity we will assume throughout this paper that $`\lambda /R=\delta `$ is constant during the afterglow phase. At each time the emitting region is then divided in two zones (see Fig 1). The first zone (region A) is the layer with dimension $`\lambda `$ just behind the shock front, where there are newly accelerated electrons and magnetic field: here electrons emit via S and SSC; the second zone (region B) of linear dimension $`R\lambda `$, is filled by electrons that have just cooled by S and SSC emissions and that now Compton scatter the radiation coming from the first layer. We suppose that the newly accelerated electrons are advected with the flow; they produce synchrotron radiation only for a short time (with a cooling break therefore only at very high energies) in region A near the shock front and then they emit most of the IC radiation while downstreaming through region B. On the other hand if the electrons diffused freely, so that they could move back into the magnetized region after having lost most of their energy via IC, then the synchrotron spectrum would have a lower-energy break. This is probably less likely therefore we concentrate on the first case. The total observed spectrum then results from the sum of region A and region B spectra.

### 2.1 region A: S and SSC emission

The consequences of $`\lambda <R`$ on the S spectrum in region A can be derived easily from what we mentioned before. The observed peak frequency $`\stackrel{~}{\nu _m}`$ remains unchanged,

$$\stackrel{~}{\nu _m}=\nu _m,$$
(1)
and

$$\stackrel{~}{\gamma _m}\stackrel{~}{\nu _m}^{1/2}=\gamma _m.$$
(2)
In our generalized model with $`\delta <1`$ we denote quantities by the symbol ’$`\stackrel{~}{q}`$’, where ’$`q`$’ is defined for $`\delta =1`$, (for the explicit definitions of quantities in the standard model we refer to Panaitescu & Kumar 2000).

The comoving peak intensity is $`I^{}\stackrel{~}{n}^{}B^{}\lambda ^{}`$, where $`\stackrel{~}{n}^{}`$, $`B^{}`$ and $`\lambda ^{}`$ are respectively the external density, the magnetic field and the emitting region linear dimension in the comoving frame; therefore the observed peak flux $`\stackrel{~}{F_p}`$ is smaller than $`F_p`$ by a factor $`\delta `$:

$$\stackrel{~}{F_p}\mathrm{\Gamma }^3(\stackrel{~}{n}^{}B^{}\lambda ^{})\frac{R^2}{\mathrm{\Gamma }^2}=\delta \times F_p,$$
(3)
where $`\mathrm{\Gamma }`$ is the bulk Lorentz factor. Since the magnetic field persists for a shorter time, a smaller fraction of electrons cool significantly, before the magnetic field disappears in that region. The Lorentz factor of the electrons, whose cooling time is equal to the typical timescale of the system, $`\stackrel{~}{\gamma _c}`$ is therefore higher:

$$\stackrel{~}{\gamma _c}\frac{1}{t_b^{}B^2(1+\stackrel{~}{Y})}=\gamma _c\times 𝒴\delta ^1,$$
(4)
where $`𝒴=(1+Y)/(1+\stackrel{~}{Y})`$. The observed cooling frequency $`\stackrel{~}{\nu _c}B^{}\stackrel{~}{\gamma _c}^2\mathrm{\Gamma }`$ is then

$$\stackrel{~}{\nu _c}=\nu _c\times 𝒴^2\delta ^2.$$
(5)
The synchrotron self-absorption optical thickness can be approximated as

$$\stackrel{~}{\tau _a}(\nu )\left(\frac{\stackrel{~}{n}\lambda }{B^{}\stackrel{~}{\gamma _p}^5}\right)\left(\nu /\stackrel{~}{\nu _p}\right)^{5/3}=\stackrel{~}{\tau _p}\left(\nu /\stackrel{~}{\nu _p}\right)^{5/3},$$
(6)
for $`\nu <\stackrel{~}{\nu _p}`$, where $`p=c`$ when the electrons are radiative ($`\stackrel{~}{\gamma _c}<\stackrel{~}{\gamma _m}`$) and $`p=m`$ when they are adiabatic ($`\stackrel{~}{\gamma _c}>\stackrel{~}{\gamma _m}`$). We will refer to the first case as fast-cooling regime and to the second case as slow-cooling regime. The absorption frequency corresponds to $`\stackrel{~}{\tau _a}(\stackrel{~}{\nu _a})=1`$ and in slow cooling regime is given by

$$\stackrel{~}{\nu _a}=\nu _a\times \delta ^{3/5},$$
(7)
where we use Eq. 1 and Eq. 6; in fast cooling regime Eqs 5, 6 and 18 ($`𝒴=1`$) give

$$\stackrel{~}{\nu _a}=\nu _a\times \delta ^{8/5}.$$
(8)
Eqs 37 and 8 are calculated in the case $`\stackrel{~}{\nu }_a<\stackrel{~}{\nu }_p`$ (see Granot & Sari 2001 for the different possible S spectra). Nevertheless $`\stackrel{~}{\nu }_a\stackrel{~}{\nu }_p`$ for high densities

$$\stackrel{~}{n}7\times 10^2\delta ^{12/7}E_{53}^{12/10}(1+\stackrel{~}{Y})^3ϵ_{b,2}^{27/10}$$
(9)
in fast cooling and in slow cooling

$$\stackrel{~}{n}8.2\times 10^5\delta ^1E_{53}^{1/2}ϵ_{b,2}^{1/2}ϵ_{e,2}^5(1+z)t_{day}^{5/2}$$
(10)
where $`E`$ is the isotropic equivalent energy and $`t_{day}`$ is the observed time in days. The inferred external density is reasonable for $`\delta =0.1`$ ($`\stackrel{~}{n}3.6\times 10^4`$cm-3 and $`\stackrel{~}{n}8.5\times 10^6t_{day}^{5/2}`$cm-3 respectively) and for $`\delta =0.01`$ ($`\stackrel{~}{n}2\times 10^6`$cm-3 and $`\stackrel{~}{n}8.5\times 10^7t_{day}^{5/2}`$cm-3) but for lower $`\delta `$ the density lower limits expressed in Eq 9 and Eq 10 are perhaps too high for the afterglow environment. For simplicity we focus thereafter on $`\stackrel{~}{\nu }_a<\stackrel{~}{\nu }_c`$ but the reader should keep in mind that these cases hold as long as Eq 9 and Eq 10 are not satisfied.

The magnetic energy $`\frac{B^2}{8\pi }n\mathrm{\Gamma }^2`$ and the Compton parameter are higher at the beginning of the afterglow evolution and so the electrons are more likely then to be radiative. Since the minimum injected $`\gamma _m\mathrm{\Gamma }`$ decreases with time while $`\gamma _c`$ increases as the efficiency of radiative loss decreases, the electrons undergo a transition to the adiabatic regime. The transition time between the fast cooling and the slow cooling regime is

$$\stackrel{~}{T}_{fs}=Tfs\times \delta ^2;$$
(11)
for

$$\stackrel{~}{n}<40\delta ^2E_{53}^1(Y_F+1)^2ϵ_{B,2}^2ϵ_{e,1}^2(1+z)^1t_{day}^1$$
(12)
the slow-cooling regime dominates the fireball evolution during times when the afterglow is observed ($`\stackrel{>}{}1`$ day).

The accelerated electrons that upscatter synchrotron photons give rise to a SSC component in the spectrum (for the SSC in the standard model see Sari & Esin 2000). During the afterglow phase the medium has a Thomson optical depth

$$\stackrel{~}{\tau }=\frac{1}{3}\stackrel{~}{n}\sigma _TR\delta $$
(13)
of the order of $`\stackrel{~}{\tau }10^5\times \delta `$ for $`\stackrel{~}{n}=10^2`$cm-3 ($`\sigma _T`$ is the Thomson cross section), therefore only a minor fraction of the emission is upscattered, and the synchrotron flux is not significantly altered. Nevertheless if $`\stackrel{~}{Y}1`$ the electron cooling is Compton dominated and SSC upscatters low energy photons, that might be detected as a Compton component in x-ray and, in the future, in higher energy band spectra (Sari, Narayan & Piran 1996). The parameter $`\stackrel{~}{Y}`$ can be defined as the mean number of scatterings times the average fractional energy change per scattering, gained by a photon traversing a finite medium (Rybicki & Lihgtman 1979). For a Thomson thin, relativistic medium $`\stackrel{~}{Y}`$ is

$$\stackrel{~}{Y}\frac{4}{3}\stackrel{~}{\tau }\stackrel{~}{\gamma }^2N(\stackrel{~}{\gamma })𝑑\stackrel{~}{\gamma }$$
(14)
where $`N(\stackrel{~}{\gamma })`$ is the normalized $`\stackrel{~}{\gamma }`$ distribution (see e.g. PK00 for explicit formulae for $`\nu _a<\nu _p`$). For $`\stackrel{~}{Y}>1`$ one may need to consider more than one scattering per photon; however the energy of a photon in the rest frame of the second scattering electron is a factor $`\gamma ^3`$ higher than the original energy, generally exceeding the Thompson regime limit (511 keV); the scattering cross section is then substantially reduced (Klein-Nishina (KN) regime) and further scatterings inhibited. As a consequence we will consider only single scattering of S photons. The KN effect limits also the electrons that mostly contribute to the SSC luminosity. They have a Lorentz factor

$$\stackrel{~}{\gamma }\stackrel{~}{\gamma }_{KN}511\mathrm{\Gamma }/h\stackrel{~}{\nu }_{Ep},$$
(15)
where $`\stackrel{~}{\nu }_{Ep}`$ is the S energy peak frequency ($`\stackrel{~}{\nu }_c`$ in slow cooling and $`\stackrel{~}{\nu }_m`$ in fast cooling). In zone A the S spectrum can have a very high $`\stackrel{~}{\nu _c}`$ (see Eq 5); consequently in slow cooling $`\stackrel{~}{\gamma }_c`$ could be greater than $`\stackrel{~}{\gamma }_{KN}`$ and the Klein-Nishina effect substantially decreases $`\stackrel{~}{Y}`$ (we will refer to this situation as “zone A in KN regime”). Moreover if $`\stackrel{~}{Y}1`$ (which is the case in most simulations) we can omit the contribution from $`\stackrel{~}{\gamma }>\stackrel{~}{\gamma }_{KN}`$. Eq. 14 for $`2<p<3`$ gives then

$$\stackrel{~}{Y_S}(1+\stackrel{~}{Y_S})^{(3p)}=Y_S(1+Y_S)^{(3p)}\times \delta ^{(p2)},$$
(16)
for $`\stackrel{~}{\gamma }_c\stackrel{~}{\gamma }_{KN}`$ and

$$\stackrel{~}{Y_S}(1+\stackrel{~}{Y_S})^{(2(p3))}=Y_S(1+Y_S)^{(2(p3))}\times \delta ^{(73p)},$$
(17)
for $`\stackrel{~}{\gamma }_c\stackrel{~}{\gamma }_{KN}`$. If the electrons are in the fast cooling regime it is very unlikely that $`\stackrel{~}{\gamma }_m>\stackrel{~}{\gamma }_{KN}`$ and Eq. 14 gives

$$\stackrel{~}{Y_F}=\frac{4}{3}\stackrel{~}{\tau }\stackrel{~}{\gamma }_m\stackrel{~}{\gamma _c}=Y_F;$$
(18)
note that $`\stackrel{~}{Y_F}`$ does not depend on $`\delta `$, because most of the electrons cool before the field disappears. Given a S seed spectrum $`\stackrel{~}{F}_S(\nu _s)`$, the resulting single scattering IC spectrum is (Rybicki & Lightman, 1979)

$$\stackrel{~}{F_\nu }^{IC}=A_0^1\stackrel{~}{F}_S\left(\frac{\nu }{4\stackrel{~}{\gamma }^2X}\right)f(X)𝑑X_{\stackrel{~}{\gamma }_1}^{\stackrel{~}{\gamma }_2}N(\stackrel{~}{\gamma })𝑑\stackrel{~}{\gamma },$$
(19)
where $`A=3\sigma _TR\delta `$ (in this model), $`X=\nu /4\stackrel{~}{\gamma }^2\nu _s`$ and $`f(X)=2X\mathrm{ln}X+X2X+1`$ includes the exact cross section angular dependence in the limit $`\stackrel{~}{\gamma }1`$ (Blumenthal and Gold,1970). For the cases we are treating it is possible to use for the first integral the analytic expressions given by Sari & Esin (2000) in their Appendix A. The IC spectral peak flux and break frequencies for $`\delta <1`$ are related to the IC spectrum for $`\delta =1`$ in the following way

$$\stackrel{~}{F_p}^{IC}\stackrel{~}{\tau }\stackrel{~}{F_p}=F_p^{IC}\times \delta ^2$$
(20)
using Eq 3 and Eq 6

$$\stackrel{~}{\nu }_m^{IC}2\stackrel{~}{\gamma _m}^2\stackrel{~}{\nu _m}=\nu _m^{IC},$$
(21)
using Eq. 21;

$$\stackrel{~}{\nu }_c^{IC}2\stackrel{~}{\gamma _c}^2\stackrel{~}{\nu _c}=\nu _c^{IC}\times \left(𝒴\delta ^1\right)^4,$$
(22)
using Eq. 45;

$$\stackrel{~}{\nu }_a^{IC}2\stackrel{~}{\gamma _m}^2\stackrel{~}{\nu _a}=\nu _a^{IC}\times \delta ^{3/5},$$
(23)
in the slow cooling regime using Eq. 27 and

$$\stackrel{~}{\nu }_a^{IC}2\stackrel{~}{\gamma _c}^2\stackrel{~}{\nu _a}=\nu _a^{IC}\times \delta ^{2/5}$$
(24)
in the fast cooling regime, using Eq. 48.

### 2.2 Region B: IC emission

Zone B extends over a distance $`(R\lambda )`$; since we are mainly interested in $`\delta 10^1`$, we approximate $`(R\lambda )R`$ throughout this paper. This region is populated by electrons that have radiatively cooled via S and SSC for a time $`t_b`$. If they were originally injected at the front shock with a power-low energy distribution $`N(\stackrel{~}{\gamma })\stackrel{~}{\gamma }^p`$ then the steady state distribution in region B depends on whether the electrons were adiabatic or radiative in region A. The $`\stackrel{~}{\gamma }_{cB}`$ for electrons that cool in a dynamical timescale $`t_{exp}`$ is in fact

$$\stackrel{~}{\gamma }_{cB}=\frac{m_ec}{\frac{4}{3}\sigma _T\frac{U_{rad}}{2}t_{dyn}},$$
(25)
where $`U_{rad}`$ is the S radiation energy density and the factor $`1/2`$ takes into account that only half of the S photons travel from region A towards region B in the electrons rest frame. The calculation of $`U_{rad}`$ in the two different regimes leads to

$$\stackrel{~}{\gamma }_{cB}\{\begin{array}{cc}(n^2e_b\mathrm{\Gamma }\stackrel{~}{\gamma }^{p1}\stackrel{~}{\gamma }_c^{3p}R^2\delta )^1\hfill & \stackrel{~}{\gamma }_m<\stackrel{~}{\gamma }_c\hfill \\ (n^2e_b\mathrm{\Gamma }\stackrel{~}{\gamma }_m\stackrel{~}{\gamma }_cR^2\delta )^{1.}\hfill & \stackrel{~}{\gamma }_m>\stackrel{~}{\gamma }_c\hfill \end{array}$$
(26)
If electrons in region A are in the slow cooling regime, the new steady state configuration is

$$N(\stackrel{~}{\gamma })\{\begin{array}{cc}\stackrel{~}{\gamma }^p\hfill & \stackrel{~}{\gamma }_m<\stackrel{~}{\gamma }<\stackrel{~}{\gamma }_{cB}\hfill \\ \stackrel{~}{\gamma }^{(p+1)}\hfill & \stackrel{~}{\gamma }_{cB}<\stackrel{~}{\gamma }<\stackrel{~}{\gamma }_c\hfill \\ \stackrel{~}{\gamma }^{(p+2)}\hfill & \stackrel{~}{\gamma }>\stackrel{~}{\gamma }_c\hfill \end{array}$$
(27)
for slow cooling ($`\stackrel{~}{\gamma }_m<\stackrel{~}{\gamma }_{cB}`$) or

$$N(\stackrel{~}{\gamma })\{\begin{array}{cc}\stackrel{~}{\gamma }^2\hfill & \stackrel{~}{\gamma }_{cB}<\stackrel{~}{\gamma }<\stackrel{~}{\gamma }_m\hfill \\ \stackrel{~}{\gamma }^{(p+1)}\hfill & \stackrel{~}{\gamma }_m<\stackrel{~}{\gamma }<\stackrel{~}{\gamma }_c\hfill \\ \stackrel{~}{\gamma }^{(p+2)}\hfill & \stackrel{~}{\gamma }>\stackrel{~}{\gamma }_c\hfill \end{array}$$
(28)
for fast cooling ($`\stackrel{~}{\gamma }_m>\stackrel{~}{\gamma }_{cB}`$). $`\stackrel{~}{\gamma }_{cB}>\stackrel{~}{\gamma }_c`$ is very unlikely, even if the electrons in region A have two sources of energy loss (S and SSC), because $`t_b<t_{exp}`$ and SSC can be inhibited further for KN effect. If electrons in region A undergo fast cooling, the most common case is $`\stackrel{~}{\gamma }_{cB}<\stackrel{~}{\gamma }_c`$ and the new steady state configuration is

$$N(\stackrel{~}{\gamma })\{\begin{array}{cc}\stackrel{~}{\gamma }^2\hfill & \stackrel{~}{\gamma }_{cB}<\stackrel{~}{\gamma }<\stackrel{~}{\gamma }_c\hfill \\ \stackrel{~}{\gamma }^3\hfill & \stackrel{~}{\gamma }_c<\stackrel{~}{\gamma }<\stackrel{~}{\gamma }_m\hfill \\ \stackrel{~}{\gamma }^{(p+2)}\hfill & \stackrel{~}{\gamma }>\stackrel{~}{\gamma }_m\hfill \end{array}$$
(29)
The “external Compton” spectrum in zone B is calculated by means of Eq 19 with $`A=\frac{1}{2}\sigma _TR`$ and the Compton parameter $`\stackrel{~}{Y}_B`$ using Eq 14 and Eq 2729 with $`\stackrel{~}{\tau }_B=\frac{\tau }{2}`$. Even if zone B has a higher number of IC emitting electrons , the Compton parameter can be lower than in region A: a higher fraction of the electrons in region B cool efficiently and they are in average less energetic then in region A, therefore increasing the density parameter $`\stackrel{~}{Y}_B`$ decreases while $`\stackrel{~}{Y}`$ increases because in zone A $`\gamma _c`$ moves towards low values and the KN effect is less and less important (until the electrons become radiative and $`\stackrel{~}{Y}_F`$ does not depend on density). In this situation the SSC radiation can actually account for most of the energy that goes into IC luminosity, (this is the case in Figs 2 and 3 right panels); anyway the SSC usually dominates the spectrum at high energy because it has a higher energy break. The peak flux is

$$\stackrel{~}{F_p}^{IC}(\nu )\stackrel{~}{\tau }_B\stackrel{~}{F_p},$$
(30)
which is a factor $`\delta ^1/2`$ higher than $`\stackrel{~}{F_p}^{IC}`$ in region A; while the spectral breaks are related to the S seed spectrum in slow cooling as follows

$$\stackrel{~}{\nu }_a^{IC}2\stackrel{~}{\gamma }_m^2\stackrel{~}{\nu }_a,$$
(31)
$$\stackrel{~}{\nu }_m^{IC}2\stackrel{~}{\gamma }_m^2\stackrel{~}{\nu }_m,$$
(32)
$$\stackrel{~}{\nu }_c^{IC}2\stackrel{~}{\gamma }_{cB}^2\stackrel{~}{\nu }_c,$$
(33)
if region B is in the slow cooling regime as well, and

$$\stackrel{~}{\nu _a}^{IC}2\stackrel{~}{\gamma }_{cB}^2\stackrel{~}{\nu }_a,$$
(34)
$$\stackrel{~}{\nu }_c^{IC}2\stackrel{~}{\gamma }_m^2\stackrel{~}{\nu }_c,$$
(35)
$$\stackrel{~}{\nu }_m^{IC}2\stackrel{~}{\gamma }_{cB}^2\stackrel{~}{\nu }_m,$$
(36)
if B is in the fast cooling regime. If S seed spectrum is in the fast cooling regime and $`\stackrel{~}{\gamma }_{cB}<\stackrel{~}{\gamma }_c<\stackrel{~}{\gamma }_m`$

$$\stackrel{~}{\nu _a}^{IC}2\stackrel{~}{\gamma }_{cB}^2\stackrel{~}{\nu }_a,$$
(37)
$$\stackrel{~}{\nu }_c^{IC}2\stackrel{~}{\gamma }_{cB}^2\stackrel{~}{\nu }_c,$$
(38)
$$\stackrel{~}{\nu }_m^{IC}2\stackrel{~}{\gamma }_c^2\stackrel{~}{\nu }_m.$$
(39)
## 3 Discussion and conclusions

In Fig 2 and Fig. 3 left panels, we show as examples how the spectrum (respectively for fast and slow cooling) modifies if a short magnetic lengthscale is taken into account. Note that S and IC peak fluxes are lower (see also Eq. 3, Eq. 20, Eq. 30), a lower energy break corresponding to the absorption frequency is present (see Eqs. 7 and 8) and the cooling frequency is at higher energies (see Eq. 5). The SSC component in the slow cooling spectra of region A (Fig 2) is depleted because the KN effect is important. In the right panels we raise the density until $`\stackrel{~}{F}_p\stackrel{~}{n}^{1/2}\delta =F_pn^{1/2}`$, thus $`\stackrel{~}{n}=\frac{n}{\delta ^2}`$. We obtain $`\stackrel{~}{\nu _c}=\nu _c`$ and the $`\delta <1`$ spectrum fails to match the standard spectrum only at radio wavelengths. In this band such high densities, in fact, overcompensate the shorter magnetic lengthscale and $`\stackrel{~}{\nu _a}`$ becomes greater then $`\nu _a`$. Another important feature is that the total (region A plus region B) $`\stackrel{~}{Y}`$ parameter is never much higher then $`Y`$, even when $`Y`$ corresponds to a density $`\delta ^2`$ lower; in fact when $`\stackrel{~}{F}_p=F_p`$, the main contribution to IC luminosity comes from region A and Eq 18 and Eq 16 hold: only $`\stackrel{~}{Y}_S`$ depends on density and $`\delta `$ but the dependence are very weak ($`\stackrel{~}{Y}_S(1+\stackrel{~}{Y}_S)^{(3p)}\stackrel{~}{n}^{(p2)/2}\times \delta ^{(p2)}`$). These illustrative examples of how the model works should now be followed by a modeling of real data. In fact the general agreement between the two spectral shapes suggests that data could also be well reproduced in this framework, with possibly different values for the break frequencies and peak flux respect to the standard model. As regard the corresponding estimated density compared to the standard model, there are two competitive facts: the peak flux is brought to observed ranges by an higher density but the absorption frequency cannot be too high (see fig. 2 and 3) (note that unlike the standard model a high density, $`\stackrel{>}{}10^2`$cm-3 does not produce an extreme inverse Compton component, which is not observed); Therefore, even if there is the possibility that a higher density is required, the suitable value can be derived within the errors only by a broadband modeling of data, with a consistent variations of all the parameters.

Similarly to the standard model there are univocal correlations between the estimated spectral parameters and the derived fireball energy, external density and shock parameters. They are

$$\stackrel{~}{E}\stackrel{~}{F_p}^{3/2}\stackrel{~}{\nu _m}^{5/6}\stackrel{~}{\nu _c}^{1/4}\stackrel{~}{\nu _a}^{5/6}\delta ^{1/2}T^{0.5}$$
(40)
$$\stackrel{~}{n}\stackrel{~}{F_p}^{3/2}\stackrel{~}{\nu _m}^{25/12}\stackrel{~}{\nu _c}^{3/4}\stackrel{~}{\nu _a}^{25/6}\delta ^{1/2}T^{7/2}$$
(41)
$$\stackrel{~}{ϵ_e}\stackrel{~}{F_p}^{1/2}\stackrel{~}{\nu _m}^{11/12}\stackrel{~}{\nu _c}^{1/4}\stackrel{~}{\nu _a}^{5/6}\delta ^{1/2}T^{3/2}$$
(42)
$$\stackrel{~}{ϵ_B}\stackrel{~}{F_p}^{1/2}\stackrel{~}{\nu _m}^{5/4}\stackrel{~}{\nu _c}^{5/4}\stackrel{~}{\nu _a}^{5/2}\delta ^{3/2}T^{5/2},$$
(43)
in slow cooling regime while for a radiative regime

$$\stackrel{~}{E}\stackrel{~}{F_p}^{3/2}\stackrel{~}{\nu _c}^{1/6}\stackrel{~}{\nu _a}^{5/6}\delta ^{1/2}T^{0.5}$$
(44)
$$\stackrel{~}{n}\stackrel{~}{F_p}^{3/2}\stackrel{~}{\nu _c}^{17/6}\stackrel{~}{\nu _a}^{25/6}\delta ^{1/2}T^{7/2}$$
(45)
$$\stackrel{~}{ϵ_e}\stackrel{~}{F_p}^{1/2}\stackrel{~}{\nu _m}^{1/2}\stackrel{~}{\nu _c}^{2/3}\stackrel{~}{\nu _a}^{5/6}\delta ^{1/2}T^{3/2}$$
(46)
$$\stackrel{~}{ϵ_B}\stackrel{~}{F_p}^{1/2}\stackrel{~}{\nu _c}^{5/4}\stackrel{~}{\nu _a}^{5/2}\delta ^{3/2}T^{5/2}.$$
(47)
Some constraints must be added to this set of equations. Because $`ϵ_e+ϵ_B+ϵ_p=1`$, where $`ϵ_p`$ is the fractional energy that goes to protons, we most likely expect $`ϵ_e+ϵ_B<0.5`$. The kinetic energy in the afterglow should not be much lower then the energy released in the $`\gamma `$-ray phase, otherwise very large efficiencies of radiation, incompatible with internal shock models, would be required (Lazzati, Ghisellini & Celotti 1999). Finally the apparent source size, the rate of expansion and the energy in each afterglow should be consistent with observations in the radio wavebands of interstellar scintillation and its quenching (Waxman, Kulkarni & Frail 1998). In principle these requirements together with observations of the spectral parameters at a particular time would set limits on $`\delta `$ and therefore on the structure of the magnetic field generated in shocks. Unfortunately the current data do not allow a model independent estimation of spectral parameters. Actually the spectral shape, and in particular the break frequencies and the peak flux, are not unequivocally constrained by data in even the best studied afterglows (e.g. GRB 970508, GRB 000418). These spectral features are derived by modeling the sparse measurement collected in various wavebands at different times: very often different versions of the standard model can fit the data, consequently finding different spectral parameters (see for example Berger et al 2000, for GRB 000418). Moreover these model dependent quantities (above all $`\stackrel{~}{\nu _a}`$) enter with higher dependences than $`\delta `$ in the estimation of the fireball parameters (Eq. 40-47); therefore the break frequencies and peak flux, derived under the assumption that $`\delta =1`$, cannot obviously be used to figure out the value of $`\delta `$. These equations and the physical limits on the fireball parameters should be used instead to constrain the overall modeling of data. The results would actually lead to very interesting conclusions.

If the modeling succeeded in reproducing the data with $`\delta <1`$, a different value for the external density, energy and shock parameters would possibly be derived. If it indicates a denser environment, this could be more compatible with the strong observational evidence that would place GRBs explosions in a star forming region. However a very short scale magnetic field (say $`\delta \stackrel{<}{}10^4`$) would be probably ruled out, since extremely high densities would be required in order to get fluxes in the observed ranges and the radio flux ($`\nu <\stackrel{~}{\nu _a}`$) would probably be too heavily absorbed: $`\stackrel{~}{F}_\nu \stackrel{~}{n}^{1/2}`$ (in adiabatic) and $`\stackrel{~}{F}_\nu \delta ^1\stackrel{~}{n}^{5/3}`$ (in radiative regime). An other implication of the possibility to fit the very same data with models based on different shock physics assumptions is that parameter estimations greatly suffer from shock theory uncertainties. In this work we have only explored the effect of a short-scale magnetic field, but other reasonable variations of the standard model (such as time dependent $`ϵ_B,ϵ_e`$ and $`p`$) may affect the parameter estimates as well. This would call for alternative model independent external density and energy estimates in order to constrain, together with the afterglow modeling, the shock physics.

On the other hand if this model can reproduce data only with $`\delta =1`$, it could be either evidence in favor of an extremely stable field or the indication that some other modifications of the theory must be taken into account. For example, as already mentioned, $`\delta `$ can vary with time: for instance if the Rayleigh-Taylor instability occurs at the contact discontinuity, the lengthscale would grow as the fireball evolves and consequently late multibands observations would be more consistent with $`\delta =1`$, while early measurements would require $`\delta 1`$. Observational tests of this possibility can be performed using the theory developed in this work; in fact we can look for a coherent changing in time of $`\delta `$, modeling instantaneous spectra at different times, with all the remaining parameters kept constant. Of, course, with a very good dataset, this approach can be extended to all the other shock parameters, including $`ϵ_B,ϵ_e`$ and $`p`$. Such a procedure would yield unbiased estimates of shock parameters variability, and the GRB afterglow would become a laboratory where shock physic can be tested. The best case up to now for constraining $`\delta `$ and its temporal behavior is probably GRB 000926, but in a year SWIFT, REM and and possibly the VLA will provide us with more simultaneous data (also at early times) in X-ray, optical, infrared and longer wavelengths and therefore more constraining fits will be performed. This work actually highlights the importance of observing afterglows at all wavelengths, from radio to hard x-ray, but emphasizes as well that current estimations from this data of the external density surrounding $`\gamma `$-ray bursts may depend on uncertain details of shock physics.

## Acknowledgments

We thank Davide Lazzati for the deep physical insight of his suggestions and Reem Sari for a very useful critic discussion. We thank also Andrei Beloborodov for encouraging discussions. ER thanks the Isaac Newton and PPARC for financial support.


text:
# Statistical mechanics of dissipative systems

## Acknowledgments

The support of the DoE through Caltech’s ASCI Center for the Simulation of the Dynamic Response of Materials is gratefully acknowledged.