###### Abstract

In this paper we review the basics of magneto-rotational properties of neutron
stars focusing on spin-up/spin-down behavior at different evolutionary stages.
The main goal is to provide equations for the spin frequency changes in
various regimes (radio pulsar, propeller, accretor, etc.). Since presently
spin behavior of neutron stars at all stages remains a subject of many
uncertainties, we review different suggestions made over the years in the
literature.

###### keywords:

neutron stars; magnetic field; radio pulsars; accretion)

1 1 0 https://doi.org/ Spin evolution of neutron stars Spin evolution of
neutron stars Pavel Abolmasov 1, Anton Biryukov 2,3 and Sergei B. Popov 2,4*
Pavel Abolmasov, Anton Biryukov and Sergei B. Popov Abolmasov, P.; Biryukov,
A.; Popov, S.B. Correspondence<EMAIL_ADDRESS>Tel.: +39-040-2240-397
(S.B.P.)

##  1 Introduction

In many respects, observational appearances of neutron stars (NSs) are
determined by their spin period and magnetic field. Thus, on the one hand, an
understanding of magneto-rotational evolution is necessary to explain observed
phenomena and to predict properties of yet undiscovered types of sources. On
the other hand, it allows us to derive parameters of NSs (e.g., magnetic
field) from observations and understand the formation and evolution of various
sources.

Magneto-rotational evolution is partly driven by the internal physics of NSs
(e.g., field decay), and partly by the interaction of an NS with the
surrounding medium. Typically, NS activity (particles wind, electromagnetic
radiation, magnetic field) works against a tendency of external matter to
penetrate down to the NS surface either due to gravity or just ram pressure.
This counterbalance of external and internal forces determines the stage of
evolution of an NS. Here we use the following classification (see e.g.,
Lipunov (1992)): ejector, propeller, accretor, and georotator. Transitions
from stage to stage often can be specified in terms of the equality of some
specific radii, see Fig. 1. Let us define the most important of them.

The equality of kinetic energy and the absolute value of potential energy of
the matter surrounding a compact object defines the gravitational capture
radius (aka Bondi radius):

$R_{\mathrm{G}}=\frac{2GM}{v^{2}}.$ (1)

Here $G$ is the Newton constant, $M$ is the compact object mass, and $v$ is
velocity relative to the medium accounting for the sound velocity
$c_{\mathrm{s}}$: $v^{2}=v_{\infty}^{2}+c_{\mathrm{s}}^{2}$ ($v_{\infty}$ is
the NS spatial velocity).

Magnetic field lines corotate with the NS. Thus, the linear velocity of a
field line grows (in the equatorial plane) as $\Omega R$, where $R$ is the
distance from the star’s center and $\Omega=2\pi/P$ is spin frequency. As this
velocity is limited by the speed of light, there is a critical distance called
the light cylinder radius:

$R_{\mathrm{l}}=c/\Omega,$ (2)

where $c$ is the velocity of light.

Plasma frozen in the magnetosphere corotates with it with the maximum velocity
$\Omega R$. If this value is larger than the local Keplerian velocity
$\sqrt{GM/R}$ then a centrifugal barrier prevents accretion down to the NS
surface. Equality of the linear and Keplerian velocity defines the corotation
radius:

$R_{\mathrm{co}}=(GM/\Omega^{2})^{1/3}.$ (3)

Recently, Lyutikov Lyutikov (2023) (see also Abolmasov and Biryukov (2020),
Sec. 2.3) reconsidered the penetration of matter in the cases of spherical and
disc accretion regarding centrifugal force. He derived a value of the
centrifugal barrier, $R_{\mathrm{br}}$, that is smaller than
$R_{\mathrm{co}}$. In particular, for the spherical accretion for an aligned
dipole in the equatorial plane: $R_{\mathrm{cb}}=(2GM/3\Omega^{2})^{1/3}$.

Magnetic pressure (energy density) of a dipolar magnetosphere behaves as
$P_{\mathrm{mag}}\propto R^{-6}$. Thus, it rapidly grows toward the compact
object. The pressure of gravitationally captured free-falling matter grows as
$P_{\mathrm{ext}}\propto R^{-5/2}$ under the standard assumption that the
sound velocity is about the free-fall velocity, i.e. $T\propto r^{-1}$, and
density is $\propto r^{-3/2}$ from the continuity relation. Then at some
distance $R_{\mathrm{m}}$ the matter can be stopped by the magnetic field (for
low fields it can happen that $R_{\mathrm{m}}$ is smaller than the NS radius
$R_{\mathrm{NS}}$, then the influence of the field can be neglected). The
magnetospheric radius $R_{\mathrm{m}}$ might be calculated differently in the
case of disc and spherical accretion. Also, the radius is calculated
differently depending on the relative value of $R_{\mathrm{m}}$ and
$R_{\mathrm{G}}$. For spherical accretion and $R_{\mathrm{m}}<R_{\mathrm{G}}$
we have:

$R_{\mathrm{m}}=\left(\frac{\mu^{2}}{2\dot{M}\sqrt{2GM}}\right)^{2/7}.$ (4)

This critical distance is also called the Alfven radius, $R_{\mathrm{A}}$.
Here $\mu=BR_{\mathrm{NS}}^{3}$ is a magnetic moment, where $B$ is the
equatorial surface dipolar field. As it is an approximate estimate, the
numerical coefficient in the denominator can be different for different
assumptions. The accretion rate can be estimated as:

$\dot{M}=\xi_{\mathrm{acc}}R_{\mathrm{G}}^{2}\rho v,$ (5)

where $\rho=m_{\mathrm{p}}\,n$ is the density of the accreted medium at
$R_{\mathrm{G}}$ ($m_{\mathrm{p}}$ is proton mass). The coefficient
$\xi_{\mathrm{acc}}$ depends on details of accretion (for example, eq. 4 is
obtained for $\xi_{\mathrm{acc}}=4\pi$), and below we discuss it in every
particular case. Note, that $\dot{M}$ can appear in equations even if
accretion is not possible at that stage; then this value just characterizes
the properties of the external medium.

Shvartsman radius $R_{\mathrm{Sh}}$ is determined by the balance between the
pressure of relativistic particles wind and external medium. For
$R_{\mathrm{Sh}}>R_{\mathrm{G}}$ and standard magneto-dipole rate of losses we
have:

$R_{\mathrm{Sh}}=\left(\frac{2\mu^{2}(GM)^{2}\Omega^{4}}{3\dot{M}v^{5}c^{4}}\right)^{1/2}.$
(6)

However, the coefficient in the numerator can differ for other models of
rotation energy losses (see Sec. 3 below).

$R_{\mathrm{Sh}}$ must be larger than the light cylinder radius as it is
assumed that the wind is generated in the vicinity of the light cylinder
boundary. If $R_{\mathrm{Sh}}<R_{\mathrm{G}}$ then the boundary between the
wind and external plasma is rapidly shrinking as the wind pressure
$P_{\mathrm{w}}\propto R^{-2}$ and the pressure of the gravitationally
captured matter grows as $P_{\mathrm{ext}}\propto R^{-5/2}$. Thus, at the
ejector stage, we have
$R_{\mathrm{Sh}}>\max{(R_{\mathrm{G}},R_{\mathrm{l}})}$.

Figure 1:  A schematic representation of relation between critical radii at
different stages of magneto-rotational evolution.

The ejector stage is over when the relativistic particle wind is not produced
anymore and the Poynting flux becomes sufficiently low. Usually, it is assumed
that this happens when external plasma can penetrate the light cylinder. This
can happen due to spin-down, magnetic field decay, or modifications of
external conditions. Inside the light cylinder, matter starts to interact with
the magnetosphere. Unless some very exotic situation is realized, immediately
after the transition from the ejector stage, $R_{\mathrm{m}}$ is just slightly
smaller than $R_{\mathrm{l}}$. This means that
$R_{\mathrm{m}}>R_{\mathrm{co}}$ as the Keplerian velocity cannot be larger
than $c$. This corresponds to the propeller regime. Finally, when
$R_{\mathrm{m}}<R_{\mathrm{co}}$ and $R_{\mathrm{m}}<R_{\mathrm{G}}$ matter
can reach the surface of the NS. So, we have accretion. Some other more exotic
situations can be realized. We discuss them in Sec. 6.

Figure 2:  A schematic representation of spin evolution of an isolated NS. The
compact object consequently passes the stages of ejection, propeller, and
accretor. Finally, it can reach the stage when the spin behavior is determined
by the accretion of turbulent angular momentum.

The basic picture of magneto-rotational evolution looks like this (see Fig.
2). An NS is born rapidly rotating ($P\lesssim 1$ s) with a sufficient
magnetic field ($B\sim 10^{12}$ G) at the ejector stage (often as a radio
pulsar). Then it slows down (and, maybe, the magnetic field decays) coming to
the propeller stage at some $t=t_{\mathrm{E}}$. The transition happens at a
critical period $P_{\mathrm{E}}$. Consequent spin-down can bring the NS to the
stage of accretion, especially in the case of an isolated object, when another
critical period, $P_{\mathrm{A}}$, is reached at $t=t_{\mathrm{P}}$ This route
can be followed by an isolated NS as well as by an NS in a binary system. At
the accretion stage, the spin period can change due to changes in the angular
momentum of the incoming matter. For example, in the case of an isolated NS,
this can be due to the interstellar turbulence (see Sec. 5.2). Turbulence
becomes important at the period $P_{\mathrm{cr}}$. At the stage of turbulent
accretion, the spin period fluctuates around the value $P_{\mathrm{turb}}$.

This simple picture can be perplexed by an initial stage of fallback, non-
trivial field evolution, variations of properties of the surrounding medium,
etc. Mainly, in this review we focus on the spin evolution. Magnetic field
evolution was recently reviewed, e.g. in Igoshev et al. (2021).

In the following sections, we discuss the main evolutionary stages and
transitions between them in different astrophysical scenarios. We describe the
spin evolution of an NS in ‘chronological’ order, starting from the early
stage of fallback of a part of supernova ejecta onto a newborn NS (Sec. 2).
Then we describe the ejector stage paying special attention to radio pulsars
(sec. 3). This stage is followed by the propeller (Sec. 4) and accretor stages
(Sec. 5). In Sec. 6 we discuss some hypothetical regimes of interaction of an
NS with surrounding medium. Our conclusions are presented in Sec. 7.

In the text we use convention that $X_{\mathrm{x}}=X/10^{x}$.

##  2 Fallback

In this section, we discuss the so-called fallback – an episode of very
intense accretion of the ejecta onto a newborn compact object. The total
amount of accreted matter can reach $\lesssim 0.1M_{\odot}$ and the accretion
rate at early stages can be from a $\sim$few up to $\sim 10^{4}M_{\odot}$
yr-1. For the first time this process was discussed more than half a century
ago in Colgate (1971); Zel’dovich et al. (1972). Modern approaches to modeling
fallback are considerably influenced by the seminal paper by Roger Chevalier
Chevalier (1989).

In Chevalier (1989) the author analyzed in detail different regimes of
fallback accounting for heating due to radioactive decay and changes in
opacity (in time and at different distances from the compact object).
Asymmetries due to the kick velocity, rotation, and pulsar emission were also
briefly discussed. The author used mainly a semianalytical approach so that
all involved physical processes could be clearly visible and understood.
Strong initial fallback is due to a reverse shock and interaction of the
ejecta with the envelope. At the initial stages, the accretion rate dependence
on time roughly follows the law $\dot{M}\propto t^{-3/2}$. Later, when the
ballistic regime is established, the standard law of fallback $\dot{M}\propto
t^{-5/3}$ is realized. Duration of a quasi-spherical fallback depends on the
properties of the progenitor and for red supergiants can reach $\sim 1$ year.

Fallback properties at early stages significantly depend on the structure of
the pre-SN. Thus, the total amount of the accreted matter also depends on the
parameters of the progenitor. Emission properties are dependent on the
trapping of photons by the envelope. Only at a low accretion rate, do they
start to diffuse out. Before this moment, the effects of the central compact
object are hidden.

Fallback can be important also to explain some properties of gamma-ray bursts
and superluminous supernovae, see e.g. Metzger et al. (2018) and references
therein. However, here we will not discuss these issues focusing on
consequences related to the magneto-rotational evolution of NSs as the
fallback can drastically modify the initial properties of these compact
objects. For example, it can reduce the surface magnetic field at the early
phases of evolution, can significantly change the spin period, and finally,
can even change the evolutionary stage of a young NS (e.g., preventing the
appearance of a radio pulsar). We discuss all these possibilities below.

Significant fallback (here the accretion rate $\dot{M}$ is more important than
the total accreted mass $\Delta M$, however, these two parameters are
interlinked as $\Delta M=\int\dot{M}(t){\rm d}t$) can strongly influence the
observational appearance of a young NS. In particular, the radio pulsar
activity can be completely switched off. This possibility was suggested and
analyzed in Muslimov and Page (1995). If the pressure of the external medium
falling onto the compact object is larger than both pressure of the
relativistic wind and magnetic field pressure for all radii larger than
$R_{\mathrm{NS}}$ then the field is submerged at depth $\sim$ a few tens or
hundred meters below the surface. Thus, the external field that determines the
radio pulsar activity and many other manifestations of an NS is orders of
magnitude lower than the initial value. For a range of parameters used in
Muslimov and Page (1995) it was found that the field diffuses out on a time
scale $\sim 10^{3}$ yrs. This value mostly depends on the depth of submergence
which is determined by the amount of accreted matter. Later this process was
studied in several papers in more realistic frameworks.

In Geppert et al. (1999) the authors modeled the field diffusion in more
detail. In particular, it was demonstrated that due to the field submergence,
the characteristic ages of young pulsars, defined as $P/2\dot{P}$, can be
significantly different (larger) in comparison with their actual ages. The
authors notice that due to compression the field in the crust is amplified by
$\sim 2$ orders of magnitude. And as the spatial scale for this field is
relatively small it is a subject of relatively rapid Ohmic decay. With the
parameters used in this paper, the field diffuses out in $\sim 10^{3}-10^{4}$
yrs. Also, the authors notice that submergence can explain the absence of
radio pulsars in some supernova remnants (SNRs).

A new interest in the fallback field submergence appeared in the 2010s. This
was related to the necessity to explain the properties of central compact
objects in SN remnants (CCOs). Wynn Ho was probably the first who used the
field submergence to explain the observational appearance of CCOs Ho (2011).
The author shows that accretion of $\Delta M\lesssim 10^{-4}M_{\odot}$ is
enough for significant submergence with the time of re-emergence $\lesssim
10^{4}$ yrs. This is sufficient to explain the properties of the most studied
CCOs.

A 2D model of submergence and re-emergence coupled with detailed calculations
of field decay and thermal evolution was presented in Viganò and Pons (2012).
The authors obtained results qualitatively in correspondence with the study Ho
(2011). Corresponding tracks in the $P-\dot{P}$ diagram can be found in Pons
et al. (2012). In 2D it was possible to analyze changes in the topology of the
field during re-emergence. In more detail, this process was studied in Igoshev
et al. (2016). The authors presented a 2D model of field re-emergence after a
fallback episode with $\Delta M=10^{-3}\,M_{\odot}$. It is shown that higher
multipoles, which are strongly suppressed during fallback, re-emerge faster
than the dipole. Thus, they can dominate at ages $\sim 10^{4}-10^{5}$yrs. This
is important for the pulsar switch-on.

Details of field submergence were calculated numerically in several papers by
Bernal et al.: Bernal et al. (2010), Bernal et al. (2013), Fraija et al.
(2018). In Bernal et al. (2013) the authors considered the interaction of the
fallback matter with a magnetic field loop anchored in the crust. They
demonstrated with 2D and 3D modeling that for sufficiently strong accretion on
a timescale $\sim 100$ ms after the reverse shock starts to interact with the
loop, the field (which is always below the accretion shock) is submerged.
Convection in the accretion envelope is crucially important for the fate of
the field. If the envelope remains convective then the field is not submerged
completely. If the reverse shock is not developed at all the field is not
modified significantly.

