# Molecular Bubble and Outflow in S Mon Revealed by Multiband Datasets

Dejian Liu  Purple Mountain Observatory, Chinese Academy of Sciences, Nanjing
210023, People’s Republic of China  School of Astronomy and Space Science,
University of Science and Technology of China, Hefei 230026, People’s Republic
of China  Ye Xu  Purple Mountain Observatory, Chinese Academy of Sciences,
Nanjing 210023, People’s Republic of China  School of Astronomy and Space
Science, University of Science and Technology of China, Hefei 230026, People’s
Republic of China  YingJie Li  Purple Mountain Observatory, Chinese Academy of
Sciences, Nanjing 210023, People’s Republic of China  Zehao Lin  Purple
Mountain Observatory, Chinese Academy of Sciences, Nanjing 210023, People’s
Republic of China  School of Astronomy and Space Science, University of
Science and Technology of China, Hefei 230026, People’s Republic of China
Chaojie Hao  Purple Mountain Observatory, Chinese Academy of Sciences, Nanjing
210023, People’s Republic of China  School of Astronomy and Space Science,
University of Science and Technology of China, Hefei 230026, People’s Republic
of China  WenJin Yang  School of Astronomy & Space Science, Nanjing
University, 163 Xianlin Avenue, Nanjing 210023, People’s Republic of China
Max-Plank-Institut f$\ddot{u}$r Radioastronomie, auf dem H$\ddot{u}$gel 69,
53121 Bonn, Germany  Jingjing Li  Purple Mountain Observatory, Chinese Academy
of Sciences, Nanjing 210023, People’s Republic of China  School of Astronomy
and Space Science, University of Science and Technology of China, Hefei
230026, People’s Republic of China  Xinrong Liu  Purple Mountain Observatory,
Chinese Academy of Sciences, Nanjing 210023, People’s Republic of China
School of Astronomy and Space Science, University of Science and Technology of
China, Hefei 230026, People’s Republic of China  Yiwei Dong  Purple Mountain
Observatory, Chinese Academy of Sciences, Nanjing 210023, People’s Republic of
China  School of Astronomy and Space Science, University of Science and
Technology of China, Hefei 230026, People’s Republic of China  Shuaibo Bian
Purple Mountain Observatory, Chinese Academy of Sciences, Nanjing 210023,
People’s Republic of China  School of Astronomy and Space Science, University
of Science and Technology of China, Hefei 230026, People’s Republic of China
Deyun Kong  School of Physics, Harbin Institute of Technology, Harbin 150001,
People’s republic of China

###### Abstract

We identify a molecular bubble, and study the star formation and its feedback
in the S Mon region, using multiple molecular lines, young stellar objects
(YSOs), and infrared data. We revisit the distance to S Mon, $\sim$ 722 $\pm$
9 pc, using _Gaia_ Data Release 3 parallaxes of the associated Class II YSOs.
The bubble may be mainly driven by a massive binary system (namely 15 Mon),
the primary of which is an O7V-type star. An outflow is detected in the shell
of the bubble, suggesting ongoing star formation activities in the vicinity of
the bubble. The total wind energy of the massive binary star is three orders
of magnitude higher than the sum of the observed turbulent energy in the
molecular gas and the kinetic energy of the bubble, indicating that stellar
winds help to maintain the turbulence in the S Mon region and drive the
bubble. We conclude that the stellar winds of massive stars have an impact on
their surrounding environment.

Interstellar molecules (849) – Young stellar objects (1834) – Stellar wind
bubbles (1635)

##  1 Introduction

The mass-loss phase is a very common phenomenon during the early evolutionary
stages of a star. Strong winds from young stars inject momentum and energy
into the surrounding environment, thus affecting the dynamics and structure of
their parent clouds (e.g., Lada, 1985; Arce et al., 2011; Frank et al., 2014;
Bally, 2016; Li et al., 2020). High-mass stars evolve very quickly and they
may still be embedded in their parent cloud when they reach the main-sequence
stage. When a high-mass star forms in the molecular cloud, its ultraviolet
(UV) radiation can ionize and heat the surrounding gas to create an HII
region. The expanding HII region can reshape the surrounding molecular gas and
drive a bubble (e.g., Churchwell et al., 2006; Zhang et al., 2015). During the
main-sequence stage, high-mass stars can drive spherical winds that blow away
the gas around the star and then form bubbles (e.g., Castor et al., 1975; Arce
et al., 2011). Bubbles and outflows are common stellar-feedback phenomena in
the star-formation process and provide information about the physical
properties of their surroundings (e.g., Lada, 1985; Arce et al., 2010, 2011;
Li et al., 2018, 2020).

There have been quite a few studies concerned with bubbles, outflows and their
feedback (e.g., see Frank et al., 2014; Dale, 2015; Bally, 2016, for reviews).
Since directly measuring the proper motion of extended structures, like
molecular clouds, is very difficult, many previous studies only analyzed the
line-of-sight dynamics in a region. Additionally, recent studies started to
investigate the 3D structure of molecular clouds in the solar neighborhood
(e.g., Green et al., 2019; Lallement et al., 2019; Guo et al., 2021), mainly
enabled by the _Gaia_ space mission (Gaia Collaboration et al., 2016). Young
stellar objects (YSOs) that only recently formed in star-forming molecular
clouds are still close to their birth sites; hence, they still, on average,
share the distance and kinematic properties of their parent clouds (e.g.,
Fűrész et al., 2008; Tobin et al., 2009; Hacar et al., 2016; Großschedl et
al., 2021). On this basis, cloud distances can be estimated by using the
distances of YSOs located within the clouds as proxies (e.g., Großschedl et
al., 2018, 2021; Zhang, 2023). The recent release of astrometric data from
_Gaia_ Data Release 3 (DR3, Gaia Collaboration et al., 2023) allows one to
obtain a large number of high-quality parallaxes and proper motions of nearby
YSOs, which also makes it possible to explore further the 3D kinematics of
molecular clouds that harbor these YSOs (e.g., Großschedl et al., 2021;
Flaccomio et al., 2023).

NGC 2264 is a nearby high-mass star-forming complex in Monoceros with a
distance of about 700–800 pc (e.g., Zucker et al., 2020; Flaccomio et al.,
2023), and it is well studied at multiple wavelengths (for a review see Dahm,
2008). Studies of the spatial and dynamical structure of NGC 2264 revealed
that it is composed of two complexes: S Mon in the northern part and the Cone
Nebula in the southern part (e.g., Sung et al., 2008; Tobin et al., 2015;
Venuti et al., 2018). The gas structure of S Mon is dominated by massive stars
at its center. For instance, 15 Mon, which consists of an O7V primary star
with a mass of 35 $M_{\odot}$ and an O9.5V secondary star with a mass of 24
$M_{\odot}$ (Dahm, 2008). 15 Mon exhibits variations in its UV resonance line
profile as well as fluctuations in its soft X-ray flux, which can be
interpreted as being induced by variations in its mass-loss rate (Snow et al.,
1981; Grady et al., 1984; Dahm, 2008).

Sung & Bessell (2010) and Venuti et al. (2018) studied the star formation
history of S Mon by analyzing the properties of young clusters, and found that
stars are still forming actively in the region. Buckle et al. (2012) found
that the CO emission surrounding S Mon shows many filaments and arcs, and
these structures are unlikely due to protostellar outflow activity. Based on
13CO emission, Tobin et al. (2015) speculated that there exists a bubble in
the S Mon region which may be driven by 15 Mon. Since studies of the S Mon
region either used only young stars or molecular gas, or did not focus on the
details of the bubble, the properties of this region, especially the bubble,
outflows, and their feedback remain unclear.

In this work, we carried out molecular line observations and collected YSO
information in the S Mon region, aiming to provide insights into the star
formation activities and their feedback. This paper is organized as follows.
In Section 2, we present the data used in this work. The correlation between
YSOs and molecular gas, the distance of the molecular gas, and details of the
bubble and outflow are presented in Section 3. We discuss the driving source
of the bubble and the energy cascade in S Mon, as well as the kinematic
features of the bubble based on the YSOs in Section 4. In Section 5, we
summarize the main results.

##  2 Data

###  2.1 Molecular Gas

Molecular line observations of S Mon were carried out between 2022 April and
November using the Purple Mountain Observatory Delingha (PMODLH) 14 m
millimeter-wavelength telescope. Four molecular lines were observed: 12CO (J =
1 - 0), 13CO (J = 1 - 0), C18O (J = 1 - 0), and HCO+ (J = 1 - 0). The three CO
lines were observed simultaneously. All lines were observed with the nine-beam
Superconducting Spectroscopic Array Receiver system in the sideband separation
mode (Shan et al., 2012), with the 13CO and C18O lines in the lower sideband
and the 12CO and HCO+ lines in the upper sideband. The on-the-fly (OTF; Sun et
al., 2018) mode was used during the observations, and the OTF raw data were
gridded in fits cube with a pixel size of $30^{\prime\prime}$ using the GILDAS
software package (Pety, 2005).111http://www.iram.fr/IRAMFR/GILDAS The typical
integration time per position was approximately 3 minutes. Each fast Fourier
transform spectrometer with a bandwidth of 1 GHz provides 16 384 channels,
yielding a spectral resolution of 61
kHz.222http://www.radioast.nsdc.cn/english/zhuangtaibaogao.php All results
presented in this work are expressed as brightness temperatures, $T_{\rm
R}^{*}=T_{A}^{*}/\eta_{\rm mb}$, where $T_{A}^{*}$ is the antenna temperature
and $\eta_{\rm mb}$ is the main beam efficiency. Table 1 lists the
observational parameters of the molecular lines.

Table 1: Observational Parameters of the Molecular Lines

Molecular line | Rest frequency | HPBW | $T_{\rm sys}$ | $\eta_{\rm mb}$ | $\delta_{v}$ | rms noise
---|---|---|---|---|---|---
(J = 1 - 0) | (GHz) | (′′) | (K) |  | (km s-1) | (K)
(1) | (2) | (3) | (4) | (5) | (6) | (7)
12CO | 115.271 | 50 | $250-300$ | 51.1% | 0.159 | 0.15
13CO | 110.201 | 52 | $150-200$ | 56.3% | 0.166 | 0.07
C18O | 109.782 | 52 | $150-200$ | 56.3% | 0.166 | 0.06
HCO+ | 89.189 | 59 | $100-150$ | 58.1% | 0.187 | 0.03

Note. — (1) Molecular line. (2) Frequency of the molecular line. (3) Half-
power beamwidth (HPBW). (4) Typical system temperature. (5) Main beam
efficiency. (6) Velocity resolution. (7) Main beam rms noise.

###  2.2 Young Stellar Object Sample

The classification of YSOs was initially proposed by Lada & Wilking (1984) to
delineate their evolutionary stages, and is defined as follows. Class I YSOs
are protostars that are still embedded in an envelope and surrounded by a
circumstellar disk, causing significant infrared-excess and a rising slope of
the spectral energy distribution (SED) in the mid- to far-infrared range;
Class II YSOs are pre-main-sequence (PMS) stars with circumstellar disks,
which show a negative SED slope, that is significantly above the slopes of
main-sequence stars (see Lada et al., 2006). Class II YSOs with inner holes in
their accretion disks are named transition disks and Class III YSOs are PMS
stars surrounded by thin disk remnants (anaemic disks, ADs) or they could be
already disk free.

To study the young cluster NGC 2264, a variety of multiwavelength datasets and
criteria have been employed to identify the evolutionary stages of YSOs. For
example, Sung et al. (2008, 2009) classified YSOs based on SED slopes of
infrared data combined with X-ray data from the Chandra X-ray Observatory
(Weisskopf et al., 2002). Broos et al. (2013) adopted a Naive Bayes approach
to classify YSOs, using data from the Massive Young star-forming Complex Study
in Infrared and X-rays (MYStIX) project (Feigelson et al., 2013). Cody et al.
(2014) classified YSOs based on the slope of SEDs at near- and mid-infrared
wavelengths. Rapson et al. (2014) utilized a color-based classification scheme
to classify YSOs with Spitzer (Fazio et al., 2004) data and Two Micron All Sky
Survey (2MASS; Skrutskie et al., 2006) photometry. Venuti et al. (2018)
classified YSOs using spectroscopic and photometric data from the _Gaia_ -ESO
Survey (GES, Gilmore et al., 2012; Randich et al., 2013) and the Coordinated
Synoptic Investigation of NGC 2264 Survey (CSI 2264, Cody et al., 2014).

