Young stellar objects are thought to commonly undergo sudden accretion events that result in a rise in bolometric luminosity. These outbursts likely coincide with the onset of planet formation and could impact the formation of planets. The reason behind this dramatic enhancement of accretion is an active area of research and the mass of the system is a critical parameter. Using the Northern Extended Millimeter Array we survey five outbursting sources three FU Ori one EX Or and one peculiar source with the primary goal of determining the systems mass using an optically thin line of CO. We estimate the mass of a central region for each object that using both continuum emission and CSUP17SUPO J  21. The CSUP17SUPO emission likely includes both disk and inner envelope material thus acts as an upper limit on the disk mass ranging from 0.33 to 3.4 MSUBSUB for our sources. These derived masses suggest that the inner 1000 au contains enough mass along the line of sight for these sources to be gravitationally unstable.
['Protoplanetary disks', 'Young stellar objects', 'FU Orionis stars', 'Submillimeter astronomy', 'Interferometry', '1300', '1834', '553', '1647', '808', 'Astrophysics - Solar and Stellar Astrophysics', 'Astrophysics - Earth and Planetary Astrophysics']
High Mass Inner Regions Found in Five Outbursting Sources
{'High Mass Inner Regions Found in Five Outbursting Sources': "Jenny K. Calahan, 1, 2 Edwin A. Bergin, 1 Merel van 't Hoff, 1 Ke Zhang, 3 Nuria Calvet, 1 and Lee Hartmann 1  \n1 University of Michigan, 323 West Hall, 1085 South University Avenue, Ann Arbor, MI 48109, USA 2 Center for Astrophysics | Harvard & Smithsonian, 60 Garden St., Cambridge, MA 02138, USA 3 Department of Astronomy, University of Wisconsin-Madison, 475 N. Charter St., Madison., WI 53706", 'ABSTRACT': "Young stellar objects are thought to commonly undergo sudden accretion events that result in a rise in bolometric luminosity. These outbursts likely coincide with the onset of planet formation, and could impact the formation of planets. The reason behind this dramatic enhancement of accretion is an active area of research, and the mass of the system is a critical parameter. Using Northern Extended Millimeter Array, we survey five outbursting sources (three FU Ori, one EX Or, one 'peculiar' source) with the primary goal of determining the system's mass using an optically thin line of CO. We estimate the mass of a central region for each object that using both continuum emission and C 17 O J=2-1. The C 17 O emission likely includes both disk and inner envelope material, thus acts as an upper limit on the disk mass, ranging from 0.33-3.4 M ⊙ for our sources. These derived masses suggest that the inner ∼ 1000 au contains enough mass along the line of sight for these sources to be gravitationally unstable.  \nKeywords: protoplanetary disk, astrochemistry", '1. INTRODUCTION': "FU Orionis-type objects (FU Ori) and EX Ori-type objects (EXOr) have been theorized to be a common evolutionary phase between the Class I and Class II stages for low-mass pre-main sequence stars (Hartmann & Kenyon 1996; Quanz et al. 2007; Vorobyov & Basu 2015). An FU Ori object is classified as such after undergoing a sudden and extreme (several magnitudes) increase in brightness in optical and near-infrared wavelengths (Herbig 1977). This outburst has been attributed to short-lived and massive accretion events, predicted to accrete up to 10 -4 - 10 -5 M ⊙ /yr (Hartmann & Kenyon 1996; Audard et al. 2014; Fischer et al. 2023) onto the star. While ExOrs also undergo sudden increases in observed magnitude, they are less extreme and outbursts happen more frequently than FU Oris. There are a number of mechanisms that could cause this sudden and fast accretion, including thermal or gravitational instabilities (Bell & Lin 1994; Armitage et al. 2001; Boley et al. 2006; Zhu et al. 2010), perturbations from a close-by massive companion (Bonnell & Bastien 1992), or sporadic infall from the circumstellar envelope (Vorobyov & Basu 2006; Vorobyov et al. 2013).  \nRobust mass constraints on the disk and inner envelope are essential in determining what mechanism(s) is/are responsible for FU Ori outbursts. Gravitational instabilities are associated with triggering a magnetorotational instability (MRI) which would cause the sudden burst of accretion, and this would require the disk to be relatively massive (M disk /M ∗ ≳ 0.1) (Hartmann & Kenyon 1996; Liu et al. 2016; Cieza et al. 2018) and thermal instabilities require regions of the disk that have high enough opacity to trap heat (Zhu et al. 2009). Regardless of the mechanism, FU Ori systems need to be massive enough to sustain the high accretion rate for many years.  \nCO is commonly used as a mass tracer for molecular clouds, protostellar envelopes, and protoplanetary disks, as it is a highly abundant molecule, is chemically stable, and easily detectable. The measured ratio of CO/H 2 = 10 -4 in the ISM (i.e. Frerking et al. 1982; Kramer et al. 1999; Parvathi et al. 2012), and is consistent with the ratio found in young Class 0/I disks (van 't Hoff et al. 2020; Zhang et al. 2020). FU Ori objects present a unique laboratory for the purpose of mass determination of young stellar systems due to their recent outbursts. The outburst heats up the surrounding disk and envelope, pushing the snowlines beyond their regular radial  \nextent (Cieza et al. 2016; van 't Hoff et al. 2018), exposing more of the disk's mass to be probed with CO. The radial location of the CO sublimation front (temperature ≈ 20 K) during an FU Ori outburst can almost quadruple or even increase tenfold (at the intermediate stage of luminosity rise, and the outburst maximum, respectively) (Molyarova et al. 2018).  \n12 CO emission is found to be optically thick towards typical molecular clouds and protoplanetary systems, thus can only provide an lower limit to the mass. Given this, we use lesser abundant isotopologues of 12 CO to probe denser regions. 13 CO and C 18 O have often been used to probe higher density regions. C 18 O emission (both the J = 1-0, and 2-1 transitions) is generally assumed to be optically thin, and has been used as the mass tracer for FU Ori objects and other young stellar objects. Recent work shows, however, that C 18 O emission may not be completely optically thin. Available observations of FU Ori objects find that C 18 O J =1-0 can have measured τ > 1 (Fehér et al. 2017). Towards those same disks, emission from the brighter and more abundantly populated C 18 O J =2-1 line is expected to have an even higher τ . Recent observations with Northern Extended Millimeter Array (NOEMA) towards three non-FU Ori Class I objects showed that C 18 O and C 17 O J =2-1 were both optically thick in the disk region, and 13 C 18 O was required to reach optically thin emission levels (Zhang et al. 2020).  \nTo determine the total mass of the inner regions of these outbursting source, an optically thin tracer (C 17 O or 13 C 18 O) is used, along with the measured temperature (derived from C 18 O) and an assumed abundance of the optically thin tracer to H 2 . We present NOEMA observations of five outbursting sources with the primary science targets of 13 CO, C 18 O, C 17 O, and 13 C 18 O. We will use the optically thin tracer to estimate the mass for each source. We focus on the central region, as defined by the peak of the continuum, covering areas with equivalent radii of hundreds of au.", '2.1. Source Selection': "For our sample, we isolate sources that have been published and observed with NOEMA or ALMA in C 18 O J =1-0 or 2-1, and are located within the observational range of NOEMA (Fehér et al. 2017; Principe et al. 2018; Ábrahám et al. 2018). There are nine total sources that fit this criteria; we chose the brightest five sources, each with C 18 O emission peaking directly on source. Within the northern hemisphere, these objects have the highest likelihood of a C 17 Oand 13 C 18 Odetection. Source properties can be found in Table 1. Our sources range from  \nclassic FU Ori objects (V1057 Cyg, V1735 Cyg, V1515 Cyg), to a peculiar source with unique spectral features but otherwise FU Ori characteristics (V1647 Ori), to an EXor/UXor type star with a higher period of luminosity variability (V2492 Cyg).  \nV1057 Cyg and V1515 Cyg, along with the namesake source FU Ori were the first three discovered and defined the class of FUOr-type objects (Herbig 1977). V1735 Cyg was discovered later, and exhibited the characteristic sharp rise in luminosity followed by a slow decay with short-period modulations (Szabó et al. 2022). V1647 Ori was at first categorized as an FU Ori type object due to its extreme rise in luminosity, however has exhibited multiple significant outburst events since the first one captured, thus it has been classified as either borderline FU Ori/EX Or, or 'peculiar' source (Audard et al. 2014; Connelley & Reipurth 2018; Fischer et al. 2023). V2492 Cyg was seen to go into an outburst stage in 2010, and since then its variability in brightness has been attributed to a combination of episodic changes in extinction and accretion (Hillenbrand et al. 2013; Giannini et al. 2018; Ibryamov & Semkov 2021).", '2.2. NOEMA Data Reduction': "These results are derived from NOEMA observations that took place July - October 2021 as a part of project S21AL. Each source was observed for ≈ 3.4 hrs with beam sizes ≈ 1. '' 58 × 1. '' 50 (V1735 Cyg, V1057 Cyg) 2. '' 45 x 1. '' 69 (V1647 Ori) 1. '' 04 × 0. '' 92 (V2492 Cyg, V1515 Cyg). We utilized Band 3, specifically in the range of 202-211 GHz (lower sideband) and 218-226 GHz (upper sideband). This setup covers our key science lines, simultaneously covering the 13 CO, C 18 O, C 17 O, and 13 C 18 O J=2-1 transitions. The initial round of data quality assurance and reduction was done with a staff member at IRAM (Institut de Radioastronomie Millimetrique). The raw data was self-calibrated and CLEAN-ed using the GILDAS package mapping 1 . The first step was to find the systematic velocity and correct the UV-tables to the accurate values, which aids in line detection. Next we isolate the continuum by zeroing out every line and other feature such as dead pixels or end-of-band wiggles. With the channels containing only continuum emission left over, we average over every ∼ 200 channels before self-calibration and imaging. That solution is then imposed on the original spectra that contains the emission lines. The self-calibrated data is then  \nTable 1. Source Properties  \nNote -Distances taken from parallaxes in Gaia Release 3 (Gaia Collaboration et al. 2020), however it is worth noting that outburst soruces are particularly difficult to determine distances to. Outburst date and luminosity taken from review Audard et al. (2014) and references therein. Stated ranges indicate the timespan where the onset of the outburst could have started.  \ncontinuum-subtracted and imaged. For the goal of determining the mass of each source, we utilized the highvelocity resolution (0.5 km/s or 62.5 kHz) spectra of 13 CO, C 18 O, C 17 O and 13 C 18 O. V1735 Cyg likely contains significant extended emission in the three brightest CO isotopologues, as its solution contains artifacts that proved to be less strong if the smallest baselines were excluded. However, to maintain consistency in the analysis of all sources, the full NOEMA configuration including all baseline was used for the final results.", '3. RESULTS': "In our NOEMA sample we have observed 13 CO, C 18 O, and C 17 O towards all five sources, and the emission is partially to fully resolved. Towards one of them, V1057 Cyg, we also have an unresolved detection of 13 C 18 O. In the moment zero maps, all sources have 13 CO and C 18 O peaked on source, coincident with the continuum. C 17 O peaks on source for V1057 Cyg and V1515 Cyg, however the other three have C 17 O surrounding or nearby the central star. It is likely that the continuum may be blocking some central C 17 O emission for these three sources. In all sources, there is extended gas emission as compared to the continuum (with V1515 Cyg being being the most compact, see Figure 1).  \nWe seek to determine the mass of the central region from both the continuum observations and an optically thin gas line. For the continuum mass derivation, we use the following equation:  \nM cont = S ν d 2 κ ν B ν ( T ) . (1)  \nHere, S ν is the flux density, d is the distance, B ν ( T ) is the Plank function for a given dust temperature T, and κ ν is the dust opacity coefficient which we assume to be 2.2 cm 2 g -1 ( κ ν = 10( ν/ 1000 GHz) ; Beckwith et al. 1990; Cieza et al. 2018). For the dust temperature, we uniformly use T=25 K for both the dust and gas mass  \ndeterminations, as that is approximately the brightness temperature for the continuum of each object. We use a gas-to-dust mass ratio of 100, and thus find inner-region masses of 0.39-0.85 M ⊙ (see Table 3).  \nOur ideal gas-mass tracer will be an optically thin transition of CO, and so we first calculate the optical depths of each of these isotopologues. The brightness temperature of a molecule is defined by:  \nT B = T ex (1 -e -τ ) (2)  \nwhere T ex is the excitation temperature and τ is the optical depth. Taking the ratio between the brightness temperatures of two isotopologues will provide a relationship between the optical depths of the two lines:  \nT B ( 13 CO ) T B ( C 18 O ) = T ex ( 13 CO )(1 -e -τ 13 ) T ex ( C 18 O )(1 -e -τ 18 ) . (3)  \nThe excitation temperatures of each isotopologue of the same rotational transition can be assumed to be equal. And if you can assume that the 13 CO is optically thick, you can simplify this relation to  \nT B ( 13 CO ) T B ( C 18 O ) ≈ 1 1 -e -τ 18 . (4)  \nOne can further state that τ 18 = τ 13 /X where X is the abundance ratio between the two isotopologues. Thus, the ratio of the optical depths is equal to the abundance ratio. In our case we use ratio values from Qi et al. (2011) which were found from models that reproduce multiple observations of CO isotopologues in a protoplanetary disk, and are consistent with the ISM (Wilson 1999): 12 C/ 13 C=67, 16 O/ 18 O=444, and 18 O/ 17 O=3.8. If R = T B ( 13 CO ) T B ( C 18 O ) then you can express the optical depth of each isotopologue line as follows:  \nτ 13 = -444 67 ln(1 -1 R ) (5)  \nFigure 1. Moment zero maps of 13 CO, C 18 O, C 17 O J=2-1 for each source observed in this survey.  \n<!-- image -->  \nTable 2. Observational Properties  \nτ 18 = -3 . 8 ln(1 -1 R ) (6)  \nτ 17 = τ 18 / 3 . 8 (7)  \nUsing the observed flux ( S ν ∝ T B ) of each isotopologue, we calculated the optical depths for each observed line. We started with moment zero maps derived from our cleaned spectral cubes using bettermoments (Teague & Foreman-Mackey 2018). All five sources show optically thin emission in the C 17 O transition, while C 18 O is marginally optically thick ( τ ≈ 1 in most sources), and all 13 CO is optically thick.  \nThe mass of each source is therefore calculated using C 17 O emission. In an environment that is optically thin and in LTE, the column density of the upper energy state of a certain molecular transition is defined as:  \nN thin u = 4 πS ν ∆ ν A ul Ωhc , (8)  \nwhere S ν ∆ ν is the integrated flux density. The frequency range, ∆ ν , is defined by the width of the C 17 O emission profile in frequency space, and we use the same frequency range for each source's mass estimate. A ul is the Einstein-A coefficient that represents the rate with which an energy level is depopulated through spontaneous emission, and this is a constant for the C 17 O J=21 transition and equal to 10 -6 . 2 s -1 (Müller et al. 2005). Ω is the solid angle of the emitting area.  \nAfter the column density of the J=2 level is calculated, an LTE assumption can easily relate it to the total C 17 O column density, N T :  \nN T = N u g u Q ( T rot )[ e -E u /kT ] -1 , (9)  \nwhere g u is the statistical weight of the J=2 level, Q ( T rot ) is the partition function, and E u is the energy  \nFigure 2. The 13 C 18 O J=2-1 detection towards V1057 Cyg, summed over the central few pixels with an equivalent area to the beam.  \n<!-- image -->  \nof the J=2 level. All of these constants are found in the JPL line catalogue (Pickett et al. 1998).  \nOnce the total C 17 O column density is found, we can use the isotope ratios to back out the total CO column and then H 2 column density. We use C 16 O/C 18 O=1687 and CO/H 2 = 10 -4 (Wilson 1999; Frerking et al. 1982; Bergin & Williams 2017). Once the H 2 column density is derived per pixel, we sum over the surface area of a set of given pixels to produce a total mass. Using the continuum data, we define a central area of the image that contains at least all of the 15 σ continuum emission (see Figure 3). Three sources (V1057 Cyg, V1515 Cyg, and V1647 Ori) only included 20 σ flux so that the central areas were comparable sizes between sources. The final calculated masses from C 17 O range between 0.253.4 M ⊙ . These regions likely contains both disk and inner envelope emission in the corresponding C 17 O image (See Figure 4), in particular, emission from the envelope along our line of sight. All moment zero maps were used by implementing the bettermoments packages (Teague & Foreman-Mackey 2018).  \nThe average C 17 Ocolumn densities, mass column densities, derived mass from C 17 O, continuum flux, the corresponding equvalent emitting radius, and the ratio between the final gas and dust masses are shown in Table 3. These derived masses from the C 17 O are likely upper estimate on the mass of the disk, while also being lower-end estimates on the total mass within this emitting area as not all of our sources have their C 17 O peak overlap with the continuum peak.  \nFor the source V1057 Cyg, we also have a detection of 13 C 18 O which can be used as another mass tracer. The 13 C 18 O flux is not located in the central pixel, but is detected after summing over the area of the beam (see Fig. 2). We find an inner-region mass of 2.6 M ⊙ using  \nTable 3. Derived Column Densities and Mass  \n13 C 18 O. This used the sum of the flux over the same area used for the mass determination from both the continuum and C 17 O, and a 12 C 16 O/ 13 C 18 O = 29,748. Compared to the mass determined from C 17 O, it is slightly lower however well within a factor of two, thus we interpret this to be consistent. Both optically thin tracers suggest V1057 Cyg to be the most massive source within our sample.", '4. DISCUSSION': 'Mass is a fundamental property when it comes to understanding the protoplanetary disk environment and FU Ori outburst events. The quickest and most often used technique to extract disk masses is using dust thermal continuum observations. Assuming optically thin dust emission and a gas-to-dust ratio of 100, a gas mass can be calculated from the dust observations. We use this technique to estimate a lower-end mass of the central region, due to the fact that this continuum emission is likely optically thick (i.e. Tobin et al. 2020; Kóspál et al. 2021). Within similar objects, the dust that is probed at ∼ 220 GHz has been shown to be optically thick, and corresponding C 17 O emission often shows a dip corresponding with the center of the continuum emission (see V1647 Ori and V1057 Cyg in Figure 4). Using optically thin observations of C 17 O, we also calculate a mass of the central region, and treat this as an upper-end estimate of the disk mass. In each of our sources, C 17 O emission is optically thin, however it comes with some caveats. The emission likely probes the mass of the in-falling envelope and disk while the continuum probes large dust grains confined to the disk. Thus, the C 17 O emission likely probes more mass over a larger area, especially along the line of sight. However, the dust emission is expected to be optically thick and blocking C 17 O emission from the backside of the disk/envelope. C 17 O will be an upper limit if a source has an envelope, however if there is no envelope the C 17 O mass could be closer to a lower limit given an optically thick dust disk. Three of our disks have total inner-region masses that have similar values when derived from the continuum  \n(with gas/dust =100) and from an optically thin gas tracer (C 17 O). V1735 Cyg, V1647 Ori, and V1515 Cyg have masses that agree within a factor of less than two, suggesting that there may be little envelope contribution. V2492 Cyg, our EXOr source, has a gas mass that is just over 2 times more massive than what would be assumed from continuum observations alone. The final gas mass along the line of sight with a radial extent of ∼ 1000 au around V1057 Cyg is 3.4 M ⊙ , over four times more massive than what would be assumed from continuum alone. V2492 Cyg and V1057 Cyg likely have a significant contribution from the surrounding envelope, perhaps feeding the disk and playing a role in the mechanism that caused an episodic accretion event for these sources.  \nWithin our small sample, there appears to be no trend in mass versus outburst classification, peak luminosity, nor most recently observed brightness (most up-to-date photometic observations for each source: Kopatskaya et al. 2013; Peneva et al. 2009; Ninan et al. 2013; Szabó et al. 2022; Ibryamov & Semkov 2021). Our two most massive sources are an FU Ori (V1057 Cyg) and EX Ori object (V2492 Cyg) which exhibit markedly different peak luminosities and time since outburst. The fact that these two sources are the two most massive in our sample agrees with previous work that measured the mass of these two sources (Fehér et al. 2017) but does not agree with a trend seen when observing the continuum probed with ALMA that EX Ori sources tend to be less massive than FU Ori sources (Cieza et al. 2016). The most thorough samples of FU Ori masses come from continuum emission alone, with ≈ 18 determined masses (Cieza et al. 2018; Kóspál et al. 2021, and references therein). Compared to previously observed sources, our five sources reside alongside the most massive that have been seen with ALMA: V883 Ori, Haro 5a IRS, V900 Mon, and L1551 IRS 5 N. Our sample was biased to target the brightest outbursting sources in the Northern Hemisphere, thus it is not particularly surprising that we find higher masses than the average source in Cieza  \nFigure 3. Continuum maps corresponding to ∼ 221 GHz, extracted from the upper sideband of our NOEMA observations. The bulk of the continuum emission resides within the black contours, which is used to calculate a continuum mass. All contours correspond to at least where 15 σ emission is found, and the corresponding areas for each source are similar spatial scales. These areas are used to calculate a gas mass from C 17 O, see Figure 4.  \n<!-- image -->  \net al. (2018); Kóspál et al. (2021). Direct comparisons of the continuum mass are not straightforward between papers, as different emitting areas are used, temperatures, and distances. However, it is insightful to determine if our sample hold true with a trend seen thus far, that the majority of FUOr sources are gravitationally unstable.  \nThe Toomre Q parameter is often used to determine if a disk is gravitationally stable or not (Toomre 1964). A Qparameter greater than one suggests the disk is stable, while less than one is unstable. We use the following expression to calculate the Toomre Q parameter:  \nQ ≡ c s Ω πG Σ . (10)  \nThe angular speed of the gas disk ( Ω ) is defined by √ GM ⋆ /R 3 disk and the surface density ( Σ ) is calculated from the final total mass and emitting radius using the formula √ M disk /πR 2 . The sound speed of the gas (c s ) is calculated using √ k B T/µm H where k B is the boltzmann constant, T is taken to be the same temperature as used for the mass calculations (25 K), µ is the mean molecular weight which we use 2.3 (i.e. Armitage 2019). Using the values from Table 3, we find Toomre Q parameters from 0.06 - 0.79, with V1057 Cyg corresponding to the most gravitationally unstable. All sources in this  \nsample have a Toomre Q parameter below one, suggesting that the inner regions of all of these objects contain enough mass to trigger a sudden outburst event. Previous mass determinations of FU Ori sorces also find that the majority of these objects tend to be gravitationally unstable (Cieza et al. 2018; Kóspál et al. 2021).  \nDue to our beam-size, we are unable to resolve disk structure, and our continuum observations are unresolved on scales of a typical disk. Additionally, moment 1 maps of C 17 O showing the intensity weighted average velocity do not show typical signatures of disk rotation. Thus we are unable to make strong constraints on the protoplanetary disk size. Given our beam sizes of ∼ 1000 au, we are also unable to see clear Keplarian rotation from the disk in a position-velocity diagram. However, with such a relatively high mass within < 1000 au diameter region around each source, our results are consistent with previous work in determine FU Ori-type object central masses. What is clear about these systems is that they are full of complex structure. The moment 0 maps of 13 CO and 18 CO for each source show extended emission, and for three sources seem to suggest asymmetric features that could be part of infalling structures or streamers (V1057 Cyg, V2492 Cyg, and  \nFigure 4. Moment zero maps of C 17 O with the white contour indicating the central area from which a mass is determined, as defined by the continuum flux.  \n<!-- image -->  \nV1735 Cyg). There is very likely larger scale emission that is resolved out due to our beam size, and follow up work using single dish data will be able to focus on the large scale environment of each source. The triggering of these outbursting sources may not only rely on the disk mass and dynamics, but also the surrounding environment and infalling material.', '5. CONCLUSION': 'Using the NOEMA interferometer, we observed five outbursting sources, three FU Ori objects and two EX Ori/Peculiar sources. C 17 O was detected in each source and was found to be optically thin, thus a mass tracer. We found masses of the centrally peaked clumps from 0.31-3.4 M ⊙ using both dust and C 17 O emission to determine a range of masses for the central region of each source (along the line of sight with a radial extent of  \n< 1000 au around the central star). While these masses seem to represent upper limits on the total disk mass, it falls in line with a general trend seen in FU Ori objects suggesting massive disks that may be gravitationally unstable.', '6. ACKNOWLEDGMENTS': 'Matplotlib (Hunter 2007), Astropy (Astropy Collaboration et al. 2013, 2018), NumPy (Harris et al. 2020)  \nJ.K.C. acknowledges support from the National Science Foundation Graduate Research Fellowship under Grant No. DGE 1256260 and the National Aeronautics and Space Administration FINESST grant, under Grant no. 80NSSC19K1534.  \nE.A.B. acknowledge support from NSF Grant#1907653 and NASA grant XRP 80NSSC20K0259', 'REFERENCES': 'Ábrahám, P., Kóspál, Á., Kun, M., et al. 2018, ApJ, 853,  \n28, doi: 10.3847/1538-4357/aaa242 Armitage, P. J. 2019, in Saas-Fee Advanced Course, Vol. 45, Saas-Fee Advanced Course, ed. M. Audard, M. R. Meyer, & Y. Alibert, 1, doi: 10.1007/978-3-662-58687-7\\_1  \nArmitage, P. J., Livio, M., & Pringle, J. E. 2001, MNRAS, 324, 705, doi: 10.1046/j.1365-8711.2001.04356.x  \nAstropy Collaboration, Robitaille, T. P., Tollerud, E. J.,  \net al. 2013, A&A, 558, A33,  \ndoi: 10.1051/0004-6361/201322068  \nAstropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M.,  \net al. 2018, AJ, 156, 123, doi: 10.3847/1538-3881/aabc4f'}