Field submergence can modify the evolutionary status of a newborn NS. In a 1D
approach (i.e., for spherical symmetry) it was analyzed in Shigeyama and
Kashiyama (2018) and Zhong et al. (2021). The authors aim to account for the
back reaction of the NS (due to a relativistic particle wind and magnetic
pressure) on the fallback flow. This allows them to explain different types of
young neutron stars (radio pulsars, CCOs, magnetars, etc.) by different
outputs of such interaction. In Shigeyama and Kashiyama (2018) the authors
presented self-similar analytical considerations and in Zhong et al. (2021) –
numerical modeling.

In Zhong et al. (2021) the model is mostly determined by two parameters:
characteristic fallback time, $t_{\mathrm{fb}}$, and the fallback mass,
$M_{\mathrm{fb}}$. The latter parameter can be calculated from the dependence
of the fallback rate on time (notice, that it is not the rate of fallback on
the NS surface as accretion can be prevented by relativistic wind or the
magnetosphere):

$\dot{M}_{\mathrm{fb}}=\dot{M}_{\mathrm{fb,ini}}\times\begin{cases}1&t\leq
t_{\mathrm{fb}}\\\ (t/t_{\mathrm{fb}})^{-l}&t>t_{\mathrm{fb}}\end{cases}.$ (7)

The fallback time was set to
$t_{\mathrm{fb,1}}=({t_{\mathrm{fb}}}/{10\,\mathrm{s}})=1$. Among other
parameters, it determines the initial encounter radius $R_{\mathrm{enc}}$ at
which the relativistic wind from the NS collides with the fallback matter:

$R_{\mathrm{enc}}\approx
R_{\mathrm{fb}}=(GMt_{\mathrm{fb}}^{2})^{1/3}=2.7\times
10^{9}\mathrm{cm}\,t_{\mathrm{fb,1}}^{2/3}.$ (8)

The luminosity of the wind, $L_{\mathrm{sd}}$, is calculated either with the
usual magneto-dipole formula, or it can be enhanced due to accretion of some
fraction of the fallback matter through a disc as it was suggested by Parfrey
et al. (2016).

The authors determine three possible regimes: significant field submergence,
disturbed magnetosphere, and matter expulsion. If
$({R_{\mathrm{enc}}L_{\mathrm{sd}}})/({GM\dot{M}_{\mathrm{fb,ini}}})<1$ then
the fallback matter reaches the surface and submerges the magnetic field or
significantly disturbs the magnetosphere (in the latter case the authors
speculate about a magnetar formation due to field enhancement).

Matter expulsion can happen due to the usual wind produced by the magnetic
dipole. This corresponds to $M_{\mathrm{fb}}<7.7\times
10^{-8}M_{\odot}\,B_{13}^{2}P_{-2}^{-4}t_{\mathrm{fb,1}}^{5/3}$. If the
fallback mass is larger but $<5.2\times
10^{-3}M_{\odot}\,B_{13}^{2}P_{-2}^{-14/3}t_{\mathrm{fb,1}}^{23/9}$ then the
fallback matter is repelled by the enhanced spin-down power.

Fallback can result in the formation of a disc around the NS if the
surrounding matter has enough angular momentum. These can be viscous accretion
discs or passive discs like the one observed around the anomalous X-ray pulsar
4U 0142+61 Wang et al. (2006). In this century, such discs have been
extensively studied by many authors, starting with the paper by Chatterjee et
al. Chatterjee et al. (2000), especially as an alternative to explain the
properties of anomalous X-ray pulsars. See, e.g. Jones (2007); Benli and Ertan
(2016) and references therein. Here we do not pretend to present an extensive
review of fallback disc studies. Instead, we focus on a few papers mostly
relevant to the subject of the present paper.

Ideally, fallback properties might be obtained from realistic calculations of
SN explosions. At the moment, numerical models of SN are far from being
complete. Still, new important results provide some optimism. In particular,
recent calculations presented in Janka et al. (2022) provide a realistic setup
for analysis of the fallback influence on magneto-rotational evolution of NSs.

In this model, the newborn NS is not fixed in the center of the grid. The
movement of the compact object allows accounting for new effects. In
particular, the authors demonstrate that an NS can be spun up by the fallback
and the obtained spin is correlated with the kick velocity. Vorticity in the
fallback matter is not related to the spin of the progenitor. Instead, it is
generated by non-radial hydrodynamic instabilities. Late fallback is
associated with the highest angular momentum. Accretion of just $\sim
10^{-2}M_{\odot}$ can result in the spin period $\sim 0.05$ s which is even
slightly shorter than a usual value for young NSs Popov and Turolla (2012).
Typically in this scenario, the angle between the spin axis and velocity
vector is $\lesssim 30^{\circ}$. Alignment is more pronounced for large kicks
as in this case, the NS can move farther from the center of the explosion, so
accretion is less symmetric. For low kicks, an NS accretes from all directions
and this results in averaging of the accreted angular momentum. According to
the proposed model, fallback discs might be widely spread among young NSs with
not very low spatial velocities.

Joint effects of fallback disc formation and field emergence on the appearance
of radio pulsars are considered in Fu and Li (2013). The authors perform
population synthesis modeling accounting for spin up/down due to the NS
interaction with the disc. In their model, the magnetic field is not submerged
by the initial strong fallback. Instead, the authors use the equation
identical to the one used for field behavior in an accreting binary system:

$B=\frac{B_{0}}{1+\Delta M/10^{-5}M_{\odot}}.$ (9)

Here $B_{0}$ is the field before accretion starts. In the population synthesis
it is assumed to have a log-normal distribution with
$\langle\log\,B_{0}/\mathrm{G}\rangle=12.35$ and
$\sigma_{\mathrm{logB_{0}}}=0.4$. And $\Delta M$ is the total amount of the
accreted matter. After accretion stops, the field re-emerges with the rate:

$\dot{B}=0.01\left(\frac{\Delta
M}{M_{\odot}}\right)^{-3}\left(1-\frac{B}{B_{0}}\right)^{2}\,\rm{G\,yr}^{-1}.$
(10)

An NS can pass stages of accretion, propeller, and finally ejector. At the
ejector stage, the NS can appear as a radio pulsar. The authors of Fu and Li
(2013), in the first place, intended to explain statistics of the braking
indices distribution. By construction, in their model braking indices might be
$\leq 3$. Thus, the results are more or less compatible with the data on
pulsars in SNRs. Still, it is very important that the authors used in their
population model complicated early spin evolution due to NS-disc interaction
and sub/re-emergence of the magnetic field.

Discovery of several long-period radio sources attracted new interest to the
model of fallback discs as the interaction between an NS and the disc can
result in a drastic spin-down of the compact object. The first new source with
a peculiar period was the 76-second pulsar PSR J0901-4046 Caleb et al. (2022).
In addition, two radio sources with periods $\sim 1000$ s were found: GLEAM-X
J1627-52 Hurley-Walker et al. (2022) and GPM J1839-10 Hurley-Walker et al.
(2023). As their properties point towards magnetar-scale magnetic field, all
three objects might be young. So, it is impossible to explain their spin
periods simply by a magneto-dipole-like spin-down. Ronchi et al. Ronchi et al.
(2022) successfully explain long spin periods of freshly discovered sources in
the framework of a fallback disc around a young strongly magnetized NS.

As specific angular momentum of the fallback matter is $\sim
10^{16}-10^{17}$cm2 s-1, see Janka et al. (2022), the circularization radius
where the accreted matter settles at a circular orbit forming a disc is
$R_{\mathrm{c}}\sim 10^{6}-10^{8}$ cm. The viscous time scale in the disc is
$t_{\nu}\approx 2.1\times 10^{3}T_{\mathrm{c,6}}^{-1}R_{\mathrm{c,8}}^{1/2}$
s. Here $T_{\mathrm{c,6}}$ is the disc’s central temperature in units $10^{6}$
K and $R_{\mathrm{c,8}}=R_{\mathrm{c}}/10^{8}$ cm. The mass flow rate in the
outer parts of the disc is:

$\dot{M}_{\mathrm{d}}(t)=\dot{M}_{\mathrm{d,0}}\left(1+\frac{t}{t_{\nu}}\right)^{-\alpha}.$
(11)

The coefficient $\alpha$ depends on the opacity and in Ronchi et al. (2022)
the authors assume it to be equal to 1.2. The accretion rate at the inner edge
of the disc is limited by the Eddington rate.

The authors do not account for the field submergence, but allow for the decay
as they are interested in times scales up to $\gtrsim 10^{5}$ yrs. The field
decay is calculated according to the analytical expression from Aguilera et
al. (2008):

$B(t)=B_{0}\frac{\exp{(-t/\tau_{\mathrm{Ohm}})}}{1+\frac{\tau_{\mathrm{Ohm}}}{\tau_{\mathrm{Hall,0}}}(1-\exp{(-t/\tau_{\mathrm{Ohm}})})}.$
(12)

Here $\tau_{\mathrm{Ohm}}$ is the Ohmic dissipation time scale which is
assumed to be $4.4\times 10^{6}$ yrs. And $\tau_{\mathrm{Hall,0}}=6.4\times
10^{4}(10^{14}\rm{G}/B_{0})$ is the initial Hall time scale which depends on
the initial magnetic field $B_{0}$.

With such a setup and realistic fallback rates, it is possible to obtain spin
periods $\sim 100$ s for $B_{0}\sim 10^{13}$ G within $10^{5}$ yrs and $\sim
1000$ s for $B_{0}\sim 10^{14}-10^{15}$ G.

##  3 Ejector stage and radio pulsar activity

With more than 3000 known rotation-powered pulsars111See the ATNF Pulsar
Catalogue Manchester et al.
(2005),https://www.atnf.csiro.au/research/pulsar/psrcat/expert.html, v.1.73,
ejectors are the most common type of Galactic NSs observed. The rotational
evolution of these sources is not affected by the surrounding interstellar
medium but is driven by the complex electrodynamical processes within their
magnetospheres.

For now, it is fairly clear, that most of the rotational energy of the
ejectors is carried away by the wind of energetic charged particles (electrons
and positrons) accelerated by strong electric fields within the magnetosphere
of the star. At the same time, the physics of this process is still not well
understood (there are detailed reviews on this topic Beskin et al. (2015);
Beskin (2018)). Typical observed ejectors have spin periods in the range
$P\sim 10^{-3}-10$ s while their surface magnetic fields are of the order
$\sim 10^{9}-10^{13}$ G. This is consistent with observed spin-down rates of
ejectors $\dot{P}\sim 10^{-20}-10^{-9}$ s/s (e.g. Lorimer and Kramer (2012)).
Therefore, their spin-down luminocities
$L_{\mathrm{sd}}=4\pi^{2}I\dot{P}P^{-3}$ are of order $10^{30}-10^{38}$ erg/s
(see Fig. 3).

Figure 3: Spin periods and their derivatives for 3473 pulsating neutron stars
collected in the ATNF Pulsar Catalogue. Middle panel: Classical single
rotation-powered pulsars are indicated by grey dots; red circles mark sources
with high-energy pulsed emission; green crosses indicate radio quietness,
while horizontal blue stripes are for pulsars in binary systems. Lines of
constant surface magnetic fields $B_{\mathrm{sd}}\propto\sqrt{P\dot{P}}$ Gs,
ages $\tau_{\mathrm{sd}}=P/2\dot{P}$ and luminosities
$L_{\mathrm{sd}}\propto\dot{P}P^{-3}$ are also drawn assuming
$R_{\mathrm{NS}}=10$ km and $I=10^{45}$ g cm2. Bottom panel: Sources with
negative $\dot{P}$ are shown in the same coordinates. Most of these pulsars
are in globular clusters and their “negative spin-down” is due to the
Shklovskii effect Shklovskii (1970). The pulsar graveyard border is plotted
according to Beskin and Istomin (2022). Upper panel: Vertical lines show the
spin period values for pulsars for which $\dot{P}$ has not yet been estimated.
Grey lines are for all such objects, green ones mark the values for radio-
quiet sources, red ones are for those exhibiting high-energy emission, while
blue lines are for pulsars in binary systems.

In classical theory, single neutron stars are viewed as rotating, highly
magnetized, conducting bodies with a primarily dipolar magnetic field. Indeed,
higher-order multipoles decrease very rapidly with the distance from the star.
However, multipolar magnetic fields have been also considered in the
literature (see e.g. Krolik (1991); Kantor and Tsygan (2004); Pétri (2015,
2017)).

The magnetosphere of the star is created by its magnetic field, electric
fields, and electric currents. If there is an electric field
$\mathbf{E_{\perp}}$ perpendicular to the local magnetic field $\mathbf{B}$, a
non-zero Poynting flux $\mathbf{S}\propto(\mathbf{E}\times\mathbf{B})$ exists,
which carries away the rotational energy of the star. This is true for both
vacuum and plasma-filled rotating magnetospheres, which are two fundamental
models for the physics of the ejector.

In the vacuum approach, the spin-down of a neutron star is due to the surface
electric currents induced by the corresponding electric fields. So the
deceleration of the star is driven by the Ampère’s force. On the other hand,
the environment of real neutron stars is not vacuum, but rather filled with
plasma Sturrock (1971); Ruderman and Sutherland (1975). In general, the
magnetosphere of a neutron star consists of two parts, containing closed and
open magnetic field lines respectively. The plasma can freely escape through
the region of open lines, giving rise to the so-called ‘pulsar wind’. Its
back-reaction will also be responsible for the loss of the rotational energy
Spitkovsky (2004). Notice, that magnetospheric currents define local electric
and magnetic fields in a self-consistent manner.

So, finally, the spin-down torque acting on the neutron star can be described
as the sum of two components: the vacuum and the wind. Hereafter we assume
that the spin-down torque is positive and that the neutron star is spherical
with a moment of inertia I. This means that for the ejector stage:

$I{\mathbf{\dot{\Omega}}}=-\mathbf{K}_{\mathrm{sd,vac}}-\mathbf{K}_{\mathrm{sd,w}}.$
(13)

For the vacuum component, the surface quadrupolar electric field $E\sim(\Omega
R/c)B$ is induced by the rotation, and so the corresponding spin-down torque
is

$K_{\mathrm{sd,vac}}=f_{\perp}\dfrac{\mu_{\perp}^{2}}{R^{3}_{\mathrm{l}}},$
(14)

where $\mu_{\perp}=\mu\sin\chi$ is the orthogonal component of the magnetic
moment, $\chi$ is the angle between the magnetic and spin axes, while
$f_{\perp}$ is a structure factor. The latter is of the order of unity to a
first approximation, although its details depend on certain assumptions. The
vacuum losses are inefficient for a perfectly aligned rotator:
$K_{\mathrm{sd,vac}}=0$ when $\chi=0$. This is similar to the behavior of the
emission from a rotating magnetic dipole. Furthermore, if $f_{\perp}=2/3$,
then the equation (14) coincides with the classical one for the
electromagnetic loss of a point-sized dipole. However, it is important to
point out that there is no real magneto-dipolar emission from a neutron star
(see the discussion in Section 2.3 of Beskin (2018)). So the above equation is
nothing more than a law for the spin-down of a magnetized, conducting,
rotating sphere in a vacuum.

As for the wind component of the ejector’s spin down, its power
$L_{\mathrm{sd,w}}$ is based on the electric potential drop $\Delta\varphi$
near the star’s surface, so that

$L_{\mathrm{sd,w}}=\mathbf{K}_{\mathrm{sd,w}}\cdot\mathbf{\Omega}=-2\pi
r_{p}^{2}c\kappa\rho_{\mathrm{GJ}}\Delta\varphi,$ (15)

where $r_{p}=R_{\mathrm{NS}}\sqrt{\Omega R_{\mathrm{NS}}/c}$ is the star’s
polar cap radius, $\rho_{\mathrm{GJ}}\approx\omega B_{\parallel}/c$ is the
Goldreich-Julian charge density Goldreich and Julian (1969), while $\kappa\sim
10^{3}$ is the coefficient describing the primary particle density Kou and
Tong (2015). Taking the maximum potential drop for a rotating dipole
$\Delta\Phi=\mu_{\parallel}/R_{\mathrm{l}}^{2}$ Ruderman and Sutherland (1975)
as an essential scale factor, the wind braking torque can be formally
expressed in a manner similar to (14) as

$K_{\mathrm{sd,w}}=f_{\parallel}\dfrac{\mu^{2}_{\parallel}}{R^{3}_{\mathrm{l}}},$
(16)

