Low Mass X-ray Binaries: The Effects of the Magnetic Braking Prescription

TL;DR

This study uses detailed stellar evolution models to analyze how different magnetic braking laws affect the evolution and observed properties of low-mass X-ray binaries, revealing that standard prescriptions underestimate mass transfer rates.

Contribution

It demonstrates that the commonly used Skumanich magnetic braking law is insufficient for modeling LMXBs and proposes modified laws that better match observations.

Findings

Standard Skumanich MB law underestimates mass transfer rates by an order of magnitude.

Modified MB laws incorporating magnetic field scaling and wind mass loss better reproduce observed LMXB populations.

The results challenge the use of traditional MB prescriptions in binary evolution modeling.

Abstract

We present a population study of low- and intermediate-mass X-ray binaries (LMXBs) with neutron star accretors, performed using the detailed 1D stellar evolution code MESA. We identify all plausible Roche-lobe overflowing binaries at the start of mass transfer, and compare our theoretical mass transfer tracks to the population of well-studied Milky Way LMXBs. The mass transferring evolution depends on the accepted magnetic braking (MB) law for angular momentum loss. The most common MB prescription ("Skumanich MB") originated from observations of the time-dependence of rotational braking of Sun-type stars, where the angular momentum loss rate depends on the donor mass $M_d$, donor radius $R_d$, and rotation rate $\Omega$, $\dot{J} \propto M_d R_d^{\gamma} \Omega^3$. The functional form of the Skumanich MB can be also obtained theoretically assuming a radial magnetic field, isotropic isothermal winds, and boosting of the magnetic field by rotation. Here we show that this simple form of the Skumanich MB law gives mass transfer rates an order of magnitude too weak to explain most observed persistent LMXBs. This failure suggests that the standard Skumanich MB law should not be employed to interpret Galactic, or extragalactic, LMXB populations, with either detailed stellar codes or rapid binary population synthesis codes. We investigate modifications for the MB law, and find that including a scaling of the magnetic field strength with the convective turnover time, and a scaling of MB with the wind mass loss rate, can reproduce persistent LMXBs, and does a better job at reproducing transient LMXBs.

Figures

Tables

Table 1. Table 1: The different scaling values used in Equation  20 for the various cases.

Case	$\beta$	$\xi$	$\alpha$
1 - Default Skumanich	0	0	0
2 - Convection Boosted	0	2	0
3 - Intermediate	0	2	1
4 - Wind Boosted	2	4	1

Table 2. Table 2: Sco X-1 Test Properties

Quantity	Observed	Boosted MB	Default MB	Ref.
Mass Ratio	0.28 - 0.51	0.42	0.28	(1,2)
Period [Hours]	18.89551	18.90	18.90	(3)
MT [$M_{\odot }{\mathrm{yr}}^{-1}$]	$2.2\times {10}^{-8}$	$2.1\times {10}^{-8}$	$1.1\times {10}^{-9}$	(3)
Effective Temperature [K]	<4800	4718	4627	(1,2)
Donor Mass [$M_{\odot }$]	0.28 - 0.70	0.56	0.50	(1,2)
NS Mass [$M_{\odot }$]		1.33	1.79	(1,2)

Table 3. Table 3: NS Binaries

Source	Flags	$M_d$	Mass Ratio	Period	Distance	Average Mass Transfer	Reference
		($M_{\odot }$)	$\left(M_d/M_a\right)$	(Hours)	(kpc)	($M_{\odot }$ yr-1)
4U 1820-303	Per, GC, UC	$-$	$-$	${0.183}_{\mathrm{S87}}$	$7.9\pm {0.4}_{\mathrm{H13}}$	$1.2\pm 0.6\times {10}_{\mathrm{H13}}^{-8}$	S87, H13
4U 0513-40	Per, GC, UC	$0.045$	${0.03}^{\ast }$	${0.283}_{\mathrm{Z09}}$	$12.1\pm {0.6}_{\mathrm{H13}}$	$1.2\pm 0.6\times {10}_{\mathrm{H13}}^{-9}$	Z09, H13
2S 0918-549	Per, UC	$0.024-{0.039}_{\mathrm{Z11}}$	$0.02-{0.03}^{\ast }$	${0.290}_{\mathrm{S01}}$	$5.4\pm {0.8}_{\mathrm{H13}}$	$2.6\pm 1.5\times {10}_{\mathrm{H13}}^{-10}$	Z11, H13
4U 1543-624	Per, UC	${0.03}_{\mathrm{W04}}$	${0.02}^{\ast }$	${0.303}_{\mathrm{W04}}$	$7_{\mathrm{H13}}$	$1.3\begin{array}{c}+1.8\\ -1.2\end{array}\times {10}_{\mathrm{H13}}^{-9}$	W04, H13
4U 1850-087	Per, GC, UC	${0.04}_{\mathrm{H96}}$	${0.03}^{\ast }$	${0.343}_{\mathrm{H96}}$	$6.9\pm {0.3}_{\mathrm{H13}}$	$2.2\pm 1.1\times {10}_{\mathrm{H13}}^{-10}$	H96, H13
M15 X-2	Per, GC, UC	$0.02-{0.03}_{\mathrm{D05b}}$	${0.02}^{\ast }$	${0.377}_{\mathrm{D05b}}$	$10.4\pm {0.5}_{\mathrm{H13}}$	$3.8\pm 1.9\times {10}_{\mathrm{H13}}^{-10}$	D05b, H13
IGR J17062-6143	T, PSR, UC	$0.0155-{0.0175}_{\mathrm{S18b}}$	${0.01}^{\ast }$	${0.633}_{\mathrm{S18b}}$	$7.3\pm {0.5}_{\mathrm{K17}}$	$2.5\times {10}_{\mathrm{K17}}^{-11}$	K17, S18b
XTE J1807-294	T, PSR, UC	$-$	$-$	${0.670}_{\mathrm{M03}}$	$8{\begin{array}{c}+4\\ -3.3\end{array}}_{\mathrm{H13}}$		M03, H13
4U 1626-67	Per, PSR, UC		${0.02}^{\ast }$	${0.700}_{\mathrm{C98a}}$	$8{\begin{array}{c}+5\\ -3\end{array}}_{\mathrm{H13}}$	$8.0\begin{array}{c}+14\\ -6\end{array}\times {10}_{\mathrm{H13}}^{-10}$	C98a, H13
XTE J1751-305	T, PSR, UC	$-$	$-$	${0.710}_{\mathrm{M02}}$	$8{\begin{array}{c}+0.5\\ -1.3\end{array}}_{\mathrm{H13}}$	$5.1\begin{array}{c}+2.6\\ -2.9\end{array}\times {10}_{\mathrm{H13}}^{-12}$	M02, H13
XTE J0929-314	T, PSR, UC	$\simeq {0.01}_{\mathrm{G02}}$	$0.007$	${0.730}_{\mathrm{G02}}$	$8{\begin{array}{c}+7\\ -3\end{array}}_{\mathrm{H13}}$		G02, H13
4U 1916-053	Per, UC	$0.064\pm {0.01}_{\mathrm{H13}}$	$0.046$	${0.833}_{\mathrm{W81}}$	$9.3\pm {1.4}_{\mathrm{H13}}$	$6.3\pm 3.7\times {10}_{\mathrm{H13}}^{-10}$	W81, H13
Swift J1756.9-2508	T, PSR, UC	$\le {0.022}_{\mathrm{K07}}$	$0.02$	${0.912}_{\mathrm{K07}}$	$8\pm 4_{\mathrm{H13}}$	$1.7\begin{array}{c}+2.3\\ -1.5\end{array}\times {10}_{\mathrm{H13}}^{-11}$	K07, H13
NGC 6440 X-2	T, PSR, GC, UC	$\simeq {0.0076}_{\mathrm{B15}}$	$0.005$	${0.955}_{\mathrm{A10}}$	$8.5\pm {0.4}_{\mathrm{H13}}$	$1.0\pm 0.5\times {10}_{\mathrm{H13}}^{-12}$	A10, B15, H13
HETE J1900.1-2455	T, PSR	$0.016-{0.07}_{\mathrm{K06a}}$	$0.01-{0.05}^{\ast }$	${1.39}_{\mathrm{W08}}$	$4.7\pm {0.6}_{\mathrm{W08}}$	$4.63\times {10}_{\mathrm{W08}}^{-11}$	K06a, W08
1A 1744-361	T	$0.07-{0.22}_{\mathrm{B06}}$	$0.07{?}^{\ast }$	$1.62\pm {0.37}_{\mathrm{W08}}$	$6\pm 3_{\mathrm{W08}}$	$2.22\times {10}_{\mathrm{W08}}^{-11}$	B06, W08
IGR J17379-3747	T, PSR	$>{0.056}_{\mathrm{S18a}}$	$>{0.04}^{\ast }$	${1.88}_{\mathrm{S18a}}$	$\sim {8.5}_{\mathrm{S18a}}$	$\sim 4\times {10}_{\mathrm{S18a}}^{-11}$	S18a, S18c
SAX J1808-3658	T, PSR	$0.04{\begin{array}{c}+0.02\\ -0.01\end{array}}_{\mathrm{W13}}$	${0.04}^{\ast }$	${2.01}_{\mathrm{C98b}}$	$3.4-{3.6}_{\mathrm{C12}}$	$1.73\times {10}_{\mathrm{C12}}^{-11}$	C98b, C12, W13
XB 1832-330	T,GC	$-$	$-$	${2.1}_{\mathrm{E12}}$	${10.0}_{\mathrm{P01}}$	$\sim 3\times {10}_{\mathrm{E12}}^{-10}$	P01, E12
IGR 00291+5934	T, PSR	$0.039-{0.16}_{\mathrm{D17}}$	$0.02-0.11$	${2.46}_{\mathrm{G06}}$	$2.6-{3.6}_{\mathrm{C12}}$	$1\times {10}_{\mathrm{D17}}^{-12}$	C12, G05, D17
4U 1822-00	Per, M	$-$	$-$	$3.20{?}_{\mathrm{S07}}$	$6.3\pm 2_{\mathrm{S07}}$	$9.0\begin{array}{c}+8.0\\ -5.0\end{array}\times {10}_{\mathrm{S07}}^{-10}$	S07
4U 1636-536	Per	$0.29-{0.48}^{\ast }$	$0.21-{0.34}_{\mathrm{W16}}$	${3.79}_{\mathrm{W08}}$	$6\pm {0.5}_{\mathrm{W08}}$	$1.25\times {10}_{\mathrm{C12}}^{-9}$	C12, W08, W16
EXO 0748-676	T	$0.1{?}_{\mathrm{D14a}}$	$0.07{?}^{\ast }$	${3.82}_{\mathrm{D14a}}$	$7.1\pm {1.2}_{\mathrm{D14a}}$		C12, D14a
4U 1254-69	Per	${0.45}_{\mathrm{C68b}}$	$0.33-{0.36}_{\mathrm{C13}}$	${3.93}_{\mathrm{W08}}$	$13\pm 3_{\mathrm{W08}}$	$1.77\times {10}_{\mathrm{W08}}^{-9}$	C86b, C13, W08
4U 1728-16 (GX 9+9)	Per	${0.4}_{\mathrm{K09b}}$	${0.29}_{\mathrm{K09b}}$	${4.20}_{\mathrm{L07}}$	$5{?}_{\mathrm{K06b}}$	$2.91\times {10}_{\mathrm{C97}}^{-9}$	C97, K06b, L07
XTE J1814-338	T	$0.19-{0.32}_{\mathrm{W17}}$	$0.123{\begin{array}{c}+0.012\\ -0.01\end{array}}_{\mathrm{W17}}$	${4.27}_{\mathrm{W08}}$	$8\pm {1.6}_{\mathrm{C12}}$		C12, W08, W17
4U 1735-444	Per	$\le 0.58$	$0.05-{0.41}^{C06}$	${4.65}_{\mathrm{W08}}$	$8.5\pm {1.3}_{\mathrm{W08}}$	$6.31\times {10}_{\mathrm{C12}}^{-9}$	C06, C12, W08
4U 1746-37	Per, GC	$--$	$-$	${5.16}_{\mathrm{B04}}$	${11.6}_{\mathrm{B04}}$	$1\times {10}_{\mathrm{B04}}^{-9}$	S01, B04
2A 1822-371	Per, M	$0.47\pm {0.04}_{\mathrm{I15}}$	${0.28}_{\mathrm{I15}}$	${5.57}_{\mathrm{B17}}$	${2.5}_{\mathrm{B17}}$	$\sim 2\times {10}_{\mathrm{B17}}^{-8}$	B17, I15
XTE J2123-058	T	$0.76\pm {0.22}_{\mathrm{C02}}$	$0.49\pm {0.1}_{\mathrm{C01}}$	${5.96}_{\mathrm{W08}}$	$9.6\pm {1.3}_{\mathrm{W08}}$		C02, C12, S03a, W08
X 1658-298	T	$0.3-{0.8}_{\mathrm{P18}}$	$0.21-{0.57}^{\ast }$	${7.12}_{\mathrm{D14a}}$	$12\pm 3_{\mathrm{W08}}$	$1\times {10}_{\mathrm{W08}}^{-9}$	D14a, P18, W08
2A 0521-720 (LMC X-2)	Per, M	$-$	$-$	${8.16}_{\mathrm{L07}}$	$50\pm 2_{\mathrm{A09}}$	$4\times {10}_{\mathrm{C12}}^{-8}$	C12, A09, L07
SAX J1748.9-2021	T, GC	$0.12-1_{\mathrm{S16}}$	$0.09-{0.71}^{\ast }$	${8.76}_{\mathrm{S16}}$	$8.5\pm {0.4}_{\mathrm{S16}}$	$\sim 7\times {10}_{\mathrm{W08}}^{-11}$	S16, W08
IGR J18245-2452	T, GC	$>{0.17}_{\mathrm{P13}}$	${0.12}^{\ast }$	${11.0}_{\mathrm{P13}}$	${5.5}_{\mathrm{P13}}$	$\lesssim 1.0\times {10}_{\mathrm{P13}}^{-10}$	P13
GRS 1747-312	T, GC	$--$	$-$	${12.36}_{\mathrm{I03}}$	${9.5}_{\mathrm{V18}}$	$1\times {10}_{\mathrm{V18}}^{-10}$	B04, V18
4U 1456-32 (Cen X-4)	T	$0.31\pm {0.27}_{\mathrm{D05a}}$	$0.18\pm {0.06}_{\mathrm{D05a}}$	${15.1}_{\mathrm{L07}}$	$1.2\pm {0.2}_{\mathrm{C12}}$	$4\times {10}_{\mathrm{C12}}^{-11}$	C12, D05a, L07
AC 211 (X2127+12)	Per, GC	$\sim {0.1}_{\mathrm{V04}}$	$\sim {0.1}_{\mathrm{V04}}$	${17.1}_{\mathrm{I93}}$	${10.4}_{\mathrm{C68b}}$	$\sim 7\times {10}_{\mathrm{I93}}^{-9}$	C86, I93, V04
H 1617-155 (Sco X-1)	Per, M	$0.28-{0.70}_{\mathrm{S15}}$	$0.28-{0.51}_{\mathrm{S15}}$	${18.9}_{\mathrm{W08}}$	$2.8\pm {0.3}_{\mathrm{S15}}$	$3\times {10}_{\mathrm{C12}}^{-8}$	C12, S15, W08
4U 1908+005 (Aql X-1)	T	$-$	$-$	${18.9}_{\mathrm{W08}}$	$4.55\pm {1.35}_{\mathrm{W08}}$	$6\times {10}_{\mathrm{C12}}^{-10}$	C12, W08
4U 1624-49	Per	$-$	$-$	${20.9}_{\mathrm{B00}}$	$15{\begin{array}{c}+2.9\\ -2.6\end{array}}_{\mathrm{X09}}$	$4.6\times {10}_{\mathrm{B09}}^{-9}$	B09, L05, X09
3A 1702-363 (GX 349+2)	Per	${0.78}^{\ast }$	$\sim {0.56}_{\mathrm{I09}}$	$21.9\pm {0.4}_{\mathrm{I09}}$	$5\pm {1.5}_{\mathrm{C12}}$	$2.37\pm 1\times {10}_{\mathrm{C12}}^{-8}$	C12, I09, W08
2A 1655+353 (Her X-1)	T	$2.03\pm 0.37$	$1.45$	${40.8}_{\mathrm{L07}}$	$6.1{\begin{array}{c}+0.9\\ -0.4\end{array}}_{\mathrm{L14}}$	$1.3\times {10}_{\mathrm{C12}}^{-8}$	C12, L07, R11
4U 2142+38 (Cyg X-2)	Per, M	$0.56\pm {0.07}_{\mathrm{P16}}$	$0.34\pm {0.01}_{\mathrm{P16}}$	${236.3}_{\mathrm{W08}}$	$10.55\pm {4.45}_{\mathrm{W08}}$	$3.0\times {10}_{\mathrm{C12}}^{-8}$	C12, M18, W08
GRO J1744-28	T, M	$0.2-{0.7}_{\mathrm{D15}}$	$0.15-{0.5}^{\ast }$	${284.0}_{\mathrm{L07}}$	$8{?}_{\mathrm{D15}}$	$\sim 1\times {10}_{\mathrm{D14b}}^{-8}$	D15, D14b, L07

Table 4. Table 4: Properties of Persistent LMXBs

System Name	${\log }_{10}\left(P\right)$	q	${\log }_{10}\left({\dot{M}}_a\right)$
Ultra-compact XRBs
4U 0513-40	[-0.57, -0.52]	[0.01, 0.06]	[-9.2, -8.6]
2S 0918-549	[-0.56, -0.51]	[0.01, 0.06]	[-9.8, -8.6]
4U 1543-624	[-0.54, -0.49]	[0.01, 0.06]	[-9.1, -8.6]
4U 1850-087	[-0.48, -0.43]	[0.01, 0.06]	[-10.0, -8.4]
M15 X-2	[-0.44, -0.39]	[0.01, 0.06]	[-9.7, -9.1]
4U 1626-67	[-0.17, -0.12]	[0.01, 0.06]	[-9.7, -8.6]
4U 1916-053	[-0.10, -0.05]	[0.03, 0.08]	[-9.6, -8.9]
4U 1636-536	[0.56, 0.61]	[0.15, 0.40]	[-9.1, -8.6]
GX 9+9	[0.60, 0.65]	[0.20, 0.33]	[-8.7, -8.2]
4U 1735-444	[0.65, 0.70]	[0.29, 0.48]	[-8.4, -7.9]
2A 1822-371	[0.73, 0.78]	[0.26, 0.36]	[-7.8, -7.3]
Sco X-1	[1.26, 1.31]	[0.15, 0.58]	[-7.8, -7.3]
GX 349+2	[1.33, 1.38]	[0.39, 0.65]	[-8.0, -7.3]
Cyg X-2	[2.35, 2.40]	[0.25, 0.53]	[-7.7, -7.2]

Table 5. Table 5: Properties of Transient LMXBs

System Name	${\log }_{10}\left(P\right)$	q	${\log }_{10}\left({\dot{M}}_a\right)$
HETE J1900.1-2455	[0.12, 0.17]	[0.01, 0.06]	[-10.5, -10.0]
1A 1744-361∗	[0.19, 0.24]	[0.05, 0.18]	[-11.7, -10.7]
SAX J1808-3658	[0.28, 0.33]	[0.02, 0.07]	[-11.0, -10.5]
IGR 00291+5394	[0.37, 0.42]	[0.02, 0.13]	[-11.8, -11.3]
EXO 0748-676∗	[0.56, 0.61]	[0.05, 0.10]	[-9.3, -8.3]
4U 1254-69	[0.57, 0.62]	[0.23, 0.38]	[-9.0, -8.5]
XTE J1814-338∗	[0.61, 0.66]	[0.10, 0.27]	[-11.2, -10.2]
XTE J2123-058∗	[0.76, 0.81]	[0.27, 0.82]	[-11.2, -10.2]
X 1658-298	[0.83, 0.88]	[0.15, 0.67]	[-9.1, -8.6]
SAX J1748.9-2021	[0.92, 0.97]	[0.06, 0.83]	[-10.3, -9.8]
IGR J18245-2452	[1.02, 1.07]	[0.09, 0.14]	[-10.2, -9.7]
Cen X-4	[1.16, 1.21]	[0.02, 0.48]	[-10.6, -10.1]
Her X-1	[1.59, 1.64]	[0.83, 2.00]	[-8.1, -7.6]
GRO J1744-28	[2.43, 2.48]	[0.25, 0.53]	[-8.2, -7.7]

Table 6. Table 6: Maximum Lifetime of an Observed LMXB and Fraction of Parameter Space

	Default		Convection Boosted		Intermediate		Wind Boost
System Name	${\tau }_{\max }$ [years]	$A_{\mathrm{sys}}/A_{\mathrm{tot}}$	${\tau }_{\max }$ [years]	$A_{\mathrm{sys}}/A_{\mathrm{tot}}$	${\tau }_{\max }$ [years]	$A_{\mathrm{sys}}/A_{\mathrm{tot}}$	${\tau }_{\max }$ [years]	$A_{\mathrm{sys}}/A_{\mathrm{tot}}$
4U 0513-40	$3.23\times {10}^6$	$4.75\times {10}^{-4}$	$1.24\times {10}^6$	$7.91\times {10}^{-4}$	$4.00\times {10}^6$	$1.90\times {10}^{-3}$	$0$	$0$
2S 0918-549	$4.09\times {10}^6$	$4.75\times {10}^{-4}$	$2.30\times {10}^6$	$7.91\times {10}^{-4}$	$4.40\times {10}^6$	$2.06\times {10}^{-3}$	$0$	$0$
4U 1543-624	$4.83\times {10}^6$	$6.33\times {10}^{-4}$	$4.80\times {10}^6$	$6.33\times {10}^{-4}$	$7.91\times {10}^6$	$2.06\times {10}^{-3}$	$0$	$0$
4U 1850-087	$8.33\times {10}^6$	$6.33\times {10}^{-4}$	$8.83\times {10}^6$	$1.11\times {10}^{-3}$	$1.28\times {10}^7$	$2.06\times {10}^{-3}$	$0$	$0$
M15 X-2	$1.21\times {10}^7$	$6.33\times {10}^{-4}$	$1.11\times {10}^7$	$1.27\times {10}^{-3}$	$1.41\times {10}^7$	$2.22\times {10}^{-3}$	$0$	$0$
4U 1626-67	$6.76\times {10}^7$	$9.50\times {10}^{-4}$	$1.06\times {10}^8$	$1.58\times {10}^{-3}$	$7.48\times {10}^7$	$6.33\times {10}^{-4}$	$0$	$0$
4U 1916-053	$1.85\times {10}^7$	$1.58\times {10}^{-4}$	$1.21\times {10}^7$	$3.17\times {10}^{-4}$	$3.76\times {10}^7$	$4.75\times {10}^{-4}$	$0$	$0$
4U 1636-536	$1.22\times {10}^8$	$4.23\times {10}^{-2}$	$6.24\times {10}^7$	$1.91\times {10}^{-2}$	$4.39\times {10}^7$	$7.74\times {10}^{-2}$	$1.22\times {10}^5$	$2.06\times {10}^{-3}$
GX 9+9	$0$	$0$	$3.38\times {10}^7$	$1.60\times {10}^{-2}$	$3.73\times {10}^7$	$5.89\times {10}^{-2}$	$1.85\times {10}^4$	$1.27\times {10}^{-3}$
4U 1735-444	$0$	$0$	$1.48\times {10}^7$	$5.70\times {10}^{-3}$	$9.83\times {10}^6$	$2.26\times {10}^{-2}$	$0$	$0$
2A 1822-371	$0$	$0$	$0$	$0$	$5.26\times {10}^6$	$3.39\times {10}^{-2}$	$0$	$0$
Sco X-1	$0$	$0$	$1.69\times {10}^7$	$3.48\times {10}^{-3}$	$1.75\times {10}^6$	$3.91\times {10}^{-2}$	$2.71\times {10}^3$	$2.37\times {10}^{-3}$
GX 349+2	$4.09\times {10}^5$	$7.91\times {10}^{-4}$	$1.50\times {10}^7$	$6.17\times {10}^{-3}$	$1.91\times {10}^5$	$2.37\times {10}^{-3}$	$0$	$0$
Cyg X-2	$1.65\times {10}^6$	$8.39\times {10}^{-3}$	$2.56\times {10}^6$	$1.74\times {10}^{-2}$	$2.00\times {10}^5$	$1.28\times {10}^{-2}$	$0$	$0$
HETE J1900.1-2455	$5.61\times {10}^8$	$6.61\times {10}^{-2}$	$3.05\times {10}^8$	$4.57\times {10}^{-2}$	$5.87\times {10}^8$	$1.01\times {10}^{-1}$	$0$	$0$
1A 1744-361∗	$1.99\times {10}^6$	$1.58\times {10}^{-4}$	$3.51\times {10}^7$	$1.58\times {10}^{-4}$	$6.80\times {10}^7$	$7.91\times {10}^{-4}$	$0$	$0$
SAX J1808-3658	$5.49\times {10}^7$	$6.33\times {10}^{-4}$	$1.05\times {10}^8$	$6.33\times {10}^{-4}$	$1.66\times {10}^8$	$1.74\times {10}^{-3}$	$0$	$0$
IGR 00291+5394	$0$	$0$	$0$	$0$	$1.88\times {10}^8$	$4.75\times {10}^{-4}$	$0$	$0$
EXO 0748-676∗	$0$	$0$	$4.09\times {10}^7$	$2.53\times {10}^{-3}$	$0$	$0$	$0$	$0$
4U 1254-69	$5.23\times {10}^7$	$1.28\times {10}^{-2}$	$1.06\times {10}^6$	$7.91\times {10}^{-4}$	$1.70\times {10}^7$	$1.65\times {10}^{-2}$	$0$	$0$
XTE J1814-338∗	$0$	$0$	$0$	$0$	$0$	$0$	$2.60\times {10}^5$	$1.27\times {10}^{-3}$
XTE J2123-058∗	$0$	$0$	$0$	$0$	$1.40\times {10}^3$	$1.58\times {10}^{-3}$	$2.25\times {10}^2$	$2.53\times {10}^{-3}$
X 1658-298	$1.19\times {10}^8$	$1.30\times {10}^{-2}$	$2.85\times {10}^7$	$8.23\times {10}^{-3}$	$1.55\times {10}^7$	$4.91\times {10}^{-3}$	$2.26\times {10}^4$	$3.17\times {10}^{-3}$
SAX J1748.9-2021	$7.55\times {10}^8$	$6.33\times {10}^{-4}$	$3.08\times {10}^8$	$1.27\times {10}^{-3}$	$4.92\times {10}^7$	$1.39\times {10}^{-2}$	$8.96\times {10}^4$	$2.69\times {10}^{-3}$
IGR J18245-2452	$8.91\times {10}^8$	$6.33\times {10}^{-4}$	$0$	$0$	$3.82\times {10}^7$	$2.37\times {10}^{-3}$	$0$	$0$
Cen X-4	$0$	$0$	$1.14\times {10}^8$	$1.58\times {10}^{-4}$	$3.53\times {10}^7$	$1.58\times {10}^{-2}$	$4.32\times {10}^3$	$1.11\times {10}^{-3}$
Her X-1	$7.38\times {10}^6$	$1.27\times {10}^{-3}$	$9.97\times {10}^6$	$1.27\times {10}^{-3}$	$5.68\times {10}^5$	$1.42\times {10}^{-3}$	$4.75\times {10}^5$	$2.69\times {10}^{-3}$
GRO J1744-28	$5.85\times {10}^6$	$8.70\times {10}^{-3}$	$5.48\times {10}^6$	$1.31\times {10}^{-2}$	$6.98\times {10}^4$	$5.70\times {10}^{-3}$	$0$	$0$

Equations

$$\begin{gathered} {\dot{M}}_{\mathrm{d}}^{\mathrm{w}}=\eta \times 4\times 10^{-13}\frac{R}{R_{\odot }}\frac{L}{L_{\odot }}\frac{M_{\odot }}{M}{M}_{\odot }\mathrm{y}{\mathrm{r}}^{-1}, \\ \eta =1 \end{gathered}$$

$$\frac{R_L}{a}=\frac{0.49q^{2/3}}{0.6q^{2/3}+\ln (1+q^{1/3})}.$$

$${\dot{M}}_{\mathrm{Edd}}=\frac{4\pi c{R}_{\mathrm{a}}}{{\kappa }_e}\approx \frac{3.4}{1+X}\times 10^{-8}{M}_{\odot }{\mathrm{yr}}^{-1}$$

$$\frac{{\dot{J}}_{\mathrm{GR}}}{J}=-\frac{32}{5}{G}^3{c}^5\frac{M_{\mathrm{d}}{M}_{\mathrm{a}}(M_{\mathrm{d}}+M_{\mathrm{a}})}{a^4}$$

$${\dot{J}}_{\mathrm{MB},\mathrm{Sk}}=-3.8\times 10^{-30}{M}_{\mathrm{d}}{R}_{\odot }^4{\left(\frac{R_{\mathrm{d}}}{R_{\odot }}\right)}^{{\gamma }_{\mathrm{mb}}}{\Omega }^3\mathrm{dyne}\mathrm{cm}$$

$$\frac{1}{2}{\rho }_{\mathrm{A}}{v}_{\mathrm{A}}^2\simeq \frac{B(r{)}^2}{8\pi }.$$

$$B(r)=B_s\frac{R_s^2}{r^2}.$$

$$\begin{array}{rl}{\dot{J}}_{\mathrm{MB}}& =-4\pi \Omega {\int }_0^{\pi /2}{\rho }_{\mathrm{A}}{v}_{\mathrm{A}}{R}_{\mathrm{A}}^2(R_{\mathrm{A}}\sin \theta {)}^2\sin \theta d\theta \\ & \simeq -\frac{8}{3}\pi \Omega {\rho }_{\mathrm{A}}{v}_{\mathrm{A}}{R}_{\mathrm{A}}^4,\\ & \simeq -\frac{2}{3}\Omega {\dot{M}}_{\mathrm{W}}{R}_{\mathrm{A}}^2.\end{array}$$

$$\begin{array}{rl}{\dot{M}}_{\mathrm{W}}& =4\pi R^2{\rho }_s{v}_s,\\ & =4\pi R_{\mathrm{A}}^2{\rho }_{\mathrm{A}}{v}_{\mathrm{A}},\end{array}$$

$$4\pi R_A^4{\rho }_A{v}_A^2=B_s^2{R}^4,$$

$$\begin{array}{rl}{\dot{M}}_{\mathrm{W}}{R}_A^2{v}_A& =B_s^2{R}^4,\\ R_A^2& =\frac{B_s^2{R}^4}{{\dot{M}}_{\mathrm{W}}{v}_A}.\end{array}$$

$${\dot{J}}_{\mathrm{MB}}\propto \Omega B_s^2{R}^4.$$

$$v_{\mathrm{A}}^2=\frac{2GM}{R},$$

$${\dot{J}}_{\mathrm{MB}}\propto {\dot{M}}_{\mathrm{W}}\Omega B_s^4{R}^4.$$

$$\begin{array}{rl}{N}_d& \approx R_0^{-2},\\ N_d& \approx {\Omega }^2{\tau }_{\mathrm{conv}}^2.\end{array}$$

$${\tau }_{\mathrm{conv}}={\int }_R^{R_s}\frac{dr}{v_{\mathrm{conv}}}.$$

$$\begin{array}{rl}\frac{B_s}{B_{s,\odot }}& =\frac{R_0}{R_{0,\odot }},\\ \frac{B_s}{B_{s,\odot }}& =\frac{\Omega }{{\Omega }_{\odot }}\frac{{\tau }_{\mathrm{conv}}}{{\tau }_{\odot \mathrm{conv}}}.\end{array}$$

$$\begin{array}{rl}{\dot{J}}_{\mathrm{MB}}& \propto \Omega B_s^2{R}^4,\\ & \propto {\Omega }^3{\tau }_{\mathrm{conv}}^2{R}^4.\end{array}$$

$$\begin{array}{rl}{\dot{J}}_{\mathrm{MB}}& \propto {\dot{M}}_{\mathrm{W}}\Omega B_s^4{R}^4\\ & \propto {\dot{M}}_{\mathrm{W}}{\Omega }^5{\tau }_{\mathrm{conv}}^4{R}^4.\end{array}$$

$${\dot{J}}_{\mathrm{MB},\mathrm{boost}}={\dot{J}}_{\mathrm{MB},\mathrm{Sk}}{\left(\frac{\Omega }{{\Omega }_{\odot }}\right)}^{\beta }{\left(\frac{{\tau }_{\mathrm{conv}}}{{\tau }_{\odot ,\mathrm{conv}}}\right)}^{\xi }{\left(\frac{{\dot{M}}_{\mathrm{W}}}{{\dot{M}}_{\odot ,W}}\right)}^{\alpha }.$$

$${\dot{M}}_{\mathrm{crit}}=k{P}_{\mathrm{hr}}^b\mathrm{g}{\mathrm{s}}^{-1}$$

$$\dot{M}=\frac{L_{\mathrm{X}}{R}_{\mathrm{a}}}{G{M}_{\mathrm{a}}}.$$

Full text

††pubyear: 2018††pagerange: Low Mass X-ray Binaries: The Effects of the Magnetic Braking Prescription–Low Mass X-ray Binaries: The Effects of the Magnetic Braking Prescription

1 Introduction

Low mass X-ray binaries (LMXBs) with a neutron star (NS) accretor are among the most well-studied binary systems in astrophysics. Over 100 of them have been observationally identified in the Milky Way over the last 50 years, while a plethora of binary parameters – mass ratios, orbital periods, mass transfer (MT) rates, and donors’ effective temperatures in some cases – have been determined for several of them. For many years, X-ray binaries have posed as an enticing problem for theorists, providing grounds to verify ways to model binary stellar evolution (Faulkner, 1971; Rappaport et al., 1983; Webbink et al., 1983; Joss & Rappaport, 1984; Podsiadlowski et al., 2002). The applications of theoretical models of LMXBs are not limited to only our Galaxy, but have also been applied to interpret the observed populations of X-ray binaries in other galaxies (e.g., Fragos et al., 2008; Tremmel et al., 2013), and even to estimate their feedback on reionization of the early universe (e.g., Fragos et al., 2013)

Understanding and interpreting LMXBs can be split into two distinct stages: the formation and life of a binary system with a NS before the start of the MT, and the evolution of the binary system during the MT. Before the onset of the MT, the evolution of the binary is expected to proceed through a a common envelope phase and a supernova explosion (Bhattacharya & van den Heuvel, 1991). Supernova natal kicks, supernova explosion mechanisms, and widely used pulsar kick distributions, while observationally derived, are not fully understood (Fryer et al., 2012; Hobbs et al., 2005). The outcomes of common envelope events are highly uncertain as well (for a review, see Ivanova et al., 2013).

In this paper, we focus on the evolution of LMXBs during MT. The driving mechanism of the MT phase is the donor’s response to the mass and angular momentum loss. In short period systems, the dominant channel of angular momentum loss is gravitational radiation, which is well understood theoretically and has been confirmed by observations (Weisberg & Huang, 2016; LIGO Scientific Collaboration & VIRGO Collaboration, 2016; LIGO Scientific Collaboration & Virgo Collaboration, 2017). In cases with longer orbital periods, magnetic braking (MB) is the dominant angular momentum loss mechanism (Rappaport et al., 1983).The strength of MB affects the evolution of LMXBs by increasing or decreasing the angular momentum loss of the binary. In cases where the angular momentum loss is greater, the binary will tend to shrink and thus undergo more aggressive MT. The default assumption for MB is to adopt the "Skumanich" empirical law where the angular momentum loss scales with donor mass $M_{d}$, donor radius $R_{d}$, and rotation rate $\Omega$, $\dot{J}\propto M_{d}R_{d}^{\gamma}\Omega^{3}$ (Skumanich, 1972). In theoretical models, this empirical law is used by employing the parameterized prescriptions for MB stated by either Verbunt & Zwaan (1981) or Rappaport et al. (1983). Numerical studies of LMXB populations have shown systematic mismatches between simulated results and observations. Observed MT rates have been found differing from the predicted rates in simulated systems by up to an order of magnitude (Podsiadlowski et al., 2002). Despite obtaining more detailed data on LMXBs, the discrepancies remained prevalent. For instance, Sco X-1 cannot be explained by the "Skumanich" law. To match the observations of Sco X-1, the rate of angular momentum loss has to be boosted, for example by taking into account the effect of stellar wind loss Pavlovskii & Ivanova (2015).

We are taking this analysis a step further by including the effects of the convective turnover time, as well as considering a non-isothermal stellar wind § 2.2. We apply different MB laws to a grid of 2136 binary systems spanning reasonable initial conditions, for a circularized binary system with a NS, where MT starts within 10 Gyr § 2.1. We describe the principal properties of the simulated MT systems in § 4 and have compiled a table with updated properties of some observed NS LMXBs (see § 3). The observational data is then used to constrain the MB laws in § 5.

2 Numerical Method

2.1 Setting up and evolving the population of binaries

We consider the set of initial binary systems as follows:

• •

The masses of donors range from $M_{\rm d}=1.0{\ M_{\odot}}$ to $M_{\rm d}=7{\ M_{\odot}}$. The grid in donor masses uses steps of $1{\ M_{\odot}}$ for donors with masses $\geq 5{\ M_{\odot}}$, $0.5{\ M_{\odot}}$ for donors between $3{\ M_{\odot}}$ and $5{\ M_{\odot}}$, 0.2 ${\ M_{\odot}}$ for donors between $2.4{\ M_{\odot}}$ and $3{\ M_{\odot}}$, and $0.1{\ M_{\odot}}$ for donors with masses $\leq 2.4{\ M_{\odot}}$. This is the mass of the donors at ZAMS. The adopted metallicity is $Z=0.02$.

• •

The initial binary orbital periods range from the periods at which the donor stars would overflow their Roche lobes at ZAMS, to the maximum orbital period at which they would start the MT while they are on the red giant branch. We define the initial period as the period at ZAMS of the donor star; which is not the same as the period that a binary would have at the start of the MT. There can be a large difference between when the initial period at the donor stars ZAMS and at the onset of RLOF. The mesh for the initial orbital periods has a fixed step of 0.05 in $\log_{10}(P)$. All orbits are circular.

• •

With the initial mesh of seed masses and periods, we ran 2136 simulations for each considered MB scheme (see § 2.2).

• •

The compact companions are NSs with mass $M_{\rm a}=1.4{\ M_{\odot}}$ and radius $R_{\rm a}=11.5$ km.

All calculations for the single stars and mass transferring binaries were performed using the one-dimensional stellar evolution code MESA 111http://mesa.sourceforge.net (Modules for Experiments in Stellar Astrophysics), revision 10398 and the August 2018 release of the MESASDK 222http://www.astro.wisc.edu/~townsend/static.php?ref=mesasdk. MESA is a modern open-source set of stellar libraries as described in Paxton et al. (2011); Paxton et al. (2013, 2015, 2018). The donors are evolved using default assumptions as in MESA. In particular, we use the mixing length $\alpha=2$, no semiconvection, and no overshooting. MESA uses mixing length theory as described in Cox & Giuli (1968). We adopt a grey atmosphere boundary condition and use the OPAL opacity tables for solar composition (Grevesse & Noels, 1993)333The inlists files can be found at the MESA marketplace website: http://cococubed.asu.edu/mesa_market/.

The donor stars are evolved using Reimer’s wind mass loss prescription (Reimers, 1975):

[TABLE]

$\eta$ is a scaling factor or efficiency of wind loss. In our calculations we use $\eta=1$.

When the donor star overfills its Roche lobe we calculate MT; we do not consider any other forms of MT. For the Roche lobe radius, we use the fitting formula for the volume equivalent one-dimensional star’s radius, as provided by Eggleton (1983):

[TABLE]

Here $R_{L}$ is the Roche lobe radius of the donor star with mass $M_{\rm d}$, $a$ is the orbital separation between the two stars, and $q=M_{\rm d}/M_{\rm a}$ is the mass ratio of the two stars. To calculate the MT via the Lagrange point between the two stars $L_{1}$, $\dot{M}_{\rm d}^{L1}$, we employ the "Ritter" mass loss scheme, as implemented in MESA (see Ritter, 1988, for details of this MT prescription). We consider that the MT may be non conservative. If $\dot{M}_{\rm d}^{L1}$ exceeds the Eddington-limited MT rate $\dot{M}_{\rm Edd}$, the excess $\dot{M}_{\rm d}^{L1}-\dot{M}_{\rm Edd}$ cannot be accreted on to the compact object, and is assumed to be lost from the system. For the Eddington limited MT rate on a NS we use,

[TABLE]

Here $\kappa_{e}$ is the opacity due to Thomson electron scattering, $\kappa_{e}=0.2(X+1)$ cm2 g*-1*, where X is the hydrogen mass fraction in the material transferred from the donor.

The angular momentum of the system is lost through gravitational radiation, or through MB, or is carried away with the mass lost from the system. The mass lost from the donor due to wind mass loss leaves with the specific angular momentum of the donor. Note the orbital evolution is calculated before the start of the MT as well. If the MT rate exceeds $\dot{M}_{\rm Edd}$, the material exceeding the Eddington limit is lost with the specific angular momentum of the accretor. The angular momentum loss due to gravitational radiation is found using the standard formula (Faulkner, 1971):

[TABLE]

2.2 Magnetic braking

We calculate the MT sequences considering several prescriptions for the angular momentum loss via MB. The first prescription considered uses the formulation from Rappaport et al. (1983):

[TABLE]

Here, $R_{\rm d}$ is the radius of the donor, $\gamma_{\rm mb}$ is a dimensionless parameter from 0 to 4 and $\Omega$ is the angular velocity of the donor. The stars are kept in corotation with the binary as the angular velocity is calculated using the binary period. With $\gamma_{\rm mb}=4$, Equation 5 describes the standard Skumanich law as derived by Verbunt & Zwaan (1981), and is the most commonly used form for MB in calculations of LMXB evolution.

The Skumanich law was scaled to describe main sequence stars similar to our Sun. In systems with donors different from the Sun, the increased rate of mass loss with the stellar wind, as well as magnetic field strength that does not scale directly with the angular velocity of the star, will play a role in the MB calculation (Mestel, 1968; Mestel & Spruit, 1987; Kawaler, 1988). To determine the effects of these additional terms, we follow similar steps as Pavlovskii & Ivanova (2015). We start with the formulation given by Mestel & Spruit (1987) which parameterizes the amount of angular momentum lost from the system through mass leaving through the Alfvèn surface. The Alfvèn surface is the surface where the ram pressure is equal to the magnetic pressure (Mestel, 1968; Mestel & Spruit, 1987):

[TABLE]

$\rho_{\rm A}$ and $v_{\rm A}$ denote the density and velocity of the wind as it crosses the Alfvèn surface. $B(r)$ is the poloidal magnetic field strength as a function of radius. This value encompasses the structure of the magnetic field of the star. In the simplest case where the magnetic field is radial, and $B_{s}$ is the surface magnetic field strength, then:

[TABLE]

Throughout this work, we will only be using a radial magnetic field, but it should be noted that other magnetic field structures are possible. In the context of magnetic braking however, the Alfvèn surface represents the maximum radius where the stellar wind is locked in corotation with the surface of the star, beyond this point the mass is assumed to be lost (Mestel, 1968; Mestel & Spruit, 1987). The equation which describes the angular momentum loss through an Alfvèn surface is:

[TABLE]

Here it is assumed that the Alfvèn surface, $R_{\rm A}$ does not depend on $\theta$, the polar angle, and that the wind coming from the star is isotropic. Should the system be rapidly rotating the scaling of the magnetic braking to the rotation rate will change as so-called dead zones may form (see discussion in Ivanova, 2006). It is important to note the Mestel & Spruit (1987) parameterization can account for “dead zones” where material is trapped in magnetic fields lines and not lost through the Alfvèn surface. The material confined within the dead zone remain in corotation within the dipole field and is not lost from the system. We don’t include the effects of dead zones in this work and as such may over estimate the amount of material lost and thus the angular momentum loss in tight binaries. The wind mass loss rate, if isotropic can be described using the following equation:

[TABLE]

Combining equations 6 and 7 gives us:

[TABLE]

Including the wind mass loss equation as given in equation 9, equation 10 becomes:

[TABLE]

Under the assumption of isothermal winds, the wind reaches a sonic wind velocity at the Alfvèn surface of $v_{\rm A}=c_{w}$ (Mestel & Spruit, 1987) where $c_{w}$ is a constant value. Combining equations 8 and 11 gives the following MB scaling equation:

[TABLE]

Equation 12 interestingly, does not contain any scaling with the stellar wind, despite the wind strength being a fundamental physical property of MB. The assumption of isothermal winds does not apply to giant stars. In the case of giant stars with cooler temperatures, winds may be accelerated by a variety of mechanisms, this requires a different self-consistent description of the wind velocity (Suzuki, 2007). In the case of a radial field and a nonthermal wind where the wind moves at speeds on order of the escape velocity,

[TABLE]

We get the following scaling:

[TABLE]

The strength of the surface magnetic field, $B_{s}$, scales with the dynamo number $N_{D}$ (Parker, 1971). The dynamo number is related to physical values in the MB through the Rossby number, $R_{0}$ (Noyes et al., 1984),

[TABLE]

Here $\tau_{\rm conv}$ is the turnover time of convective eddies,

[TABLE]

$R$ and $R_{s}$ are the bottom and the top of the outer convective zone respectively, while $v_{\rm conv}$ is the local convective velocity. We follow a simple approximation made by Ivanova (2006) where $B_{s}\propto N_{d}^{1/2}$. This allows us to adopt the scaling relations:

[TABLE]

Work by Aurière et al. (2015) has shown a general correlation between the semi-empirical Rossby number and the observed magnetic field strength of a star. Aurière et al. (2015) also noted that dwarf stars have a steeper relation between these properties than giant stars. In general, dwarf stars have been found to have shorter rotation periods with strong magnetic fields. As such, the relation used likely underestimates the strength of the magnetic field in dwarf stars. Rewriting the magnetic field scaling in Equation 12 gives us the following:

[TABLE]

The radial isothermal approximation results in the Skumanich scaling with $\Omega^{3}R^{4}$ if we ignore the convective turnover time $\tau_{\rm conv}$. The radial non-thermal approximation from Equation 14 becomes

[TABLE]

Rewriting the Skumanich law to include the additional terms for wind $\dot{M}_{\rm W}$, convective turnover time $\tau_{\rm conv}$ and rotation rate $\Omega$, the general MB equation we use will be

[TABLE]

The value used to normalize the convective turnover time, $\tau_{\odot,\rm conv}=2.8\times 10^{6}\ \rm s$, was found by evolving a $1{\ M_{\odot}}$ star at $Z$=0.020 to 4.6 Gyrs. Similarly, the solar wind value $\dot{M}_{\odot,\rm W}=2.5\times 10^{-14}M_{\odot}\rm\ yr^{-1}$ (Carroll & Ostlie, 2006), and $\Omega_{\odot}\approx 3\times 10^{-6}\ \rm s^{-1}$ is the angular frequency of the Sun using an orbital period of 24 days.

The power $\xi$ can vary, where $\xi=0$ describes the same simplified assumptions for which the Skumanich law is valid with $\alpha=0$ and $\beta=0$ (i.e., radial magnetic field and isothermal winds). $\xi=2$ is the case described in equation 18 which results in the convection boosted Skumanich case. $\xi$ may be as high as 4 for the case of winds from giants where the velocity of the wind grows linearly with distance; we note that in this case, the dependence on the angular velocity will also have to be modified to $\Omega^{5}$, vs. the Skumanich law’s factor of $\Omega^{3}$. Therefore in this case, $\xi=4$, $\beta=2$ and $\alpha=1$.

We will use these additional scaling terms and define different MB cases for the tested grid of binaries:

“Default”: We use the default MB scheme described by Rappaport et al. (1983), without the additions mentioned in Pavlovskii & Ivanova (2015). $\gamma_{\rm mb}=4$ in this case, and all subsequent cases. 2. 2.

“Convection-boosted”: We adopt the scaling found in Equation 18, which is the Skumanich law, scaled by the convective turnover time $(\tau_{\rm conv})^{\xi}$. The value of $\xi=2$ will be used in this prescription. 3. 3.

“Intermediate”: We use the convection-boosted MB scheme and apply an additional wind scaling term, linear in wind mass loss rate ($\alpha=1$). 4. 4.

“Wind-boosted”: This MB scheme uses the scaling values from Equation 19. This prescription includes all three scaling terms shown in Equation 20 with $\beta=2$, $\xi=4$ and $\alpha=1$.

These systems are evolved to 10 Gyrs, or until the donor star loses its envelope and detaches. If the simulation encounters dynamically unstable MT, which MESA was not designed to adequately model, the system will likely encounter numerical issues and stop. We do not consider irradiation effects on the companion star during its evolution (see §6).

2.3 Verification against a previous study

Of particular importance is the binary system Scorpius X-1 as this system allows us to compare our results to the work of Pavlovskii & Ivanova (2016). In Pavlovskii & Ivanova (2016), the authors used a modified MB prescription where $\alpha=1$, $\beta=0$, and $\xi=0$, with an MB gamma $\gamma_{\rm mb}$ of 3. In this work, we tested all MB prescriptions described in §2.1 and used an MB gamma $\gamma_{\rm mb}$ of 4. To ensure the changes to MB were correctly implemented, the comparisons between models used in Pavlovskii & Ivanova (2016) were rerun.

Sco X-1 is an LMXB, observed to have a mass ratio in the range from 0.28 - 0.51, favouring a value of $\approx$ 0.30 (Mata Sánchez et al., 2015; Steeghs & Casares, 2002). The NS is constrained to have a mass of $<$ 1.73 ${\ M_{\odot}}$ (Mata Sánchez et al., 2015). The period of Sco X-1 is 18.8951 hours, and the MT rate is estimated to be at least $\sim 2.2\times 10^{-8}~{}{\ M_{\odot}}yr^{-1}$ (Watts et al., 2008; Pavlovskii & Ivanova, 2016). Observations of this system provided upper limits on the spectral class of the donor of K4 or later, with the luminosity class IV, and the implied effective temperature less than 4800 K.

The models tested in Pavlovskii & Ivanova (2016) were composed of a donor at ZAMS with masses ranging from 0.9 to 1.8${\ M_{\odot}}$, and a NS varying from 1.24 to 1.6${\ M_{\odot}}$. They assumed solar metallicity and Reimer’s wind prescription. That study, to find the mass transferred the $L_{1}$ Lagrange point, $\rm\dot{M}_{\rm d}^{\it L1}$, used the method described in Pavlovskii & Ivanova (2015), while we use the "Ritter" prescription. The method of Pavlovskii & Ivanova (2015) is important for determining the initial MT stability in systems with a very high mass ratio, and the use of the "Ritter" prescription should not play a role in the test of Sco X-1, or in finding long-lived LMXBs. Pavlovskii & Ivanova (2016) considered the case of standard MB and wind-boosted MB, adopting $\gamma_{\rm mb}=3$ for both. This is the default value for $\rm\gamma_{mb}$ in MESA444We note this value was chosen to be default by the MESA core developers groups, to make a test comparison to results published in the past, and is not motivated by physics; we remain convinced that the standard Skumanich law should be used with $\gamma_{\rm mb}=4$.. Pavlovskii & Ivanova (2016) find the default prescription of MB gives insufficient mass transfer to reproduce the observed mass transfer rate of Sco X-1, by at least an order of magnitude. To produce the observed properties of Sco X-1 (within $\approx 10\%$ estimated uncertainty), a modified MB law must instead be used.

In this test only, for comparison purposes, we have run similar MT sequences with $\gamma_{\rm mb}=3$. We considered one of the sets of binary companion masses presented in Pavlovskii & Ivanova (2016): a 1.3${\ M_{\odot}}$ NS and a 1.0 ${\ M_{\odot}}$ donor. The test run for the default MB was done using an initial period of 2.7 days, while for the test of the wind-boosted MB we used an initial period of 7.6 days. Both initial periods were taken from systems in Pavlovskii & Ivanova (2016). The binary systems were evolved until they had similar orbital periods as Sco X-1. The MT tracks of the two systems are shown in Figure 1. They are similar to those shown in Figures 1 and 2 of Pavlovskii & Ivanova (2016). The results of the simulations with the two MB prescriptions as well as the observed values are listed in Table 2.

As can be seen, we confirm that the modified MB prescription better reproduces the observed value of the MT rate in Sco X-1. We will return to the case of Sco X-1 in section 5 to review which MB prescriptions can reproduce Sco X-1.

3 Observational Data for NS LXMBs

To compare our results to observations, we have compiled a list of up-to-date properties of some Milky Way NS LMXBs. See Table 3 for a list of systems, relevant data, and references.

The systems with periods shorter than 80 minutes are ultra-compact X-ray binaries (UCXBs), where the donor star must be partly or completely degenerate. Our binary evolution method may not necessarily be the dominant method to produce these systems, but we keep them in the consideration. Within the table there are also systems labelled as GC systems, which are systems found in globular clusters. Binaries formed within a globular cluster are not likely to be produced from primordial binaries, but instead are more effectively produced via dynamical encounters between binaries, as well as via physical collisions between NSs and subgiants (Verbunt, 1987; Ivanova et al., 2005, 2008).

The other flags shown in the second column denote if the source is a persistent or transient system. Persistent systems are those where bright ($L_{X}>10^{35}$ erg/s) X-ray emissions have been consistently seen whenever X-ray monitoring missions have observed these systems, over a 40-60 year timespan. Transient systems have large changes in their X-ray emission, typically exceeding $L_{X}=10^{36}$ erg/s at some points (outbursts) and declining below $L_{X}=10^{35}$ erg/s at other times (quiescence). There are several systems currently thought to be persistent which could instead be in a long-term outburst state (lasting $>$50 years), and may be reclassified as transient systems in the future.

Meyer & Meyer-Hofmeister (1981) predicted the existence of a critical MT rate separating the persistent and transient systems. The disc instability model (DIM) predicts under what circumstances an accretion disk will experience instabilities (see Lasota, 2001, for a review of DIM). Stability in the context of the DIM means an accretion disc does not experience outbursts; a stable disc that remains hot will produce a persistent X-ray binary. The criterion for stability is given by Coriat et al. (2012):

[TABLE]

Here, $\rm P_{\rm hr}$ is the period of the system in hours. For the non-irradiated case, $k=(2.6\pm 0.9)\times 10^{16}$, and $b=1.76$. If there is irradiation of the accretion disc, and the accretor is a NS, $b=1.59$ and $k=(2.9\pm 0.9)\times 10^{15}$. Systems that lie above the DIM line are expected to be persistent while systems below the line are transient.

Here, we only consider the systems where the MT rate is known. If the cited papers do not provide the MT rate but instead provide the X-ray luminosity, $L_{x}$, we find $\dot{M}$ using

[TABLE]

Observationally, an upper limit for NS mass has been found to be $2.01\pm 0.04M_{\odot}$ (Antoniadis et al., 2013), our calculations will be done assuming the mass of the NS is $M_{\rm a}=1.4~{}M_{\odot}$ with a radius of $R_{\rm a}=11.5$ km (Özel & Freire, 2016). The key properties of interest for this work are the mass ratio $q=M_{d}/M_{a}$, the period $p$ and the average mass transfer rate $\dot{M}$.

4 Principal Results

4.1 Evolutionary Tracks

Here we present the results for the 2136 binary models of each MB prescription by plotting the evolutionary tracks grouped by donor mass on subplots in figures 2 - 5. These figures show the donor mass and period evolution of each simulated binary over the course of its lifetime. Mass in these figures is meant as a proxy for time as the donor stars all decrease in mass over the course of the binary evolution. The binary population simulations each required on order of hours to days to finish, totalling approximately thirty core years of simulation time. For comparison purposes, we overlay the observed data points from Table 3 on the plots. Sco X-1 is denoted with a star in the subplots as it was our test case from 2.3.

We can see that our choice of initial conditions affects the resulting evolutionary track of the system. In general, the more massive the donor star, the higher the MT rate. The choice of initial period results in changes that are less monotonic in comparison. Many systems starting with short periods may initiate Roche lobe overflow (RLOF) immediately, while longer period systems may lose significant mass through winds before this.

The separation between systems which evolve to longer orbital periods during their evolution, and those in which the orbital periods shrink, is known as the bifurcation period (Podsiadlowski et al., 2002; van der Sluys et al., 2005). As the initial donor mass increases, the bifurcation period also increases until we encounter systems undergoing dynamically unstable MT. We will refer to a binary with the orbital period near the end of its evolution larger than the period at the start of RLOF as growing in period, while systems with a shorter period near the end of their evolution are shrinking. One can see that most of the observed LMXBs are located in the region of shrinking systems. Specifically, in the case using the default Skumanich law, several of them are in the parameter space where the shrinking and growing systems bifurcate.

4.1.1 Default

The default Skumanich MB results are shown in Figure 2. This MB prescription produces a clear bifurcation in periods through the low and intermediate mass donors. The bifurcation period at RLOF is $\sim 1$ day and shows a distinct split in periods as no binaries pass through the region $q<0.4$, $1\lesssim\log_{10}(P)\lesssim 2.5$. The sharp transition seen in the low mass high period systems is a result of the binary reaching a stopping conditions given in Section 2.1.

4.1.2 Convection Boosted

From Figure 3, we cannot see many differences between the default MB and the “convection boosted" MB. One difference between the default and convection boosted MB schemes is that evolutionary tracks run through the region $q<0.4$, $1\lesssim\log_{10}(P)\lesssim 2.5$. The additional systems passing through this region cause some ambiguity in determining bifurcation periods in systems with initial donor masses $M\lesssim 1.5M_{\odot}$, since the binaries near the bifurcation period show very little change in period over their evolution.

4.1.3 Intermediate

The “Intermediate" case includes the additional scaling factor which accounts for the effects of wind mass loss as seen in equation 20. Figure 4 shows the behaviour of the binary systems with this MB prescription. The additional wind scaling plays a significant role in wider binary systems. The stronger MB scheme brings the binary systems together on a shorter timescale. In these systems, the MB and total angular momentum loss are consistently an order of magnitude higher than the default case. Therefore, the system loses enough angular momentum due to MB that gravitational radiation comes into play once MB stops.

4.1.4 Wind Boosted

The “wind boost" case shown in figure 5, includes the effects of the convective turnover time and accounts for the rotation rate of the star. In this case, the individual values of magnetic field strength as calculated using equation 17, the turnover time, and the wind mass loss rates are all within reasonable ranges. The magnetic field, which reaches a maximum of 100G shown in Figure 6, is within the range expected for giant stars (Aurière et al., 2015). The convective turnover time is also similar to those calculated by Pavlovskii & Ivanova (2015) for systems that are predicted to reproduce Sco X-1. It appears that the individual properties are all within reasonable ranges. However, the combination of all these boost factors produces MB that is too strong, resulting in MT that consistently exceeds $1M_{\odot}\ \rm yr^{-1}$. MT at these rates results in evolution on a dynamical timescale, and as such the results from this highly boosted MB should not be trusted. It is likely that this prescription has reached and exceeded a saturation limit that is not accounted for in this work (Mestel & Spruit, 1987). It has been shown that in cases where the rotation rate is high, additional magnetic field structure effects must be included to dampen the angular momentum loss (Ivanova & Taam, 2003). For completeness, we include the results from the wind-boosted case, but since the MT rate is so high, it is unlikely these simulations accurately describe reality.

4.2 Binary Properties

Figure 7 shows evolutionary tracks in the donor mass-period plane with the colour of each point representing the neutron star mass at that point. If a system experiences MT below the Eddington limit, the accretor mass can grow significantly. The significant increase in mass is a common outcome for default Skumanich MB, and the convection boosted cases, where most of the LMXBs evolve to contain a NS more massive than $2M_{\odot}$. A “stronger MB”, such as our intermediate case, produces fewer systems where NS masses have increased significantly, due to the MT exceeding the Eddington limit for portions of the evolution. The lack of substantial accretion of material onto the NS is even more apparent in the wind-boosted case where there is no significant increase in mass despite the large decrease in donor mass. However, the wind-boosted systems generally die too quickly – most simulations don’t reach the relevant M, P range for the majority of observed binaries. From observations there is a lack of NS detected near $2.0M_{\odot}$ with the most massive detected at $2.01\pm 0.04M_{\odot}$ (Antoniadis et al., 2013). The results shown in the default and convection boosted case suggest that all short period binaries contain a NS with a mass exceeding $2.0M_{\odot}$, whereas observations have found that these systems contain NS accretors in the range of $1.4M_{\odot}$. This preliminary result strongly supports the stronger MB prescriptions over the weaker convection boost and default cases.

Measurement of the surface chemical composition is possible in select binary systems. In binary systems involving white dwarf accretors, the surface chemical abundances can constrain the possible formation channels, with high $\rm N/C>10$ or $\rm N/O\gtrsim 10$ implying a helium donor (Nelemans et al., 2010). Uncertainties in model spectra for UCXBs result in unreliable abundance ratios (Werner et al., 2006). The existence of strong C and O lines but weak He and N lines imply a helium donor star. However, in many cases more detailed observations are necessary to classify the possible donors (Nelemans et al., 2010). Figures 8 and 9 show the C/N and N/O surface ratios from the simulated systems, respectively. In our donors, the initial C:N:O ratio is 0.37:0.108:1.0.

This change of abundance indicates where the donor star is in its evolution. In systems with a lower-mass donor, $M\sim 1M_{\odot}$, the donors show high N abundances as the donor star is stripped of C and O. In binaries with more massive donors, $M\gtrsim 1.5M_{\odot}$, the CN and CNO cycles can result in substantial changes in chemical composition. These significant changes can be observed in the material that is being transferred in these compact binaries.

4.3 Relative Densities

Figure 10 shows the data plotted in the period-MT plane with the colour bar representing a normalized frequency. Here we can see which systems are more or less likely to appear in each respective bin. We calculate frequency as follows:

$\tau^{mn}_{\rm tot}$ is the total evolution time of a binary system, as defined by the initial $m$-period and $n$-mass. 2. 2.

$\tau_{ij}^{mn}$ is the time that a system defined by the initial $m$-period and $n$-mass spends in an $i$-period bin and a $j$-MT rate bin. 3. 3.

$f_{ij}^{mn}=\tau_{ij}^{mn}/\tau^{mn}_{\rm tot}$ is the fraction of time (or frequency) that a particular system (given by an initial $mn$) appears in a bin defined by a particular $i$-period and $j$-MT rate. 4. 4.

$f_{ij}=$ max $(f_{ij}^{mn})$ the frequency plotted in Fig. 10, is the maximum frequency of all systems to appear in a particular $i$,$j$ bin.

The frequency gives an indication of where the evolutionary tracks spend the maximum amount of time during the evolution in the period-MT plane. This frequency however, does not necessarily represent a likelihood of detecting a binary in said bin as we have neglected effects of how likely a binary is to form and the number of systems which may cross through a bin. This is done to avoid equating regions where many systems cross into a bin for short periods of time to regions where one system spends a large fraction of its lifetime. For example, using this method, a bin where one binary spend $10^{7}$ years out of its total lifetime within the parameter space results in a much higher frequency than a bin where 1000 binaries each spend $10^{4}$ years.

We see the densest regions of the figure are at short periods, with low MT rates. These systems correspond to the binaries that have shrinking periods over the course of their evolution, including ultra-compact binaries. In general, the higher the MT rate, the lower the frequency; this is not surprising as it is difficult to maintain high MT rates. In the default MB, there is a gap in the simulated density between $0.5\lesssim\log_{10}(\rm P/hr)\lesssim 1.0$ and $\log_{10}(\dot{M}/{\ M_{\odot}}\rm yr^{-1})\lesssim-11$. While there are no observed systems in this range, this gap begins to get populated once convection is accounted for and is filled in with the "intermediate" prescription. In general, as the MB boost is increased, the MT rate at lower periods is also increased. The "convection boost", "intermediate" and "wind boost" cases appear to more effectively reproduce the systems with periods $\log_{10}(\rm P/hr)\sim 1$ and MTs near the Eddington limit. These systems include GX 9+9, 4U 1735+444, and 2A 1822-371, which do not overlap with any simulated systems in figure 10. Sco X-1, on the other hand, appears to be reproducible in figure 10, as there is significant overlap with the simulated systems in the period-MT plane. We argue, however, that to reproduce an observed LMXB, in addition to the period and MT rate, we must also match the mass ratio of the system.

5 Comparison to the observed population of LMXBs

We describe a set of binary systems that have an observationally determined MT rate, orbital period, and mass ratio (see Table  3, not all observed systems can be used). We bin each binary system within a range of period, MT and mass ratio, with the observed values used as the central bin values. These bins are then used to analyze the systems, and the adopted range for $P$, $\dot{M}$ and $q$ are described in Tables 4 and 5. The “observational” bin sizes are large enough to accommodate the anticipated observational errors.

The likelihood of a given MB scheme being correct depends on how effectively it can reproduce the observed binary systems. To find this value, we check if a simulated MT system passes through any of the “observational” bins. If a simulated system passes through one of the observational bins, we can find the total time that the system spent $\tau_{ij}^{n}$ in that “observational” bin. Here, $n$ stands for the number of the “observational” bin, and $ij$ describes the initial orbital period and the initial mass of the donor. We then find the size of the initial parameter space in the initial orbital periods and the donor masses, the parameter space from which the systems could evolve through the particular “observational” bin $n$.

In Table 6 we provide the maximum amount of time, $\tau_{\rm max}^{n}$ that any of the simulated systems can spend in the $n$ observational bin of interest, and the fraction of the initial parameter space that can produce the observed systems. These results are separated by MB prescription, with the wind-boosted case producing clearly ineffective results. With the high MT rate, the simulated parameter space overlaps with very few observed data points in 5 and 10. This high boost rate is likely invalid, and including a dampening factor should yield more realistic results.

