Thermonuclear Burst Observations for Model Comparisons: A Reference Sample

Abstract

We present observations of thermonuclear (type-I) X-ray bursts, selected for comparison with numerical models. We provide examples of four distinct cases of ignition: He-ignition in mixed H/He fuel (case 1); He-ignition in pure He fuel, following exhaustion of accreted H by steady burning (case 2); ignition in (almost) pure He accumulated from an evolved donor in an ultracompact system; and a superburst, thought to arise from ignition of carbon fuel produced as a by-product of H/He bursts. For regular bursts, we measured the recurrence time and calculated averaged burst profiles from RXTE observations. We also estimated the recurrence time for pairs of bursts, including those observed during a transient outburst, modelled using a numerical ignition code. For each example we list the burst properties including recurrence time, fluence, peak flux, the persistent flux level (and inferred accretion rate), and the ratio of persistent flux to fluence. In the accompanying material, we provide a bolometric lightcurve for each burst, determined from time-resolved spectral analysis. Along with the inferred or adopted parameters for each burst system, including distance, surface gravity, and redshift, these data are suggested as suitable test cases for ignition models.

1 INTRODUCTION

Thermonuclear (type-I) bursts arise from unstable ignition of accreted fuel (typically mixed hydrogen and helium) on the surface of accreting neutron stars in low-mass binary systems. Such events are of high priority for observers, as they provide information about the fuel composition, accretion rate, and even neutron star spin, mass, and radius (e.g. Heger et al. 2007; Strohmayer & Bildsten 2006; Steiner, Lattimer, & Brown 2013). A key component contributing to our understanding of the burning physics is numerical modelling of the complex series of nuclear reactions which trigger and power the bursts (e.g. Fujimoto et al. 1981). To date, the degree to which predictions of numerical models have been compared in detail to observations is limited; this paper is part of a wider effort to address this situation.

The physical conditions and processes which broadly influence the burst behaviour are relatively well understood. The primary determinant is the local accretion rate $\dot{m}$ , frequently expressed as a fraction of the Eddington rate $\dot{m}_{\rm Edd}$ ¹ , at which the outwards force due to radiation pressure (assuming spherical symmetry) equals the surface gravity. The Eddington rate is also the predicted threshold at which the burning is expected to stabilise (e.g. Narayan & Heyl 2003), although observationally bursting behaviour appears to cease for many sources at a substantially lower level (Cornelisse et al. 2003; Galloway et al. 2008).

For systems accreting a mix of hydrogen and helium at typical accretion rates ( $\sim 0.1\dot{m}_{\rm Edd}$ ), the fuel layer is hot enough that H will burn stably via the hot-CNO cycle. This process is limited by β-decay, and will exhaust the H at the base of the layer in a time

(1) $$\begin{equation} t_{\rm CNO} = 9.8 \left(\frac{X_0}{0.7}\right)\left(\frac{Z_{\rm cno}}{0.02}\right)^{-1}\ {\rm h} \end{equation}$$

(e.g. Lampe, Heger, & Galloway 2016) which depends only on the accreted H-fraction X ₀ and the metallicity Z _cno. Thus, the composition of the fuel layer at ignition may vary, and a number of categories have been identified by numerical studies. If the recurrence time is short, the burst will ignite in a mixed H/He environment, leading to a ‘case 1’ (Fujimoto et al. 1981) or ‘prompt mixed’ (Narayan & Heyl 2003) burst. At lower accretion rates, the recurrence time may exceed t _CNO, leading to ignition in a pure-He layer, i.e. ‘case 2’ or ‘delayed helium’ bursts. At even lower rates, it is thought that H-burning too will become unstable, leading to ignition via that mechanism; but no unambiguous examples of such bursts have been observed.

An additional ignition case is provided by neutron stars in ultracompact systems, which likely accrete (almost?) pure He. Here, the ignition will always be in a H-poor environment, largely independently of the accretion rate.

The nuclear energy generation during the burst is related to the average H-fraction in the fuel < X > , with Q _nuc = 1.6 + 4 < X >MeVnucleon⁻¹. The composition of the fuel may be inferred from the α parameter, the ratio of burst energy to accretion energy over the burst recurrence time. For a neutron star of mass M, radius R,

(2) $$\begin{equation} \alpha = \frac{Q_{\rm grav}}{Q_{\rm nuc}}(1+z) \approx \frac{GM}{R} Q_{\rm nuc}^{-1} (1+z), \end{equation}$$

where z is the gravitational redshift at the neutron star surface, given by 1 + z = (1 − 2GM/Rc ²)^−1/2. Provided the accretion rate is proportional to the persistent flux F _p, α may be measured (up to a factor corresponding to the ratio of the anisotropy of burst and persistent emission) as

(3) $$\begin{equation} \alpha _{\rm obs} = \frac{F_{\rm p} c_{\rm bol} \Delta t}{E_{\rm b}}, \end{equation}$$

where c _bol is the bolometric correction to the (band-limited) flux F _p, Δt is the (regular) burst recurrence time, and E _b the burst bolometric fluence (total energy).

