Constraining the Source of the High Velocity Ejecta in the Type Ia SN 2019ein

We present multi-wavelength photometric and spectroscopic observations of SN 2019ein, a High Velocity Type Ia supernova discovered in the nearby galaxy NGC 5353 with a two day non-detection limit. SN 2019ein exhibited some of the highest measured expansion velocities of any Type Ia supernovae, with a Si II absorption minimum blueshifted by 24,000 km s$^{-1}$ at 14 days before peak brightness. More unusually, we observed the emission components of the P Cygni profiles to be blueshifted upwards of 10,000 km s$^{-1}$ with respect to the host galaxy redshift before B-band maximum light. This blueshift, among the highest in a sample of other Type Ia supernovae, is greatest at our earliest spectroscopic epoch and subsequently decreases toward maximum light. We discuss models of progenitor systems and explosion mechanisms that could explain these extreme absorption and emission velocities. Radio observations beginning 14 days before B-band maximum light yield non-detections at the position of SN 2019ein, which rules out symbiotic progenitor systems, most models of white dwarfs powering fast optically-thick accretion winds, and optically thin shells of mass $\lesssim 10^{-6}$ M$_\odot$ at radii $<100$ AU. Comparing our spectra to models and observations of other High Velocity Type Ia supernovae, we find that the observational properties of SN 2019ein match those of a delayed detonation explosion. We propose that the high velocity ejecta and blueshifted emission may be the result of abundance enhancements due to ejecta mixing in an asymmetric explosion, or optical depth effects in the photosphere of the ejecta at early times. These findings may provide evidence for common explosion mechanisms and ejecta geometries among High Velocity Type Ia supernovae.


INTRODUCTION
Type Ia supernovae (SNe Ia) are thermonuclear explosions involving at least one white dwarf (WD) progenitor star (Bloom et al. 2012). A unique characteristic of SNe Ia is that they show a relationship between their peak luminosity and the width of their light curve, known as the Phillips relation (Phillips 1993). This correlation al-lows the calibration of absolute brightness by light curve shape, which enables the determination of distances on cosmological scales. As standardizable candles, observations of SNe Ia have revealed the existence of dark energy (e.g. Riess et al. 1998;Schmidt et al. 1998;Perlmutter et al. 1999) and allow for a low-redshift measurement of the Hubble constant (e.g. Riess et al. 2019). A better understanding of their progenitor systems, explosion mechanisms, and observational characteristics is important to mitigate systematic uncertainties in order to use these objects for cosmological measurements.
Over the last several decades, sky surveys and deep imaging have led to the discovery of thousands of SNe Ia (e.g. Guy et al. 2010;Macaulay et al. 2019). Large samples have shown that significant diversity exists within the population of SNe Ia (e.g. Parrent, Friesen, & Parthasarathy 2014). Obtaining detailed observations of SNe Ia is important to understanding the sources of this diversity. While the majority of SNe Ia are "normal" and obey the Phillips relation, a sizeable minority are either overluminous, 91T-like SNe (Phillips et al. 1992;Filippenko et al. 1992a) or underluminous, 91bg-like SNe (Filippenko et al. 1992b) compared to this relation. These peculiar objects tend to show varied spectral evolution around peak brightness, suggesting that fundamental differences beyond luminosity exist in the population of SNe Ia.
Models of progenitor systems and explosion mechanisms have attempted to explain the observed photometric and spectroscopic heterogeneity. Most SN Ia progenitor systems are modeled by accretion onto a degenerate WD from a non-degenerate companion (the singledegenerate scenario, e.g. Whelan & Iben 1973) or by the accretion or merger of two degenerate WDs (the double-degenerate scenario, e.g. Iben & Tutukov 1984). In addition, a variety of theoretical explosion models have been able to reproduce observed characteristics of SNe Ia. One such model is a delayed detonation explosion, where a (subsonic) deflagration flame transitions to a (supersonic) detonation at some transition density (Nomoto et al. 2013;Iwamoto et al. 1999). Delayed detonation simulations are able to reproduce a wide variety of light curve widths, 56 Ni masses, ejecta compositions, and ejecta velocities in Chandrasekhar-mass progenitor WDs (Khokhlov 1991;Seitenzahl et al. 2013).
Several observational classification schemes have been proposed that may indirectly probe these different physical models. Branch et al. (2006) propose one such scheme, in which SNe Ia are classified by the strength of their Si II absorption features at maximum light. Additionally, Wang et al. (2009) sort SNe Ia by their Si II velocities, measured from the minimum of the absorption trough at B -band maximum light, into two classes: High Velocity (HV), with v Si II 12,000 km s −1 , and Normal, with v Si II 12,000 km s −1 . This diversity may arise from different distributions of Si in the outer layers of the ejecta, which in turn depend on the explosion mechanism. For instance, Mazzali et al. (2005) studied the Si II and Ca II absorption features in the Type Ia SN 1999ee and found that two separate components, separated by over 7,000 km s −1 , were visible in the spectra before B -band maximum light. The authors described these as high velocity features (HVFs) and photospheric velocity features (PVFs) and suggested they could be the result of additional mass at high velocities. Other studies have attributed HVFs, particularly of the Ca II NIR feature, to interactions between the SN shock wave and a shell of circumstellar material formed from the SN progenitor system (e.g. Gerardy et al. 2004;Mulligan, & Wheeler 2018;Mulligan & Wheeler 2017).
One distinguishing feature between explosion mechanisms is the symmetry of the ejecta. Kasen & Plewa (2007) modeled spectroscopically normal SNe Ia and found that, in asymmetric explosions, the color evolution and Si II 6355Å velocity evolution exhibit significant viewing-angle dependence. Additionally, Maund et al. (2010) found an empirical relation between the Si II 6355Å velocity gradient, as originally defined by Benetti et al. (2005), and the polarization across the same line, which traces the degree of the Si asymmetry in the ejecta (see e.g. Wang & Wheeler 2008 for a review). Therefore, a better understanding of the spectroscopic differences of SNe Ia is crucial to constraining their explosion mechanisms and ejecta geometries.
In this paper we present observations of SN 2019ein, an extreme HV SN Ia. Our early-time observations, beginning two weeks before maximum light, make SN 2019ein one of the best-studied HV SN Ia. The earliest spectral data at 14 days before B -band maximum light reveal some of the highest ejecta velocities ever measured. Perhaps more interestingly, the emission features in the P Cygni profile of SN 2019ein are blueshifted with respect to the redshift of its host galaxy. This systematic offset is greatest at very early times (several days after explosion) and gradually decreases as the SN evolves. Such a large emission shift sets SN 2019ein apart from other SNe Ia and hints at a puzzling explosion mechanism and ejecta geometry.
In Section 2, we detail our data acquisition, reduction, and analysis procedure. In Section 3, we present comprehensive early-time light curves from the nearultraviolet (NUV) to the near-infrared (NIR), along with model fits and fitted parameters. In Section 4, we present spectra, measure velocities of absorption features, and compare observations with a delayed detonation explosion model. In Section 5, we place limits on the nature of the progenitor system and the source of high velocity ejecta using early-time radio observations. In Section 6 we offer several possible explanations for the high velocity ejecta and blueshifted emission fea-  (Tonry et al. 2019). The last nondetection of the transient was by the Zwicky Transient Facility (Bellm et al. 2019) in r -band at a limit of 19.72 mag on 58602.27 (2019 April 29.27), implying that SN 2019ein was discovered within two days after explosion. SN 2019ein exploded in the outskirts of NGC 5353, a lenticular galaxy in a nearby galaxy group. An image of NGC 5353 along with SN 2019ein is shown in Figure  1. Using surface brightness fluctuation measurements, Jensen et al. (in prep) measure the distance to NGC 5353 to be 32.96 ± 1.68 Mpc.
We began daily photometric and spectroscopic followup observations of SN 2019ein starting on 2019 May 2 with the Global Supernova Project using Las Cumbres Observatory (LCO; Brown et al. 2013). Our first spectrum, obtained with the FLOYDS spectrograph on the 2m telescope at Haleakalā, allowed LCO to classify SN 2019ein as a young SN Ia (Burke et al. 2019). UB-Vgri -band data were obtained with the SBIG, Sinistro, and Spectral cameras on 0.4m, 1m, and 2m telescopes, respectively. Using lcogtsnpipe (Valenti et al. 2016), a PyRAF-based photometric reduction pipeline, PSF fitting was performed (Stetson 1987). UBV -band photometry were calibrated to Vega magnitudes using Landolt standard fields (Landolt 1992), while gri -band photometry were calibrated to AB magnitudes using the Sloan Digital Sky Survey (SDSS; SDSS Collaboration et al. 2017). Because the SN was offset from the host galaxy, image subtraction was not necessary.
We obtained four epochs of NIR photometry in JHK s filters with the 2MASS camera on the Minnesota-60 telescope on Mt. Lemmon, AZ, as part of the Arizona Transient Exploration and Characterization (AZTEC) program. The data were reduced and stacked with the IRAF 1 -xdimsum package. Aperture photometry was obtained with IRAF and calibrated to 20 local standards from the 2MASS catalog (Skrutskie et al. 2006).
We requested observations from the Ultraviolet and Optical Telescope (UVOT) on the Neil Gehrels Swift Observatory (Gehrels et al. 2004) after the detection of SN 2019ein. The first epoch of Swift data was obtained on MJD 58605.404 (2019 May 2.4), coincident with the first LCO photometric and spectroscopic epochs. Data were obtained in uvw2, uvm2, uvw1, u, b, and v filters and reduced using the data reduction pipeline for the Swift Optical/Ultraviolet Supernova Archive (SOUSA; Brown et al. 2014), including applying aperture corrections and zeropoints from Breeveld et al. (2011). Galaxy subtraction was not performed.
Spectroscopic observations are detailed in Table 1. 14 LCO spectra were obtained using the FLOYDS instruments on LCO 2m telescopes at Siding Springs and Haleakalā between -14 days to 60 days with respect to B -band maximum. Our spectra cover approximately the entire optical range from 3500Å to 10000Å at resolution R ≈ 300 to 600. Data were reduced using the floydsspec custom pipeline which performs flux and wavelength calibration, cosmic ray removal, and spectrum extraction 2 . In addition, we obtained several spectra in the optical and NIR using the B&C Spectrograph on the Bok 90 telescope, the Blue Channel Spectrograph on the MMT at the Fred Lawrence Whipple Observatory, and the SpeX spectrograph (Rayner et al. 2003) in PRISM mode with a 0.5 × 15" slit on the NASA Infrared Telescope Facility, which was obtained and reduced following the methods in Hsiao et al. (2019). These data are presented in Section 4.

