A SPECTROSCOPIC STUDY OF THE EXTREME BLACK WIDOW PSR J1311–3430

, , and

Published 2015 May 11 © 2015. The American Astronomical Society. All rights reserved.
, , Citation Roger W. Romani et al 2015 ApJ 804 115 DOI 10.1088/0004-637X/804/2/115

0004-637X/804/2/115

ABSTRACT

We report on a series of spectroscopic observations of PSR J1311–3430, an extreme black-widow gamma-ray pulsar with a helium-star companion. In a previous study we estimated the neutron star mass as ${{M}_{{\rm NS}}}=2.68\pm 0.14\;{{M}_{\odot }}$ (statistical error), based on limited spectroscopy and a basic (direct heating) light-curve model; however, much larger model-dependent systematics dominate the mass uncertainty. Our new spectroscopy reveals a range of complex source behavior. The variable He i companion wind emission lines can dominate broadband photometry, especially in red filters or near minimum brightness, and the wind flux should complete companion evaporation in a spin-down time. The heated companion face also undergoes dramatic flares, reaching ∼40,000 K over ∼20% of the star; this is likely powered by a magnetic field generated in the companion. The companion center-of-light radial velocity is now well measured with ${{K}_{{\rm CoL}}}=615.4\pm 5.1$ km s−1. We detect non-sinusoidal velocity components due to the heated face flux distribution. Using our spectra to excise flares and wind lines, we generate substantially improved light curves for companion continuum fitting. We show that the inferred inclination and neutron star mass, however, remain sensitive to the poorly constrained heating pattern. The neutron star's mass, ${{M}_{{\rm NS}}}$, is likely less than the direct heating value and could range as low as 1.8 M$_{\odot }$ for extreme equatorial heating concentration. While we cannot yet pin down ${{M}_{{\rm NS}}}$, our data imply that an intrabinary shock reprocesses the pulsar emission and heats the companion. Improved spectra and, especially, models that include such shock heating are needed for precise parameter measurement.

Export citation and abstract BibTeX RIS

1. INTRODUCTION

PSR J1311–3430 (hereafter J1311) is a "black widow" millisecond pulsar (MSP; Pletsch et al. 2012) with spin-down luminosity $\dot{E}=5\times {{10}^{34}}\;{{I}_{45}}$ erg s−1 (for a neutron star moment of inertia ${{10}^{45}}\;{{I}_{45}}$ g cm$^{2}$) in a PB = 93.8 minute orbit (Romani 2012), the shortest of any confirmed spin-powered pulsar.5 J1311 is particularly interesting, as photometric and spectral studies (Romani et al. 2012, hereafter R12) show that it has a strong, violently variable evaporative wind driven from an ${{M}_{2}}=0.01\;{{M}_{\odot }}$ He companion and that the neutron star is particularly massive. In R12 we found that with standard assumptions for heating of the companion by the pulsar, the light-curve modeling gave ${{M}_{1}}=2.68\pm 0.14\;{{M}_{\odot }}$ (statistical errors), although with different assumptions for asymmetric heating one could obtain masses ∼2.1–2.9 M$_{\odot }$ (systematics range). This system promises to reveal much about close binary evolution (van Haften et al. 2012; Benvenuto et al. 2013), the efficacy of the companion evaporation channel for producing isolated MSPs, and (possibly) fundamental constraints on the equation of state of matter at supernuclear densities.

Measurement of the orbital motion, heating, and wind stripping of the companion is the key to such studies. To this end, we have made optical spectroscopic observations at several epochs. Here we summarize the measurements, describe the new phenomena discovered, introduce new models of the system, fit for the component masses, and comment on the results.

2. SPECTROSCOPIC CAMPAIGNS

In R12, we discussed photometric evidence for substantial optical flaring activity of J1311. These flares appear to have a red spectral energy distribution (SED), but further study of this behavior requires spectroscopic characterization of the variable component. Also, multi-orbit photometry shows variable asymmetry about the optical maximum. This limits the reliability of model fitting (with simple symmetric pulsar heating) in determining the system inclination i and the correction factor ${{K}_{{\rm cor}}}={{K}_{{\rm CoM}}}/{{K}_{{\rm CoL}}}$, where ${{K}_{{\rm CoL}}}$ is the observed radial-velocity amplitude from spectral features on the heated face of the companion (the center of light; CoL) and ${{K}_{{\rm CoM}}}$ is the true radial-velocity amplitude of the companion center of mass (CoM). These factors are critical, since for a given observed radial-velocity amplitude the pulsar mass is ${{M}_{{\rm PSR}}}\propto {{({{K}_{{\rm cor}}}/{\rm sin} \;i)}^{3}}$. The detailed behavior of the photospheric absorption lines can probe the surface-temperature distribution, and emission-line components can probe the evaporating wind.

Thus, intensive spectroscopic campaigns offer the best hope for further understanding this important system. We report here on our efforts to date with the Keck I/Low Resolution Imaging Spectrometer (LRIS; Oke et al. 1995) and the Gemini/GMOS-S (Hook et al. 2005) systems. For both, we sought long-slit multi-orbit coverage with high resolution and sensitivity. By placing the slit at a position angle (PA) of 10°, we were able to simultaneously observe the pulsar and an ${{r}^{\prime }}=18.05$ mag G-type comparison star 15'' from J1311 (Figure 1). The latter allowed us to stably cross-calibrate the flux and wavelength, permitting coherent analysis of both data sets. With Keck I we suffered from poor weather, with only half a night of the two scheduled nights delivering J1311 spectra, under marginal conditions. Nevertheless, we obtained three orbits of continuous coverage. The GMOS-S observations were queue scheduled and so all data were usable, with three visits covering slightly more than one orbit each. We describe these observations in detail, along with several short supplementary spectral sequences that give context to the variable behavior, especially near binary optical maximum.

Figure 1.

Figure 1. 120 s ${{g}^{\prime }}$ SOAR/SOI image of PSR J1311–3430 at maximum ${{\phi }_{B}}\approx 0.75$, showing the pulsar, the comparison star, and the 1''-wide slit position.

Standard image High-resolution image

2.1. 2013 March LRIS Spectroscopy