The composition of the fuel may also be deduced from the shape of the lightcurve. He burns during the burst via the triple-α reaction, which proceeds on a much faster timescale than the hot-CNO, rp, and (α, p) reactions which burn H. The higher the He mass fraction, the faster the burning will proceed during the burst, reflecting in the burst rise and duration, as well as the timescale τ (the ratio of the fluence to the peak flux).

The numerical codes that have been used to compare to observations fall into three broad classes. The first class determines the ignition conditions given the accretion rate and fuel composition, and predicts the burst recurrence time and energetics. While this model has been developed to compare with the inferred atmospheric expansion during a burst as suggested by measurements of burst oscillations (Cumming & Bildsten 2000), it has also been used to compare with observations of bursts at low accretion rates (Galloway & Cumming 2006). The primary disadvantage of such models is that they do not follow the time-dependent compositional structure of the atmosphere, which is significantly modified by the thermonuclear burning.

The next class of models are one-zone time-dependent codes that simulate some fraction of the nuclear reaction network. Such codes have been used to demonstrate the extent of the rp-process that powers mixed H/He bursts (Schatz et al. 2001), as well as serving to probe the sensitivity of burst lightcurves to individual reaction rates (Cyburt et al. 2016).

The current (practical) state of the art for burst modelling is in 1-D multi-zone models that track the full extent of the nuclear reaction networks. Codes in this class include kepler (Woosley et al. 2004); the general relativistic hydrodynamics code AGILE, coupled with a nuclear reaction network (Fisker et al. 2006), and MESA (Paxton et al. 2015).

1-D multi-zone models have been used to make the most detailed comparisons of observations to models, most notably for GS 1826 − 24 (Heger et al. 2007; Zamfir, Cumming, & Galloway 2012). However, these comparisons are generally limited in the extent of the sources and the observational data utilised. These limitations, despite the extensive modelling capabilities (and accumulated observational data) motivate the present study.

Here, we present a sample of observations of thermonuclear bursts, intended for comparison with numerical models. The objectives of the sample assembly are two-fold: first, as test cases to quantify variations between different codes (and hence intrinsic model uncertainties); and second, as examples that may be used to refine the system parameters by direct comparison with individual codes.

We deliberately take a pragmatic approach which incorporates constraints on system parameters obtained by comparison with particular numerical models (of varying degrees of fidelity). While this approach introduces dependence of the system parameters on the specific models chosen for these analyses, this dependence is unimportant for either of the two objectives described above. For the code comparisons, the specific choice of input parameters does not matter, provided they are plausible; and the parameters we have adopted are the best current parameter estimates for the systems studied. As for improving the system parameters via more detailed comparisons with individual models, the adopted system parameters provided here should be viewed as a suggested starting point.

This paper is presented as follows. In Section 2, we describe the sources selected to make up the sample. In Section 3, we describe the observational data from which the sample lightcurves and other parameters are drawn. In section Section 4, we briefly describe the inferred properties of the bursting sources, and compare the lightcurves in the sample. Finally, in Section 5, we suggest how the assembled data may be applied and used to test numerical models.

2 SOURCE SELECTION

We selected four well-studied burst sources, which span the range of theoretical ignition cases identified observationally (Table 1). We also selected a fourth source to serve as an example of a superburst, 4U 1636 − 536, which is notable as one of only two such events observed at high sensitivity with RXTE (e.g. Keek et al. 2014).

Table 1. Target thermonuclear burst source properties.

Values in italics indicate there are no constraints specific for that source.

References: (1) Heger et al. (2007); (2) this work; (3) Kuulkers et al. (2003); (4) Cumming (2003); (5) Özel et al. (2016).

We list in Table 1 those quantities required to simulate thermonuclear bursts and compare to the observations, including the distance, mass fraction of hydrogen X ₀ and of CNO nuclei, Z _CNO, and the adopted redshift 1 + z and surface gravity g. We caution that in most cases these parameters are not precisely known; however, for the purposes of deciding upon suitable values for simulation tests, we adopt values which are generally consistent with the observed burst properties and/or the neutron star population.

The surface gravity and fuel composition parameters are key input to the models, while the redshift is used to transform the burst lightcurves (and recurrence times) from the neutron-star frame to the observer’s frame (see Section 5). The source distance d is required to convert model predictions to observed quantities, and is subject to specific uncertainties depending upon the method by which it is estimated, as described further in appendix A.

Below, we describe the basic observational properties of each source.

2.1. GS 1826 − 24

This source, discovered as a new transient by Ginga in 1988 September (Tanaka 1989), has been X-ray bright and exhibiting bursts ever since. It is known as the ‘clocked’ or ‘textbook’ burster, due to its unusually regular and consistent bursts (Ubertini et al. 1999; Galloway et al. 2004). The inferred accretion rate varies little on short timescales, and is within the range 5–13% $\dot{m}_{\rm Edd}$ (Galloway et al. 2008; Chenevez et al. 2016). This pattern of regular bursting behaviour has been broken only recently, with the transition to weaker, irregular bursts accompanied by a soft spectral state, in 2014 June (Chenevez et al. 2016). By virtue of the typically consistent bursting behaviour, GS 1826 − 24 is also one of the few sources for which detailed efforts have been made to deduce the composition of the accreted fuel. However, these efforts have not yet resulted in unambiguous results.

