# ALMA Survey of Lupus Class III Stars: Early Planetesimal Belt Formation and
Rapid Disk Dispersal

J. B. Lovell1, M. C. Wyatt1, M. Ansdell2, M. Kama${}^{3},1$, G. M. Kennedy4,5,
C. F. Manara6, S. Marino1,7, L. Matrà8, G. Rosotti9, M. Tazzari1, L.
Testi6,10, J. P. Williams11
1Institute of Astronomy, University of Cambridge, Madingley Road, Cambridge,
CB3 0HA, UK
2National Aeronautics and Space Administration Headquarters, 300 E Street SW,
Washington DC 20546, USA
3Tartu Observatory, University of Tartu, 61602 Tőravere, Estonia
4Department of Physics, University of Warwick, Coventry, CV4 7AL, UK
5Centre for Exoplanets and Habitability, University of Warwick, Gibbet Hill
Road, Coventry CV4 7AL, UK
6European Southern Observatory, Karl-Schwarzschild-Strasse 2, 85748, Garching
bei München, Germany
7Max Planck Institute for Astronomy, Königstuhl 17, 69117 Heidelberg, Germany
8School of Physics, National University of Ireland Galway, University Road,
Galway, Ireland
9Leiden Observatory, Leiden University, P.O. Box 9513, NL-2300 RA Leiden, the
Netherlands
10INAF - Osservatorio Astrofisico di Arcetri, L.go E. Fermi 5, 50125, Firenze,
Italy
11Institute for Astronomy, University of Hawai’i at Mänoa, Honolulu, HI, USA
E-mail<EMAIL_ADDRESS>
(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

Class III stars are those in star forming regions without large non-
photospheric infrared emission, suggesting recent dispersal of their
protoplanetary disks. We observed 30 class III stars in the 1-3 Myr Lupus
region with ALMA at ${\sim}856\mu$m, resulting in 4 detections that we
attribute to circumstellar dust. Inferred dust masses are
$0.036-0.093M_{\oplus}$, ${\sim}1$ order of magnitude lower than any previous
measurements; one disk is resolved with radius ${\sim}80$ au. Two class II
sources in the field of view were also detected, and 11 other sources,
consistent with sub-mm galaxy number counts. Stacking non-detections yields a
marginal detection with mean dust mass ${\sim}0.0048M_{\oplus}$. We searched
for gas emission from the CO J=3-2 line, and present its detection to NO Lup
inferring a gas mass $(4.9\pm 1.1)~{}\times~{}10^{-5}~{}M_{\oplus}$ and gas-
to-dust ratio $1.0\pm 0.4$. Combining our survey with class II sources shows a
gap in the disk mass distribution from $0.09-2M_{\oplus}$ for $>0.7M_{\odot}$
Lupus stars, evidence of rapid dispersal of mm-sized dust from protoplanetary
disks. The class III disk mass distribution is consistent with a population
model of planetesimal belts that go on to replenish the debris disks seen
around main sequence stars. This suggests that planetesimal belt formation
does not require long-lived protoplanetary disks, i.e., planetesimals form
within ${\sim}2$ Myr. While all 4 class III disks are consistent with
collisional replenishment, for two the gas and/or mid-IR emission could
indicate primordial circumstellar material in the final stages of
protoplanetary disk dispersal. Two class III stars without sub-mm detections
exhibit hot emission that could arise from ongoing planet formation processes
inside ${\sim}1$ au.

###### keywords:

circumstellar matter - planetary systems - planets and satellites: dynamical
evolution and stability - techniques: interferometric - submillimetre:
planetary systems.

††pubyear: 2020††pagerange: ALMA Survey of Lupus Class III Stars: Early
Planetesimal Belt Formation and Rapid Disk Dispersal–C

##  1 Introduction

Significant progress over the past ${\sim}$decade has been made in our
understanding of the statistics of the planetary systems of nearby stars. We
now know ${\sim}50\%$ of Sun-like stars have tightly packed systems of low
mass planets at <1 au, ${\sim}5\%$ have Jupiter-like planets, and ${\sim}1\%$
have Hot Jupiters (Winn & Fabrycky, 2015). We also know ${\sim}20\%$ of stars
have detectable debris disks (i.e., belts of planetesimals at 20-100au) that
are analogous to, but much more massive than, the Solar System’s Kuiper belt
(Wyatt, 2008; Sibthorpe et al., 2018; Hughes et al., 2018). The formation of
these planets and planetesimals must occur in a star’s protoplanetary disk,
which survives up to ${\sim}10$ Myr, however these can disperse on timescales
$<1\rm$ Myr (Williams & Cieza, 2011). The finer details of planetesimal and
planet formation processes and their timescales are uncertain.

The growth of initially $\rm{sub-\mu\rm{m}}$ dust into cm-sized pebbles, as
well as their radial drift due to gas drag, is well understood (Birnstiel et
al., 2010). However, growth beyond cm-sized pebbles into planetesimals (that
then either grow into planets, or feed the debris disks seen later on) is not.
Growth beyond these pebbles is problematic due to radial drift and the
bouncing barrier (Windmark et al., 2012; Pinte & Laibe, 2014; Booth et al.,
2018). Additionally such large objects cannot be observed directly. The way to
form planetesimals may be to concentrate pebbles in regions sufficiently dense
for gravitational instabilities to form them directly (Johansen & Youdin,
2007). Such concentrations could arise in the earliest phases of evolution
($<0.1\rm{Myr}$) while disks are gravitationally unstable (Booth & Clarke,
2016; Nixon et al., 2018), or at the end of a protoplanetary disk’s life as it
disperses (e.g., between 3-10 Myr, see Carrera et al., 2017), or indeed any
time between these (i.e., $0.1-10\rm{Myr}$) if the disk is radially stratified
(Pinilla et al., 2016; Carrera et al., 2020). Theoretically therefore,
planetesimal formation could occur at any point within this broad range of
ages.

The collisional growth of planetesimals into planets is however well
understood, it just takes time. This is not a problem close to the star where
orbital timescales are short. Indeed the formation of close-in super-Earths
must be largely complete by ${\sim}3\rm{Myr}$ as ${\sim}50\%$ of stars have
such planets and ${\sim}50\%$ of stars have lost their protoplanetary disks by
this age (e.g., Haisch et al., 2001). However, at 10s of au, the long
timescales required for planet growth may be more prohibitive (e.g., Kenyon &
Bromley, 2010). Whilst pebble accretion may account for the growth of planets
at these distances (Bitsch et al., 2015), these require planetary seeds that
are $\sim$ Ceres sized, so this does not fully resolve this time scale issue.
Nevertheless, annular gaps are seen at 10s of au in $0.5-1$ Myr protoplanetary
disks (ALMA Partnership et al., 2015; Sheehan & Eisner, 2018; Andrews et al.,
2018), and these could be carved by planets (Dipierro et al., 2015; Rosotti et
al., 2016), suggesting planet formation processes in outer disks could be well
advanced by ${\sim}1\rm{Myr}$.

While planetesimals are not detectable, the dust that is created in their
collisions is. Such “debris disks” are seen around ${\sim}20\%$ of
($>10\rm{Myr}$) main sequence stars (Wyatt, 2008) showing planetesimal
formation is common at 10s of au. It is not possible however to say by looking
at main sequence stars when their planetesimals formed; i.e., these could have
had either short-lived or long-lived protoplanetary disks. Neither can
observations of protoplanetary disks answer this question, as second
generation debris disk dust would be fainter than the primordial
protoplanetary disk it is embedded within (e.g., see Fig.3 Wyatt, 2008).
Observations of young stars $<10$ Myr that have recently dispersed their
protoplanetary disks would shed light on disk dispersal processes, and provide
valuable constraints on where and when planetesimals form.

Young Stellar Objects (YSOs) in star forming regions can be classified based
on their mid-infrared emission (e.g., $24$ $\mu\rm{m}$). Such emission
originates from dust inside 10s of au of the star and is indicative of the
presence of a protoplanetary disk (Adams et al., 1987). The different classes
can be considered as an evolutionary sequence, with newly born stars thought
to evolve from class I for which a protoplanetary disk and envelope is
present, to class II for which the protoplanetary disk emission dominates, and
on to class III which are often considered to be diskless since no (or
relatively low levels of) mid-infrared emission is present above the stellar
photosphere. Such evolution can be readily seen in the spectral energy
distributions of YSOs in these stages, although we note in between the class
II and III stage also exist transition disks, with SEDs that can resemble a
further intermediary evolutionary stage between these (Espaillat et al.,
2014). Based on the amount of dust available and analysis of grain growth,
evidence exists that planet core formation may have begun prior to the class
II stage (Andrews et al., 2018; Manara et al., 2018b), and potentially as
early as the class 0/I phases (Greaves & Rice, 2011; Tychoniec et al., 2020).

A related classification scheme is based on evidence for ongoing accretion,
placing YSOs into classical T-Tauri stars (CTTS) and Weak-line T-Tauri stars
(WTTS), which are broadly similar to the class II and III except that some
WTTS have significant infrared emission and so massive disks (Manara et al.,
2013; Manara et al., 2017; Alcalá et al., 2019). Thus it is studies of class
III stars, which is those with a mid-IR spectral index $\alpha_{\rm{IR}}<-1.6$
(Lada, 1987; Williams & Cieza, 2011), which provide information on the nature
of systems with recently dispersed protoplanetary disks (although note that
care must be taken when compiling YSO lists based on this classification,
given that a range of wavelengths have been used to assess the mid-IR spectral
index, see Espaillat et al., 2014).

Despite this, class III stars have received little attention due to an
inherent (apparent) lack of infrared emission. However, constraints on their
circumstellar dust levels are usually weak, since star forming regions are
$>100\rm{pc}$. Moreover, the absence of $24$ $\mu\rm{m}$ emission (i.e., at
the longest wavelength for which large scale surveys exist) does not rule out
the presence of cold dust. Far-IR ($70-160\mu\rm{m}$) measurements are more
sensitive to cold dust, but are not available for all stars, at least not at
levels deep enough to detect dust emission except at protoplanetary disk
levels. Nowadays, sub-mm observations with ALMA provide the best sensitivity
to disks at 10s of au, and there have been several sub-mm surveys of class II
stars (e.g., Lupus, Taurus, Upper Sco, Ophiuchus; Ansdell et al., 2016;
Andrews et al., 2013; Barenfeld et al., 2016; Cieza et al., 2019). However,
only two sub-mm surveys included class III stars (Hardy et al., 2015; Cieza et
al., 2019), and these were only sensitive enough to detect the most massive
disks (at $M_{\rm{dust}}>0.2M_{\oplus}$). Thus, class III stars could host
debris disks comparable in mass to those around nearby stars (e.g., such as
Fomalhaut with $\sim 0.03M_{\oplus}$, see Holland et al., 2017), a possibility
also not ruled out by far-IR surveys of star forming regions.

Figure 1: Gaia data for the 69 candidate class III stars and 85 class II stars
(for which Gaia data is available): (left) histogram of parallaxes, (middle)
proper motions, (right) G band absolute magnitude vs colour. The 30 class III
sources deemed to be background objects are shown in black, while the final
sample of 30 class IIIs are shown with filled orange markers, the remaining 9
class III sources (not observed) with empty orange markers, and the observed
class II sources shown with blue markers. Isochrones were included on the
magnitude - colour plot to show the sample consistency with the expected age
of Lupus.

In this work we report on the first dedicated ALMA survey of class III stars
in the young ($1-3\rm{Myr}$) and nearby ($125-200\rm{pc}$) Lupus star forming
region (Comerón, 2008). Previous surveys have observed all class II stars in
this region at both $890$ $\mu\rm{m}$ and 1.3 mm down to dust masses of
${\sim}0.2M_{\oplus}$ resulting in 71/96 detections in Ansdell et al. (2016,
2018) with a further 3/6 detections of class II sources in Sanchis et al.
(2020). However, at this age ${\sim}25\%$ of stars have lost their
protoplanetary disks (e.g., Fig.1a of Dunham et al., 2015), so these class II
stars are the ${\sim}75\%$ with mid-IR indications of a protoplanetary disk.
The class III stars we present here are those without protoplanetary disks,
and therefore correspond to those within the population that for whatever
reason had the shortest lived protoplanetary disks. The ALMA survey discussed
here is deeper than any previous observations of class III YSOs by a factor of
5-10. At such levels we are able to detect emission down to debris disk levels
(i.e., $0.024M_{\oplus}$).

This work is the first of two papers reporting on these class III
observations. Here we discuss our survey of class III sources alongside a
wider sample of Lupus Young Stellar Objects, while Lovell, submitted. (2020)
presents a more detailed analysis of NO Lup, a class III source found in this
survey to have a CO gas counterpart. In $\S$2 we discuss how we selected the
sources for this survey, while further context for the photometric
classification of the full sample of Lupus YSOs considered in this work is
given in $\S$3. The observational set up of the ALMA measurements is described
in $\S$4. We analyse the continuum measurements in $\S$5, the CO measurements
in $\S$6, provide a detailed discussion of our results and their possible
interpretations in $\S$7, and summarise our key findings and conclusions in
$\S$8.

##  2 Source Selection

We compiled a sample of 82 candidate class III sources from Hughes et al.
(1994); Padgett et al. (2006); Merín et al. (2008); Mortier et al. (2011);
Alcalá et al. (2014). We used Gaia astrometry to distinguish candidate class
III young stellar objects from contaminant background objects from this
sample. Indeed, it has been shown that many candidate class III young stellar
objects in the Lupus star forming region are in reality background objects
(Manara et al., 2018a). Of the 82 class III candidates, 13 did not have Gaia
astrometric solutions (Gaia Collaboration et al., 2018) and were removed from
the sample (although we note may still be class III Lupus YSOs). Of the
remaining 69, 30 were deemed to be background sources due to Gaia parallax
smaller than 2 mas, corresponding to distances much greater than $200$ pc
whereas the distance to the majority of Lupus targets is ${\sim}160$ pc
(Manara et al., 2018a; Alcalá et al., 2019). In Fig. 1 the 30 background
sources are seen with empty black markers, where it is also shown that the 39
class III sources (orange, either filled/unfilled) have proper motions
consistent with Lupus YSOs, $\mu_{\alpha}$ and $\mu_{\delta}=(-16\pm 2,-21\pm
2)$ $\rm{mas}\,\rm{yr}^{-1}$ (Galli et al., 2013) and lie close to the 3 Myr
isochrone in magnitude-colour space (see right plot, except for one class III
faint G-band magnitude source). The final sample of 30 was then chosen from
these sources (excluding the faint G-band source, and preferentially choosing
the closer sources whilst maintaining a spread of distances), shown as filled
orange markers in Fig. 1.

The final sample of 30 Lupus class III’s only included sources from Lupus
Clouds I and III, and include 1 A-type, 14 K-type, and 15 M-types stars. Six
of the sample of 30 are in known binary systems, with the survey star in each
of these binaries being either an K-type or M-type star (Zurlo, submitted,
2020), see Table 2 and $\S$7.1.5 for more details. Stellar masses were
available from the literature for 16/30 of the sources (Manara et al., 2013;
Galli et al., 2015; Péricaud et al., 2017). We estimate the stellar masses of
the remaining 14 class III’s in $\S$5.4. We note that the average distance to
sources in the sample is $141$ ${\rm{pc}}$, slightly closer than the average
distance for the class II sources of ${\sim}160\rm{pc}$.

For the class II’s also being considered as part of a larger sample of Lupus
YSOs, we include the 95 ALMA observed sources discussed in Ansdell et al.
(2016, 2018), the two class II sources J16104536-3854547 and J16121120-3832197
(both of which were ruled out by Ansdell et al., 2018), and the additional 5
brown dwarf class II sources discussed in Sanchis et al. (2020) in Fig. 1 with
blue markers. Although 17 of these 102 class II sources were not found to have
Gaia astrometric solutions (and are not plotted) we still include these in our
sample for it to be as complete as possible. In the left plot of Fig 1, 6
sources have Gaia parallaxes smaller than 3 mas, and 5 of those 6 also have
low proper motions inconsistent with Lupus (Galli et al., 2013). Although we
cannot rule these out as having erroneous Gaia data, we excluded these 6
sources from our analysis; these are J15430227-3444059, J16104536-3854547, Sz
123B, J16120445-3809589 and J16121120-3832197 from Ansdell et al. (2016) and
Sz 102 from Ansdell et al. (2018); Sanchis et al. (2020), i.e., we also rule
out 2 of the same sources as Ansdell et al. (2018). We therefore consider a
total of 96 class II and 30 class III YSOs to form a representative sample of
126 Lupus class II and III stars. Further detailed discussion of the class II
sources is included in Ansdell et al. (2016, 2018); Sanchis et al. (2020),
which we refer to throughout.

##  3 Photometric Classification

###  3.1 Mid-IR Photometry

The sample of class III stars described in §2 used the classification given in
the literature (Merín et al., 2008; Evans et al., 2009; Cutri & et al., 2013;
Alcalá et al., 2014). Fig. 2 considers how the infrared properties of the
class III sub-sample compare with those of the class II sub-sample, by showing
their [K]-[12] and [K]-[24] colours. The grey dashed lines indicate the
regions for which $d\log{\nu F_{\nu}}/d\log{\lambda}>-1.6$ between K-band
(2MASS) and either 12$\mu$m (WISE) or 24$\mu$m (Spitzer) fluxes (i.e., the
photometrically defined boundary between class II and III, for which all
sources below both are defined as being in the lower left quadrant).

Figure 2: Lupus stars in our sample are shown in orange if class III, or else
if they are class II, in blue or purple (if these were sourced from Ansdell et
al. (2016) or Sanchis et al. (2020) respectively). Upper limits for those
without 24 $\mu$m detections are shown with downward pointing triangles, those
marked as blue "+" crosses were deemed as obscured (see Ansdell et al., 2016)
for the class II sources, and those marked as orange "x" crosses are deemed
confused class III sources in this work. The grey dashed lines show that
corresponding to $\alpha_{\rm{K-12}}$ and $\alpha_{\rm{K-24}}$ = -1.6 (with
the class III region shaded in green), and the green dot-dashed line
represents the cut to which all YSOs in the Ophiuchus survey (Cieza et al.,
2019) conformed (discussed further in $\S$7). We note that 24$\mu$m photometry
for J16071400-385238 was not available, therefore this data point is plotted
using its 22$\mu$m WISE 4 upper limit flux instead.

It can be seen that 29/30 of the class III stars are inside this lower-left
quadrant, alongside two class II sources. We note that both of these class II
sources are defined in Ansdell et al. (2016) as being obscured, and although
form part of the sub-sample of class II stars, are removed from subsequent
analysis for this reason. One class III source J16071400-385238 is outside of
this region (in the upper left quadrant), however further inspection of this
finds this source to be confused at 12$\mu$m (i.e., it does not appear to be a
point source in WISE W3). This source also lacks MIPS 24$\mu$m photometry
(i.e., it was too faint to be detected by Spitzer), and has its [24] magnitude
estimated from its WISE W4 (22$\mu$m) photometry, which was an upper limit,
and also happens to be the faintest in the sample (in the K band). One class
II source identified in Sanchis et al. (2020), AKC2006-18, falls below the
$\alpha_{\rm{K-24}}$ threshold, but not the $\alpha_{\rm{K-12}}$ threshold,
placing it in the lower-right quadrant (i.e., although it may have 24$\mu$m
photometry consistent with our class III definition, this is not the case for
its 12$\mu$m photometry). Nevertheless, this source remains distinct from the
population of class III YSOs. Therefore except for the three confused/obscured
sources, the complete class II and III sample are seen to be well separated,
and are shown to correspond as expected from $\alpha_{\rm{IR}}$ between the
2.2$\mu$m K band and 24$\mu$m mid-IR fluxes.

In addition, we find that the majority of our class III sample sit in a
cluster in the lower left corner of Fig. 2, which is where stars exhibiting
purely photospheric emission would be expected. However, falling outside or
inside this cluster is no guarantee of a star having discernable excess
emission or not. For that a more comprehensive analysis of the full SED is
required, which we discuss in $\S$5.4. For example, 4 of the sources in the
lower-left class III cluster are found in such an analysis to have infrared
excesses (MU Lup, J155526-333823, MV Lup, CD-31 12522), while 2 of the non-
confused sources external to this cluster are found to have no significant
excesses (J1609085-390343 and SONYC Lup 3-5).

Considering the class III sources that lie outside the lower-left cluster in
more detail, both NO Lup and Sz 107 are later confirmed to have mid-IR
excesses. Two class III sources are found outside the cluster but with
24$\mu$m upper limits that, while consistent with $\alpha_{\rm{K-24}}<-1.6$,
could also be compatible with an excess being present (these are also broadly
consistent with having no 12$\mu$m excess). There are two class III stars
outside the cluster for which we consider the 12$\mu$m emission to be
confused, and so their position on the figure should not be taken as evidence
for excess emission (Sz 109 and Sz 108). Inspecting Sz 109 with the IRSA
Finder Chart, it is found with significant levels of diffuse emission nearby
and does not well resemble a point source in WISE W3. This is one of the two
faintest YSOs in the sample (along with J16071400-385238) and is therefore
deemed as having confused photometry. Sz 108 is nearby (${\sim}4\arcsec$) to
the class II star, Sz 108B. This star could have contributed to Sz 108’s $12$
$\mu\rm{m}$ flux, and we note that the WISE photometric flag suggests this is
variable. Accounting for these confused sources, we therefore note that whilst
the existence of a photospheric excess is not ruled-out by being part of the
lower-left cluster, there appears to be weak positive correlation between the
sources outside/inside of this either having/not having a real mid-IR excess.
We consider in more detail the significance of these $12$ $\mu\rm{m}$ and $24$
$\mu\rm{m}$ excesses in $\S$5.4.

Six of the class II sources are found to have $\alpha_{\rm{IR}}<-1.6$ if 12
$\mu$m is used instead of 24 $\mu$m. This suggests that whilst these may have
large 24 $\mu$m fluxes, their warmer 12 $\mu$m emission is less excessive in
comparison to their photospheres. This might be expected if their disks have
an inner cavity, i.e., the morphology of the class of protoplanetary disk
known as a transition disk. From left-right (in [K]-[12]) these are
J16083070-3828268, Sz 111, Sz 84, Sz 115, J16000060-4221567, and
J16101857-3836125, and indeed the first three have previously been identified
as transition disks (Merín et al., 2008, 2010; Ansdell et al., 2016; van der
Marel et al., 2018). The left-hand side of this region (i.e., [K]-[12]<1.5,
[K]-[24]>3.3) may be a reasonable probe for transition disk objects. We note
that whilst Sz 107 sits below the [K]-[24] magnitude threshold that results in
classification as a class III, it is closer to the class II’s Sz 84 and Sz 111
than any of the other class III stars in our sample, which we will consider
further in $\S$7.

Figure 3: ALMA survey image showing all 30 observations in scale of
significance (with respect to each image rms). Significant emission is noted
with either an "A", "B" or "C", and numerically counted if more than one
source type in a single image is detected. The source numbers in the bottom-
left of each panel correspond to those in Table 2. In all images, North is up,
East is left, the beam and 200 au scale bar are shown in the upper and lower
right respectively.

###  3.2 H-alpha Equivalent Widths

We briefly consider here whether or not the sample of 30 class III YSOs shows
evidence for accretion based on their H$\alpha$ equivalent widths (EW). Using
the classification scheme of White & Basri (2003) (i.e., based on the width of
the H$\alpha$ line with respect to stellar type), we note here that 26/30 of
our class III stars showed no evidence of accretion (Krautter et al., 1997;
Cieza et al., 2007; Wahhaj et al., 2010; Manara et al., 2013; Galli et al.,
2015). In the case of 4/30 sources (J16071400-385238, J160758.9-392435,
J160908.5-390343 and J160917.1-392710) data was not available in the
literature to check this (i.e., these appear to be missing spectra).
Therefore, although in the case of 4 YSOs evidence for ongoing accretion
cannot be determined, none of the class III stars for which data was available
showed evidence of accretion. We will revisit this in $\S$7 when we consider
the evidence for other processes in these YSOs.

###  3.3 Summary of Photometric Classification

In summary, the classification in the literature translates well into the cuts
shown in colour-colour space. Accounting for confusion and flux upper limits,
we have shown that our sample of class III YSOs is photometrically distinct to
the class II sample in having low $12\mu\rm{m}$ and $24\mu\rm{m}$ spectral
indices, defined as $[\rm{K}]-[24]<3.3$ and $[\rm{K}]-[12]<2.2$. We will
discuss further the presence of mid-IR excesses in $\S$5.4. Additionally, we
showed that for all class III sources with data available, none showed
evidence for accretion.

##  4 Observations

ALMA was set up to observe the 30 class III YSOs in 21 separate scheduling
blocks, each comprising of two scans of each star as part of programme
2018.1.00437.S. In total each source was observed 42 times (i.e., for a total
of ${\sim}33.6$ minutes on each). All measurements were made in ALMA cycle 8
over the course of the first two months of 2019, using 43 antennae and
baselines ranging between 15-360m. The observing routine and sources used as
bandpass, flux and phase calibrators are summarised in Table 1. The single
phase calibrator J1610-3958 was used for all observations, and as both flux
and bandpass calibrators, either J1517-2422 or J1427-4206 were used (as
detailed in the table), and the precipitable water vapour (PWV) during these
measurements ranged between 0.4-1.7 mm.

ALMA was used in Band 7, having correlators set up with 3 spectral windows
centred on 343.240 GHz, 355.240 GHz, and 357.135 GHz each with a bandwidth of
2.0 GHz and channel widths of 15.625 MHz (for a total of 128 channels each)
set up for continuum observations. A fourth spectral window was centred on
345.135 GHz with a bandwidth of 1.875 GHz and channel widths of 488.28 kHz
(for a total of 3840 channels) set up for $\rm{CO}$ J=3-2 spectral line
observations. From these frequencies and bandwidths, we calculate a bandwidth
weighted wavelength of $\lambda=856\mu$m and quote this as our ALMA Band 7
wavelength herein.

The Band 7 visibility data sets were calibrated using the CASA software
version 5.1.1-5 with the standard pipeline provided by the ALMA Observatory.
No data flagging was applied. Data for each source was concatenated using the
CASA package concat, and then the concatenated measurement sets (for the
spectral windows conducting continuum measurements) were time averaged to
10-second widths and averaged into 4 channels (i.e., per spectral window) for
each of the field observations. Separate measurement sets were produced for
the analysis of J=3-2 CO emission, which were not channel or time-averaged,
(i.e., these retained the full 3840 channels). The CASA statwt task was run on
all measurement sets with which we conducted our continuum analysis, as
outlined in $\S$5.

Table 1: Observation data for the ALMA observing routine for the 42 scans on each source. Note that the Scan Time column represents the amount of time spent for each number of scans in a given schedule block. D.M.Y | $N_{\rm{scans}}$ | Scan Time | Flux and | Phase |
---|---|---|---|---|---
|  | [mins] | Bandpass |  |
23.01.19 | 2 | 48 | J1517-2422 | J1610-3958 |
09.03.19 | 2 | 48 | J1517-2422 | J1610-3958 |
11.03.19 | 4 | 48 | J1517-2422 | J1610-3958 |
12.03.19 | 2 | 48 | J1427-4206 | J1610-3958 |
12.03.19 | 2 | 48 | J1517-2422 | J1610-3958 |
13.03.19 | 4 | 48 | J1517-2422 | J1610-3958 |
14.03.19 | 2 | 48 | J1517-2422 | J1610-3958 |
15.03.19 | 2 | 48 | J1517-2422 | J1610-3958 |
16.03.19 | 2 | 48 | J1517-2422 | J1610-3958 |
17.03.19 | 2 | 48 | J1427-4206 | J1610-3958 |
17.03.19 | 4 | 48 | J1517-2422 | J1610-3958 |
18.03.19 | 6 | 48 | J1517-2422 | J1610-3958 |
19.03.19 | 2 | 48 | J1427-4206 | J1610-3958 |
19.03.19 | 6 | 48 | J1517-2422 | J1610-3958 |

##  5 Survey Continuum Analysis

This section presents analysis of the sub-mm continuum observations of the 30
class III stars. Detection types are described in §5.1. In §5.2 constraints
are presented for the the sizes and masses of the detections, and upper limits
are discussed in §5.3 for those where no such detection was made. SEDs of the
class III sample are considered in §5.4 to further constrain their dust
distribution, and we consider the potential for extra-galactic confusion in
§5.5.

###  5.1 Continuum Detections

Non-interactive t-clean was used to image each measurement set, over 50,000
iterations with a threshold of ${\sim}122$ $\mu$Jy (${\sim}$3 times the rms,
see later) and auto-multithresh masking applied, shown in Fig. 3. The field of
view of our images are set to show the entire primary beam (${\sim}12\arcsec$
radius). These images are scaled with respect to the image rms, which has a
mean level across all images of $36.7\pm 1.3$ $\mu\rm{Jy}\,\rm{beam}^{-1}$.
The average synthesised beam has a position angle $84.7^{\circ}$ and size
($0.75\times 0.62$)$\arcsec$ when imaged with natural weighting. The range of
major and minor axis extents, and position angles were $<0.005\arcsec$,
$<0.017\arcsec$ and $<0.9^{\circ}$ respectively, demonstrating strong
consistency in the synthesised beams for each source as expected given the
observing sequence.

The Target ID and image number in Fig. 3 correspond to the "No." column in
Table 2. We identified potential sources as those where the images exhibited
brightnesses above 3$\sigma$. To confirm these as sources we required a
$>5\sigma$ detection of emission within an aperture centred on the peak pixel
of the source. Three aperture sizes were used ${\sim}(1.1\times 0.9)\arcsec$,
${\sim}(1.5\times 1.2)\arcsec$, or ${\sim}(2.3\times 1.9)\arcsec$, to account
for the possibility that the sources are unresolved, partially resolved, or
well resolved, respectively. If the unresolved/partially resolved aperture
sizes measured consistent flux values, then we report the measurement based on
the smaller aperture (i.e., the one with the lower error). As expected, the
unresolved aperture flux measurement is consistent with that of the source’s
peak pixel value in units of Jy beam-1. In the case of the resolved source
J155526.2-333823 (the only source requiring the largest aperture), the peak
pixel is not well-centred on the source, and we centred this aperture on the
ALMA phase centre (which is expected to correspond to the centre of the disk).

In total, we detected 17 sources in 11 separate images, which we classified as
either type "A", "B" or "C" depending on their proximity to known class II and
III stars. Detections within $1.0\arcsec$ of the image centre were given an ID
of "A" and considered associated with the class III star. Detections
$>1.0\arcsec$ from the image centre consistent with the position of a known
class II star were given an ID of "B". Significant detections $>1.0\arcsec$
from the image centre, and inconsistent with the position of known class II
stars, were given an ID of "C". The separation between the image centre and
the centroid of the detected sources is given for all significant source
detections in the "Sep." column of Table 2. This analysis finds 4 A sources
(1A, 2A, 3A, and 4A), 2 B sources (5B and 6B), and 11 C sources (4C, 5C, 8C1,
8C2, 11C1, 11C2, 11C3, 12C1, 12C2, 19C, and 22C), and Table 2 is divided into
two parts based on this classification (i.e., depending on whether the images
contain a type A or B detection in the upper section, or contain either a type
C detection or no detection at all in the lower section). The following
apertures were used for the source fluxes given in Table 2: the large aperture
(2A), the partially-resolved aperture (5B, 5C, 8C2, 12C1, 12C2, 19C and 22C)
and the unresolved aperture (4C, 6B, 8C1, 11C1, 11C2 and 11C3). For the 26
class III stars where no type A detections were made, we determined the
$3\sigma$ upper limit on their fluxes from the rms within apertures with an
area equal to the size of the beam (i.e., $3\times 36.7\mu\rm{Jy}$, which
assumes disk radii $\leq\rm{FWHM_{beam}}\times d/2{\sim}50$ au). Given that
the expected sub-mm stellar flux for class III and class II sources is
significantly lower than the detection limits, we interpret these 4 A
detections and 2 B detections as due to circumstellar dust.

Table 2: Class III sample including stellar information and results from ALMA continuum imaging. Notes a: source partially lost at edge of image; total flux uncertain, b: two sources overlap; flux may be confused, c: Merín et al. (2008), d: Zurlo, submitted (2020), e: Hardy et al. (2015) f: Padgett et al. (2006), g: Köhler et al. (2000), h: Torres et al. (2006), i: Gaia Collaboration et al. (2018), j: Manara et al. (2017), k: Krautter et al. (1997), l: Manara et al. (2013), m: Galli et al. (2015), n: this work. No. | Target | Other Name | RA | Dec | Lupus | SpT | $M_{\star}$ | d | Sep. | ID | $F_{\rm{856\mu m}}$ | $M_{\rm{dust}}$ | Notes
---|---|---|---|---|---|---|---|---|---|---|---|---|---
|  |  | [J2000] | [J2000] | Cloud |  |  [$M_{\odot}$] | [pc] |  [$\arcsec$] |  | [mJy] |  [$M_{\oplus}$] |
1 | J154041.2-375619 | MU Lup |  15:40:41.170i |  -37:56:18.540i | I |  K6.0j | $0.89\pm 0.15^{\rm{m}}$ | $135.3\pm 0.8$ | $0.12$ | A | $0.26\pm 0.05$ | $0.045\pm 0.008$ | -
2 | J155526.2-333823 | - | 15:55:26.220 | -33:38:23.24 | I |  K5.0g | $0.91\pm 0.15^{\rm{m}}$ | $119.2\pm 0.7$ | $0.0$ | A | $0.68\pm 0.07$ | $0.093\pm 0.010$ | -
3 | J160311.8-323920 | NO Lup | 16:03:11.812 | -32:39:20.31 | I |  K7.0e,f |  0.7e | $133.7\pm 0.7$ | $0.0$ | A | $0.21\pm 0.04$ | $0.036\pm 0.007$ | -
4 | J160841.8-390137 | Sz 107 | 16:08:41.79 | -39:01:37.0 | III |  M5.5j,l | $0.16\pm 0.10^{\rm{l}}$ | $152.3\pm 2.5$ | $0.38$ | A | $0.22\pm 0.04$ | $0.050\pm 0.010$ | -
| - | - | - | - | - | - | - | - | $12.3$ | C | $1.45\pm 0.20^{\rm{a}}$ | - | -
5 | J160842.7-390618 | Sz 108 | 16:08:42.73 | -39:06:18.3 | III |  M1.0c |  0.43n | $151.3\pm 1.2$ | - | - | <0.11 | <0.024 |  Binaryd
|  | Sz 108B | 16:08:42.86 | -39:06:15.04 | III |  |  | $169\pm 3$ | $4.1$ | B | $23.5\pm 2.4$ | $6.5\pm 0.6$ | -
| - | - | - | - | - | - | - | - | $5.4$ | C | $1.20\pm 0.13$ | - | -
6 | J160957.1-385948 | Sz 119 | 16:09:57.07 | -38:59:47.9 | III |  M4.0c | $0.30\pm 0.05^{\rm{m}}$ | $110\pm 4$ | - | - | <0.11 | <0.013 | -
|  | Lup 818s | 16:09:56.28 | -38:59:51.90 | III |  |  | $157.0\pm 2.7$ | $9.8$ | B | $6.8\pm 0.7$ | $1.63\pm 0.21$ | -
7 | J153802.7-380723 | MT Lup | 15:38:02.724 | -38:07:22.17 | I |  K5.0k |  0.90n | $133.0\pm 0.7$ | - | - | <0.11 | <0.019 | -
8 | J154038.3-342137 | Sz 67 | 15:40:38.27 | -34:21:36.7 | I |  M4.0c |  0.26n | $110.1\pm 1.5$ | - | - | <0.11 | <0.013 | -
| - | - | - | - | - | - | - | - | $1.8$ | C1 | $1.00\pm 0.11$ | - | -
| - | - | - | - | - | - | - | - | $4.1$ | C2 | $0.41\pm 0.06$ | - | -
9 | J154306.3-392020 | - |  15:43:06.245i |  -39:20:19.482i | I |  K6.0j |  0.91n | $166\pm 7$ | - | - | <0.11 | <0.030 |  Binaryd
10 | J154641.2-361848 | CD-35 10498 | 15:46:41.199 | -36:18:47.44 | I |  K1.0h | $1.38\pm 0.46^{\rm{m}}$ | $149.5\pm 1.7$ | - | - | <0.11 | <0.024 |  Binaryd
11 | J154959.2-362958 | MV Lup | 15:49:59.219 | -36:29:56.59 | I |  K2.0f | $1.55\pm 0.55^{\rm{m}}$ | $139.3\pm 0.7$ | - | - | <0.11 | <0.021 | -
| - | - | - | - | - | - | - | - | $9.5$ | C1 | $2.75\pm 0.29$ | - | -
| - | - | - | - | - | - | - | - | $3.6$ | C2 | $2.77\pm 0.28^{\rm{b}}$ | - | -
| - | - | - | - | - | - | - | - | $2.9$ | C3 | $1.48\pm 0.15^{\rm{b}}$ | - | -
12 | J155219.5-381932 | MW Lup | 15:52:19.522 | -38:19:31.36 | I |  K7.0f | $0.80\pm 0.12^{\rm{m}}$ | $130.3\pm 0.8$ | - | - | <0.11 | <0.018 | -
| - | - | - | - | - | - | - | - | $10.2$ | C1 | $2.73\pm 0.29$ | - | -
| - | - | - | - | - | - | - | - | $11.2$ | C2 | $1.39\pm 0.19$ | - | -
13 | J155533.8-370941 | MX Lup | 15:55:33.795 | -37:09:41.08 | I |  K6.0f |  0.85n | $129.6\pm 0.8$ | - | - | <0.11 | <0.018 | -
14 | J160159.2-361256 | NN Lup | 16:01:59.272 | -36:12:53.89 | I |  K3.0f | $1.20\pm 0.29^{\rm{m}}$ | $136.3\pm 0.7$ | - | - | <0.11 | <0.020 | -
15 | J160352.5-393902 | CD-39 10292 | 16:03:52.502 | -39:39:01.28 | III |  K3.0h | $1.01\pm 0.25^{\rm{m}}$ | $130\pm 4$ | - | - | <0.11 | <0.019 |  Binaryd
16 | J160430.6-320729 | CD-31 12522 | 16:04:30.564 | -32:07:28.81 | III |  K2.0h |  1.17n | $130.6\pm 0.8$ | - | - | <0.11 | <0.018 | -
17 | J160713.7-383924 | - | 16:07:13.698 | -38:39:23.88 | III |  K7.0f | $0.75\pm 0.11^{\rm{m}}$ | $158.9\pm 2.3$ | - | - | <0.11 | <0.027 | -
18 | J160714.0-385238 | Lup 605 | 16:07:14.00 | -38:52:37.9 | III |  M6.5c |  0.03n | $124\pm 4$ | - | - | <0.11 | <0.017 | -
19 | J160749.6-390429 | Sz 94 | 16:07:49.60 | -39:04:29.0 | III |  M4.0j,l |  0.28l | $116.2\pm 0.9$ | - | - | <0.11 | <0.015 | -
| - | - | - | - | - | - | - | - | $10.4$ | C | $1.67\pm 0.17$ | - | -
20 | J160758.9-392435 | Lup 714 | 16:07:58.90 | -39:24:34.9 | III |  M5.0c |  0.05n | $162\pm 3$ | - | - | <0.11 | <0.028 | -
21 | J160816.0-390304 | SONYC Lup3-5 | 16:08:16.03 | -39:03:04.2 | III |  M6.5j,l |  0.05n | $164\pm 5$ | - | - | <0.11 | <0.029 | -
22 | J160834.6-390534 | V1027 Sco | 16:08:34.60 | -39:05:34.0 | III |  A3.0c |  3.16n | $161.5\pm 2.2$ | - | - | <0.11 | <0.028 | -
| - | - | - | - | - | - | - | - | $5.5$ | C | $1.09\pm 0.15$ | - | -
23 | J160848.2-390419 | Sz 109 | 16:08:48.16 | -39:04:19.2 | III |  M5.0c |  0.05n | $184\pm 5$ | - | - | <0.11 | <0.036 | -
24 | J160908.5-390343 | Lup 608s | 16:09:08.50 | -39:03:43.1 | III |  M5.0c |  0.08n | $160.3\pm 2.4$ | - | - | <0.11 | <0.028 | -
25 | J160917.1-392710 | Lup 710 | 16:09:17.13 | -39:27:09.7 | III |  M5.0c |  0.03n | $131\pm 3$ | - | - | <0.11 | <0.019 | -
26 | J160942.6-391941 | Sz 116 | 16:09:42.57 | -39:19:40.8 | III |  M1.5c | $0.43\pm 0.10^{\rm{m}}$ | $146.5\pm 0.8$ | - | - | <0.11 | <0.023 |  Binaryd
27 | J161012.2-392118 | Sz 121 | 16:10:12.21 | -39:21:18.3 | III |  M3.0j,l | $0.31\pm 0.08^{\rm{m}}$ | $158.0\pm 1.7$ | - | - | <0.11 | <0.028 | -
28 | J161016.4-390805 | Sz 122 | 16:10:16.44 | -39:08:05.4 | III |  M2.0j,l | $0.39\pm 0.05^{\rm{m}}$ | $138.9\pm 0.6$ | - | - | <0.11 | <0.021 | -
29 | J161153.4-390216 | Sz 124 | 16:11:53.35 | -39:02:16.1 | III |  M0.0c |  0.75n | $125.4\pm 0.7$ | - | - | <0.11 | <0.017 | -
30 | J161302.4-400433 | V1097 Sco | 16:13:02.416 | -40:04:32.93 | III |  K7.0f | $0.80\pm 0.13^{\rm{m}}$ | $139.9\pm 2.5$ | - | - | <0.11 | <0.021 |  Binaryd
Table 3: Modelled parameters and constraints for the unresolved class III detections of MU Lup, NO Lup and Sz 107 in the upper section, and the modelled parameters for the resolved detection of J155526.2-333823 in the lower section. Source | $I_{0}$ | $\sigma$ | $R_{95}$ | $R_{\rm{c}}$ | $\gamma_{1}$ | $\gamma_{2}$ | Inc | PA | $\Delta\rm{RA}$ | $\Delta\rm{Dec}$
---|---|---|---|---|---|---|---|---|---|---
|  | [arcsec] | [au] | [au] |  |  | [deg] | [deg] | [mas] | [mas]
MU Lup | $8.2\pm 0.3$ | $0.14\pm 0.07$ | $<53$ | - | - | - | $50\pm 30$ | $90\pm 50$ | $-50\pm 50$ | $-90\pm 30$
NO Lup | $8.0\pm 0.4$ | $0.19\pm 0.07$ | $<56$ | - | - | - | $50\pm 30$ | $120\pm 50$ | $20\pm 40$ | $30\pm 40$
Sz 107 | - | - | $<60$ | - | - | - | - | - | - | -
J155526.2-333823 | $7.58\pm 0.10$ | - | - | $80\pm 15$ | $3.3\pm 1.5$ | $3.8\pm 0.8$ | $31\pm 12$ | $138\pm 21$ | $-30\pm 60$ | $-90\pm 50$

In image 5 we did not detect emission from the class III star Sz 108, but
emission consistent with Sz 108’s class II binary companion, Sz 108B. Given
the $>4.0\arcsec$ angular separation these sources are not confused. This is
similarly true in image 6, in which the class III star Sz 119 was not
detected, yet the nearby (${\sim}9.8\arcsec$ separated) class II source was
(i.e., J16095628-3859518, which we refer to as Lup 818s herein, see López
Martí et al., 2005). In the five other binary systems noted, neither the
source nor its companion were detected. Further discussion of stellar
multiplicity in Lupus for these 6 binaries is provided in Zurlo, submitted
(2020), and the effect of binarity on circumstellar disk evolution in our
discussion section $\S$7.1.5.

###  5.2 Continuum Source Analysis

Figure 4: Model fits to the ALMA image of J155526.2-333823. (Left) ALMA
cleaned image. The contours represent emission at $\pm 2$ and $\pm 3$ times
the rms. (Centre) Best-fit model image based on values in Table 3 with
contours set to demonstrate the ring-like nature of the model. (Right)
Residual image (data minus model), with contours set at $\pm 2$ and $\pm 3$
times the rms, and the dashed lines showing the model contours for reference.
In all images North is up and East is left, and the ellipse in the lower left
corner represents the size of the synthesised beam (${\sim}90\times
80\rm{au}$).

####  5.2.1 Continuum Source Modelling

The class III (type A) detections were unresolved for MU Lup, NO Lup and Sz
107, but resolved for J155526-333823. To find limits on the radii of the
unresolved sources, and disk parameters for the resolved source, we used
GALARIO (Tazzari et al., 2018) to model the emission profiles of the 4 sub-mm
detections.

For the unresolved observations the upper limit radii of the sources was
determined using a simple parametrised Gaussian model defined as

$I(r)=I_{0}\exp{\Big{[}-\Big{(}\frac{r}{2\sigma}\Big{)}^{2}\Big{]}},$ (1)

where $I_{0}$ is a normalisation factor, r is the deprojected sky radius, and
$\sigma$ the Gaussian width of the profile. These parameters were fitted along
with the positional offsets in RA and Dec, inclination (where $i=0^{\circ}$
corresponds to a face-on disk) and position angle (with positive values taken
as anti-clockwise from North). By fitting with such a model, we find the
average radius within which $95\%$ of the source flux is bounded, $R_{95}$,
which we define as the radius limit for our unresolved sources (in a manner
consistent with Sanchis et al., 2020). For the resolved observation of
J155526.2-333823, we followed the same procedure, but with a modified self-
similar disk model with free parameters of the disk radius ($R_{\rm{c}}$) and
inner and outer exponents ($\gamma_{1}$ and $\gamma_{2}$), defined as

$I(r)=I_{0}\Big{(}\frac{r}{R_{\rm{c}}}\Big{)}^{\gamma_{1}}\exp{\Big{[}-\Big{(}\frac{r}{R_{\rm{c}}}\Big{)}}^{\gamma_{2}}\Big{]}.$
(2)

All the values determined from our emcee modelling (for which 100,000 steps
were used and either 36 or 48 walkers for the unresolved/resolved source
modelling respectively Foreman-Mackey et al., 2013) are shown in Table 3. For
Sz 107 the nearby bright type C source (likely background) means that the disk
parameters found by GALARIO are much less accurate, so we use the results for
NO Lup to set an upper limit on the disk radius given their similar flux
measurements (taking account of their different distances).

The constraints on the disk sizes for MU Lup and NO Lup (and estimate for Sz
107) show that this emission arises from dust within ${\sim}60$ au from their
stars (and is entirely consistent with the size scales found by considering
the FWHM of the source emission and the beam). Although the centroids of MU
Lup and NO Lup are both well fitted by GALARIO, their inclinations and
position angles are essentially unconstrained (as expected with these being
unresolved). On the other hand, the disk of J155526.2-333823 is sufficiently
large and bright that it was resolved and is well constrained by GALARIO,
putting relatively tight bounds on its parameters (i.e., the radius, PA and
inclination). This disk is the largest detected, with a radius
${\sim}80\rm{au}$, significantly greater than the blackbody radii obtained
from SED fitting (see Table 4), where the positive $\gamma_{1}$ determines
this to be a ring with an inner gap. Fig. 4 shows the ALMA image, the disk
model that was subtracted from this, and the resulting residual image, within
which the lack of residuals $>3\sigma$ demonstrates the reasonable fit.

We note here that the aim of the survey was to maximise the S/N of the
detections. As such, the ${\sim}0.8\arcsec$ resolution was chosen that the
majority should be unresolved, explaining why only 1 of our 4 detections was
resolved. Follow-up high resolution measurements of the three unresolved
sources (and additionally the resolved source) would be needed to gather much
tighter constraints on the radii, widths, and other morphological features of
these disks.

####  5.2.2 Continuum Dust Masses

Based on the assumption that the type A and B detections are from optically
thin, isothermal circumstellar dust, we estimated their dust masses (in
$M_{\oplus}$) as

$M_{\rm{dust}}[M_{\oplus}]=\frac{F_{\nu}d^{2}}{\kappa_{\nu}B_{\nu}(T_{\rm{dust}})}\approx
0.22\Big{(}\frac{d}{150}\Big{)}^{2}F_{856\mu\rm{m}},$ (3)

where $d$ is distance in $\rm{pc}$ and $F_{856\mu\rm{m}}$ is the Band 7 flux
in mJy, for an assumed dust temperature of $T_{\rm{dust}}{\sim}20\rm{K}$
(e.g., the median for Taurus disks, see Andrews & Williams, 2005) and a dust
grain opacity $\kappa_{\nu}$ of $10\rm{cm}^{2}\,\rm{g}^{-1}$ at 1000GHz, with
an opacity power-law index $\beta=1$ (e.g., Hildebrand, 1983; Beckwith et al.,
1990). These dust masses are noted in the Table 2. For the class III disks
these masses range between $0.036-0.093M_{\oplus}$, and for the class II disks
they are $6.5$ and $1.63M_{\oplus}$ respectively for sources Sz 108B (5B) and
Lup 818s (6B). Given Ansdell et al. (2016) assumed these 2 class II sources to
be at $200\rm{pc}$, whereas here we have calculated their masses using their
Gaia DR2 parallax distances, our dust mass measurements are smaller than those
previously derived. After accounting for the different distances, these are
consistent. For the type A sources with $3\sigma$ upper limits on their
fluxes, eq. 3 is used to provide the $3\sigma$ upper limits on their masses.
Given the type C detections do not coincide with known objects, these are most
likely emission from background sub-mm extragalactic (SMG) sources, common in
deep sub-mm observations (discussed further in $\S$5.5), rather than stellar
sources within Lupus. We therefore do not provide dust mass estimates for type
C detections.

Figure 5: Stacked image of all 26 sources in scale of significance (with respect to the rms of the combined image stack, ${\sim}7.2\mu\rm{Jy}\,\rm{beam}^{-1}$). Within the inner ${\sim}0.68\arcsec$ radius circle there is a $2.8\sigma$ significant detection. Note that all other sources with significant emission are the non-type A detections as discussed in $\S$5. In this stacked image, North is up, East is left. |
---|---
|
|
Figure 6: SEDs for MU Lup (SpT: K6.0): top-left, J155526.2-333823 (SpT: K5.0):
top-right, NO Lup (SpT: K7.0): middle-left, Sz 107 (SpT: M5.5): middle-right,
MV Lup (SpT: K2.0): bottom-left, and CD-31 12522 (SpT: K2.0): bottom-right.
Each plot shows the fit to the stellar photosphere and the circumstellar
emission as a single modified blackbody (solid line). In the case of MU Lup a
single black body fit is plotted (red dashed line), and in the cases of
J155526.2-333823, NO Lup and Sz 107, two-component modified blackbody fits are
plotted (red dashed lines). The dots and triangles show the available
photometry and upper limits.

###  5.3 Survey Non-Detections

For the class III sources for which sub-mm emission was not detected, we
produced a stacked image of their continuum images to constrain the average
dust mass. The 26 images for which this was true combined to form the image
shown in Fig. 5. Noticeable in this image are the B and C type detections that
fall within the fields of view, however we are only concerned with the inner
${\sim}0.68\arcsec$ that is coincident with the class III stars (within one
beam), as highlighted in the expanded central region. By integrating the
emission within this inner region we found a total of $20.2\pm 7.2\mu\rm{Jy}$,
i.e., a tentative ${\sim}2.8\sigma$ detection coincident with all of the class
III stars which individually have non-detections.

The mean distance to the stars in our survey is $141{\rm{pc}}$, and so
equation 3 shows that this tentative stacked detection (if from circumstellar
dust) corresponds to a dust mass of $M_{\rm{dust}}=0.0048\pm
0.0017M_{\oplus}$. This implies that whilst we only detected significant
emission to 4/30 class III stars (i.e., type A sources), we may have detected
further low intensity emission from the other sources. In the context of known
dust masses, this would be a limit far below any previous measurements for
stars of this age, and nearly an order of magnitude below the lowest mass type
A detection in $\S$5.2.2.

###  5.4 SED Analysis

To investigate the nature of the infrared emission seen in the class III
stars, Spectral Energy Distributions (SEDs) for all sources were produced
using the methodology of Yelverton et al. (2019), i.e., we fit a stellar
photosphere to determine $L_{\star}$ and $T_{\rm{eff}}$, plus a blackbody disk
where merited, to the available photometry. We include mid- and far-IR data
from WISE, Spitzer and Herschel, but also data at shorter wavelengths that
constrain the stellar photospheric emission. All 30 sources had photometry in
J-, H- and K-bands (2MASS; Cutri & et al., 2013), Gaia G, $B_{P}$ and $R_{P}$
bands (Gaia Collaboration et al., 2018), WISE bands W1, W2, W3 and W4 (Wright
et al., 2010), and the following number in each instrument and band: 28/30 in
both IRAC 3.6$\mu$m and IRAC 4.5$\mu$m and 29/30 in both IRAC 5.8$\mu$m and
IRAC 8.0$\mu$m, 29/30 in MIPS 24$\mu$m (Werner et al., 2004; Rieke et al.,
2004; Merin et al., 2009), 2/30 in PACS 70$\mu$m (Pilbratt et al., 2010;
Poglitsch et al., 2010; Cieza et al., 2013), alongside the ALMA Band 7
(${\sim}856\mu$m) measurements presented in this work.

We found from the analysis of the full SEDs, that 24 of the sources are well
explained with purely photospheric emission; i.e., there is no evidence for an
additional circumstellar component in addition to their stellar photospheres.
An additional emission component is required for the remaining 6 sources (MU
Lup, J155526.2-333823, NO Lup, Sz 107, MV Lup and CD-31 12522), which includes
the 4 for which sub-mm continuum emission was detected, which we consider to
be circumstellar emission. The SEDs of the 6 with excess emission are shown in
Fig. 6, and the 24 purely photospheric SEDs in Appendix A.

The SED modelling process resulted in a fit to the stellar photospheric
emission, predominantly constrained by the shorter wavelength data. In
addition to being essential for predicting the stellar contribution at longer
wavelengths, we used the resulting stellar parameters $L_{\star}$ and
$T_{\rm{eff}}$ to estimate the stellar masses shown in Table 2 for the 14/30
class III stars for which values were not available in the literature (i.e.,
those denoted with superscript "n"). To do so, we applied the models of Siess
et al. (2000, for V1027 Sco) and Baraffe et al. (2015) for the 13 others,
except Sz 109 and J160816.0-390304 for which estimates based on their spectral
types were required (i.e., since these were poorly fitted by stellar models).
In all cases where the stellar mass uncertainty is not stated in Table 2
(i.e., in the literature or as derived here), we assume this to be 10%.

We first consider the $12\mu\rm{m}$ and $24\mu\rm{m}$ excesses, given the
discussion of [K]-[12] and [K]-[24] in $\S$3 and the colour-colour diagram in
Fig. 2. For the 6 sources with excess emission, these are all found to have
significant excesses at 24$\mu$m. Significant excesses were found at 12$\mu$m
for 4 sources (MU Lup, NO Lup, MV Lup and CD-31 12522), while both
J155526.2-333823 and Sz 107 were found to have tentative 12$\mu$m excesses
($2.6\sigma$ and $2.2\sigma$ respectively). Given that MU Lup,
J155526.2-333823, MV Lup and CD-31 12522 were not deemed as outliers in Fig. 2
(i.e., these sit with the majority of class III sources in the lower-left
cluster), this illustrates that low [K]-[12] or [K]-[24] indices do not
necessarily rule out the presence of excess emission. It is also notable that
4/6 of the class III sources with a 24 $\mu$m excess were found to have
detectable sub-mm emission in the ALMA survey. This suggests that the presence
of a mid-IR excess may raise the sub-mm detectability for class III YSOs. We
discuss the implications of this warm emission further in $\S$7.

The excess emission for these 6 sources was first fitted with a single-
component modified blackbody, parameterised by a temperature, T, fractional
luminosity, $f=L_{\rm{disk}}$/$L_{\star}$, $\lambda_{0}$ and $\beta$ (the
latter two parameters describe the location and gradient beyond which the
modification to the blackbody fit is implemented). This temperature is readily
converted to a black body radius, Rbb, by assuming this emission is from
optically thin dust in thermal equilibrium with the stellar radiation of
luminosity $L_{\star}$. Whilst this radius provides an indication of the size
of the dust disk, it is expected that this radius will be systematically
underestimated (e.g., since the dust is likely to be small enough to emit
inefficiently and so be hotter than a black body, see Pawellek et al., 2014;
Pawellek & Krivov, 2015). These single-component blackbody SED fits are shown
as brown solid lines ($F_{\rm{\nu,single}}$) on Fig. 6, and the model fits for
the temperature, associated blackbody radius and fractional luminosity are
shown in the upper part of Table 4. We have not included the $\lambda_{0}$ and
$\beta$ values as these are poorly constrained and highly degenerate. In all 6
cases $\beta$ can be 0, for which $\lambda_{0}$ is then meaningless and
unconstrained. However, the joint 2D parameter space is somewhat constrained
for NO Lup and Sz 107 because the dust temperatures are well constrained by
their IRS spectra (Lebouteiller et al., 2011). For MU Lup, J155526.2-333823,
MV Lup and CD-31 12522 the SED of the excess emission is poorly sampled which
means that a wide range of single-component modified black body fits can model
the data by choosing different $\lambda_{0}$ and $\beta$. In the case of MU
Lup, the solid line is plotted with $\lambda_{0}{\sim}120\mu$m and
$\beta=1.2$, whereas the dashed line is a pure black body (for which the same
degeneracy exists in the cases of J155526.2-333823, MV Lup and CD-31 12522).
However, despite this allowed range of model fits, a pure black body model
(i.e., $\beta=0$) provides the warmest temperature fit (i.e., this corresponds
to the lowest radius). For MU Lup, such a pure black body radius is at 10 au.
Therefore, accounting for the range of allowed models, this suggests that the
emission measured for MU Lup is likely from a component at tens of au.

Table 4: SED blackbody model parameters for all class III stars determined to
have excess emission above the stellar photosphere (note J155526.2-333823 is
written here as J155526). The upper part of the table includes the 4 objects
with sub-mm detections and 2 objects undetected in the sub-mm but with
archival mid-IR observations found evidence for a warm component, and their
respective modified blackbody fits. The lower part of the table includes the
three sources determined as well-fitted by a two-component modified blackbody,
though we note here the broad parameter space that can be fitted by these
modified black bodies, and therefore whilst errors are defined, these are only
weakly constrained.

Single Component Modelling
Source T Rbb $f$ [K] [au] $\times\,0.01$ MU Lup $50\pm 10$ $16\pm 9$ $0.1\pm
0.5$ J155526 $52\pm 8$ $14\pm 5$ $0.3\pm 0.3$ NO Lup $117\pm 2$ $3.1\pm 0.1$
$0.34\pm 0.01$ Sz 107 $144\pm 3$ $1.16\pm 0.05$ $1.39\pm 0.03$ MV Lup $600\pm
200$ $0.2\pm 0.1$ $0.09\pm 0.06$ CD-31 12522 $750\pm 40$ $0.13\pm 0.01$
$0.13\pm 0.05$
Two Component Modelling
Source $T_{1}$ $R_{\rm{bb1}}$ $f_{1}$ $T_{2}$ $R_{\rm{bb2}}$ $f_{2}$ [K] [au]
$\times\,0.01$ [K] [au] $\times\,0.01$ J155526 >50 <15 0.01-1 20-60 10-100
0.01-1 NO Lup 117 3.1 0.34 25-120 2-56 0.00002-0.4 Sz 107 144 1.16 1.39 20-150
1-60 0.00005-4

Although their single-component black body radii are poorly constrained, the
16 au radius found for MU Lup is consistent with this source being unresolved
in the sub-mm (i.e., <53 au in $\S$5.2.1), and at $\leq 0.2$ au for both MV
Lup and CD-31 12522 it is not surprising that these were not detected in the
sub-mm. However the single modified black body fits for J155526.2-333823, NO
Lup and Sz 107 are problematic. For J155526.2-333823, whilst this is found
with a single black body radius of ${\sim}14$ au, this source was resolved
with a sub-mm radius of ${\sim}80$ au (i.e., more than a factor of 5 larger),
and both NO Lup and Sz 107 (as accurately measured with IRS spectra, see
Lebouteiller et al., 2011) are found to have well-constrained temperatures
that result in blackbody radii less than ${\sim}3$ au. In the case of Sz 107,
it can also be seen that the single component modified black body is not as
well fitted to both the sub-mm measurement and the warm emission. For these
reasons, we therefore modelled J155526.2-333823, NO Lup and Sz 107 as having
two separate modified black body components, representing both warm and cool
emission. While such an additional component is not required to fit the SED,
this is required to explain the excess emission seen towards many older main
sequence stars (e.g., Kennedy & Wyatt, 2014) and so it is not unreasonable
that these stars might exhibit two-component SEDs. Thus it is important to
consider how such an additional component might affect our estimate of the
radii of these disks.

The two-component models for J155526.2-333823, NO Lup and Sz 107 are under-
constrained in that their SEDs can be fitted by a range of two-component
blackbody models, as illustrated by the shaded regions on the plots in
Appendix B, where further discussion of their interpretation is provided. In
the lower part of Table 4 we provide ranges of allowed temperatures,
associated blackbody radii and fractional luminosities. Fig. 6 shows examples
for each system of one such possible two-component model with red dashed
lines, with the warm component modelled to overlap with emission
$\lambda<100\mu$m and the cool component modelled to overlap with the sub-mm
detection. Both NO Lup and Sz 107’s warm components are tightly constrained,
so match the results of the single modified blackbody fits, while the sub-mm
emission originates in a cooler component that is consistent with being at 10s
of au from the star (though must both still remain unresolved in our ALMA
observations at $<56$ au and $<60$ au respectively).

In summary, despite the limitations of modelling the SEDs for the 6 sources
shown in Fig. 6, we have shown that reasonable single modified black body
models are found for all sources, alongside alternative 2-component modified
black body models for J155526.2-333823, NO Lup and Sz 107. This modelling has
shown that disks with black body radii of tens of au remain consistent with
all four sub-mm detections (i.e., the circumstellar material detected by ALMA
may be >10s of au away from the star).

Figure 7: Sky map of the far-IR galactic dust extinction in the Lupus star
forming region, where the locations of the class III stars in this survey have
all been noted. Stars represented by "o" markers indicate that they are
located in Lupus cloud I, whereas markers represented by a filled circle
indicate the star is in Lupus cloud III. White markers indicate that
background CO emission was detected in the channel maps, whereas green markers
indicate that no background CO was detected, discussed further in $\S$ 6.
Figure 8: ALMA CO moment 0 images for all 30 class III stars. The colour scale
is in units of significance with respect to each image’s rms. The source
numbers in the bottom-left of each panel correspond to those in Table 2. Note
that the plots in which rms is higher were those which suffered from CO cloud
contamination, further discussed in $\S$ 6. In all images, North is up, East
is left, and the beam and 200 au scale bar are shown in the upper and lower
right respectively.

###  5.5 Background Sources

In our sample we find 11 type C detections, 10 of which have fluxes $>1.0$
mJy, and 1 with a flux ${\sim}0.41$ mJy. We previously inferred these to be
sub-mm galaxies (SMGs) in $\S$5.1, given that their locations are inconsistent
with other known stellar sources in Lupus. This detection rate can be compared
with that expected from surveys designed to determine SMG number counts at
${\sim}850\mu$m (Simpson et al., 2015; Popping et al., 2020). To do so, we
have to account for the primary beam correction, which means that while we
could detect sources brighter than 5 times $36.7$ $\mu$Jy (i.e., $\approx 0.2$
mJy) towards the centre of the image, sources would need to be ${\sim}2$ times
brighter than this at a radius of $8.4\arcsec$, and ${\sim}5$ times brighter
than this at $12.5\arcsec$ (i.e., at $\sim 0.4$ mJy and $\sim 1.0$ mJy
respectively).

The Popping et al. (2020) number counts predict that we should detect
${\sim}8$ sources above 1 mJy inside $12.5\arcsec$ and ${\sim}8$ sources above
0.4 mJy inside $8.4\arcsec$ (all of which are scaled to cover the total survey
area, i.e., inside an area equal to 30 circles with these three radii).
Compared with the number of type C sources we measured in each of these
regions, 10 and 6 respectively, these numbers therefore appear strongly
consistent. At lower flux levels it may be notable that we only detected 1 SMG
in the range $0.4-1$ mJy inside $8.4\arcsec$, whereas the Popping et al.
(2020) number counts predict we should have detected ${\sim}4$. Given that
these numbers are both low, we do not ascribe great significance to this
difference, but note that a significant difference might be evidence for
cosmic variance or a turn-over in SMG number counts at low flux levels.
Indeed, our survey could be considered as an independent measurement of
856$\mu$m SMG number counts in this region of sky down to 1 mJy and indeed
below, even if its interpretation is complicated by the possibility of
Galactic sources in the field of view.

We further note that the 4 sub-mm detections towards class III stars (i.e.,
type A) were all detected within $0.4\arcsec$ of their corresponding ALMA
phase centres, with fluxes exceeding 0.2 mJy. Inside this region over the
entire survey, the Popping et al. (2020) number counts correspond to an SMG
detection probability of just ${\sim}3\%$. The contamination rate may be even
lower than this, given that our survey measured a lower number of $>0.4$ mJy
sources than predicted by Popping et al. (2020). Thus, whilst it is possible
that one of these class III (type A) detections is an SMG, this also appears
highly unlikely.

##  6 Survey CO Analysis

In this section we use the ALMA observations to search for CO emission
coincident with the class III stars, and with the two class II sources that
fall in the field of view.

###  6.1 Survey Detections

To produce image cubes for each measurement set, we first transformed this
data to the barycentric reference frame, and subtracted a linear fit of the
continuum that was fitted over line-free channels (i.e., those external to the
region $\pm 15$ km s-1 from the CO J=3-2 spectral line) using the CASA task
$\rm{uvcontsub}$. These were then all imaged with the CASA task $\rm{tclean}$
non-interactively, for only the initial major cycle to produce data cubes,
with a rest frequency set to the CO J=3-2 spectral line frequency,
$f_{\rm{CO}}=345.79599$ $\rm{GHz}$. The FOV was likewise chosen to show the
entire primary beam.

####  6.1.1 Cloud Contamination

Channel maps were produced from the sub-cubes of all survey sources and in
16/30 of these high intensity, large scale CO emission was observed in the
field of view in the channels in the range $f_{\rm{CO}}\pm 10$ km
$\rm{s}^{-1}$. In none of these cubes were emission signals found further away
than a few channels from this frequency at this intensity. For all 30 sub-
cubes, velocities from +13 to -16 km s-1 are shown in Appendix C in Figures 19
and 20, with 16 of these demonstrating large scale CO J=3-2 emission unrelated
to any detected continuum emission (in images 4, 5, 6, 7, 12, 15, 19, 20, 21,
22, 23, 24, 26, 27, 29, 30). The full cubes for these 16 sources were
investigated and the intensity and scale of the observed CO emission was found
to be inconsistent with being due to flagged data, missing baselines, or
visibility values, and found irrespective of the $\rm{uvcontsub}$ task being
applied.

Table 5: Shown are the total CO fluxes, estimated CO gas masses for an optically thin, LTE (50K) assumption, and whether or not CO cloud emission was present in the data. No. | Name | Flux (CO) | $M_{\rm{CO,LTE}}$ | Cloud
---|---|---|---|---
|  |  [$\rm{Jy}\,\rm{km}\,\rm{s}^{-1}$] |  [$10^{-5}\times M_{\oplus}$] | Gas
1 | MU Lup | <0.12 | <2.1 | -
2 | J155526.2-333823 | <0.12 | <1.6 | -
3 | NO Lup | $0.29\pm 0.07$ | $4.9\pm 1.1$ | -
4 | Sz 107 | <0.19 | <4.2 | Y
5 | Sz 108 | <0.19 | <4.2 | Y
| Sz 108B | $0.56\pm 0.06$ | $12.5\pm 1.4$ | Y
6 | Sz 119 | <0.19 | <2.1 | Y
| Lup 818s | $0.037\pm 0.012$ | $0.43\pm 0.14$ | Y
7 | MT Lup | <0.19 | <3.2 | Y
8 | Sz 67 | <0.12 | <1.4 | -
9 | J154306.3-392020 | <0.12 | <3.2 | -
10 | CD-35 10498 | <0.12 | <2.6 | -
11 | MV Lup | <0.12 | <2.3 | -
12 | MW Lup | <0.19 | <3.1 | Y
13 | MX Lup | <0.12 | <1.9 | -
14 | NN Lup | <0.12 | <2.2 | -
15 | CD-39 10292 | <0.19 | <3.1 | Y
16 | CD-31 12522 | <0.12 | <2.0 | -
17 | J160713.7-383924 | <0.12 | <2.9 | -
18 | J160714.0-385238 | <0.12 | <1.7 | -
19 | Sz 94 | <0.19 | <2.4 | Y
20 | J160758.9-392435 | <0.19 | <4.7 | Y
21 | J160816.0-390304 | <0.19 | <4.7 | Y
22 | V1027 Sco | <0.19 | <4.7 | Y
23 | Sz 109 | <0.19 | <6.0 | Y
24 | J160908.5-390343 | <0.19 | <4.6 | Y
25 | J160917.1-392710 | <0.12 | <2.0 | -
26 | Sz 116 | <0.19 | <3.9 | Y
27 | Sz 121 | <0.19 | <4.4 | Y
28 | Sz 122 | <0.12 | <2.2 | -
29 | Sz 124 | <0.19 | <2.9 | Y
30 | V1097 Sco | <0.19 | <3.5 | Y

We plotted the positions of all 30 class III stars on Fig. 7, overlaid onto an
image of far-IR Milky Way galactic dust extinction, obtained from the "SFD"
dustmap as described in Schlegel et al. (1998) and updated by Schlafly &
Finkbeiner (2011). This can be used to indicate whether galactic dust or CO
gas might have contaminated our observations near the CO J=3-2 transition
frequency (discussed later in $\S$7), and showed a correlation between the
dust extinction coefficient and the presence of CO emission.

We therefore conclude that such large scale CO emission is highly likely from
galactic CO gas, present in large abundance in the Lupus star forming region,
and detected by our deep ALMA imaging (although we note that the intensity and
extent of the gas emission may not be well traced by our images). All 16 of
these sources are noted in Table 5 in the "Cloud Emission" column with a "Y"
if the channel maps were affected. We note that cloud contamination was
reported in the class II Lupus survey of Ansdell et al. (2018), however the
depth of the observations in this class III survey has resulted in cloud
emission being much brighter than in any previous ALMA imaging of Lupus.

Given the size scale of the emission, in the 16 measurement sets where CO
cloud contamination was detected we implemented a uv-baseline cut during the
imaging process to remove all visibility data with $\rm{uv-
distance}<100\rm{k}\lambda$, significantly reducing the structure from
emission on large scale sizes ($>2\arcsec$, i.e., larger than any CO emission
that would be circumstellar). Moment-0 maps for all 30 sources were then
produced from emission within 10 km s-1 of the CO spectral line, shown in Fig.
8. These CO moment-0 maps show significant source detections in positions
coincident with continuum emission sources 3A (NO Lup, i.e., one of the class
III sources) and both class II sources, Sz 108B (5B) and Lup 818s (6B). No CO
emission was detected within $1.0\arcsec$ of any other class III sources or
type-C detections from our continuum maps.

Figure 9: CO spectra of the class III star NO Lup (left), and the class II
stars Sz 108B (middle) and Lup 818s (right). The red, blue and green dashed
lines show the per-channel $1\sigma$, $3\sigma$ and $5\sigma$ noise levels,
respectively, calculated from the rms signal of the integrated flux in
channels distant from the CO emission.

###  6.2 Class II and III CO Analysis

Of the class III stars for which sub-mm continuum emission was detected (i.e.,
those designated type A), only 1 source was affected by the CO cloud
contamination, Sz 107, for which we can only set an upper limit on the CO gas
flux/mass (see later in this section). MU Lup and J155526-333823 showed no
significant CO emission coincident with the star, and we can only derive upper
limit fluxes/masses for these two sources also. NO Lup, however, demonstrates
highly significant CO emission, in fact at higher significance than its
continuum detection. Cloud contamination affected some channels for both of
the class II sources (Sz 108B and Lup 818s), however the brightness of these
two sources meant that their CO flux was still detectable above the background
emission. We note here that it is also possible that cloud contamination has
masked low level detections of CO gas in sources 19, 23, 24, and 30, given
that these appear most affected at the class III stellar position. Further
high resolution imaging of these sources may reveal if CO gas is present.

Spatially integrated spectra were produced for NO Lup, Sz 108B and Lup 818s,
and are shown in Fig. 9. For the NO Lup spectrum, emission was integrated
within an elliptical mask with an aperture of $(2.5\times 2.1)\arcsec$ (i.e.,
such that the mask fully enveloped all significant CO emission). For the class
II spectra, emission was integrated within elliptical masks centred on the two
class II sources with apertures of $(1.6\times 1.4)\arcsec$ and $(1.2\times
1.0)\arcsec$ for Sz 108B and Lup 818s, respectively. CO fluxes were determined
by integrating the area under these spectra between the points at which these
cross their respective $1\sigma$ thresholds; -11.0 to +8.1
$\rm{km}\,\rm{s}^{-1}$ for NO Lup, -9.0 to +2.0 $\rm{km}\,\rm{s}^{-1}$ for Sz
108B, and -5.4 to +1.0 $\rm{km}\,\rm{s}^{-1}$ for Lup 818s. The median
velocity of the detected CO emission is consistent with that expected from the
radial velocities of Lupus stars, $\rm{RV_{Lup}}=2.8\pm 4.2\rm{kms^{-1}}$ (see
Frasca et al., 2017). For NO Lup, emission is measured over a broad range of
velocities $\sim{19}$ $\rm{km}\,\rm{s}^{-1}$. The spectrum may be consistent
with the double-peaked structure typical of spectrally resolved, inclined
rotating disks. Therefore, whilst the inclination of the NO Lup disk was not
well constrained from its continuum emission (see Table 3), the CO spectrum
suggests that this source may be inclined with respect to the plane of the
sky. In both Sz 108B and Lup 818s, the spectral line shows only a single sharp
peak; the lower negative peak of Sz 108B close to $-2\rm{km}\,\rm{s}^{-1}$ is
likely due to cloud contamination absorbing the red-shifted side of the disk.

To determine the uncertainty on the CO line fluxes, we found for all non-
contaminated sources, for an aperture with a size of $(1.2\times 1.0)\arcsec$
located at the image centre, the integrated emission uncertainty is
$\rm{rms_{CO}=3.0\,mJy}$ per 0.42 km s-1 channel (based on the emission
excluding channels within $\pm 20$ km s-1 of the spectral line). For the
sources affected by cloud contamination the rms is higher once uv-baseline
cuts have been applied, and for these we instead find the integrated emission
uncertainty $\rm{rms_{CO}=4.7\,mJy}$ per 0.42 km s-1 channel (i.e., the
uncertainty increased by ${\sim}50\%$). Although the channel width is $\Delta
v=0.42$ km s-1, the effective spectral bandwidth is ${\sim}2.667$ larger than
this, since adjacent ALMA channels are not fully independent 111For a complete
discussion, see
https://safe.nrao.edu/wiki/pub/Main/ALMAWindowFunctions/Note_on_Spectral_Response.pdf.
With these considerations, and by including an additional 10% flux calibration
error in quadrature, we find the integrated fluxes and uncertainties for NO
Lup, Sz 108B and Lup 818s as shown in Table 5. Further, for an assumed
spectral line width of ${\sim}5$ $\rm{km}\,\rm{s}^{-1}$, we place $3\sigma$
upper bounds on the CO flux of all 29 other class III survey sources as either
$<0.19$ or $<0.12\,\rm{Jy}\,\rm{km}\,\rm{s}^{-1}$ (i.e., depending on whether
these sources were affected by cloud contamination or not respectively).

By further assuming that the CO gas is in Local Thermodynamic Equilibrium
(LTE), has a low optical depth (i.e., $\tau{<}1$), and a temperature of
$50\rm{K}$, we derived CO gas masses for NO Lup, Sz 108B and Lup 818s, and
upper limit masses on the 29 of our survey target sources. These calculations
used equations 2 and 8 (for the J=3-2 transition) from Matrà et al. (2015),
for which all derived $M_{\rm{CO,LTE}}$ masses are shown in Table 5. If the CO
abundance (i.e., (x${}_{\rm{CO}}=\rm{[CO]}$/$\rm{[H_{2}]}$) has a
protoplanetary disk-like abundance (i.e., $10^{-4}$, as per Williams & Best,
2014), we find the total gas-to-dust mass ratio for NO Lup of $1.0\pm 0.4$.
For the class III sources with dust measurements but no CO counterpart, the
$3\sigma$ upper limit gas-to-dust ratios are <0.25, <0.16 and <0.60 for MU
Lup, J155526.2-333823 and Sz 107, respectively. A caveat to the above
calculations is that if the emission comes from a sufficiently compact region,
it may be optically thick, which would imply that the values stated in Table 5
could all be larger. Class II protoplanetary disks are usually optically thick
in CO, but this is not necessarily the case for the low mass class III disks
examined here, with dust masses comparable to debris disks (see later, Fig.
10). Reinforcing this point, we measure class III CO fluxes that are either
upper limits, or in the case of NO Lup, a factor of ${\sim}$5 times lower than
the average FCO of the Lupus class II sources (see Ansdell et al., 2016), and
similar in line luminosity to nearby, optically thin CO-bearing debris disks
(see Kral et al., 2017). Therefore, while higher angular resolution
observations are needed to better estimate disk sizes and optical depth, we
present optically thin masses as a reasonable estimate of the CO contents in
these low FCO class III sources, but note that the estimates for the gas-rich
class II sources (e.g., Sz 108B) are likely underestimated. A second caveat is
that the assumption of LTE may not apply if H2 is not abundant in these disks,
e.g. if CO gas is secondary, that is, released from planetesimals instead of
being primordial (e.g., Dent et al., 2014; Marino et al., 2016; Moór et al.,
2017; Matrà et al., 2017). If the CO gas is not in LTE, then CO gas mass
limits are typically much higher as shown by Matrà et al. (2015, 2018), but
this does not mean that the total gas-to-dust ratio is necessarily higher due
to the absence of H2. We discuss the implications of this further in
$\S$7.1.2.

###  6.3 Summary of CO Analysis

By inspecting the CO emission from the 30 sources in our survey, whilst we
found in 16 images large scale emission from presumed CO in the Lupus cloud,
we found that one of the class III stars with a continuum detection also has a
CO counterpart, NO Lup, that will be discussed further in (Lovell, submitted.,
2020). After accounting for the likely galactic CO emission from the Lupus
clouds, we also detected CO towards both class II stars for which we also
detected continuum emission, Sz 108B and Lup 818s. None of the other class III
stars or background/other (type C) detections contained CO emission, and upper
limit CO fluxes/masses were derived for all 29 class III sources. From these
measurements we have constrained the CO emission coincident with all class III
and II sources present in this survey.

##  7 Discussion

###  7.1 Summary of Circumstellar Material around Class III Stars

In $\S$3 we showed the full set of YSOs that comprise our Lupus sample, and
how the [K]-[12] and [K]-[24] photometry compares for the class II and III
sources. We showed that the class II and III sources could be distinguished
photometrically consistent with the Lada & Wilking (1984) classification,
since the class II’s all have either $\alpha_{\rm{K-24}}>-1.6$ or
$\alpha_{\rm{K-12}}>-1.6$, which is not the case for the class III’s (although
some of these only have upper limits to their mid-IR emission). For the
majority the former criterion using the 24 $\mu$m photometry is sufficient to
classify the sources (i.e., all except AKC2006-18). Furthermore, whilst the
class II sources showed a large range in the colour-colour values, the class
III sources were mostly confined to near-stellar colours. Nevertheless,
detailed consideration of their SEDs showed that low level non-photospheric
emission is present for 6/30 of the class III sources. Our analysis led us to
different conclusions about the nature of the circumstellar material for
different stars in the 30 class III YSOs.

####  7.1.1 Likely Debris Disks: MU Lup and J155526.2-333823

These both have mm-sized dust which is inferred to be at 10s of au, with one
of these resolved as a ring at 80 au. Neither demonstrates any CO J=3-2
emission. Upper limits on the CO gas masses imply a $3\sigma$ upper limit on
the gas-to-dust ratios of <0.25 and <0.16 (for MU Lup and J155526.2-333823
respectively) if these had abundances similar to protoplanetary disks. We
therefore consider these to be inconsistent with gas dominated primordial
disks. In addition, while these also fell inside the $\alpha_{\rm{IR}}{<}-1.6$
classification for class III YSOs at both $12\mu\rm{m}$ and $24\mu\rm{m}$, low
level excess emission is present at 24 $\mu$m for both and at 12 $\mu$m for MU
Lup. Neither shows evidence of ongoing accretion on to their stars. In
conjunction, we find these two sources consistent with being young debris
disks.

Figure 10: Disk mass as a function of age for known debris disks and
protoplanetary disk populations, alongside the 30 class III Lupus sources
outlined in this work. The vertical dotted line at 10 Myr is only an
approximate division between where protoplanetary and debris disks might be
expected, and the horizontal line is the approximate division between sources
in star forming regions that are considered class II and class III (noting
that counter examples exist on either side of these lines). We plot as a green
dashed line an illustrative path that a system may follow as it evolves from
being a protoplanetary disk to a debris disk, via the class III stage explored
in this work.

####  7.1.2 Uncertain if Debris Disks or Protoplanetary Disks: NO Lup and Sz
107

NO Lup is found to have both dust and gas. The SED of the dust component can
be fitted as a single temperature component which would put it at ${\sim}3$
au, but it is equally likely that there is an additional component at 10s of
au (or a broad disk extending that far) that dominates the sub-mm continuum
emission. Most surprising is the bright gas which has a high velocity width
($\sim 19\,\rm{km}\,\rm{s}^{-1}$) with a gas-to-dust mass ratio measured as
$1.0\pm 0.4$, based on the assumption of the CO gas being in LTE and having
x${}_{\rm{CO}}=0.0001$. Whilst such assumptions are reasonable for
protoplanetary disks, they may not be for NO Lup. For example, if the LTE
assumption is invalid (i.e., in the case that the observed CO gas is secondary
and instead replenished by planetesimals), so too may be the assumption that
x${}_{\rm{CO}}=0.0001$ (i.e., since primordial H2 gas may have already
dispersed, and in such a system xCO could be much higher). Therefore, whilst
non-LTE conditions can result in higher CO gas masses (as previously discussed
in $\S$6.2), the physics involved which violate such conditions (i.e., the
dispersal of H2 gas) can produce lower total gas masses. If this is the case
for NO Lup (and indeed the other 3 sub-mm detected sources), then this total
gas-to-dust mass ratio could be much lower than calculated, and the disk would
instead be dust dominated. To assess how much second generation gas we may
expect around NO Lup, we can compare it with the $\sim 10$ Myr-old M-dwarf,
TWA 7, where second-generation CO gas was recently detected (Matrà et al.,
2019). If planetesimals in the NO Lup disk have the same volatile content as
those around TWA 7, we can use the properties of its dust disk (setting the
rate at which planetesimals are destroyed) to determine the level of second
generation CO that may be produced. Accounting for the differences in stellar
mass and luminosity, fractional luminosity, belt radius and width (assuming
$r{\sim}30$ au and a dr/r of 0.5 for NO Lup), we would expect ${\sim}18$ times
more CO gas mass around NO Lup than TWA 7. Thus, by taking the same LTE,
optically-thin excitation conditions for TWA 7, we would then predict
${\sim}2.4\times 10^{-5}M_{\oplus}$ of second-generation CO to be present
around NO Lup. This is only a factor ${\sim}2$ lower than the CO gas mass
measured (see Table 5), and thus these estimates show broad consistency
between the measured CO gas mass and a scenario in which this may be
continuously replenished from exo-cometary ice (i.e., one in which the CO gas
may no longer be primordial). This will be explored further in Lovell,
submitted. (2020).

Sz 107 has mid-IR colours that puts it close to other class II transition
disks, but its sub-mm dust mass is orders of magnitude lower (the dust masses
of nearby [K]-[12] and [K]-[24] class II YSOs Sz 84 and Sz 111 are ${\sim}7.7$
and ${\sim}75M_{\oplus}$, respectively). This suggests that whilst the mm-
grains in this system have dispersed, modest levels of mid-IR emission persist
(i.e., at similar levels to much more massive transition disks). Sz 107
therefore may provide a useful constraint on the rate at which this rapid mass
dispersal takes place (i.e., since this is likely to have happened relatively
recently), and also on the order, since it reinforces that the mm-sized grains
in the outer disk seem to disappear before the mid-IR emission from the inner
regions is finally cleared (e.g., Wyatt et al., 2015). Like NO Lup its dust
emission could be fitted by a single temperature component within a few au,
but there could be an additional component at 10s of au (or a broad disk).
There is no evidence for gas, but the limits are not as strong due to
confusion from cloud emission. Nevertheless, we could only place a $3\sigma$
upper limit on the total gas-to-dust ratio of <0.6 (i.e., whilst it is likely
this system is dust dominated, it may still have a large reservoir of gas
comparable in mass to its dust). Neither of these sources showed evidence of
ongoing accretion on to their stars.

We find these two sources therefore cannot be as clearly defined, being debris
disk-like in terms of their dust mass and accretion properties, but having
[K]-[24] excesses that are only slightly lower than known class II disks or
transition disks. Moreover, in the case of NO Lup, the gas detection is
consistent with a protoplanetary disk composition and a gas-to-dust ratio of
unity, but is also with the gas being secondary. These two systems are
therefore prime candidates for further deep ALMA imaging to constrain their
physical sizes, CO gas masses and origins in further detail.

####  7.1.3 Hot Dust from Asteroid Belts or Planet Formation: MV Lup and CD-31
12522

Two stars have 24 $\mu$m-only excesses, and neither showed evidence of ongoing
accretion on to their stars. The SED fitted temperatures puts the dust in
these systems within 1 au, which could be from a close-in asteroid belt
(Wyatt, 2008; Chen et al., 2009; Morales et al., 2009; Kennedy & Wyatt, 2014),
or dust created during collisions between forming planets (e.g., the shake-
down of Kepler planetary systems Kenyon & Bromley, 2004; Genda et al., 2015;
Wyatt & Jackson, 2016; Su et al., 2020). We therefore place these in a
separate category to the other class III YSOs, since while the mid-IR excess
could be an indicator that planetesimal formation has taken place (and/or that
planet formation is ongoing), it is not clear that this has taken place beyond
a few au (i.e., at debris disk radii) given the limits set by the non-
detection of sub-mm dust.

####  7.1.4 Photospheres: The Remaining Class III Stars

The remaining 24 class III’s appear to be bare photospheres as they do not
exhibit infrared excesses or CO gas. For 20/24 of these class IIIs there is no
evidence of ongoing accretion on to their stars, (the remaining 4 of which we
were unable to confirm, see $\S$3.2). In all 24 cases however, the SED
analysis concluded that no significant excess emission is present and so we
believe that these are appropriately classed, but lack evidence for
planetesimal formation, either in the mid-IR or in the sub-mm.

####  7.1.5 The Effect of Stellar Binarity?

Zurlo, submitted (2020) observed 29/30 of the same class III stars in our
survey for binarity (J160714.0-385238 was too faint to observe), within a
separation range of 20-1200 au, and found that our sample contains 6 binaries.
Three of these have wide (>100 au) separations: Sz 108B at ${\sim}$683 au, Sz
116 at ${\sim}$245 au, and CD-35 10498 at ${\sim}$172 au, and 3 of these have
intermediate (20-100 au) separations: J154306.3-392020 at ${\sim}$24 au, CD-39
10292 at ${\sim}$47 au and V1097 Sco at ${\sim}$21 au. Stellar binaries are
known to impact the detection rate of disks, and thus may have influenced
those observed in our class III survey, however for circumstellar emission
with the extent probed by our observations, only those with intermediate
separations are likely to significantly influence circumstellar disk
properties (Jensen et al., 1996; Alexander, 2012; Rosotti & Clarke, 2018;
Manara et al., 2019; Long et al., 2019). Indeed, Yelverton et al. (2019)
showed that the detection rate of older debris disks was lower for binary
separations between 25-135 au, thus if the circumstellar emission measured in
our class III survey is dust from the progenitors of the debris disks observed
at 10s of au, then a similar relationship may exist. While no circumstellar
emission was detected towards the 3 class III stars with intermediate
separation binaries, a detection rate of 0/3 is consistent with being the same
as that of the class III stars without confirmed intermediate separation
binaries (6/26). Thus there is no evidence in our sample that stellar binarity
influences the detection fraction of circumstellar emission, although we
cannot rule out that the binary companions of J154306.3-392020, CD-39 10292
and V1097 Sco accelerated the dispersal of their disks, and we cannot comment
on the influence of binaries with separations outside of the 20-1200 au range
probed by Zurlo, submitted (2020). There is also no evidence that intermediate
separation binaries are over-represented in the class III sample compared with
that of class II’s, as might be expected if such binaries disperse their
protoplanetary disks more quickly; e.g., Zurlo, submitted (2020) found
intermediate separation binaries towards 4/76 of their Lupus class II sample,
which is consistent with the incidence rate of 3/29 towards class III’s.
Larger sample sizes would help to ascertain the influence of binarity on disk
detection rates.

Figure 11: Disk masses in Lupus as a function of stellar mass for the class II
and class III surveys. The dashed lines show the mean upper limit mass from
the class II data in blue (Ansdell et al., 2016), purple (Sanchis et al.,
2020) and for the class III sources (this work) in orange. The blue solid line
shows the relationship for class II YSOs between their stellar mass and dust
mass from Ansdell et al. (2016).

###  7.2 The Lupus Disk Mass Distribution

This survey of 30 class III stars was sensitive to disk masses down to
$0.024M_{\oplus}$, which for context we include on a plot of disk mass vs age
in Fig. 10, alongside the survey populations of Upper Sco (Barenfeld et al.,
2016), Lupus (class II) (Ansdell et al., 2016) and Ophiuchus (Cieza et al.,
2019; Williams et al., 2019), and for the debris disk population in Holland et
al. (2017). We note here that in comparing the masses of disks in different
populations it is most important that the value of $\kappa_{\nu}$ is used
consistently. Our value of 3.5 cm2 g-1 (which does not account for the gas
mass, since we are comparing their dust masses) is consistent with that used
to estimate the dust masses for the Lupus class II, Upper Sco and Ophiuchus
disks. However, the $\kappa_{\nu}$ used to estimate the debris disk masses in
Holland et al. (2017) and dust masses in Sanchis et al. (2020) were both lower
by ${\sim}60\%$, and thus we re-scaled these accordingly. For illustrative
purposes we randomly assigned the ages of stars in Upper Sco, Lupus and
Ophiuchus between 5-11 Myr, 1-3 Myr and 0.75-1.0 Myr respectively (consistent
with the ages of stars in these regions). This plot demonstrates that whilst
the masses of protoplanetary disks are significantly higher than those of the
much older debris disks, the dust masses measured for the Lupus class IIIs are
consistent with those of the much older debris disk population. It also shows
that there has been a significant improvement in sensitivity compared with
previous surveys which has allowed us to detect the lowest mass disks at 1-3
Myr. Whilst it may appear from the plot that an increase in dust mass is seen
from those of the class III’s we detected at ${\sim}2$ Myr and the brightest
debris disks at ${\sim}$10-100 Myr, similar to that discussed in Turrini et
al. (2019), this may be due to a spectral type bias; i.e., the class III
sample is dominated by K and M stars, whereas the debris disk sample is
dominated by G, F and A type stars which we would expect to host more massive
disks. The same caveat applies when comparing the detections on this plot
between different regions; e.g., the large number of $0.1-0.3M_{\oplus}$ disk
detections in Upper Sco could result from the different stellar mass, disk
mass and disk radius distribution compared with Lupus (Ansdell et al., 2016;
Hendler et al., 2020), but also from other features of the sample selection as
noted later in this section.

Here we combine these class III results with the previous surveys of class II
stars (as outlined in $\S$2), which were sensitive to disk masses of
${\sim}0.2M_{\oplus}$, to consider the disk mass distribution in Lupus as a
whole. The previous class II surveys had showed that there is a dependence of
disk mass on stellar mass, and Fig. 11 shows disk mass vs stellar mass for the
combined sample (excluding class II sources without stellar masses, see
Ansdell et al., 2016, 2018; Sanchis et al., 2020, for the complete list). Note
that for the class II star V1094 Sco, we re-derived the disk mass from the
flux in van Terwisga et al. (2018) using equation 1 of Ansdell et al. (2016)
and the Gaia DR2 distance to this source. This recovers the aforementioned
correlation for the class II population, with higher mass disks found only
around higher mass stars, and lower mass disks only around lower mass stars
(for which there was also a greater proportion of non-detections). Fitting to
the detected disks in the class II population, Ansdell et al. (2016) found a
positive relationship between $\log M_{\rm dust}$ and $\log M_{\star}$ with a
slope $1.8\pm 0.4$ (which we re-plot on our Fig. 11).

The detected class III disks all have a similar dust mass, and fall below the
detection threshold of the class II survey. Whilst the class II and III
detections cover a range of stellar masses, when non-detections are included
it may be clearer that neither this survey nor the Ansdell et al. (2016)
survey were particularly sensitive to disks around the lowest mass stars
${<}0.30M_{\odot}$. The Ansdell et al. (2016) class II survey deliberately
excluded Brown Dwarf YSOs (i.e., those with masses ${<}0.1M_{\odot}$), but we
include those analysed by Sanchis et al. (2020) which observed sources down to
YSO stellar masses of ${\sim}0.02M_{\odot}$ (of which 4/9 with
${<}0.1M_{\odot}$ are found with upper limit disk masses). For the sub-sample
of 31 stars in Lupus that are ${>}0.7M_{\odot}$ (including both class II and
III stars), there were no non-detections of class II stars, and so this sub-
sample can be considered to have been searched down to a uniform depth of
$0.024M_{\oplus}$. Considering such stars it is clear from Fig. 11 that there
is a large gap in disk mass between the lowest mass class II disk at
${\sim}2M_{\oplus}$ and the highest mass class III disks at
${\sim}0.09M_{\oplus}$. The same is true when considering the 54 stars with
masses $M_{\star}>0.3M_{\odot}$, but the gap appears to get smaller because of
the lower disk masses of the class II’s. For even lower stellar masses, there
is a large fraction of class II’s with non-detections, and the correlation
would suggest that many of these simply have disks just below the detection
threshold of that survey.

Figure 12: Cumulative fraction of disks with a mass greater than those in the
binned range for the combined Lupus class II and III surveys, for all stars
with either $M_{\star}>0.3M_{\odot}$ (blue dashed) or those with
$M_{\star}>0.7M_{\odot}$ (orange dotted). The shaded regions represent
$1\sigma$ confidence intervals.

Rather than comparing the disk mass distribution of class II’s with class
III’s, we first consider in Fig. 12 the disk mass distribution for the
combined sample of class II’s and III’s, but cutting this in stellar mass.
This was calculated by comparing the number of stars with disks detected at a
given mass with the number of stars for which such detections were possible
(which is usually the full sample until the survey sensitivity is reached,
though this is variable due to the different distances to the stars), with
asymmetric error bars given by Poisson statistics (Gehrels, 1986, see also
$\S$7.3). This further quantifies the gap mentioned above. Also evident from
Fig. 12 is that a larger fraction of higher mass stars have lost their
protoplanetary disks at this age (50% have $M_{\rm dust}>0.3M_{\oplus}$ for
$M_{\star}>0.7M_{\odot}$) compared with lower mass stars (60% for
$M_{\star}>0.3M_{\odot}$).

The presence of a gap in disk masses (see Figs. 11 and 12) suggests that the
process responsible for dispersing proto-planetary disks acts quickly, perhaps
being triggered once the disk mass has decreased below a certain threshold.
This is reminiscent of what has been called "two-timescale behaviour" in the
observational literature about disk dispersal (see e.g., the seminal papers of
Strom et al., 1989; Skrutskie et al., 1990), on the basis of infrared
excesses. Because infrared excesses trace only the warm inner regions, our
results imply that the disk dispersal process removes also the bulk of the
disk’s dust mass that is concentrated at 10s of au on roughly the same
timescale, depleting this by at least an order of magnitude. Theoretically,
the best explanation for the two-timescale behaviour is the process of photo-
evaporation (see e.g., Clarke et al., 2001; Alexander et al., 2014; Ercolano &
Pascucci, 2017), i.e., a disk wind driven by the gas thermal energy. Photo-
evaporation disperses the inner part (<1 au) of the disk when the mass
accretion rate onto the star becomes lower than the mass-loss rate of the
wind. Because there is currently a large uncertainty about photo-evaporation
mass-loss rates depending on which energy bands of the radiation from the
central star provide the energy to drive mass-loss, it is interesting to
estimate the mass-loss rates implied by our observations. A lower proto-
planetary disk mass threshold of ${\sim}3M_{\rm{Jup}}$ (such as that found in
Lupus for solar mass stars, assuming a dust-to-gas ratio of 100), for a
typical viscous time-scale of 1 Myr (Rosotti et al., 2017; Lodato et al.,
2017) would imply a photo-evaporative mass-loss rate of ${\sim}3\times
10^{-9}M_{\odot}$ yr-1. This value is in the range of what can be explained by
current models (Gorti & Hollenbach, 2009; Owen et al., 2011; Picogna et al.,
2019) including far ultra-violet and X-ray photons, while being too high to be
explained by earlier models based only on extreme ultra-violet radiation (see
e.g., Alexander et al., 2006).

As an alternative, over the last few years it has been proposed that disk
evolution and dispersal is controlled by disk winds driven by the disk
magnetic field, rather than a thermal wind as in photo-evaporation. Because
these models are still in their infancy, depending on how the magnetic field
flux evolves over the disk lifetime they can predict a wide variety of
outcomes (see e.g., Armitage et al., 2013; Suzuki et al., 2010; Bai, 2016),
ranging from a slow "fading" of the disk to a behaviour similar to that
predicted by photo-evaporation. Our observational results constrain that only
the latter alternative is compatible, at least for ${>}0.3M_{\odot}$ stars
with disks that disperse within ${\sim}2$ Myr.

An important caveat is that our observations are not sensitive to gas mass,
rather to dust mass which is dominated by mm-cm-sized grains. Photo-
evaporation, or other types of disk winds, remove the gas, and while these can
also entrain grains up to ${\sim}1\,\mu$m in size (Takeuchi et al., 2005) they
cannot be directly responsible for the rapid removal of larger dust grains.
While Alexander & Armitage (2007) showed that the interplay with dust dynamics
can remove the dust in the inner part of the disk (<1 au), dust at larger
distances still remains in the disc. The problem of where the mm-sized dust
ends up was highlighted by Wyatt et al. (2015) where it was argued that this
must be removed in collisions, either grinding it into small dust or growing
it into larger bodies. This feasibility of the former possibility was
demonstrated by Owen & Kollmeier (2019), who showed that collisional
destruction can grind the grains into micron-sized dust that is removed by
radiation pressure. The latter possibility was also demonstrated by Carrera et
al. (2017) which showed that the dust can grow into planetesimals, either
through collisions or via the streaming instability. What our observations
show is that whichever of these the dispersal mechanisms is acting, its effect
on the mm-sized dust must be both rapid and complete, in removing over 90% of
the dust mass. Although depletion may not be so substantial for lower mass
stars (i.e., since the gap is not so pronounced), the transition disk of Sz 84
has a disk >100 times more massive than the disk of the class III star, Sz
107, despite these having approximately the same stellar mass (i.e., 0.18 and
0.16$M_{\odot}$ respectively), and mid-IR excess (see Fig. 2). Deeper
observations of the lower mass class II sources (i.e., to at least the same
sensitivity as this work) would therefore allow us to quantify if this gap
persists for the lowest mass YSOs between the class II and class III stages.

###  7.3 Comparison with Class IIIs in other Star Forming Regions

The disk mass distribution is shown in Fig. 13 in the format shown commonly in
other surveys of star forming regions, wherein the samples are split into
their different classes. The complete sample of Lupus YSOs is included,
excluding only sources without confirmed stellar masses. This was produced by
binning all such sources in their disk masses, and evaluating the uncertainty
on the fraction in each bin using binomial statistics (i.e., for the total
number of sources detected, and the total number of sources in each bin, see
Gehrels, 1986). This reinforces that the class III’s in Lupus are orders of
magnitude lower in mass than their class II counterparts. This class III
distribution will be compared in $\S$7.4 with that expected for debris disks
at this age, which shows excellent agreement.

First we use Fig. 13 to compare the Lupus class III disk mass distribution
with that of the only other two surveys that observed class III stars in
Ophiuchus and Upper Sco (Williams et al., 2019; Barenfeld et al., 2016,
respectively). However, it is important to note that there are two reasons why
these three class III samples may not be equivalent. First, Ansdell et al.
(2016) demonstrated that there is a different stellar mass distribution for
Lupus compared with Upper Sco YSOs from Barenfeld et al. (2016), which may be
important given the trend of class II disk mass with stellar mass in Fig. 11,
and a similar trend for debris disk masses (Greaves & Wyatt, 2003). Second,
and perhaps most importantly, the class III sample selection criteria are
different for the different surveys. While all have the same requirement for
class III stars not having large 12 and 24 $\mu$m excesses (i.e., they are all
samples of class III stars that lie in the green shaded region of Fig. 2), the
Ophiuchus and Upper Sco surveys set additional constraints to exclude class
III’s in the bottom left corner for which there was no evidence for mid-IR
excess emission. The ODISEA survey in Ophiuchus excluded class III stars with
[K]-[24]<1 (equivalent to excluding those in our survey in the cluster of near
photospheric fluxes in the bottom corner of Fig. 2; see the horizontal green
line on that figure) and the survey of Upper Sco excluded all class III stars
with bare photospheres (see the discussion in both Luhman & Mamajek, 2012;
Barenfeld et al., 2016). This complicates the comparison between these
different associations, and to acknowledge the different sample selection
criteria we add the suffix X1 and X2 to class III in reference to the
Ophiuchus and Upper Sco samples in Fig. 13.

Figure 13: Cumulative fraction of disks with a mass greater than those in the
binned range for Lupus sources (blue), where class III are dashed and class II
are dotted, and similarly for Ophiuchus (in orange) and class III sources in
Upper Sco (in red), noting that for both Ophiuchus and Upper Sco the class III
sources were selected based on different constraints to our Lupus survey, and
are hence labelled as Class III X1 or Class III X2 here. The shaded regions
represent $1\sigma$ confidence intervals.

We first compare Lupus with Ophiuchus. The class II’s in Lupus and Ophiuchus
have comparable mass distributions, although the masses are slightly lower for
Ophiuchus. This is different to expectations of disks decreasing in mass with
age, given that Ophiuchus is slightly younger than Lupus at ${\sim}1$ Myr.
However, this may simply reflect a stellar mass distribution with more lower
mass stars given the correlation in Fig. 11, or it may be due to differences
in the local environment (see further discussion in Williams et al., 2019).
The possibility of overlapping populations (and therefore ages) in the YSOs of
Ophiuchus (Luhman & Esplin, 2020), and presumably Lupus, should also be
considered. Nevertheless, it could also be the case that this mass
distribution reflects an increase in dust production (e.g., from planetesimal
collisions) at the slightly older age of Lupus in comparison to Ophiuchus
(e.g., as indicated by the toy model of Gerbig et al., 2019). The Ophiuchus
class III X1 survey was not as deep as that performed here, instead reaching a
sensitivity comparable with the class II surveys of ${\sim}0.2M_{\oplus}$.
Thus it appears surprising that the rate of detection of class III’s is
comparable for the two surveys, which translates into the disk mass
distribution for Ophiuchus class III’s looking like it is higher than those of
Lupus by around half to one order of magnitude (see Fig. 13). However, as
noted above, sources in the Ophiuchus class III X1 sample all have large
[K]-[24] excesses, which we have shown implies a larger disk mass, given that
mid-IR excess is likely present. A deep survey of all Ophiuchus class III’s
would be required for a proper comparison. We can however, consider what a
Lupus class III X1 survey would have found, noting that 4 of our class III’s
meet the class III X1 K-24>1 colour threshold (excluding those that are
confused or which have mid-IR upper limits). Two of these had sub-mm
detections, both with disk masses almost an order of magnitude lower than
those in Williams et al. (2019). While the highest class III disk masses could
be rarer in Lupus than in Ophiuchus, this is more likely a result of small
number statistics. For example, we should not be surprised by the fact that we
did not detect disks more massive than ${\sim}0.1M_{\oplus}$ in our class III
sample, since the model in § 7.4 predicts these should occur at a rate $\sim
3$% (i.e., about 1 in 30). Such high disk masses would be expected to occur
more frequently in a class III X1 sample.

In the case of Upper Sco, we are also not comparing with a consistently
selected sample nor observing depth. In the Barenfeld et al. (2016) survey all
sources were measured to a sensitivity of ${\sim}0.1M_{\oplus}$, and only
sources named "Debris/Ev. Trans" overlap with the same constraints we set for
our class III sub-sample. As shown by Luhman & Mamajek (2012), there is a
large population of low [K]-[12] vs [K]-[24] sources (i.e., bare photospheres)
which were not observed, whereas in our class III sample such bare
photospheres were included (i.e., the detection fractions cannot be compared
as directly as implied by Fig. 13). Nevertheless, we found in our sample 6
stars not deemed to be bare photospheres (i.e., 6 stars which would have met
the Luhman & Mamajek, 2012, criteria and fall in what we are calling a class
III X2 sample). We detected sub-mm circumstellar dust in 4 of these, and all
of these were less massive than the 5/31 Upper Sco detections. Similar to the
case in Ophiuchus, the lack of high mass class III disks could be down to
small number statistics; e.g., if the rates of $>0.1M_{\oplus}$ class III X2
disks is the same in the two regions then we would expect to have detected
$\sim 1$ (i.e., $6\times 5/31$) in Lupus. However, it remains possible that
the two disk mass distributions are different either due to evolution, since
Upper Sco is older (i.e., 5-11 Myr, see Mathews et al., 2012), or because of
the different stellar mass distributions (i.e., for which in comparison to our
Lupus sample, there are a higher fraction of M-types and a lower fraction of
K-types in Upper Sco, along with a number of G-type stars which did not
feature in our sample). To consider these possibilities, comparisons between
samples must be done in a consistent manner, and to fully consider the
complete disk mass distribution of class III stars, sources without excesses
(including bare photospheres) should be included.

###  7.4 Population Model

The population of debris disks seen around older main sequence stars has been
well characterised by far-IR surveys of both nearby field stars of co-eval
populations in moving groups. These show a decay in disk brightness with
stellar age that can be understood within the context of a population model in
which all stars are born with a planetesimal belt that depletes due to steady
state collisional erosion on a different timescale according to its radius and
initial mass (Wyatt et al., 2007). We use such a population model,
specifically that with parameters that fit observations of nearby Sun-like
stars from Sibthorpe et al. (2018), to predict what this debris disk
population would have looked like at the age of Lupus, for which we used the
stellar masses (see Table 2), temperatures and luminosities (derived for the
sources as modelled by our SEDs in $\S$5.4), and $\kappa_{\nu}{\sim}3.5$ cm2
g-1 to be representative of the Lupus class III sample and estimated
consistently with our previous mass calculations. Fig. 14 shows the model
predictions for the dust masses and black body radii (for a distribution
ranging between 1-1000 au) of a population of 10,000 debris disks after just
${\sim}2\rm{Myr}$ of evolution, and in Fig. 13 the cumulative disk mass
distribution is shown alongside the measured rate from this survey.

There are reasons why extrapolating this population model back to ${\sim}2$
Myr could result in erroneous predictions. For example, the model is
constrained by the 20% of nearby stars for which their debris disks can be
detected. This means that the disks of the remaining 80% of the population,
other than being too faint to detect at old ages, are unconstrained. The model
has accounted for these by assuming they have small radii $\ll 5$ au, since
this means they have collisionally depleted below the detection threshold by
Gyr ages, whereas this might not be the case for Lupus ages. Indeed, Fig. 14
shows that the model population includes a large fraction of small bright
disks. However, where the population model is expected to make robust
predictions is for bright disks with radii $>5$ au, since it is those that
evolve into the disks that remain detectable after Gyrs of evolution.

This population model was used to predict the detection rate for this survey,
and finds at a mass sensitivity of ${\sim}0.024M_{\oplus}$ (i.e., one
consistent with our survey) a rate of ${\sim}7\%$. Thus it was anticipated
that, if class III stars host debris disks that are representative of the
population seen around older stars, we should achieve detections of debris
disks for $N_{\rm{det}}=2\pm\sqrt{2}$ of the sample, consistent with the
observed detection rate of 4/30. The debris disk mass distribution predicted
in Fig. 13 is also close to that observed for the class III sample.
Furthermore, the model predicted that the population of non-detectable disks
would have a mean disk flux of $0.0136$ mJy. Using Eq.3 at the $141$ pc
average class III source distance, we find this flux to correspond to a mass
of $0.0026M_{\oplus}$, consistent with (but slightly lower than) the measured
stacked non-detection mass of $0.0048\pm 0.0017M_{\oplus}$. This demonstrates
that the model has both accurately predicted the detection rate of Lupus class
III disks, and the mean disk mass of those below our detection threshold.

Figure 14: Top: Predicted population dust mass and radii for a sample of
10,000 disks based on the models of Wyatt et al. (2007); Sibthorpe et al.
(2018), and the expected detection thresholds at $24\mu\rm{m}$, $70\mu\rm{m}$
(i.e., with Spitzer), $100\mu\rm{m}$, $160\mu\rm{m}$ (i.e., with Herschel) and
$850\mu\rm{m}$ (i.e., with ALMA). Shown are the positions of the 4 sub-mm
detections (lilac) and the 2 hot IR-excesses (pink). Note that the plotted
radii are approximated by the black body SED fits (based on temperature and
fractional luminosity which have a large range of values see Table 4) but that
these span a wide range of values consistent with their detections/non-
detections with either Spitzer (at 24$\mu$m) or ALMA (at 850$\mu$m).

While for two of class III detections there remains some uncertainty about
whether primordial gas or dust is present (i.e., NO Lup and Sz 107), the two
remaining detections (i.e., of MU Lup and J155526-333823) are compatible with
the predicted debris disk population. The good agreement of the model
predictions with the disks observed around class III stars has a number of
implications. First, it would seem that whatever processes have caused the 30
class III stars at this age to have lost their protoplanetary disks faster
than the sample that are class II, those processes have not adversely affected
the formation of planetesimals at 10s of au. This suggests that planetesimal
formation, or more specifically the formation of the parent bodies of the
debris disks seen around older stars, occurs earlier than 2 Myr, also in
agreement with studies of young (class II) Brown Dwarfs in Ophiuchus and Lupus
(e.g., Testi et al., 2016; Sanchis et al., 2020, respectively).

Having inferred that many (if not all) the detected class III disks are
examples of young debris disks, given these are already bright this further
suggests that debris disks are already stirred at this age (i.e., collisions
are grinding planetesimals into dust, rather than collisions leading to
growth). We confirm this for the 4 sub-mm detected disks by calculating the
collisional lifetime of mm-cm sized grains (i.e., the timescale over which
bodies of this size collide in a collisional cascade, see Wyatt, 2008,
equation 15), which for dr/r=0.5, e=0.1, $\rm{Q_{D}^{\star}}=10^{3}$ J kg-1,
the stellar masses in Table 2 and the disk radii in Table 3 (which for MU Lup,
NO Lup and Sz 107 are upper limits), we find collision lifetimes in the ranges
${\sim}1-10$ kyr, ${\sim}3-30$ kyr, ${\sim}3-30$ kyr and ${\sim}20-200$ kyr
for MU Lup, J155526-333823, NO Lup and Sz 107 respectively (i.e., all much
less than the lower-bound age of Lupus of 2 Myr, although only ${\sim}10\%$ in
the case of cm-sized grain collisional lifetimes for the disk around Sz 107).
Accounting for this range of ages, this means that even for recent
protoplanetary disk dispersal (i.e., within 100-200 kyr), the collisional
lifetime of mm-cm grains in at least 3 of these systems is much shorter than
the system age (i.e., the dust must be replenished via the break-up of larger
bodies).

Further, we can estimate the lower limit on the maximum planetesimal size
required for this dust replenishment to occur based on these collisional ages,
and find a range of sizes between $10-300$ m (i.e., 10s to 100s of metre-sized
bodies are required to replenish the mm and cm sized grains observed in the
sub-mm). This is in good agreement with the debris disk population model, and
measurements of aged debris disks, which both require the presence of ${\sim}$
km-sized planetesimals. Although there are models in which debris disks are
stirred by the growth of large Pluto-sized bodies within them (so-called self-
stirring models, Kenyon & Bromley, 2010), at tens of au on timescales of 1-3
Myr these would be unable to stir debris disks (i.e., in order to produce the
observed dust). For example, in the resolved 80 au belt of J155526.2-333823,
this would take >1 Gyr, far in excess of the age of Lupus. Self-stirring of
1-3 Myr old disks at tens of au by Pluto-sized bodies is therefore
incompatible with our observations.If these belts are not self-stirred, an
alternative scenario may be one in which the stirring is forced by secular
interactions between the planetesimals and unseen planets. Equation 15 of
Mustill & Wyatt (2009) requires any such perturbing planets to be close to the
belt and massive (i.e., given the timescale of 2 Myr). For example, in the
case of J155526.2-333823, a planet 10 au inside the belt with a relatively low
eccentricity of e=0.1 would need to have more mass than 2MJ in order to stir
the disk within the the age of Lupus. On the other hand, disks may simply be
born stirred (as discussed by Booth & Clarke, 2016), or stirred during the
process of protoplanetary disk dispersal.

##  8 Conclusions

We report the ALMA observations of 30 class III stars in the $1-3$ Myr old
Lupus Star Forming Region, having measured these to the deepest sensitivity
thus far, to probe the disk masses present in each system. We detected
significant continuum emission from 4 class III sources, 2 nearby class II
sources (within the observation field of view), and marginally detect emission
from the 26 class III sources without individual detections after stacking. We
also report the first CO gas detection around a class III star, NO Lup.

We attribute this continuum emission to circumstellar dust present in these
Young Stellar Objects, and determine the range of dust masses from these 4
detections as $0.036-0.093M_{\oplus}$. Through stacking our sub-sample of non-
detections, we tentatively report dust mass present at the level
${\sim}0.0048M_{\oplus}$. By an order of magnitude, these disks are less
massive than those previously detected in other star forming regions. We found
one of the class III detections to be spatially resolved, for which we
constrained the size of the disk to be $80\pm 15\rm{au}$, and for the three
unresolved class III detections, we set upper limits on their radii (all being
$<60\rm{au}$).

The detection of 4 class III disks is in accordance with steady-state debris
disk model predictions, which suggest that the planetesimal belts observed
around older systems (i.e., the ${\sim}10$s of au debris disks) are already in
place by ${\sim}2\rm{Myr}$. By comparing the frequency of detections as a
function of disk mass, for the first time, we probe the disk mass distribution
of YSOs down to ${\sim}0.024M_{\oplus}$ levels, and show that the evolution of
dust mass between the class II and III stages is both rapid and complete,
removing over 99% of dust. We show how this changes when considering different
cuts in stellar mass, demonstrating that this may not be as substantial for
lower mass stars. By also comparing this disk mass distribution with the
debris disk evolution model, we show very good agreement between the predicted
detection rate across the complete range of masses in our class III sample.

Considering the complete sample of class II and III Lupus YSOs outlined here
(comprising 30/126 class III and 96/126 class II), we show that the processes
which remove the protoplanetary disks around these stars have not affected the
formation of planetesimals at tens of au (i.e., the planetesimals that form
the parent bodies of debris disks). This further implies that such young disks
are already stirred, such that they are bright enough to be detected.

## Acknowledgements

We thank the anonymous reviewer for their comments which improved the quality
of this work. We also thank Alice Zurlo for useful discussions in relation to
stellar multiplicity in Lupus. JBL is supported by an STFC postgraduate
studentship. GMK is supported by the Royal Society as a Royal Society
University Research Fellow. GR acknowledges support from the Netherlands
Organisation for Scientific Research (NWO, program number 016.Veni.192.233).
M.T. has been supported by the UK Science and Technology research Council
(STFC), and by the European Union’s Horizon 2020 research and innovation
programme under the Marie Sklodowska-Curie grant agreement No. 823823 (RISE
DUSTBUSTERS project). MK gratefully acknowledges funding by the University of
Tartu ASTRA project 2014-2020.4.01.16-0029 KOMEET. This project has received
funding from the European Union’s Horizon 2020 research and innovation
programme under the Marie Sklodowska-Curie grant agreement No 823823
(DUSTBUSTERS). This work was partly supported by the Deutsche Forschungs-
Gemeinschaft (DFG, German Research Foundation) - Ref no. FOR 2634/1 TE
1024/1-1. This work was partly supported by the Italian Ministero
dell’Istruzione, Università e Ricerca (MIUR) through the grantProgettiPremiali
2012 iALMA (CUP C52I13000140001) and by the DFG cluster of excellence Origin
and Structure of the Universe (www.universe-cluster.de).

## Data Availability

This work makes use of the following ALMA data: $ADS/$ $JAO.ALMA$
$\\#2018.1.00437.S$. ALMA is a partnership of ESO (representing its member
states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST and
ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic
of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ.

## References

  * ALMA Partnership et al. (2015) ALMA Partnership et al., 2015, ApJ, 808, L3
  * Adams et al. (1987) Adams F. C., Lada C. J., Shu F. H., 1987, ApJ, 312, 788
  * Alcalá et al. (2014) Alcalá J. M., et al., 2014, Astronomy and Astrophysics, 561, A2
  * Alcalá et al. (2019) Alcalá J. M., Manara C. F., France K., Schneider C. P., Arulanantham N., Miotello A., Günther H. M., Brown A., 2019, A&A, 629, A108
  * Alexander (2012) Alexander R., 2012, ApJ, 757, L29
  * Alexander & Armitage (2007) Alexander R. D., Armitage P. J., 2007, MNRAS, 375, 500
  * Alexander et al. (2006) Alexander R. D., Clarke C. J., Pringle J. E., 2006, MNRAS, 369, 216
  * Alexander et al. (2014) Alexander R., Pascucci I., Andrews S., Armitage P., Cieza L., 2014, in Beuther H., Klessen R. S., Dullemond C. P., Henning T., eds, Protostars and Planets VI. p. 475 (arXiv:1311.1819), doi:10.2458/azu_uapress_9780816531240-ch021
  * Andrews & Williams (2005) Andrews S. M., Williams J. P., 2005, ApJ, 631, 1134
  * Andrews et al. (2013) Andrews S. M., Rosenfeld K. A., Kraus A. L., Wilner D. J., 2013, ApJ, 771, 129
  * Andrews et al. (2018) Andrews S. M., et al., 2018, ApJ, 869, L41
  * Ansdell et al. (2016) Ansdell M., et al., 2016, ApJ, 828, 46
  * Ansdell et al. (2018) Ansdell M., et al., 2018, ApJ, 859, 21
  * Armitage et al. (2013) Armitage P. J., Simon J. B., Martin R. G., 2013, ApJ, 778, L14
  * Bai (2016) Bai X.-N., 2016, ApJ, 821, 80
  * Baraffe et al. (2015) Baraffe I., Homeier D., Allard F., Chabrier G., 2015, A&A, 577, A42
  * Barenfeld et al. (2016) Barenfeld S. A., Carpenter J. M., Ricci L., Isella A., 2016, ApJ, 827, 142
  * Beckwith et al. (1990) Beckwith S. V. W., Sargent A. I., Chini R. S., Guesten R., 1990, AJ, 99, 924
  * Birnstiel et al. (2010) Birnstiel T., Dullemond C. P., Brauer F., 2010, A&A, 513, A79
  * Bitsch et al. (2015) Bitsch B., Lambrechts M., Johansen A., 2015, A&A, 582, A112
  * Booth & Clarke (2016) Booth R. A., Clarke C. J., 2016, MNRAS, 458, 2676
  * Booth et al. (2018) Booth R. A., Meru F., Lee M. H., Clarke C. J., 2018, MNRAS, 475, 167
  * Carrera et al. (2017) Carrera D., Gorti U., Johansen A., Davies M. B., 2017, ApJ, 839, 16
  * Carrera et al. (2020) Carrera D., Simon J. B., Li R., Kretke K. A., Klahr H., 2020, arXiv e-prints, p. arXiv:2008.01727
  * Chen et al. (2009) Chen C. H., Sheehan P., Watson D. M., Manoj P., Najita J. R., 2009, ApJ, 701, 1367
  * Cieza et al. (2007) Cieza L., et al., 2007, ApJ, 667, 308
  * Cieza et al. (2013) Cieza L. A., et al., 2013, ApJ, 762, 100
  * Cieza et al. (2019) Cieza L. A., et al., 2019, MNRAS, 482, 698
  * Clarke et al. (2001) Clarke C. J., Gendrin A., Sotomayor M., 2001, MNRAS, 328, 485
  * Comerón (2008) Comerón F., 2008, The Lupus Clouds. p. 295
  * Cutri & et al. (2013) Cutri R. M., et al. 2013, VizieR Online Data Catalog, p. II/328
  * Dent et al. (2014) Dent W. R. F., et al., 2014, Science, 343, 1490
  * Dipierro et al. (2015) Dipierro G., Price D., Laibe G., Hirsh K., Cerioli A., Lodato G., 2015, MNRAS, 453, L73
  * Dunham et al. (2015) Dunham M. M., et al., 2015, ApJS, 220, 11
  * Ercolano & Pascucci (2017) Ercolano B., Pascucci I., 2017, Royal Society Open Science, 4, 170114
  * Espaillat et al. (2014) Espaillat C., et al., 2014, in Beuther H., Klessen R. S., Dullemond C. P., Henning T., eds, Protostars and Planets VI. p. 497 (arXiv:1402.7103), doi:10.2458/azu_uapress_9780816531240-ch022
  * Evans et al. (2009) Evans N. J., et al., 2009, VizieR Online Data Catalog, p. J/ApJS/181/321
  * Foreman-Mackey et al. (2013) Foreman-Mackey D., Hogg D. W., Lang D., Goodman J., 2013, Publications of the Astronomical Society of the Pacific, 125, 306
  * Frasca et al. (2017) Frasca A., Biazzo K., Alcala J. M., Manara C. F., Stelzer B., Covino E., Antoniucci S., 2017, VizieR Online Data Catalog, pp J/A+A/602/A33
  * Gaia Collaboration et al. (2018) Gaia Collaboration et al., 2018, Astronomy and Astrophysics, 616, A1
  * Galli et al. (2013) Galli P. A. B., Bertout C., Teixeira R., Ducourant C., 2013, A&A, 558, A77
  * Galli et al. (2015) Galli P. A. B., Bertout C., Teixeira R., Ducourant C., 2015, A&A, 580, A26
  * Gehrels (1986) Gehrels N., 1986, ApJ, 303, 336
  * Genda et al. (2015) Genda H., Kobayashi H., Kokubo E., 2015, ApJ, 810, 136
  * Gerbig et al. (2019) Gerbig K., Lenz C. T., Klahr H., 2019, A&A, 629, A116
  * Gorti & Hollenbach (2009) Gorti U., Hollenbach D., 2009, ApJ, 690, 1539
  * Greaves & Rice (2011) Greaves J. S., Rice W. K. M., 2011, MNRAS, 412, L88
  * Greaves & Wyatt (2003) Greaves J. S., Wyatt M. C., 2003, MNRAS, 345, 1212
  * Haisch et al. (2001) Haisch Karl E. J., Lada E. A., Lada C. J., 2001, ApJ, 553, L153
  * Hardy et al. (2015) Hardy A., et al., 2015, A&A, 583, A66
  * Hendler et al. (2020) Hendler N., Pascucci I., Pinilla P., Tazzari M., Carpenter J., Malhotra R., Testi L., 2020, ApJ, 895, 126
  * Hildebrand (1983) Hildebrand R. H., 1983, QJRAS, 24, 267
  * Holland et al. (2017) Holland W. S., et al., 2017, MNRAS, 470, 3606
  * Hughes et al. (1994) Hughes J., Hartigan P., Krautter J., Kelemen J., 1994, AJ, 108, 1071
  * Hughes et al. (2018) Hughes A. M., Duchêne G., Matthews B. C., 2018, ARA&A, 56, 541
  * Jensen et al. (1996) Jensen E. L. N., Mathieu R. D., Fuller G. A., 1996, ApJ, 458, 312
  * Johansen & Youdin (2007) Johansen A., Youdin A., 2007, ApJ, 662, 627
  * Kennedy & Wyatt (2014) Kennedy G. M., Wyatt M. C., 2014, Monthly Notices of the Royal Astronomical Society, 444, 3164
  * Kenyon & Bromley (2004) Kenyon S. J., Bromley B. C., 2004, ApJ, 602, L133
  * Kenyon & Bromley (2010) Kenyon S. J., Bromley B. C., 2010, ApJS, 188, 242
  * Köhler et al. (2000) Köhler R., Kunkel M., Leinert C., Zinnecker H., 2000, Astronomy and Astrophysics, 356, 541
  * Kral et al. (2017) Kral Q., Matrà L., Wyatt M. C., Kennedy G. M., 2017, Monthly Notices of the Astronomical Society, 469, 521
  * Krautter et al. (1997) Krautter J., Wichmann R., Schmitt J. H. M. M., Alcala J. M., Neuhauser R., Terranegra L., 1997, Astronomy and Astrophysics Supplement Series, 123, 329
  * Lada (1987) Lada C. J., 1987, in Peimbert M., Jugaku J., eds, IAU Symposium Vol. 115, Star Forming Regions. p. 1
  * Lada & Wilking (1984) Lada C. J., Wilking B. A., 1984, Astrophysical Journal, 287, 610
  * Lebouteiller et al. (2011) Lebouteiller V., Barry D. J., Spoon H. W. W., Bernard-Salas J., Sloan G. C., Houck J. R., Weedman D. W., 2011, ApJS, 196, 8
  * Lodato et al. (2017) Lodato G., Scardoni C. E., Manara C. F., Testi L., 2017, MNRAS, 472, 4700
  * Long et al. (2019) Long F., et al., 2019, ApJ, 882, 49
  * López Martí et al. (2005) López Martí B., Eislöffel J., Mundt R., 2005, A&A, 440, 139
  * Lovell, submitted. (2020) Lovell, submitted. J. B., 2020, Lupus Class III Disks with ALMA, Submitted