text:
# Mass loading of pulsar winds

## 1 Introduction

In this paper we explore the dynamics of magnetized relativistic pulsar winds subject to mass loading. We are motivated by several possible applications: (i) slowing down of the high-$`\sigma `$ flows ($`\sigma `$ is conventionally defined as a ratio of Poynting to particle fluxes) inside the pressure confined pulsar wind nebulae (referred below to as a “$`\sigma `$ paradox”), (ii) low X-ray efficiencies and a conspicuous lack of radio PWN around many energetic pulsars and (iii) the structure of the ram pressure confined pulsar wind nebulae.

(i). The $`\sigma `$ paradox is a long standing problem in pulsar physics (Rees & Gunn 1974, Kennel & Coroniti 1984). Models of the pulsar magnetosphere (Goldreich & Julian 1969, Arons & Scharlemann 1979, Ruderman & Sathurland 1975) predict that near the light cylinder most of the spin-down luminosity of a pulsar should be in a form of Poynting flux, $`\sigma 1`$. On the other hand, modeling of the dynamics of the Crab nebula gives a low value of $`\sigma `$ at what is commonly believed to be a reverse shock - strongly magnetized flows cannot match boundary conditions (Kennel & Coroniti 1984). Several possibilities have been proposed to resolve the $`\sigma `$ paradox. One possibility is that $`\sigma `$ may change between the light cylinder and the reverse shock (Michel 1994, Coroniti 1993, Melatos & Merose 1996, critique by Lyubarsky & Kirk 2001). Alternatively, $`\sigma `$ may change between the reverse shock and the contact discontinuity. Begelman (1998) argued that the post-shock flow may be distorted by kink instabilities which would annihilate the magnetic flux and allow strongly magnetized flows to match the boundary conditions at the edge of the nebula (see also critique by Arons 1998). Arons (2001) has suggested that interaction of the wind with the ejecta filaments may slow down the wind. In this paper we investigate the latter possibility.

(ii). Most PWNs, with the notorious exception of the Crab, have very low efficiencies for the conversion of the pulsar spin-down luminosity $`\dot{E}_0`$ into X-rays, $`L_X10^3\dot{E}_0`$ (Seward & Wang 1988, Becker & Truemper 1997). No simple combination of the pulsar parameters (like period and period derivative) can be fitted to the X-ray luminosity, which suggests that environmental effects may be important (Possenti et al. 2001). Similarly, sensitive radio searches have failed to detect PWNs towards selected energetic pulsars (Gaensler et al. 2000). We show that mass loading of the pulsar wind is very effective in slowing down the wind and thus weakening the termination shock. This may be responsible for low X-ray and radio efficiencies of the PWNs.

(iii) In the past several years a new well defined class of objects - ram pressure confined pulsar wind nebulae - has been established. Apart from a simple model by Wang et al. (1993) no realistic model has been proposed to account for the often unusual structure structure of such PWNs. Bucciantini & Bandiera (2001) have argued that classical bow-shock models cannot be applied to the ram pressure confined PWNs and suggested that mass loading may be important.

The importance of mass loading on the pulsar wind may be estimated using the Crab nebula as an example. It has long been known that that some thermal plasma does evaporate off the trapped ejecta (Wilson 1974, Michel et al. 1991). It is only natural to expect that the coupling of a “light” pair plasma by the “heavy” particles from the ejecta would strongly affects the dynamics of the light fluid. Recent reevaluation of the mass of Crab filaments (Fesen 1997) gives $`4.6\pm 1.8M_{}`$, while the total mass ejected into the Crab nebula by the pulsar over its lifetime is tiny $`10^{12}M_{}`$ (at a rate $`\dot{N}_0=10^{38}`$ pairs/sec (e.g., Kennel & Coroniti 1984); ions may have a two order of magnitude larger mass flux; Galant et. al 1992). Below we show that mass loading starts to strongly affect the dynamics of sub-sonic flows when the pick-up mass flux, $`\dot{M}`$ becomes of the order of the total energy flux $`\dot{E}_0`$ divided by $`c^2`$ (eq. (38)). The total pick-up mass flux is $`\dot{M}M_{\mathrm{picked}}v/r`$ (with velocity given by the expansion velocity $`v2000`$ km/sec of the Crab nebula with the radius $`r2`$ pc; Kennel & Coroniti 1984). The mass flux required to slow down the flow to match the boundary $`M_{\mathrm{picked}}\dot{E}_0r/(c^2v)=3\times 10^{27}g`$, a tiny fraction of the available mass. Thus even a weak coupling of filaments to the flow may strongly affect its evolution. This may be contrasted with the theoretical estimates of Balbus (1981) who have argued that filament evaporation by conduction may occur on a very short time scale, of the order of the dynamical time scale.

In the ram pressure confined PWNs, we should substitute for $`v`$ the velocity of the pulsar $`v200`$ km/sec and for $`r`$ the stand-off distance, which is typically a fraction of a parsec. The total pick-up mass flux, $`M_{\mathrm{picked}}10^{24}`$ g, is nevertheless smaller since typical pulsars have $`\dot{E}_010^{35}`$ erg/sec. This is approximately the amount of the ISM material contained in the head part of the ram pressure confined PWN.

## 2 Model assumptions

We assume that an interaction between the trapped ejecta and the pulsar winds occurs through mass loading of the wind. Filament material may be coupled to the pulsar wind by thermal evaporation, hydrodynamic ablation and by photo- and electro-evaporation by high energy photons and relativistic electrons. Interaction of the ISM clumps with non-relativistic winds have been studied extensively (Borkowski et al. 1990, Cowie & McKee 1977, McKee & Cowie 1977, Hartquist et al. 1986). No treatment of the interaction of the ISM material with relativistic magnetized flow exist; estimates based on the models of suprathermal evaporation (Balbus & McKee 1982) predict a very high rate (Balbus 1981). Here we are interested in the dynamical effects of the mass loading on the relativistic flow; hence we do not consider here the microphysics of the flow-filament interaction and instead we parameterize it by a capture rate, (eq. (3)), taken to be a prescribed function of a distance from a pulsar. The newly acquired particles are assumed to be immediately coupled to the flow. 111Relative streaming of the newly ionized component with respect to the bulk flow may result in ion-cyclotron instabilities and turbulence generation (e.g., Galeev et al. 1996) We neglect the energy losses required for the possible ionization and the energy associated with the thermal motion of donor particles. Mass loading inside a shock may also be important (Zank & Oughton 1991). Mass loading of relativistic shocks is briefly considered in the Appendix A.