We observed J1311 with Keck I/LRIS on 2013 March 11 (MJD 56,362.429–56,362.626), using the 1''-wide slit, a dichroic splitter at 5600 Å, the 600/4000 grism in the blue camera, and the 400/8500 grating in the red camera. Table 1 lists the wavelength coverage and resolution for all observations. Our goal was to obtain high signal-to-noise ratio (S/N) spectra at all orbital phases. With a high southern declination, observations from Mauna Kea were perforce at moderate airmass (${\rm sec} z$ = 1.7–2.5). The seeing was poor (mostly $\gt 1\buildrel{\prime\prime}\over{.} 3$) and the transparency variable. Each exposure was 300 s in duration; this was a compromise between improving the S/N but not having too much orbital smearing.

Table 1.  PSR J1311–3430 Spectroscopic Campaigns

Date Tel./Instr. Range Res Exp. Orb. Phase $\langle {\rm He}\;{\rm I}\;6678\rangle $
(MJD)   (Å) (Å) ($N\times s$) (${{\phi }_{B}}$) (${{10}^{-18}}\;{\rm erg}\;{\rm c}{{{\rm m}}^{-2}}\;{{{\rm s}}^{-1}}$)
56064a Keck I/LRIS 3100–10,500 4,7b 6 × 300 0.56–0.93 <3
56362 Keck I/LRIS 3100–10,500 4,7b 49 × 300 −0.18–2.85 200.
56390 Keck II/DEIMOS 4450–9630 4.6 3 × 300 0.77–0.90 390.
56418 Gemini/GMOS-S 3950–5100 0.9 38 × 175 0.58–1.85 28.c
56446 Gemini/GMOS-S 3950–5100 0.9 38 × 175 −0.16–1.10 <10c
56448 Gemini/GMOS-S 3950–5100 0.9 38 × 175 −0.17–1.09 <8c
56422 Keck I/LRIS 3100–10,500 4,7b 3 × 300 0.01–0.14 9.0

Notes.

aPreviously described by Romani et al. (2012); all other observations presented here for the first time. bBlue-side resolution, red-side resolution. cNo red-side coverage; estimate from He i λ4471.

Download table as:  ASCIITypeset image

Standard processing and optimal extraction were applied; fluxes were calibrated against the red-side and blue-side standard stars BD + 26°2606 and BD + 28°4211, respectively. The comparison-star flux varied by a factor of 1.7 during this observation sequence, largely owing to variable seeing-induced slit losses. After normalizing the flux, as monitored by the comparison star, we show (Figure 2) "trailed spectra" for the central portions of the blue and red sides.

Figure 2.

Figure 2. Trailed spectra from the 2013 March Keck/LRIS campaign. Fluxes have been normalized to remove variations in the simultaneously observed comparison star. Upper panel: blue side, showing Doppler-shifting absorption features from the companion photosphere. Lower panel: red side, showing variable He i emission features from the companion wind. The scale is inverted: darker indicates more flux.

Standard image High-resolution image

The 49 spectra cover just over three orbits (${{\phi }_{B}}=-0.18$ to 2.85), with the midpoint of the barycentric arrival time of each exposure referenced to the ephemeris of Pletsch et al. (2012). Here we set ${{\phi }_{B}}=0$ at the first pulsar ascending node after the start of observations. For convenience, we will follow Breton et al. (2013) in referring to phase near optical maximum (pulsar superior conjunction, ${{\phi }_{B}}=0.75$) as "Day" and phases viewing the companion's unheated face (${{\phi }_{B}}=0.25$) as "Night." "Dawn" and "Dusk" are at phases ${{\phi }_{B}}=0.5$ and ${{\phi }_{B}}=0.0$, respectively.

Several features are evident in these trailed negative (darker means more flux) spectra. First, we see the continuum fading and brightening from spectrum 1 (bottom) to 49 (top). The continuum modulation is stronger in the blue than in the red, as expected from the high temperature of the companion's heated Day face. Superposed on this continuum are Doppler-shifted absorption features (R12). We also see broad and variable emission lines, strongest in the Dusk and Night phases. While little emission is visible during the first orbital minimum (spectra 5–10), the emission lines are very strong in the third Night (spectra 35–40) and indeed have such large equivalent width that they dominate in the red. Since the companion is very H-poor (R12), both the photosphere and the wind trail driven from it are dominated by He i.

2.2. 2013 Gemini Spectroscopy

The Keck observations have good sensitivity and broad wavelength coverage, but lack the resolution to detect the line-shape distortions imposed by nonuniform companion heating. To search for such effects, we observed J1311 with Gemini, using the GMOS-S spectrograph in a novel mode to maximize resolution in the blue part of the spectrum, which displays many photospheric He i and metal lines. These observations used the R831 grating in second order through a 0farcs 75 slit, which covered 3950–5100 Å at $R\equiv \lambda /\delta \lambda \approx 5000$, $\sim 60\;{\rm km}\;{{{\rm s}}^{-1}}$ resolution. This configuration required an order-blocking filter. The best available filter, ${{g}^{\prime }}$, unfortunately has a 4000 Å cutoff that decreases sensitivity to Ca H; moreover, the protected silver coating on the primary mirror greatly limits its effective area below 4000 Å. To minimize the orbital velocity smearing, we restricted exposures to 175 s (i.e., maximum smearing of $60\;{\rm km}\;{{{\rm s}}^{-1}}$ at quadrature), obtaining 38 consecutive exposures during each of three visits to the source. Each visit started near an optical maximum and covered ${\Delta }{{\phi }_{B}}=1.27$ in phase. The spectrum orbital phase is set by the exposure midpoint at the solar system barycenter.

The data were processed using the Gemini IRAF package. Fluxes were established using observations of the A0 standard CD 32-9927. Comparison-lamp exposures were only taken once per visit, but, as for the LRIS spectroscopy, we monitored velocity stability and slit losses using the nearby comparison star. This also facilitated comparison between observing runs. The first visit, on 2013 May 6 (MJD 56,418.15328–56,418.23601), featured excellent ∼0farcs5 seeing and high throughput, while the next two visits, 2013 June 3 and 5 (MJD 56,446.04081–56,446.12281 and 56,448.05900–56,448.14101), had more typical ∼0farcs8 image quality and lower S/N.