We describe in detail the YSO data collection, classification criteria, and
astrometric parameters acquisition in Appendix A, and briefly introduce these
procedures here. We collected YSO candidates in the S Mon region from Sung et
al. (2009), Broos et al. (2013), Cody et al. (2014) and Rapson et al. (2014),
within an area of $202.8^{\circ}<l<203.1^{\circ}$, and
$2.0^{\circ}<b<2.3^{\circ}$. Duplicates were eliminated through coordinate
cross-matching with a radius of 1$\arcsec$, yielding in an initial sample of 1
023 YSO candidates. As these candidates were classified by different methods,
we reclassified them through a homogeneous classification scheme based on
infrared data. We adopted the SED slope ($\alpha={\rm d}\log\lambda
F_{\lambda}/{\rm d}\log\lambda$, Lada, 1987) and the classification scheme of
Greene et al. (1994) to classify the evolutionary stage of the YSOs. Flat-
spectrum sources are classified as Class II YSOs in this study. In line with
Sung et al. (2009), we use $\alpha\geq 0.3$ to classify objects as Class I and
$-1.8$ to 0.3 for Class II. When using selection criteria based on infrared
photometry, only the part of Class III YSOs with infrared-excess (i.e., ADs)
can be classified. Therefore, X-ray data collected by Flaccomio et al. (2023)
were also employed to classify additional Class III YSOs. The Class III YSOs
in this work are comprised of sources with ADs and disk-less PMS stars. For
the former, the lower limit for a source with infrared-excess is given by Lada
et al. (2006) with $\alpha>$ -2.56. This slope is typical for stellar
photospheres when there is no additional infrared-excess from a disk present,
and when using the Spitzer IRAC bands to calculate the infrared SED (see e.g.,
Lada et al., 2006; Teixeira et al., 2012). Given the uncertainties in
photometry, however, we took -2.3 as a rough limit to exclude possible
contaminants from field stars. For disk-less PMSs, they have a SED slope of
below -2.3 and X-ray emission. The SED fitting was performed without reddening
corrections following Cody et al. (2014). In this case, we gain a sample of
seven Class I, 160 Class II, and 245 Class III YSOs (101 ADs + 144 disk-less
PMSs). The parallaxes and proper motions of these YSOs were taken from _Gaia_
DR3 (Gaia Collaboration et al., 2023). The radial velocities (RVs) of the YSOs
were collected from the Sloan Digital Sky Survey (SDSS) DR17 APOGEE-2
(Abdurro’uf et al., 2022) and the catalog of Tobin et al. (2015).

To summarize, we end up with a sample of one Class I, 145 Class II, and 231
Class III YSOs that had matched with the Gaia DR3 catalog, and 47 Class II and
91 Class III YSOs with available RV information.

To compare the gas and stellar RVs, the heliocentric RVs ($V_{\rm Helio}$) of
the YSOs were converted to the local standard of rest velocity ($V_{\rm
LSR}$). We followed the conversion method described by Reid et al. (2009) via:

$V_{\rm LSR}=V_{\rm Helio}+(V_{\odot}\sin l+U_{\odot}\cos l)\cos
b+W_{\odot}\sin b,$ (1)

where $l$ and $b$ are the Galactic longitude and latitude, respectively, and
$U_{\odot},V_{\odot}$, and $W_{\odot}$ denote the solar motion components in
the direction of the Galactic Center, along Galactic rotation, and toward the
north Galactic pole, respectively. The used values for $U_{\odot},V_{\odot},$
and $W_{\odot}$ are 10.3, 15.3, and 7.7 km s-1, respectively (Kerr & Lynden-
Bell, 1986).

##  3 Results

In the following we present the relationship between the molecular gas and
YSOs in Section 3.1. Then, the distance to the S Mon region is estimated in
Section 3.2 using YSOs that are closely correlated to the gas. In Section 3.3,
we present the bubble observed in the S Mon region. An outflow is detected in
one of clumps, as described in Section 3.4.

###  3.1 Relationship between the Molecular Gas and YSOs

Figure 1 shows a position–velocity (P-V) diagram (a) and the normalized
$V_{\rm LSR}$ distribution (b) for YSOs and molecular gas traced by 12CO in S
Mon. The 12CO emission shows two velocity components in this region, with one
ranging from 0 to 7 km s-1 and the other from 7 to 15 km s-1. The presence of
a gap or minimum in the CO emission at $V_{\rm LSR}\sim 7$ km s-1 is in
agreement with the results of Tobin et al. (2015). In panel (b) of Figure 1,
the $V_{\rm LSR}$ of the majority (62%) of the Class II YSOs peak at about 10
km s-1, which is consistent with the 12CO emission component ranging from 7 to
15 km s-1. Most (71%) Class III YSOs have LSR velocities between 0 and 7 km
s-1, with a peak at 5 km s-1.

We investigate the two velocity components of molecular gas in S Mon
separately. Figure 2 shows integrated maps of 12CO and HCO+ within the
velocity range from 0 to 7 km s-1 (panels (a) and (b)) and from 7 to 15 km s-1
(panels (c) and (d)). The integrated map of molecular gas is superposed on a
Herschel (Pilbratt et al., 2010) ${\rm H_{2}}$ column density map of this
region, with a resolution of 18$\arcsec$ (Nony et al., 2021). The 12CO and
HCO+ emission between 0 and 7 km s-1 appears to be diffuse and is not
consistent with the Herschel ${\rm H_{2}}$ column density map. This component
seems to dominate in the easternmost area of the region. In particular in the
12CO map (Figure 2 (a)) there appears to be a somewhat stronger component
starting from the left of the image, beyond the shown sky cutout. Our
investigation does not include the easternmost region. The 12CO and HCO+
emission between 7 and 15 km s-1 shows an arc-like structure, possibly
indicative of a “bubble” (centered at (0, 0) offset, see more details in
Section 3.3). We note that the distribution of the molecular gas within this
velocity range is in agreement with the distribution of the Herschel H2 column
density map.

Figure 3 depicts the projected distribution of the YSOs in combination with
12CO emission and Herschel H2 column density map. Both the Class II and Class
III YSOs, which have LSR velocities between 7 and 15 km s-1 (yellow pluses),
generally coincide with the denser regions of the molecular cloud. More
members associated with the molecular gas are found belong to the Class II
YSOs.

Compared with the Class III YSOs, the Class II YSOs are better correlated with
the molecular gas in both spatial distribution and kinematic properties.
Considering the association between the molecular gas and the Class II YSOs,
we are able to use the Class II YSOs to characterize the properties of the
molecular gas, such as its distance and kinematics. The proper motions of the
YSOs are used to derive the kinematics of the molecular gas in the plane of
the sky, since they are the missing part, and cannot be measured directly for
clouds.

Figure 1: (a) P-V diagram showing 12CO (background), Class II (red pluses),
and Class III (magenta crosses) YSOs in S Mon. (b) Normalized histogram of
$V_{\rm LSR}$ for 12CO (gray line), Class II (magenta line), and Class III
(red line) YSOs.

(a) 12CO, 0–7 km s-1 (b) HCO+, 0–7 km s-1 (c) 12CO, 7–15 km s-1 (d) HCO+, 7–15
km s-1

Figure 2: Integrated maps of 12CO and HCO+ over the velocity ranges of 0 to 7
km s-1 and 7 to 15 km s-1. The contour levels are shown from 30% to 90% with
steps of 10% of the peak integrated intensity, and the black open circle in
the lower left corner of each panel represents the beam size of PMODLH. The
gray background shows the Herschel H2 column density map (Nony et al., 2021).
In all panels, the (0, 0) offsets corresponds to
$(l,b)=(202.95^{\circ},2.15^{\circ})$.

(a) Class II (b) Class III

Figure 3: Projected distribution of Class II (left) and Class III (right) YSOs
with different velocity ranges. The blue pluses (crosses) indicate Class II
(III) YSOs with $V_{\rm LSR}$ ranging from 0 to 7 km s-1. The yellow pluses
(crosses) indicate Class II (III) YSOs with $V_{\rm LSR}$ ranging from 7 to 15
km s-1. The magenta pluses (crosses) represent Class II (III) YSOs with
unknown $V_{\rm LSR}$. The red diamonds represent Class I YSOs. The gray
contours and background are the same as in Figure 2.

###  3.2 Distance to S Mon

Accurate distances to molecular clouds are essential to study their physical
properties. We employ the _Gaia_ DR3 parallaxes of the Class II YSOs to
estimate the distance to S Mon. In this work, only Class II YSOs with parallax
accuracies better than 10% are selected for distance estimation. First,
outliers of parallaxes beyond $3\sigma$ from the median value are excluded.
Class II YSOs, which do not coincide with the gas distribution (CO at 7–15 km
s-1, magenta pluses in Figure 4 (a)), are also excluded from the distance
estimation. Finally, a total 35 YSOs that meet these criteria are selected to
estimate the distances shown (yellow pluses in Figure 4 (a)).

The distance is determined from _Gaia_ parallaxes with $d({\rm
pc})=1\,000/\varpi$, where $\varpi$ is the parallax in units of mas. Lindegren
et al. (2021) provided a method to correct the parallax zero points, however
for red sources with $m_{G_{\rm BP}}-m_{G_{\rm RP}}>1.6$, the correction could
lead to a large parallax bias of about 40 $\mu$as. Since the photometric data
of the used Class II YSOs have generally redder colors, we do not apply the
parallax zero-point correction in this work.

We estimate the distance to the cloud directly using the YSOs associated with
the cloud and their Gaia DR3 parallaxes. The cumulative distribution function
(CDF) of _Gaia_ distances is shown in panel (b) of Figure 4, where the median
distance estimate is of 722 $\pm$ 9 pc. The uncertainty is derived by
$\sigma/\sqrt{N}$, where $\sigma$ represents the dispersion of the distance,
and $N$ represents the number of stars.

Flaccomio et al. (2023) conducted a study on the NGC 2264 star-forming region
using X-ray data obtained with XMM-Newton telescope. The S Mon core region in
their study appears to be close to the region we studied, although we note
that they did not give the specific coordinates of the region they
investigated. Although only half of the samples we used for distance
estimation aligned with their sample counterparts, our distance estimation is
similar to their results (719 $\pm$ 9 pc). However, the distance we derive is
smaller than the distance estimated ($759_{-26}^{+10}\pm 37$ pc) by Zucker et
al. (2020), who inferred the distance by fitting a line-of-sight dust model
with Gaia DR2 data. We conducted a parallax difference test between Gaia DR2
and DR3 using a selected set of Class II YSOs in our region. These YSOs have
LSR velocities ranging from 0 to 15 km s-1, and both DR2 and DR3 exhibit
parallax accuracies better than 10%. The difference between these two groups
is approximately 25 pc, meaning the two estimates agree within their given
uncertainties. Although we did not test the method of Zucker et al. (2020)
explicitly, variations in methods could also introduce some differences in
distance estimations. Therefore, the difference in distance between this work
and that of Zucker et al. (2020) could be caused by the difference in both
method and data. In the subsequent analysis, we adopt 722 $\pm$ 9 pc as the
distance to the S Mon cloud, which has a velocity ranging between 7 and 15 km
s-1.

(a) Distribution of Class II YSOs (b) CDF of _Gaia_ distance

Figure 4: (a) Projected distribution of Class II YSOs with parallax
uncertainties smaller than 10%. The yellow pluses are the Class II YSOs used
for the distance estimation. The blue pluses show YSOs with LSR velocities
smaller than 7 km s-1. The magenta pluses show YSOs which do not coincide with
the gas distribution. The red triangle represents the projected position of 15
Mon. (b) CDF of _Gaia_ distances. The red vertical line indicates the median
value of 722 pc, with a standard error of 9 pc denoted by the gray shadow. The
upper and lower percentiles within 1$\sigma$ around the median are labeled by
the dashed gray lines.