The most striking result is that the default Skumanich MB law cannot reproduce most of the persistent systems - in fact, no observed persistent LMXBs with orbital periods between about 4.5 and 23 hours can be produced once we account for mass ratio. The main reason is that MB is not strong enough to drive the observed MT rate. In principle, in addition to convection or wind-boosted MB laws, one can invoke also another alternative MB law to explain systems like Sco X-1 - as was done, e.g., by Chen (2017), who considered initial donors of $1.5-2.5M_{\odot}$ star and applied the MB law derived from Justham et al. (2006). The best fitting progenitor systems from Chen (2017) have $1.6-1.8M_{\odot}$ donors with a 300G fossil magnetic field. It would, however, be rather intriguing if most of the observed persistent systems must be descendants of low-mass A stars with magnetic fields about an order of magnitude weaker than that of Ap stars, but 100 times larger than that of regular stars.

A significant result is that while we can create UCXBs, the initial parameter space to form these systems, and the lifetimes of the systems in those data bins, suggest that binary systems with a NS accretor and a non-perturbed donor are unlikely to be the main progenitors. They are instead likely to be produced either in globular clusters via physical collisions of a NS with a red giant Ivanova et al. (2005), or as a result of common evolution in the field, where the MT can start either from a cooled-down stripped core or, conversely, from a hot stripped core (Heinke et al., 2013). While the UCXBs have low MT rates, the short time spent in the observed bin is due to the simulated binary having an MT rate near the boundary of the bin.

As the MB strength increases, the number of persistent systems that can be reproduced increases. The binaries that could not be reproduced by the default MB prescription, GX 9+9, 4U 1735-444 and Sco X-1, are those that Podsiadlowski et al. (2002) found had MT rates much higher than their simulations reached. The convection-boosted MB is necessary to reproduce these persistent systems. Similarly, the available evidence indicates that the MT rate of 2A 1822-371 is super-Eddington (Bak Nielsen et al., 2017). For 2A 1822-371, the convection-boosted case is still insufficient to reach the high MT observed. The wind-boosted case pushes the MT rate high enough to reproduce this system at an appropriate period. The mass ratio however, does not match with the observed binary.

The reproducibility of a transient binary is affected by the choice of MB prescription. Unfortunately, unlike persistent systems, a clear trend isn’t apparent from one prescription to another. Instead, different systems are reproduced by different MB prescriptions, and with the uncertainties in determining an average MT in these systems, we cannot use the transient binaries to reliably draw any conclusions without a clear trend. One result that can be seen from the reproducibility of the transient systems in Table 6 is that the systems that are most difficult to reproduce are those where the MT rate is only constrained by an upper limit. The intermediate case reproduces the largest number of observed transient systems with only EXO 0748-676 and XTE J1814-338 not being reproduced.

The wind-boosted MB simulations cannot reproduce the majority of the observed binaries regardless of whether they are persistent or transient. The wind-boost prescription gives very short lifetimes for all reproduced binaries (which makes their detection unlikely), and cannot reproduce the UCXBs. This suggests that the simulated wind-boosted case is exceeding some saturation point for MB and the systems are losing too much angular momentum too quickly.

5.1 Other effects

Irradiation may play a significant role in driving winds from the donor star (Ruderman et al., 1989), and in causing the donor star to expand to a larger radius than expected for its orbital period (Podsiadlowski, 1991). Such an increase in donor radius due to irradiation may cause cycles of increased MT rates (up to a factor of $\sim$30), followed by a decrease to below-average MT rates (e.g. Hameury et al., 1993). The detailed physics of such irradiation-driven MT cycles have not yet been established, but current work suggests that these cycles should require small convective timescales in the donor, and thus may operate on systems with periods between 4 and 15 hours (Büning & Ritter, 2004). Irradiation-induced MT cycles could potentially produce the large observed MT rates in some of our transient and persistent systems in this period range, but should not be relevant for longer-period systems such as Sco X-1. Whether irradiation-induced MT can play a significant role depends on as-yet-undetermined details of the heating efficiency of the irradiating flux, and the fraction of time at the increased MT rate.

We have difficulty reproducing two transient systems (EXO 0748-676 and XTE J1814-338), both with orbital period near 4 hours, with any of our MB schemes. Our best MB scheme, the "intermediate" case, predicts MT rates higher than observed for these 2 systems. It is possible that irradiation-driven MT cycles might alter these systems’ evolution enough to match their observed characteristics (although such MT cycles are thought to cut off around 4 hours, Büning & Ritter 2004). An alternative possibility is that these systems may turn on as millisecond radio pulsars intermittently, during which they eject all mass transferred from their companion (Burderi et al., 2002; Burderi et al., 2003), as the transitional millisecond pulsars appear to do (Archibald et al., 2009; Papitto et al., 2013). If so, the time-averaged MT rate onto their NSs would be lower than we calculate. There is indeed evidence that these two systems may be transitional millisecond pulsars. XTE J1814-338 is known to show accretion-induced X-ray pulsations (Markwardt & Swank, 2003). Its donor star also shows evidence of irradiation by an unknown energy source, which may be spin-down energy from a radio pulsar (Baglio et al., 2013; Wang et al., 2017). EXO 0748-676 has not shown detected X-ray pulsations (despite sensitive RXTE X-ray observations). However, a careful study during quiescence showed that no accretion disk was present, which may indicate that a transitional millisecond pulsar had turned on, and is ejecting transferred mass (Ratti et al., 2012).

6 Conclusions

In this work, we have examined how different MB prescriptions affect the evolution of LMXBs. The observational data to which we compare our simulations is given in section 3. By systematically studying the parameter space of interest, we cover a range of possible seed masses and periods for these binary systems. The results of comparing our simulations to the observations are given in sections 4 and 5. The key results of this work are:

• •

Using "weaker" MB schemes such as the default Skumanich prescription, and even the "convection boosted" case described in this work, results in an overabundance of highly massive NS accretors.

• •

The highest density region in our parameter space, as seen in figure 10, is the region of short periods and low MT rates. This high-density region is found in all MB prescriptions.

• •

In the default, convection boosted and intermedaite MB prescriptions, all UXCBs of interest can be reproduced. Although, $\tau_{\rm max}$ is small suggesting these systems are difficult to form using this method.

• •

The "default" MB scheme reproduces results similar to Podsiadlowski et al. (2002). This weak MB scheme cannot reproduce some observed persistent systems in our simulations; these simulated binaries differ from observed binaries by up to an order of magnitude in MT rate.

• •

The "convection-boosted" prescription reproduces persistent systems much better than the default scheme, as it successfully simulates the properties of GX 9+9, 4U 1735-444 and Sco X-1. It cannot, however, reproduce the suspected super-Eddington system 2A 1822-371 (Bak Nielsen et al., 2017).

• •

Once we account for wind in the MB scheme, we can reach high enough MT rates to reproduce 2A 1822-371. Super-Eddington MT rates are achieved in the "intermediate" MB prescription.

• •

The intermediate prescription produces the largest number of observed transient LMXBs. Only EXO 0748-676 and XTE J1814-338 cannot be reproduced. These two systems only have an upper limit for mass transfer rate which my be the reason why these systems are difficult to reproduce.

• •

Including the effects of a non-thermal wind in our "wind boost" case results in very high MB. The high angular momentum loss results in MT rates that exceed $1M_{\odot}\rm\ yr^{-1}$. The "wind boost" case likely has reached and exceeded a saturation point with MB, and additional effects must be considered to dampen the angular momentum loss for this scheme.

The systematic mismatches between observed and predicted NS LMXB properties seen in previous work such as Podsiadlowski et al. (2002) are again found in our work when the default MB prescription is used. These discrepancies between observations and simulations begin to disappear, however, once we include the effects of convective turnover time and non-isothermal winds. With these changes, the MT rates approach those seen in observed systems, and our simulations more effectively reproduce the samples of persistent, and transient, binaries. With these results in mind, numerical studies of LMXBs that begin their calculations at long periods should no longer use the Skumanich prescription for MB as it does not adequately reproduce observed systems. Instead, the studies need to include additional effects in their MB schemes.

A clear extension of this work would be to include saturation effects in the MB prescription (Mestel & Spruit, 1987). One source of decreasing the MB strength that is not accounted for in this work, is the change in the magnetic field structure as the period of the system changes. It has been shown that in short period binaries a so-called ’dead zone’ is produced, trapping wind material (Mestel & Spruit, 1987; Ivanova & Taam, 2003). This trapped material cannot escape the system, reducing the angular momentum loss through MB. Additionally, the inclusion of irradiation-induced wind and more complex magnetic field structures, such as a dipolar field similar to what is done in Justham et al. (2006), is possible.

Possibilities for future analysis with the simulations produced for this work include the possibility of determining viable progenitors of observed LMXBs. Using the reproducibility search for observed systems in section 5, we can find a rough parameter space that produces progenitors for each of our LMXBs of interest for a given MB prescription. Finding the possible progenitors will act as a "reverse population synthesis" method, where instead of providing initial conditions, we use observed binaries and their progenitors to infer what the initial conditions may have been. The reverse population synthesis method cannot, however, be used with the simple Skumanich law, where many persistent systems are not reproduced.

Acknowledgements

We would like the thank the anonymous referee for helpful comments. NI acknowledges support from the CRC program, funding from an NSERC Discovery grant, and acknowledges that a part of this work was performed at the KITP, which is supported in part by the NSF under grant No. NSF PHY17-48958. CH acknowledges support from an NSERC Discovery Grant, and an Accelerator Supplement. This research was enabled by the use of computing resources provided by WestGrid and Compute/Calcul Canada.