In Figure 3 we show the average spectrum from maximum light (${{\phi }_{B}}=0.65$–0.85, Day, all GMOS-S orbits) and near minimum light (${{\phi }_{B}}=0.15$–0.35, Night, all GMOS-S orbits), Doppler corrected to rest using the radial-velocity solution (see below). For comparison, the upper trace shows an He-dominated model atmosphere with ${{T}_{{\rm eff}}}=12,000$ K and ${\rm log} \;g=4.5$ (Jeffery et al. 2001). In the spectrum from J1311 maximum, many narrow absorption features of He i and low-excitation metal lines are seen. In contrast, the average spectrum at minimum shows broad emission features, with He i and Mg i most prominent. Unlike the red LRIS spectra, these do not strongly dominate the 4000–5000 Å continuum, so we cannot independently measure the wind emission line profile and variation. In the trailed spectra this line emission is most obvious around ${{\phi }_{B}}\approx 0$ (centered on "Dusk" phases) during the first visit. It is weaker but present in the rest of the GMOS-S data.

Figure 3.

Figure 3. Average spectra near maximum light (${{\phi }_{B}}=0.65$–0.85, black curve) and minimum light (${{\phi }_{B}}=0.15$–0.35, green curve) from the GMOS-S data. A comparison He-dominated model atmosphere (magenta curve) shows that many strong, narrow absorption features are present in the maximum-light "Day" spectrum; some of the strongest lines are marked. Broad emission is apparent at minimum light.

Standard image High-resolution image

2.3. Other Spectroscopy

Since the wind features sometimes show dramatic changes from one orbit to the next, we also describe several other observations of J1311, designed to check the system's status. In R12, we discussed our initial six consecutive 300 s LRIS exposures covering ${{\phi }_{B}}=0.56$–0.93 on 2012 May 17 (MJD 56,064.306–56,064.331). These data, obtained with the same configuration as on 2013 March 11, established the companion He dominance and provided an estimate ${{K}_{{\rm obs}}}=609.5\pm 7.5$ km s−1 for a sinusoid fitting the CoL radial velocity. At this epoch J1311 appeared to be in quiescence, with no emission lines.

We also observed J1311 on 2013 April 13 (MJD 56,390.424–56,390.433) for 3 × 300 s near optical maximum (${{\phi }_{B}}=0.75$) using the DEep Imaging Multi-Object Spectrograph (DEIMOS; Faber et al. 2003) on Keck II. DEIMOS provided coverage of 4450–9630 Å, with the dichroic set at 6990 Å. Standard extraction, flux calibration, and cleaning were applied. Observations were made at the parallactic angle (Filippenko 1982) and did not include the comparison star.

The companion was in a particularly bright state with a dramatic flare during the second spectrum (centered at ${{\phi }_{B}}=0.835$). During this flare the continuum flux increased by a factor of ∼2.5. He ii absorption appears quite strong, particularly in the second spectrum. Figure 4 shows the dramatic continuum variability over ∼300 s at this epoch. Note that the difference spectrum at bottom shows that the emission-line variations, if any, are quite small. We thus infer that this heating is applied below the photosphere, allowing relatively narrow absorption features to appear. The broader superimposed emission lines arise from a different location.

Figure 4.

Figure 4. Sequential spectra from the blue side of the Keck II/DEIMOS observations. He i emission from the wind is present in all spectra. The flux doubles during the second exposure (green curve), and the residual spectrum (cyan, at bottom) shows the absorption features after subtracting the continuum power law: He ii absorption is strong. This flare component is only weakly present in the first spectrum (red curve), and no He ii absorption is apparent during the third spectrum (blue curve) 600 s later.

Standard image High-resolution image

Finally, we returned to J1311 with Keck I/LRIS on 2013 May 10 (MJD 56,422.347–56,422.355) for three additional 300 s exposures, this time covering phases ${{\phi }_{B}}=0.012$–0.139. The source appears to have returned to quiescence with faint line emission.

In Table 1 we summarize these various spectroscopic campaigns. To give some idea of the state of the wind line emission, we list the flux of the He i λ6678 line (we measure the He i λ4471 line for GMOS-S spectra), after averaging over the full visit. Of course, as apparent from Figure 2, the line strength is highly variable over even a single orbit, so these fluxes should be considered only a rough characterization of the wind activity. For example, although the line emission was not detected during the average of the second and third GMOS-S visits, we can detect He i emission during the Night $({{\phi }_{B}}\approx 0.25$) portion of the orbit (Figure 3).

Figure 5 compares the $\lambda {{F}_{\lambda }}$ spectra near ${{\phi }_{B}}=0.75$ from these five data sets. The GMOS observation and the LRIS/56,064 spectrum (here and elsewhere, the label refers to MJD) have fluxes similar to that of LRIS/56,362, but have been offset for readability. We also offset the LRIS/56,422 spectrum by −0.75; note that, unlike the other spectra in this figure, this is at quadrature (${{\phi }_{B}}=0.01$) and thus is intrinsically fainter. Two DEIMOS/56390 spectra near maximum light are also shown: phase ${{\phi }_{B}}=0.77$ before the main flare (only $\lambda \gt 5000$ Å plotted to avoid overlap) and the flare peak at ${{\phi }_{B}}=0.84$. By ${{\phi }_{B}}=0.90$, J1311 returns to the quiescence level. All quiescent spectra have similar continuum fluxes, but the He i lines vary greatly: absent on 56,064, strong on 56,362, and a factor of ∼2 stronger on 56,390.

Figure 5.

Figure 5. Spectra near optical maximum (${{\phi }_{B}}=0.65$–0.85) from our five observing campaigns (grouping the three GMOS-S visits). The GMOS-S (cyan) and LRIS/56,064 (blue) spectra have very similar flux to that of our long LRIS/56,362 campaign (black) and are offset for readability. The LRIS/56,422 (magenta) is also offset by −0.75, but is at "Dusk" phase and so is intrinsically fainter. The other spectra have flux as plotted. Note the dramatic continuum increase during the DEIMOS/56,390 flare (green), during which He ii absorption is prominent. Also note the strong and variable He i emission in the red.

Standard image High-resolution image

3. THE WIND KINEMATICS

The broad He i emission, which must arise in the companion wind, is clearly episodic and variable. Our best measurement of this emission is in the 2013 March Keck/LRIS data (MJD 56,362). Since the lines appear optically thin with relatively constant line ratios, we have formed a wind-velocity profile by taking the weighted, continuum-subtracted average of the three strongest He i lines in the red: 5875.61, 6678.1, and 7065.18 Å. The result (Figure 6) shows the flaring, varying emission-line profile. Scaled versions of this line profile provide good matches to the emission-line structure in the blue, as well.

Figure 6.