Radio Observations
Radio observations of SN 2019ein were obtained with the Karl G. Jansky Very Large Array (VLA) on 2019  The observations were obtained in wide-band continuum mode, yielding 4 GHz of bandwidth sampled by 32 spectral windows, each 64 MHz wide sampled by 2 MHz-wide channels. We used 3C286 as our flux and bandpass calibrator, and J1419+3821 as our phase calibrator. Table  2 contains details of the observations. We obtained the datasets processed by the VLA CASA calibration pipeline, run on CASA version 5.4.1. 3 The pipeline consists of a collection of algorithms that automatically loads the raw data into a CASA measurement set (MS) format, flags corrupted data (e.g. due to antenna shadowing, channel edges, radio frequency interference or RFI), applies various corrections (e.g. antenna position, atmospheric opacity) and derives delay, flux-scale, bandpass and phase calibrations which are applied to the data. 3 https://science.nrao.edu/facilities/vla/data-processing/ pipeline For each epoch, the C-band data was split into 4-6 GHz and 6-8 GHz datasets, and each one was imaged using the CASA routine tclean. We use Briggs weighting of the data with a robust=0.7 to provide reasonable balance of angular resolution and source sensitivity. We used multi-term, multi-frequency synthesis as our deconvolution algorithm (set with deconvolver='mtmfs' in tclean), which performs deconvolution on a Taylorseries expansion of the wide-band spectral data in order to minimize frequency-dependent artifacts (Rau & Cornwell 2011). We set nterms=2 which uses the first two Taylor terms to create images of intensity and spectral index. Multiple bright radio sources appear offcenter in the 8.4 field of view, so we use "w-projection" (applied with gridder='wproject' in tclean) to account for non-coplanar effects when deconvolving these sources (Cornwell et al. 2008). The radio nucleus of the host galaxy is the brightest radio source in the field (peak flux ∼ 25 mJy) and forms artifacts near the site of the SN, so we performed a phase-only self-calibration with a solution interval of 2 mins to further clean and reduce the RMS noise in the image. The cleaned and self-calibrated 6-8 GHz image was then convolved to the resolution of the 4-6 GHz image using CONVL in AIPS. Both images were then combined using COMB in AIPS, No radio source was detected at the site of SN 2019ein in any of the cleaned deconvolved images down to 3σ limits of ∼ 18 µJy in the first image, and 25 and 23 µJy in the subsequent images. We discuss the constraints on progenitor models set by these limits in Section 5.