The dynamics of the mass loaded flows turns out to be non trivial and somewhat counterintuitive. The overall dynamics of the mass loaded flows resembles a nozzle-type flow where the rate of mass-loading contributes a term which may be considered as a negative pressure (eq. (20)). Hence the effects of mass loading are opposite to that of pressure. For example, contrary to the naive guess, the loading of strongly magnetized pulsar wind leads not to slowing down but to acceleration  of the flow. Acceleration formally proceeds to arbitrary larger velocities which, naturally, cannot be realized in reality. We argue that such acceleration of the wind will result in development of internal instabilities (which we have not addressed yet) that would destroy the toroidal magnetic flux and, in a way similar to suggestion by Begelman (1998), would eventually allow a wind to slow down.

The pulsar wind model that we adopt is similar to the Kennel & Coroniti model: spherical outflow of a polytropic gas with a toroidal magnetic field. In addition, the wind is assumed to be loaded by ejecta particles which are initially at rest. By neglecting the radial component of the magnetic field we limit the applicability of this approach to asymptotically large distances, far from the acceleration region of the flow (Weber & Davis 1967, Goldreich & Julian 1970, Kennel et al. 1983). Typically the launching of the pulsar wind is thought to occur close to the light cylinder, with the flow reaching super fastmagnetosonic velocities by passing successfully through the slow, Alfvén and fast sonic points (Goldreich & Julian 1970). At distances large compared with the light cylinder radius, the flow usually can be well approximated as a radial motion with toroidal magnetic field, launched with some given parameters: mass, energy and magnetic fluxes loss rates.

Pulsar winds are launched at super-fastmagnetosonic velocities and have to decelerate in a shock transition in order to be matched onto the non-relativisticly moving external media (the SNR ejecta or ISM). Under the hydrodynamic approximation, a contact discontinuity separates the shocked wind material from the external medium. In some circumstances, the external material may get inside the contact discontinuity and even inside the reverse shock propagating in the wind. This may be due to either incomplete ISM or ejecta ionization (neutral particles may not be well coupled to the ionized component on time scales of interest to us) or by evaporation of the ISM or ejecta clumps trapped by the expanding wind. Coupling of the newly added material to the flow may considerable change its evolution.

The evolution of non-relativistic mass-loaded flows has been investigated in a number of papers, starting with the Solar wind-comet interaction (Biermann et al. 1967). Lately Hartquist et al. (1986), Arthur et al. (1993), Smith (1996), Williams et al. (1999), Toniazzo (2001) have investigated effects of mass loading on stellar outflows. Here we show that the presence of a magnetic field makes it formally impossible for the mass loaded flow to reach infinity. In order to be able to reach infinity, the flow must either become non-stationary or destroy its magnetic flux, presumably by developing internal instabilities which result in reconnection.

## 3  Governing Equations

The formal treatment of the problem starts with the set of relativistic magnetohydrodynamic equations which can be written in terms of conservation laws (e.g. Landau & Lifshitz 1975):

$`T_{,i}^{ij}=S_i,`$ (1)
$`F_{,i}^{ij}=0,`$ (2)
$`(\rho u^i)_{,i}=R`$ (3)
where

$$T^{ij}=(w+b^2)u^iu^j+(p+\frac{1}{2}b^2)g^{ij}b^ib^j$$
(4)
is the stress-energy tensor, $`w=\rho +\frac{\mathrm{\Gamma }}{\mathrm{\Gamma }1}p`$ is the plasma proper enthalpy, $`\rho `$ is proper plasma density and $`p`$ is pressure, $`b^2=b_ib^j`$ is the plasma proper magnetic energy density times $`4\pi `$, $`p`$ is pressure, $`u^i=(\gamma ,\gamma \beta )`$ are the plasma four-velocity, Lorentz-factor and three-velocity, $`g^{ij}`$ is the metric tensor, $`b_i=\frac{1}{2}\eta _{ijkl}u^jF^{kl}`$ are the four-vector of magnetic field, Levy-Chevita tensor and electro-magnetic field tensor, $`S_i`$ is a four-vector representing energy and momentum sources and $`R`$ is the density source (dimensions of R are $`\mathrm{g}\mathrm{cm}^3\mathrm{sec}^1`$). We assume that the adiabatic index $`\mathrm{\Gamma }`$ is constant for algebraic simplicity. Since we are interested in transonic transitions for relativistic fluid motion, we expect that a relativisticly hot flow with $`p\rho `$ will have an adiabatic index of $`4/3`$ while strongly subsonic (non-relativistic) flows will asymptote to $`\mathrm{\Gamma }=5/3`$.

We assume that loading is due to the medium at rest in the laboratory frame, $`𝐒=0`$, $`S_0=R`$. 222 The filaments are moving with velocities up to 300 km/sec (Trimble 1968), much smaller than the velocity of the pulsar wind near the termination shock $`v_{\mathrm{shock}}c`$ or near the edge of the remnant $`v2000`$ km/sec (Kennel & Coroniti 1984).  Writing out eqns. (1-3) in coordinate form and assuming a stationary, spherically symmetric outflow with toroidal magnetic field, we find

$`{\displaystyle \frac{1}{r^2}}_r\left[r^2(w+b^2)\beta \gamma ^2\right]=R`$ (5)
$`{\displaystyle \frac{1}{r^2}}_r\left[r^2\left((w+b^2)\beta ^2\gamma ^2+(p+b^2/2)\right)\right]{\displaystyle \frac{2p}{r}}=0`$ (6)
$`{\displaystyle \frac{1}{r}}_r\left[rb\beta \gamma \right]=0`$ (7)
$`{\displaystyle \frac{1}{r^2}}_r\left[r^2\rho \beta \gamma \right]=R`$ (8)
The above relations can be simplified if one introduces energy, mass and magnetic flux loss rates,

$`=\beta \gamma ^2\left(b^2+{\displaystyle \frac{\mathrm{\Gamma }}{\mathrm{\Gamma }1}}p+\rho \right)`$ (9)
$`=\beta \gamma \rho `$ (10)
$`𝒦=\beta \gamma b`$ (11)
which evolve according to equations

$`_r={\displaystyle \frac{2}{r}}+R`$ (12)
$`_r={\displaystyle \frac{2}{r}}+R`$ (13)
$`_r𝒦={\displaystyle \frac{2𝒦}{r}}`$ (14)
which have solutions