Figure 6. Velocity profile of the He i emission from the companion wind on MJD 56,362. For reference, the red curve shows the radial velocity of the simple sinusoidal fit to the absorption features, marking the CoL velocity of the heated surface of the companion star. The emission lines do not track this motion, but a likely sinusoidal component is shown by the dotted line. Over three periods are shown, with dark colors marking bright line emission.

Standard image High-resolution image

In Figure 6, although the emission is patchy, there is approximate orbital periodicity, with one component following the dotted sinusoid. Note that the emission other than this sinusoidal component seems to generally lie redward of the companion CoM velocity. We conclude that some preferred vectors of emission in the rotating system exist, but that the emission is intermittent along these vectors.

Figure 7 show a schematic model illustrating a possible kinematic origin for the wind line structure. Here we compute the cometary wind trail for matter stripped at low velocity (in the corotating frame) from the companion and then accelerated by pulsar radiation pressure; the only free parameter in this model is ${{a}_{{\rm PSR}}}/{{g}_{{\rm PSR}}}$, the ratio of wind acceleration to gravity. A similar model with the addition of a finite launch speed of 2/3 ${{v}_{{\rm orb}}}$ does a good job of defining the radio eclipse envelope for the original black widow PSR B1957+20 (Rasio et al. 1989). Phase ${{\phi }_{B}}=0.3$ is shown in the upper panel, with the arrows and histogram indicating the possible radial velocities along this wind trail; a patch of He-containing gas will provide an emission line within this envelope. The lower panel illustrates the range of velocities as a function of phase. This figure can be compared with Figure 6; we reproduce the two curves from that figure to guide the comparison. The wind trail plotted fades after ${\Delta }{{\phi }_{B}}=0.7$; this boundary matches well in both amplitude and phase the quasi-sinusoid traced in Figure 6 by the dotted curve. If the wind emission is relatively strong at lag ${\Delta }{{\phi }_{B}}=0.7$, patchy at smaller ${\Delta }{{\phi }_{B}}$, but stronger at ${{\phi }_{B}}\approx 0.75$–0.25, this pattern provides a reasonable match to the data. Many orbits would need to be monitored to determine if this represents a dominant persistent emission pattern.

Figure 7.

Figure 7. Simple ballistic model for the cometary trail emission. Mass is lost from the companion and accelerated radially at $3/4\;{{g}_{{\rm PSR}}}$, via the pulsar spin-down flux. One phase, with the trail (out to a lag ${\Delta }{{\phi }_{B}}=0.7$), along with the line profile (histogram) for uniform illumination, is shown at top. The red arrows show the companion velocity in space and at the red edge of the line profile, for ${{\phi }_{B}}=0.3$. Two periods of the trail velocity range, along with the companion velocity and the sinusoidal component curves (curves exactly as in Figure 6), are shown in the lower panel.

Standard image High-resolution image

This kinematic model is attractive, but the present data do not allow us to determine the physical origin of the emission-line velocities; of course, any simple sinusoid can be reproduced by a rotating vector fixed in the binary frame, for example, at the companion star. The dotted curve shown can be produced by a vector (i.e., a jet) from the companion with ${{v}_{j}}\approx 1.3{{v}_{{\rm orb}}}$ in the star's rest frame, directed in the orbital plane ∼60° ahead of the line of centers. The physical nature of the pulsar flux driving the trail is also unclear. With $\dot{E}=5\times {{10}^{34}}\;{{I}_{45}}\;{\rm erg}\;{{{\rm s}}^{-1}}$, ${{a}_{{\rm PSR}}}/{{g}_{{\rm PSR}}}=3/4$ requires a mean wind cross section of $\langle \sigma /m\rangle \approx 1.5\times {{10}^{3}}\;{\rm c}{{{\rm m}}^{2}}\;{{{\rm g}}^{-1}}$. This is about right for a ionized pair plasma, but high for a baryonic wind. It is also interesting to consider what this wind tells us about pulsar evaporation. The momentum density of the pulsar wind is $40\;{{I}_{45}}\;{\rm g}\;{\rm c}{{{\rm m}}^{-1}}{{{\rm s}}^{-2}}$ at the orbital radius. To ram-pressure balance this wind along the line of centers requires a companion wind density ${{\rho }_{{\rm W}}}=1.0\times {{10}^{-14}}\;{{I}_{45}}v_{600}^{-2}\;{\rm g}\;{\rm c}{{{\rm m}}^{-3}}$ and a total mass flux ${{\dot{m}}_{{\rm W}}}=\pi R_{*}^{2}{{\rho }_{{\rm W}}}{{v}_{{\rm W}}}\approx 6\times {{10}^{13}}\;{{I}_{45}}v_{600}^{-1}\;{\rm g}\;{{{\rm s}}^{-1}}$, for a wind speed $600\;{{v}_{600}}\;{\rm km}\;{{{\rm s}}^{-1}}$, comparable to the orbital speed and wind only from the heated face of the companion. Note that this gives an evaporation timescale $\tau \approx {{m}_{2}}/{{\dot{m}}_{{\rm W}}}\approx 10\;{{m}_{-2}}{{v}_{600}}/{{I}_{45}}$ Gyr, for ${{m}_{2}}={{10}^{-2}}\;{{m}_{-2}}\;{{M}_{\odot }}$. Although m2 was larger when accretion stopped and evaporation began, $\dot{E}$ was higher in the past so the companion may be completely ablated in less than a spin-down time, especially for a slow wind or ${{I}_{45}}\gt 1$.

4. THE PHOTOSPHERE FLARE VARIABILITY

Figures 4 and 5 show that in addition to the emission-line variability, there can be epochs of true continuum increase. The dramatic ${{i}^{\prime }}$ flare shown in Figure 1 of Romani (2012), which erupts to 6 times the normal flux at ${{\phi }_{B}}\approx 0.8$, is likely a similar event. Figure 4 gives an important clue to the nature of these variations. The difference spectrum shows that the flare does not affect the emission-line fluxes, but produces strong absorption lines. The He ii absorption dramatically increases (including a very strong, unusual detection of He ii λ6560, generally confused with the absent Hα line). We infer that these flares must be caused by impulsive deep heating well below the photosphere. However, since the photometric flares and spectral flares have decay times of ∼300 s, if we assume a sound speed $\sim 8\;T_{4.5}^{1/2}\;{\rm km}\;{{{\rm s}}^{-1}}$ in the ∼30,000 K He gas seen during these bursts, we infer a characteristic depth of $\sim 2\times {{10}^{8}}$ cm ($\sim 0.03\;{{R}_{*}}$).