where $\mu_{\parallel}=\mu\cos\chi$ and
$f_{\parallel}=2\kappa\Delta\varphi/\Delta\Phi$ is another structure factor
that depends ultimately on the microphysics of the acceleration gap. The wind
torque calculated as shown above disappears for a perfectly orthogonal rotator
($\chi=90^{\circ}$ ) but is present for a perfectly aligned one.

At the turn of the millennium Harding et al. Harding et al. (1999) and then Xu
& Qiao Xu and Qiao (2001) suggested that in fact, both spin-down mechanisms
have a right to exist and that the total spin-down is driven by their
combination. However, as the pulsar spins down, the potential drop may not
sustain the need to accelerate particles Ruderman and Sutherland (1975). Thus,
the maximum period $P_{\mathrm{death}}$ exists where wind losses are
efficient. This is the so-called ‘death line’ (or ‘death period’) of
radiopulsars, beyond which these objects cannot be observed as radio-loud
sources (since plasma acceleration and secondary particle production are
thought to be related to pulsar radio emission). After crossing the death
line, the wind component become negligible, so the pulsar’s spin-down is
expected to follow the vacuum approximation Contopoulos and Spitkovsky (2006).

In this way, a complete torque that is responsible for the loss of rotational
energy of a single ejector can be written as follows:

$K_{\mathrm{sd}}=\left\\{\begin{array}[]{lr}K_{\mathrm{sd,vac}}+K_{\mathrm{sd,w}}=\dfrac{\mu^{2}}{R_{\mathrm{l}}^{3}}(f_{\perp}\sin^{2}\chi+f_{\parallel}\cos^{2}\chi),&\text{if
}P<P_{\mathrm{death}}\\\
K_{\mathrm{sd,vac}}=f_{\perp}\dfrac{\mu^{2}}{R_{\mathrm{l}}^{3}}\sin^{2}\chi,&\text{if
}P\geq P_{\mathrm{death}}\end{array}\right.$ (17)

The values of the these factors and their dependence on the parameters of the
star are of particular interest for an accurate understanding of the rotation
history of ejectors. Both quantities depend on the details of the physics of
neutron star magnetospheres – the topology of the field and currents, the
interaction of the currents with the stellar surface, and the physical
structure and composition of the surface itself. And, in principle, there are
two approaches to unravelling these factors theoretically: analytical and
numerical ones. The purely analytical approach works well when considering a
neutron star with a vacuum magnetosphere (and hence in determining
$f_{\perp}$, see the section 3.2), whereas the complex structure of the filled
magnetosphere (taking into account both $f_{\perp}$ and $f_{\parallel}$) can
only be revealed numerically as described in the section 3.3. But first, a
brief review of the pulsar ‘death condition’ must be given.

###  3.1 Pulsar ‘death line’

There has been much discussion regarding the ‘death period’
$P_{\mathrm{death}}$ since the discovery of active radio pulsars. It is widely
believed that pulsar radio emission is generated by the secondary plasma near
the polar regions of neutron stars Sturrock (1971); Ruderman and Sutherland
(1975). The ‘death’ of pulsars is therefore closely linked to the possibility
of the formation of electron-positron pairs. As shown by Ruderman and
Sutherland Ruderman and Sutherland (1975), the critical condition for this
process is the equality of the height of the vacuum gap $H\propto
B^{-4/7}R_{\mathrm{l}}^{3/7}\rho_{\mathrm{c}}^{2/7}$ (where
$\rho_{\mathrm{c}}$ is the radius of curvature of the local magnetic fields)
to the polar cap radius $r_{\mathrm{p}}$. Corresponding ‘death’ condition can
be written as

$P_{\mathrm{death}}=a\dot{P}_{-15}^{b}$ (18)

where $a\approx 0.6$ s, $b=4/11$ and $\dot{P}_{-15}=\dot{P}/(10^{-15}\mbox{
s/s})$. This equation represents a monotonous curve (in particular, a straight
line) on the $P-\dot{P}$ diagram for neutron stars when its plotted in
logarithmic scales. This is why the pulsar ‘death condition’ is often referred
to as the ‘death line’.

However, equation (18) only partially explains the properties of the observed
population of radio pulsars since there are many sources with
$P>P_{\mathrm{death}}$ if the ‘death period’ is calculated using eq. (18).
Then Chen and Ruderman Chen and Ruderman (1993) considered several variants of
the magnetic field structure near the polar cap and obtained refined
parameters for the death condition. In particular, they found $(a,b)=$(3.2 s,
4/9) and (11 s, 1/2) as possible solutions. Accounting for relativistic frame
dragging affects the condition of the electron-positron pairs production due
to additional induced electric fields Muslimov and Tsygan (1992). Thus, Zhang
et al. Zhang et al. (2000) have considered two polar cap acceleration models
and found $P_{\mathrm{death}}\approx 0.67\dot{P}^{2/5}_{-15}$ s for a
realistic death line in the case of dipolar and $P_{\mathrm{death}}\approx
5.75\dot{P}^{1/2}\rho_{c,6}^{-1/2}$ s for multipolar field configurations,
respectively.

In the influential paper by Faucher-Giguère and Kaspi Faucher-Giguère and
Kaspi (2006) describing the population synthesis of isolated radiopulsars, the
death line condition was used in the form

$\dfrac{B}{P^{2}}=0.17\times 10^{12}\mbox{ G/s}^{2}$ (19)

which corresponds to $P_{\mathrm{death}}\approx 3.28\dot{P}^{1/3}_{-15}$ s.
This relation was earlier established empirically by Rawley et al. Rawley et
al. (1986) from observations of a newly discovered long-period pulsar. In
spite of its purely empirical nature, it was successful in reproducing the
observed population of active pulsars. However, many authors (like Gullón et
al. Gullón et al. (2014) and Graber et al. Graber et al. (2023)) have not
assumed the existence of a death line at all, since pulsars become
undetectable earlier due to small luminosities and beaming fractions.

At the same time, Beskin et al. recently presented a more sophisticated
analysis of the processes near the polar cap of active pulsars Beskin and
Litvinov (2022); Beskin and Istomin (2022) and found that the classical
equation (18) with $b=4/11$ is still a good approximation to the death line
condition but either $a=8.27$ or $a=4.6$ s depending on the spin-down model.
In fact, they have considered both a numerical MHD model based on Spitkovsky’s
work Spitkovsky (2006) and their own analytical model Beskin et al. (1993).
They argued, that the analytical model (with $a=8.27$ s) fits the observations
much better since it allows many long-period pulsars to exist. However, the
pulsars with the longest known periods (J0250+5854 with $P\approx 23.5$ s,
J0901-4046 with $P\approx 76$ s) are still contradict this model. So the
physics of pulsar ‘death’ is far from complete.

In addition, it is worth mentioning that the ‘death line’ is likely to depend
not only on its magnetic field (or spin-down rate), but also on the magnetic
angle (see e.g. Beskin (2018)). This fact is often ignored in the analysis. A
reasonable form for the ‘death’ condition can then be approximately expressed
as

$\cos\chi<AP^{15/7}B_{12}^{-8/7},$ (20)

where $A\approx 0.4-0.5$ from the analysis of more than 150 pulsars with known
$\chi$ Novoselov et al. (2020) (see also Kantor and Tsygan (2004)),
$B_{12}=B/(10^{12}\ G)$ – is the equatorial magnetic field of a pulsar and
period $P$ is in seconds. This means that nearly orthogonal pulsars shut down
much earlier than aligned ones, which affects their statistics.

###  3.2 Spin-down of a star with vacuum magnetosphere

The braking torque resulting from the non-zero Poynting flux of a magnetized
star rotating in a vacuum can be calculated analytically and in great detail.

As early as November 1967 (a few months before the official discovery of
radiopulsars by Hewish and Bell Hewish et al. (1968)) Pacini published a note
briefly discussing the energy budget of a rotating neutron star with a dipolar
magnetic field Pacini (1967). He suggested that such a star could lose its
rotational energy with luminosity
$L_{\mathrm{sd}}\sim\mu^{2}\Omega/R^{3}_{\mathrm{l}}$. But more than a decade
earlier, Deutsch had derived Deutsch (1955) the first vacuum solution for the
electromagnetic field of a perfectly conducting, rigidly rotating spherical
star.222See also an interesting pedagogical note by Satherley and Gordon
Satherley and Gordon (2022) regarding Deutsch’s classic work. In general, the
braking torque vector $\mathbf{K}_{\mathrm{sd,vac}}$ acting on such a star
consists of three orthogonal terms. If $\mathbf{e}_{\Omega}$ and
$\mathbf{e}_{\mu}$ are unit vectors collinear with the spin and magnetic axes
of the star, respectively, then

$\dfrac{\mathbf{K}_{\mathrm{sd,vac}}}{K_{\mathrm{sd,0}}}=\eta_{1}(r_{\mathrm{NS}})\cos\chi(\mathbf{e}_{\mu}\times\mathbf{e}_{\omega})+\eta_{2}(r_{\mathrm{NS}})\cos\chi(\mathbf{e}_{\mu}\times\mathbf{e}_{\omega})\times\mathbf{e}_{\omega}+\eta_{3}(r_{\mathrm{NS}})\sin^{2}\chi\mathbf{e}_{\omega}.$
(21)

Here $K_{\mathrm{sd,0}}~{}=~{}\mu^{2}/R^{3}_{\mathrm{l}}$, while
$\eta_{1,2,3}(r_{\mathrm{NS}})$ are dimensionless functions, depending on the
assumptions made about the structure of the fields inside and outside the star
Ostriker and Gunn (1969); Davis and Goldstein (1970); Good and Ng (1985);
Melatos (2000). These functions depend on the small parameter
$r_{\mathrm{NS}}~{}=~{}R_{\mathrm{NS}}/R_{l}$.

In terms of the spin-down equation (17) one then has
$f_{\perp}~{}\equiv~{}-\eta_{3}$. As the third term of (21), the second term
represents the braking torque arising due to radiation in the far (wave) zone.
It is responsible for the evolution of the magnetic angle $\chi$. In the
vacuum solution $\eta_{2}~{}>~{}0$, therefore the magnetic angle decreases
systematically ($\dot{\chi}~{}<~{}0$) on the timescale of the pulsar spin-
down. Moreover, typically $\eta_{3}~{}=~{}-\eta_{2}$.

The Euler equations describing the evolution of a spherical star with moment
of inertia $I$, are therefore as follows:

$P\dot{P}=4\pi^{2}\eta_{3}(r_{\mathrm{NS}})\dfrac{\mu^{2}}{Ic^{3}}\sin^{2}\chi$
(22)

for the evolution of the spin period and

$P^{-2}\dot{\chi}=-4\pi^{2}\eta_{2}(r_{\mathrm{NS}})\dfrac{\mu^{2}}{Ic^{3}}\sin\chi\cos\chi$
(23)

for the evolution of the magnetic angle. As for the first term in (21),
proportional to $\eta_{1}(r_{\mathrm{NS}})$, it does not directly affect the
observed spin evolution of the star. Since it is always perpendicular to both
the rotational and magnetic axes, it causes a precession around the magnetic
axis with angular velocity

$w_{1}=\eta_{1}(r_{\mathrm{NS}})\dfrac{\mu^{2}}{R_{\mathrm{l}}^{2}}\dfrac{\cos\chi}{Ic},$
(24)

or, in terms of period, $T_{1}~{}\approx~{}6\times
10^{4}PB^{-2}_{12}~{}\cos^{-1}\chi$ years. Here we assumed $I=10^{45}$ g cm2,
$R_{\mathrm{NS}}=10^{6}$ cm, and the spin period $P$ is taken in seconds Good
and Ng (1985); Biryukov et al. (2012). The value of $T_{1}$ is so short
relative to the spin-down time scale at the ejector stage

$t_{\mathrm{ej}}\sim\dfrac{I\Omega
R_{\mathrm{l}}^{3}}{\eta_{3}(r_{\mathrm{NS}})\mu^{2}}\sim 3\times
10^{7}P^{2}B_{12}^{-2}\mbox{ years,}$ (25)

because $\eta_{3}(r_{\mathrm{NS}})\sim 1$ to a first approximation, while
$\eta_{1}\propto 1/r_{\mathrm{NS}}$.

Indeed, the first term in (21) can be interpreted either as a result of the
inertia of the magnetic dipole itself Goldreich (1970), or as a near-field
radiation torque Melatos (2000), or as the appearance of an internal field
momentum Beskin (2018). This torque is often referred to as the ‘anomalous
braking torque’ (see e.g. Beskin and Zheltoukhov (2014)). Melatos Melatos
(2000) has given an accurate calculation of $\eta_{1,2,3}$, assuming the
internal magnetic field to be a point dipole field:

$\eta_{1}(r_{\mathrm{NS}})=\dfrac{6(r_{\mathrm{NS}}^{2}+6)}{5r_{\mathrm{NS}}(r_{\mathrm{NS}}^{6}-3r_{\mathrm{NS}}^{4}+36)}+\dfrac{6-4r_{\mathrm{NS}}^{2}}{15r_{\mathrm{NS}}(r_{\mathrm{NS}}^{2}+1)}\approx\dfrac{3}{5}r_{\mathrm{NS}}^{-1}-\dfrac{38}{17}r_{\mathrm{NS}}+...$
(26)

and

$\eta_{2}(r_{\mathrm{NS}})=-\eta_{3}(r_{\mathrm{NS}})=\dfrac{2r_{\mathrm{NS}}^{4}}{5(r_{\mathrm{NS}}^{6}-3r_{\mathrm{NS}}^{4}+36)}+\dfrac{2}{3(r_{\mathrm{NS}}^{2}+1)}\approx\dfrac{2}{3}-\dfrac{2}{3}r_{\mathrm{NS}}^{2}+...$
(27)

respectively. Thus it is clear that if only leading order terms are taken the
account, then $\eta_{2,3}~{}=~{}2/3$ is constant but
$\eta_{1}~{}\propto~{}r_{\mathrm{NS}}^{-1}$. As a result, unlike the
‘conventional torques’ which are proportional to $1/R_{\mathrm{l}}^{3}$, the
‘anomalous’ one scales as $\propto 1/R_{\mathrm{l}}^{2}$, making it up to four
orders of magnitude larger, since $r_{\mathrm{NS}}\sim 10^{-4}-10^{-2}$ for
most real pulsars. Extensive analysis made in Melatos (2000) has shown that
this torque could be of high importance for a non-spherical neutron star,
causing significant quasi-periodic variations in the spin-down rate and
magnetic angle of a star with a typical timescale $T_{1}$. The first term
($\propto r_{\mathrm{NS}}^{-1}$) of the equation (26) has also been considered
by Beskin & Zheltoukhov in Beskin and Zheltoukhov (2014). They have shown,
that if the angular momentum of the internal electric field is also taken into
account, then the anomalous torque depends on the structure of the magnetic
field inside the star. Thus, if there is a uniform magnetic field up to radius
$R_{\mathrm{in}}<R_{\mathrm{NS}}$, while the dipole configuration is relevant
for $R>R_{\mathrm{in}}$, then

$\eta_{1}(r_{\mathrm{NS}})\approx\left(\dfrac{8}{15}-\dfrac{1}{5}\dfrac{R_{\mathrm{NS}}}{R_{\mathrm{in}}}\right)r_{\mathrm{NS}}^{-1}.$
(28)

It is interesting, that $R_{\mathrm{in}}=0$ leads to the divergence of the
$\eta_{1}$ and, therefore the divergence of the anomalous braking torque.
However, this situation seems to be unphysical and would never be realized.

###  3.3 Spin-down of a star with plasma-filled magnetosphere

Much work has been done in the last two decades to obtain self-consistent
solutions of pulsar magnetospheres in both axisymmetric Contopoulos et al.
(1999); Gruzinov (2005); McKinney (2006); Timokhin (2006); Parfrey et al.
(2012) and oblique Spitkovsky (2006); Kalapotharakos and Contopoulos (2009);
Tchekhovskoy et al. (2013); Pétri (2012); Philippov et al. (2014, 2015)
configurations. Significant progress has therefore been made in this area.
Most models have assumed the force-free approximation, but radiative
magnetospheres have also been considered (e.g. Pétri (2020)). It has been
shown that a magnetized spherical rotator with a plasma-filled magnetosphere
is indeed affected by a torque proportional to $\mu^{2}/R^{3}_{\mathrm{l}}$.
However, the corresponding dimensionless coefficient is slightly different
from that found in the vacuum approximation.

In his seminal paper Spitkovsky obtained a general numerical solution for the
spin-down torque of an oblique rotator Spitkovsky (2006):

$K_{\mathrm{sd}}=\dfrac{\mu^{2}}{R^{3}_{\mathrm{l}}}(k_{0}+k_{1}\sin^{2}\chi),$
(29)

where $k_{0}\sim k_{1}\sim 1$. Later Philippov et al Philippov et al. (2014)
refined the parameters of the spin-down law:

$k_{0}\approx 1\mbox{ and }k_{1}\approx 1.4-3.1r_{\mathrm{NS}}^{2}.$ (30)

That is, in terms of (17) we get $f_{\perp}=k_{0}+k_{1}$ and
$f_{\parallel}=k_{0}$. Since for most real pulsars $r_{\mathrm{NS}}\sim
10^{-4}-10^{-2}$ is small, a good approximation is to keep $k_{0}$ and $k_{1}$
constant. Further PIC modelling Philippov et al. (2015) also produced to a
result that is consistent with the results of MHD simulations. Namely, it was
found that $k_{0}=1.0\pm 0.1$ and $k_{1}=1.1\pm 0.1$. If the radiative
magnetosphere is taken into account then $k_{0}=1.36-1.42$ and
$k_{1}=1.73-1.76$ Pétri (2020).

The solution by Philippov et al. also contains the component responsible for
the magnetic angle evolution. It is with the second term of (21) (or with
scalar equation (23)) assuming $\eta_{2}\equiv 1$. So, numerical consideration
of pulsar magnetospheres predicts magnetic alignment of radio pulsars on spin-
down timescales.

This result is the most complete and accurate theoretical description to date
of the spin evolution of a single neutron star at the ejection stage. If one
does not care about the value of the magnetic angle and its evolution, then
spin-down losses can be estimated quantitatively simply by assuming an
isotropic distribution of magnetic angles such that
$\langle\sin^{2}\chi\rangle=2/3$. Thus, assuming (30), one finally obtains

$\langle L_{\mathrm{sd}}\rangle\approx\dfrac{\xi\mu^{2}\Omega^{4}}{c^{3}},$
(31)

where $\xi\approx 1.93$. For actual numerical analysis, $\xi=2$ remains a good
approximation.

###  3.4 Observational verification

In general, one has to accept that the theory of pulsar spin-down described
above allows one to understand and predict the properties of the actual
ejectors in much detail. On the other hand, the large number of neutron stars
observed as active radio pulsars and the high precision of measurements of
their timing parameters ($P$, $\dot{P}$, $\ddot{P}$, etc.) suggest that this
theory can indeed be well verified in observations. But the reality is much
more complicated. Despite many indirect arguments, there is no direct
observational evidence to support the theory.

In this subsection, without pretending to provide a complete review, we
briefly describe three relevant issues that remain unresolved and deserve to
be discussed.

####  3.4.1 Magnetic angle evolution

The evolution of the magnetic angle $\chi$ during the spin-down of an isolated
neutron star is an important part of its evolution. In a more general
approach, the Euler equations (22) and (23) can be rewritten in terms of the
so-called symmetric and anti-symmetric components of the magnetospheric
currents normalized to the Goldreich-Julian current:
$i_{\mathrm{s}}=j_{\mathrm{s}}/j_{\mathrm{GJ}}$ and
$i_{\mathrm{a}}=j_{\mathrm{a}}/j_{\mathrm{GJ}}$ respectively (see e.g. Beskin
(2018)). In this approach one gets Beskin et al. (1993):

$P\dot{P}\propto i_{\mathrm{s}}+(i_{\mathrm{a}}-i_{\mathrm{s}})\sin^{2}\chi$
(32)

and

$P^{-2}\dot{\chi}\propto(i_{\mathrm{a}}-i_{\mathrm{s}})\sin\chi\cos\chi.$ (33)

It is clear from these equations, that the magnetic angle tends to evolve in
such a way that the loss of rotational energy decreases with time. This means
that if a perfectly aligned rotator ($\chi=0$) will spin down less quickly
relative to the perfectly orthogonal one ($\chi=90^{\circ}$), then one should
expect a magnetic alignment ($\dot{\chi}<0$). On the other hand, if the
physics of pulsar spin-down suggests that an orthogonal rotator loses
rotational energy less efficiently, then magnetic orthogonalisation should
occur ($\dot{\chi}>0$).

As shown in Section 3.2, spin down is suppressed in vacuum magnetospheres when
$\chi=0$, and hence $\dot{\chi}<0$ there. The same result has also been
obtained in numerical simulations of plasma-filled magnetospheres. However, in
the 1980s a theory was developed by Beskin, Gurevich and Istomin Beskin et al.
(1984, 1993); Beskin and Nokhrina (2007); Beskin et al. (2013) which assumes a
non-vacuum magnetosphere of a neutron star, but predicts magnetic counter-
alignment. The MHD-based theory presented in the paper by Philippov Philippov
et al. (2014) and the analytical description by Beskin et al. (BGI theory)
give clear predictions about the evolution of the magnetic angle from
fundamental physics. To our knowledge, these are the only two such
considerations published to date. And their predictions are not consistent
with each other. On the other hand, there is still no observational evidence
in favor of one of any of these approaches, although many analyses have been
carried out Rankin (1990); Tauris and Manchester (1998); Faucher-Giguère and
Kaspi (2006); Weltevrede and Johnston (2008); Young et al. (2010); Lyne et al.
(2013); Arzamasskiy et al. (2015); Novoselov et al. (2020). So this part of
the rotational evolution of an isolated NS is not entirely clear.

####  3.4.2 Pulsar timing irregularities and braking indices

The spin-down equation (17) assumes that an isolated ejector loses its
rotational energy monotonically. However, the actual observed spin evolution
of such objects is not like this. At least about 6 per cent of known radio
pulsars exhibit sudden and rapid spin-up events, known as glitches Espinoza et
al. (2011); Haskell and Melatos (2015); Zhou et al. (2022). However, the spin-
down rate ${\rm d}P/{\rm d}t$ of almost every pulsar appears to be
significantly contaminated by additional irregular, almost stochastic
variations on short time scales of months to years, typically referred to as
‘timing noise’ Boynton et al. (1972); Cordes and Helfand (1980); D’Alessandro
et al. (1995); Urama et al. (2006); Hobbs et al. (2010); Nice et al. (2013).

The strength of this unmodelled component of the pulsar spin-down can be
characterized numerically in various ways. The standard approach is to use
either the second derivative of the pulsar spin period Arzoumanian et al.
(1994) or the variance of the timing residuals (e.g. Cordes and Helfand
(1980)). These two measures are correlated. There is also a widely used
dimensionless combination:

$n_{\mathrm{br}}=\dfrac{\ddot{\Omega}\Omega}{\dot{\Omega}^{2}}=2-\dfrac{\ddot{P}P}{\dot{P}^{2}}$
(34)

– the so-called ‘braking index’. It has a simple and clear physical meaning.
In particular, if the loss of rotational energy follows the expression
$L_{\mathrm{sd}}=-K\Omega^{n+1}$ with $K=const$, then the combination (34) is
nothing other than $n$.

Since the spin-down law for a spherical star is
$I\dot{\Omega}=-f(\chi)\mu^{2}\Omega^{3}/c^{3}$, by direct differentiation it
is clear that

$n_{\mathrm{br}}=3-2\left[2\dfrac{1}{\mu}\dfrac{{\rm d}\mu}{{\rm
d}t}-\dfrac{1}{I}\dfrac{{\rm d}I}{{\rm d}t}+\dfrac{1}{f}\dfrac{{\rm d}f}{{\rm
d}\chi}\dot{\chi}\right]\tau_{\mathrm{ch}},$ (35)

where the parameters of a neutron star are assumed to be variable and

$\tau_{\mathrm{ch}}=-\dfrac{\Omega}{2\dot{\Omega}}=\dfrac{P}{2\dot{P}}$ (36)

characteristic age of a pulsar. Classical vacuum losses ($f_{\perp}=2/3$,
$f_{\parallel}=0$) with constant magnetic fields and moment of inertia give
$n_{\mathrm{br}}=3+2\cot^{2}\chi\lesssim 10^{3}$ for real pulsars. Only the
third term in the brackets is non-zero in this case due to $\dot{\chi}<0$. The
theory by Beskin et al. Beskin et al. (1984) gives a similar result with
$n_{\mathrm{br}}\approx 1.93+1.5\tan^{2}\chi$. On the other hand, Spitkovsky’s
model Spitkovsky (2006) predicts a narrower interval for this quantity under
the same assumptions (see also Ekşi et al. (2016)). Namely,

$n_{\mathrm{br}}=3+2\dfrac{\sin^{2}\chi\cos^{2}\chi}{(1+\sin^{2}\chi)^{2}}\in
3-3.25.$ (37)

But a big problem is that the estimated values of $n_{\mathrm{br}}$ for
hundreds of pulsars are surprisingly far from all these predictions. The
values have been found to range from $\sim-10^{6}$ to $\sim 10^{6}$, and are
negative for about half of the objects (e.g. Hobbs et al. (2004); Biryukov et
al. (2012); Zhang and Xie (2012)). They are therefore unlikely to represent
the secular spin evolution of radio pulsars. In addition, the braking indices
of most pulsars do not seem to be stable from observation to
observation333However, there are a few dozen ‘low noise’ pulsars whose
$n_{\mathrm{br}}$ are constant within spans of 10-30 years, but are still
extremely anomalous Biryukov et al. (2007). There are only about a dozen
sources that are accepted as having meaningful and stable $n_{\mathrm{br}}\in
0.03$-$3.15$ Archibald et al. (2016); Marshall et al. (2016).

Although the physics of pulsar timing irregularities remains generally
unclear, the proposed solutions to the ‘anomalous braking index’ problem can
be qualitatively divided into two categories. On the one hand, there are
theoretical models that assume the relatively slow variability of the NS
parameters affecting equation (35). The authors of these ideas assume either
the variability of the magnetic moment Pons et al. (2012); Zhang and Xie
(2012); Ou et al. (2016); Igoshev and Popov (2020) or the magnetic angle
Melatos (2000); Lyne et al. (2013); Arzamasskiy et al. (2015), or the
effective moment of inertia Janssen and Stappers (2006); Tsang and
Gourgouliatos (2013); Melatos and Link (2014); Hamil et al. (2015, 2016) of
the star. On the other hand, some researchers suggest the existence of an
additional either quasi-periodic or purely stochastic external component in
the spin-down torque $\dot{\Omega}_{\mathrm{ext}}$. The underlying physics has
been proposed to be related either to magnetospheric perturbations Cheng
(1987); Kramer et al. (2006); Contopoulos (2007); Lyne et al. (2010), or the
imprint of the anomalous vacuum torque discussed above Biryukov et al. (2007);
Barsukov and Tsygan (2010); Biryukov et al. (2012).

Biryukov et al. Biryukov et al. (2012) attempted to estimate
$\dot{\Omega}_{\mathrm{ext}}$ from the statistics of the observed second
derivatives of pulsar spin frequencies. Assuming that anomalous braking
indices are due to long-term variations in pulsar spin-down (with typical
timescales of hundreds and thousands of years), they found
$0.5<\varepsilon<0.8$, where
$\varepsilon=|\dot{\Omega}_{\mathrm{ext}}/\dot{\Omega}|$ – is the relative
strength of the effect.

On the other hand, it has been shown in many papers that the observable
parameters of Galactic pulsars can be well reproduced within a population
synthesis adopting (17) with various forms of $f_{\perp}$ and $f_{\parallel}$
Faucher-Giguère and Kaspi (2006); Ridley and Lorimer (2010); Gullón et al.
(2014) but neglecting the corrections due to $\varepsilon\neq 0$. Moreover,
the model-independent analysis of pulsar kinetics in the $P-\dot{P}$ diagram
even led to a reasonable value of their birth rate of $\sim 1$ per century,
simply assuming $f_{\perp}=f_{\parallel}\equiv 1$ and $\varepsilon\equiv 0$
(Keane and Kramer, 2008; Vranešević and Melrose, 2011). This ultimately means
that the simple equation for the isolated ejector spin-down torque
$K_{\mathrm{sd}}=\xi\mu^{2}\Omega^{3}/c^{3}$, where $\xi\sim 1$ is a
‘universal constant’, is still meaningful within a numerical analysis.

####  3.4.3 Isolated neutron stars birthrate

Knowledge of the birth rate $\dot{N}_{\mathrm{INS}}$ of Galactic isolated
neutron stars (INSs) is also important in the context of neutron star
population analysis and, in particular, details of their spin down. In
general, several types of objects contribute to the total value of
$\dot{N}_{\mathrm{INS}}$: (i) ordinary rotation-powered pulsars (RPPs), (ii)
rotating radio transients (RRATs), and (iii) young radio-quiet neutron stars
(XDINS or the “Magnificent Seven”)444Although few other objects have also been
attributed to XDINS (e.g. (Pires et al., 2009; Rigoselli et al., 2019)) in
addition to seven classical sources Turolla (2009)., (iv) central compact
objects in supernova remnants (CCOs) and (v) magnetars.

The birthrate of each of them, however, is not simply proportional to the
number of observed sources but must be revealed assuming models of their spin-
down, activity, magnetic and thermal evolution. For instance, despite the
overwhelming number of ordinary pulsars, a significant fraction of neutron
stars appear to be born as magnetars Gotthelf and Vasisht (2000); Beniamini et
al. (2019).

However, the total INSs birthrate is not simply the sum of birthrates of their
different types. Keane and Kramer have shown that being simply summarized they
give the value $\gtrsim 5-10$ century-1 Keane and Kramer (2008) which is in
tension with a proposed rate of Galactic core-collapse supernova explosions

$\dot{N}_{\rm CCSN}=1-3\mbox{ century}^{-1}$ (38)

as it has been constrained by different methods Tammann et al. (1994); Diehl
et al. (2006); Keane and Kramer (2008); Schmidt et al. (2014).

Within these calculations, the birthrate of RPPs, as the most numerous class
of observed sources, can be estimated with formally small uncertainties. And
understanding of their spin-down evolution is crucial, since
$\dot{N}_{\mathrm{RPP}}\sim N_{\mathrm{RPP}}/t_{\mathrm{ej}}$.

So the analysis of pulsar kinetic current by Vranešević et al. Vranesevic et
al. (2004) initially gave $\dot{N}_{\mathrm{RPP}}\sim 1-2$ century-1. Later
Lorimer et al. Lorimer et al. (2006) used the same approach and estimated
$\dot{N}_{\mathrm{RPP}}$ to be $1.4\pm 0.2$ century-1. And finally, Vranešević
& Merlose Vranešević and Melrose (2011) perform an even more accurate
consideration, taking into account the pulsar spin-down law in the more
general form $P^{n-2}\dot{P}=$const. They found $\dot{N}_{\mathrm{RPP}}=0.8\pm
0.3$ century-1 for $n=2.9\pm 0.1$. All these results are more or less
consistent with each other.

Another way of measuring the birth rate of RPPs – a population synthesis
approach – tends to give even higher values (probably due to selection effects
being taken into account). For example, Faucher-Giguère and Kaspi Faucher-
Giguère and Kaspi (2006) reported $\dot{N}_{\mathrm{RPP}}=2.8\pm 0.5$, Ridley
and Lorimer found $\dot{N}_{\mathrm{RPP}}$ from $\sim 2$ to $5.3$ century-1
depending on the details of the assumed spin-down model Ridley and Lorimer
(2010). Although a spin-down law of the form
$P\dot{P}\propto\mu^{2}(k_{0}+k_{1}\sin^{2}\chi)$ was included in the
calculations in the latter paper, the exponential decay of a pulsar obliquity
was also assumed. which is inconsistent with the MHD model of pulsar
rotational evolution. Finally, Gullón et al. found
$\dot{N}_{\mathrm{RPP}}=2.5\pm 0.1$ century-1 including pulsar magneto-thermal
evolution in their simulations Gullón et al. (2014). Therefore, the pulsar
population synthesis generally predicts values of $\dot{N}_{\mathrm{RPP}}$
within 2-3 century-1 – about one pulsar per century more than the kinetic
current approach and which is close to the rate of core-collapse supernovae
themselves.

The obvious solution to the paradox of the overestimated birthrates of
Galactic neutron stars is to assume that their populations are connected.
Their observable properties depend on their current spin period, age, magnetic
field strength and surface temperature, while the transition between the
populations is driven by the magneto-rotational evolution Popov et al. (2010);
Kaspi (2010); Igoshev and Popov (2014); Gullón et al. (2015). The latter is
therefore crucial for understanding the properties and fate of INS. They need
to be considered as a connected population of sources, unified by the common
origin, rather than a set of separate classes555Kaspi proposed the term ’Grand
unification’ for this Kaspi (2010).

##  4 Propeller stage

The propeller stage is, probably, the least understood one. On one hand, from
the observational point of view, it is not easy to study this phase of
evolution. This is so, partly because an NS luminosity can be rather low at
this stage. In addition, the stage can be relatively short in many cases. On
the other hand, the theoretical analysis of propellers is complicated and
results of modeling often are not conclusive. We start by presenting different
theoretical approaches to model the propeller stage. And then, observational
results are reviewed and discussed.

###  4.1 Theory

The idea that a rapidly rotating magnetosphere can prevent penetration of
matter down to the surface of an NS was initially proposed by Shvartzman
Shvartsman (1970). However, as the paper by Shvartsman has been published in a
not-so-well-known journal, this idea became popular only later, mainly due to
the paper by Illarionov and Sunyaev Illarionov and Sunyaev (1975) who also
proposed the name for the evolutionary stage: propeller. The existence of the
centrifugal barrier for accretion onto NSs was also mentioned by Pringle and
Rees (1972), and discussed in more detail by Davidson and Ostriker (1973), who
also provided an estimate for the spin-down torque:
$K_{\mathrm{sd}}={\mu^{2}}/{R_{\mathrm{m}}^{3}}.$ Later on, a more detailed
study of the propeller stage was presented by Davies et al. (1979), Davies and
Pringle (1981).

Realistic 3D modeling of the interaction between the magnetosphere and
surrounding matter is a very challenging task. Not many attempts have been
made. That’s why, early studies were analytical and typically limited by some
simplified but physically motivated approaches.

At the propeller stage, a rapidly rotating magnetosphere transfers rotational
energy and angular momentum to the medium. Thus, the matter is heated and
moves away carrying some angular momentum, so that the NS spins down.
Depending on the relative efficiency of heating and cooling, the properties of
the shell around the magnetosphere can be different. External pressure at the
magnetospheric boundary determines its radius, $R_{\mathrm{m}}$, which can
significantly deviate from the Alfven radius, $R_{\mathrm{A}}$, see eq. (4).
Without a disc formation, typically $R_{\mathrm{m}}>R_{\mathrm{A}}$.

In the simplest approach, e.g. Shakura (1975), it is assumed that the external
matter is gravitationally captured and cools effectively. Thus,

$\rho(R_{\mathrm{m}})=\rho_{\infty}\left(\frac{R_{\mathrm{m}}}{R_{\mathrm{G}}}\right)^{-3/2}.$

This provides a large density at the magnetospheric boundary. The
magnetospheric boundary is determined from

$\frac{B(R_{\mathrm{m}})^{2}}{8\pi}=\frac{1}{2}\rho(R_{\mathrm{m}})(\Omega
R_{\mathrm{m}})^{2}.$

So, $R_{\mathrm{m}}\propto B^{4/13}P^{4/13}\rho_{\infty}^{-2/13}v^{6/13}$.
Then, it is assumed that the matter is expelled with the rotational velocity
of the magnetosphere: $\Omega R_{\mathrm{m}}$. Thus, it carries away a huge
amount of angular momentum. Altogether, this results in a very effective spin-
down:

$K_{\mathrm{sd}}=\dot{M}\Omega R_{\mathrm{m}}^{2}.$ (39)

So, the spin period increases rapidly.

Both assumptions (effective cooling and large velocity of the expelled matter)
in this approach were strongly criticized. At first, due to significant
heating, the density at $R_{\mathrm{m}}$ might be much lower than

$\rho(R_{\mathrm{m}})=\rho_{\infty}\left(\frac{R_{\mathrm{m}}}{R_{\mathrm{G}}}\right)^{-3/2}.$

Then, the velocity $\Omega R_{\mathrm{m}}$ is too high for the outflow.

Davidson and Ostriker Davidson and Ostriker (1973) assumed that the angular
momentum is carried away with the free-fall velocity:

$K_{\mathrm{sd}}=\dot{M}v_{\mathrm{ff}}(R_{\mathrm{m}})R_{\mathrm{m}}.$ (40)

In the seminal paper Illarionov and Sunyaev (1975) the authors define the
magnetospheric boundary as the usual Alfven radius and it is assumed that
rotational energy loss is $(1/2)\dot{M}v_{\mathrm{ff}}(R_{\mathrm{m}})^{2}$.
Then the spin-down torque is:

$K_{\mathrm{sd}}=\dot{M}v_{\mathrm{ff}}(R_{\mathrm{m}})R_{\mathrm{m}}\frac{v_{\mathrm{ff}}(R_{\mathrm{m}})}{\Omega
R_{\mathrm{m}}}.$ (41)

Fabian Fabian (1975) obtained the radius of the magnetospheric boundary by
equating the magnetic pressure and stellar wind pressure. To calculate the
spin-down torque he used energetic considerations. In this approach, the
rotational energy losses $\frac{1}{2}I(\Omega_{1}^{2}-\Omega_{2}^{2})$ are
equal to the kinetic energy of escaping matter

$\frac{1}{2}\dot{M}\frac{2GM}{R_{\mathrm{m}}}.$

In this model, the spin-down rate is identical to the one by Illarionov,
Sunyayev Illarionov and Sunyaev (1975). It can be written in an alternative
form:

$K_{\mathrm{sd}}=\dot{M}\frac{GM}{R_{\mathrm{m}}}\frac{1}{\Omega}.$ (42)

Note, that $R_{\mathrm{m}}$ can be defined in different ways. What is
important, is that $K_{\mathrm{sd}}\propto\Omega^{-1}$. Then the
characteristic time scale of spin-down $\tau_{\mathrm{sd}}\equiv
P/\dot{P}\propto(P_{1}^{-2}-P_{2}^{-2})$. This means the time scale is defined
by the initial period.

In their very influential paper Davies and Pringle Davies and Pringle (1981)
discussed several regimes of the propeller. In their notation, the most
standard case, typically discussed in the literature, is dubbed “supersonic
propeller”. In this regime, $R_{\mathrm{m}}>R_{\mathrm{co}}$ and also
typically $R_{\mathrm{m}}\ll R_{\mathrm{l}}$. At the subsonic propeller stage
$R_{\mathrm{m}}<R_{\mathrm{co}}$, accretion is still not possible as matter
cannot cool efficiently. Finally, “a very rapid rotator” regime can be
realized in the situation when $R_{\mathrm{m}}\lesssim R_{\mathrm{l}}$, i.e.
when the magnetospheric radius is just slightly below the light cylinder
radius.

One of the main points, used and discussed by Davies and Pringle is related to
the thermal balance in the envelope around the magnetosphere. The existence of
a hot envelope results in a different pressure outside the magnetosphere with
respect to the case of an efficiently cooled inflow assumed for calculation of
$R_{\mathrm{A}}$.

At the supersonic propeller stage, the magnetospheric radius is:

$R_{\mathrm{m}}=R_{\mathrm{A}}\left(\frac{R_{\mathrm{G}}}{R_{\mathrm{A}}}\right)^{2/9}.$
(43)

It is assumed that at this stage the energy released at the magnetospheric
boundary is carried away by convection or turbulent motions. The turbulent
velocity at any radius is assumed to be equal to the sound velocity. These
assumptions result in the polytropic index $n=1/2$ and $\rho\propto r^{-1/2}$.

Thus, we obtain:

$K_{\mathrm{sd}}=\frac{1}{8\sqrt{2}}\left(\frac{R_{\mathrm{co}}}{R_{\mathrm{m}}}\right)^{3/2}\frac{\mu^{2}}{R_{\mathrm{m}}^{3}}.$
(44)

Matter is expelled from the magnetospheric boundary with the velocity
significantly lower than the free-fall velocity.

At the phase of a very rapid rotator, the linear velocity of the rotating
magnetosphere, $\Omega R_{\mathrm{m}}$, is very high. This heats up the
surrounding matter significantly. Effectively, this means that pressure is
constant in the envelope up to “infinity”, i.e. up to $R_{\mathrm{G}}$ where
properties are determined by the external medium (e.g., ISM in case of an
isolated NS). This stage can be realized at low accretion rates or with high
magnetic fields when the outer boundary of the heated envelope extends beyond
$R_{\mathrm{G}}$. Otherwise, the pressure gradient becomes important and so,
from the ejector stage the NS immediately come to the supersonic propeller
stage.

At the very rapid rotator stage matter is expelled with the rotational
velocity of the magnetosphere, but the density is low, so the spin-down rate
is low, too:

$K_{\mathrm{sd}}\approx\frac{\mu^{2}}{R_{\mathrm{m}}^{3}}.$ (45)

In addition to the propeller stage for which $R_{\mathrm{m}}>R_{\mathrm{co}}$,
Davies and Pringle Davies and Pringle (1981) discussed also a “subsonic
propeller” stage with $R_{\mathrm{m}}<R_{\mathrm{co}}$ (the stage was later
re-considered e.g., in Ikhsanov (2001)). Here, despite the centrifugal barrier
does not exist, matter cannot effectively accrete due to its high temperature.
For low accretion rates, the subsonic stage can be quite long, and accretion
can be postponed. However, instabilities at the magnetospheric boundary might
allow the matter to accrete. This situation was recently analyzed in the
framework of the so-called settling accretion model Shakura et al. (2012).

Now, let us say a few words about the numerical modeling of the propeller
stage. Wang and Robertson Wang and Robertson (1985) used numerical simulations
to model spin-down at the propeller stage. Specifically, they accounted for
instabilities at the magnetospheric boundary. They obtained a very rapid spin-
down:

$K_{\mathrm{sd}}=\dot{M}\frac{2\zeta_{\mathrm{WR}}}{\pi\beta}\Omega
R_{\mathrm{m}}^{2}.$ (46)

The pre-factor $({2\zeta_{\mathrm{WR}}})/({\pi\beta})$ is expected to be $>1$.
Notice the similarity with the equation by Shakura Shakura (1975), see eq.
(39).

What is also important, in their model $R_{\mathrm{m}}$ is a function of the
spin period, which is quite natural at the propeller stage, but this is
typically not taken into account explicitly in many analytical models.

Numerical modeling of the propeller stage in the case of an isolated NS
interacting with the interstellar medium was performed in Romanova et al.
(2003), see also their earlier study in Toropin et al. (1999). Fitting of the
numerical results provided the following scaling: $K_{\mathrm{sd}}\propto
B^{0.8}\Omega^{1.3}$. This is a very rapid spin-down, similar to the one
proposed in Shakura (1975), where $K_{\mathrm{sd}}\propto\Omega$. Curiously,
the dependence $\propto B^{4/5}$ can also be obtained analytically, as
discussed in Romanova et al. (2003). At first, the authors define the radius
of the magnetospheric boundary from the condition

$\frac{B^{2}}{8\pi}\frac{1}{R_{\mathrm{m}}^{6}}=\rho(R_{\mathrm{m}})(\Omega
R_{\mathrm{m}})^{2}.$

So, that

$R_{\mathrm{m}}=\left(\frac{f\mu^{2}}{2\dot{M}\Omega}\right)^{1/5}.$

And then, the spin-down torques is defined as

$I\frac{{\rm d}\Omega}{{\rm d}t}=-\frac{f\mu^{2}}{2R_{\mathrm{m}}^{3}}.$

Note, that here again $R_{\mathrm{m}}$ depends on the spin period.

Francischelli and Wijers Francischelli and Wijers (2002) compared several
variants of the outflow velocity at the propeller stage and application to
fossil discs around young NSs and to the source SMC X-1. They concluded that
the highest velocity – $\Omega R_{\mathrm{m}}$, – is not compatible with
observations of SMC X-1.

###  4.2 Observations

Observational data about NSs at the propeller stage are not very abundant.
Mainly, they are limited to compact objects in low-mass (LMXBs) or high-mass
(HMXBs) X-ray binaries with a disc formation. Unfortunately, measurements of
the spin-down rate at the propeller stage are not done, yet. Still, it is
possible to determine the critical luminosity at which the transition to this
stage occurs.

The first identification (not very secure at that time) of the propeller stage
was probably made for the X-ray pulsar V0332+53 in Stella et al. (1986). Some
other candidates have been also identified in this paper, including 4U 0115+63
which was later studied by Campana et al. (2001). This source demonstrated a
decrease of the X-ray flux by more than two orders of magnitude in nearly one-
half of a day. Another example – Aql X-1, – was presented in Campana et al.
(1998). This object is of special interest as it has a millisecond spin
period.

Several other sources demonstrate behavior which can be explained as on-set of
the propeller stage: SAX J1808.4-3658 Campana et al. (2008), ULX M82 X-2
Tsygankov et al. (2016), GRO J1744-28, and GX 1+4 Cui (1997).

In Tsygankov et al. (2017) the authors present observations of GRO J1008-57.
The initial goal was to detect the transition to the propeller stage.
Surprisingly, at some value (larger than the expected critical number) the
flux decrease stopped. Thus, the authors conclude that the source instead of
reaching the propeller stage, changed the mode of accretion. In this case, the
disc at $R_{\mathrm{m}}$ has temperature $T<6500$ K. Interestingly, this
regime is valid only for pulsars with relatively long spin periods:
$P>36.6\,B_{12}^{0.49}$ s. For smaller periods a usual propeller regime
exists.

In Campana et al. (2018) the authors consider various systems (LMXBs, HMXBs,
cataclysmic variables, and young stars) for which transitions to/from the
propeller stage have been proposed. Among the systems with NSs the authors
discuss three LMXBs – SAX J1808.4–3658, IGR J18245–2452, and GRO J1744–28, –
and three HMXBs – 4U 0115+63, V 0332+53, and SMC X-2. In addition, three
cataclysmic variables with accreting white dwarf and one young stellar object
have been considered.

The authors tested the basic equation for the critical luminosity at which the
transition to the propeller stage happens:
$L_{\mathrm{lim}}\propto\mu^{2}P^{-7/3}R^{-1}$. The analyzed systems provide
the possibility to check this equation for a wide range of magnetic moments,
spin periods, and radii of accretors. Surprisingly, the observed data
perfectly confirm the dependences. In particular, the authors parametrized the
equation as: $L_{\mathrm{lim}}\propto\mu^{\alpha}P^{-\beta}R^{-\gamma}$. And
the fit gives: $\alpha=1.9^{+0.6}_{-0.3}$, $\beta=2.3^{+0.8}_{-0.4}$,
$\gamma=0.9^{+0.3}_{-03}$.

More recently, the source GRO J1750–27 was proposed as a system with
accretion/propeller transition Lutovinov et al. (2019). The authors observed a
drastic drop in the X-ray flux accompanied by a decrease in the spin-up rate
of the X-ray pulsar. Interpretation of this decrease in luminosity as
transition to the propeller stage allows (for the known distance)
determination of the magnetic field of the NS: $\sim 4\times 10^{12}$ G.

An X-ray binary SWIFT J0850.8-4219 with a red supergiant donor is proposed to
be at the propeller stage due to its low luminosity in comparison with
expectations based on typical stellar wind mass losses from red supergiants De
et al. (2023). Still, this hypothesis deserves a more detailed analysis as the
present-day data are not sufficient to make a firm conclusion.

Future instruments, like Athena Barcons et al. (2012), with large collecting
area might allow detecting pulsations at low fluxes corresponding to the
propeller stage. This will give an opportunity to test directly theoretical
models of spin-down.

##  5 Accretion

###  5.1 General properties of accretion flows onto magnetized NSs

NSs may accrete mass under various circumstances and from different sources,
including ISM, gravitationally bound supernova ejecta (the case of fallback,
see section 2), or a donor star. Interaction with the accreting matter is
determined by the set of characteristic radii introduced in section 1: if both
Shvartsman and magnetospheric radii are smaller than $R_{\rm G}$, mass is
gravitationally captured. The flows inside the gravitational capture radius
depend on the amount of linear and angular momentum in the matter captured by
the star. If both are negligibly small, the star accretes spherically in the
regime of Bondi accretion (Bondi, 1952). Such a regime would be realized if
the star is immersed in a static uniform gas. The gravity of the star creates
a radial transonic inflow of matter that in its inner parts is modified by the
presence of the solid surface and magnetosphere. Though unrealistic, Bondi
accretion model allows to estimate the mass accretion rate in the case of
subsonic motion through a uniform ambient medium according to formula (5).

Accretion from a _supersonically_ moving uniform media may be thought of as
gravitational focussing of the flow: gravity changes the direction of a
streamline with a shooting parameter $b$ by an angle $\sim R_{\rm G}/b$,
resulting in a shock wave of a roughly conical shape that dissipates some of
the kinetic energy of the flow and allows for mass accretion onto the star
from the downstream direction. This is an idealized case of Bondi-Hoyle-
Lyttleton accretion (Edgar, 2004).

An important factor that modifies the two scenarios mentioned above is the net
angular momentum of the captured matter. For the net angular momentum of $j$
the size of the disc is determined by the circularization radius

$R_{\rm c}=\frac{j^{2}}{GM}.$ (47)

If $R_{\rm c}\gg R_{\rm m}$, the dynamics of the magnetosphere are determined
by its interaction with the disc. In the case of Roche lobe overflow,
$j\sim\sqrt{GM_{\rm tot}a}$, where $a$ is the binary separation in the system,
and $M_{\rm tot}$ is the total mass of the system, hence $R_{\rm c}\sim a$
(though usually several times smaller). If the NS captures the stellar wind of
its companion, $j$ is much smaller, and the accretion flow is much better
described by the Bondi-Hoyle-Lyttleton scenario (see however section 5.3.6).

All the accretion flows described above exist only outside the magnetosphere.
Within the magnetosphere, magnetic stresses redirect the plasma along the
field lines, creating _funnel flows_ that start at $R\sim R_{\rm m}$ in the
regions where the field lines interact with the inflow. Making assumptions
about the shape of the magnetic field lines allows us to predict the shape of
the funnel flow and the regions of the NS surface where the accreting matter
finally falls and releases its gravitational energy. For a dipolar magnetic
field, these accretion sites are associated with the polar caps or rings or
arcs around them. If a dipolar magnetosphere acquires mass in its equatorial
regions at a distance of $R_{\rm m}$, the matter, following a field line, hits
the NS surface at a polar angle $\sin\theta\sim\sqrt{R_{\rm NS}/R_{\rm m}}\sim
0.01$ (for the values typical for an accreting strongly magnetized NS, see
section 5.3.1).

Accretion onto magnetized stars is a subject of many numerical studies, that
are often designed for young stellar objects but are also applicable for NS
accretion with an identical set of dimensionless parameters. To describe the
accretion flow geometry and torques acting on an NS with a dipolar magnetic
field accreting from a thin disc, the required dimensionless parameters are
the so-called fastness parameter $\omega=\Omega/\Omega_{\rm K}(R_{\rm m})$
(where $\Omega_{\rm K}(R)=GM/R^{3}$ is the local Keplerian frequency), $R_{\rm
m}/R_{\rm NS}$, and magnetic angle $\chi$. The rotation of the NS is assumed
aligned with the rotation in the disc. We will talk about such simulations in
more detail later in Sec. 5.3.4. Simulations of Bondi and Bondi-Hoyle-
Lyttleton accretion onto magnetized stars are much less numerous (see however
the following section 5.2).

###  5.2 Isolated neutron stars

Accretion onto isolated NSs from the ISM was proposed more than half a century
ago by Shvartsman Shvartsman (1971, 1970) and Ostriker, Rees, and Silk
Ostriker et al. (1970) (the papers by Shvartsman were submitted slightly
earlier). All of these authors correctly noticed that due to gradual spin-down
or/and significant magnetic field decay an NS finally can start accreting the
interstellar gas. In Ostriker et al. (1970) it is estimated that for a typical
pulsar magnetic field accretion onto a low-velocity NS can start in a few
hundred million years due to spin-down and corresponding decrease of the
pulsar wind power which prevents accretion at earlier times.

These authors applied Bondi accretion to estimate the luminosities of
accreting isolated neutron stars (AINSs), see eq. (5). As they mainly
considered low-velocity NSs (in 1969-1970 it was not known that compact
objects could obtain large kicks), the estimated luminosity appeared to be
$\sim 10^{31}$ \- $10^{32}$ erg s-1. Shvartsman also accounted for a ‘negative
feedback’ related to the heating of the external medium by the NS emission.
This somehow lowered the expected luminosity and shifted the spectrum from
X-rays to UV.

Interestingly, in the original papers Shvartsman (1971); Ostriker et al.
(1970) the authors did not discuss in detail the possibility of detecting
AINSs as individual X-rays sources but mainly focused on the influence
(heating and ionization) of the accretion luminosity on the surrounding
medium. Later in the 1990s different authors intensively discussed
observations of individual accreting isolated sources with space X-ray
observatories.

Blaes et al. analyzed the observability of AINSs in several papers Blaes and
Rajagopal (1991); Blaes and Madau (1993); Madau and Blaes (1994); Blaes et al.
(1995). In particular, in Blaes and Madau (1993) the authors calculated the
number of AINSs potentially detectable by ROSAT, and obtained that about
several thousand such sources can be found. This very optimistic prediction is
mainly due to the fact that the initial velocities of NSs were underestimated
at that time. Despite the authors discussed some evolutionary aspects,
finally, for their estimates they assumed that all old NSs accrete at the
Bondi rate. This is, of course, an oversimplification. Reasonable accounting
for spin evolution for realistic velocities of NSs has been done only in Popov
et al. (2000) and the realistic magnetic field distribution was added later by
Boldin and Popov (2010). In Blaes et al. (1995) the authors considered the
influence of the AINS emission onto the surrounding medium and ‘negative
feedback’ on a deeper level than it has been done in early studies. Naturally,
this process reduces the number of potentially observable AINSs.

3D hydrodynamical simulations of accretion onto an INS were presented in
Blondin and Raymer (2012). Magnetic field was not explicitly included in
calculations but the authors modeled accretors of different sizes comparable
with magnetospheric radius. The accretion rate was found to be reduced by a
factor $\sim 2$ in comparison with the Bondi-Hoyle-Littleton rate (eq. 5 with
$\xi_{\mathrm{acc}}=\pi$).

Accretion onto isolated magnetized NSs was studied numerically in a series of
papers by Romanova and her colleagues Toropin et al. (1999); Toropina et al.
(2001); Romanova et al. (2003); Toropina et al. (2012). The main result
presented in Toropina et al. (2012) is related to the dependence of the
accretion rate on the magnetic field of the NS. For magnetic field $\sim
10^{12}$ G the authors obtain a reduction in the accretion rate by a factor
$\sim 2$ with approximate dependence $\dot{M}\propto\mu^{-0.4}$.

An important question related to accretion onto INSs is about the spin
properties of these objects. Naively, one can expect that an AINS is just
spinning down, and the braking torque is
$\propto{\mu^{2}}/{R_{\mathrm{co}}^{3}}$ as the external medium does not
bring, on average, angular momentum. However, the picture is more complicated
due to turbulence. Accreting matter brings fluctuating angular momentum which
is equal to zero only after averaging over a relatively long time $\sim
R_{\mathrm{G}}/v$. Up to our knowledge, for the first time, this situation was
briefly considered in Lipunov and Popov (1995). In this paper the authors
proposed that an AINS does not spin down continuously, but can reach some
quasi-equilibrium period $P_{\mathrm{turb}}\propto\mu^{2/3}$ around which the
spin fluctuates. The situation was considered in more detail in a brief note
Popov et al. (2001).

The evolutionary equation for the spin frequency can be written as:

$\frac{{\rm d}\mathbf{\Omega}}{{\rm
d}t}=\mathbf{F}+\mathbf{\Phi},\,F=-\frac{k_{\mathrm{t}}\mu^{2}}{I\,R_{\mathrm{co}}^{3}},\,\Phi\propto\frac{\dot{M}j}{I}.$
(48)

Here $\mathbf{\Phi}$ is the turbulent torque and
$\left<\mathbf{\Phi}\right>=0$. Specific angular momentum $j$ is determined as
min$[v_{t}(R_{\mathrm{G}})R_{\mathrm{G}},v_{\mathrm{K}}(R_{\mathrm{A}})R_{\mathrm{A}}]$,
where $v_{\mathrm{K}}$ is Keplerian velocity and $v_{t}$ is turbulent velocity
taken at the gravitational capture radius. The coefficient $k_{\mathrm{t}}$ is
usually taken to be equal to 1/3.

The evolution of the spin frequency can be represented as diffusion. Then for
$\Omega_{\mathrm{turb}}=2\pi/P_{\mathrm{turb}}$ we have:

$\Omega_{\mathrm{turb}}^{2}=\int_{0}^{\infty}\Omega^{4}\exp{[-V(\Omega)/D]}{\rm
d}\Omega/\int_{0}^{\infty}\Omega^{2}\exp{[-V(\Omega)/D]}{\rm d}\Omega.$ (49)

Here

$D=\frac{1}{6}\left(\frac{\dot{M}j}{I}\right)^{2}\frac{R_{\mathrm{G}}}{v}$

and

$V(\Omega)=\frac{\mu^{2}}{3GMI}|\Omega|^{3}.$

For $j=v_{\mathrm{t(}}R_{\mathrm{G}})R_{\mathrm{G}}$ we obtain:

$P_{\mathrm{turb}}=3.9\times
10^{8}\,\mu_{30}^{2/3}n_{1}^{-2/3}I_{45}^{1/3}\left(\frac{M}{1.4\,M_{\odot}}\right)^{-26/9}\left(\frac{v}{100\,\mathrm{km}\,{\mathrm{s}}^{-1}}\right)^{43/9}\left(\frac{R_{t}}{2\times
10^{20}\mathrm{cm}}\right)^{2/9}\rm s.$ (50)

Here $R_{\mathrm{t}}$ is the characteristic scale of interstellar turbulence
for which $v_{\mathrm{t}}(R_{\mathrm{t}})=10$ km s-1 and the turbulent
velocity depends on the spatial scale as $v_{\mathrm{t}}\propto r^{1/3}$.
Normalized units are $I_{45}=I/10^{45}\,\mathrm{g}\,\mathrm{cm}^{2}$ and
$\mu_{30}=\mu/10^{30}\mathrm{G}\,\mathrm{cm}^{3}$.

Numerically, the spin evolution of AINSs was studied in Prokhorov et al.
(2002). The authors studied the spin evolution at the stage of accretion of a
population of isolated NSs with realistic velocity and magnetic field
distributions, but for a constant ISM density. After the onset of accretion,
an NS is spinning down: ${\rm d}\Omega/{\rm
d}t\approx\mu^{2}/R_{\mathrm{co}}^{3}$. The influence of turbulence starts to
be significant when the spin period reaches the value

$P_{\mathrm{cr}}=2\pi\frac{\mu}{\sqrt{GM\dot{M}j}}.$ (51)

Typically this happens in $\sim 10^{5}$ \- $10^{7}$ yrs. After that, the
period on average continues to increase, finally reaching $P_{\mathrm{turb}}$
on a time scale $\gtrsim 10^{9}$ yrs, and fluctuating around this value. A
typical time scale for $\dot{P}$ fluctuations is $R_{\mathrm{G}}/v.$ The spin
period distribution of AINSs appears to be very wide with typical values in
the range $\sim 10^{5}$ \- $10^{7}$ s. Note, that the authors did not discuss
the influence of shocks around AINSs on the angular momentum transfer by
turbulent accreted matter.

As accretion onto INSs proceeds at low rates, the settling accretion regime
might be applicable. We describe this mode of accretion in Sec. 5.3.6. To the
problem of AINSs this model was applied in Popov et al. (2015). It was
demonstrated that in the settling accretion regime, low-velocity NSs are
expected to be transient sources with a typical duration of bursts of about a
few hours and luminosity $\sim 10^{31}$ erg s-1. Persistent luminosities are
expected to be three orders of magnitude lower, and so out of the reach of the
present-day X-ray surveys. In addition, the authors proposed that inside
convective envelopes AINSs can spin down much faster, and so they can reach
quasi-equilibrium in a shorter time than it has been discussed in Prokhorov et
al. (2002).

In the early 1990s, high hopes to detect AINSs were related to ROSAT
observations Treves and Colpi (1991). Step by step, it became clear that due
to different effects, the number of observable AINSs cannot be very high (see
e.g., Manning et al. (1996) and references therein). In the late 1990s, early
enthusiasm was changed by pessimism related to the non-detection of any AINS
by ROSAT (see a review in Treves et al. Treves et al. (2000)). Finally,
calculations presented in Popov et al. (2000, 2000) demonstrated that it is
impossible to expect even a few AINSs in the ROSAT data.

One of the most detailed searches for INSs in the ROSAT data was presented in
Turner et al. (2010). Up to now there are no known AINSs or good candidates.
Still, the discovery of such sources will be a major breakthrough in
understanding low-rate accretion and long-term evolution of isolated NSs.

###  5.3 Accretion in binaries

####  5.3.1 Variety of NSs in binary systems

Among the massive stars producing NSs, most are components of binary or
multiple systems (Sana et al., 2012). Existence of a binary companion may
affect the evolution of the NS in multiple ways, from tidal torques adjusting
the rotation of the progenitor to mass exchange. Observationally, NSs in
binaries manifest themselves as X-ray binaries (XRBs) and as millisecond
pulsars (MSP). In both cases, the presence of a binary companion is crucial
for the rotational properties of the NS. MSPs are radiopulsars thought to have
accretion episodes in their past (Manchester, 2017). The connection between
the two classes of systems is confirmed by the existence of accreting
millisecond X-ray pulsars (AMXPs, Archibald et al. (2009)). We will cover
these classes of objects later in section 5.3.7.

Classification of XRBs may be done in different ways, using such criteria as
mass transfer type (Roche-lobe overflow or wind capture) or magnetic field of
the accretor, but the traditionally adopted distinction is between high-mass
and low-mass XRBs (HMXBs and LMXBs, respectively), with the mass of the donor
being the crucial parameter. Higher donor mass (several Solar masses or
larger) implies a relatively young age of the NS, a large magnetic field, and
often mass transfer through the wind. Low donor mass makes the formation of
the binary an unlikely scenario, as progenitor systems should be disrupted by
supernova explosions. Thus, many LMXBs, especially in dense stellar
environments like globular clusters, should have been formed by exotic
processes such as binary companion exchange (Bahramian and Degenaar, 2023), or
require a finely tuned natal kick. In terms of the number of active systems,
the low formation probability is compensated by the abundance of low-mass
stars and their longevity. The NSs in such systems are on average much older,
which implies low magnetic field strengths. Also, accretion from a low-mass
donor is usually associated with Roche-lobe-overflow mass transfer.

Low magnetic fields of LMXBs make it difficult to study their rotational
evolution. With the exception of AMXPs, their spin periods are inferred
indirectly from the other variability patterns such as post-burst oscillations
in X-ray bursters (Strohmayer et al., 1996) or quasi-periodic oscillations in
the sources lacking coherent periodic pulsations. The latter includes the so-
called ’Z‘ and ’Atoll‘ sources introduced by Hasinger and van der Klis (1989),
containing most of the known LMXBs (see Bahramian and Degenaar (2023) for a
review.) It is however unclear if the QPO frequencies are reliable indicators
of the spins of these objects, see Méndez and Belloni (2007)).