The measured α-values for the source, at ≈ 40, strongly suggest a high fraction of H in the burst fuel, and so we adopt the solar value, i.e. X ₀ = 0.7, for the accreted H fraction. The lack of variation in burst fluence (i.e. total burst energy) over a range in accretion rate was taken as evidence by Galloway et al. (2004) that the fuel must be deficient in CNO nuclei. Otherwise, steady (hot-CNO) H-burning between the bursts would significantly change the fuel composition at ignition, as the burst recurrence time decreased from around 6 to 4 h. A subsequent comparison of burst properties, including the detailed shape of the lightcurve, with time-dependent 1-D model predictions by kepler (Woosley et al. 2004) led instead to the conclusion that the accreted fuel must be approximately solar in composition, i.e. Z _CNO ≈ 0.02, at odds with the earlier analysis (Heger et al. 2007). Since kepler includes substantially more physics than the code used in the earlier analysis, we adopt solar composition following the most recent study.

Despite the uncertainty in the fuel composition, the properties of the bursts from GS 1826 − 24 strongly suggest burning of mixed H/He fuel by unstable He ignition, i.e. ‘case 1’ of Fujimoto et al. (1981) [or alternatively, ‘prompt mixed’ bursts in regime 3 of Narayan & Heyl (2003)].

2.2. SAX J1808.4 − 3658

This source was discovered as a new transient in 1996 with BeppoSAX (in ’t Zand et al. 1998), and has been active every 2–3 yrs ever since. Observations made during the 1998 outburst by RXTE revealed 401 Hz X-ray pulsations, making it the first accretion-powered millisecond pulsar (Wijnands & van der Klis 1998). Bright, radius-expansion bursts have been observed in almost every outburst, and the detection of burst oscillations also at 401 Hz in bursts during the outburst of 2002 October confirmed the link between oscillation frequency and neutron-star spin (Chakrabarty et al. 2003).

The transient outbursts typically reach a maximum inferred accretion rate of only $0.05\ \dot{m}_{\rm Edd}$ , and the bursts have fast ( ≲ 0.5 s) rises indicative of He-rich fuel. The wide (2.01 h) orbit can accommodate a (H-rich) companion (Chakrabarty & Morgan 1998), and it is thought that the accreted H is burned away by β-limited hot CNO burning prior to ignition. The ≈ 20–30 h recurrence times, given the low accretion rates, are sufficient provided the metallicity, Z _cno is not far below the solar value. This type of burst corresponds to ‘case 2’ of Fujimoto et al. (1981), or ‘delayed helium’ (regime 4) of Narayan & Heyl (2003).

Comparisons of the burst properties with the predictions of numerical ignition models indicate a range of possible compositions, with Z _cno∝X (Galloway & Cumming 2006). The lowest plausible H-fraction in the accreted fuel was X ₀ ≈ 0.35.

2.3. 3A 1820 − 303

Located in the globular cluster NGC 6624, 4U 1820 − 303 (also known as Sgr X-4) was first detected as a bursting source with the Astronomical Netherlands Satellite (ANS; Grindlay et al. 1976). When in the bursting state, the source exhibits quasi-regular radius-expansion bursts with recurrence times in the range 2–4 h (e.g. Cumming 2003). The source intensity is modulated on a ≈ 176-d periodicity (Priedhorsky & Terrell 1984), and it only exhibits bursts in the (relatively short) low-intensity phase of the cycle, when the inferred accretion rate is ≈ 0.2–0.95 $\dot{M}_{\rm Edd}$ (e.g. Galloway et al. 2008).

The orbital period of 685 s is known from periodic modulation of both persistent X-rays and the optical counterpart (King & Watson 1986; Stella, White, & Priedhorsky 1987), and indicates an ‘ultracompact’ binary which is too close for a Roche-lobe filling main-sequence star. The mass donor is thus assumed to be a very low-mass He white dwarf, and thus accreting pure He (i.e. X ₀ = 0). Some evolutionary models predict a small (X ~ 0.1) amount of H in the surface layers, and this cannot be ruled out from comparing the observed bursts with ignition models (Cumming 2003).

The mass and radius of the neutron star in this system have been estimated by equating the blackbody normalisation in the cooling tail and the peak flux of PRE bursts, and adopting the inferred cluster distance (Özel et al. 2016). Although unresolved systematic errors may yet affect such measurements (e.g. Steiner et al. 2013), these values represent the best current estimates for the source, and are adopted in Table 1.

2.4. 4U 1636 − 536

4U 1636 − 536 is a persistently accreting burst source in a 3.8 h orbit (van Paradijs et al. 1990). The system is one of the most prolific type-I (thermonuclear) bursters known, and has been studied extensively with most major observatories (e.g. Lewin, van Paradijs, & Taam 1993; Galloway et al. 2008). The neutron star spin has been measured from burst oscillations and transient pulsations during a superburst at 579.3 Hz (Zhang et al. 1997; Strohmayer & Markwardt 2002). The relatively wide orbit, spectral features from hydrogen (Augusteijn et al. 1998), and the (at times) long profiles of the bursts strongly suggest the mass donor is a main-sequence star which accretes hydrogen-rich fuel.