By measuring the He ii/He i line-intensity ratios in the flare spectrum, compared with the ratios in model spectra of Jeffery et al. (2001), we find an effective flare temperature of ${{T}_{{\rm Fl}}}$ = 39,000 ± 1000 K. The continuum flux increase from quiescence is a factor of 2.2 at 5500 Å assuming Planckian continua in this range. This implies a flare-emitting area 0.27 times that of the quiescent star. Thus, we infer that a large fraction of the star brightened by a factor of 8.3 during the eruption. At the estimated (from the dispersion measure) distance of 1.4 kpc, the equivalent radius of the flare area is 0.05 R$_{\odot }$, which is $\sim 0.5\;{{R}_{*}}$, in good agreement with the spectral estimate.

The average thermal luminosity of this burst is $\sim 5\times {{10}^{33}}\;{\rm erg}\;{{{\rm s}}^{-1}}$. This is ∼80 times the pulsar flux hitting the companion, so we infer that these outbursts represent stored energy released deep below the photosphere. The natural candidate energy source is a magnetic field, with characteristic magnitude $\langle B\rangle \approx 10$ kG required to power the flare. The origin of such a field is unclear. Typically stars with ${{T}_{{\rm eff}}}$ = 12,000 K are radiative throughout. However, here we expect that high-energy particles from the pulsar are deeply heating the photosphere, with (for example) TeV ${{e}^{\pm }}$ depositing their energy at ${\Sigma }\approx 300\;{\rm g}\;{\rm c}{{{\rm m}}^{-2}}$. While this does not in itself induce convection, we expect latitudinal flow because of the large temperature gradients from the differential heating. If convective, given the short effective spin period for the tidally locked companion, we have the ingredients needed for a robust dynamo, large magnetic fields, and energetic flare events.

5. IMPROVED QUIESCENCE LIGHT CURVES FROM SPECTRA

The large equivalent width of the variable He i lines can certainly affect the broad band fluxes, especially at binary minimum and for red colors. This is likely the origin of the red "flaring events" noted by R12 in simultaneous multi-band GROND light curves. We can use our spectra to partly correct for this effect, synthesizing light curves from the continuum. To do this, we integrate our normalized spectra over the SDSS ${{g}^{\prime }}{{r}^{\prime }}{{i}^{\prime }}{{z}^{\prime }}$ filter passbands, after excising ±1500 km s−1 wide bands around the strong wind lines (mostly He i, but also Ca and Mg) seen at binary minimum (Figure 6). To normalize these filter bandpasses to standard SDSS magnitudes, we integrated the comparison-star spectra over the same bands. Carefully calibrated GROND observations give magnitudes for this star of ${{g}^{\prime }}=18.58$, ${{r}^{\prime }}=18.05$, ${{i}^{\prime }}=17.79$, and ${{z}^{\prime }}=17.70$ mag (A. Rau, private communication). This determines an offset to the measured magnitudes for each spectrum. Comparing with synthesized magnitudes for the full SDSS passbands, we do indeed find decreased source fluxes, by as much as a magnitude when the wind lines are strong.

However, even this line excision does not completely remove the light-curve fluctuations. To best estimate a "quiescence" light curve, we use our spectra from MJD 56,362 and compare the synthesized magnitude estimates in ${\Delta }\phi =0.05$ phase bins. We adopt the largest (faintest) magnitude in each band in each bin as the quiescence value and plot the resulting light curves in Figure 8. Statistical errors are obtained from the fluctuations in the spectroscopic fluxes integrated to derive the band magnitudes. We further estimate the remnant light-curve fluctuations associated with faint wind lines or continuum flares, by taking the difference to the second faintest magnitude observed for a given bin. These error estimates are added in quadrature to the statistical errors and plotted on the bin points in Figure 8.

Figure 8.

Figure 8. Synthesized band magnitudes from the Keck MJD 56,362 spectra. These ${{g}^{\prime }}{{r}^{\prime }}{{i}^{\prime }}{{z}^{\prime }}$ magnitude estimates are formed by integrating over the SDSS passbands after excising strong (mostly He i) emission lines. Magnitudes are normalized by matching the synthesized comparison-star magnitudes to the direct imaging values. Two periods are shown with estimated errors during the second cycle (see text).

Standard image High-resolution image

Several features of the derived light curves deserve comment. First, our synthesized magnitudes, normalized to the comparison star, are roughly 0.6 mag fainter than those reported from direct GROND imaging in R12. It is possible that this apparent 40% flux decrease is caused by differential slit losses between the pulsar and comparison star. However, the GROND magnitudes at MJD 56118 may have been well above quiescence, since these data show larger fluctuations (even at ${{\phi }_{B}}\approx 0.75$) than seen in these synthesized magnitudes. Of course, even these new light curves will not represent complete quiescence—fluctuations are still visible across the maximum—but they provide a substantially improved target for quiescence modeling. Overall the light curves follow the form of the direct imaging SOAR/GROND curves reported by R12, although we appear to have better captured the red colors at binary minimum. It should be noted that as in R12, the light curves remain slightly asymmetric, brighter at "Dusk" (${{\phi }_{B}}\approx 1.0$) by ∼0.1 mag compared to Dawn (${{\phi }_{B}}\approx 0.5$).

6. THE COMPANION RADIAL VELOCITY

The dramatic heating offsets the CoL in a poorly understood, temperature-dependent way. Also the variable wind emission in precisely the He i lines dominating the companion photosphere can shift the absorption-line centroids. These effects make measurement of the companion radial velocity challenging. We focus here on the 2013 March Keck/LRIS and the 2013 May–June Gemini/GMOS-S data sets. For the Keck data we start by removing a scaled version of the emission-line velocity profiles (Figure 6) from each detected line feature in the blue channel. We then use an He-dominated model atmosphere template from the atlas of Jeffery et al. (2001) (${{T}_{{\rm eff}}}$ = 12,000 ± 1000 K, log $g=4.6\pm 0.2$) as a template for velocity cross-correlation to obtain the best radial velocity of each individual spectrum. Cross-correlation is performed using the IRAF RVSAO package (Kurtz & Mink 1998). The quality of the fit is monitored by the correlation statistic R (the ratio of the correlation peak height to random peak height), which varies from ∼12 on the day side, where velocity uncertainties are as small as ${{\sigma }_{v}}=15$ km s−1, to an insignificant $R\lt 1$ on the Night side. Measurement of the comparison star with a G-star template showed a slow drift in radial velocity from +30 to −30 km s−1, with typical uncertainties of 6 km s−1. We accounted for this residual calibration drift and also corrected to heliocentric velocities.