Because the matter accreted by a non-magnetized NS through an accretion disc
has the net angular momentum equal to Keplerian at the inner disc edge $R_{\rm
in}\simeq 10$km (determined either by the surface of the star or by its
gravity field), angular momentum is gained at the time scales only several
times smaller than the time scales of mass growth. For accretion at the
Eddington limit ($\dot{M}\simeq 10^{18}\rm g\,s^{-1}$), spinning up to a
several millisecond period would take millions of years.

Strong magnetic field makes the spin-up process much more efficient, as a
large magnetosphere intercepts accreting matter with a larger
($\propto\sqrt{R_{\rm m}}$, if accretion runs through a disc) net angular
momentum. On the other hand, rapid rotation in combination with strong
magnetic fields suggests pulsar torques and pulsar wind formation. Even more
importantly, rapid rotation of the NS may lead to propeller regime and thus
quenching of mass transfer (see section 4). Hence, the spin evolution in HMXBs
is usually thought of as deviations from an equilibrium between spin-up and
spin-down torques. The equilibrium point shifts with the changing mass
accretion rate and some other parameters, likely related to the geometry of
the accretion flows.

Many HMXBs are observed as X-ray pulsars (XRPs, i. e., sources of coherent
periodic variability in the X-ray range, clearly related to rotation).
Studying X-ray pulsations allows to probe their rotational evolution in real
time and on long archival data series. One of the most important projects of
this kind is the monitoring program held out with Fermi/GBM in the hard X-ray
range (Malacaria et al., 2020). The series obtained by Fermi/GBM span more
than 10 years of observations on 39 X-ray pulsars.