###  3.3 Bubble

We identified the bubble in S Mon based on four criteria summarized from
previous studies (e.g., Arce et al., 2011; Li et al., 2015; Feddersen et al.,
2018): (1) the shell shows a circular or arc-like structure in integrated CO
or HCO+ maps; (2) the P-V diagram of the shell shows an expansion signature
(i.e., a “U” or “V” shape); (3) the CO or HCO+ shell matches the morphology of
infrared data in at least one band; and (4) the shell contains a candidate
driving source. An expanding bubble model requires three parameters: bubble
radius ($R_{\rm bub}$), expansion velocity ($V_{\rm exp}$), and central
velocity ($V_{\rm cnt}$). We followed the method of Arce et al. (2011) to
estimate these parameters. The radius is determined by a Gaussian fit to the
azimuthally averaged profile of the integrated intensity map. The central and
expansion velocities can be roughly read from a P-V diagram, where the upper
part of the “U” or “V” structure corresponds to the central velocity and the
difference between the central velocity and the lowest point of “U” or “V”
structure indicates the expansion velocity. Then, the central and expansion
velocities are corrected by visual fitting based on the channel maps. As
indicated by Arce et al. (2011), the uncertainty in the central velocity of a
bubble is about $\pm 0.5$ km s-1, and since the bubble is not detected over
its entire velocity range, the expansion velocity is obtained as a lower
limit.

An arc-like structure is clearly visible in the channel maps of 12CO and HCO+
(see the details in Figures 5–6), spanning a broad $V_{\rm LSR}$ range from
7.5 to 12.5 km s-1. This structure is also evident in the integrated maps of
12CO and HCO+ from 7.5 to 12.5 km s-1 and the Herschel H2 column density map
(see panels (a) and (d) of Figure 7, Nony et al., 2021). The P-V diagrams of
12CO and HCO+ (panels (b) and (e) of Figure 7) both show an upside down “U” or
“V” shape. Based on the P-V diagram of 12CO (HCO+), $V_{\rm exp}$ and $V_{\rm
cnt}$ are derived as $\sim 3.1$ km s-1($\sim 3.4$ km s-1) and $\sim 8.2$ km
s-1($\sim 8.4$ km s-1), respectively. The $V_{\rm exp}$ values estimated from
12CO and HCO+ are similar, thus an average of 3.3 km s-1 is adopted as a lower
limit for the expansion velocity. The center of the bubble is set to be (0,
0). The radius and thickness of the bubble traced by the 12CO emission (see
panel (c) of Figure 7) are $\sim 5.4$ and $\sim 4.2$ arcmin, respectively,
corresponding to 1.13 and 0.88 pc at a distance of 722 pc. These values are
compared with the bubble traced by HCO+ emission (see panel (f) of Figure 7),
i.e., a radius of $\sim 5.3$ arcmin ($\sim 1.11$ pc) and a thickness of $\sim
4.4$ arcmin ($\sim 0.92$ pc). Based on both the 12CO and HCO+ emission and
considering the uncertainties from the distance and the Gaussian fits, we
derive the radius and thickness of the bubble to be $1.1\pm 0.1$ pc and
$0.9\pm 0.1$ pc, respectively. The mean column density of the bubble is $\sim
2.5\times 10^{20}~{}{\rm cm}^{-2}$ and the mass ($M_{\rm shell}$) is $\sim
1\,400\pm 400~{}M_{\odot}$ (see the calculation and estimation of
uncertainties in Appendix B).

Another arc-like structure appears to be located in the southwest region of
the Herschel H2 column density map (see the boxed region of Figure 8 (a)). The
low-sensitivity CO data, observed using PMODLH, were employed to investigate
the structure. The rms noise of the 12CO, 13CO and C18O spectra is 0.37, 0.20,
and 0.20 K, respectively. A small arc-like structure is only evident in the
13CO integrated map (see panel (b) of Figure 8), so we use the 13CO data to
produce a P-V diagram for the arc-like structure. However, there is no clear
“U” or “V” shape in the P-V diagram (see panel (c) of Figure 8). In addition,
we do not find a potential driving source for this arc-like structure.
Moreover, using high-spatial-resolution CO maps (20$\arcsec$), Tauber et al.
(1993) found that the properties of the northern region of the arc-like
structure ($\Delta b>-7\arcmin$) are consistent with the ionizing radiation
from 15 Mon, while those of the southern region are not. In conclusion, the
arc-like structure shown is unlikely to be a bubble. Given the wind-blown
appearance, this feature might be caused by feedback from the northeast side,
coincident with the general location of the massive O- and B-type stars in the
region.

There are eight clumps, labeled as A–H in Figure 7, identified from the
integrated map of C18O emission. Figure 9 plots the four molecular lines we
observed toward the C18O peak position of each clump. Table 2 lists the fits
for the C18O peak emission for each clump. There are three main components
among these eight clumps, i.e., $3~{}{\rm km~{}s^{-1}}<V_{\rm LSR}<7~{}{\rm
km~{}s^{-1}}$, $7~{}{\rm km~{}s^{-1}}<V_{\rm LSR}<10~{}{\rm km~{}s^{-1}}$, and
$10~{}{\rm km~{}s^{-1}}<V_{\rm LSR}<13~{}{\rm km~{}s^{-1}}$. As discussed in
Section 3.1, the bubble is not associated with the component of $3~{}{\rm
km~{}s^{-1}}<V_{\rm LSR}<7~{}{\rm km~{}s^{-1}}$, but instead with the other
two components. Clumps A to D correlate with the velocity component of
$7~{}{\rm km~{}s^{-1}}<V_{\rm LSR}<10~{}{\rm km~{}s^{-1}}$ and centered at
$\sim$ 8.5 km s-1, and clumps E to H correlate with the velocity component of
$10~{}{\rm km~{}s^{-1}}<V_{\rm LSR}<13~{}{\rm km~{}s^{-1}}$ and centered at
$\sim$10.5 km s-1. Figure 10 shows integrated maps of C18O with different
velocity components, which are in agreement with the parameters given in Table
2. These two groups of clumps correspond to different lobes of the bubble, so
the different kinematic characteristics may be caused by the expansion of the
bubble. Furthermore, the velocity dispersion around the bubble is larger (see
Figure 11), indicating possible interaction between the bubble and surrounding
gas. Table 3 lists the H2 column density, mass, and turbulent energy of each
clump dervied using the optically thin 13CO line at the peak of the clump (see
the calculations in Appendix B).

Figure 5: Channel maps of 12CO emission. The red triangle represents the star
15 Mon. The blue dashed circle shows the bubble obtained by fitting the gas,
which has a radius and thickness of $\sim 5.4$ and $\sim 4.2$ arcmin,
respectively. The magenta circle indicates the beam size of PMODLH. The
velocity range for each map is labeled in the upper right of the map.

Figure 6: Channel maps of HCO+. The blue dashed circle shows the bubble
obtained by fitting the gas, which has a radius and thickness of $\sim 5.3$
and $\sim 4.4$ arcmin, respectively. The other labels and colors are the same
as in Figure 5.

(a) Intensity map of 12CO (b) P-V diagram of 12CO (c) Radial intensity profile
of 12CO (d) Intensity map of HCO+ (e) P-V diagram of HCO+ (f) Radial intensity
profile of HCO+

Figure 7: (a) Integrated map of 12CO emission(blue contours) superposed on the
Herschel H2 column density map (gray background, Nony et al., 2021), and the
green contour lines are the emission of C18O. The integrated velocity range is
$7.5-12.5$ km s-1 for both 12CO and C18O, and the contour levels are shown
from 30% to 90% with steps of 10% of the peak integrated intensity of the
corresponding emission. The orange dashed circle shows the bubble. Clumps
identified based on the C18O integrated map are labeled A to H. (b) P-V
diagram of 12CO along the arrow in panel (a), where the contour levels are
shown from 10% to 90% with steps of 10% of the peak value. The red dashed
lines show the expansion velocity range from the visual inspection, i.e.,
$8.2-11.3$ km s-1. (c) Azimuthally averaged radial intensity profile of the
12CO bubble (black curve), where a Gaussian fit (red curve) to the intensity
profile was performed to estimate the radius (peak) and thickness (Full Width
at Half Maximum, FWHM). (d) Integrated map of HCO+ (magenta contours). The
other labels are the same as in panel (a). (e) P-V diagram of HCO+ along the
arrow in panel (d), and the expansion velocity range is $8.4-11.8$ km s-1. The
other features are the same as in panel (b). (f) Azimuthally averaged radial
intensity profile of the HCO+ bubble, where the other labels and symbols are
the same as in panel (c).

(a) Herschel H2 column density map (b) Intensity map of 13CO (c) P-V diagram
of 13CO

Figure 8: (a) Herschel H2 column density map, showing a larger area to the
southwest of the investigated region (Nony et al., 2021). The red circle
depicts the bubble we identified, and the black box shows the region housing
the small arc-like structure. (b) Integrated map of 13CO emission (green
contours) in the boxed region in panel (a). The integrated velocity range is
$9.5-11.0$ km s-1, and the contour levels are shown from 30% to 90% with steps
of 10% of the peak integrated intensity of 13CO emission. (c) P-V diagram of
13CO along the arrow in panel (b), where the contour levels are shown from 10%
to 90% with steps of 10% of the peak integrated intensity of the corresponding
emission.

(a) clump A (b) clump B (c) clump C (d) clump D (e) clump E (f) clump F (g)
clump G (h) clump H

Figure 9: Spectra of 12CO, 13CO, C18O, and HCO+ of the C18O peak position of
clumps A to H. The red lines represent the sum of the multiple Gaussian
fitting components.

Figure 10: Integrated map of C18O superposed on the Herschel column density
map. The velocity range for the left panel is 7–10 km s-1, and for the right
panel is 10–13 km s-1. The other features are the same as in Figure 7 (a).
Figure 11: Distribution of the velocity dispersion ($\sigma$). The green
circles indicate the bubble region. Table 2: Observed Parameters of Each Clump

Clump |  12CO($J$ = 1-0) |  |  13CO($J$ = 1-0) |  | C18O($J$ = 1-0) |  | HCO+($J$ = 1-0)
---|---|---|---|---|---|---|---
$T_{\rm R}^{*}$ | FWHM | $V_{\rm LSR}$ |  | $T_{\rm R}^{*}$ | FWHM | $V_{\rm LSR}$ |  | $T_{\rm R}^{*}$ | FWHM | $V_{\rm LSR}$ |  | $T_{\rm R}^{*}$ | FWHM | $V_{\rm LSR}$
(K) | (km s-1) | (km s-1) |  | (K) | (km s-1) | (km s-1) |  | (K) | (km s-1) | (km s-1) |  | (K) | (km s-1) | (km s-1)
(1) | (2) | (3) | (4) |  | (5) | (6) | (7) |  | (8) | (9) | (10) |  | (11) | (12) | (13)
A | 19.1 | 2.5 | 9.1 |  | 6.5 | 1.5 | 9.2 |  | 2.1 | 0.4 | 9.2 |  | 0.4 | 2.6 | 9.2
B | 17.9 | 3.3 | 8.7 |  | 6.0 | 2.3 | 8.6 |  | 2.1 | 0.8 | 8.7 |  | 0.6 | 2.1 | 8.5
C | 19.5 | 3.3 | 8.5 |  | 8.1 | 1.7 | 8.2 |  | 2.7 | 0.6 | 8.3 |  | 0.6 | 1.8 | 8.3
D | 26.1 | 3.0 | 8.8 |  | 6.9 | 1.9 | 9.0 |  | 2.3 | 0.4 | 9.2 |  | 0.7 | 2.9 | 9.0
E | 29.7 | 2.8 | 10.6 |  | 10.8 | 1.7 | 10.6 |  | 2.7 | 0.6 | 10.6 |  | 1.9 | 2.5 | 10.8
F | 25.3 | 2.0 | 10.5 |  | 7.6 | 1.7 | 10.3 |  | 1.9 | 0.6 | 10.2 |  | 2.3 | 1.9 | 10.4
G | 26.5 | 2.6 | 10.8 |  | 7.7 | 2.0 | 10.8 |  | 1.6 | 0.7 | 10.8 |  | 1.7 | 2.4 | 11.1
H | 23.1 | 2.9 | 9.8 |  | 6.7 | 1.9 | 9.9 |  | 1.5 | 0.5 | 10.1 |  | 1.8 | 1.9 | 10.3