The system has exhibited long-term variations in its persistent intensity. Between 1996 and 2001, the X-ray flux was in the range 3–6 × 10⁻⁹ergcm⁻²s⁻¹, indicating an accretion rate of 11–21% $\dot{m}_{\rm Edd}$ (cf. with Galloway et al. 2008). Since then, the flux has been consistently approximately 50% lower. At these accretion rates, it would be expected that bursts would consistently exhibit long profiles arising from ignition of mixed H/He fuel, although this is not always the case. 4U 1636 − 536 is also well known for short-recurrence time bursts occurring in trains of up to 4 (Keek et al. 2010).

The first superbursts were detected from RXTE/ASM data (Wijnands 2001), and two additional events were subsequently reported, one also detected by the PCA (Kuulkers et al. 2004; Kuulkers 2009). All the superbursts were detected while the system was in its higher flux range, prior to 2001. The 2001 event, on February 22 (MJD 51962.7), was observed with the PCA instrument, and has offered the highest signal-to-noise data during such an event, comparable only to the event seen from 4U 1820 − 303 (Strohmayer & Brown 2002). Preliminary spectral fitting of the 4U 1636 − 536 burst indicated that the spectrum could not be adequately fitted with blackbody models (Kuulkers et al. 2004), and subsequently it was shown that the burst spectrum included substantial contributions arising from reflection of the burst emission from the accretion disk (Keek et al. 2014).

3 OBSERVATIONS AND ANALYSIS

We used preliminary results from the Multi-INstrument Burst ARchive (MINBAR ² ), which includes bursts observed by the Rossi X-ray Timing Explorer Proportional Counter Array (RXTE/PCA), BeppoSAX Wide-Field Camera (WFC; Jager et al. 1997; in ’t Zand et al. 2004) and INTEGRAL Joint European X-Ray Monitor (JEM-X; Lund et al. 2003).

We calculated burst lightcurves from time-resolved spectroscopic analysis of selected bursts observed by RXTE/PCA (Jahoda et al. 1996). The MINBAR analysis follows the approach of Galloway et al. (2008), incorporating the latest PCA responses and the effects of deadtime, as for the measured peak PRE burst fluxes (see appendix A). Earlier analyses of samples including these data have demonstrated evidence for an increased contribution of persistent flux during bursts, possibly arising from the effects of Poynting–Robertson drag on the inner accretion disk (Worpel, Galloway, & Price 2013, 2015). While the fractional increase in the persistent flux can be ~ 10, this contribution is a small fraction of the burst flux, and so any correction required to the burst flux is also small (of order a few per cent, similar to the intrinsic uncertainty of the time-resolved flux measurements). Thus, we neglect this effect for the purposes of determining our mean profiles.

3.1. Selecting burst sequences

For GS 1826 − 24 and 4U 1820 − 303, we identified burst sequences in MINBAR for which we could reliably infer the burst recurrence time, even when the observations were interrupted. Such interruptions arise for instruments in low-Earth orbit, due to Earth occultations and passages through the South Atlantic Anomaly. An example burst train is shown in Figure 1.

Figure 1. Example lightcurve showing bursts observed in 2000 September by BeppoSAX/WFC and the Rossi X-ray Timing Explorer/PCA from GS 1826 − 24. The top panel shows the lightcurves from the two instruments, rebinned to 10-s resolution, with the detected bursts marked with open symbols. The WFC lightcurve is shown in black, with open diamonds indicating the bursts; for PCA, the lightcurve is red, with open squares marking the bursts. The vertical lines above the symbols indicate the predicted times of bursts according to a constant recurrence time model; note that the unobserved bursts in the model train consistently fall in the data gaps. In the lower panel, we show the residual to the constant recurrence time model. Here, the best-fit recurrence time is 4.177 ± 0.010 h, with an RMS error of 9.53 min.

We then selected and averaged bursts observed by the RXTE/PCA within each burst train. MINBAR also includes bursts observed by BeppoSAX/WFC and INTEGRAL/JEM-X, which are useful for establishing the regularity of the bursts, but the sensitivity of these instruments is much lower than RXTE, and so are unsuitable for producing high-quality burst lightcurves.

For GS 1826 − 24, we augmented our burst sample with optical bursts observed during RXTE observations in 1998 June, as also analysed by Thompson et al. (2008).

We measured the average recurrence time of each train by assigning a trial integer value to each burst, and then performing a linear least-squares fit to the times. The trial values were determined by dividing the time since the first burst, by the minimum burst separation of the train; we subsequently adjusted the sequence numbers to minimise the RMS value of the fit. We typically obtained an RMS value of ≈ 20 min or better.

From the set of bursts in each train N _burst, we then selected the RXTE bursts (N _av, Table 2), and created average lightcurves from each train. We also calculated averaged burst fluences, peak fluxes and α = ΔtF _p c _bol/E _b, where F _p is the persistent (accretion) flux, c _bol the bolometric correction (see Section 3.2), and E _b the burst fluence.

Table 2. Properties of thermonuclear bursts observed from target sources by the Rossi X-ray Timing Explorer.

Burst index number in the catalogue of Galloway et al. (2008).

^b Persistent flux in the energy range 3–25 keV.