The majority of X-ray pulsars are members of binaries with Be-class donors
(Be/X-ray systems, Raguzova and Popov (2005); Reig (2011)). The population of
such systems is diverse, with the orbital periods in the range $20-200$d and
NS spins $1-10^{3}$s. Orbital and spin periods are known to correlate with
each other (Corbet, 1984), supporting the idea that that their rotation is
close to equilibrium between different torques acting on the star (see section
5.3.4). The exact origin of this correlation is however unclear. Mass transfer
rate in Be/X-ray systems is strongly modulated by orbital eccentricity and
other factors such as the mass accumulation in the decretion disc of the star.
This results in a complex phenomenology of periodic and aperiodic flares,
during which the NS is spun up by accreting matter, and quiescent states where
the rotation slows down. This picture is not universal and allows for objects
that are, for example, spun up even during the quiescent state, suggesting
that the equilibrium is reached beyond the time span of observations. For a
given object, however, $\dot{P}$ is well correlated with the luminosity
(Malacaria et al., 2020).

Another important category of XRPs are objects accreting from the winds of
early-type supergiants. The prototype and the best studied object of this
class is Vela X-1 (Kretschmar et al., 2021). Accretion disc may form in such
systems but not necessarily, as the circularization radius is not necessarily
larger than the radius of the magnetosphere.