Similar measurements were made of the GMOS-S spectra. Here we reach $R\approx 7$ during the Day phases. With the high resolution, this delivered velocity uncertainties as small as ${{\sigma }_{v}}\approx 5$ km s−1. Night-time cross-correlations were often insignificant, but occasionally delivered apparently significant $R\approx 1$–2; some of these measurements were close to the orbital curve, while others had large departures. The comparison-star radial velocity was monitored to typically ${{\sigma }_{v}}\approx 4$–6 km s−1, although uncertainties grew to as large as ∼10 km s−1 when conditions were poor. No strong velocity trend was evident in these data.

The corrected cross-correlation velocities for these two data sets are shown in Figure 9. Two periods are illustrated, phased to the orbital solution of Pletsch et al. (2012). During the first, all points with $R\gt 1$ are plotted, with their ${{\sigma }_{v}}$ error flags. The radial-velocity curve during the Day phase is very well measured. During the Night phase a number of GMOS-S points appear to follow the expected curve, but many are displaced by ∼+30–60 km s−1; the few plotted LRIS points in this phase range seem largely random. Recalling that the GMOS-S data measure only the blue portion of the spectrum, it may be that the radial velocity is tracking flux from the bright lune of the heated side, visible at the companion limb for $i\lt 90{}^\circ $.

Figure 9.

Figure 9. J1311 radial-velocity measurements from LRIS and GMOS-S data. The first cycle shows all $R\gt 1$ correlation points, with velocity error flags. During the second cycle the data used in the simple sinusoidal fit (solid line; dashed line marks the systemic velocity Γ) are in green and cyan; red points are excluded. The bottom panel shows the residuals to this best-fit sinusoid. The curve plots residuals to the ELC velocity curve ("ELC-Sp" model) after subtracting the simple sinusoid, demonstrating that the model's non-sinusoidal radial-velocity terms are present in the data (see text).

Standard image High-resolution image

In addition to the low-statistics Night-time points, it is clear that there are large variations from the radial-velocity curve at the Dawn and, especially, Dusk phases. At times the velocities follow the simple sinusoid; at others the velocities lie up to ∼100 km s−1 "inward" of the expected curve, to the blue ($-v$) near ${{\phi }_{B}}=0.5$ and to the red ($+v$) near ${{\phi }_{B}}=0.0$. These offsets are highly significant and do correlate among spectra. However, we are not able to associate these with, say, periods of strong wind emission lines. One interpretation is that these offsets represent radial outflow near the companion terminator; an emission component to the red (blue) at Dawn (Dusk) would impose a variable blueshift (redshift) on the absorption features measured against the cross-correlation template. More high-quality data are needed to pin down the physical origin of this effect.

One consequence of the terminator shifts is a difficulty in fitting a simple radial-velocity curve. If we fit only to correlation measurements with $R\gt 3$ and restrict to ${{\phi }_{B}}=0.5$–1.0, we obtain an observed radial-velocity amplitude ${{K}_{{\rm CoL}}}=611.4\pm 5.4$ km s−1 and systemic velocity ${\Gamma }=39.4\pm 3.7$ km s−1 (statistical errors only). The fit departures remain substantial, with a ${{\chi }^{2}}$ per degrees of freedom (dof) of 9.8 for 74 dof. If we restrict further to the Day side ($0.55\lt {{\phi }_{B}}\lt 0.95$), we obtain ${{K}_{{\rm CoL}}}=615.4\pm 5.1$ km s−1, with ${{\chi }^{2}}$/dof = 4.0. In contrast, fitting all points with $R\gt 3$ includes many discrepant Dusk points and gives ${{K}_{{\rm CoL}}}=590.7\pm 7.2$ km s−1, with ${{\chi }^{2}}$/dof = 15.6. It is evident that non-sinusoidal terms are present. The bottom of Figure 9 shows the residuals to the simple sinusoidal fit. Despite the scatter, some systematic trend is visible, especially during the Day side. The curve plotted, for a best-fit ELC model with an L1 cool region (see below), shows a similar residual to the fit sinusoid, and the fit to this curve decreases ${{\chi }^{2}}$ by 53.7, i.e., ${{\chi }^{2}}$/dof decreases from 9.8 to 9.1.

7. MODELS OF THE HEATED PHOTOSPHERE

There are a number of codes used to model light curves and radial velocities in interacting binaries. We find that none are fully adequate to cover the extremes of temperature, surface composition, and pulsar heating found in J1311. Nevertheless, such modeling fails in instructive ways and provides a guide toward future, more detailed solutions. We report here on modeling with the ELC code (Orosz & Hauschildt 2000) and the ICARUS code (Breton et al. 2012), both kindly made available to us by their authors. For both we can run in "MSP" mode, with the period PB = 0.065115 days and pulsar orbit ${{x}_{1}}={{a}_{1}}{\rm sin} i=0.010581$ lt-s held fixed. Both codes model the effect of pulsar heating of the companion as illumination by a point (X-ray) source.

For the ICARUS code we ran with color tables generated from the BTSettl atmospheres (Allard et al. 2011). In all cases the best-fit models preferred large Roche-lobe fill factors; we set the effective fill factor to 0.99. In fitting our quiescent ${{g}^{\prime }}{{r}^{\prime }}{{i}^{\prime }}{{z}^{\prime }}$ light curves, the main parameters are the orbital inclination i, the Night side (unheated) star temperature TN, and the heated side (Day) temperature TD. We report this last parameter as an effective heating luminosity ${{L}_{H}}=4\pi {{a}^{2}}\sigma (T_{D}^{4}-T_{N}^{4})$, where the orbital separation is $a\approx 0.95\;{{R}_{\odot }}$ and the heating luminosity can be compared with the pulsar spin-down power $\dot{E}=5\times {{10}^{34}}\;{{I}_{45}}\;{\rm erg}\;{{{\rm s}}^{-1}}$. ICARUS also fits the absolute fluxes and so returns a distance estimate and an effective extinction. If left free, the extinction fits to an insignificantly small value; to be conservative, we set the value to the low net Galactic extinction AV = 0.173 mag in this direction (Schlafly & Finkbeiner 2011). The fit values are given in Table 2, with uncertainties as reported by the ICARUS fitter. The distance is surprisingly large at 2.6 kpc, and the heating flux slightly exceeds the standard pulsar spin-down power. Most importantly, the fit inclination is so small that with the observed radial-velocity amplitude ${{K}_{{\rm CoL}}}=615$ km s−1, the minimum pulsar mass would be $2.77\;{{M}_{\odot }}$. Since the true kinematic amplitude is larger (typically with ${{K}_{{\rm cor}}}\gt 1.03$, Table 2), then ${{K}_{{\rm CoM}}}\gt 633$ km s−1 and the inferred mass would be unphysically large at $\gt 3\;{{M}_{\odot }}$. Also, the fit residuals (Figure 10) demonstrate that the maximum is substantially wider than the model light curve, and that the colors are redder at peak than predicted by the model. Clearly, the present ICARUS model is inadequate for a detailed fit.