^c Inferred from ignition model comparisons.

^d Calculated based on linear interpolation between observations falling within the burst interval.

^e Only two bursts detected, so may not have been part of a regular train.

^f Recurrence time upper limit of 1.75 yr based on the next earliest event, detected by the RXTE/ASM (Kuulkers et al. 2004).

^g Values taken from the ‘optimal’ fits of Keek et al. (2014).

We estimated the accretion rate $\dot{m}$ as a fraction of the Eddington rate $\dot{m}_{\rm Edd}$ as

(4) $$\begin{equation} \dot{m}=4\pi d^2 F_{\rm p}c_{\rm bol}/\dot{m}_{\rm Edd}, \end{equation}$$

where F _p is the average persistent flux over the burst train, and d the estimated source distance from Table 1. This estimate neglects the possible anisotropy of the persistent emission, so must be taken as a guide only.

The identification of suitable burst trains for 4U 1820 − 303 was made difficult by the relatively small number of bursts available in the MINBAR sample. Just 67 bursts were detected by BeppoSAX/WFC, INTEGRAL/JEM-X, or RXTE/PCA, with only 16 events observed with RXTE. We ultimately identified eight burst trains with three or more bursts and evidence for regular bursting behaviour; the range of inferred recurrence time was 2.7–4.4 h.

Only two of these sequences included bursts detected with RXTE, each with just three bursts. For each sequence, we observed a pair of bursts within a few hours of each other, and a third burst some days earlier or later. We ultimately rejected one of these sequences, in 2009 May (MJD 54980; see also in’t Zand et al. 2012) because there were significant variations in the persistent flux (and spectral shape) over the interval, and it was not possible to be confident about the recurrence time for the close pair of bursts, due to data gaps between them. This was not the case for the remaining train, in 1997 May (MJD 50572), and so we retained this set of bursts. With only one burst observed by RXTE in this sequence, we simply provide the time-resolved spectroscopy results from that burst, as was done for SAX J1808.4 − 3658 (see Section 3.4).

Due to the sparse data for 4U 1820 − 303, we augmented the sample with a pair of bursts observed by RXTE/PCA on 2009 June 12, separated by only 1.89 h. The nearest observations in time preceding (following) that observation were 13 d before (3 d after), so any possible test for periodicity over that time range would likely be uninformative given the likely phase drift of the burst sequence. Assuming that the separation of 1.89 h does reflect a regular recurrence time, this pair would be the most frequent bursts ever observed from this system (cf. with Clark et al. 1977).

The constraints on recurrence times for the superbursts from 4U 1636 − 536 is poor due to the low duty cycle of observations between the detected events. The superburst in our sample was preceded by another 1.75 yrs prior. However, the inferred accretion column is approximately twice the ignition column, suggesting that the actual recurrence time is closer to 0.9 yrs (Keek et al. 2015). In addition to the properties of the burst itself reported in Table 2, we constrained the ‘quench time’ corresponding to the period following the superburst during which normal thermonuclear burning is suppressed. The first thermonuclear (H/He) event detected from 4U 1636 − 536 following the superburst was on MJD 51985.5 (BeppoSAX observation 10898), giving t _quench ⩽ 22.8 d. The observation during which that event was detected commenced 17.5 h before the burst was detected, and had a duty cycle of ≈ 70%, with interruptions due to the low-Earth orbit. The coverages suggests a ≈ 30% chance that an earlier burst could have been missed, so that the quench time could have been shorter by up to 17.5 h. However, no intervening observations were made to rule out additional bursts missed earlier than that observation.

3.2. Bolometric corrections

We measured the persistent flux in the 3–25 keV band from RXTE and BeppoSAX observations by fitting commonly used phenomenological models, and converted the integrated flux to a bolometric value by applying a bolometric correction for each observation. The correction was estimated by fitting to RXTE PCA and HEXTE data a broadband model typically comprising an absorbed Comptonisation component (compTT in xspec; Titarchuk 1994), usually with a Gaussian component representing Fe Kα emission in the range 6.4–6.7 keV. We measured the 3–25 keV flux from each fit, and then created an ideal response (using dummyrsp in xspec) covering the energy range 0.1–1 000 keV, and integrated the flux over this range, excluding the effects of neutral absorption. We estimated the bolometric correction as the ratio of the unabsorbed 0.1–1 000 keV flux to the absorbed 3–25 keV flux; for the observations here, the correction was typically in the range 1.4–2.

This approach is reliable provided that the contribution to the bolometric flux outside the range to which PCA and HEXTE are sensitive (typically 3–100 keV) is small (or at least consistent). However, it is known that for GS 1826 − 24, that additional low-energy contributions may arise at times, so that the persistent flux measured by RXTE in the 3–25 keV band is not always a reliable estimator of the accretion rate (e.g. Thompson et al. 2008). In some observations, there is evidence for unusually low α-values, suggesting that a significant fraction of the X-ray flux is emitted outside the RXTE band; contemporaneous Chandra or XMM-Newton observations suggest that this flux may be emitted as a low-temperature ( ≲ 1 keV) thermal component. Thus, we avoided from our data selection bursts from such epochs identified by Thompson et al. (2008), including those from 2003 April, and also the observations from 1998 June and 1997 November. An additional epoch of unusually low α ≈ 26 is from observations in 2006 August, and was also excluded. The remaining set of three burst trains have broadly consistent α-values, but we caution that there may yet be additional systematic errors contributing to the determination of the accretion rate in this (and possibly other) systems.

3.3. Lightcurve fitting

Here, we describe the approach by which the values of the neutron star radius, redshift, and surface gravity for GS 1826 − 24 were chosen. We replicated the analysis of Heger et al. (2007) with updated analysis results and an improved method to determine the best-fit system parameters. We selected the train of bursts observed in 2000 September, with a recurrence time of 4.177 h, similar to that used previously. The seven bursts observed with RXTE are augmented by the observation of 25 additional bursts with BeppoSAX over the same interval, with the complete train spanning 12 d. However, the entire sequence of bursts is not consistent with a steady recurrence time, likely due to small variations in accretion rate. Thus, we restricted the sequence to a shorter interval covering the RXTE bursts, beginning on MJD 51810 and spanning almost 4 d (as plotted in Figure 1). We then compared the averaged lightcurve (calculated from the RXTE bursts only) with examples predicted by the kepler code, tabulated in Lampe et al. (2016).

Heger et al. (2007) found excellent agreement between the observed lightcurve and a model referred to as ‘A3’, with solar composition, an accretion rate of 1.58 × 10⁻⁹M_⊙yr⁻¹. In the sample of Lampe et al. (2016), the accretion rate is quantified in terms of the Eddington value, defined as 1.75 × 10⁻⁸M_⊙yr⁻¹. Thus, for model A3, the accretion rate was $\dot{M}=0.09\dot{M}_{\rm Edd}$ , and this model corresponds to ‘a05d’ in the Lampe et al. (2016) sample.

In our revised comparison, we simultaneously matched the lightcurve and the recurrence time, using a Markov-Chain Monte Carlo (MCMC) code emcee (Foreman-Mackey et al. 2013) to marginalise over the three parameters of interest: the redshift 1 + z (by which the lightcurve will be stretched when transforming into the observer’s frame), the relative intensity scale, including the distance and redshift:

(5) $$\begin{equation} F_X = \frac{L_X}{4\pi d^2\xi _{\rm b}(1+z)}, \end{equation}$$

where L_X is the burst luminosity predicted by kepler, d is the distance, ξ_b takes into account the possible anisotropy of the burst emission (see also appendix A, and F_X the observed burst flux). A third parameter is related to the offset in time between the model and observation, and has no physical meaning.

Our initial runs with this code did not match the recurrence time, and to ensure that the model-predicted burst interval (when scaled by 1 + z) matched that observed, we weighted this parameter in the likelihood function by a factor of 50. This ensured in the subsequent comparisons that the redshifted model recurrence time matched that observed to within a few percent.

The comparison with model a05d, as used by Heger et al. (2007), gave good agreement with the observed lightcurves, although requiring a redshift of 1 + z = 1.392. Following the approach of Lampe et al. (2016), we transformed the model results (calculated in a Newtonian frame) to the observer’s frame assuming that the gravity, and the Newtonian and the general-relativistic neutron star masses are identical. This implies that the radius required for the Newtonian and general-relativistic gravity will be different, and the corresponding radius (for a neutron star mass of 1.4M_⊙) implied by our comparison of model a05d was 8.5 km, well below the lower limit of 10.4 km inferred for a wide range of neutron stars (e.g. Steiner et al. 2013).

Thus, we made a limited exploration of alternative models, focussing on those with solar composition. The comparison with model a028 gave a best-fit redshift of 1.234, implying a radius of 12.1 km, well within the expected range. The accretion rate (in units of the Eddington rate) for this model is 0.082, which is somewhat higher than that inferred for the 2000 Sep epoch, at 0.0692; but may be consistent taking into account the persistent emission anistropy. Heger et al. (2007) estimated the ratio of persistent to burst anisotropies at ξ_p/ξ_b = 1.55, implying ξ_p > 1 (e.g. Fujimoto 1988); in that case, the luminosity (and hence accretion rate) inferred from the persistent flux would be higher, bringing the observation estimate roughly in line with the value assumed for run a028. The corresponding distance is dξ^1/2 _b = 6.1 kpc. We do not quote errors on these parameters, as the lightcurve comparison (Figure 2) is not of sufficient quality to have confidence in the posterior distributions (the estimated probability densities for the model parameters). However, the comparison is of comparable quality as in Heger et al. (2007), and so for the present purposes, we adopt the values for the redshift and distance.

Figure 2. Best-fitting comparison of averaged, observed burst lightcurves from GS 1826 − 24, with kepler model a028. The symbols with errors show the average lightcurve observed by RXTE in 2000 September, at which time the recurrence time was 4.177 h The model curve, rescaled on both axis based on the best-fit distance and redshift, is the red solid curve. The best-fitting distance and redshift, determined from an MCMC analysis, were dξ^1/2 = 6.1 kpc and 1 + z = 1.23, respectively.

3.4. Outburst fitting