In a homogeneous wind, the net angular momentum is (see for instance
Illarionov and Sunyaev (1975)) $j_{\rm 0}\sim\Omega_{\rm orb}R_{\rm G}^{2}$.
This angular momentum is aligned with the rotation of the binary. Strong
inhomogeneity of the wind would create an additional randomly oriented net
angular momentum component $j_{\rm w}\sim R_{\rm G}v_{\rm w}\gg j_{0}$, where
$v_{\rm w}\gg R_{\rm G}\Omega_{\rm orb}$ is the velocity of the wind. Hence,
for such systems, the angular momentum of the accreting matter is likely to be
randomly aligned.

####  5.3.2 Spin-up by an aligned disc

First, let us consider a planar thin conducting accretion disc interacting
with the dipolar magnetic field of a star. Even in its simplified version, the
problem of disc-magnetosphere interaction does not have a general analytical
solution. A thorough review of the existing models and approaches to disc-
magnetosphere interaction was made, for instance, by Lai (2014).

The general picture is dictated by the existence of a magnetosphere and the
strong radial decrease in magnetic tensions. One might expect magnetic field
lines threading the disc within some depth $\Delta R$ that is determined by
intermediate-scale processes such as turbulence and reconnections between the
stellar magnetic fields and the fields frozen into the accreting matter. Large
$\Delta R\gtrsim R$ models are usually referred to as ‘magnetically threaded
disc’, or MTD (Ghosh and Lamb, 1979; Wang, 1995) models. Within the
penetration region, the internal angular momentum flow in the disc, normally
governed by effective viscosity forces (Shakura and Sunyaev, 1973), is
contributed by the torque from the distorted magnetospheric lines. Purely
poloidal magnetic field of an aligned dipole does not produce any torque, but
the relative motion of the magnetosphere with respect to the disc creates
toroidal field components $B_{\varphi}^{\pm}$, having different signs on the
different surfaces of the disc. Toroidal fields should also change sign at the
corotation radius and are likely of the same order as poloidal fields.
Magnetic field torque acting on the NS has two contributions that are not
clearly distinguishable: the angular momentum coming with the accreting matter
itself and the torque produced by Maxwell stresses
$B_{z}B_{\varphi}^{\pm}/4\pi$ (where $B_{p}$ is the poloidal field). As the
flow inside the magnetosphere is bound to corotate with the NS, its angular
momentum is also primarily transmitted by magnetic stresses. While the net
angular momentum at the magnetospheric boundary is $\Omega R_{\rm m}^{2}$, its
value near the surface of the NS is $\sim\Omega R_{\rm NS}^{2}$, that is
normally several orders of magnitude smaller.

The torque acting on the surface of the NS from the disc may be calculated as
the star-disc interaction (SDI) torque, containing both spin-up and spin-down
parts

$K_{\rm SDI}=-\int RB_{z}B_{\varphi}R{\rm d}R{\rm d}\varphi,$ (52)

where integration is performed over the surfaces of the disc threaded by
magnetic field lines (usually, upper and lower surfaces). Straightforward
integration of this expression requires knowledge about both field components.
In general, the SDI torque may contain a spin-up part created by the field
lines threading the disc inside the corotation radius, and a spin-down part
from the lines extending outside (their tensions marked with green and red
arrows in Fig. 4). Magnetic field opening and screening (see below) tend to
decrease the contribution of the outer parts of the disc, making SDI torque
essentially a spin-up torque.

The MTD approach assumes that the poloidal field interacting with the disc is
equal or at least close to the original dipolar field in a broad range of
radii $\Delta R\gtrsim R$, and the toroidal fields are generated by the
difference in the rotation velocity between the star and the disc. The
expression used for the toroidal fields has the form $B^{\pm}_{\varphi}(R)=\pm
B_{z}(R)f(\Omega/\Omega_{\rm K})$, where the function $f(x)$ should meet a
couple of requirements: $f(1)=0$, $f(x<1)<0$, $f(x>1)>0$, and $|f(x)|\lesssim
1$. The last condition is related to the stability of predominantly toroidal
force-free magnetic field configurations. In these assumptions, the torque
related to disc-magnetosphere interaction is convenient to normalize with the
accretion torque spinning the NS up,

$K_{\rm SDI}=\dot{M}\sqrt{GMR_{\rm m}}n(\omega).$ (53)

This kind of normalization may also be used in a more general case when many
different aligned torques affect the rotation of the star at the same time.
The simplified picture where all the matter accreting through the disc gives
its angular momentum at $R_{\rm m}$ to the star corresponds to $n(\omega)=1$.
In general, the dimensionless factor $n(\omega)$ accounts for all the
processes of angular momentum exchange between the disc and the star (Wang,
1995; Shi et al., 2015).

Figure 4:  Three-dimensional sketch illustrating the MTD model of field
penetration into an aligned disc. The blue circle in the disc plane marks the
corotation radius. The field lines have dipolar and toroidal components, the
latter related to the motion of the disc and limited by the absolute value of
the seed poloidal field. Arrows show the rotation of the disc and the tension
along the field lines.

The classic MTD model assumes that the magnetic field retains its poloidal
component well outside the radius of the magnetosphere, as one would expect if
magnetic diffusivity in the disc were high. Rotation of the disc creates
additional azimuthal components scaling with the poloidal dipolar field
components $B_{\varphi}\sim B_{\rm p}$. Integrating the torques over the
surface of the disc, one yields for the torque (assuming $B_{\varphi}=\pm
B_{z}$)

$K_{\rm SDI}=\frac{1}{3}\left(\frac{1}{R_{\rm m}^{3}}-\frac{2}{R_{\rm
co}^{3}}\right)\mu^{2}.$ (54)

Such a configuration allows for a rotational equilibrium, where part of the
disc (inside the corotation radius) spins the star up, while the outer regions
spin the star down (see Figure 4). At a constant mass accretion rate, this
kind of equilibrium keeps the corotation radius close to the magnetospheric
radius.

Decreasing the penetration depth restricts the integration range in (52) to a
narrow strip where the rotation frequency of the disc is likely to be faster
than $\Omega$. The spin-up torque should still be close to the accretion
torque (53), because of the angular momentum advected with the matter.
However, as we show later in section 5.3.4, this spin-up is opposed by other
forces, not directly related to star-disc interaction.

####  5.3.3 Spin-up of an inclined dipole

The very existence of an XRP phenomenon requires an NS with misaligned
rotational and magnetic axes. Observational estimates of magnetic angles
$\chi$ (see section 3) in X-ray binaries cover a broad range, from nearly
aligned objects towards nearly orthogonal (see in particular the recent
results from IXPE Doroshenko et al. (2022)). The observational estimates for
magnetic angles in MSPs are strongly model-dependent (Johnson et al., 2014)
but also favour a broad distribution. For all these reasons, it is important
to be able to relax the alignment condition for the dipolar magnetic field of
the star. At the same time, its rotational axis is likely to be aligned with
the spin-up torque, with a couple of important exceptions we consider later in
section 5.3.5. Relaxing the magnetic obliquity leads to several consequences.
First, the condition for the magnetospheric boundary is now satisfied at
different radii for different azimuthal angles. According to Bozzo et al.
(2018), this leads to a correction factor related to the magnetic field
strength dependence on the polar angle. Even more importantly, magnetic
inclination breaks the axial symmetry of the flow, and the funnel flow becomes
concentrated to a couple of streams restricted in the azimuthal direction.

An MTD approach allows calculating the spin-up torque, assuming that the
magnetic field interaction with the plasma of the disc creates a non-dipolar
component (generalization of $B_{\varphi}^{\pm}$ in the aligned case) and
consequently angular momentum exchange between the disc and the star. Accurate
calculation in the MTD approach was performed by (Wang, 1997) and involved
three surfaces: upper, lower, and inner faces of the disc. The unperturbed
magnetic field line of an inclined dipole was complemented by an additional
field component directed along the relative velocity difference between the
disc and the magnetosphere. Integration over all the surfaces leads to a
torque that changes sign at some critical fastness parameter values close to
unity.

####  5.3.4 Spin-down torques and equilibrium

The MTD picture described above has a major caveat related to the opening of
the field lines due to asynchronous rotation between the disc and the star
(Aly and Kuijpers, 1990; Lai, 2014). A more realistic picture involves field
lines loaded with mass and extending to infinity, which creates possibilities
for mass ejections. A lot of inspiration here comes from the studies initially
aimed at magnetized young stellar objects that evolve under a similar
combination of torques. In Matt and Pudritz (2005), the existence of
magnetized winds draining the star of angular momentum and carrying away some
part of the accreting matter was proposed as an explanation for the spin
period distribution in T Tau stars. Many detailed 3D MHD simulations (Romanova
et al., 2013; Zanni and Ferreira, 2013; Ireland et al., 2022) designed for
young star magnetospheres are scalable for the smaller and rapidly rotating
magnetospheres of NSs with identical dimensionless parameters (fastness
$\omega$, magnetic angle, and disc thickness), that allows to make
implications about XRPs as long as relativistic effects are ignored. All the
simulations reproduce a moderate penetration depth of the order disc
thickness, field line opening and magnetized ejections. The estimated torque
related to the ejections is of the order $K_{\rm sd}\sim\Omega R_{\rm m}^{2}$,
similar to the braking torque given by the propeller theory $K_{\rm
sd}\sim\mu^{2}/R_{\rm co}^{3}$. The exact value contains correction factors
depending on the fastness parameter (Ireland et al., 2022). Fully relativistic
MHD simulations (Parfrey and Tchekhovskoy, 2023; Das and Porth, 2023) confirm
the general picture, but are so far limited only to very compact
magnetospheres.

In an accreting magnetosphere, some of the field lines remain open and free of
the matter supplied from the disc. As these lines tend to corotate with the
star and extend beyond the light cylinder, they create a braking torque
similar to the pulsar torque considered in section 3. Particle creation due to
pair cascades is in this case probably suppressed, but the contribution from
the vacuum losses remains as long as there are open magnetic field lines. On
the other hand, interaction with the disc changes the field shape outside the
radius of the magnetosphere, and the resulting electromagnetic torque is
likely to be higher than in the case of an ejector. In the aligned case,
Parfrey et al. (2016) estimate the value of the torque as

$K_{\rm Parfrey}\simeq\zeta^{2}\left(\frac{R_{\rm l}}{R_{\rm
m}}\right)^{2}\frac{\mu^{2}}{R_{\mathrm{l}}^{3}}$ (55)

where $\zeta\sim 1$ is a dimensionless parameter, and $K_{\rm sd}$
corresponding to the pulsar braking torque given by equation (29). Three-
dimensional MHD simulations confirm the picture drawn by Parfrey et al. (2016)
and also show that the open field lines of the accreting magnetosphere behave
like Poynting-dominated jets in accreting black hole systems powered by the
rotation of the compact object (Das and Porth, 2023; Parfrey and Tchekhovskoy,
2023).

It is unclear how applicable these results are to a more general geometry. In
particular, is this spin-down torque aligned with the rotation of the NS or
with the disc that opens up the field lines? How does this particular spin-
down torque depend on the magnetic inclination of the rotator? Answering these
questions would require full-scale simulations with non-zero magnetic and disc
inclination angles.

Combined, two or more torques of different signs, dependent on mass accretion
rate and the spin on the NS, result in a rotational evolution equation

$\frac{d}{dt}\left(I\Omega\right)=K_{\rm su}-K_{\rm sd}=n(\omega)K_{\rm su},$
(56)

which has a single equilibrium point at $n=0$. Different models using
different assumptions about the spin-down torques differ in the dimensionless
factor $n(\omega)$. The solution $n(\omega)=0$ is usually close to the
propeller limit $\omega=1$, or, in terms of spin period

$P_{\rm eq}\simeq 2\pi(GM)^{5/7}\mu^{6/7}\dot{M}^{-3/7}\simeq
5\left(\frac{\mu}{10^{30}\rm
G\,cm^{3}}\right)^{6/7}\left(\frac{\dot{M}}{10^{16}\rm
g\,s^{-1}}\right)^{-3/7}\rm s.$ (57)

Applying this expression to an XRP with a known spin period and a reliable
estimate of the mass accretion rate (that is proportional to the bolometric
accretion luminosity) allows to estimate the magnetic moment in the assumption
of equilibrium. In a strong flare, when the spin-up torque dominates over the
spin-down one, the magnetic moment may be probed by the observed spin-up rate.
Several methods of magnetic field estimation using the spin-up and spin-down
torques are described in Chashkina and Popov (2012).

Depending on the variability of the mass supply rate and potentially other
parameters, some objects might oscillate between accretor and propeller
regimes, while others proceed to accrete in the low-accretion-rate state.
There is observational evidence that some Be/X-ray systems enter the propeller
stage (see section 4), while others retain a small persistent mass accretion
rate during quiescence (Tsygankov et al., 2017). The reason for such a
dichotomy is unclear and may be connected to additional parameters such as
magnetic angle, or with the mass accretion rate variations in the disc.

The time required to reach the equilibrium state may be estimated as

$t_{\rm su}\simeq\frac{I\Omega_{\rm eq}}{K_{\rm su}}\sim\frac{I\Omega_{\rm
eq}^{4/3}}{\dot{M}(GM)^{2/3}}\sim 10^{3}\frac{I}{10^{45}\rm
g\,cm^{2}}\left(\frac{\dot{M}}{10^{18}\rm g\,s^{-1}}\right)^{-3/7}\rm yr,$
(58)

where the radius of the magnetosphere was assumed $R_{\rm m}\sim R_{\rm co}$
(that allows to exclude the value of $\Omega_{\rm eq}=2\pi/P_{\rm eq}$ from
the equation), that is a reasonable assumption near the equilibrium. Cases,
when the spin period changes considerably during the observational period, are
scarce and related to unique and probably catastrophic events. A bona fide
example of such a case is the super-Eddington XRP NGC300 ULX1 (Vasilopoulos et
al., 2019) which was spun up by about a factor of five during the period of
observations.