Note. — (1) Source name. (2), (5), (8), (11) Brightness temperature(s) of
12CO, 13CO, C18O, and HCO+, respectively, at the emission peak of C18O. (3),
(6), (9), (12) FWHM of 12CO, 13CO, C18O, and HCO+, respectively. (4), (7),
(10), (13) Central velocity of 12CO, 13CO, C18O, and HCO+, respectively.

Table 3: Physical Parameters of the Clumps

Clump | $T_{\rm ex}$ | $d_{\rm clump}$ | $\int T_{\rm mb,13}dv$ | $N_{\rm H_{2},clump}$ | $M_{\rm clump}$ | $E_{\rm turb,clump}$ | $E_{\rm grav,clump}$
---|---|---|---|---|---|---|---
(K) | (pc) | (K km s-1) | ($10^{20}$ cm-2) | ($M_{\odot}$) | (1044 erg) | (1043 erg)
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8)
A | 22.6 | $0.2\pm 0.1$ | 15.4 | 2.7 | $8.1\pm 2.1$ | $1.3\pm 0.3$ | $5.6\pm 0.8$
B | 21.3 | $0.2\pm 0.1$ | 16.5 | 2.4 | $5.8\pm 2.0$ | $2.2\pm 0.8$ | $2.9\pm 0.7$
C | 23.0 | $0.3\pm 0.1$ | 20.5 | 3.7 | $14.3\pm 3.6$ | $3.0\pm 0.8$ | $11.7\pm 2.2$
D | 29.6 | $0.3\pm 0.1$ | 16.2 | 3.4 | $10.4\pm 2.7$ | $2.7\pm 0.7$ | $6.2\pm 1.2$
E | 33.2 | $0.5\pm 0.1$ | 16.8 | 6.1 | $24.8\pm 5.1$ | $5.2\pm 1.1$ | $21.0\pm 4.4$
F | 28.8 | $0.4\pm 0.1$ | 16.1 | 3.7 | $11.5\pm 3.0$ | $2.4\pm 0.6$ | $5.7\pm 1.5$
G | 30.0 | $0.5\pm 0.1$ | 22.6 | 3.9 | $16.7\pm 4.2$ | $4.9\pm 1.2$ | $9.5\pm 3.0$
H | 26.6 | $0.6\pm 0.1$ | 10.5 | 3.1 | $7.7\pm 1.7$ | $2.0\pm 0.4$ | $1.7\pm 0.5$

Note. — (1) Source name. (2) Excitation temperature of the clump, which is
obtained by the 12CO spectral line from the C18O peak emission of each clump.
(3) Diameter of the clump. (4) Integrated intensity of 13CO emission. (5) H2
density of the clump. (6) Mass of the clump. (7) Turbulent energy of the
clump. (8) Gravitational energy of the clump.

###  3.4 The outflow of clump E

In this section we investigate clump E in more detail, since the velocity
profile of this clump alone displays typical outflow characteristics and
because the spectra of the other clumps are too complex to identify reliable
outflows. Outflows are a direct signature of ongoing star formation, and they
can be identified by examining their line profiles, integrated intensity maps,
and P-V diagrams. The central velocity and position of the driving source are
estimated through the C18O line. The emission peaks of C18O can be seen in
Figure 7. The initial velocity ranges of the blue and red wings are determined
from the line wings of the spectral line where the C18O emission reaches the
1$\sigma$ noise level. According to morphologies of the P-V diagram and the
integrated map of the line wings of the 13CO and HCO+ emission, then visually
adjusted the velocity range of the line wings. To improve the signal-to-noise
ratio, we smooth the CO and HCO+ lines to a velocity resolution of 0.25 km
s-1.

Figures 12 and 13 show the identified outflow in clump E based on 13CO and
HCO+ emission, respectively. The spatial distributions of the 13CO and HCO+
outflows are consistent with each other and are aligned along the east–west
direction. The blue and red lobes are symmetrically distributed with respect
to the C18O emission peak, which are similar in size. We note that the C18O
and HCO+ emission at $V_{\rm LSR}\sim 8.6$ km s-1, and the offset of $\sim
1.5$ arcmin in the P-V diagram, are not associated with the blue lobe of the
outflow, thus we did not cut the edge of the blue wing to the 1$\sigma$ noise
level.

The method presented in Liu et al. (2021) was employed to estimate various
properties of the outflows, such as the mass, momentum, and energy. Therein,
they derived the H2 column density of the outflow lobes assuming that the gas
is in local thermodynamic equilibrium (LTE) and the excitation temperature is
30 K (15 K) for 13CO (HCO+). The mass of an outflow lobe can be then estimated
based on its size and the H2 column density. The other physical properties of
the outflows are obtained based on the mass, velocity, and size of the lobes
(see Appendix A of Liu et al., 2021). In our work, we use the ratios of column
densities of $N_{\rm H_{2}}/N_{\rm{}^{13}CO}\approx 5\times 10^{5}$ and
$N_{\rm H_{2}}/N_{\rm HCO^{+}}\approx 10^{8}$ (Turner et al., 1997; Simon et
al., 2001). We also use the inclination of the outflow, i.e., the angle
between the long axis of the outflow and the line of sight, to correct the
parameters of the outflow. Since the inclination of an outflow cannot be
reliably determined, we adopt 57$\fdg$3 as a proxy in this work, which is the
mean value for a random distribution (see details in Bontemps et al., 1996; Li
et al., 2018). We note that the calculation of the outflow’s physical
parameters (see Table 4), involve many assumptions, such as the excitation
temperature and the inclination. Therefore, we only consider the uncertainties
in the outflow lobe region and distance when estimating the errors of the
physical parameters. The mass, momentum, and dynamical timescale of the
outflow are a few $M_{\odot}$, a few tens of $M_{\odot}$ km s-1, and $\sim
10^{5}$ yr, respectively. Based on Yang et al. (2018), the typical mass of a
low-mass outflow is about $0.1\sim 1~{}M_{\odot}$ and the typical mass loss
rate is about $10^{-7}\sim 10^{-6}~{}M_{\odot}~{}{\rm yr^{-1}}$. These values
are much smaller than those found for the outflow in clump E. Moreover, all
physical parameters of the outflow in clump E align with the typical values of
high-mass outflows as summarized in table 1 of Yang et al. (2018). Therefore,
it is unlikely that the outflow in clump E is driven by a single low-mass YSO.
We do not find any high-mass YSOs in this region based on the current data,
yet there are several low-mass YSOs embedded in the vicinity of the outflow
region (see panel (a) of Figures 12 and 13). Therefore, the outflow might be
driven by a deeply embedded protostar that has not been detected yet.

(a) Distribution of the outflow (b) P-V diagram (c) Spectra of the blue side
of the lobe (d) Spectra of the red side of the lobe

Figure 12: Outflow in clump E identified by 13CO emission. (a) Integrated C18O
map with 13CO blue and red lobe contours. The contour levels are shown from
40% to 90% with steps of 10% of the peak intensity of each lobe. The black
diamonds and cyan pluses represent Class I and Class II YSOs, respectively.
The magenta, blue, and red stars are the emission peaks of C18O, and the blue
lobe and red lobe, respectively. (b) P-V diagram along the white arrow in
panel (a). The black contour levels are shown from 10% to 90% with steps of
10% of the peak value. The blue and red dashed lines are the velocity ranges
of the blue and red lobes, respectively. (c) The blue spectrum presents 13CO
at the peak position of the blue emission of the 13CO outflow. The black
spectrum is C18O at the emission peak of C18O. The blue shading of the
spectrum indicates the blue line wing velocity of 13CO. The feature at 8.6 km
s-1 is likely unrelated to the outflow, so emission with a velocity below 8.6
km s-1 is not included in the blue wing of the outflow profile. (d) The red
spectrum is 13CO at the peak position of the red emission of the 13CO outflow.
The black spectrum is C18O at the emission peak of C18O. The red shading of
the spectrum indicates the red line wing velocity of 13CO.

(a) Distribution of the outflow (b) P-V diagram (c) Blue spectra (d) Red
spectra

Figure 13: Outflow in clump E identified by HCO+ emission. The description of
each panel is the same as that in Figure 12, except the contour levels in
panel (a) are shown from 50% to 90% with steps of 10% of the peak intensity of
each lobe.  Table 4: Physical Properties of the Outflow of Clump E

Line | Lobe | $\Delta v$ | $\left\langle\Delta v_{\rm lobe}\right\rangle$ | $N_{\rm H_{2}}$ | $M_{\rm lobe}$ | $l_{\rm lobe}$ | $P_{\rm lobe}$ | $E_{\rm lobe}$ | $t_{\rm lobe}$ | $L_{\rm m~{}(lobe)}$
---|---|---|---|---|---|---|---|---|---|---
(km s-1) | (km s-1) | ($10^{21}$ cm-2) | ($M_{\odot}$) | (pc) | ($M_{\odot}$ km s-1) | ($10^{44}$ erg) | ($10^{5}$ yr) | ($10^{-1}~{}L_{\odot}$)
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8) | (9) | (10) | (11)
13CO | blue | (9, 9.6) | 2.1 | 1.1 | $6.2\pm 1.2$ | $0.8\pm 0.2$ | $24.6\pm 4.7$ | $9.1\pm 1.8$ | $2.2\pm 0.3$ | $1.8\pm 0.2$
13CO | red | (11.5, 12.6) | 2.7 | 1.1 | $9.6\pm 2.8$ | $0.9\pm 0.2$ | $49.8\pm 13.9$ | $24.3\pm 6.4$ | $2.0\pm 0.3$ | $5.3\pm 0.6$
HCO+ | blue | (9, 9.6) | 2.1 | 0.4 | $2.5\pm 0.9$ | $0.8\pm 0.1$ | $10.1\pm 3.4$ | $3.7\pm 1.3$ | $2.3\pm 0.4$ | $0.7\pm 0.1$
HCO+ | red | (12, 13) | 2.7 | 0.5 | $2.9\pm 0.6$ | $0.8\pm 0.1$ | $14.7\pm 3.3$ | $7.1\pm 1.6$ | $1.8\pm 0.2$ | $1.8\pm 0.2$

Note. — (1) Molecular line. (2) Blue/red lobe. (3) Outflow velocity range for
the blue/red lobe. (4) Velocity of the blue/red lobe. (5) H2 density of the
blue/red lobe. (6) Mass of the blue/red lobe. (7) Length of the blue/red lobe.
(8) Momentum of the blue/red lobe. (9) Kinetic energy of the blue/red lobe.
(10) Dynamical timescale of the blue/red lobe. (11) Mechanical luminosity of
the blue/red lobe.

##  4 Discussion

In this section, we discuss the physical properties of the bubble (Section
4.1), and the energy cascade of S Mon region (Section 4.2). In Section 4.3, we
use the Class II YSOs to investigate the kinematic characteristics of the
bubble.

###  4.1 Physical Properties of the Bubble

In the following, we discuss the physical properties of the bubble. Using the
shell mass ($M_{\rm shell}$) and expansion velocity ($V_{\rm exp}$) of the
bubble, we can estimate the shell’s momentum ($P_{\rm shell}$), kinetic energy
($E_{\rm shell}$), and kinetic timescale ($t_{\rm kinetic}$). These quantities
are determined as $P_{\rm shell}=M_{\rm shell}V_{\rm exp}$, $E_{\rm
shell}=0.5M_{\rm shell}V^{2}_{\rm exp}$, and $t_{\rm kinetic}=R_{\rm
bub}/V_{\rm exp}$, respectively. For the bubble we find that $P_{\rm shell}$
is 4 600 $\pm$ 1 300 $M_{\odot}$ km s-1, $E_{\rm shell}$ is $(1.5\pm
0.4)\times 10^{47}$ erg, and $t_{\rm kinetic}$ is $(3.3\pm 0.3)\times 10^{5}$
yr. Considering that the expansion velocity serves as a lower limit, the
uncertainties of these parameters only come from the mass.