Light Curves of SN 2019ein
Swift uvw2, uvm2, uvw1, Johnson-Cousins UBV, SDSS gri, and 2MASS JHK s light curves are shown in Figure 2, along with SALT2 (Guy et al. 2007) fits to BVgri data between -14 days and 40 days with respect to B -band maximum light. Our high cadence observations make the rise of this light curve extremely well-sampled. Since SN 2019ein was discovered quite early, we are able to tightly constrain the rise time and explosion time. Given the SALT2 fits to the light curve, we find the B -band maximum occurred at MJD 58619.45 ± 0.03 (2019 May 16.5) which implies a rise time of ≈ 14 days since the beginning of observations and a maximum of 17 days since explosion (all phases hereafter are given in terms of B -band maximum light). Using an expanding fireball model Kawabata et al. (2019) estimated the explosion time of SN 2019ein to be MJD 58602.87 ± 0.55, giving a B -band rise time of ≈ 16.5 days, which is consistent with our estimates. This fast rise supports the suggestion that HV SNe tend to have shorter B -band rise times (Ganeshalingam et al. 2011).
SALT2 fitted parameters are given in Table 3 along with calculated values of absolute magnitude M, ∆m 15 (B), Milky Way E(B − V ), and distance modulus µ. Additionally, we correct for host galaxy red-dening by adopting the value presented in Kawabata et al. (2019), who estimate the host extinction as E(B − V ) host = 0.09 ± 0.02 mag. Overall, the fitted parameters show that SN 2019ein is a photometrically normal SN Ia, albeit with a slightly lower absolute magnitude at peak brightness (M Bmax = −18.81 ± 0.059) than expected. For a decline rate of ∆m 15 (B) = 1.40 ± 0.004, SNe Ia have on average M Bmax ≈ −19 (Hamuy et al. 1996). Therefore SN 2019ein falls slightly below the average, even with the modest reddening correction. We find good agreement between our estimated parameters and those derived in Kawabata et al. (2019), although our peak B -band absolute magnitude is fainter than their estimates, perhaps due to our use of a different distance modulus.

Color
The B-V color evolution of SN 2019ein is plotted in Figure 3, along with the color curve of the delayed detonation explosion model of Blondin, Dessart, & Hillier (2015), hereafter the B15 model. The model broadly matches the data at all phases, particularly around Bband maximum, although it tends to predict a bluer color at later phases. Similar trends can be seen in comparisons of other HV SNe with both the B15 model (Gutierrez et al. 2016) and NV SNe (Wang et al. 2009). At early times, the B-V color evolution matches the red group of Stritzinger et al. (2018), although among this sample SN 2019ein has a unique Branch classification (Branch et al. 2006), as described in Section 4.1. After correcting for host reddening, around B -band maximum the B-V color of SN 2019ein is 0.08 ± 0.04. This value falls in the overlap between the Normal and HV subsamples of Foley & Kasen (2011).
Additionally we measure the NUV-optical colors using our Swift photometry. The uvw1-v and u-v colors one day after B -band maximum are 1.58 ± 0.08 and 0.25 ± 0.06, respectively. These colors place SN 2019ein in the NUVR group of Milne et al. (2013), consistent with results that show HV and HVG SNe are all members of the NUVR group (Milne et al. 2013;Brown et al. 2018). Given our velocity measurements discussed in Section 4, the colors of SN 2019ein fit well with those of other HV SNe.

Bolometric Luminosity and 56 Ni Mass
Using our maximum light photometry, we estimate the bolometric maximum luminosity and the corresponding 56 Ni mass. We follow the methods outlined in Howell et al. (2009): first, the measured flux at B -band maximum is calculated by integrating the magnitudes in U, g, r, and i filters to ensure optimal wavelength coverage of the optical region. In order to calculate the flux across the rest of the spectrum we adopted the synthetic spectrum produced by B15 model. As we show in Section 4.3, this model produces close fits to the spectrum of SN 2019ein at maximum light. We scale the B15 spectrum flux to match the distance of SN 2019ein. Next we scale the flux of the synthetic spectrum within our filter wavelength ranges to match the observed flux. We then divide this "warped" flux by the ratio of the measured flux to the total flux in the synthetic spectrum, and define this quantity to be the bolometric flux.
Following this procedure we find a maximum bolometric luminosity of L ∼ 7.28 × 10 42 erg s −1 and a corresponding 56 Ni mass of ∼ 0.33 M , according to Arnett's rule (Arnett 1979) assuming a rise time of ≈ 16.5 days. This mass is on the low end for SNe Ia (Stritzinger et al. 2006). However, since SN 2019ein is a relatively fast decliner with ∆m 15 (B) = 1.40 and is slightly sublumi- nous (M B = -18.81), we conclude that this 56 Ni mass is a reasonable estimate.