Detailed timing studies of XRPs allow testing of the models of magnetospheric
accretion flows. In particular, the flares of Be/X-ray systems are often
thought of as the case when spin-up torque dominates over all the spin-down
mechanisms, hence $\dot{\Omega}\propto K_{\rm su}\propto\dot{M}\sqrt{R_{\rm
m}}$. Classical scaling $R_{\rm m}\propto R_{\rm A}\propto\dot{M}^{-2/7}$,
implying $\dot{\Omega}\propto\dot{M}^{6/7}$, is generally supported by
observations (Bildsten et al., 1997). There are indications for a steeper
dependence at higher luminosities (Parmar et al., 1989) that might be a
signature for a more elaborate disc-magnetosphere interaction (Kluźniak and
Rappaport, 2007; Tsygankov et al., 2017).

####  5.3.5 Misaligned spin-up torques

If the spin-up torque does not change its direction (in particular, if its
direction is fixed by the orbital plane of the system), one would expect the
rotation axis of the NS to align with the torque within the spin-up time (58).
But on some of the evolutionary stages and in certain accretion regimes in
X-ray binaries the rotational axis of the star may deviate strongly from the
external spin-up torque. First, before the mass exchange, the accreting NS is
likely to have a rotational axis misaligned with the orbital plane due to the
kick received during the supernova explosion (Wang et al., 2006) or a three-
body process responsible for the creation of the system.

Another important case is accretion from an inhomogeneous wind when the
density and/or velocity gradients within the wind make the direction of the
torque essentially decoupled from the orbital plane (Wang, 1981). Numerical
simulations of such an accretion regime were done, for instance, by El Mellah
et al. (2018). The spin of the NS in this case is lower than the equilibrium
value, and its direction is involved in a random walk adjusting to the rapidly
changing external torque. The observed spin derivative in this case should be
weakly correlated with the flux, because of the relative orientation of the
torque and the angular momentum. Such a behavior is observed for Vela X-1
(Kretschmar et al., 2021) and other XRPs fed by supergiant winds.

Misaligned torques are also important as a possible trigger for magnetic angle
evolution. Adding the spin-up torque to Euler equations of rotational dynamics
(Biryukov and Abolmasov, 2021) reveals an additional term proportional to
$\sin\chi\cos\chi\sin\alpha\cos\alpha$, where $\alpha$ is the misalignment
angle between the rotational axis of the NS and the external spin-up torque,
that increases the magnetic angle during the alignment stage. As a result, the
alignment of the NS rotational axis increases the magnetic angle, increasing
the possibility of observing the NS as an XRP. Wind accretion with a randomly
aligned spin-up torque results in a random walk evolution of $\chi$,
increasing the number of strongly oblique pulsars.

####  5.3.6 Quasi-spherical accretion regime

If the circularization radius becomes smaller than $R_{\rm m}$, the Keplerian
disc can not form, and the torques acting on the NS should be revised. The net
angular momentum of the matter interacting with the magnetosphere is now much
smaller than Keplerian at the magnetosphere radius, meaning that the spin-up
torque is severely diminished. For a homogeneous wind, $j\sim\Omega_{\rm
orb}R_{\rm G}^{2}$ (Illarionov and Sunyaev, 1975).

The matter near the magnetospheric boundary is heated by shock waves to
temperatures close to virial, which implies a quasi-spherical envelope around
the magnetosphere. The rotation of the NS adjusts to the rotation of the
envelope via magnetic stresses at the edge of the magnetosphere.

The structure of the envelope depends on the presence of cooling and angular
momentum transfer within it. In the case of weak cooling, an elaborate
theoretical model was developed by Shakura et al. (2012). As the main process
responsible for angular momentum transfer, the authors considered an
anisotropic turbulent velocity field. Anisotropy allows to change the
direction of the angular momentum transfer (Kato and Yoshizawa, 1997). In
particular, predominantly radial motions lead to a flat distribution in
angular momentum, meaning that the rotation frequency in the envelope
$\Omega_{\rm env}\propto R^{-2}$. Shakura et al. (2012) adopt a more general
power-law scaling for the rotation frequency.

Equilibrium spin periods predicted by the quasi-spherical model are much
longer than for disc accretors, roughly corresponding to the equal net angular
momenta of the outer magnetosphere rim and the captured matter $\Omega R_{\rm
m}^{2}\simeq\Omega_{\rm orb}R_{\rm G}^{2}$, implying, up to a coefficient of
the order unity,

$P_{\rm eq,\,sph}\sim 10^{3}\left(\frac{v}{10^{8}\rm
cm\,s^{-1}}\right)^{4}\left(\frac{\mu}{10^{30}\rm
G\,cm^{3}}\right)^{8/7}\left(\frac{\dot{M}}{10^{16}\rm
g\,s^{-1}}\right)^{-4/7}\rm s.$ (59)

For identical magnetic moments and mass accretion rates, this value is about
two orders of magnitude longer than the equilibrium spin for disc accretion
(compare equation 57).

Such an accretion scenario may be applied to wind-fed HMXBs like Vela X-1, to
symbiotic wind-fed LMXBs like GX 1+4 (González-Galán et al., 2012), and
potentially to some of the Be/X-ray systems having very long spin periods that
require suspiciously large magnetic fields if interpreted as disc accretors
(Chashkina and Popov, 2012).

####  5.3.7 Recycled and transitional pulsars

About 10 percent of the RPP population belongs to millisecond pulsars (MSPs),
a category of short-period, weakly magnetized NSs with the estimated ages of
about Gyrs, well beyond the expected period of pulsar lifetime (Manchester,
2017). This population is understood as recycled, i. e. spun up by accretion
in a binary system to a period short enough to support a pulsar cascade. The
majority of MSPs have binary companions, and some of them have masses in
excess of $1.8$$M_{\odot}$.

Among the NSs with reliably measured rotation periods, the fastest rotator is
PSR J1728-2446ad (Hessels et al., 2006) having a rotation frequency of more
than 700Hz. Several other MSPs are rotating at similar rates $\sim 500-600$Hz.
There are indications that some non-pulsating NSs in LMXBs have frequencies
around 500-600Hz (Patruno et al., 2017). These values are several times
smaller than the Keplerian rotation rate ($\sim 2.5$kHz for a 2$M_{\odot}$star
with a radius of 10km), which is sometimes interpreted as a signature of an
additional braking mechanism. If the excess $\sim 0.3-0.5M_{\odot}$ of mass
observed in MSPs was gained via disc accretion, the expected rotation
frequencies should be pretty close to critical. For $I\simeq 10^{45}\rm
g\,cm^{2}$ and total angular momentum gained via accretion $\Delta
M\sqrt{GMR_{\rm NS}}\sim 10^{49}\rm g\,cm^{2}s^{-1}$, this implies a spin
period of $P_{\rm s}\lesssim 0.6$ms, or a frequency $\sim 1.5$kHz. One factor
limiting neutron star rotation frequency is gravitational waves (Bildsten,
1998; Owen et al., 1998). The upper spin limit may be also related to the
electromagnetic braking term (see Parfrey et al. (2016) and sec. 5.3.4), or
simply to the decrease in spin-up efficiency with the shrinkage of the
magnetosphere (Patruno et al., 2012). Population synthesis simulations run by
Tauris (2012) shows that in certain assumptions the observed spin limit may be
interpreted as a consequence of the braking torques slowing the NS down on the
propeller stage.

Our understanding of MSPs as a product of binary evolution is supported by the
existence of LMXBs containing NSs with properties similar to MSPs (accretion-
powered millisecond pulsars, AMXPs, see Patruno and Watts (2021)) and, most
importantly, transitional MSPs (TMSPs) that combine the radio emission typical
for RPPs with an evidence for an accretion disc (Papitto and de Martino,
2022). Observations show TMSPs switch between different activity states. The
prototypical object PSR J1023+0038, in particular, is observed either as a
radiopulsar with an X-ray luminosity $L_{\rm X}<10^{33}\rm erg\,s^{-1}$, or in
an X-ray bright state during which it has pulsations in the X-ray range, but
is not observed in radio (Archibald et al., 2015). X-ray pulsations are
observed in both phases, but the pulse shapes are different. In the AMXP
state, the pulse profiles are dominated by the first overtone. An obvious
interpretation of this behavior is the transition between ejector and accretor
states (see sec. 1). Pulsating X-ray emission is interpreted as a signature of
a funnel flow within the magnetosphere. An alternative explanation proposed by
Veledina et al. (2019) considers that the source is always rotation-powered,
and the X-ray emission observed during the bright states is produced by the
disc existing outside the light cylinder and reprocessing part of the pulsar
wind power. Anisotropy of the pulsar wind leads to a modulation of the X-ray
emission with a doubled spin frequency of the pulsar.

####  5.3.8 Magnetic field burial by accretion

Ohmic magnetic field decay given by equation (12) predicts a late-time
exponential decay on millions of years, but is definitely inapplicable to
objects as old as hundreds of millions of years, such as MSPs and AMXPs. The
situation is unclear, and most studies assume that Ohmic losses are negligibly
small during late evolution (see discussion in Faucher-Giguère and Kaspi
(2006)).

The proposed long accretion histories and large accreted masses in MSPs imply
that most of the crust material in these objects is renewed. This should also
affect the magnetic field strength and structure. The idea of the field being
screened by a layer of accreted matter was proposed already in Bisnovatyi-
Kogan and Komberg (1974). It is known to be an important effect in other
circumstances, such as fallback accretion (see sec. 2).

Due to the presence of the magnetosphere, an XRP is primarily accreting onto a
relatively small fraction of its surface, about the size of the two polar caps
($\sim 2\pi R_{\rm NS}^{4}/R_{\rm m}^{2}$) or smaller. Magnetic field tensions
are able to support only a small layer with a surface density of

$\Sigma\simeq\frac{B^{2}}{4\pi g}\sim 3\times 10^{8}B_{12}^{2}\,{\rm
g\,cm^{-2}},$ (60)

where $g$ is surface gravity. This led to the suggestion that the newly
accreted matter spreads over the surface of the NS from the magnetic poles
toward the equator, screening the field and distorting the shapes of the field
lines (Melatos and Phinney, 2001). In an idealized axisymmetric approach, the
dipolar component of a buried magnetic field decays approximately
quadratically with the accreted mass. For instance, Payne and Melatos (2004)
approximate their results with a power-law relation

$\mu\simeq\mu_{\rm i}\left(\frac{\Delta M}{4.6\times
10^{-5}M_{\odot}}\right)^{-2.25\pm 0.22}.$

Accretion of only about $10^{-4}M_{\odot}$ would lead to a magnetic moment
decrease by a factor of several. Local magnetic fields near the magnetic
equator in this scenario are not buried and are even expected to grow due to
flux conservation, forming a “magnetic tutu”. The increase in field gradients,
both vertical and latitudinal, increases Ohmic losses, which was also supposed
to contribute to the long-term evolution of the field (Konar and Bhattacharya,
1999).

It should be noted that the spreading of the accreted matter over the surface
of the NS is accompanied by a number of instabilities acting on different
scales and caused by different processes (Manchester, 2017). Local interchange
instabilities (Litwin et al., 2001; Kulsrud and Sunyaev, 2020) operate on very
small time scales (several local dynamics) and effectively result in magnetic
field diffusion in the favorable parts of the flow. Also, the strongly
magnetized layer buried under a less-magnetized layer of accreted matter is
susceptible to Parker instability (Payne and Melatos, 2007). The resulting
magnetic fields are suppressed, decayed, and have a much more complex topology
than the initial dipole.

##  6 Exotic and hypothetical stages

In this section, we briefly discuss several stages of magneto-rotational
evolution that remain hypothetical, i.e. no direct evidence in favor of their
existence is known.

###  6.1 Georotator

The name of this stage is related to the fact that the shape of the
magnetosphere and some processes are similar to the case of the Earth
magnetosphere in the Solar wind. The stage with
$R_{\mathrm{m}}<R_{\mathrm{co}}$ but $R_{\mathrm{m}}>R_{\mathrm{G}}$ was
mentioned already in Illarionov and Sunyaev (1975). The authors just note that
in this case matter flows around the magnetosphere as the gravitational
influence of the compact object is too small. Such a situation can be realized
if the velocity of matter relative to the NS is too high (e.g., due to very
fast stellar wind in a wide binary, or due to a large spatial velocity of the
compact object in the case of an INS), if the magnetic field is too large, or
if matter density is too low (this can happen e.g., for an INS far above the
Galactic plane).

It is worth to note, that the condition $R_{\mathrm{m}}<R_{\mathrm{co}}$ is
important. Correspondingly, the georotator stage is an alternative to the
stage of accretion, not to the propeller stage. Thus, one can consider
evolution as Ejector$\rightarrow$Propeller$\rightarrow$Georotator, or
switching from accretor to georotator. The latter scenario can be realized
e.g., while an INS moves in the Galactic potential rising high above the
Galactic plane and then coming back.

Also, we note that the magnetospheric radius is defined by different equations
in the case of accretion, at the propellers stage, and in the case when
$R_{\mathrm{m}}>R_{\mathrm{G}}$. That’s why the critical condition for
transition to/from the georotator stage might be carefully calculated
accounting for necessary effects. If we consider a typical INS with
$\mu\approx 10^{30}$ G cm3 and the ISM density $n\approx 1$ cm-3, then the
critical velocity for transition to the georotator stage is $\sim 400-500$ km
s-1. Note, that for such parameters duration of the ejector plus propeller
stages is longer than the Galactic age. Still, e.g. INSs with large magnetic
fields can become georotators.

To our knowledge, spin-down at the georotator stage was never analyzed in
detail. From general considerations, it can have some similarities with the
subsonic propeller stage Davies and Pringle (1981) or/and with the settling
accretion Shakura et al. (2012). When a long spin period is reached, the
influence of turbulence on the rotational frequency might become important
Prokhorov et al. (2002); Popov et al. (2001), similar to the case of accreting
INSs.

The case $R_{\mathrm{m}}>R_{\mathrm{G}}$ was studied numerically for non-
rotating INSs in Toropina et al. (2001). The authors dub this stage a
‘magnetic plow’. In this study, just the case of aligned magnetic axis and
spatial velocity vector was analyzed. It is demonstrated that a tiny amount of
matter still can accrete onto the NS surface, but the energy release is small
and hardly can be detected. Instead, the authors focus on the energy
dissipation in the bow shock and in the magnetic tail.

In the bow shock, the energy release can be estimated as:

$\dot{E}_{\mathrm{shock}}\approx\frac{\pi}{2}R_{\mathrm{m}}^{2}\rho v^{3}\sim
10^{21}B_{12}^{2/3}n^{2/3}\frac{v}{200\,\mathrm{km}\,\mathrm{s}^{-1}}\,\mathrm{erg}\,\mathrm{s}^{-1}.$
(61)

It is expected that the energy is emitted in X-ray and, maybe, in softer
bands.

The case of energy dissipation in the magnetic tail is more interesting as
different regimes can be realized. If the tail is not expanding, i.e. even far
from the NS its cross-section has the size $\sim R_{\mathrm{m}}$, then the
total energy is:

$E_{\mathrm{tot}}\sim
10^{27}B_{12}n^{1/2}\frac{v}{200\,\mathrm{km}\,\mathrm{s}^{-1}}\frac{l}{100\,R_{\mathrm{m}}}\,\mathrm{erg}.$
(62)

Here $l$ is the length of the tail. This energy can be dissipated due to
reconnection on the dynamical time scale

$t_{\mathrm{dyn}}=l/v\sim
10^{6}B_{12}^{1/3}n^{-1/6}\frac{v}{200\,\mathrm{km}\,\mathrm{s}^{-1}}\frac{l}{100\,R_{\mathrm{m}}}\,\mathrm{s}.$

Then

$\dot{E}_{\mathrm{rec}}=E_{\mathrm{tot}}/t_{\mathrm{dyn}}\sim
10^{21}B_{12}^{2/3}n^{2/3}\left(\frac{v}{200\,\mathrm{km}\,\mathrm{s}^{-1}}\right)^{7/3}\,\mathrm{erg}\,\mathrm{s}^{-1}.$
(63)

Energy can be also dissipated in bursts in a volume $\sim\pi
R_{\mathrm{m}}^{3}$, so the emitted amount is: