Weather on Other Worlds. VI. Optical Spectrophotometry of Luhman 16B Reveals Large-amplitude Variations in the Alkali Lines

The 500K–2000K effective temperatures characteristic of L and T dwarfs enable the formation of condensate clouds of types not observable in our solar system — including clouds of silicates and iron vapor. Such clouds have a profound influence on the spectra of brown dwarfs, and indeed there is a range of spectral types near the L-T transition (roughly L8–T5; Kirkpatrick et al., 2000) where the effective temperature is believed to be nearly constant and the spectral sequence is due to decreasing thickness of cloud cover as we move from the L dwarfs to the T dwarfs (Ackerman & Marley, 2001).

If clouds thick enough to affect the spectrum are distributed inhomogeneously in longitude over the atmosphere of a brown dwarf, rotationally modulated photometric and spectroscopic variability should result (e.g., Artigau et al., 2009). Studying such variability in detail can reveal the properties of the clouds. In particular, the variation in photometric amplitude as a function of wavelength can be diagnostic of the height of the cloud tops. In spectral regions with very low gas opacity (i.e., the bright continuum away from strong spectral absorption lines) regions with thick clouds extending to high altitudes are expected to be much darker than regions with thin, low-altitude clouds, because in the latter areas, the low gas opacity allows us to observe strong flux from deep, warm layers of the atmosphere. Hence, the less cloudy regions appear very bright in contrast to the dark, cloudy regions, and high-amplitude rotationally modulated variations are the result. By contrast, in spectral regions of high gas opacity (i.e., strong absorption lines), we cannot observe flux from deep in the atmosphere regardless of cloud conditions — so the contrast between thicker and thinner clouds is expected to be much less and the photometric amplitude is expected to be lower. This is in some respects an oversimplified picture — for example, the vertical temperature profile of the atmosphere (as well as the profiles of pressure and gas-phase elemental or molecular abundances) could be affected by the thickness of the clouds. Nevertheless, several observational results confirm the expectation of reduced photometric amplitude in regions of strong gas opacity — in particular the near-infrared H₂O absorption from 1.35–1.50$\mu$m (e.g., Apai et al., 2013; Buenzli et al., 2015a).

The binary brown dwarf system Luhman 16 (WISE J1049; Luhman, 2013) is composed of two objects with spectral types L7.5 and T0.5, respectively (Burgasser et al., 2013; Kniazev et al., 2013). With a distance of only 1.996 pc (Bedin et al., 2017), Luhman 16 is the third closest system to the Sun (after $\alpha$ Cen and Barnard’s star) and its two components are the nearest known brown dwarfs. The system has an orbital period of $27.5\pm 0.4$ years and a semimajor axis of $3.56\pm 0.03$ AU, with dynamically determined masses of $33.5\pm 0.3$ Jupiter masses for Luhman 16A and $28.6\pm 0.3$ for Luhman 16B (Lazorenko & Sahlmann, 2018).

Luhman 16B, the T0.5 secondary, is not only the closest T dwarf but also one of the most highly variable brown dwarfs known. Gillon et al. (2013) found it to vary with a period of about five hours and an amplitude exceeding 20% in the optical on some nights. Among periodically varying brown dwarfs, this amplitude is exceeded only by the 26% $J$-band variability detected by Radigan et al. (2012) in 2MASS J21392676+0220226, which was one of only two objects showing periodic variations with amplitudes greater than 10% among more than 150 brown dwarfs monitored in surveys by Koen (2013), Radigan et al. (2014), and ourselves (Metchev et al., 2015; Heinze et al., 2015).

The very high photometric amplitudes observed for Luhman 16B suggest that its clouds are unusually large, dark, and inhomogeneous. They also mean Luhman 16B offers an opportunity to make measurements of wavelength-dependent amplitude variations with unusually high signal-to-noise ratio. Many studies have been made of the photometric and spectroscopic variability in this remarkable system: Gillon et al. (2013); Biller et al. (2013); Crossfield et al. (2014); Buenzli et al. (2015a, b); Apai et al. (2021), to list just a few. There is broad agreement that the rotation period of Luhman 16B is about $5\pm 0.2$ hours, but finding a more precise value is difficult because the amplitude and light curve shape change dramatically on a timescale of just a few rotations (e.g., Gillon et al., 2013; Apai et al., 2021). The L7.5 primary Luhman 16A varies with an amplitude that is generally much smaller. This makes its rotation period more difficult to discern, but it has been estimated at about eight hours (Mancini, 2015; Apai et al., 2021).