Figure 10.

Figure 10. Residuals to the model fits. Top: photometric (${{g}^{\prime }}{{r}^{\prime }}{{i}^{\prime }}$) residuals with (left to right) the ICARUS fit, the ELC baseline fit, and the ELC fit with a $0.6\times {{L}_{1}}$ spot. Bottom left: radial-velocity points used in the ELC fitting. Middle and right panels show the residuals to the ELC baseline and L1 spot fits, respectively.

Standard image High-resolution image

Table 2.  J1311–3430 Light Curve/Spectral Fits

Parameter Units ICARUS ELC ELC-Sp ELC-Sp/iq
Incl. Deg. 57.1 ± 1.5 64.0 ± 3.5 $81_{-6}^{+4{\rm b}}$ $77_{-4}^{+3}$
TN K 4240 ± 150 4800 ± 400 5000 ± 400 4500c
${\rm Log}({{L}_{H}})$ erg s−1 34.80 ± 0.03 $35.5_{-0.3}^{+0.15}$ 35.6 ± 0.2 35.0c
q $184.5\pm 1$ $179.5_{-0.5}^{+1}$ $180\pm 0.5$
${{K}_{{\rm cor}}}$ 1.062 1.033 1.036
$M_{1}^{{\rm a}}$ M$_{\odot }$ 2.63$_{-0.2}^{+0.3}$ $1.82_{-0.06}^{+0.15}$ $1.88_{-0.06}^{+0.09}$
$M_{2}^{{\rm a}}$ ${{10}^{-2}}\;{{M}_{\odot }}$ 1.42 ± 0.15 $1.01_{-0.02}^{+0.09}$ $1.04_{-0.03}^{+0.05}$
${{\chi }^{2}}$/dof 148 3.74 2.60 2.66

Notes.

aValue inferred from fit—projected uncertainties, errors inflated by ${{\chi }^{2}}$/dof. bBounded by lack of an eclipse. cFixed in fit, with spot depth $0.6\times $.

Download table as:  ASCIITypeset image

The ELC code was used to fit the light curves and spectra simultaneously. For these fits we used an atmosphere model table computed for SDSS filters by J. Orosz, extracted from the "NextGen" atmosphere grid with an extension to ${{T}_{{\rm eff}}}\;\gt $ 10,000 K using Castelli & Kurucz (2004) models. As before, the fits preferred large Roche-lobe fill factors and this value was fixed at 0.99. We again fit for the inclination, TN, and heating, but additionally constrained the mass ratio $q={{M}_{1}}/{{M}_{2}}$ for direct estimates of the component masses. The results, which include radial-velocity constraints, prefer larger inclinations and lower masses than the ICARUS solution. The TN temperature is slightly higher, but most importantly the fit LH is much higher and is indeed substantially larger than the total spin-down power. Again, the observed maximum is appreciably wider than the heating model suggests, and the mid-Day (ϕ = 0.25) temperature is cooler than predicted. The fit is not statistically satisfactory and implies a large pulsar mass ${{m}_{1}}=2.6_{-0.2}^{+0.3}$. Here, and for all ELC fits, we quote "$1\sigma $" errors as the full multi-dimensional projection over all fit parameters of the region about the fit minimum with a fit statistic increase of $\chi _{{\rm min} }^{2}/{\rm dof}=1$; in this sense our uncertainties are conservative, with marginalization over all other parameters and error inflation for imperfections of the model.

The ELC code has the option of applying artificially heated or cooled star spots to the secondary. As in R12, we can substantially improve the fit by cooling the inner Lagrange point. We apply a spot that decreases the local unperturbed surface temperature with a Gaussian profile of radius 45°. The best fits invoke a ∼40% temperature decrease and do give a substantially better model for the light-curve maximum. The result is that the companion is brightest in a ring pulsar-ward of the terminator, with the CoL much closer to the CoM than for direct heating (small ${{K}_{{\rm cor}}}$). In addition, the heating pattern broadens the light-curve maximum. This drives the inclination i to larger values, to maintain the depth of the light-curve minimum. Both of these effects serve to decrease the inferred pulsar mass.

While there are physical effects (see below) that can accommodate a larger effective heating power, ${{L}_{H}}={{10}^{35.6}}\;{\rm erg}\;{{{\rm s}}^{-1}}$ is much larger than the pulsar can supply. Also, the large TN gives a poor match to the relatively red colors and large magnitude at ${{\phi }_{B}}=0.25$. The drive to these large values is a consequence of the spot approximation for the surface-temperature distribution. If we fix ${{L}_{H}}={{10}^{35}}\;{\rm erg}\;{{{\rm s}}^{-1}}$ and TN = 4500 K, the largest values consistent with the observations, we find that the inferred inclination (and thus neutron-star mass) remains a very strong function of the spot temperature decrement (Figure 11), which is rather poorly determined. The best fit, with a 0.6 × spot (40% spot decrement), implies ${{M}_{1}}=1.9\;{{M}_{\odot }}$. Weaker spots lead to large inferred masses; stronger spots require edge-on orbits in conflict with the lack of an X-ray eclipse (Romani 2012). Thus, our conclusions depend sensitively on details of the heating model. Notice in Figure 10 that the radial-velocity residuals vary only weakly between quite different heating patterns. The main effect is an overall increase in amplitude, which is reflected in the changed ${{K}_{{\rm cor}}}$; the temperature (color) is more sensitive to the heating model.

Figure 11.

