Convective dynamos of black widow companions

If evaporation by the pulsar’s radiation is too weak to significantly reduce the companion’s mass (Ginzburg & Quataert, 2020), then the companion must instead lose mass through Roche-lobe overflow in order to reach the low masses observed in black widow systems. We follow GQ21 and assume that the loss of orbital angular momentum always keeps the companion close to filling its Roche lobe and dictates stable mass loss (we calculate the rate in Section 2.3). Modelling of black widow optical light curves seems to support this picture for a large portion of the systems (Draghis et al., 2019; Kandel & Romani, 2023; Mata Sánchez et al., 2023).

The radius $r$ of a Roche-lobe filling companion is given by

where $M\gg m$ is the pulsar’s mass and $a$ is the binary separation (Eggleton, 1983). The orbital period $P$ of the Keplerian orbit is therefore

where $G$ is the gravitational constant. Our $M\gg m$ approximation (typically $M\gtrsim 1.4\,{\rm M}_{\odot}$; see Kandel & Romani, 2023) and Roche-lobe filling assumption allow us to evolve the companion in mesa as a single star, with the orbital period given at any moment by its mean density $\sim mr^{-3}$ using equation (2).

Irradiation of the companion by the pulsar’s high-energy photons is implemented in mesa using the $F_{\star}$–$\Sigma_{\star}$ method (Paxton et al., 2013) by depositing a flux

at an outer mass column density $\Sigma_{\star}$, assuming that heat is efficiently redistributed between the tidally locked companion’s day and night hemispheres (this assumption holds for high-energy photons that are deposited deep enough in the atmosphere; see GQ21, ). We choose $\Sigma_{\star}=180\textrm{ g cm}^{-2}$, given by the pair production cross section, which dominates the interaction of the pulsar’s energetic GeV photons with the companion’s atmosphere (Bethe & Heitler, 1934).

The pulsar’s high-energy photon luminosity $L_{\rm irr}$ is one of the free parameters of our model. Fermi measured $L_{\rm irr}$ for a subset of spider pulsars, finding that it is a significant – but not constant – fraction of the pulsar’s total spin-down power (Abdo et al., 2013). This subset has been used to test black widow evolutionary models (GQ21; Swihart et al. 2022).

Gravitational waves carry orbital angular momentum from the binary on a time-scale of

where $c$ is the speed of light (Landau & Lifshitz, 1971). We take $M=1.4\,{\rm M}_{\odot}$ for consistency with the ATNF pulsar catalogue (Manchester et al., 2005), even though many spider pulsars reach masses $M\approx 2\,{\rm M}_{\odot}$ (Linares, 2020; Kandel & Romani, 2023).

Gravitational waves carry significant angular momentum only from the shortest period black widows (e.g. fig. 1 in GQ21, ). Another sink of angular momentum is the companion’s wind, which is forced to corotate with its magnetic field up to a large radius – such that even a low-mass wind can carry significant angular momentum and spin down the star (Weber & Davis, 1967). Black widow companions are close enough to their host pulsars for their spin to be tidally locked to the orbit. Therefore, this magnetic braking mechanism ultimately decreases the binary’s orbital angular momentum.

Previous studies incorporated an empirical magnetic braking prescription that was calibrated to solitary main-sequence stars (Verbunt & Zwaan, 1981; Rappaport et al., 1983). El-Badry et al. (2022) recently demonstrated that this prescription cannot be extrapolated to short-period binaries such as CVs. More crucially, magnetic braking in black widows and redbacks is driven by pulsar ablation, which dominates over the companion’s spontaneously emitted wind. For this reason, we adopt a model that was developed from first principles to black widows (Ginzburg & Quataert, 2020). Using a hydrodynamical calculation to evaluate the strength of the ablated wind, the time-scale on which magnetic braking removes orbital angular momentum from a Roche-lobe filling companion is given by