SPECTROSCOPIC ANALYSIS
Figure 4 (left) shows the spectral evolution of SN 2019ein, from -14 to 60 rest-frame days with respect to B -band maximum light. In our earliest spectrum  ( Figure 4, top right) the most striking features include the broad absorption trough centered at approximately 7500Å, which is most likely the result of blended Ca II and O I absorption at high velocities (> 30,000 km s −1 ), as well as the broad Si II absorption centered at a wavelength less than 6000Å. The Si II absorption minimum corresponds to a velocity of approximately 24,000 km s −1 , which is one of the highest velocities ever measured in a SN Ia (Gutierrez et al. 2016). Additionally, the Ca II H&K absorption feature is not well defined in this spectrum. This could be due to blending with other absorption lines, or it may be that the line is blueshifted outside of the sensitivity of our spectrograph, although such a blueshift would correspond to a seemingly-unphysical velocity of ∼45,000 km s −1 . At this phase, the entire spectrum of SN 2019ein is noticeably blueshifted with respect to that of SN 2011fe. Before maximum light the blueshifts of the emission peaks remain prominent; the shifts are greatest in our first epoch, where both the Si II 6355Å and Ca II NIR components are displaced with velocities upwards of 10,000 km s −1 (Figure 4, bottom right). Also seen in the earliest spectrum is a small absorption "notch" at approximately 6250Å. This feature is most likely C II 6580Å at ≈ 20,000 km s −1 . Unburnt carbon in early-time spectra is not unusual (e.g. Parrent et al. 2011;Blondin et al. 2012;Folatelli et al. 2012;Maguire et al. 2014); however, few SNe Ia show Si II 6355Å absorption velocities higher than C II 6580Å absorption velocities at early times (Parrent et al. 2011;Folatelli et al. 2012;, with a notable exception being SN 2011fe (Parrent et al. 2012;Pereira et al. 2013).
The fact that SN 2019ein shows the opposite trend at this phase may shed light on the explosion mechanism and ejecta geometry. Parrent et al. (2011) note that v C II /v Si II < 1 if the C II feature comes from an asymmetric ejecta distribution viewed at an angle with respect to the observer's line of sight.
Close to B -band maximum we notice a possible HVF in the Ca II H&K line, as a weaker, lower velocity component becomes visible at roughly -4 days. However, HVFs usually develop at earlier phases, and this feature we observe is equally well fit by Si II absorption, making identification of HVFs at this phase difficult. Besides this exception we do not find evidence for HVFs in the spectra of SN 2019ein, and all the velocities we report here are measured from the center of the dominant absorption feature for each line. The lack of two distinct absorption components sets SN 2019ein apart from most other HV SNe Ia. We discuss possible reasons for this difference in Section 6.  Table 1. The phase with respect to B -band maximum is shown to the right of each spectrum. All fluxes are plotted on a linear scale. Top right: The spectrum of SN 2019ein at -14 days (red) compared with that of SN 2011fe at the same phase (black). Bottom right: The emission components of the Si II 6355Å (left) and Ca II NIR (right) P Cygni profiles before B -band maximum light. The rest wavelengths of these features are shown with dashed lines.
Using the MMT Observatory we obtained a medium resolution (R ≈ 3900) spectrum centered on the Si II 6355Å absorption feature at +18 days with respect to B -band maximum ( Figure 5). At this phase the feature takes on an unusual asymmetric appearance. In particular there appears to be multiple overlapping absorption troughs, each with a different line strength and Doppler shift. This may be caused by significant Si II mixing at this epoch, in which different distributions of Si II are moving at different velocities. This possibility is explored further in Section 6.
By approximately three weeks after maximum light the Si II feature begins to blend with iron group element (IGE) lines which dominate the spectrum. These IGE features are most easily seen in the NIR spectrum obtained 32 days after B -band maximum, shown in Figure  6. Line blanketing from IGEs are seen in between the two telluric regions and at wavelengths greater than 2.0 µm. At this later phase most C I and intermediate mass element (IME) lines, usually seen around maximum light (e.g. Hsiao et al. 2013Hsiao et al. , 2015, have disappeared from the NIR spectrum. Branch et al. (2006) showed that the ratio of the pseudo-equivalent width (pEW) of the Si II 6355Å absorption line to that of the Si II 5972Å line can be used as a spectroscopic classification of SNe Ia. Here we 5600 5800 6000 6200 6400 6600 6800 7000

Branch Classification
Rest classify SN 2019ein in the same way. We measure the pEWs with the following procedure: first, the spectrum is smoothed with a Savitzky-Golay filter to reduce the effects of noise. Next, the absorption feature of interest is defined and maxima blueward and redward of the absorption minimum along the continuum are found. We define the pseudo-continuum as simply the linear curve connecting the two maxima, so long as the curve does not intersect the spectral feature. Finally, the pEW is calculated using the formula (e.g. Garavini et al. 2007) where f (λ i ) is the measured flux, f c (λ i ) is the flux of the pseudo-continuum, and ∆λ i = λ i+1 − λ i is the size of the wavelength bin at each wavelength interval λ i . At maximum light, we find that the pEW of Si II 6355 is 125 ± 2.1Å and the pEW of Si II 5972 is 22.5 ± 2.8Å. The corresponding Branch diagram is plotted in Figure  7.

Absorption Velocities
SN 2019ein shows some of the highest expansion velocities of any SNe Ia in its early-time spectra. Velocities were calculated following the method outlined in Childress et al. (2014): first we select the absorption line of interest and define a pseudo-continuum by fitting a linear curve to the continuum maxima on both sides of the absorption trough. We normalize the flux with respect to this pseudo-continuum before fitting a Gaussian to the normalized absorption line. The minimum of the Gaussian is taken to be the Doppler-shifted observed wavelength, and the expansion velocity is calculated by comparing this measured absorption minimum to the known rest value of the line. Figure 8 compares the Si II 6355Å and Ca II H&K absorption velocity evolution of SN 2019ein to several other HV SNe Ia from Gutierrez et al. (2016). These objects all show similar spectral features (Figure 9), including strong Ca II NIR absorption at early times, broad Si II 6355Å absorption at maximum light, and high velocity Si II and Ca II before maximum. We do not report a Ca II H&K velocity -14 days with respect to B -band maximum light, as at this epoch no clear absorption minimum is identified within the wavelength range of our spectrograph.
The velocity evolution of all lines is rapid. The first epoch, corresponding to 14 days before maximum light and at most three days after explosion, shows the highest Si II velocity in this sample. By maximum light the ejecta velocity remains high, yet falls within the range of the other HV SNe. After maximum light, we measure a velocity gradient to the Si II velocity following the example of Blondin et al. (2012), who found that measuring the change in Si II velocity between maximum light and 10 ± 2 days after maximum gives the most consistent result. Using this method, we calculate a Si II velocity gradientv = 122 ± 25 km s −1 day −1 , placing SN 2019ein in the high velocity gradient (HVG) class ). Dessart et al. (2014) found that delayed detonation explosions best model BL, HVG SNe Ia. Figure 10 compares the spectra of SN 2019ein at various phases to delayed detonation model spectra produced by Blondin, Dessart, & Hillier (2015). The B15 model simulates the spherically symmetric delayed detonation of a Chandrasekhar mass WD, imposed with radial mixing to match abundance stratifications observed in SNe ejecta, particularly those of IMEs and IGEs. Synthetic spectra from explosion to nearly 100 days after maximum light are produced.

Comparison to a Delayed Detonation Explosion Model
Beginning 10 days before B -band maximum the synthetic spectra match the strengths and velocities of most of the absorption features, including the Si II 6355Å and Ca II NIR and H&K troughs. However, at our earliest epoch of -14 days the spectrum of SN 2019ein significantly deviates from the B15 model spectrum at the same phase. The model fails to reproduce the ex-tremely high absorption velocities, the broad mix of O I and Ca II NIR absorption, and the overall blueshift of the emission features with respect to rest. The authors found similar discrepancies when comparing the earliest model spectra to the early-time spectra of SN 2002bo, and suggested this may be due to underestimated outward mixing or a more complicated explosion than their one-dimensional, spherically symmetric model. Observational evidence for this enhancement of IMEs in the outer layers of the ejecta, possibly due to an extended burning front or significant mixing, was also found for SN 2002bo (Benetti et al. 2004) and the HV SN 2004dt (Altavilla et al. 2007).
In order to investigate the cause of this discrepancy we now explore possible sources of the HV ejecta in SN 2019ein, including interaction with a circumstellar shell of material from the progenitor system or mixing and optical depth effects in the ejecta.

PROGENITOR CONSTRAINTS FROM RADIO OBSERVATIONS
Radio emission is a sensitive probe of the progenitor environment (which we will refer to as circumstellar medium or CSM). The CSM is modified by mass-loss from the progenitor in the pre-SN stage, and interaction of the SN ejecta with this CSM accelerates electrons to relativistic energies and amplifies the ambient magnetic field, producing synchrotron radio emission (Chevalier 1982(Chevalier , 1984(Chevalier , 1998   provided constraints on the CSM environment and progenitor properties for both core-collapse (e.g. Ryder et al. 2004;Soderberg et al. 2006;Chevalier & Fransson 2006;Weiler et al. 2007;Salas et al. 2013) and SNe Ia (Panagia et al. 2006;Chomiuk et al. 2016). Radio emission is yet to be detected from a SN Ia, but non-detections have provided stringent constraints on progenitor scenarios (Chomiuk et al. 2016), particularly for nearby events like SN 2011fe (Horesh et al. 2012;Chomiuk et al. 2012) and SN 2014J (Pérez-Torres et al. 2014). We can similarly interpret possible progenitor scenarios of SN 2019ein by comparing our VLA observations with models of radio emission from circumstellar interaction.

Wind Model (∝ r −2 )
For single-degenerate progenitors a fraction of the mass, transferred via accretion from a non-degenerate companion, is expected to be lost in the form of a wind. Chevalier (1982)  ] Shells with f = 0.2 10 15 cm, 2 × 10 19 g/cc 10 15 cm, 7 × 10 19 g/cc 2 × 10 15 cm, 2 × 10 19 g/cc 2 × 10 15 cm, 7 × 10 19 g/cc 3 × 10 15 cm, 2 × 10 19 g/cc 3 × 10 15 cm, 7 × 10 19 g/cc 6 GHz Radio (VLA) such a wind, characterized by a constant mass-loss rate (Ṁ ) and wind velocity (v w ), which leads to a CSM whose density (ρ) varies with radius (r) as The synchrotron radio light curve from a shock propagating through such a CSM is described in Chevalier (1982) and Chevalier (1998). In this work, we follow the formalism of Chomiuk et al. (2016) (hereafter C16), which adopted the self-similar solutions of Chevalier (1982) for radio observations of SNe Ia. We assume a Chandrasekhar mass WD progenitor that exploded with 10 51 ergs of kinetic energy, consistent with our optical observations, and a steep outer ejecta profile of ρ ej ∼ v −10 ej interacting with the above CSM. Electrons are accelerated to a power-law spectra (∼ E −p , with p = 3). The average fraction of the shock energy shared by the cosmic ray electrons and the amplified magnetic field in the shock vicinity is parameterized as e and b respectively. As in Chomiuk et al. (2012) and Chomiuk et al. (2016), we set e = 0.1 and b = [0.1, 0.01], consistent with values expected in Type Ib/c SNe (Chevalier & Fransson 2006;Sironi & Spitkovsky 2011;Soderberg et al. 2012).
The light curve models for different values of the free parametersṀ and v w is shown in Figure 11(a). The rising part of the light curves corresponds to the regime where the ejecta are still optically thick to synchrotron self-absorption at 5 GHz. Once the ejecta are optically thin, the light curve declines. Higher ratios ofṀ /v w correspond to denser outflows, which leads to brighter light curves and a delayed transition to the optically thin stage, which explains why the peaks are shifted to later epochs. Figure 11(b) shows our constraints on theṀ /v w parameter space from the VLA upper limits in Section 2.1. We are able to rule out the parameter space forṀ /v w > 1.9 × 10 −10 M yr −1 (km/s) −1 for e = b = 0.1. For e = 0.1 and b = 0.01, we find thatṀ /v w > 9.5×10 −10 M yr −1 (km/s) −1 .
The above constraints onṀ /v w can provide some insight into possible single-degenerate progenitor models for SN 2019ein by comparing with typical values oḟ M and v w expected in these models as compiled in Chomiuk et al. (2012). Our observations are sensitive enough to rule out symbiotic progenitors, i.e. a WD that accretes from the wind of a giant companion, which is generally characterized byṀ > 10 −8 M yr −1 and v w ≈ 30 km s −1 (Seaquist & Taylor 1990;Patat et al. 2011;Chen et al. 2011). A symbiotic channel was also deemed unlikely for the nearest events SN 2011fe and SN 2014J, and was found to contribute no more than 16% of a sample of 85 SNe Ia with available radio observations studied by (Chomiuk et al. 2016). For e = 0.1 and b = 0.01, our increased upper limit oḟ M /v w > 9.5 × 10 −10 M yr −1 (km/s) −1 still excludes the majority of symbiotic progenitors observed in the Galaxy (Seaquist et al. 1993;Chomiuk et al. 2016).
WDs can also be in single-degenerate systems with a main-sequence or a slightly evolved companion undergoing mass-transfer via Roche lobe overflow. For mass accretion rates 3 × 10 −7 M yr −1 , steady nuclear burning occurs on the surface of the WD, and about ∼ 1% of the mass is lost from the outer Lagrangian point with velocities ∼ few 100 km s −1 (Shen & Bildsten 2007). The expectedṀ /v w in such a scenario falls within our VLA limits, and therefore such a progenitor channel cannot be ruled out for SN 2019ein from just our radio observations. With increasing accretion rate, however, the nuclear burning shell will drive fast, opticallythick winds with v w ≈ few ×1000 km s −1 (Hachisu et al. 1999), and some part of this parameter space is ruled out by our VLA upper limits. For accretions rates ≈ (1 − 3) × 10 −7 M yr −1 , the steady burning will be interrupted by recurrent nova flashes. Novae with short recurrence time will likely create a series of dense shells that the SN shock will interact with, and such shells typically have values ofṀ and v w shown in Figure 11(b). For longer recurrence times, the SN shock is more likely to interact with CSM created with a steady wind witḣ M ≈ 10 −9 − 10 −8 M yr −1 in between distant novae shells (Wood-Vasey, & Sokoloski 2006). Both these cases are allowed within our upper limits in the context of the r −2 model, but we will also analyze the presence of nova shells with a more appropriate shell-interaction model in Section 5.2.
We note here the importance of radio observations taken soon after explosion or discovery for SNe Ia. The first observation, which was triggered < 2 days of discovery and 4 days of explosion, provided almost a factor of five deeper constraint onṀ /v w than the observation a week later (Figure 11(a)). This is because lowerṀ /v w shifts the peak of the radio light curve to earlier times, as seen in Figure 11(a). The prompt observation resulted in more stringent constraints on Type Ia progenitor models involving symbiotic systems and optically thick winds.

Shell Interaction Model
Interaction and acceleration of material in solarcomposition CSM shells has been proposed as a way to explain HVFs in Type Ia spectra (Gerardy et al. 2004;Mulligan, & Wheeler 2018;Mulligan & Wheeler 2017). Such shells can be expected in WD progenitors that undergo nova outbursts. Shells can consist of recently ejected material or swept-up material from previous outbursts. Shock interaction with such a shell can produce detectable radio emission, and radio light curves for such a CSM created by discrete mass-loss events cannot be appropriately described by a continuous mass-loss model.
We therefore use the models described in Harris et al. (2016), hereafter H16, for radio emission from a CSM shell interacting with SN ejecta. H16 performed hydrodynamical simulations of a ρ ej ∼ v −10 ej ejecta profile interacting with a single, solar-metallicity, fully-ionized shell defined by an inner radius R, fractional width f = ∆R/R, and constant shell density ρ sh . Interaction creates a forward shock in the shell and reverse shock in the ejecta, but the dynamics do not reach selfsimilarity unlike the Chevalier (1982) case. The forward shock subsequently accelerates the CSM and sets it in free-expansion.
In the optically thin approximation, the H16 light model can be analytically expressed in terms of f, ρ sh and R (see Eqs. 5-13 in H16). Figure 12 shows example light curves from the H16 model. The light curves are characterized by a rapid brightening at the beginning of the interaction, reaching a peak luminosity when the forward shock reaches the outer edge of the CSM shell, and a steep decline once the shock breaks out. 4 For larger f , the light curves peak at later times as the shock takes . The left panel corresponds to a thin shell while the right corresponds to a thick shell. The gray shaded region is the parameter space where the H16 light curves are inconsistent with the constraints described in Section 5.2. The red shaded region is where optical depth due to synchrotron self-absorption is > 1 and the optically thin ejecta assumption of H16 is no longer valid.
longer to reach the CSM outer edge. For larger R, the light curves begin at a later time, and larger ρ sh produces brighter light curves. Similar to the analysis in Cendes et al. (2020), we explore the parameter space of R-ρ sh for a given f that produce light curves within our VLA upper limits at the observed epochs. We explore two cases of shells: a "thin" shell (f = 0.2) characteristic of shells expected in nova eruptions, and a "thick" shell (f = 1) to show the effects of increasing shell width. Similar to the wind model, the shell models assume a standard Chandrasekhar mass WD explosion with 10 51 ergs of kinetic energy, and e = b = 0.1. We also use an additional constraint: the peak of the light curve must occur before the first epoch (i.e. at 3.87 days after explosion). This is because any shell interaction leading to highvelocity absorption features must have occurred before the first spectral observation (i.e. after the shell has been accelerated by the forward shock, Gerardy et al. 2004). Figure 13 shows the result of applying the H16 shell models to our VLA observations. For our fiducial model parameters in both the thin and thick shell cases, the VLA limits only allow CSM shells 10 −6 M within radii < 100 AU. In comparison, CSM masses ∼ 10 −3 − 10 −2 M are generally required to explain HVFs observed in SNe Ia (Gerardy et al. 2004).
A caveat however is that the H16 model approximates the radio emission as optically thin, whereas at densi-ties > 10 5 cm −3 for the radii explored here, effects such as synchrotron self-absorption, free-free absorption and radiative cooling of the shock gas will become important. Light curve models in this optically-thick regime would require a formal solution of the radiative transfer equation, which will be explored in an upcoming paper (Harris et al., in prep), and will help provide more accurate constraints on the presence of dense and massive CSM shells.

MIXING AND OPTICAL DEPTH EFFECTS
The high emission velocities before maximum light make SN 2019ein unusual, even among HV SNe Ia. More specifically, although P Cygni emission blueshifts have been theoretically predicted and observed in Type II SNe (Dessart & Hillier 2005) and were discussed by Blondin et al. (2006) in a sample of low-redshift SNe Ia, the emission velocities seen in the spectra of SN 2019ein are the highest ever measured. Figure 14 shows the evolution of the emission velocities for four lines in the spectra of SN 2019ein, compared to the sample from Blondin et al. (2006). It is clear that at early times the emission peaks in SN 2019ein are substantially more blueshifted than in any of the objects in the comparison sample. This extreme behavior is most clearly seen in the plots of the Ca II H&K and S II emission velocities, where the emission velocities at -14 days with respect to B -band maximum are ≈ 15,000 km s −1 . At the same phase, the Si II emission component of the P Cygni profile is blueshifted by ≈ 10,000 km s −1 . Even around maximum light the velocities of these lines are among the highest ever measured. After maximum the emission peaks are either no longer resolvable or become distorted due to line overlap, possibly of multiple Doppler-shifted emission features (see Figure 5). Here we only present emission velocities up to B -band maximum, where we trust our measurements have not been biased.
In order to investigate whether specific ejecta compositions or abundance enhancements could cause both the high absorption and emission velocities at early times, we compare SYN++ (Thomas, Nugent, & Meza 2011) model spectra to our spectrum of SN 2019ein at -10 days. In particular, we focus on the Si II 6355Å feature and test whether multiple components of the ejecta, such as a HV component with velocity above the photospheric velocity (PV), can reproduce the measured Doppler shifts.
Our synthetic spectrum is shown in Figure 15 compared to our spectrum of SN 2019ein at -10 days. We find that an ejecta with only a HV Si II component offset from the PV by several thousand km s −1 provides the best fit to our data. This matches the lack of separate HVFs and PVFs, particularly in the Si II and Ca II lines, at early times in our observed spectra, as well as the analysis of our radio observations which places stringent constraints on the mass of a CSM shell, which has been proposed to produce HVFs. Altogether, this indicates that only a HV Si II component is present in the ejecta of SN 2019ein. We attempt to explain the existence of this HV component as being due to ejecta mixing from an asymmetric explosion or being caused by optical depth effects in the outer layers of the ejecta.

Evidence for Asymmetries
In Section 4.3 we discuss how objects similar to SN 2019ein exhibit significant mixing of their IMEs to higher velocities. This may be evidence of an aspherical ejecta distribution due to an asymmetric explosion, in which "clumps" of IMEs are mixed to higher velocities along the observer's line of sight, producing high velocity absorption and emission features. Similar "clumps" are produced in models of off-center delayed detonation explosions (Seitenzahl et al. 2013). The connection between mixing and asymmetries has observational support. Polarization measurements show that HVFs in SNe ejecta are more polarized than PVFs, indicating that HV ejecta have more asymmetric distributions (Maund et al. 2013;Bulla et al. 2016). Additionally, Nagao et al. (2019) observed high polarization alongside a blueshift of the Hα line during the photospheric phase of the Type II SN 2017gmr.
There is evidence for such an aspherical ejecta distribution in the spectra of SN 2019ein; in the early-time spectra of SN 2019ein we see Si II with a higher absorption velocity than C II at -14 days, whereas the opposite relation is true at this phase for the SNe studied by Parrent et al. (2011). In a spherically symmetric model of a Type Ia explosion the inner burning regions are surrounded by a shell of unburnt material comprised mostly of C and O. However, an asymmetric explosion with strong outward mixing could force IMEs produced in the nuclear burning to higher velocities. In addition, our medium-resolution spectrum obtained 18 days after B -band maximum light reveals multiple overlapping Si II absorption features, with each absorption minimum offset by several thousand km s −1 . This may be evidence of significant mixing of the Si II ejecta to lower and higher velocities, rather than a stratified shell-like structure proposed in studies of other SNe Ia (Cain et al. 2018).
Studying the nebular phase spectroscopy of SNe Ia, Maeda et al. (2010) found a relationship between the Doppler shift of the nebular phase emission features and the Si II velocity gradient at early times, whereby HVG SNe exhibit redshifted nebular phase emission lines and LVG SNe show blueshifted nebular emission lines. This correlation is suggested to detail information about the symmetry of the explosion, as the nebular phase emission features trace the deflagration ash in the core of the progenitor WD. In this model, HVG SNe are viewed from the direction opposite the initial deflagration.
Finally, Maund et al. (2010) and Cikota et al. (2019) found correlations between the Si II line polarization and velocity evolution around maximum, with more polarized SNe belonging to the HV and HVG classes. Maund et al. (2010) argue this relationship implies the existence of global asymmetries in the ejecta which, along with the correlation between velocity evolution and nebular phase velocity shifts (Maeda et al. 2010), connects early and late time velocity behavior to the three-dimensional geometry of the explosion. It is possible that the high absorption and emission velocities seen in the spectra of SN 2019ein are signatures of IMEs in the ejecta that were outwardly mixed in an off-center explosion.

Optical Depth Effects
Blueshifted emission features have been suggested to be caused by optical depth effects in Type II SNe (e.g. Dessart & Hillier 2005Anderson et al. 2014). These features arise from steep density profiles in the expanding ejecta (see Figure 16 in Dessart & Hillier 2005). Blondin et al. (2006) were the first to model these features in SNe Ia. Using the radiative transfer code CMFGEN (Hillier & Miller 1998), the authors found that differences in optical depths of Si II and S II lines resulted in overall blueshifted emission peaks when the flux was integrated over a range of impact parameters.
In this picture, contours of constant optical depth in the photosphere trace out a variable amount of the emitting ejecta, where the amount of emission from the ejecta above the photosphere depends on the optical depth of the line being considered. In the classical P Cygni profile, the flux from the ejecta moving perpendicular to the observer's line of sight makes up the majority of the emission, resulting in an emission peak centered on the line's rest wavelength. However, when lines with low optical depth are considered, there is little to no emission from the ejecta at large impact parameters due to the steep density gradient in the outer ejecta layers. Instead, the flux is dominated by ejecta moving toward the observer even if the ejecta is distributed more or less spherically. The result is an overall blueshift of the emission peak, proportional to the ejecta velocity. Because the authors modeled individual lines, the emission blueshifts cannot be the result of line overlap.
As seen in Figure 14, SN 2019ein has some of the largest emission peak velocities at early times compared to the sample from Blondin et al. (2006). This extreme behavior can be understood in the context of the above explanation; since blueshifted emission is dominated by flux from the ejecta moving toward the observer, one would expect that the emission velocity would be corre-lated with the absorption velocity. Blondin et al. (2006) found that this trend exists, with the ratio of v peak to v abs approaching 0.6 around -10 days relative to B -band maximum for the S II 5454Å line. Over time the photosphere quickly recedes from the low density material and the emitting region becomes more spherical, causing this ratio to approach 0 around B -band maximum light. However, at early times the photosphere is very far out in the ejecta, so the emitting ejecta we observe must be at high velocities, leading to a greater emission peak blueshift. Our SYN++ spectrum supports this picture; the blueshift of the Si II 6355Å emission peak is proportional to the Si II ejecta velocity. Since SN 2019ein has some of the highest absorption velocities measured at its earliest phases, the emission velocities of those lines are among the highest as well.

Discussion
We find that the high absorption and emission velocities at early times can be explained by a HV only ejecta component, possibly due to mixing in an asymmetric explosion or optical depth effects in the outer layers of the ejecta. It is possible that both effects are at play; the models of Maeda et al. (2010) predict that the outer regions of the SN ejecta on the side opposite from an off-center ignition are less dense and produce HVG SNe. It could be that in this lower density environment the ejecta is optically thinner, leading to a majority of the flux stemming from material moving along the observer's line of sight and producing blueshifted emission features.
Another nearby SN Ia with a well-studied density structure is SN 2012fr (Childress et al. 2013;Maund et al. 2013;Contreras et al. 2018). Cain et al. (2018) found that SN 2012fr showed signs of a shell-like density enhancement at low velocities, which could explain the unusual Si II velocity evolution as well as the presence of separate HV and PV features. However, SN 2012fr and similar SNe Ia tend to be slow decliners (∆m 15 1), HV yet LVG, and CN or SS on the Branch diagram (Contreras et al. 2018). These classifications are at odds with those we present for SN 2019ein. Therefore, we suggest that SN 2019ein most likely has a different density enhancement than the shell-like structure of SN 2012fr.
One potential bias in our measurements is line blanketing. Line blanketing can warp the shape of the emission peaks, potentially biasing measurements of the peak wavelength. As noted by Branch et al. (2007), synthetic spectra rarely if ever produce the significant emission peak Doppler shifts around maximum light observed by Blondin et al. (2006) and again here in SN 2019ein. However, synthetic spectra also rarely produce only HV ejecta components as we have done in our SYN++ model. Furthermore the observed emission shifts are not seen in just one line but globally, and seem to follow a similar evolution over time. This can be seen in the comparison between the B15 model spectra and the real data in the earliest epoch ( Figure 10). Therefore, we conclude that line blanketing is unlikely able to reproduce the peculiar emission blueshifts at all wavelengths and phases.
Future observations will be necessary to provide more conclusive results on the geometry of the ejecta. For example, measuring a Doppler shift of the nebular phase emission peaks could support the argument that SN 2019ein has signatures of an aspherical explosion. In addition, early-time polarimetry data would provide an additional measurement of the asymmetries in both specific spectral features, such as the Si II absorption line, and globally via the continuum polarization. However, as discussed by Dessart & Hillier (2011) and Kasen & Plewa (2007), low continuum polarization does not necessarily imply that the explosion was spherically symmetric; both authors find that even in models with significant asphericity the line and continuum polarization could be low due to density and ionization effects. In the case of some geometries presented in Dessart & Hillier (2011) both an emission blueshift and a low polarization signal are produced, regardless of the underlying symmetry of the ejecta.
It is reasonable to question why SN 2019ein is so unusual, even among a sample of other BL, HV SNe. One possible explanation is early time observations; it is possible that SN 2019ein was first observed mere hours after explosion, allowing us to see the extremely high absorption and emission velocities at an earlier phase than other HV SNe. Another explanation is that SN 2019ein was observed from a rare viewing angle, as would be the case if the explosion were strongly asymmetric. Either way, the analysis of early-time photometry and spectroscopy presented here demonstrates the importance of finding and observing SNe Ia quickly after explosion.

SUMMARY
We have presented photometric and spectroscopic observations of SN 2019ein, a SN Ia with some of the highest early-time ejecta velocities ever measured. We observe a Si II 6355Å absorption velocity of 24,000 km s −1 14 days before B -band maximum light. In addition, the early time emission components of the P Cygni profiles appear blueshifted with respect to the host galaxy redshift, with emission peaks of Si II, Ca II, and S II moving at velocities up to or above 10,000 km s −1 . This emission blueshift is also among the highest ever mea-sured, making SN 2019ein an outlier even among other High Velocity SNe.
Radio observations taken as early as <4 days after explosion provide insight into the progenitor system of SN 2019ein as well as the source of the high velocity ejecta. Our 3σ VLA upper limits of 18, 25, and 23 µJy at 3.87, 11.57, and 17.58 days after explosion are sensitive enough to rule out symbiotic progenitors for SN 2019ein. We also rule out part of the parameter space of a single-degenerate model involving accretion from a main-sequence or slightly evolved companion at accretion rates > 3 × 10 −7 M yr −1 , since the resulting fast optically thick winds would have likely created detectable circumstellar material. Such progenitor scenarios were also ruled out for the nearest and most wellstudied SNe Ia 2011fe and 2014J. Our upper limits cannot rule out models of a white dwarf accreting at lower rates ∼ (1 − 3) × 10 −7 M from a main-sequence or slightly evolved companion via winds that are sometimes interrupted by recurrent nova flashes. With our shell interaction model (Harris et al. 2016) we can rule out the presence of optically thin shells, which have been theoretically predicted to source high velocity ejecta, of masses > 10 −6 M at distances < 100 AU from the progenitor. However, denser or more massive shells in the optically thick regime cannot be ruled out by the current model, and will be revisited in the future with a more sophisticated shell model taking into account synchrotron self-absorption and radiative losses.
We find that SN 2019ein is well fit by a delayed detonation explosion mechanism except at early times, where our measured ejecta velocities are even higher than those predicted. By modeling the early spectra of SN 2019ein we find that both the high absorption and emission velocities may be due to a high velocity component of the ejecta that is detached from the photosphere. This detached component of the ejecta may be evidence of an aspherical distribution of intermediate mass elements, perhaps due to mixing in an asymmetric explosion (Seitenzahl et al. 2013). Additionally, optical depth effects in the very outer layers of the ejecta may lead to an overall blueshift in the spectrum, as the majority of the flux observed comes from material moving along the observer's line of sight. These results highlight the need for more detailed modeling of SN ejecta, especially at early times.
By studying a larger sample of High Velocity SNe Ia we can begin to probe the overlap between explosion models, asymmetries, and ejecta velocities. Results from such a sample would have implications on theories of Type Ia progenitor systems and explosion mechanisms. It is possible that a united picture will emerge, one in which the ejecta geometry and the viewing angle to a SN affect observables such as color, velocity, and light curve width. A similar intrinsic difference has already been noted in the colors and host environments of High Velocity and Normal SNe (Wang et al. 2013;Zheng, Kelly, & Filippenko 2018), and may be used to reduce uncertainties in Type Ia distances, improving the precision of cosmological measurements. C.P. thanks G. Mirek Brandt for helpful comments and discussions, Claudia Gutiérrez for generously sharing the velocity measurements and spectra used in Figures 8 and 9, Peter Iláš for his efforts in creating Figure  1, and Chelsea Harris for guidance with her shell model and its application to the radio observations. C.P. is supported by NSF grant 1911225 and NASA Swift GI 1518168. S.K.S. and L.C. are supported by NSF grants AST-1412549, AST-1412980, and AST-1907790. D.J.S. is supported by NSF grants AST-1821967, 1821987, 1813708, 1813466, and 1908972. Research by S.V. is supported by NSF grant AST-1813176. M. S. is a Visiting Astronomer at the Infrared Telescope Facility, which is operated by the University of Hawaii under contract NNH14CK55B with the National Aeronautics and Space Administration.
The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. This research made use of observations from the LCO network, as well as the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with NASA.