$`={\displaystyle \frac{\dot{M}_0+\dot{M}}{4\pi r^2}}`$
$`={\displaystyle \frac{\dot{E}_0+\dot{M}}{4\pi r^2}}`$
$`𝒦={\displaystyle \frac{}{2\sqrt{\pi }r}}`$ (15)
where $`\dot{M}_0`$ and $`\dot{E}_0`$ are the central source’s mass and energy loss rates, $``$ is the electromotive force and we have introduced

$$\dot{M}=4\pi Rr^2𝑑r$$
(16)
for the acquired mass flux.

It is convenient to introduce two other parameters of the flow: magnetization $`\sigma `$ and a fast magnetosonic wave phase velocity $`\beta _f`$ and Lorentz factor $`\gamma _f`$

$`\sigma ={\displaystyle \frac{b^2}{w}}={\displaystyle \frac{𝒦^2}{\beta 𝒦^2}}={\displaystyle \frac{^2}{\dot{E}_0\beta ^2}}`$
$`\beta _f^2={\displaystyle \frac{\sigma }{1+\sigma }}+{\displaystyle \frac{\mathrm{\Gamma }p}{(1+\sigma )w}}`$
$`\gamma _f^2{\displaystyle \frac{1}{1\beta _f^2}}={\displaystyle \frac{3(1+\sigma )}{2+\gamma (1+\sigma )}}`$ (17)
In particular for a relativisticly hot plasma, $`\rho p`$, we have $`w=\mathrm{\Gamma }p/(\mathrm{\Gamma }1)`$ and

$$\beta _f^2=\frac{\sigma +\mathrm{\Gamma }1}{1+\sigma }=\frac{1+3\sigma }{3(1+\sigma )}\text{ for }\mathrm{\Gamma }=4/3\text{.}$$
(18)
We will need also an expression for pressure in terms of the fluxes:

$$p=(\mathrm{\Gamma }1)\frac{\gamma (1+\sigma )}{\mathrm{\Gamma }\beta \gamma ^2(1+\sigma )}=\frac{\left(\mathrm{\Gamma }1\right)\left(\gamma \gamma _{f}^{}{}_{}{}^{2}\right)}{\beta \gamma ^2\left(2\mathrm{\Gamma }\right)\mathrm{\Gamma }\gamma _{f}^{}{}_{}{}^{2}}$$
(19)
Since $`p>0`$ it follows that $`>\gamma (1+\sigma )`$ and $`>\gamma \gamma _f^2`$.

Eliminating $`𝒦`$ in favor of $`\gamma _f`$ we get a particularly transparent form for the evolution of Lorentz factor

$$\left(\frac{\left(2\mathrm{\Gamma }\right)}{\beta ^2\gamma ^3\left(\mathrm{\Gamma }1\right)}\right)\left(\gamma ^2\gamma _{f}^{}{}_{}{}^{2}\right)_r\gamma =\frac{2\left(\gamma \gamma _{f}^{}{}_{}{}^{2}\right)}{r}\frac{\left(\gamma 1\right)\left(2\mathrm{\Gamma }\right)\left(1+\gamma \mathrm{\Gamma }\right)\gamma _{f}^{}{}_{}{}^{2}}{\left(\mathrm{\Gamma }1\right)}R$$
(20)
Equation (20) has a form of nozzle-type hydrodynamical flows (e.g., Landau & Lifshitz 1999); in astrophysical context it is best known for Parker’s solutions of the solar wind (Parker 1960). The lhs of eq. (20) contains a familiar special point at the sonic transition $`\gamma =\gamma _f`$. The positively defined first term on the rhs describes the evolution of Lorentz factors due to pressure effects, and the negatively defined second term is due to mass loading. Thus, the rate of mass loading may be considered as a negative pressure.

By neglecting the radial magnetic field far from the acceleration region we have lost the magnetic sling-shot effect often evoked for the acceleration of plasma (Michel 1969, Goldreich & Julian 1970, Kennel et. al 1983). The upshot of those works is that the flow is accelerated to supersonic velocities with a terminal Lorentz factor $`\gamma \sqrt{\sigma }`$. Since it is believed that the plasma near the light cylinder has $`\sigma _010^310^6`$, but the inferred pulsar wind Lorentz factor is $`10^6\sqrt{\sigma _0}`$, an additional acceleration is required.

Without mass loading the rhs of eq. (20) is always positive so that super-fastmagnetosonic flows accelerate while sub-fastmagnetosonic flows decelerate. With mass loading the situation is more complicated. Mass loading enters in the eq. (20) in two ways: explicitly via the rate of mass loading $`R`$, given by the negatively defined second term on the rhs, and implicitly via the accumulated mass and energy fluxes in $``$ and $``$. The negatively defined second term could in principle become larger than the pressure term. Below we show that this indeed happens for magnetized flows.

A quick examination of eq. (20) reveals an important fact: for strongly relativistic flows the mass loading term is enhanced by a factor $`\gamma ^2`$. For typical pulsar winds $`\gamma 10^6`$, thus mass loading is extremely efficient in slowing down the wind. We defer a detail examination of the eq (20) until Section 5 and review beforehand the evolution of unloaded flows.

## 4 Dynamics of relativistic flows without mass loading

To guide us through the effects of mass-loading we first derive the relations governing the evolution of unloaded relativistic magnetized winds. The properties of these solutions have been discussed extensively (e.g. Kennel & Coroniti 1984), but the exact analytical form given below has not been written, according to our knowledge. First, we introduce two new parameters (instead of $`\dot{M}_0`$ and $``$) - the unloaded wind terminal Lorentz factor $`\gamma _0`$ and the initial magnetization parameter $`\sigma _0`$ defined by the following reactions

$`\gamma _0={\displaystyle \frac{\dot{E}_0}{\dot{M}_0(1+\sigma _0)}},`$
$`\sigma _0={\displaystyle \frac{^2}{\gamma _0\beta _0\dot{M}_0}}.`$ (21)
### 4.1 Unmagnetized unloaded relativistic flows

In the absence of a magnetic field the fast magnetosound velocity is

$$\beta _f^2=\left(\mathrm{\Gamma }1\right)\frac{\left(\gamma \right)}{}.$$
(22)
In this case the equation for the evolution of the wind Lorentz factor (20) can be integrated exactly:

$$r\frac{1}{\sqrt{\beta }}(\gamma _0\gamma )^{1/(2(\mathrm{\Gamma }1))}\gamma ^{(2\mathrm{\Gamma })/(2(\mathrm{\Gamma }1))}.$$
(23)
This shows that there are two branches of solutions: supersonic and subsonic. The supersonic branch after initial acceleration with $`\gamma r^{2(\mathrm{\Gamma }1)/(2\mathrm{\Gamma })}`$ ($`\gamma r`$ for $`\mathrm{\Gamma }=4/3`$) reaches a terminal Lorentz factor $`\gamma _0=\dot{E}_0/\dot{M}_0`$ (with $`\gamma _0\gamma r^{2(\mathrm{\Gamma }1)}`$). In the accelerating part the flow has $`p\rho `$ \- it is a pressure driven acceleration. Since pressure falls off with distance faster than density (for $`\mathrm{\Gamma }>1`$), a transition to the coasting phase with $`\gamma \gamma _0`$ occurs at $`p\rho `$. The subsonic branch decelerates to zero velocity at infinity with $`\beta 1/r^2`$, keeping pressure and density almost constant, determined by the mass conservation ratio.

### 4.2  Relativistic Magnetized Unloaded Flows

Magnetized flows also have two branches: subsonic and supersonic. The terminal velocity is determined from the condition $`_r\gamma =0=p`$, which, using eq. (17), can be written

$$\frac{\gamma \gamma _0}{\gamma _0}=\frac{\beta \beta _0}{\beta }\sigma _0.$$
(24)
Eq. (24) generally has two solutions: a supersonic one

$$\gamma =\gamma _0$$
(25)
(it is highly supersonic for $`\gamma _0\sqrt{\sigma _0}`$) and a subsonic one, which in the limit $`\gamma _0\sqrt{\sigma _0}`$ gives

$$\beta =\frac{\sigma _0}{1+\sigma _0}.$$
(26)
No physical solutions exist for $`\beta <\sigma _0/(1+\sigma _0)`$ \- such flows cannot transport enough magnetic flux.

The supersonic solution behaves similarly to the unmagnetized case. The flow is accelerated by pressure effects as long as $`p\rho `$, reaching a coasting phase with $`\gamma \gamma _0`$ when $`p\rho `$. Behavior of the subsonic branch is qualitatively different from the unmagnetized case. To study the behavior of the subsonic branch we can use a simplifying assumption $`\gamma _0\mathrm{}`$ since for subsonic flows the Lorentz factor is usually only weakly relativistic, $`\gamma \sqrt{(1+\sigma _0)/(2\mathrm{\Gamma })}`$ while $`\gamma _0\sqrt{\sigma _0}`$. The limit $`\gamma _0\mathrm{}`$ is equivalent to neglecting the mass loss rate of the central source if compared with the energy loss rate. The evolution of the flow is then given by

$$r\left(\beta \gamma \right)^{(2\mathrm{\Gamma })/(2(\mathrm{\Gamma }1))}\left(\beta \sigma _0/(1+\sigma _0)\right)^{1/(2(\mathrm{\Gamma }1))}=\frac{\beta \gamma }{(\beta \sigma _0/(1+\sigma _0))^{3/2}}\text{for }\mathrm{\Gamma }=4/3.$$
(27)
This shows that subsonic flows reach a minimum velocity given by eq. (26). The fact that magnetized flows cannot slow down to zero velocity since they have to transport magnetic flux is crucial in determining the asymptotic dynamics of the mass loaded flows.

## 5 Evolution of mass loaded pulsar winds

A pulsar produces a wind with properties determined by three parameters: energy flux $`\dot{E}_0`$, terminal Lorentz factor $`\gamma _0`$ and magnetization $`\sigma _0`$. We also assume that near the acceleration region (located presumably close to the light cylinder) the mass loading is unimportant, so that the wind has reached a strongly super-fastmagnetosonic terminal Lorentz factor $`\gamma _0\sqrt{1+\sigma _0}`$ (this requires $`\dot{M}_0^2/\sigma _0,\dot{E}_0\beta _0/(1+\sigma _0)`$) and its internal pressure has dropped to 0.

The ratio of the mass flux due to loading to the initial energy flux is

$$\zeta =\frac{4\pi Rr^2𝑑r}{\dot{E}_0}=\frac{4\pi Rr^3}{(3n)\dot{E}_0}$$
(28)
for a power-law dependence of mass loading on radius, $`Rr^n`$. Parameter $`\zeta `$ will often be used instead of the radius $`r`$.

### 5.1  Weakly loaded flow

Here we show that mass loading of relativistic flows is extremely efficient. Consider evolution of a strongly super-fastmagnetosonic cold flow when mass loading may be approximated as a small perturbation. In the limit $`r0`$ and $`\gamma ,\gamma _01`$, eq. (20) gives

$$\frac{\mathrm{ln}\gamma }{\mathrm{ln}r}=\frac{2(\gamma _0\gamma )}{\gamma +(2\mathrm{\Gamma })\gamma _0/(\mathrm{\Gamma }1)}\frac{4\pi r^2R}{\dot{E}_0}\frac{\mathrm{\Gamma }\gamma ^2\gamma _0(1+\sigma _0)}{\gamma (\mathrm{\Gamma }1)+(2\mathrm{\Gamma })\gamma _0}.$$
(29)
If initially $`\gamma =\gamma _0`$, then

$$\frac{\mathrm{ln}\gamma }{\mathrm{ln}r}=\frac{4\pi r^3R}{\dot{E}_0}\mathrm{\Gamma }\gamma _0^2(1+\sigma _0)=\frac{4\pi r^3R}{\dot{M}_0}\mathrm{\Gamma }\gamma _0$$
(30)
Thus, the flow starts to decelerate quickly when the accumulated mass flux becomes comparable to $`\dot{M}_0/\gamma _0`$ \- an increase in efficiency of loading by a factor of $`\gamma _0`$ ($`10^6`$!). The typical scale for the deceleration of the flow is

$$r_d\left(\frac{\dot{E}_0}{4\pi R\mathrm{\Gamma }\gamma _0^2(1+\sigma ))}\right)^{1/3}=\left(\frac{\dot{M}_0}{4\pi R\mathrm{\Gamma }\gamma _0}\right)^{1/3}$$
(31)
which is $`\gamma _0^{1/3}`$ times smaller than in the non-relativistic case.

Efficient deceleration of supersonic flows by mass loading come from the requirement that the newly acquired particles, which were initially at rest, have to be accelerated to large Lorentz factors. This costs a lot of energy and momentum.

### 5.2 Special points in the flow

Special point in the flow occur when the rhs of the eq. (20) becomes 0. It is shown below that in unmagnetized flows, the only special point in the flow occurs where $`\beta =\beta _f`$, in order to keep $`_r\beta `$ finite, while mass loading of magnetized flows introduces another special point where $`\beta \beta _f`$ and, thus, $`_r\beta =\mathrm{}`$.