where $B$ is the companion’s magnetic field (Ginzburg & Quataert 2020; GQ21). As explained in Ginzburg & Quataert (2020), the ablated wind is launched primarily by photons close to the electron’s rest mass ($\sim$ MeV). These photons deposit energy through Compton scattering high enough in the companion’s atmosphere, where cooling is inefficient and heating accelerates the material beyond the escape velocity. We denote the pulsar’s MeV luminosity by $L_{\rm MeV}$, which is typically unconstrained by observations (Fermi was sensitive to more energetic GeV $\gamma$ rays from pulsars; see Abdo et al., 2013). We limit this free parameter to $L_{\rm MeV}\lesssim L_{\rm irr}$, such that it is capped by the broader spectrum total irradiation.

In stable Roche-lobe overflow, the donor star loses its mass on the same time-scale that it loses its angular momentum, up to an uncertain order-unity coefficient which we omit (this depends on the details of the mass transfer; see Rappaport et al., 1982). The companion’s mass-loss rate in our model is therefore given by

with $t_{\rm GW}$ and $t_{\rm mag}$ calculated in equations (4) and (5).

To complete our computational model, we need to evaluate the companion’s magnetic field $B$, which together with the irradiation parameters $L_{\rm irr}$ and $L_{\rm MeV}$ uniquely determines the system’s evolution. The main improvement of this work compared to GQ21 is a more physically motivated model for $B$, which we describe in this section.

Convective dynamos are thought to generate magnetic fields by a combination of rotation and convection in a conducting fluid. However, the scaling laws that determine the strength of these fields are still in question (Christensen, 2010; Augustson et al., 2019). We follow GQ21 and assume that the magnetic energy density is in equipartition with the kinetic energy density of convective eddies

where $\rho$ is the fluid density and $v_{\rm conv}$ is the convective velocity. This scaling is thought to describe the fields of fast rotating convective dynamos with eddy turnover times that are much longer than the spin period (Christensen et al., 2009); this is the regime of tidally locked spider companion stars with periods of several hours. For simplicity, and due to the uncertainties of the dynamo theory, we omit potential order-unity coefficients from equation (7). Equation (5) indicates that our magnetic braking prescription is sensitive only to the combination $L_{\rm MeV}B^{3}$, such that these order-unity coefficients are absorbed by the parameter $L_{\rm MeV}$. The convection carries an energy flux of

such that the magnetic field scales as (Christensen et al., 2009)

In fully convective single stars and brown dwarfs, the convective flux can be related to the luminosity that has to be transported outwards $L\sim r^{2}F_{\rm conv}$, yielding a simple expression $B\sim m^{1/6}r^{-7/6}L^{1/3}$, where we used the mean density $\rho\sim mr^{-3}$ (Reiners et al., 2009). The luminosity of black widow companions, on the other hand, is dominated by atmospheric re-emission of the pulsar’s irradiating flux $L\approx F_{\star}\pi r^{2}$, with only a much smaller internal luminosity $L_{\rm int}\ll L$ delivered by convection from the companion’s interior. GQ21 accounted for this effect by scaling $B\sim m^{1/6}r^{-7/6}L_{\rm int}^{1/3}$.

Here we take a more direct approach. Instead of scaling $B$ with global quantities ($m$, $r$, and $L_{\rm int}$), we evaluate $B$ locally at each point using equation (7) with the local $\rho$ and $v_{\rm conv}$ (as computed by mesa). We then average $B$ over the convective envelope

assuming that convection effectively transports $B$. Averaging by radius rather than by mass does not change the results significantly; in any case the field is dominated by the bottom of the convection zone (perhaps it is more physically motivated to average the energy density by volume, yet all averages yield similar results in our case).

Our method has two major advantages compared to GQ21. First of all, at high values of $L_{\rm irr}$, the internal luminosity $L_{\rm int}$ is several orders of magnitude lower than the total luminosity $L_{\rm int}\ll L\approx L_{\rm irr}(r/2a)^{2}$ (see GQ21), making it numerically challenging to compute accurately. Equation (10) bypasses this technical difficulty by using $v_{\rm conv}$, which is directly computed by mesa from the net flux through each cell. More crucially, as we show below, strong irradiation inhibits convection and pushes the tachocline outwards, thinning the convective envelope. This effect was missed by GQ21 who assumed that convection always penetrates deep inside the companion star, leading to an overestimate of $B$ for high values of $L_{\rm irr}$. Equation (10) corrects this by integrating only over the convective zone. The new method will thus enable us to adequately describe the dynamos of various black widow and redback companions, which range from fully convective to almost fully radiative structures, depending on $L_{\rm irr}$.

Fig. 1 shows evolutionary tracks of initially $1\,{\rm M}_{\odot}$ main-sequence companions that are evolved according to our model (Section 2) with the mesa version r22.11.1. The black widow and redback observations are taken from the ATNF Pulsar Catalogue http://www.atnf.csiro.au/research/pulsar/psrcat (Manchester et al., 2005), version 1.70 (May 2023). Specifically, we plot all pulsars with spin periods below 30 ms and companions on orbits shorter than 1 d that are classified as either ‘ultra light’ or ‘main sequence’ in the catalogue (this removes white dwarf and neutron star companions from the sample).¹¹1We added to our sample the black widows J1518+0204C (Pallanca et al., 2014) and J1653$-$0158 (Nieder et al., 2020), which were miss-classified as helium white dwarfs in the catalogue. The circular grey markers indicate median companion masses and the error bars encompass 90 per cent of the possible inclinations (i.e. a minimal inclination angle of $25.8^{\circ}$).

Unsurprisingly, we find the same range of pulsar irradiation luminosities $L_{\rm irr}\approx 0.1-3\,{\rm L}_{\odot}$ as GQ21. These luminosities are required in order to lengthen the Kelvin–Helmholtz cooling and contraction time to several Gyr such that companions remain inflated and fill their Roche lobes at the observed orbital periods. Our fitted values of $L_{\rm MeV}$, on the other hand, are somewhat different from GQ21 because of our more careful calculation of the magnetic field $B$. We nominally fit a constant fraction $f\equiv L_{\rm MeV}/L_{\rm irr}=0.2$ in order for the evolutionary tracks to reach the observed mass cutoff at $m\approx 10^{-2}\,\rm{M}_{\odot}$ in about 10 Gyr; we interpret this cutoff as due to a maximum system age (the two extremely low-mass companions below the cutoff, J2322$-$2650 and J1719$-$1438 with $m\approx 10^{-3}\,{\rm M}_{\odot}$, are considered outliers here). However, at the highest values of $L_{\rm irr}$, the magnetic field weakens as we explain below, forcing us to fit a higher $f=0.6$ for $L_{\rm irr}=2\,{\rm L}_{\odot}$ in order to reach the observed mass cutoff. For $L_{\rm irr}=3\,{\rm L}_{\odot}$ even the maximal $f=1$ is not enough – the star evolves too slow – unable to reproduce the longest period observed black widows.

Fig. 2 shows our estimate for the companion’s magnetic field $B$, assuming equipartition between the magnetic and kinetic energies in the convective envelope (see Section 2.4). For most of the tracks, our results are in reasonable agreement with GQ21 – fields of order $B\sim 10^{2}-10^{3}\textrm{ G}$ that weaken with rising $L_{\rm irr}$, because it reduces the internal luminosity $L_{\rm int}$ and with it the convective flux. Both here and in GQ21 the fields converge to similar values at low companion masses (which comprise the majority of spider pulsars). Interestingly, Wadiasingh et al. (2018) inferred $\sim\,$kG fields from a pressure balance between the companion’s magnetosphere and the pulsar wind, which is required to explain the observed intrabinary shock (see also Polzin et al., 2018; Miao et al., 2023, who derive somewhat weaker fields).

For extremely strong pulsar irradiation $L_{\rm irr}\gtrsim 3\,{\rm L}_{\odot}$, however, our magnetic fields fall by an order of magnitude below GQ21. The reason for this can be found in Fig. 3, where we show the extent (by mass $\Delta m$ and by radius $\Delta r$) of the companion’s convective envelope. The internal luminosity $L_{\rm int}$ is so low in this case that it can be evacuated from the star almost without convection – the star becomes almost fully radiative, with only a thin convective envelope. The density $\rho$ of this outer shell is much lower than the stellar interior density, leading to a reduction in $B$ for a given convective flux $F_{\rm conv}$, as seen in equation (9). This effect was missed by GQ21 who assumed that $\rho$ is close to the companion’s mean density and thus overestimated $B$. Such high values of $L_{\rm irr}$ have been directly measured for several spider pulsars (Abdo et al., 2013; Swihart et al., 2022), necessitating our improved treatment of the convective dynamo.

In addition to its lower density, a geometrically thin convective shell presumably generates a smaller-scale magnetic field, with dominant multipole orders $\sim r/\Delta r\gg 1$, rather than a dipole. It has been argued that such small-scale fields are much less efficient in removing angular momentum (Garraffo et al., 2015, 2016; Réville et al., 2015; Garraffo et al., 2018), potentially weakening magnetic braking for high $L_{\rm irr}$ even beyond our model. We do not take this additional effect into account here; our highest $L_{\rm irr}$ models evolve too slow to become black widows in any case.²²2Sarkar & Tout (2022) also weaken magnetic braking for thin envelopes using their $\gamma$ parameter, but the physical motivation is different.

Fig. 3 demonstrates that spider pulsars can be used to benchmark magnetic dynamo and magnetic braking theories. The structure of their companions changes non-trivially as they lose mass, depending on the value of $L_{\rm irr}$ (which can be mapped to their orbital period, see Fig. 1). Short-period spider companions ($L_{\rm irr}\lesssim 0.5\,{\rm L}_{\odot}$) resemble main-sequence stars and become fully convective at $m\approx 0.3\,{\rm M}_{\odot}$. Companions orbiting at intermediate periods ($L_{\rm irr}\approx{\rm L}_{\odot}$) emit less internal heat because of the irradiation boundary condition and retain their radiative cores down to much lower masses $m<0.1\,{\rm M}_{\odot}$, as already noticed by GQ21. At the longest orbital periods ($L_{\rm irr}\gtrsim 2\,{\rm L}_{\odot}$), the radiative core even penetrates outwards as companions lose mass (at least initially), leaving only a thin convective envelope.

This variety of spider companion structures, ranging from fully convective to almost fully radiative, depending on their mass and orbital period, may be utilized to test the relative roles of the convective dynamo and of the tachocline in generating magnetic fields (e.g. Sarkar & Tout, 2022). Specifically, the companion mass distribution is a direct indicator of the magnetic braking rate, which dictates the mass-loss timescale. In the future, we plan to constrain and calibrate different magnetic braking prescriptions by testing them against this observed distribution. Here, as a first step, we limited ourselves to fitting the specific model of Christensen et al. (2009) and GQ21 in order to demonstrate the potential of these exotic binary systems.

Fig. 1 shows that there are several observations that cannot be explained by Roche-lobe overflow of main-sequence stars. In this section we demonstrate that most of these observations can be reproduced by the well-known bifurcation of the evolutionary trajectories of evolved donor stars (Pylyser & Savonije, 1988, 1989; Podsiadlowski et al., 2002). Specifically, for stars that overflow their Roche lobe at a time $t_{0}$ that is close to the end of their main-sequence lifetime $t_{\rm ms}$ there is an interesting interplay between the mass loss and the nuclear burning. As $t_{0}/t_{\rm ms}$ increases, donors enrich their interiors in helium before losing their envelopes. These dense helium-rich companions evolve to short orbital periods. Above some critical $t_{0}\approx t_{\rm ms}$, however, the donor develops a degenerate helium core, while the envelope expands to become a red giant. These low-density giants evolve to long orbital periods before their envelopes are lost, leaving behind a detached white dwarf (i.e. the degenerate helium core) orbiting the pulsar.