The approach adopted for SAX J1808.4 − 3658 was necessarily different, as this source has not exhibited regular bursts in previous observations. The persistent flux (and hence accretion rate) varied significantly over the time at which the bursts were observed in 2002 October, and so the burst separations (and hence, presumably, the recurrence times) vary significantly. Furthermore, the burst separations are all longer than 12 h, so that many data gaps are present between each burst, introduced by the interruptions in the visibility of the source due to the observing satellite’s ≈ 90-min low-Earth orbit.

Galloway & Cumming (2006) compared the observed burst times with the predictions of a numerical ignition model, to constrain the system parameters. For the purposes of this paper, we have revised the analysis of those authors, to reconcile the changes in the RXTE/PCA calibration since that work was completed, to correctly treat the persistent and burst anisotropy factors, and to use a more rigorous parameter space exploration to determine the best-fit system parameters and uncertainties.

We matched the burst times, fluences, and α-measurements to the settle code developed by Cumming & Bildsten (2000) using MCMC and marginalising over the accreted hydrogen fraction X ₀, CNO metallicity Z _CNO, base flux Q_b , and scaling factors incorporating the source distance and emission anisotropy. The settle model predicts burst parameters in the neutron star frame, given an inferred accretion rate. Necessarily, we adopt the average inferred accretion rate (based on linear interpolation between each observation) for each burst interval, as settle does not simulate varying $\dot{m}$ . The redshift 1 + z and (Newtonian) surface gravity g are fixed in the settle code, at the values listed in Table 1 so the results are implicitly for those values.

The best-fit parameters we report are dependent on the fidelity of the model, which does not take into account all of the underlying physics, such as the effects of thermal and compositional inertia. More detailed models, such as kepler (Woosley et al. 2004) can provide more rigorous modelling of the underlying physics, however are computationally expensive and so impractical for the Monte Carlo approach adopted here. However, kepler simulations to reproduce the bursts observed from SAX J1808.4 − 3658 during the 2002 outburst have been performed recently, demonstrating that the assumption of a constant accretion rate between bursts results in the recurrence time being underestimated, when compared to the scenario in which the accretion rate decreases between bursts (Johnston et al., in preparation). To account for this bias, we corrected the recurrence time predictions of settle using correction factors determined by kepler. The correction factor depends on the $\dot{m}$ gradient, and thus is different for each burst. We adopt correction factors for each of the seven predicted bursts of SAX J1808.4 − 3658, which vary from ≈ 1 for the first burst to ≈ 1.2 for the final burst.

We infer the persistent and burst anisotropy factors (ξ_p and ξ_b, respectively) from the scaling factors and find ξ_p = 1.05 ± 0.01 and ξ_b = 0.91^+0.04 _−0.03. This implies an inclination of 62 ±1○. Note that the error on the inclination is statistical only, and does not include the contribution from uncertainties in the mass and radius, due to the fixed values adopted in the settle model.