Figure 11. Confidence regions for two-parameter ELC model fits, holding LH and TN fixed near their maximum values and varying the assumed spot deficit from $1.0\times $ (no spot) to 0.5× (deep spot, excluded by the lack of an eclipse). The best-fit point and the ${{\chi }^{2}}$/dof are marked for each deficit. The uncertainty ellipses are for $\chi _{{\rm min} }^{2}+({{\chi }^{2}}$/dof) and $\chi _{{\rm min} }^{2}+3\times ({{\chi }^{2}}$/dof). Lines of constant M1 (solid, dashed) and M2 (dotted) are drawn. The inclination and mass ratio, and hence the inferred pulsar mass, depend dramatically on the heating pattern.

Standard image High-resolution image

8. CONCLUSIONS: INDIRECT HEATING, NEUTRON-STAR MASS

Our goal of measuring the orbital inclination and the ${{K}_{{\rm cor}}}$ factor with light-curve fitting has been frustrated by the poor fit to the light-curve maximum. A strong clue to the nature of the difficulty is that simple direct-heating fits require a heating power larger than can be supplied by the $5\times {{10}^{34}}\;{{I}_{45}}\;{\rm erg}\;{{{\rm s}}^{-1}}$ spin-down power. For the ICARUS model, LX exceeds the nominal spin-down power by a factor of 1.3. For the basic ELC model, the excess is a factor of 6, and for the ELC spot model, the heating power exceeds that supplied by the pulsar by a factor of 8. This is especially worrisome as the typical heating efficiency estimated by Breton et al. (2013) is $\eta =0.15$ of the pulsar spin-down.

This discrepancy is so large that a number of factors are likely relevant. First, the effective moment of inertia and resulting spin-down power may be larger than assumed, especially if M1 is $\gt 2\;{{M}_{\odot }}$ where ${{I}_{45}}\approx 2$–4, depending on the equation of state (Lattimer & Schutz 2005). Also, the present modeling assumes normal composition in the atmosphere; He-dominated spectra can be bluer, allowing a lower fit ${{T}_{{\rm eff}}}$. Most importantly, the assumption of direct, radiative heating in the ICARUS and ELC models is likely inadequate for J1311. Instead, we expect that the pulsar wind terminates at an intrabinary shock wrapped around the L1 point, reprocessing the pulsar wind into high-energy particles that deliver heat below the companion photosphere. This reprocessing can both increase the effective heating flux over direct radiation, as the intrabinary shock subtends a larger angle at the pulsar, and preferentially heat the companion terminator. Further, Coriolis forces break the shock symmetry in the rotating frame and may thus provide a mechanism for explaining the asymmetric light-curve maximum. The presence of strong companion fields at the L1 point (indicated by the large flares) may also preferentially direct intrabinary shock particles away from L1. Such effects are dramatic for J1311, but may be present at lower levels in other evaporating pulsar binaries, where high-precision photometry shows asymmetric light curves (e.g., Schroeder & Halpern 2014). This may bias mass estimates inferred from fitting light curves with simple direct-heating models.

Our spectroscopic campaign has revealed a number of new aspects of this remarkable system. First, J1311 undergoes violent flaring in the optical, which suggests strong magnetic activity on the heated face of the companion. Second, strong (but intermittent) line emission is generated in the wind driven off the companion. This outflow spectrum is dominated by neutral helium lines and appears to follow preferred paths corotating with the binary system. Since the line flux along these paths is sporadic, further study is needed to determine the persistent patterns, but emission along a pulsar-wind-driven spiral outflow provides a viable model. Selecting against companion flares and spectrally excluding the outflow line flux, we have made improved measurements of the companion photosphere's heated light curve and color variation. These data strengthen the conclusion that simple radiative heating models are not adequate to describe strongly interacting black-widow binaries like J1311. Instead, we posit that the heating flux is reprocessed in an intrabinary shock whose illumination of the companion surface preferentially heats away from L1, toward the Day–Night terminator. If this deep heating is mediated by shock-accelerated particles, coherent magnetic fields near L1 can also help redirect the flux and decrease radiation from the L1 point.

The net effect of asymmetric heating is to prefer smaller ${{K}_{{\rm cor}}}$ and larger sin i, but the precise values of these critical parameters are very sensitive to the heating pattern. There are two important paths to improving their determination. (1) Improved, more physical modeling can better exploit existing measurements. However, statistically acceptable ${{\chi }^{2}}$ will probably require excellent exclusion of flare perturbations, even if the heating model is physically correct. (2) Observationally, one would like to directly measure the heating pattern across the companion Day side. Improved spectroscopy remains the key. While well-measured colors, especially near minimum (${{\phi }_{B}}=0.75$), can help constrain the fits, we have seen in this study that these can best be obtained by synthesizing colors from spectra, after excluding flares and the intermittent wind emission lines. Of course, with sufficient S/N, the spectra themselves will help determine the ${{T}_{{\rm eff}}}$ distribution as a function of phase, by measuring relative absorption-line strengths for species with differing excitation. Even more directly, higher spectral resolution with good S/N can measure rotational broadening of the companion's absorption lines and its variation with orbital phase, thereby providing perhaps the cleanest constraints on sin i.

Until such improvements are made, we cannot exclude any mass in the range 1.8–2.7 M$_{\odot }$ for PSR J1311–3430. This range of values still allows a wide range of system evolutions and may, or may not, constrain the equation of state at high densities. Thus, more work is needed to fully exploit the promise of this unusual system.

We thank Kelsey Clubb, Ori Fox, and Melissa Graham for assistance with some of the Keck observations, as well as German Gimeno for help with the Gemini campaign. We also thank the referee for many detailed comments. This work was partially supported by NASA grants NNX11AO44G and NNX12A068G. A.V.F. and S.B.C. were supported by Gary and Cynthia Bengier, the Richard and Rhoda Goldman Fund, the Christopher R. Redlich Fund, the TABASGO Foundation, and NSF grant AST-1211916. Some of the data presented herein were obtained at the W. M. Keck Observatory, which is operated as a scientific partnership among the California Institute of Technology, the University of California, and NASA; the Observatory was made possible by the generous financial support of the W. M. Keck Foundation.

Footnotes

  • J1653-0158 (Romani et al. 2014) may be a similar system with an even smaller Pb = 74.8 minute orbit, but it has not yet had a pulse detection.

Please wait… references are loading.
10.1088/0004-637X/804/2/115