We find 13 B-type stars and a massive binary system, 15 Mon, located in the S
Mon region using SIMBAD111https://simbad.u-strasbg.fr/simbad. Most B-type
stars are located in the shell of the bubble rather than at its center (see
Figure 19 in Appendix C), and only two B-type stars are located at the center
of the bubble. We employ Equation (2) from Arce et al. (2011) to estimate the
mass-loss rate ($\dot{m}_{\rm w}$), which can be calculated by:

$\dot{m}_{\rm w}=\frac{P_{\rm shell}}{v_{\rm w}\tau_{\rm w}},$ (2)

where $v_{\rm w}$ is the wind velocity, and $\tau_{\rm w}$ is the wind
timescale (i.e., the duration of wind activity). Here, the kinetic timescale
($t_{\rm kinetic}$) is employed to estimate the mass-loss rate of the driving
source roughly. Given the typical wind velocity of O- or B-type stars is about
1 000–2 000 km s-1 (Chen et al., 2013), the mass-loss rate would be
$10^{-6}-10^{-5}~{}{M_{\odot}}~{}{\rm yr^{-1}}$. The value is much larger than
the typical mass-loss rate of B-type stars
($10^{-8}-10^{-11}~{}M_{\odot}~{}{\rm yr^{-1}}$) but close to that of O-type
stars ($10^{-6}-10^{-7}~{}M_{\odot}~{}{\rm yr^{-1}}$, Chen et al., 2013).
Therefore, it is improbable that the B-type stars are the main driving source
of the bubble. Although 15 Mon is not precisely at the center of the bubble,
it stands out as the closest high-mass star to the bubble’s center, apart from
the two B-type stars. Therefore, 15 Mon might be the main driving source of
the bubble.

The dynamic age of the HII region can be used to judge whether the bubble is
triggered by 15 Mon. Assuming that the HII region expands in a homogeneous
medium, its dynamical age, $t_{\rm dyn}$ is (see the model described by Dyson
& Williams, 1980):

$t_{\rm dyn}=\frac{4R_{\rm s}}{7c_{s}}\left[\left(\frac{R_{\rm H_{II}}}{R_{\rm
s}}\right)^{7/4}-1\right],$ (3)

where $R_{\rm H_{II}}$ is the radius of the HII region (Xu et al., 2017), for
which we adopt an inner radius of 0.7 $\pm$ 0.1 pc. $c_{s}$ is the isothermal
sound speed of the ionized gas, and a value of 10 km s-1 is adopted. The
radius of the associated Strömgren sphere ($R_{\rm s}$) reads:

$R_{\rm s}=\left(\frac{3Q_{\rm Ly}}{4\pi n_{0}^{2}\alpha_{B}}\right)^{1/3},$
(4)

where $\alpha_{B}=2.6\times 10^{-13}~{}{\rm cm^{3}~{}s^{-1}}$ is the hydrogen
recombination coefficient to all levels above the ground level, $n_{0}$ is the
initial number density of the gas, and $Q_{\rm Ly}$ is the ionizing
luminosity. We only consider the number of atoms that can be ionized, so
$n_{0}=\frac{M_{\rm shell}}{4/3\pi R_{\rm bub}^{3}}\approx$ (2.6–6.3) $\times
10^{4}~{}{\rm cm^{-3}}$. Since the actual ionized medium within the bubble is
likely less massive as compared to the whole shell mass determined from the
broad CO ring, the number of atoms that can be ionized might be overestimated.
For an O7V star, the ionizing luminosity ($Q_{\rm Ly}$) is $10^{47.62}~{}{\rm
s^{-1}}$ (Martins et al., 2005). Based on these calculations, the obtained
dynamical age is (3.5–8.9) $\times 10^{5}$ yr. Even if the secondary member of
15 Mon, an O9.5V star, is also taken into account, the dynamical timescale is
only reduced by less than 0.1%. We conclude that the kinetic timescale of the
bubble is comparable to the dynamical age of the HII region.

Montillaud et al. (2019) estimated the age of 15 Mon using two methods. The
first method, based on SED fitting, led to a median age of $1.5-2$ Myr, with a
dispersion between 0.2 and 3 Myr. The second method was based on isochrone
fitting in color–magnitude diagrams and led to a median age of 3 Myr, with a
dispersion between 1 and 6 Myr. Our estimation on the dynamical timescale of
the HII region may be less than the age of 15 Mon. Hence, the bubble is likely
to be driven by the expanding HII region of 15 Mon. As the projected position
of 15 Mon is not located at the center of the bubble, we presume that 15 Mon
may have moved away from its once central position in the bubble.

###  4.2 Energy Cascade in S Mon

Stellar winds inject energy and momentum into the ambient clouds that may help
sustain turbulence and disrupt their surroundings (Nakamura & Li, 2007; Arce
et al., 2011). To assess whether winds can sustain the turbulence in S Mon, we
compare the wind energy injection rate into the cloud, $\dot{E}_{\rm w}$, with
the cloud turbulent dissipation rate, $L_{\rm turb,cloud}$. Based on the 13CO
emission, we can derive the 3D turbulent velocity dispersion ($\sigma_{\rm
3d}$) with an approximation of $\sigma_{\rm 3d}=\sqrt{3}\sigma_{\rm 1d}$ (Li
et al., 2018), where $\sigma_{\rm 1d}$ is the 1D turbulent velocity dispersion
along the line of sight and we consider the mass ($M_{\rm cloud}$) of S Mon.
The obtained $\sigma_{\rm 3d}$ is 1.7 $\pm$ 0.2 km s-1 and $M_{\rm cloud}$ is
2 200 $\pm$ 500 $M_{\odot}$. The cloud turbulent energy, $E_{\rm turb,cloud}$,
is about (6.3 $\pm$ 2.5) $\times 10^{46}$ erg (see the calculation in Appendix
B). The cloud turbulent dissipation rate can be calculated as:

$L_{\rm turb,cloud}=\frac{E_{\rm turb,cloud}}{t_{\rm diss}},$ (5)

where $t_{\rm diss}$ is the cloud turbulent dissipation time. We roughly
estimate $t_{\rm diss}$ via (McKee & Ostriker, 2007):

$t_{\rm diss}\sim 0.5\frac{d_{\rm cloud}}{\sigma_{\rm 1d}},$ (6)

where the cloud diameter is $d_{\rm cloud}\sim$ (4.0 $\pm$ 0.4) pc and
$\sigma_{\rm 1d}\sim$ (1.0 $\pm$ 0.1) km s-1. Based on Equations (5)–(6), we
obtain a cloud turbulent dissipation time of $t_{\rm diss}\sim$(2.0 $\pm$ 0.4)
$\times 10^{6}$ yr, and a cloud turbulent dissipation rate of $\sim$ (0.5–1.7)
$\times 10^{33}~{}{\rm erg~{}s^{-1}}$.

Since the bubble might be driven by the expanding HII region of 15 Mon, the
dynamical age ($t_{\rm dyn}$) of the bubble is more suitable than the
kinemtaic timescale ($t_{\rm kinetic}$) to estimate the mass-loss rate
($\dot{m}_{\rm w}$). Therefore, we reestimated this value, which is
approximately (1.9–8.6)$\times 10^{-6}~{}{M_{\odot}~{}{\rm yr}^{-1}}$, and the
wind energy injection rate can be estimated using the following equation from
McKee (1989):

$\dot{E}_{\rm w}=\frac{1}{2}(\dot{m}_{\rm w}v_{\rm w})\sigma_{\rm 3d}.$ (7)

The derived wind injection rate is about (4.5–8.1)$\times 10^{33}~{}{\rm
erg~{}s^{-1}}$, which is larger than the cloud turbulent dissipation rate.
This indicates that strong winds serve as a significant energy origin for S
Mon, consisting of sustaining the observed turbulence. This result agrees with
the findings of other studies (Arce et al., 2011; Li et al., 2015; Feddersen
et al., 2018).

In the following, we investigate the turbulent energy in this region. The
total wind energy input into the interstellar medium is about $E_{\rm w}=\int
0.5\dot{m}_{\rm w}(t)v_{\rm w}^{2}dt\simeq\Delta Mv_{\rm w}^{2}/2\simeq$
(0.7–1.7) $\times 10^{50}$ erg (Lamers & Cassinelli, 1999). Both the wind
energy of 15 Mon and the kinetic energy of the bubble ($E_{\rm shell}\sim$
(1.5 $\pm$ 0.4) $\times 10^{47}$ erg) are higher than the turbulent energy of
the cloud ($E_{\rm turb,cloud}\sim$ (6.3 $\pm$ 2.5) $\times 10^{46}$ erg),
thus both of them are help to sustain the turbulence in the S Mon region.
Moreover, the wind energy probably plays a crucial role in driving the bubble,
since it is about three orders of magnitude higher than the kinetic energy of
the bubble. Despite the substantial energy carried away by ionizing photons
from the star throughout its lifetime, only 0.01–0.1 percent of this energy is
converted into kinetic energy in the shell, with the majority being lost as
radiation (Geen et al., 2015). Walch et al. (2012) conducted hydrodynamic
simulations to explore the dynamical effects of a single O7 star. They
determined that the total radiative energy injected by the ionizing source
reaches $7\times 10^{51}$ erg within approximately 1 Myr. A portion of this
injected energy is then transformed into kinetic energy, amounting to around
$3\times 10^{48}$ erg, which surpasses the kinetic energy we estimated in this
work ($E_{\rm shell}\sim$ (1.5 $\pm$ 0.4) $\times 10^{47}$ erg). Thus,
radiative energy contributes to driving the bubble.

Next, we studied the turbulent energy of the outflow of clump E. The total
kinetic energy of the outflow in clump E ($\sum{E_{\rm lobe}}\sim 4.4\times
10^{45}$ erg) is higher than its turbulent energy ($E_{\rm turb,clump}\sim$
(5.2 $\pm$ 1.1) $\times 10^{44}$ erg), and it accounts for $\sum{E_{\rm
lobe}}/E_{\rm turb,cloud}\sim 4-10\%$ of the turbulent energy of the cloud.
This suggests that the outflow in clump E helps to maintain the turbulence in
this clump and introduces additional energy to surrounding gas, thereby
affecting the turbulence of the local environment (Arce et al., 2010; Li et
al., 2018). This result is also consistent with the conclusion made by Li et
al. (2020) that the outflow is sufficient to maintain turbulence on a scale of
0.1–0.4 pc.

Moreover, we also estimate the gravitational binding energy of the S Mon
region. The detailed calculation is presented in Appendix B. We find that the
total wind energy of 15 Mon is over two orders of magnitude higher than the
gravitational binding energy of the cloud, i.e., (2.1 $\pm$ 1.0)$\times
10^{47}$ erg, indicating that the wind can disrupt the cloud. Based on the
above calculations, we find that half of the gravitational binding energy of
the cloud is approximately equal to the kinetic energy of the bubble, and much
higher than the total kinetic energy of the outflow in clump E. According to
the virial theorem (Nakamura & Li, 2014), the bubble has the potential to
disrupt the cloud, while the outflow is unable to disrupt the cloud.

We also investigated the gravitational binding energy of clump E ($E_{\rm
grav,clump}\sim$(2.1 $\pm$ 0.2)$\times 10^{44}$ erg), which is lower than the
total kinetic energy of the outflow, indicating that the outflow contributes
to disrupting the clump. The results for the outflow found here are consistent
with the conclusions made by Li et al. (2020) in that the outflow activity
could potentially disperse material away from the parent cloud on a scale of
$0.2-1.0$ pc.