The rhs of eq. (20) is equal to 0 only at one point, where the lhs of the eq. (20) is also $`0`$. This is a well studied critical sonic point of the flow where $`\beta =\beta _f`$. For a magnetized flow the rhs of eq. (20) become equal to 0 at two points: a critical sonic point and another point located for a given flow characteristic at radii larger than the sonic point. For even larger radii a sub-fastmagnetosonic mass-loaded flow starts to accelerate \- this behavior is opposite to the unmagnetized flow. Acceleration of the flow proceeds very quickly; the flow reaches a point $`\beta =\beta _f`$ and $`_r\beta =\mathrm{}`$. A steady state flow cannot exist beyond this limiting radius. These conclusion hold for both relativistic and non-relativistic flows (considered in the appendix B. The phase portrait of the mass loaded relativistic magnetized flows is given in Fig. 2, and non-relativistic in Fig. 3.

### 5.3 Sonic point

Sonic point occur when both side rhs of the eq. (20) are equal to zero. It is possible to find the location of the point $`\beta =\beta _f`$ on the phase diagram $`\beta \zeta `$. The general relations are complicated, yet simple approximations may be obtained in the limiting case of small and larger $`\sigma _0`$. In the limit $`\sigma _01`$ the sonic point is located at a Lorentz factor satisfying

$$3\left(2\mathrm{\Gamma }\right)\mathrm{\Gamma }\gamma _{f}^{}{}_{}{}^{3}+\left(4\mathrm{\Gamma }\right)\gamma _{f}^{}{}_{}{}^{2}\left(2+\mathrm{\Gamma }\right)\gamma _f3=0.$$
(32)
The position of the sonic point in radius is given by $`\zeta _f=f(\mathrm{\Gamma })=O(1)`$, where $`f(\mathrm{\Gamma })`$ is a complicated function of the order of unity. For example for relativistic flows, $`\mathrm{\Gamma }=4/3`$, $`\gamma _f=1.09`$, $`\beta _f=0.40`$ and $`\zeta _f=0.91`$ while for non-relativistic flows, $`\mathrm{\Gamma }=5/3`$, $`\gamma _f=1.31`$, $`\beta _f=0.65`$ and $`\zeta _f=0.39`$.

In the opposite limit $`\sigma _01`$, we find

$$\gamma _f^2\frac{4\mathrm{\Gamma }}{2(2\mathrm{\Gamma })}\sigma _0$$
(33)
and $`\zeta _f=\stackrel{~}{f}(\mathrm{\Gamma })/\sigma _0^{5/2}`$ where again $`\stackrel{~}{f}`$ is a complicated function of the order of unity. This gives

$`\gamma _f^2=2\sigma _0,\zeta _f=1.19\times \sigma _0^{5/2}`$ for $`\mathrm{\Gamma }=4/3`$ (34)
$`\gamma _f^2={\displaystyle \frac{7}{2}}\sigma _0,\zeta _f=0.23\times \sigma _0^{5/2}`$ for $`\mathrm{\Gamma }=5/3`$
Thus, low $`\sigma `$ flow experiences a shock transition when the swept up mass flux becomes of the order of the luminosity of the central source, while high $`\sigma `$ flow experience shock transition at a distance $`\sigma _0^{5/6}`$ closer to the source. Location of the sonic points for two choices of $`\sigma _0`$ are shown in Fig. 2.

Since in the hyper-sonic regime $`\gamma _0`$ falls out of the equations, flows with different $`\gamma _0`$ but the same $`\sigma _0`$ experience a shock transition at the same location in $`\zeta `$ coordinate (Fig. 1).

### 5.4  Heavily loaded unmagnetized winds

For unmagnetized flows the effects of mass-loading may completely dominate the dynamics of the flow at larger distances. This occurs when the loaded mass flux becomes larger than the initial mass flux. In this limit $`\dot{M}\dot{E}_0,\dot{M}_0`$ and $`=0`$ it follows from the eq. (20) that only non-relativistic strongly loaded flows can extend to infinity. Asymptotically

$$\beta r^{(5\mathrm{\Gamma }+1)/(\mathrm{\Gamma }+1)}r^{7/2}$$
(35)
for $`\mathrm{\Gamma }=5/3`$.

### 5.5 Evolution of the magnetization parameter $`\sigma `$

Finally we consider evolution of the magnetization parameter $`\sigma `$. Restricting to $`\mathrm{\Gamma }=4/3`$ we find

$$\frac{\gamma ^2\gamma _f^2}{\sigma (1+\sigma )}_r\sigma =\frac{\gamma \gamma _f^2}{r}+\frac{\gamma (\gamma 1+2\beta _f^2\gamma )\gamma _f^2}{2}R$$
(36)
with $`\beta _f`$ and $`\gamma _f`$ given by eq. (17). Thus effects of mass loading are opposite to that of pressure. In the absence of mass loading, pressure effects lead to a $`\sigma `$ decrease for supersonic flows and a $`\sigma `$ increase for subsonic flows. Mass loading contributes to an increase in $`\sigma `$ for supersonic flows and decrease for subsonic flows.

In mass loaded supersonic flows, the pressure term always dominates over the mass loading term in the eq. (36), so $`\sigma `$ always increases, formally diverging at the point $`\beta =\beta _f`$.

## 6 Discussion

In this paper we have investigated the dynamics of relativistic magnetized mass loaded flows. We found that (i) super-relativistic flows are effectively slowed down by mass loading; (ii) subsonic magnetized flows subject to mass loading have to accelerate beyond some radius and formally cannot reach infinity. We suggest that internal instabilities would develop in the accelerating flow, destroying the toroidal magnetic flux. Alternatively, the flow may become non-stationary.

To estimate whether mass loading is important in slowing down pulsar winds we need to know the loading rate $`R`$ which depends on the microphysics of the flow-ejecta coupling. To make a simple estimates we assume that a considerable fraction of the ISM neutrals (or trapped ejecta particles) with density $`(1\xi )n_0`$, ionization fraction $`\xi `$ and mass $`m_p`$ are captured by the flow on a characteristic time scale $`\tau `$: $`Rn_0m_p/\tau `$. Since in the Crab we do see filaments (i.e., they are not evaporated on a dynamical time scale), for static PWNs the age of a pulsar $`T`$ may be considered as a lower limit on $`\tau T`$ and thus an upper estimate of the rate $`R`$. From the eq. (31) it follows that the flow’s initial Lorentz factor decreases on a scale

$$r10^{14}\left(\frac{\dot{E}_0}{10^{36}\mathrm{erg}}\right)^{1/3}\left(\frac{(1\xi )n_0}{0.1\mathrm{cm}^3}\right)^{1/3}\left(\frac{\tau }{10^4\mathrm{yrs}}\right)^{1/3}\left(\frac{\gamma _0}{10^6}\right)^{2/3}\left(\frac{1+\sigma _0}{1001}\right)^{1/3}\mathrm{cm}$$
(37)
We would like to stress that this is a scale for a change of the initial Lorentz factor and not the scale for the mass-loading-induced shock transition, which is much larger (see eq. (38). Distance (37) has a weak dependence on the characteristic time scale $`\tau `$, so that the deceleration of relativistic winds due to pick-up ions may be an important factor affecting the evolution of pulsar winds.

To estimate the position where a mass-loaded shock would occur we assume that pick-up rate is independent of radius. Then the shock transition occurs at

$$r\left(\frac{3\dot{E}_0\tau }{4\pi (1\xi )n_0m_p}\right)^{1/3}\times \begin{array}{cc}1& ,\text{ for }\sigma _01\\ \sigma _0^{5/6}& ,\text{ for }\sigma _01\end{array}10^{17}\times \frac{1}{1+\sigma _0^{5/6}}\mathrm{cm}$$
(38)
(for the same parameters as in eq. (37)). Radius (38) also gives a typical location where a flow becomes dominated by mass loading. Radius (38) is of the order of the shock distance in the Crab nebula. Since it is only a lower bound (since $`\tau T`$) we conclude that a mass loaded shock transition is unlikely to occur in the static PWNs. The shock observed in the Crab nebula nebula is likely an ordinary magnetohydrodynamical shock. On the other hand radius (38) is much smaller than the radius of the Crab nebula. Thus mass loading of the subsonic pulsar wind should be very efficient.

Loading of supersonic flows greatly reduces the strength of the termination shock (the change of the Lorentz factor of the flow), but not the kinetic energy flux through it. This may affect the efficiency of particle acceleration and PWN luminosity both in radio and X-rays. Indeed, Williams et al. (1999) considered the strength of termination shocks in mass loaded isothermal flows and found that shocks become weaker. They have also suggested that weakening of the shocks may contribute to the radio quietness of some wind-blown bubbles. Here we extend this possibility to PWNs. Unfortunately, at this point our understanding of relativistic shock acceleration is not good enough to argue whether weaker shocks (in terms of a change of a Lorentz factor) are less efficient at particle acceleration.

It is unlikely on numerical grounds that loading of supersonic flow will slow it down to the critical point. If it still happens, then, since for a broad variety of the loading profiles the critical point of the flow is a focus and for a focus a smooth transition from the super- to sub-fastmagnetosonic flow is not possible, a flow would shock at a Mach number $`>1`$. For example, Galeev & Khabibarkhanov (1990) have argues that weakly magnetized non-relativistic mass loaded flows shock at a Mach number $`M=2`$. Similar analysis for relativistic flows needs to be done as well.

We have shown that, contrary to naive expectations, loading of the relativistic shocked pulsar wind by particles from the filaments initially would result not in a slowing down of the flow, but in acceleration. To see the dynamic reason for this strange behavior we first note that even if a subsonic flow is weakly magnetized at the source, it becomes strongly magnetized as it reaches its terminal velocity. The total energy flux then consist of the Poynting flux $`P=4\pi \gamma ^2\beta b^2`$ and a particle flux $`(\dot{M}_0+\dot{M})\beta ^2/2`$ (using the nonrelativistic expression for simplicity). The Poynting flux, subject to the requirement to transport magnetic flux, $`\beta \gamma b=/2\sqrt{\pi }r`$, is inversely proportional to the velocity: $`P=^2/\beta `$. Mass loading increases particle energy flux, so that Poynting flux should decrease and velocity should increase.

Acceleration of the wind would result in development of instabilities that would try to destroy the reason of the acceleration: the need to transport the magnetic flux. After the instabilities have developed and the reconnection destroyed magnetic flux the flow will be allowed to decelerate, in line with the idea of Begelman (1998). This deceleration will be very quick, $`\beta r^{7/2}`$, allowing it to match to the boundary of the PWN.

The other possible application of the present work relates to the structure of the ram pressure confined PWNs. The pulsar motion through ISM constantly brings a new neutral material in the wind. The typical scale in this case is the stand-off distance for the bow shock

$$r_b=\sqrt{\frac{\dot{E}_0}{4\pi c\xi n_0m_pV_{NS}^2}}$$
(39)
and a typical time is $`r_b/V_{NS}`$ where $`V_{NS}`$ is the neutron star velocity through ISM. Combining with (38) we find that the ratio of the mass loaded induced shock transition to the stand-off distance

$$\frac{r}{r_b}=\left(\frac{3\xi }{1\xi }\right)^{1/3}\frac{1}{(1+\sigma _0)}\left(\frac{V_{NS}}{c}\right)^{1/3}$$
(40)
which is typically smaller than unity. Thus, mass loading is extremely important for the structure of the ram pressure confined PWNs (see also Bucciantini & Bandiera 2001).

I would like to thank Mikhail Medvedev, Vicky Kaspi, Elena Amato and Steve Balbus for interesting and stimulating discussions.

## Appendix A Mass loading of relativistic magnetized shocks

In this appendix we briefly consider mass loading of relativistic shocks. Assume that at the shock extra material with density $`\mathrm{\Delta }\rho `$ in the rest frame of the shocked plasma has been added to the flow. Then the energy, momentum, mass conservation and induction equations are

$`\left[\left({\displaystyle \frac{\mathrm{\Gamma }}{\mathrm{\Gamma }1}}p+\rho +b^2\right)\beta \gamma ^2\right]=\alpha `$
$`\left[\left({\displaystyle \frac{\mathrm{\Gamma }}{\mathrm{\Gamma }1}}p+\rho +b^2\right)\beta ^2\gamma ^2+p+b^2/2\right]=0`$
$`\left[\rho \beta \gamma \right]=\alpha `$
$`\left[b\beta \gamma \right]=0`$ (A1)
where $`\left[\right]`$ implies a difference between shocked quantities (denoted with a subscript 2) and unshocked quantities (denoted below with a subscript 1). Introducing the ratio of the two four-velocities $`𝒩=\beta _1\gamma _1/\beta _2\gamma _2`$, the ratio of the three-velocities $`=\beta _1/\beta _2`$, upstream magnetization parameter $`\sigma `$ as a ratio of the magnetic energy density to particle energy density $`h_1^2=\sigma \left(\rho _1+\frac{\mathrm{\Gamma }}{\mathrm{\Gamma }1}p_1\right)`$, and normalizing we find

$`\rho _1\{𝒩^2({\displaystyle \frac{1}{}}(\beta _2^2\gamma _2^2(1+\sigma )\sigma /\mathrm{\Gamma }){\displaystyle \frac{(1+\sigma )\mathrm{\Gamma }1}{\mathrm{\Gamma }}}+{\displaystyle \frac{\sigma }{2}})+`$
$`𝒩{\displaystyle \frac{\mathrm{\Gamma }1}{\mathrm{\Gamma }}}++{\displaystyle \frac{\sigma }{2}}\}+`$
$`p_1\left\{1+{\displaystyle \frac{\mathrm{\Gamma }\sigma }{2(\mathrm{\Gamma }1)}}+𝒩^2\left(\sigma {\displaystyle \frac{\mathrm{\Gamma }\sigma }{2(\mathrm{\Gamma }1)}}{\displaystyle \frac{1+\sigma }{}}+{\displaystyle \frac{(1)\mathrm{\Gamma }}{\mathrm{\Gamma }}}(1+\sigma )\beta _2^2\gamma _2^2\right)\right\}+`$
$`{\displaystyle \frac{(\mathrm{\Gamma }1)(\gamma _21)\beta _2^2\gamma _2^2\mathrm{\Gamma }}{\beta _2\gamma _2^2\mathrm{\Gamma }}}\alpha `$ (A2)
The mass loading term is positive definite.

For $`\alpha =0`$ relation (A2) reduces to the known jump conditions for the relativistic magnetized shock of arbitrary strength. For example for a strong ($`𝒩1`$, $`p_1\rho _1`$) magnetized shock this gives

$$\beta _2=\frac{1}{4(1+\sigma )}\left(2(\mathrm{\Gamma }1)+\mathrm{\Gamma }\sigma \pm \sqrt{(2(\mathrm{\Gamma }1)+\mathrm{\Gamma }\sigma )^2+8(2\mathrm{\Gamma })\sigma (1+\sigma )}\right)$$
(A3)
(this gives $`\beta _2=\mathrm{\Gamma }1`$ for $`\sigma =0`$).

Relation (A2) shows that mass loading starts to affect the properties of a strong shock when $`\alpha 𝒩^2\rho _1`$. A full investigation of the relation (A2) is beyond the scope of this paper.

## Appendix B Non-relativistic mass loaded flows

Keeping in mind possible applications of the above results to stellar outflows, below we briefly consider non-relativistic mass loaded magnetized winds. The governing equations

$`{\displaystyle \frac{1}{r^2}}_r\left[r^2\left(\rho v^2/2+{\displaystyle \frac{\mathrm{\Gamma }}{\mathrm{\Gamma }1}}p+h^2\right)v\right]=0`$
$`{\displaystyle \frac{1}{r^2}}_r\left[r^2\left(\rho v^2+h^2/2\right)\right]+_rp=0`$
$`{\displaystyle \frac{1}{r}}_r\left[rhv\right]=0`$
$`{\displaystyle \frac{1}{r^2}}_r\left[r^2\rho v\right]=R,`$ (B1)
may be simplified if one introduces energy, mass and magnetic fluxes,

$$=v\left(h^2+\frac{\mathrm{\Gamma }}{\mathrm{\Gamma }1}p+\rho v^2/2\right),=v\rho ,𝒦=vh$$
(B2)
which obey the equations

$$_r=\frac{2}{r},_r=\frac{2}{r}+R,_r𝒦=\frac{2𝒦}{r}$$
(B3)
with solutions

$$=\frac{\dot{M}_0}{4\pi r^2}+\frac{1}{r^2}Rr^2𝑑r,=\frac{\dot{E}_0}{4\pi r^2},𝒦=\frac{}{2\sqrt{\pi }r}$$
(B4)
Introducing a fast magnetosonic wave phase velocity,

$$v_f^2=\frac{h^2}{\rho }+\frac{\mathrm{\Gamma }p}{\rho },$$
(B5)
and eliminating $`𝒦`$ in favor of $`v_f`$ we get the equation for the evolution of velocity:

$$\left(v^2v_f^2\right)\frac{}{v}_rv=\frac{\mathrm{\Gamma }1}{2\mathrm{\Gamma }}\frac{2(v^2+2v_f^2)}{r}\frac{\mathrm{\Gamma }1}{2}Rv^2$$
(B6)
The first term on the rhs represents the effects of pressure on the velocity; it is always positive since $`p=(\mathrm{\Gamma }1)/(2\mathrm{\Gamma }v(2\mathrm{\Gamma }))(2(v^2+2v_f^2))>0`$. The negatively defined second term is proportional to the rate of mass loading.

When $`R=0`$ the evolution of the magnetized outflows can be integrated analytically. Treating velocity as an independent variable we find from the eq. (B6):

$$rv^{(2\mathrm{\Gamma })/(2(\mathrm{\Gamma }1))}\left(2\dot{E}_0v\dot{M}_0v^32^2\right)^{1/(2(\mathrm{\Gamma }1))}$$
(B7)
Which shows that a terminal velocity at $`r\mathrm{}`$ is determined by the third order equation for $`v`$ which involves the energy, mass and magnetic fluxes at the source.

It is convenient to introduce velocity $`v_0`$ (the terminal velocity of the super-fastmagnetosonic flow) and $`\sigma _0`$ (magnetization parameter) by the following relations

$`\dot{E}_0={\displaystyle \frac{\dot{M}_0v_0^2}{2(1\sigma _0)}}{\displaystyle \frac{\dot{M}_0v_0^2}{2}},`$
$`^2={\displaystyle \frac{\dot{M}_0v_0^3\sigma _0}{2(1\sigma _0)}}{\displaystyle \frac{\dot{M}_0v_0^3\sigma _0}{2}}`$ (B8)
When mass loading is not important the fast sound speed expressed in terms of $`v`$ and $`\sigma _0`$ is

$$v_f^2=\frac{\mathrm{\Gamma }1}{2}\left(v_0^2v^2\right)+\frac{2\mathrm{\Gamma }}{2}\frac{v_0^3\sigma _0}{v}$$
(B9)
For weak magnetization $`\sigma _01`$ the fast sonic flow with $`v=v_f`$ is located at

$$v=\sqrt{(\mathrm{\Gamma }1)/(\mathrm{\Gamma }+1)}v_0=v_0/2\text{ for }\mathrm{\Gamma }=5/3.$$
(B10)
From the eq. (B7) it follows that unmagnetized supersonic flow reaches a terminal velocity $`v_0`$ (and zero pressure), while subsonic flow slows down as $`v1/r^20`$ (and constant pressure). Magnetized flows at $`r=\mathrm{}`$ behave somewhat different: at $`r=\mathrm{}`$ they satisfy

$$v^3vv_0^2+v_0^3\sigma _0=0$$
(B11)
Eq. (B11) has two real solutions for $`\sigma _01`$ corresponding to super- and sub-fastmagnetosonic branches. For a small but finite magnetization parameter, the super-fastmagnetosonic branch remains almost unchanged, reaching $`vv_0`$ at $`r=\mathrm{}`$, while the sub-fastmagnetosonic branch asymptotes to a constant velocity determined by total magnetic and energy fluxes: $`vv_0\sigma _0^2/\dot{E}_0`$ (for $`\sigma _01`$). Thus, magnetized flows cannot slow down to zero velocity since they have to transport magnetic flux.

Both sub- and super-fastmagnetosonic branches at large radii have $`\rho r^2,pr^{2\mathrm{\Gamma }}`$ and $`hr^1`$, so that the plasma $`\beta `$ parameter (inverse of magnetization) decreases with radius $`\beta =2p/h^2r^{2(\mathrm{\Gamma }1)}`$. Both super- and sub-fastmagnetosonic branches expand to zero pressure at infinity. Flows which were weakly magnetized at the source (large $`\beta `$) become strongly magnetized as they expand. In both cases magnetized flow have zero internal pressure at infinity, so that the terminal fast sound speed is equal the shear Alfvén wave velocity:

$$v_A^2=\frac{B^2}{4\pi \rho }=\frac{v_0^3\sigma _0}{2v}\{\begin{array}{cc}\frac{v_0^2\sigma _0}{2}& \text{at the fast branch }vv_0\\ \frac{v_0^2}{2}& \text{at the slow branch }vv_0\sigma _0\end{array}$$
(B12)
Next we consider special points of the flow. The rhs of the eq. (B6) becomes zero at

$$\left(\frac{\mathrm{\Gamma }+1}{2(\mathrm{\Gamma }1)}Rr+\frac{1}{r^2}Rr^2𝑑r\right)\frac{4\pi r^2}{\dot{M}_0}=1+\frac{v_0^2}{v^2}\frac{v_0^3\sigma _0}{v^3}$$
(B13)
where parameterization (B8) was used. This equation has solution if the maximum on the rhs (which is reached at $`v=3v_0\sigma _0/2`$) is larger than 0, which limits magnetization parameter to $`\sigma _0<2/(3\sqrt{3})0.384`$.

The condition $`v=v_f`$, given by

$$(1+\mathrm{\Gamma })\left(1+\frac{4\pi Rr^2𝑑r}{\dot{M}_0}\right)=(\mathrm{\Gamma }1)\frac{v_0^2}{v^2}+(\mathrm{\Gamma }2)\frac{v_0^3\sigma _0}{v^3},$$
(B14)
can be satisfied at a given radius only for one particular value of $`v`$.

To proceed further we assume that mass loading has a power law dependence on radius: $`Rr^n`$. Eqns (B13) and (B14) then become

$`{\displaystyle \frac{(5+\mathrm{\Gamma }n(\mathrm{\Gamma }+1))}{2(\mathrm{\Gamma }1)}}\zeta =1+{\displaystyle \frac{v_0^2}{v^2}}{\displaystyle \frac{v_0^3\sigma _0}{v}}`$
$`(1+\mathrm{\Gamma })(1+\zeta )\left({\displaystyle \frac{v}{v_0}}\right)^3(\mathrm{\Gamma }1){\displaystyle \frac{v}{v_0}}(2\mathrm{\Gamma })\sigma _0=0`$ (B15)
where a ratio of the mass flux due to loading to the initial mass flux

$$\zeta =\frac{4\pi Rr^2𝑑r}{\dot{M}_0}=\frac{4\pi Rr^3}{(3n)\dot{M}_0}$$
(B16)
has been introduced. Below we assume that $`\zeta `$ is an increasing function of radius and that $`n<(5+\mathrm{\Gamma })/(\mathrm{\Gamma }+1)`$.

The system (B15) determines the location $`\zeta `$ and velocity at the fast sonic point. It can be resolve as an implicit function of $`v`$:

$`\sigma _0={\displaystyle \frac{v\left(\left(3n\right)v^2\left(1+\mathrm{\Gamma }\right)^2+\left(1+\mathrm{\Gamma }\right)\left(1+n\left(3n\right)\mathrm{\Gamma }\right)v_{0}^{}{}_{}{}^{2}\right)}{\left(3\left(3\mathrm{\Gamma }\right)\mathrm{\Gamma }n\left(2\mathrm{\Gamma }\right)\left(1+\mathrm{\Gamma }\right)\right)v_{0}^{}{}_{}{}^{3}}}`$
$`\zeta ={\displaystyle \frac{2\left(1+\mathrm{\Gamma }\right)\left(3v^2+v_{0}^{}{}_{}{}^{2}\right)}{v^2\left(3\left(3\mathrm{\Gamma }\right)\mathrm{\Gamma }n\left(2\mathrm{\Gamma }\right)\left(1+\mathrm{\Gamma }\right)\right)}}`$ (B17)
see Fig 4. In the absence of magnetic field, $`\sigma _0=0`$, relations (B17) can be resolved explicitly for the location and velocity at the fast sonic point:

$`\zeta ={\displaystyle \frac{4}{\mathrm{\Gamma }(3n)n1}}`$
$`v=\sqrt{{\displaystyle \frac{(\mathrm{\Gamma }1)(\mathrm{\Gamma }(3n)n1)}{(1+\mathrm{\Gamma })^2(3n)}}}v_0`$ (B18)
For non-vanishing $`\sigma _0`$ the location of the sound point is pushed to smaller radii (and larger $`v`$), reaching $`\zeta =0`$ at $`\sigma _0=2/3\sqrt{3}`$ and $`v=v_0/\sqrt{3}`$ independent of $`\mathrm{\Gamma }`$ and $`n`$. In particular, for the simple case of a homogeneous loading, $`n=0`$, and $`\mathrm{\Gamma }=5/3`$ the fast sonic point is at $`\zeta =v_0^23v^2/(5v^2),\sigma _0=2v/5v_0\left(1+8v^2/v_0^2\right)`$.