Spectrally resolved studies of Luhman 16B’s variability have found significantly reduced photometric amplitude in the near-infrared water absorption band from 1.35–1.50$\mu$m (Buenzli et al., 2015a; Kellogg et al., 2017), consistent with the expectation discussed in Section 1.1 above. At the same time, Buenzli et al. (2015b) did not find evidence of an amplitude change in the iron hydride absorption band at 0.99 $\mu$m, which they interpreted as indicating that there are no completely unclouded regions in the atmosphere of Luhman 16B: the photometric variations are produced by variations in cloud thickness rather than by the complete absence of clouds in some areas. Gradations in cloud thickness as opposed to a stark cloudy vs. cloudless dichotomy would also be consistent with the results of Crossfield et al. (2014), who used Doppler imaging based on high-resolution spectra in the 2.288–2.345$\mu$m regime to construct remarkable, disk-resolved maps of Luhman 16B. The maps show a total variation in photospheric surface brightness of only about 10%. Large, coherent features would be required to sustain several-percent photometric variability with such a limited photospheric contrast range, and the Crossfield et al. (2014) maps show such features — in particular, a kidney-shaped dark spot spanning about 70^∘ in latitude and only slightly less in longitude.

We present a new spectrally resolved analysis of Luhman 16B, using Gemini/GMOS observations to probe shorter wavelengths than previous studies. In particular, our data cover the extremely broad, strong lines of the neutral alkali metals: the potassium doublet centered near 7682Å and the sodium doublet at 5893Å. In Section 2 we describe our observations and data reduction, and in Section 3 we validate our GMOS spectrophotometry by comparison with simultaneously obtained broadband imaging photometry. In Section 4, we follow Apai et al. (2013); Buenzli et al. (2015a, b); and Kellogg et al. (2017) in exploring spectrally resolved variability by constructing ratios of spectra obtained in bright vs. faint photometric states. Section 5 describes the dramatic amplitude-reversal these spectral ratios reveal in the alkali lines, and in Section 6 we construct photometry in customized spectral bandpasses probing the time-resolved photometric behavior in these lines. We compare our results with expectations for simple cloud-driven variability in Section 7, and consider but reject explanations based on lightning or aurorae in Section 8. We discuss more complex and realistic cloud-based explanations in Section 9, and offer our conclusions in Section 10.

We used the GMOS spectrograph (Hook et al., 2004) on the 8-meter Gemini South telescope in its multi-slit configuration to observe Luhman 16 simultaneously with 10 similar-brightness field stars. To avoid slit losses and obtain high-precision spectrophotometry, we used a very large, 5-arcsec width for our slits. To compensate for the significant sky background admitted through such wide slits, we employed the nod-and-shuffle mode of GMOS, a remarkable capability that enables very accurate sky subtraction even in the presence of strong fringing. In order to distinguish real spectral features from any column defects in the CCDs, we alternated pairs of images between two ‘spectral beams’ by changing the grating tilt to move the spectra 50 pixels along the direction of dispersion. We will refer to these two grating tilt settings as beam 1 and beam 2.

Due to the complexity of the nod-and-shuffle observations, the elapsed time for the acquisition of each GMOS science spectrum was about 480 seconds, but only 180 seconds consisted of actual integration on Luhman 16. This reduction in observing efficiency was more than compensated by the excellent sky subtraction enabled by the nod-and-shuffle mode. On February 24, we obtained twenty-six good quality science spectra spanning 7.1 hours of time, and on February 25 we obtained thirty-six good spectra spanning 8.2 hours. We used the OG515_G0330 order-blocking filter on February 24, which blocks light shortward of 5150Å. On February 25 we switched to the RG610_G0331, which blocks light shortward of 6100Å. This change was made out of concern that the spectra of comparison stars in the 10000–11000Å  region would be contaminated by their second order 5500–6000Å  spectra in the February 24 data, though such contamination would not directly affect the Luhman 16 spectra due to the object’s extremely faint flux at blue wavelengths. Ultimately, the blue sensitivity of the February 24 spectra turned out to be more interesting than the theoretically more precise measurements in the 10000–11000Å  regime from the following night, because Luhman 16 did not exhibit any unusual behavior in the longer-wavelength regime.

The seeing was much better (0.69 arcsec) on Feb 24 than on Feb 25 (1.46 arcsec), but through deconvolution we were able to obtain resolved spectrophotometry of Luhman 16 A and B on both nights. The good seeing on Feb 24 enabled us to measure the separation between Luhman 16 A and B (perpendicular to the dispersion direction) at 7.738 pixels, or 1.131 arcsec. This precise measurement was important for deconvolving the spectral traces of the binary on both nights, as discussed below. Concurrently with our GMOS observations, we obtained unresolved photometry in the Sloan $z^{\prime}$ filter using the TRAPPIST telescope (see, e.g., Jehin et al., 2011; Gillon et al., 2013). The TRAPPIST data provide an important context for the GMOS spectrophotometry.

Since our images were taken using the nod-and-shuffle mode of GMOS, the science frames contain two spectra for each slit, offset vertically in the image by means of electronic ‘shuffling’ on the CCD during the exposure. One contains the target object and the other a sky background spectrum that is an extremely accurate match to the background that underlies the target spectrum. Thus, sky subtraction is a simple matter of shifting the sky spectrum vertically by an integer number of pixels and subtracting it from the science spectrum. We find that the nod-and-shuffle mode is very effective, and the complex, fringed sky background vanishes from the science spectrum essentially completely (Figure 1).