In conclusion, we find the bubble produced by 15 Mon helps to maintain the
turbulence in S Mon and is potentially disrupting the cloud, and while the
outflow from clump E only adds a small fraction of the overall turbulence and
energy of the system, it still may be able to perturb the local gas.

###  4.3 Bubble from the YSO Perspective

Since some YSOs may be too young to have moved away from their birth sites,
they often have the same kinematic properties as their parent clouds (Tobin et
al., 2015; Hacar et al., 2016; Großschedl et al., 2021). Thus, we can utilize
the motions of YSOs associated with a cloud to represent the motions of the
cloud. According to the correlation between YSOs and gas (see Section 3.1 for
details), the Class II YSOs, whose parallax accuracies are better than 10% and
projected locations are approximately aligned with the cloud were initially
selected to investigate the kinematic characteristics of the bubble.

A Gaussian function was used to estimate the average proper motion of the YSO
sample and to remove outliers (beyond a 3$\sigma$ threshold). A total of 38
YSOs were ultimately selected to investigate the kinematic features of the
bubble on the plane of the sky. Figure 14 (a) shows the distribution of proper
motion of the selected YSOs. The mean proper motion of these sources is
$(\bar{\mu_{l}},\bar{\mu_{b}})$ = (2.65, -3.17) ${\rm mas~{}yr^{-1}}$, and the
$1\sigma$ uncertainty range of the Gaussian fitting is (0.14, 0.17) ${\rm
mas~{}yr^{-1}}$.

The intrinsic movement of the Class II YSOs on the plane of the sky was
calculated by subtracting the mean proper motions. Then the tangential
velocities in $l$ and $b$ directions (in units of km s-1) are calculated as
follows:

$\displaystyle\upsilon_{l}=\kappa\frac{1}{\varpi}*(\mu_{l}-\bar{\mu_{l}}),$
(8)
$\displaystyle\upsilon_{b}=\kappa\frac{1}{\varpi}*(\mu_{b}-\bar{\mu_{b}}),$

where $\kappa=4.74~{}{\rm km~{}yr~{}s^{-1}}$, $\varpi$ is the parallax in
units of mas, and $\mu_{l}$ and $\mu_{b}$ are the proper motion along the $l$
and $b$ directions, respectively, in units of mas yr-1.

Figure 14 (b) shows the projected distributions and tangential velocities of
the Class II YSOs. The Class II YSOs are mainly located in the molecular gas
surrounding the bubble, and their motions are relatively random on the plane
of the sky. Figure 15 shows P-V diagrams of the selected YSOs. There is no
clear evidence of expansion in the figure. In order to determine whether the
YSOs are expanding, we regard the expansion velocity as the radial outward
component of the intrinsic velocity from the center of the bubble. The median
expansion velocity of the Class II YSOs is found to be 0.0 $\pm$ 0.1 km s-1.
Therefore, the YSOs associated with the gas may not exhibit expansion on the
plane of the sky.

Figure 14: (a) Distribution of the proper motions of the YSOs, where the
average proper motion of these sources is (2.65, -3.17) ${\rm mas~{}yr^{-1}}$,
and the $1\sigma$ uncertainty range of the Gaussian fitting is (0.14, 0.17)
${\rm mas~{}yr^{-1}}$. (b) Projected distribution and tangential velocity of
the YSOs. The magenta arrows depict the velocity for each YSO with respect to
the mean velocity of the YSO sample. The red asterisk (0, 0) is the center of
the bubble.

Figure 15: P-V diagrams of the Class II YSOs used to investigate the kinematic
features.

##  5 Summary

We observed four molecular lines (12CO, 13CO, C18O, and HCO+ $J$ = 1–0), and
collected YSOs in S Mon to investigate the bubble, an outflow identified in
one of the clumps, as well as the energy cascade in the S Mon region. This
includes the dynamical timescale, possible driving sources, and feedback of
the outflow and bubble. The main results are summarized as follows.

  * 1

The Class II YSOs are closely associated with the molecular gas in terms of
spatial distribution and kinematic properties. According to the _Gaia_ DR3
parallaxes of the Class II YSOs, we estimate the distance to S Mon to be
$\sim$ 722 $\pm$ 9 pc.

  * 2

We discover a molecular bubble in the S Mon region which is consistent with
the Herschel H2 column density maps. The molecular bubble has a radius of
$\sim 1.1\pm 0.1$ pc, a mass of $\sim 1\,400\pm 400~{}M_{\odot}$, a momentum
of $\sim$ 4 600 $\pm$ 1 300 $M_{\odot}$ km s-1, and kinetic energy of
$\sim(1.5\pm 0.4)\times 10^{47}$ erg. The dynamical timescale of the HII
region of 15 Mon is (3.5–8.9) $\times 10^{5}$ yr, which is likely to be the
main driving source of the bubble.

  * 3

We detect a molecular outflow in the shell of the bubble, which indicates that
star formation activity is ongoing around the bubble. The outflow has a mass
of a few $M_{\odot}$, a momentum of a few tens of $M_{\odot}~{}{\rm
km~{}s^{-1}}$, kinetic energy of $\sim 10^{45}$ erg, a dynamical timescale of
$\sim 10^{5}$ yr, and a mechanical luminosity of a few times
$10^{-1}~{}L_{\odot}$.

  * 4

The wind energy helps to sustain the turbulence in the S Mon region and drive
the bubble. The bubble can also provide sufficient energy input to sustain the
cloud turbulence. The measured outflow in clump E is enough to maintain the
turbulence of the clump E. Outflows (like the one detected in clump E) could
generally contribute to the turbulence of a cloud, while in S Mon further
outflows are not yet confirmed.

  * 5

The Class II YSOs correlated to the molecular bubble do not show expansion on
the plane of the sky.