As with the analysis of Galloway & Cumming (2006), our maximum-likelihood solution predicts an additional pair of bursts between the first pair of observed bursts, on MJD 52562.41363 and 52564.30515 (#1 and 2 from Galloway et al. 2008). These events would have fallen in data gaps, which explains why they were not observed. We thus report the burst separations for the last two bursts observed during the 2002 October outburst as the recurrence time; for the second burst, we estimate the recurrence time inferred from the burst train modelling based on the simulated time of the preceding event. We note that the inferred recurrence time of 16.55 h for the second burst observed with RXTE is comfortably within the overall range for the source, with the minimum separation of 12.6 h set by a pair of bursts observed by BeppoSAX (in ’t Zand et al. 1998).

4 RESULTS

The derived (or adopted) properties for the sources of interest are listed in Table 1. The properties of the bursts selected for our sample are summarised in Table 2. We plot the burst lightcurves for each of the sequences in Figure 3. The accompanying material for this paper includes lightcurves (averaged in cases where there are more than one RXTE burst in the sequence) for comparison with numerical model results.

Figure 3. Representative burst lightcurves from 4U 1820 − 303, SAX J1808.4 − 3658, GS 1826 − 24, and 4U 1636 − 536 observed by the Rossi X-ray Timing Explorer. In each panel, we plot the lightcurve at the indicated recurrence time, with flux converted to luminosity at the distances indicated in Table 1. The bursts contributing to each profile are listed in Table 2; where more than one profile is observed by RXTE in the train, we show the average profile. Note the change in x- and y-axis between the top panels (4U 1820 − 303 and SAX J1808.4 − 3658) and the lower two panels.

The burst lightcurves show clear differences related to the inferred composition and accretion rate. For 4U 1820 − 303, the burst lightcurves are short, with durations (defined as the interval over which the luminosity exceeds 10% of the maximum) of 15 s or so. Each of these bursts exhibits strong photospheric radius expansion, and the quality of the spectral fits during the radius expansion episodes are relatively poor (reduced-χ² values in the range 2–8). As $\dot{m}$ increases, the burst timescale increases slightly (τ = 6.24 to 6.55 s), which qualitatively supports the presence of a small amount of H in the accreted fuel. With smaller recurrence time, there is less time to burn the accreted H, and the mass fraction at ignition will be larger.

The bursts from SAX J1808.4 − 3658, although also thought to be powered primarily by He, have a notably different morphology, with an extended radius-expansion phase at roughly constant luminosity. The timescale variation for SAX J1808.4 − 3658 is in the opposite direction, with the burst timescale and duration becoming shorter at higher accretion rates. Here, the recurrence time is long enough that any accreted H has already been exhausted at the base, so that lower accretion rate allows a larger pile of He to accumulate. The resulting increased fluence thus requires a longer interval emitting at the Eddington luminosity, with some of the additional energy likely exported as kinetic energy of outflowing material.

The bursts from GS 1826 − 24 are different again, showing the much slower ( ≈ 5 s) rises and long rp-process powered tails characteristic of these mixed H/He bursts. The peak luminosity is significantly smaller than for the He-rich bursts, and durations are approximately a factor of 5 longer. The burst profiles show only very subtle variations from epoch to epoch, with the peak flux decreasing slightly, while the timescale increases, as the accretion rate increases.

The superburst profile is principally distinguished from the shorter, H or H/He events by the orders-of-magnitude longer timescale, with a duration of more than 2 h.

5 DISCUSSION

We have assembled a sample of burst observations intended for comparisons with, and between, numerical models. Comparisons of different model codes for the same test cases will allow estimates of the intrinsic uncertainty that affects burst simulations. These uncertainties may also depend on the particular type of burst being simulated; that is, the variation in model predictions for some types of bursts may be larger than for others. For example, Galloway & Cumming (2006) suggest that the simple analytic model adopted for the simulations of bursts from SAX J1808.4 − 3658 offer comparable fidelity to more sophisticated models, including kepler. This is less true for simulations of H-rich bursts from GS 1826 − 24 observed at higher accretion rates, as the simple models do not include the effects of compositional or thermal inertia.

It is also anticipated that detailed comparisons of the burst measurements here may allow more stringent constraints to be placed on the neutron star parameters, as has been demonstrated as a proof-of-principle in Sections 3.3 and 3.4.

The anticipated procedure to replicate each of the burst measurements in Table 2 is as follows:

1. 1. Set the surface gravity g and neutron star radius R to the values appropriate for the source of interest in Table 1.

2. 2. Set the composition in the accreted material, X ₀ and Z _CNO.

3. 3. Perform a model run with a constant accretion rate chosen from one of the entries for that source from Table 2, and extract the burst luminosity as a function of time.

4. 4. Re-scale the predicted luminosity appropriately for the redshift parameter 1 + z for this source, listed in Table 2, and also the burst anisotropy parameter ξ_b, where known. For the second and third bursts from SAX J1808.4 − 3658, the predicted recurrence time and burst fluence may be adjusted by the factor described in Section 3.4 to account for the effects of the declining accretion rate during each burst.

5. 5. Calculate from the simulated burst train the recurrence time Δt _pred, and from the burst lightcurve, the burst fluence, peak flux, α-value, and timescale.

6. 6. Compare these predictions with the observed lightcurve and the measured values in Table 2.

In the event that the model recurrence time is significantly longer (shorter) than the observation, the cause may be that the true accretion rate is different from the estimate due to the unknown degree of anisotropy ξ_p of the persistent emission. In that case, the adopted accretion rate might be increased (decreased), by roughly the ratio of the predicted and observed recurrence times, and the procedure repeated from step 3 onwards.

One additional parameter that is commonly included in burst simulations is the heat flux from below the simulation region, commonly parameterised as a ‘base flux’ Q _b. This quantity is not well constrained by observations, and furthermore depends on the accretion rate (e.g. Figure 18 of Cumming et al. 2006); most studies to date have adopted values in the range 0.1–0.15 MeV nucleon⁻¹. For the purposes of model-to-model comparisons, we recommend adopting a value of 0.1, but suggest that for the purposes of matching the observations, this parameter can be considered free within the range 0.1–0.8 MeV nucleon⁻¹, depending upon the choice of core neutrino emissivity model.

The approach to replicate the superburst observation is necessarily different. Modelling a superburst by accreting hydrogen/helium comes at great computational expense, because hundreds of hydrogen/helium flashes must be simulated to create the carbon fuel. Furthermore, it remains challenging to reproduce the observationally inferred carbon mass fraction (e.g. Stevens et al. 2014; Keek & Heger 2016). An alternative approach is to directly accrete a mixture of carbon and heavier elements (e.g. Fe or Ru; Cumming & Bildsten 2001; Schatz, Bildsten, & Cumming 2003). For the 2001 superburst from 4U 1636 − 536, we recommend a fuel composition of 26% ¹²C and 74% ⁵⁶Fe (Keek et al. 2015). Such simulations are able to describe the observed superburst light curves in detail, but generally predict larger ignition depths than observed (Keek & Heger 2011).

The ultimate goal of this study is to constrain the rates of particular nuclear reactions, or masses of particular isotopes, based on the burst observations. It has been shown that several reactions are likely to affect the shape of the burst lightcurve (e.g. Cyburt et al. 2016). Similarly, nuclear masses have also been shown to have an effect (Schatz & Ong 2016). It is hoped that the carefully assembled and calibrated data presented here offer the first step on the road towards precision nuclear physics from thermonuclear bursts.