We calibrate the wavelength solution using arc spectra. Because of our wide, 5-arcsec slits, lines in the arc spectra appear as wide blocks (Figure 2). We deconvolve these spectra using an edge-detection methodology, which begins with summing along image columns within the band of rows corresponding to each slit, and creating two 1-D vectors containing the positive and negative derivatives of the column-sums as a function of x pixel position. We discard all negative values of the derivatives, thus the positive derivative is nonzero only at the left-hand edges of ‘slit blocks’ and the negative derivative is nonzero only on the right-hand edges of the blocks. We find the cross-correlation peak of the two derivative vectors, which corresponds to the true slit width (5 arcseconds, or about 37 pixels). Finally, we align the two derivative vectors by shifting each by an amount equal to half the slit width, and average them to obtain a final deconvolved arc spectrum as shown in the bottom right panel of Figure 2.

We obtain a spectral wavelength solution by fitting the x-centroids of lines in the deconvolved spectrum vector as a quadratic function of the known wavelengths. We reject lines that show large residuals due to cosmic rays, bad columns, or proximity to the gaps in the detectors. As 9–15 lines survive unrejected for each slit, the quadratic fit to wavelength remains well-constrained. The flux from Luhman 16 is bright enough for variability analysis from 6000–10000Å, and the spectral pixel scale varies from 3.5 to 3.6Å per pixel over this range. In our Feb 24 data, when we used an order-blocking filter with a bluer cutoff wavelength of 5150Å, Luhman 16’s extremely faint ‘blue’ flux is detected all the way down to the cutoff. The actual spectral resolution element depends on the seeing, but is about 25 Å under typical conditions: thus, our spectral resolution varies from $\mathrm{R}=200$ to $\mathrm{R}=400$.

We extract both unresolved and resolved spectra from our images. The extraction methods are quite different, so in some cases the two types of spectra can offer quasi-independent confirmation that observed spectrophotometric behavior is intrinsic to Luhman 16 and not an artifact of our data processing.

The unresolved spectra contain light from both components of the binary, while in the resolved extraction we use a novel form of deconvolution to obtain distinct spectra for Luhman 16A and Luhman 16B. In both cases, our goal is to obtain precise, spectrally resolved relative photometry. We are not primarily interested in producing a physically calibrated spectral energy distribution, although we do attempt such calibration in an approximate sense in Section 2.7.

In processing the unresolved spectrophotometry, our reference for removing photometric effects imposed by Earth’s atmosphere is an average over the simultaneously-observed reference stars that occupy the other slits in our MOS mask (after confirming the reference stars have no significant astrophysical variability). In the resolved photometry, we target Luhman 16B, and our reference for removing photometric effects from Earth’s atmosphere is Luhman 16A — hence, we make the assumption (which can later be checked) that the latter object is not significantly varying.

The advantage of the unresolved spectrophotometry is its simplicity. Its disadvantages are that the variations of Luhman 16B are diluted by the light from its less-variable companion, and that obtaining accurate spectrophotometry is problematic at wavelengths longer than 10000Å due to the faintness of the reference stars’ spectra at these wavelengths and (in the February 24 data) their contamination by short wavelength flux from the stars’ second-order spectra.

The advantages of the resolved spectrophotometry are the opportunity to study the two objects independently; the exquisite removal of telluric effects by probing two objects of identical color observed in the very same slit; and the ability to probe all the way out to the red limit of the GMOS CCDs at 11000Å. Its disadvantages are the challenge of reliably deconvolving the blended spectra, and the possibility that the reference object Luhman 16A could itself be slightly variable — although at most wavelengths we ultimately avoid the latter problem by solving algebraically for the photometric variations of both objects (Equation 1).

We extract unresolved spectra by simply summing along columns in each slit. Figure 3 shows spectra from 34 different images taken on UT Feb 25. The narrow width of the bands formed by the overplotted spectra illustrates the relatively small change in raw total flux through the whole night: atmospheric conditions were consistently photometric.

We reject cosmic rays from our raw spectra by first creating a mean spectrum for each object by normalizing all its spectra and then combining them using a trim mean to create a clean master spectrum. We identify bad data in each individual 1-D spectrum using a 5$\sigma$ threshold, where $\sigma$ as a function of wavelength is determined simultaneously with the trim mean. We correct these bad spectral pixels by replacing them with a scaled value from the normalized mean spectrum for the object in question. The left-hand panel in Figure 4 illustrates the result.