We appreciate the anonymous referee for the instructive comments that helped
us to improve the paper. This work was funded by the NSFC Grands 11933011,
National SKA Program of China (Grant No. 2022SKA0120103) and the Key
Laboratory for Radio Astronomy. L.Y.J. thanks the support of the NSFC grant
No. 12203104, the Natural Science Foundation of Jiangsu Province (grant No.
BK20210999), the Entrepreneurship and Innovation Program of Jiangsu Province.
This work has made use of data from the European Space Agency (ESA) mission
_Gaia_ (https://www.cosmos.esa.int/gaia) processed by the _Gaia_ Data
Processing and Analysis Consortium (DPAC,
https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has
been provided by national institutions, in particular the institutions
participating in the _Gaia_ Multilateral Agreement.

## References

  * Abdurro’uf et al. (2022) Abdurro’uf, Accetta, K., Aerts, C., et al. 2022, ApJS, 259, 35, doi: 10.3847/1538-4365/ac4414
  * Arce et al. (2011) Arce, H. G., Borkin, M. A., Goodman, A. A., Pineda, J. E., & Beaumont, C. N. 2011, ApJ, 742, 105, doi: 10.1088/0004-637X/742/2/105
  * Arce et al. (2010) Arce, H. G., Borkin, M. A., Goodman, A. A., Pineda, J. E., & Halle, M. W. 2010, ApJ, 715, 1170, doi: 10.1088/0004-637X/715/2/1170
  * Bally (2016) Bally, J. 2016, ARA&A, 54, 491, doi: 10.1146/annurev-astro-081915-023341
  * Bontemps et al. (1996) Bontemps, S., Andre, P., Terebey, S., & Cabrit, S. 1996, A&A, 311, 858
  * Broos et al. (2013) Broos, P. S., Getman, K. V., Povich, M. S., et al. 2013, ApJS, 209, 32, doi: 10.1088/0067-0049/209/2/32
  * Buckle et al. (2012) Buckle, J. V., Richer, J. S., & Davis, C. J. 2012, MNRAS, 423, 1127, doi: 10.1111/j.1365-2966.2012.20941.x
  * Castor et al. (1975) Castor, J., McCray, R., & Weaver, R. 1975, ApJ, 200, L107, doi: 10.1086/181908
  * Chen et al. (2013) Chen, Y., Zhou, P., & Chu, Y.-H. 2013, ApJ, 769, L16, doi: 10.1088/2041-8205/769/1/L16
  * Churchwell et al. (2006) Churchwell, E., Povich, M. S., Allen, D., et al. 2006, ApJ, 649, 759, doi: 10.1086/507015
  * Cody et al. (2014) Cody, A. M., Stauffer, J., Baglin, A., et al. 2014, AJ, 147, 82, doi: 10.1088/0004-6256/147/4/82
  * Dahm (2008) Dahm, S. E. 2008, in Handbook of Star Forming Regions, Volume I, ed. B. Reipurth, Vol. 4, 966
  * Dale (2015) Dale, J. E. 2015, New A Rev., 68, 1, doi: 10.1016/j.newar.2015.06.001
  * Dyson & Williams (1980) Dyson, J. E., & Williams, D. A. 1980, Physics of the interstellar medium
  * Fazio et al. (2004) Fazio, G. G., Hora, J. L., Allen, L. E., et al. 2004, ApJS, 154, 10, doi: 10.1086/422843
  * Feddersen et al. (2018) Feddersen, J. R., Arce, H. G., Kong, S., et al. 2018, ApJ, 862, 121, doi: 10.3847/1538-4357/aacaf2
  * Feigelson et al. (2013) Feigelson, E. D., Townsley, L. K., Broos, P. S., et al. 2013, ApJS, 209, 26, doi: 10.1088/0067-0049/209/2/26
  * Fűrész et al. (2008) Fűrész, G., Hartmann, L. W., Megeath, S. T., Szentgyorgyi, A. H., & Hamden, E. T. 2008, ApJ, 676, 1109, doi: 10.1086/525844
  * Flaccomio et al. (2023) Flaccomio, E., Micela, G., Peres, G., et al. 2023, A&A, 670, A37, doi: 10.1051/0004-6361/202244872
  * Frank et al. (2014) Frank, A., Ray, T. P., Cabrit, S., et al. 2014, in Protostars and Planets VI, ed. H. Beuther, R. S. Klessen, C. P. Dullemond, & T. Henning, 451, doi: 10.2458/azu_uapress_9780816531240-ch020
  * Gaia Collaboration et al. (2016) Gaia Collaboration, Prusti, T., de Bruijne, J. H. J., et al. 2016, A&A, 595, A1, doi: 10.1051/0004-6361/201629272
  * Gaia Collaboration et al. (2023) Gaia Collaboration, Vallenari, A., Brown, A. G. A., et al. 2023, A&A, 674, A1, doi: 10.1051/0004-6361/202243940
  * Garden et al. (1991) Garden, R. P., Hayashi, M., Gatley, I., Hasegawa, T., & Kaifu, N. 1991, ApJ, 374, 540, doi: 10.1086/170143
  * Geen et al. (2015) Geen, S., Rosdahl, J., Blaizot, J., Devriendt, J., & Slyz, A. 2015, MNRAS, 448, 3248, doi: 10.1093/mnras/stv251
  * Gilmore et al. (2012) Gilmore, G., Randich, S., Asplund, M., et al. 2012, The Messenger, 147, 25
  * Grady et al. (1984) Grady, C. A., Snow, T. P., & Cash, W. C. 1984, ApJ, 283, 218, doi: 10.1086/162296
  * Green et al. (2019) Green, G. M., Schlafly, E., Zucker, C., Speagle, J. S., & Finkbeiner, D. 2019, ApJ, 887, 93, doi: 10.3847/1538-4357/ab5362
  * Greene et al. (1994) Greene, T. P., Wilking, B. A., Andre, P., Young, E. T., & Lada, C. J. 1994, ApJ, 434, 614, doi: 10.1086/174763
  * Großschedl et al. (2021) Großschedl, J. E., Alves, J., Meingast, S., & Herbst-Kiss, G. 2021, A&A, 647, A91, doi: 10.1051/0004-6361/202038913
  * Großschedl et al. (2018) Großschedl, J. E., Alves, J., Meingast, S., et al. 2018, A&A, 619, A106, doi: 10.1051/0004-6361/201833901
  * Guo et al. (2021) Guo, H. L., Chen, B. Q., Yuan, H. B., et al. 2021, ApJ, 906, 47, doi: 10.3847/1538-4357/abc68a
  * Hacar et al. (2016) Hacar, A., Alves, J., Forbrich, J., et al. 2016, A&A, 589, A80, doi: 10.1051/0004-6361/201527805
  * Kerr & Lynden-Bell (1986) Kerr, F. J., & Lynden-Bell, D. 1986, MNRAS, 221, 1023, doi: 10.1093/mnras/221.4.1023
  * Lada (1985) Lada, C. J. 1985, ARA&A, 23, 267, doi: 10.1146/annurev.aa.23.090185.001411
  * Lada (1987) Lada, C. J. 1987, in Star Forming Regions, ed. M. Peimbert & J. Jugaku, Vol. 115, 1
  * Lada & Wilking (1984) Lada, C. J., & Wilking, B. A. 1984, ApJ, 287, 610, doi: 10.1086/162719
  * Lada et al. (2006) Lada, C. J., Muench, A. A., Luhman, K. L., et al. 2006, AJ, 131, 1574, doi: 10.1086/499808
  * Lallement et al. (2019) Lallement, R., Babusiaux, C., Vergely, J. L., et al. 2019, A&A, 625, A135, doi: 10.1051/0004-6361/201834695
  * Lamers & Cassinelli (1999) Lamers, H. J. G. L. M., & Cassinelli, J. P. 1999, Introduction to Stellar Winds
  * Li et al. (2015) Li, H., Li, D., Qian, L., et al. 2015, ApJS, 219, 20, doi: 10.1088/0067-0049/219/2/20
  * Li et al. (2018) Li, Y., Li, F.-C., Xu, Y., et al. 2018, ApJS, 235, 15, doi: 10.3847/1538-4365/aaab67
  * Li et al. (2020) Li, Y., Xu, Y., Sun, Y., & Yang, J. 2020, ApJS, 251, 26, doi: 10.3847/1538-4365/abc34b
  * Lindegren et al. (2021) Lindegren, L., Bastian, U., Biermann, M., et al. 2021, A&A, 649, A4, doi: 10.1051/0004-6361/202039653
  * Liu et al. (2021) Liu, D.-J., Xu, Y., Li, Y.-J., et al. 2021, ApJS, 253, 15, doi: 10.3847/1538-4365/abcece
  * Martins et al. (2005) Martins, F., Schaerer, D., & Hillier, D. J. 2005, A&A, 436, 1049, doi: 10.1051/0004-6361:20042386
  * McKee (1989) McKee, C. F. 1989, ApJ, 345, 782, doi: 10.1086/167950
  * McKee & Ostriker (2007) McKee, C. F., & Ostriker, E. C. 2007, ARA&A, 45, 565, doi: 10.1146/annurev.astro.45.051806.110602
  * Montillaud et al. (2019) Montillaud, J., Juvela, M., Vastel, C., et al. 2019, A&A, 631, A3, doi: 10.1051/0004-6361/201834903
  * Nakamura & Li (2007) Nakamura, F., & Li, Z.-Y. 2007, The Astrophysical Journal, 662, 395, doi: 10.1086/517515
  * Nakamura & Li (2014) Nakamura, F., & Li, Z.-Y. 2014, ApJ, 783, 115, doi: 10.1088/0004-637X/783/2/115
  * Nony et al. (2021) Nony, T., Robitaille, J. F., Motte, F., et al. 2021, A&A, 645, A94, doi: 10.1051/0004-6361/202039353
  * Perryman et al. (1997) Perryman, M. A. C., Lindegren, L., Kovalevsky, J., et al. 1997, A&A, 323, L49
  * Pety (2005) Pety, J. 2005, in SF2A-2005: Semaine de l’Astrophysique Francaise, ed. F. Casoli, T. Contini, J. M. Hameury, & L. Pagani, 721
  * Pilbratt et al. (2010) Pilbratt, G. L., Riedinger, J. R., Passvogel, T., et al. 2010, A&A, 518, L1, doi: 10.1051/0004-6361/201014759
  * Randich et al. (2013) Randich, S., Gilmore, G., & Gaia-ESO Consortium. 2013, The Messenger, 154, 47
  * Rapson et al. (2014) Rapson, V. A., Pipher, J. L., Gutermuth, R. A., et al. 2014, ApJ, 794, 124, doi: 10.1088/0004-637X/794/2/124
  * Reid et al. (2009) Reid, M. J., Menten, K. M., Zheng, X. W., et al. 2009, ApJ, 700, 137, doi: 10.1088/0004-637X/700/1/137
  * Rieke et al. (2004) Rieke, G. H., Young, E. T., Engelbracht, C. W., et al. 2004, ApJS, 154, 25, doi: 10.1086/422717
  * Shan et al. (2012) Shan, W., Yang, J., Shi, S., et al. 2012, IEEE Transactions on Terahertz Science and Technology, 2, 593, doi: 10.1109/TTHZ.2012.2213818
  * Simon et al. (2001) Simon, R., Jackson, J. M., Clemens, D. P., Bania, T. M., & Heyer, M. H. 2001, ApJ, 551, 747, doi: 10.1086/320230
  * Skrutskie et al. (2006) Skrutskie, M. F., Cutri, R. M., Stiening, R., et al. 2006, AJ, 131, 1163, doi: 10.1086/498708
  * Snow et al. (1981) Snow, T. P., J., Cash, W., & Grady, C. A. 1981, ApJ, 244, L19, doi: 10.1086/183470
  * Sun et al. (2018) Sun, J. X., Lu, D. R., Yang, J., et al. 2018, Acta Astronomica Sinica, 59, 3
  * Sung & Bessell (2010) Sung, H., & Bessell, M. S. 2010, AJ, 140, 2070, doi: 10.1088/0004-6256/140/6/2070
  * Sung et al. (2008) Sung, H., Bessell, M. S., Chun, M.-Y., Karimov, R., & Ibrahimov, M. 2008, AJ, 135, 441, doi: 10.1088/0004-6256/135/2/441
  * Sung et al. (2009) Sung, H., Stauffer, J. R., & Bessell, M. S. 2009, AJ, 138, 1116, doi: 10.1088/0004-6256/138/4/1116
  * Tauber et al. (1993) Tauber, J. A., Lis, D. C., & Goldsmith, P. F. 1993, ApJ, 403, 202, doi: 10.1086/172194
  * Teixeira et al. (2012) Teixeira, P. S., Lada, C. J., Marengo, M., & Lada, E. A. 2012, A&A, 540, A83, doi: 10.1051/0004-6361/201015326
  * Tobin et al. (2015) Tobin, J. J., Hartmann, L., Fűrész, G., Hsu, W.-H., & Mateo, M. 2015, AJ, 149, 119, doi: 10.1088/0004-6256/149/4/119
  * Tobin et al. (2009) Tobin, J. J., Hartmann, L., Furesz, G., Mateo, M., & Megeath, S. T. 2009, ApJ, 697, 1103, doi: 10.1088/0004-637X/697/2/1103
  * Turner et al. (1997) Turner, B. E., Pirogov, L., & Minh, Y. C. 1997, ApJ, 483, 235, doi: 10.1086/304228
  * Venuti et al. (2018) Venuti, L., Prisinzano, L., Sacco, G. G., et al. 2018, A&A, 609, A10, doi: 10.1051/0004-6361/201731103
  * Walch et al. (2012) Walch, S. K., Whitworth, A. P., Bisbas, T., Wünsch, R., & Hubber, D. 2012, MNRAS, 427, 625, doi: 10.1111/j.1365-2966.2012.21767.x
  * Weisskopf et al. (2002) Weisskopf, M. C., Brinkman, B., Canizares, C., et al. 2002, PASP, 114, 1, doi: 10.1086/338108
  * Xu et al. (2017) Xu, J.-L., Xu, Y., Yu, N., et al. 2017, ApJ, 849, 140, doi: 10.3847/1538-4357/aa8ee0
  * Yang et al. (2018) Yang, A. Y., Thompson, M. A., Urquhart, J. S., & Tian, W. W. 2018, ApJS, 235, 3, doi: 10.3847/1538-4365/aaa297
  * Zhang (2023) Zhang, M. 2023, ApJS, 265, 115, doi: 10.1088/0004-637X/783/2/115
  * Zhang et al. (2015) Zhang, Q., Wang, K., Lu, X., & Jiménez-Serra, I. 2015, ApJ, 804, 141, doi: 10.1088/0004-637X/804/2/141
  * Zucker et al. (2020) Zucker, C., Speagle, J. S., Schlafly, E. F., et al. 2020, A&A, 633, A51, doi: 10.1051/0004-6361/201936145

##  Appendix A YSO Classification Scheme

We describe the classification scheme for YSOs in this section. The YSO
candidates in S Mon are collected from the catalogs of Sung et al. (2009),
Broos et al. (2013), Cody et al. (2014), and Rapson et al. (2014), within a
field of $202.8^{\circ}<l<203.1^{\circ},2.0^{\circ}<b<2.3^{\circ}$. We
consider two stars separated by less than 1$\arcsec$ to be the same star, and
remove duplicated YSO candidates. A total of 1 023 YSO candidates were
initially selected. However, this initial catalog not only included numerous
samples that were not YSOs, but also the YSOs in the catalog were obtained
through different classification schemes. To maintain the homogeneity of the
YSO catalog, we classified the YSO candidates based on a homogeneous method.

We applied the spectral index of the SED to classify the evolutionary stage of
the YSOs, where $\alpha={\rm d}\log\lambda F_{\lambda}/{\rm d}\log\lambda$
(Greene et al., 1994). To build this scheme, we use photometric data from the
2MASS $K_{s}$ band (Skrutskie et al., 2006), four IRAC bands (3.6, 4.2, 5.8,
and 8.0 $\mu$m) of Spitzer (Fazio et al., 2004), and the MIPS 24 $\mu$m band
(Rieke et al., 2004). Note that not all six bands are required simultaneously
for the classification; however, at least three IRAC bands are necessary. Only
sources with photometric uncertainties $\sigma<0.2$ mag were selected for
classification. Finally, we compiled a sample of seven Class I, 160 Class II,
and 245 Class III YSOs. Figure 16 shows color–color diagrams of the classified
YSOs in the S Mon region.

To obtain the parallaxes and proper motions of the classified YSOs, we cross
matched them with stars in Gaia DR3 within 1$\arcsec$. One Class I, 145 Class
II and 231 Class III YSOs were assigned astrometric parameters from Gaia DR3.

In the next step we added RV data to our YSO sample. First, we checked the
availability of _G_ aia DR3 RVs. There were only 16 YSOs with measured Gaia
DR3 RVs, of which most fell outside the velocity range of interest (about 63%
are beyond the range of 0–15 km s-1), and/or they had large RV uncertainties.
Therefore, we did not use Gaia RVs in our analysis. The RVs of the classified
YSOs were collected from APOGEE-2 DR17 (Abdurro’uf et al., 2022) and the
catalog of Tobin et al. (2015), who measured the RVs of 695 stars in NGC 2264
(using multi-fiber echelle spectroscopy at the 6.5 m MMT and Magellan
telescopes). We cross matched our YSO candidates (1 023) with the APOGEE
catalog using the source_id in the Gaia DR3 catalog, and with the catalog of
Tobin et al. (2015) within 1$\arcsec$, yielding 73 YSOs with RVs from APOGEE
and 148 YSOs with RVs from Tobin et al. (2015). A total of 44 YSO candidates
were matched in both catalogs. We made a comparison of the $V_{\rm LSR}$
(converted from heliocentric RVs to RVs relative to the local standard of
rest, LSR) collected from the two catalogs for the YSO candidates, as shown in
Figure 17. The average difference between them is 0.7 km s-1, and the
$1\sigma$ dispersion is 1.0 km s-1. The Pearson’s correlation coefficient
between them is 0.99. We therefore added the RVs from Tobin et al. (2015) to
our sample and applied a systematic shift of 0.7 km s-1, applied because the
RV errors in the APOGEE catalog are consistently smaller than those in the
catalog of Tobin et al. (2015). Finally, we obtained 47 Class II, and 91 Class
III YSOs with RVs.

To sum up, we obtained seven Class I, 160 Class II, and 245 Class III YSOs in
the S Mon region. Among them, one Class I, 145 Class II, and 231 Class III
YSOs have counterparts in _Gaia_ DR3, where 47 Class II and 91 Class III YSOs
have RVs. Table 5 describes the YSO samples.

Figure 16: Color–color diagrams showing the classified YSOs. The green
diamonds, red triangles, and magenta crosses present Class I, Class II, Class
III YSOs, respectively. Figure 17: Comparison of the $V_{\rm LSR}$ values in
the APOGEE and Tobin et al. (2015) samples. The black line represents where
$V_{\rm LSR}$ is the same in both catalogs. The red dashed line indicates the
mean value of the difference between the APOGEE and Tobin et al. (2015)
$V_{\rm LSR}$ values, and the gray filled region depicts the corresponding
$1\sigma$ region. Table 5: Description of the YSO Samples

Column | Format | Unit | Description
---|---|---|---
ID | A4 | — | ID
Gaia | I19 | — | Gaia DR3 identifier
RAdeg | F8.4 | deg | Right ascension (J2000)
DEdeg | F8.4 | deg | Declination (J2000)
plx | F7.3 | mas | Gaia DR3 parallax
e_plx | F7.3 | mas | Uncertainty in the parallax
pmRA | F7.3 | mas/yr | Proper motion in the RA* direction
e_pmRA | F7.3 | mas/yr | Uncertainty in pmRA
pmDE | F7.3 | mas/yr | Proper motion in the DEC direction
e_pmDE | F7.3 | mas/yr | Uncertainty in pmDE
Vlsr | F7.1 | km/s | LSR velocity
e_Vlsr | F7.1 | km/s | Uncertainty in the LSR velocity
Jmag | F5.2 | mag | 2MASS J band magnitude
e_Jmag | F5.2 | mag | Uncertainty in Jmag
Hmag | F5.2 | mag | 2MASS H band magnitude
e_Hmag | F5.2 | mag | Uncertainty in Hmag
Ksmag | F5.2 | mag | 2MASS Ks band magnitude
e_Ksmag | F5.2 | mag | Uncertainty in Ksmag
3.6mag | F5.2 | mag | Spitzer/IRAC 3.6 micron band magnitude
e_3.6mag | F5.2 | mag | Uncertainty in 3.6mag
4.5mag | F5.2 | mag | Spitzer/IRAC 4.5 micron band magnitude
e_4.5mag | F5.2 | mag | Uncertainty in 4.5mag
5.8mag | F5.2 | mag | Spitzer/IRAC 5.8 micron band magnitude
e_5.8mag | F5.2 | mag | Uncertainty in 5.8mag
8.0mag | F5.2 | mag | Spitzer/IRAC 8.0 micron band magnitude
e_8.0mag | F5.2 | mag | Uncertainty in 8.0mag
24mag | F5.2 | mag | Spitzer/MIPS 24 micron band magnitude
e_24mag | F5.2 | mag | Uncertainty in 24mag
Class | A10 | — | Source classification
OClass | A10 | — | Other source classification
Ref | A15 | — | Reference
Ref_RV | A15 | — | Reference for the LSR velocity

Note. — This table is available in its entirety in machine-readable form.

##  Appendix B Calculation of the Physical Parameters

In this study, the properties of the molecular gas are estimated under the
assumption that all of the molecular gas is in LTE. Optically thin 13CO lines
are employed to estimate the column density and mass of the molecular gas.
According to Garden et al. (1991), the 13CO column density can be estimated
via:

$N({\rm{}^{13}CO})=4.71\times 10^{13}\frac{T_{\rm ex}+0.88}{\exp{(-5.29/T_{\rm
ex}})}\times\frac{\tau}{1-\exp{(-\tau)}}\int T_{\rm mb,13}dv,$ (B1)

where $\tau$ is the optical depth, $T_{\rm ex}$ is the mean excitation
temperature of the molecular gas, and $T_{\rm mb,13}$ is the main-beam
temperature of 13CO.

We assume that the 12CO emission is optically thick. The excitation
temperature, $T_{\rm ex}$, can be estimated by (Garden et al., 1991):

$T_{\rm ex}=\frac{5.53}{\ln{[1+5.53/(T_{\rm mb,12}+0.82)]}},$ (B2)

where $T_{\rm mb,12}$ is the peak main-beam temperature of 12CO.

The optical depth can be derived using the following equation (Garden et al.,
1991):

$\tau({\rm{}^{13}CO})=-\ln{\left[1-\frac{T_{\rm
mb,13}}{5.29/[\exp{(5.29/T_{\rm ex})}-1]-0.89}\right]}.$ (B3)

The relation $N_{\rm H_{2}}/N_{\rm{}^{13}CO}\approx 5\times 10^{5}$ (Simon et
al., 2001) is used to estimate the H2 column density. The mass, $M$, can be
determined by:

$M=\mu m_{\rm H}N({\rm H_{2}})S,$ (B4)

where $\mu=2.72$ is the mean molecular weight, $m_{\rm H}$ is the mass of the
hydrogen atom (Garden et al., 1991), and $S$ is the projected 2D area.

The turbulent energy, $E_{\rm turb}$, can be given approximately by:

$E_{\rm turb}=\frac{1}{2}M\sigma_{\rm 3d}^{2},$ (B5)

where $\sigma_{\rm 3d}$ is the 3D turbulent velocity dispersion, which can be
calculated by:

$\sigma_{\rm 3d}=\sqrt{3}\sigma_{\rm
1d}=\frac{\sqrt{3}}{2\sqrt{2\ln{2}}}\Delta V_{\rm FWHM},$ (B6)

where $\Delta V_{\rm FWHM}$ is the 1D FWHM velocity dispersion based on
typical 13CO.

The gravitational binding energy, $E_{\rm grav}$, can be calculated by:

$E_{\rm grav}=\frac{GM^{2}}{R},$ (B7)

where $R$ is the radius.

The excitation temperatures of the eight clumps are between 21.3 and 33.2 K.
Since these clumps are regularly located within the bubble, a mean excitation
temperature of about 26.9 K is adopted as that of the bubble. The mean column
density of the bubble is $\sim 2.5\times 10^{20}~{}{\rm cm^{-2}}$ and its mass
is $\sim 1\,400~{}M_{\odot}$. In this work, we consider the uncertainties in
gas temperature, shell region, and distance for the shell mass estimation. The
total uncertainty caused by these factors is about $\pm 400~{}M_{\odot}$.

Figure 18 shows the projected area of the cloud. The red lines denoted its
boundary, determined by pixels whose main-beam brightness temperature is
larger than $3\times{\rm rms}$ in at least three successive channels. The
velocity range adopted here is between 7 and 15 km s-1 and the mass of the
cloud is derived as 2 200 $\pm$ 500 $M_{\odot}$. The uncertainties in the
excitation temperature and distance are considered in the estimation of the
cloud mass. $\sigma_{\rm 1d}$ is obtained from Gaussian fitting of the average
13CO line, which is about 1.0 $\pm$ 0.1 km s-1. The turbulent energy and
gravitational binding energy of the cloud are estimated by Equation (B5) and
Equation (B7), respectively, which correspond to (6.3 $\pm$ 2.5)$\times
10^{46}$ erg and (2.1 $\pm$ 1.0)$\times 10^{47}$ erg.

The clump’s column density ($N_{\rm H_{2}}$), mass ($M_{\rm clump}$),
turbulent energy ($E_{\rm turb,clump}$), and gravitational binding energy
($E_{\rm grav,clump}$) were calculated using the same equations as applied to
the cloud. The uncertainties in the distance and clump region were considered
in the calculation.

Figure 18: Intensity map of 13CO for the velocity range 7–15 km s-1. The red
lines depict the boundary of the cloud.

##  Appendix C Distribution of OB stars

Figure 19 shows the projected distributions and relative proper motions of the
O- and B-type stars in the S Mon region. Table 6 lists the parameters of these
stars, which are collected from SIMBAD. The astrometry parameters of 15 Mon
are from Hipparcos catalog (Perryman et al., 1997), and those of other sources
are from _G_ aia DR3 (Gaia Collaboration et al., 2023).

Figure 19: Projected distributions and relative proper motions of the O-(red
triangle) and B-type (blue diamonds) stars. The average proper motion utilized
being (2.65, -3.17) mas yr-1. Table 6: Summary of the Parameters of the O- and
B-type Stars

Identifier | SpType | $\alpha$ | $\delta$ | $\varpi$ | $\sigma_{\varpi}$ | $\mu_{\alpha}\cos\delta$ | $\sigma_{\mu_{\alpha}\cos\delta}$ | $\mu_{\delta}$ | $\sigma_{\mu_{\delta}}$ | Gaia DR3 ID
---|---|---|---|---|---|---|---|---|---|---
(deg) | (deg) | (mas) | (mas) | (mas yr-1) | (mas yr-1) | (mas yr-1) | (mas yr-1)
(1) | (2) | (3) | (4) | (5) | (6) | (7) | (8) | (9) | (10) | (11)
15 Mon† | O7V+B1.5/2V | 100.2444 | 9.8958 | 3.550 | 0.500 | -2.610 | 0.560 | -1.610 | 0.390 | 3326740865570982400
15 Mon B | B1.5/2V | 100.2440 | 9.8951 | 1.402 | 0.098 | -1.971 | 0.113 | -4.225 | 0.082 | 3326740865571935232
V* V641 Mon | B1.5IV+(B2V) | 100.1191 | 9.8179 | 1.389 | 0.045 | -1.943 | 0.051 | -3.775 | 0.040 | 3326716813754146176
LS VI +09 13 | B1V | 100.1552 | 9.7916 | … | … | … | … | … | … | 3326716470156761472
V* V684 Mon | B2.5V | 100.1599 | 9.7878 | 1.540 | 0.041 | -1.482 | -0.050 | -3.960 | -0.041 | 3326715714242517248
HD 47961 | B2.5V | 100.3638 | 9.8540 | 1.248 | 0.079 | -1.976 | 0.091 | -4.034 | 0.077 | 3326737120359493504
HD 261810 | B2.5Vn | 100.1801 | 9.7671 | 1.424 | 0.033 | -2.172 | 0.042 | -3.942 | 0.033 | 3326715439364610816
HD 261938 | B3V | 100.2578 | 9.8800 | 1.540 | 0.059 | -1.258 | 0.056 | -3.247 | 0.041 | 3326740693772848896
EM* LkHA 25 | B4Ve | 100.1860 | 9.8006 | 1.414 | 0.111 | -1.076 | 0.140 | -3.826 | 0.125 | 3326715851681467904
HD 261878 | B6V | 100.2148 | 9.8637 | 1.331 | 0.037 | -1.504 | 0.041 | -3.747 | 0.030 | 3326717260430731648
CSI+09-06383 | B7 | 100.2750 | 9.8833 | … | … | … | … | … | … | …
HD 261841 | B8IV-Ve | 100.2037 | 9.8623 | 1.408 | 0.021 | -1.722 | 0.023 | -3.568 | 0.018 | 3326717226070995328
HD 261969 | B9IV | 100.2933 | 9.8838 | 1.351 | 0.036 | -1.744 | 0.032 | -3.916 | 0.028 | 3326740006577519360
NGC 2264 181 | B9Vn | 100.2968 | 9.8821 | 1.346 | 0.029 | -1.576 | -0.032 | -3.880 | 0.027 | 3326740002281694592

Note. — (1) Identifier. (2) Spectral type. (3) Right ascension (J2000). (4)
Declination (J2000). (5) Parallax. (6) Uncertainty in parallax; (7) Proper
motion in RA* direction. (8) Uncertainty in $\mu_{\alpha}\cos\delta$. (9)
Proper motion in DEC direction. (10) Uncertainty in $\mu_{\delta}$. (11) Gaia
DR3 identifier.

†: The astrometry parameters of 15 Mon (HIP 31978) are from Hipparcos catalog,
and those of other sources are from Gaia DR3.