Interference fringing within the thickness of the silicon detector becomes significant at a wavelength of 8200Å and continues out to the red limit of GMOS sensitivity. These fringes are not from sky emission, which is cleanly subtracted thanks to the nod-and-shuffle mode of our observations. Instead, the locally monochromatic light from Luhman 16 is interfering with itself in the thickness of the silicon detector, at wavelengths where silicon is no longer sufficiently opaque to absorb all the light on one passage through the detector. This optical interference affects the fraction of light that is actually absorbed by the CCD, and depends in complex ways not only on the local properties (e.g., thickness variations) of the CCD, but also on the spatial illumination pattern of the incoming light. Hence, though fringes with similar morphology appear at the same locations in our spectral flatfield images, the flatfields cannot be used to correct the science data because, due to their more even illumination, the fringes they show have different and typically much lower amplitudes relative to those in the science images. The dependence on local properties of the CCD also means that, at a given wavelength, the fringing pattern is quite different between our two spectral beams (Figure 4, right panel). Because of this, we perform cosmic ray reduction independently on the spectra from each beam. We find that if instead we attempt to clean all the spectra together, high-amplitude fringes are incorrectly rejected as if they were cosmic rays and the spectrophotometry is adversely affected. Note that fringing ceases to be a problem at wavelengths shorter than about 8200Å, and in particular has no effect in the vicinity of the neutral potassium doublet centered near 7682Å.

The very different spectra of Luhman 16 versus the reference stars makes it desirable to obtain relative spectrophotometry by taking the ratio of Luhman 16B’s flux to that of Luhman 16A, without using the reference stars at all. Previous observations have generally indicated that the variations of Luhman 16A, if any, are much smaller in amplitude than those of Luhman 16B: thus, Luhman 16A can serve as a nearby, spectrally matched photometric standard star for studying its highly variable companion, Luhman 16B.

In order to study the variations of Luhman 16B in this way, we must first deconvolve the traces, which are heavily blended especially in our data from Feb 25, when the seeing was less good than on Feb 24. The requirement of high photometric precision necessitates a very careful approach to the deconvolution. In particular, we find that the centroids of the spectral traces must be accurately determined as a function of wavelength (equivalently, x pixel position) across the images. A preliminary fitting of the reference star traces turns out to be essential to deconvolving the blended traces of Luhman 16A and B.

The cross sections of the spectral traces of the reference stars are well-fit by Moffat functions. To perform these fits, we first smooth the spectral traces of the reference stars using a trim mean in a boxcar of width 101 pixels (about 350Å) in the dispersion (x) coordinate. In each column, we perform a least squares fit over four parameters: y centroid, Moffat $r$, Moffat $\beta$, and spectral flux. As expected, the best-fit Moffat profile is sharpest at the longest wavelengths and becomes blurrier as we move toward the blue, where the seeing is less good.

Since these four-parameters fits to the reference stars produce good results, we could in principle perform analogous six-parameter fits to the blended traces of Luhman 16A and B: we would simply require y centroids and spectral flux values for two traces rather than just one. In such fits, however, the y centroids of the traces for the two objects are highly unstable, making accurate photometry impossible. Hence, we must forward-model the y centroids of the Luhman 16 traces.

The first step is to precisely measure the separations of the two traces on our images. To do this, we fit the Moffat $r$ and $\beta$ parameters derived for Star 08 (chosen because it is closest to Luhman 16 on the sky) as cubic functions of wavelength. We adopt these cubic fits as defining the Moffat parameters to be used to fit the blended trace of Luhman 16A and B, and fit just the two trace centroids and spectral fluxes. Using a sequence of sharp images from February 24, and constraining the fit to the spectral region where both Luhman 16A and B are bright, we find that the traces are separated by 7.738 pixels (1.131 arcsec), with an RMS variation of only 0.012 pixels (0.002 arcsec).

This accurate value for the separation removes one free parameter from the Luhman 16 fit, since trace centroids for Luhman 16A and Luhman 16B are no longer independent. However, the extracted A/B flux ratio is extremely sensitive to variations in the y coordinate of the trace centroid, and even with the separation fixed, our fits to the centroid are not sufficiently stable as a function of wavelength. Hence, we find it necessary to forward-model the trace centroid as well as fixing the separation. This requires taking into account the effect of differential chromatic refraction (DCR) in Earth’s atmosphere, as well as accounting for any optical distortion or misalignment.

To achieve this, we use the traces of all reference stars except the contaminated Star 03. We register them by subtracting the mean y coordinate over the range 7000–7500Å, and then we average the stellar traces to reduce the noise. The resulting average stellar trace changes from image to image, but in a very predictable way: it can be modeled extremely well by a small, time-invariant linear slant plus the time-dependent contribution from DCR. Hence, we can model the Luhman 16 trace centroid as a function of wavelength (which we have already mapped to x pixel position), even as the trace curves over time due to DCR. All that is needed is a precise value for a constant offset in the y coordinate. A slightly different offset must be used for each image, presumably due to guiding errors, gravitational flexure in the telescope, etc. We find that this constant offset cannot be forward-modeled with sufficient precision from the reference stars, so it is the sole positional parameter that we determine by fitting the Luhman 16 trace directly. We can obtain it with an accuracy of about 0.01 pixels by minimizing the RMS residuals of the double Moffat fit to the Luhman 16 trace. Thus, even though the trace y coordinate for Luhman 16 is a complicated function involving a linear trend added to a continuously changing nonlinear DCR component, we are able to forward-model every part of it except the tiny constant offset.

In our Moffat fits deconvolving the spectra of Luhman 16 A and B, we require $r$ and $\beta$ to be the same for both objects at any given wavelength¹¹1For the final fits, we cannot take Moffat parameters from the reference stars because the latter are too faint at long wavelengths.. This in itself makes spectrophotometry based on the ratio of the two objects fairly insensitive to the Moffat parameters, but we additionally fit each parameter as a global cubic function of wavelength on each given image — thus preventing unphysical local ‘glitches’ in the Moffat parameters that might otherwise produce spectrophotometric anomalies.

We note that except for the wavelength calibration, our resolved spectral extraction has almost nothing in common with our unresolved extraction. By design, the unresolved extraction is extremely simple: just a sum over a constant set of rows in the image. As described in this section, the resolved extraction is a sophisticated, customized deconvolution. The large methodological differences between the resolved and unresolved extractions mean that considerably increased confidence may be placed in any results that are confirmed in both.

Extracting physically calibrated spectra is not our primary focus, but since we have obtained observations of the A0V star HD 98072 as a telluric standard, we can obtain an approximate physical calibration for our data, with an estimated absolute uncertainty of about 10%. To do this, we compare our spectrum of HD 98072 with a resampled and blurred model spectrum of Vega (downloaded from http://kurucz.harvard.edu/stars.html; see, e.g., Kurucz (2011)). We first resample the Kurucz spectrum at 1Å, and then blur it with a Gaussian of $\sigma=$13Å, determined empirically to produce the best resolution match to our spectrum of HD 98072. We divide the blurred model Vega spectrum by the observed raw spectrum of HD 98072 to produce a calibration vector that naturally includes both telluric absorption and the spectral response curve of the instrument. As the Kurucz spectrum gives $F_{\lambda}$, and (we assume) the spectrum of Vega is a sufficiently good match to HD 98072, multiplying by our calibration vector converts an observed spectrum into $F_{\lambda}$ times a constant factor related to the ratio of Vega’s flux to that from HD 98072. To calculate this constant of proportionality, we use the Vega calibration from Hayes et al. (1985), which gives a monochromatic flux of $3.44\times 10^{-9}~{}\mathrm{erg}~{}\mathrm{cm}^{-2}~{}\mathrm{sec}^{-1}\mathrm{\AA}^{-1}$ for Vega at 5556Å, and the $V$ magnitude of $9.40\pm 0.02$ for HD 98072 from Hog et al. (2000). To physically calibrate Luhman 16, we apply an additional scaling to correct from the 5 second exposures of HD 98072 to the 183 second effective exposure time on Luhman 16.

For our input spectra, we use the average of six HD 98072 spectra taken throughout the night of February 24, and all 26 good science spectra of Luhman 16. Our uncertainty of 10% is estimated based on the extinction/reddening of HD 98072 (which we have not attempted to correct); possible differences between HD 98072 and Vega; and the imperfect (but fairly close) mean airmass match between our averaged spectra of HD 98072 and Luhman 16. A more precise calibration is beyond the scope of this paper focused on variability analysis.

In Figure 5 we show the resulting spectra in two forms: first, in a linear scale from 6400–9800Å; and second, in a log-scale showing the faint blue flux down to 5250Å. The spectra are dominated, first, by the very red spectral slope of these low-temperature objects, and second by the enormously strong lines of the neutral alkali metals: the sodium D doublet near 5893Å, and the potassium doublet at 7682Ų²2The sodium doublet consists of two lines at 5890Å and 5896Å, and the potassium doublet of lines at 7665Å and 7699Å, but herein we consistently quote the average wavelength for the two lines in each doublet, because each doublet appears in L and T dwarf spectra as a single, enormously broadened line.. The H₂O absorption band extending redward from 9270Å is lost in the interference fringing in Figure 5. This band is known to exist both in telluric absorption and in the intrinsic spectra of L and T dwarfs (Burgasser et al., 2008), and it seems to be detectable in our raw spectra (Figures 3 and 4), but not in the physically calibrated spectra of Figure 5, where telluric absorption has been corrected. The well-known flux reversal between L and T dwarfs is evident in Figure 5: the T dwarf Luhman 16B is brighter than Luhman 16A at wavelengths longer than 9000Å. This persists out to about 1.3 $\mu$m (far beyond the coverage of our data), after which the intrinsically more luminous L dwarf again shows a brighter spectrum (Burgasser et al., 2013).

Lightning is one possible source of bright emission with a spectrum very different from the photosphere. It has long been considered likely on brown dwarfs, due to their dynamic, high-energy clouds and the ubiquity of lightning in Solar System atmospheres (Helling et al., 2011a, b; Hodosán et al., 2016). Lighting could explain our Feb 24 data by producing faint blue continuum emission correlated with the dark clouds. Where the dominant red continuum of the brown dwarf photosphere is bright, this blue continuum would have negligible effect on the amplitude and phase. In regions where the red continuum is faint — i.e., the alkali lines — the blue continuum from the lightning would overwhelm the variations of the underlying photosphere and produce the observed amplitude reversal. Our Feb 25 results could be explained by supposing that only a subset of clouds is able to produce lightning, and then in the course of the general atmospheric disorganization that appeared to take place on Feb 25, this subset of clouds fell out of longitudinal alignment with the mean continuum obscuration, which perhaps by then was dominated by clouds unable to produce lightning.

This hypothesis implies extremely energetic lightning on Luhman 16B. It would remain a tiny fraction of the bolometric flux, but would still be orders of magnitude brighter than the lightning luminosities estimated by Hodosán et al. (2016) for the extrasolar planet HD 189733b, which is a cloudy object of comparable temperature and radius to Luhman 16B (albeit with much lower mass). However, the best estimates of extrasolar lightning still rely on scaling from Solar System bodies (e.g., Jupiter), and are highly uncertain when extrapolated to hotter and more massive objects such as Luhman 16B. It is not impossible that lighting could be more common in the violently evolving silicate clouds of brown dwarfs than on any object in the Solar System.

The lightning hypothesis makes a further prediction that can be tested using our data. It predicts that the flux of Luhman 16B should extend into the extreme blue where no measurable photospheric flux is expected, and that the amplitude of variation there should be extremely large, approaching 100%. Shortward of 6400Å, the trace from Luhman 16 is too faint for us to deconvolve the separate objects, but we have already seen (Figure 12) high-amplitude variations in unresolved photometry from 6000-6400Å that are correlated with those in the potassium line. Shortward of 6000Å the flux becomes too faint for time-series photometry, but by stacking all of our spectra we have managed to extract meaningful unresolved spectra at wavelengths down to 5250Å (Figure 5).

There is less flux in the extreme blue than the lightning hypothesis would seem to predict. Our unresolved synthetic photometry (Figure 12) shows a variability amplitude of 7% in the potassium channel (7530–7830Å). The average specific flux in this regime is $8.19\times 10^{-20}\mathrm{W}\mathrm{m}^{-2}\AA^{-1}$. If the amplitude reversal is caused by lightning, the specific flux from the lightning must be at least 7% of this value, or $5.73\times 10^{-21}\mathrm{W}\mathrm{m}^{-2}\AA^{-1}$, a flux level indicated by the dashed green line in Figure 5. If the lightning is producing a blue continuum, it should have similar or greater specific flux out to 5000Å  and shorter wavelengths. However, Figure 5 shows that the total flux drops below this value from 5400-6100Å. This means the total measured flux for Luhman 16 is less than we would predict for the flux from lightning alone on Luhman 16B alone. Making matters even worse, the very faint spectral trace in the 5250-5400Å  regime on the stacked images clearly aligns with Luhman 16A rather than Luhman 16B: there is no evidence that any of the flux at wavelengths shorter than 6000Å  comes from Luhman 16B.

This lack of blue flux does not explain away the huge variations we have detected in the 6000-6400Å band. These variations seem (based on their correlation with the potassium band) to come from Luhman 16B. Perhaps the lightning simply has an unexpectedly red spectrum, with little flux blueward of 6000Å. To evaluate this possibility, we offer a brief analysis of the energetics. In unresolved synthetic photometry, we measure an amplitude reversal in both the potassium channel (7530–7830Å) and the extreme blue (6000–6400Å). The integrated fluxes in these regimes are $2.46\times 10^{-17}\mathrm{W}\mathrm{m}^{-2}$ and $5.27\times 10^{-18}\mathrm{W}\mathrm{m}^{-2}$, respectively, and the amplitudes are 7% and 15%. Multiplying the integrated fluxes by the amplitudes and adding, we set a spectrum-independent lower limit of $2.52\times 10^{-18}\mathrm{W}\mathrm{m}^{-2}$ on the lightning flux from Luhman 16B, noting that in the case of realistic spectrum that would not have all its flux in the two regions where we can measure it, the total integrated flux would be larger. This spectrum-independent lower limit on the flux corresponds to a lightning luminosity of $1.20\times 10^{17}\mathrm{W}$ for Luhman 16B. For comparison, Hodosán et al. (2016) scale from measurements made on Jupiter and Saturn to estimate the total lightning energy released on the giant extrasolar planet HD 189733b in the course of a 1.89 hr transit to be $10^{17}$ J, of which a maximum of 10% could be released in optical light. This produces a luminosity of $1.5\times 10^{12}\mathrm{W}$, which is 80,000 times less than the absolute minimum lightning flux required to explain our Luhman 16 results in terms of lightning. Adopting any realistic spectrum for the lightning would raise this number by a factor of several.

Finally, the amplitude reversal in the potassium line itself is probably inconsistent with lightning. Lightning would be down in the clouds, with plenty of neutral alkali atoms still above it to wipe out the flux. The alternative would be something analogous to the ‘sprites’ observed on Earth, which extend very far above the clouds — but presumably such emission would involve only a tiny fraction of the total energy released in each lightning bolt, making the luminosity requirements even harder to meet.

The lack of blue flux from Luhman 16B; the amplitude reversal in the potassium line where lightning would probably be absorbed; and the implausible considerations in terms of total energy all suggest that lightning is not the explanation for our observations.

Aurora-like magnetic activity has been observed for some brown dwarfs (Hallinan et al., 2007; Berger et al., 2009), and would naturally appear above the alkali absorption in the atmosphere of Luhman 16B. It might explain our data by filling in spectral regions where the photosphere is faint. However, aurorae are not usually expected to produce continuum emission, and our results do not appear consistent with line emission. High-altitude fluorescence in the potassium line, for example, would create a highly localized signature at the core of the line, which would not match the broad effects we see. If aurora did create a hot continuum, it would be inconsistent with the lack of flux in the extreme blue that we have already invoked to rule out lightning. Finally, Osten et al. (2015) have made extremely sensitive radio observations of Luhman 16, ruling out auroral activity except at a very low level.

On February 24, we observed an amplitude reversal in the 7682Å potassium doublet of Luhman 16B relative to its red continuum, while on the same night in infrared lines of the same element, Kellogg et al. (2017) observed an amplitude increase that remained in phase with the red continuum. This apparently surprising discrepancy can likely be resolved by recognizing that because the optical lines of the alkali metals are far stronger than the near-infrared lines, the two wavelength regimes do not probe the same layers of the atmosphere. The optical lines we have observed probe far higher altitudes than the infrared lines measured by Kellogg et al. (2017), and the neutral alkali abundances could be different. In considering this, we depart from the simple conceptual model of Section 7, in which longitudinal variations in cloud cover were the only important cause of variability. Now, we picture a more complex and realistic atmosphere in which clouds, gas composition, and temperature profile are mutually interrelated; all three are longitudinally inhomogeneous; and the dominant cause of photometric variability could change as a function of wavelength. Specifically, we suggest that longitudinal inhomogeneity of obscuring clouds dominates the variability in the red continuum, but that inhomogeneous abundance of the neutral alkalis (influenced by clouds) becomes an important source of variability in the alkali lines.

Kellogg et al. (2017) note that their observations could be explained by a greater abundance of potassium above the dark clouds in Luhman 16B relative to the clear regions (though they do not ultimately favor this explanation). Our data could similarly be explained by a lower abundance of potassium over the dark clouds – but at a much higher altitude. Hence, both data sets could be explained if clouds in Luhman 16B affect the temperature and chemistry of the atmosphere above them in a way that increases the abundance of neutral alkali metals at low altitude (just above the clouds), while reducing it at high altitude. Proposing a detailed mechanism to produce such effects is beyond the scope of this work, but we note that in T dwarf atmospheres, the neutral alkali metals are in a chemical equilibrium with their chlorides (e.g., Burrows et al., 2002; Kellogg et al., 2017), and hence their abundances are sensitive to local conditions that could foster either the formation or dissociation of NaCl and KCl. We speculate that such conditions could involve not only cloud-driven changes in the vertical temperature profile, but also horizontal or vertical transport of gas with temporarily non-equilibrium chemical abundances — and additionally that solid particles in the clouds might act as catalysts for either the formation or dissociation of the alkali chlorides.

We can also interpret the Feb 25 results under this basic framework. On Feb 25, the decrease in overall amplitude indicates that the atmosphere of Luhman 16B was becoming less organized. As part of this disorganization, the alkali-poor masses of the high atmosphere might have fallen out of longitudinal alignment with the mean obscuration of the red continuum, causing the observed phase shift. We would then predict that at any time when the atmosphere of Luhman 16 is strongly organized into cloudy and less-cloudy regions, and hence the amplitude of variation is large, the well-organized clouds will again deplete the alkali gas at high altitudes above them (whatever the physical mechanism by which they do this), and an amplitude-reversal will again be observed in the Na D and potassium 7682Å  lines. It appears that such conditions were already observed by Biller et al. (2013), though the inverted amplitudes were diluted by their broadband $r^{\prime}$ and $i^{\prime}$ filters. Our prediction could be further tested by photometrically monitoring Luhman 16 using customized filters to deliberately probe the sodium and potassium lines, and comparing the results with simultaneous monitoring in the red continuum. Such monitoring could be an interesting tool for studying not only Luhman 16 but also other optically variable brown dwarfs.

Several previous investigations have found significant wavelength-dependent differences in simultaneously-measured brown dwarf light curves, and have argued convincingly that these are due to the differing pressure levels probed by the various wavelengths (e.g., Buenzli et al., 2012; Yang et al., 2016; Biller et al., 2018). Our results clearly fall into the same broad category. However, we do not take the further step, made in some previous works, of characterizing our results as indicating a pressure-dependent phase shift — even though it is tempting to interpret our Feb 25 data this way. We avoid doing this because we believe the concept of phase shifts smoothly varying with pressure is not astrophysically useful. The pressure-dependent phase shifts that have been suggested in some past works cannot be reconciled with any coherent physical picture of the cloud structures on the brown dwarf because the range of phases is too large. For example, a cloud band with a longitude dependent altitude (i.e., a cloud configured like a ramp) would produce a pressure-dependent phase shift, but the full range of the phase variation could not be more than a few degrees unless the ramp was very shallow. Even if it were an extremely shallow ramp that encircled the entire brown dwarf, the range of phases would have to be less than 180 degrees — unless there were a contrast reversal, which would produce a 180-degree phase shift (i.e., an amplitude inversion).

Such amplitude inversions at different wavelengths (i.e., pressure levels) have been observed by, e.g., Buenzli et al. (2012), Biller et al. (2018), and ourselves — but where a 180 degree phase shift exists between wavelengths that probe different pressure levels, there is no physical reason to expect any range of pressure levels to fill in a smooth variation of phase shift from zero to 180 degrees. If anything, the quantity that might vary smoothly with pressure in such a case would be the amplitude. Our data appear to show something like this for Luhman 16B on February 24, though on February 25 we see something more like a phase shift. However, the lack of a coherent physical interpretation for (wide-ranging) pressure dependent phase shifts leads us to regard the observations as simply showing large pressure-dependent changes in the light curves observed by ourselves and others, not simple phase-shifts. This broader view is better able to accommodate the observations reported by Yang et al. (2016), which in several cases show large pressure-dependent changes in light curve shape that are very different from pure phase shifts — as well as our own observations of pressure-dependent changes that manifest primarily as amplitude variations on February 24 and phase variations on February 25, with additional shape changes that cannot be simply captured with either amplitude or phase.

In Section 9.1 we put forth our best attempt to explain the observed variability of Luhman 16B on February 24 and 25 in terms of effects connected with a single dominant cloud system that was strongly organized on Feb 24 and became less organized on Feb 25. Following the principle of Occam’s Razor, we avoid proposing multiple major cloud systems where a single system might explain the data. Nevertheless it is worth noting that on Jupiter — the best Solar System analog to a brown dwarf — multiple cloud systems do exist (i.e., in different ranges of latitude), and very large pressure-dependent light curve changes are seen in part because different cloud systems dominate the rotational light curve at different wavelengths (Ge et al., 2019). For example, in the 8890Å  methane band, Ge et al. (2019) find Jupiter’s rotational light curve to be dominated by the Great Red Spot, while in the 5$\mu$m atmospheric window it is dominated by the North Equatorial Belt. The difference between the variations seen in the two wavelengths cannot be well characterized either by an amplitude inversion or by a phase shift: the 8890Å  light curve is double-peaked (though asymmetrical), while the 5$\mu$m light curve is single-peaked (though not sinusoidal). At yet a third wavelength (8.80$\mu$m), Ge et al. (2019) find a double-peaked light curve so symmetrical that if it were the only information we had about Jupiter, we would certainly claim its rotational period to be half the true value. Since there is no reason to think brown dwarf atmospheres are simpler than that of Jupiter, complex wavelength-dependent changes not fitting any simple model should perhaps be our default expectation.

Both from our data and many other results (e.g., Gillon et al., 2013; Apai et al., 2021), the light curve of Luhman 16B is known to evolve very rapidly with time, changing shape and amplitude (by a factor of two or more) on the timescale of one or two rotations. This suggests violent winds driving global changes in cloudcover, but it is interesting to note that the optical light curve of Neptune (Simon et al., 2016) exhibits similarly large fractional changes on timescales comparable to its rotation period. Neptune is far colder, less massive, and slower-rotating ($P\sim 17$ hr) relative to Luhman 16B, and its typical photometric amplitude ($\sim$1%) is about ten times smaller. Hubble Space Telescope observations suggest its variations are mostly due to relatively small, high-contrast clouds (Simon et al., 2016). A light curve produced by cloud features much smaller than a hemisphere could change dramatically without extremely fast winds or global atmospheric upheaval. On the other hand, unless the cloud features responsible for the variations of Luhman 16B have extremely high contrast, its much larger photometric amplitude implies the cloud patches must be larger relative to the object’s radius than those of Neptune — a scenario supported by the Crossfield et al. (2014) finding of large dark and bright regions with only moderate brightness contrast. The large size of these regions relative to the clouds of Neptune suggests that more violent weather is required to explain rapid light curve evolution in Luhman 16B. Such transformations could involve either vertical or horizontal transport of large masses of gas, perhaps faster than the alkali-chloride chemical equilibrium could adjust to the changing conditions. Transport of this type might explain the apparent transformation of the alkali effect from an amplitude inversion to a phase shift between Feb 24 and Feb 25.