FOREST Unbiased Galactic plane Imaging survey with the Nobeyama 45 m telescope (FUGIN): Molecular clouds toward W 33; possible evidence for a cloud–cloud collision triggering O star formation

Abstract

We observed molecular clouds in the W 33 high-mass star-forming region associated with compact and extended H ii regions using the NANTEN2 telescope as well as the Nobeyama 45 m telescope in the J = 1–0 transitions of ¹²CO, ¹³CO, and C¹⁸O as part of the FOREST Unbiased Galactic plane Imaging survey with the Nobeyama 45 m telescope (FUGIN) legacy survey. We detected three velocity components at 35 km s⁻¹, 45 km s⁻¹, and 58 km s⁻¹. The 35 km s⁻¹ and 58 km s⁻¹ clouds are likely to be physically associated with W 33 because of the enhanced ¹²CO J = 3–2 to J = 1–0 intensity ratio as |$R_{\rm 3\mbox{-}2/1\mbox{-}0} > 1.0$| due to the ultraviolet irradiation by OB stars, and morphological correspondence between the distributions of molecular gas and the infrared and radio continuum emissions excited by high-mass stars. The two clouds show complementary distributions around W 33. The velocity separation is too large to be gravitationally bound, and yet not explained by expanding motion by stellar feedback. Therefore, we discuss whether a cloud–cloud collision scenario likely explains the high-mass star formation in W 33.

1 Introduction

1.1 High-mass star formation

High-mass stars have a huge influence on the interstellar medium (ISM) and galactic evolution via stellar feedback and supernova explosions. Feedback from high-mass stars is supposed to induce the formation of next-generation stars as a sequential star formation (Elmegreen & Lada 1977; Lada 1987). Supernova explosions scatter heavy elements in the ISM that drives chemical evolution of galaxies. However, the formation mechanism of high-mass stars is not clearly understood because it is difficult to achieve high mass accretion rate (∼10⁻³ M_⊙ yr⁻¹; Wolfire & Cassinelli 1987). It is, therefore, an important issue to investigate what the necessary condition is for high-mass star formation in astrophysics.

Based on theoretical studies, two scenarios for high-mass star formation are proposed: core accretion (monolithic collapse) and competitive accretion (see Zinnecker & Yorke 2007; Tan et al. 2014 for reviews). In the core accretion (monolithic collapse) model, high-mass stars are formed by the collapse of isolated gravitationally bound massive cores. It is a similar process to low-mass star formation with more massive aggregation (e.g., Nakano et al. 2000; Yorke & Sonnhalter 2002; McKee & Tan 2003; Krumholz et al. 2007, 2009; Hosokawa & Omukai 2009). On the other hand, for the competitive accretion model, high-mass stars are formed by growth of low-mass protostellar seeds by mass accretion from surrounding gas (e.g., Bonnell et al. 1997, 2001, 2004). One of the differences between the two models is the initial condition of the natal cloud. It is therefore important to observe high-mass star-forming regions at a very early stage of evolution which holds the initial condition of high-mass star formation in order to verify the theories; such a comparison between theories and observations has not been well made so far, and we do not have compelling evidence for either theory (Tan et al. 2014). Difficulties lie in the rareness of high-mass star-forming regions in the solar neighborhood, and the short timescale of the feedback processes that are heavily mixed up with star formation processes.

1.2 Cloud–cloud collisions as a trigger of high-mass star formation

During the past ten years, it has become increasingly probable that cloud–cloud collisions play an important role in high-mass star formation. In observational studies, a cloud–cloud collision was first reported in the star-forming region NGC 1333 by Loren (1976). The Sagittarius B2 star-forming region in the Galactic center was also suggested to have starburst triggered by a cloud–cloud collision based on the complementary distributions between the cloud of two velocity components (Hasegawa et al. 1994; Sato et al. 2000). In 2009, molecular observations with the NANTEN2 telescope showed two molecular clouds with different radial velocities toward a super star cluster Westerlund 2, and their distributions are interpreted as that the cluster formation was triggered by the collision of the two clouds (Furukawa et al. 2009). Further evidence for the physical association of the two molecular clouds with the cluster reinforced the interpretation (Ohama et al. 2010). Subsequently, three other super star clusters including 10–20 O stars were found to be associated with two molecular clouds with different velocities, and formation of O stars triggered by a cloud–cloud collision is likely to be a common process for clusters having more than 10 O stars (NGC 3603: Fukui et al. 2014; RCW 38: Fukui et al. 2016; DBS[2003]179: S. Kuwahara et al. in preparation). Among the eight super star clusters listed in the review article of Portegies Zwart, McMillan, and Gieles (2010), only three are known to be associated with localized nebulosities, indicating that the three are young and still associated with the remnant of natal molecular gas without heavy ionization. Additional possible cases of multiple O star formation triggered by a cloud–cloud collision are reported for M 42 (Fukui et al. 2017a), NGC 6334-NGC 6357 (Fukui et al. 2018a), M 17 (Nishimura et al. 2018), W 49 A (Miyawaki et al. 1986, 2009; Buckley & Ward-Thompson 1996), W 51 (Okumura et al. 2001; Kang et al. 2010; Fujita et al. 2017), and R136 (Fukui et al. 2017b). In addition to the above, many star-forming regions and dense clumps in the Milky Way have been suggested to be triggered by a cloud–cloud collision: LkHα198 (Loren 1977); IRAS 19550+3248 (Koo et al. 1994); IRAS 2306+1451 (Vallee 1995); M 20 (Torii et al. 2011, 2017a); RCW 120 (Torii et al. 2015); N37 (Baug et al. 2016); GM 24 (Fukui et al. 2018b); M 16 (Nishimura et al. 2017); RCW 34 (Hayashi et al. 2018); RCW 36 (Sano et al. 2018); RCW 166 (Ohama et al. 2018b); S116, S117, S118 (Fukui et al. 2018c); RCW 79 (Ohama et al. 2018a); Sh2-48 (Torii et al. 2017b); NGC 2024 (Ohama et al. 2017); NGC 2068, NGC 2071 (Tsutsumi et al. 2017); NGC 2359 (Sano et al. 2017); S87 (Xue & Wu 2008); S87E, S88B, AFGL 5142, AFGL 5180 (Higuchi et al. 2010); G0.253+0.016 (Higuchi et al. 2014); Circinus-E cloud (Shimoikura & Dobashi 2011); Sh2-252 (Shimoikura et al. 2013); L 1641-N (Nakamura et al. 2012); Serpens Main Cluster (Duarte-Cabral et al. 2011); Serpens South (Nakamura et al. 2014); L 1004E in the Cygnus OB 7 (Dobashi et al. 2014); G35.20−0.74 (Dewangan 2017); S235 (Dewangan & Ojha 2017); L 1188 (Gong et al. 2017); the Galactic Center 50 km s⁻¹ molecular cloud (Tsuboi et al. 2015); and N 159 West (Fukui et al. 2015) and N 159 East (Saigo et al. 2017) in the Large Magellanic Cloud.

In theoretical studies, hydrodynamical numerical simulation of a cloud–cloud collision was first carried out by Stone (1970a, 1970b). Habe and Ohta (1992) made numerical simulations of head-on collisions for two clouds of different sizes. They showed that gravitationally unstable cores are created by compression between the two clouds (see also Anathpindika 2010; Takahira et al. 2014, 2018; Shima et al. 2018). Balfour et al. (2015) and Balfour, Whitworth, and Hubber (2017) also presented numerical simulation of head-on and non-head-on cloud–cloud collisions. In addition, magneto-hydrodynamical (MHD) simulations of giant molecular cloud (GMC) collision were carried out by several authors (e.g., Wu et al. 2015, 2017a, 2017b; Christie et al. 2017; Bisbas et al. 2017; Li et al. 2018). They showed that GMC collisions enhanced star formation rate and efficiency (Wu et al. 2017b). Inoue and Fukui (2013) studied the interface layer of colliding clouds in three-dimensional MHD simulations, and showed formation of massive molecular cores, which likely lead to forming high-mass protostars gaining high mass accretion rate helped by amplified turbulence and magnetic fields via supersonic collision (see also Inoue et al. 2018). Kobayashi et al. (2017, 2018) discussed the evolution of GMC mass functions including cloud–cloud collisions. In global scale numerical simulations, cloud–cloud collisions are an important mechanism of star formation in the Galaxy (e.g., Tan 2000; Tasker & Tan 2009; Fujimoto et al. 2014a, 2014b; Dobbs et al. 2015; Li 2017).

These observational and theoretical studies suggest that a cloud–cloud collision is a promising mechanism for massive star formation, whereas there still remain large numbers of massive star-forming regions where cloud–cloud collisions have not been investigated well.

1.3 High-mass star-forming region W 33

W 33 is a high-mass star-forming region first cataloged by the 1390 MHz radio continuum survey (Westerhout 1958), extending for 10 pc × 10 pc centered on (l, b) ∼ (⁠|${12.\!\!^{\circ}8}$|⁠, |${-0.\!\!^{\circ}2}$|⁠). The parallactic distance of W 33 was measured as 2.4 kpc based on water maser observations by Immer et al. (2013), indicating that W 33 is located in the Scutum spiral arm in the Milky Way. Figure 1 shows a three-color composite image of the Spitzer space telescope observations (GLIMPSE: Benjamin et al. 2003; MIPSGAL: Carey et al. 2009), where blue, green, and red correspond to the 3.6 μm, 8 μm, and 24 μm emissions, respectively. The contours in figure 1 indicate the MAGPIS 90 cm radio continuum emission (Helfand et al. 2006).

W 33 harbors many star-forming clumps, OB stars, and H ii regions. There are six dust clumps (W 33 Main, W 33 A, W 33 B, W 33 Main1, W 33 A1, and W 33 B1) as shown in the pink contours in figure 1, obtained by the Atacama Pathfinder Experiment (APEX) Telescope Large Area Survey of the GALaxy (ATLASGAL) 870 μm survey (Schuller et al. 2009; Contreras et al. 2013; Urquhart et al. 2014). The sizes and masses of the dust clumps are listed, together with their physical properties, in table 1. These six dust clumps are suggested to be on different evolutional stages, indicated by their spectral properties and associations with radio continuum sources (Immer et al. 2014). W 33 Main includes a (compact) H ii region detected by radio continuum observations (Haschick & Ho 1983; Ho et al. 1986), while W 33 A and W 33 B harbor hot cores by the chemical line survey. W 33 Main1, W 33 A1, and W 33 B1 are identified as high-mass protostellar objects (Immer et al. 2014). Radio observations by Haschick and Ho (1983) revealed the presence of an obscured (proto)cluster which includes a number of high-mass stars with spectral types from O7.5 to B1.5 by assuming a distance of 4 kpc. Associations of the maser sources also support the nature of the young massive stellar objects, i.e., water and methanol masers in W 33 Main, W 33 A, and W 33 B, and OH masers in W 33 A and W 33 B (Immer et al. 2013; Menten et al. 1986; Caswell 1998; Colom et al. 2015). The Submillimeter Array (SMA) high-resolution observations showed that W 33 Main harbors multiple ultra-compact H ii regions and three high-density clumps (W 33 Main-A, W 33 Main-B1, W 33 Main-B2) embedded in a dense gas envelope detected in C₂H (Jiang et al. 2015). W 33 A has been studied based on observations at various wavelengths (Gibb et al. 2000; Roueff et al. 2006; Davies et al. 2010; de Wit et al. 2007, 2010). Galván-Madrid et al. (2010) detected high-velocity gas associated with an outflow in the CO J = 2–1 transition by SMA 230 GHz band observations in W 33 A. They suggested that star formation activity in W 33 A was triggered by filamentary convergent gas flows from two different velocity components. Recently, Maud et al. (2017) carried out Atacama Large Millimeter/Submillimeter Array (ALMA) observations in Band 6 and Band 7 with 500 au-scale resolution toward W 33 A, finding spiral and filamentary structures around the central massive young stellar object in W 33 A.

H ii regions in W 33 are shown in the radio continuum emissions shown in figure 1 as white contours, where the names of the H ii regions cataloged by the radio recombination line survey (G012.745−00.153: Downes et al. 1980; Lockman 1989) and Wide-Field Infrared Survey Explorer (WISE) satellite (Anderson et al. 2014, 2015), G012.692−00.251, G012.820−00.238, G012.884−00.237, and G012.907−00.277, are labeled.

Messineo et al. (2015) identified many high-mass stars in W 33 based on the near-infrared K-band observations using Spectrograph for INtegral Field Observations in the Near Infrared (SINFONI) on the Very Large Telescope (VLT). The positions of the OB stars identified by Messineo et al. (2015) are indicated as circles on the three-color image of the Spitzer observations in figure 1, where the circles colored in orange indicate the O-type stars, while those in white indicate early B-type stars. Two O4-6 are distributed in G012.745−00.153, and Messineo et al. (2015) discussed that the two O stars have evolved to be (super-)giant. The authors discussed the ages of these two O stars as less than 6 Myr based on a stellar model by Ekström et al. (2012) or 2–4 Myr by comparing the O4-6 supergiants in the Arches cluster in the Galactic center (Martins et al. 2008), which show similar K-band spectra to the present two O stars in W 33.

CO rotational transition line observations at millimeter wavelengths covering the entire W 33 region were carried out using the Five College Radio Observatory (FCRAO) 14 m telescope and Millimeter Wave Observatory (MWO) 5 m telescope by Goldsmith and Mao (1983), revealing that W 33 has complicated velocity structures, whereas detailed spatial and velocity distributions of molecular gas remained unclear. Observations of the H₂CO absorption and radio recombination lines indicate that W 33 Main has a velocity of ∼35 km s⁻¹, while W 33 B’s is ∼60 km s⁻¹ (Gardner & Whiteoak 1972; Bieging et al. 1978). Immer et al. (2013) detected water maser emission at a velocity of ∼35 km s⁻¹ toward W 33 Main and W 33 A, and ∼60 km s⁻¹ toward W 33 B, finding that these two velocity components are at the same annual parallaxial distance of 2.4 kpc. The authors concluded that W 33 is a single star-forming region in spite of the large velocity difference between the two velocity components, ∼25 km s⁻¹, although the origin of the two velocity components is still ambiguous.

In order to reveal the detailed spatial and velocity distributions of the molecular gas in W 33, we carried out high-resolution observations in the J = 1–0 transition of ¹²CO, ¹³CO, and C¹⁸O using the NANTEN2 and Nobeyama 45 m telescopes. This paper is organized as follows: section 2 describes the observations, section 3 gives the observational results and comparison with archive datasets, and section 4 discusses a cloud–cloud collision scenario for W 33. Section 5 concludes the paper.

2 Dataset

2.1 NANTEN2 ¹²CO J = 1–0 observations

The NANTEN2 4 m millimeter/submillimeter telescope of Nagoya University situated in Chile was used to observe a large area of W 33 in ¹²CO J = 1–0 emission with the on-the-fly (OTF) mode from 2012 May to 2012 December. The half-power beam width (HPBW) is |${2{^{\prime}_{.}}7}$| at 115 GHz. This corresponds to 1.9 pc at the distance of 2.4 kpc. A 4 K cooled superconductor–insulator–superconductor (SIS) mixer receiver provided a typical system temperature of ∼250 K in double side band (DSB), and a 16384-channel digital spectrometer with a bandwidth and resolution of 1 GHz and 61 kHz respectively, which corresponds to 2600 km s⁻¹ and 0.16 km s⁻¹, respectively, at 115 GHz. We smoothed the obtained data to a velocity resolution of 1 km s⁻¹ and angular resolution of |${200^{\prime \prime }}$|⁠. The pointing accuracy was confirmed to be better than |${15^{\prime \prime }}$| with daily observations of the Sun and IRC + 10216. We used the chopper wheel method to calibrate the antenna temperature |$T_{\rm a}^*$| (Penzias & Burrus 1973; Ulich & Haas 1976; Kutner & Ulich 1981). The absolute intensity fluctuation was calibrated by daily observations of IRAS 16293−2422 [|$\alpha _{\rm J2000.0} = {16^{\rm h}32^{\rm m}23{^{\rm s}_{.}}3} , \delta _{\rm J2000.0} ={-24^{\circ}28^{\prime }39{^{\prime\prime}_{.}}2}$|], and the intensity scale was converted into the T_mb scale by assuming its peak to be T_mb = 18 K (Ridge et al. 2006). The intensity uncertainty of NANTEN2 datasets is <20%. The typical root-mean-square (rms) noise level is ∼0.6 K with |${200^{\prime \prime }}$| smoothing data and a velocity resolution of 1.0 km s⁻¹.

2.2 Nobeyama 45 m telescope ¹²CO J = 1–0, ¹³CO J = 1–0, C¹⁸O J = 1–0 observations

Detailed CO J = 1–0 data around W 33 were obtained by using the Nobeyama 45 m telescope at the Nobeyama Radio Observatory (NRO). The HPBW is |${14^{\prime \prime }}$| at 115 GHz, and |${15^{\prime \prime }}$| at 110 GHz. This corresponds to 0.2 pc at the distance of 2.4 kpc. We simultaneously observed in ¹²CO, ¹³CO, and C¹⁸O as part of the FUGIN (FOREST Unbiased Galactic plane Imaging survey with the Nobeyama 45 m telescope; Minamidani et al. 2015; Umemoto et al. 2017) legacy survey with the OTF mode (Sawada et al. 2008) from 2014 March to 2015 May. FOREST (FOur-beam REceiver System on the 45 m Telescope) is a four-beam, dual-polarization, two-sideband (2SB) receiver providing a typical system temperature of ∼150 K in ¹³CO J = 1–0 and ∼250 K in ¹²CO J = 1–0 (Minamidani et al. 2016). We used SAM45 (Spectral Analysis Machine for the 45 m telescope: Kuno et al. 2011), an FX-type digital spectrometer the same as the ALMA ACA Correlator (Kamazaki et al. 2012). It has 4096 channels with a bandwidth and resolution of 1 GHz and 244.14 kHz respectively, which corresponds to 2600 km s⁻¹ and 1.3 km s⁻¹, respectively, at 115 GHz. We smoothed the obtained data to the angular resolution of |${30^{\prime \prime }}$|⁠. The typical pointing accuracy was confirmed to be better than |${3^{\prime \prime }}$| observing SiO maser sources, such as V VX Sgr [|$\alpha _{\rm B1950} ={18^{\rm h}05^{\rm m}02{^{\rm s}_{.}}959} , \delta _{\rm B1950} = {-22^{\circ}13^{\prime }55{^{\prime\prime}_{.}}58}$|], for the observing run toward W 33 every hour with the 40 GHz HEMT receiver named H40. The intensity variation was calibrated with daily observations of W51D [|$\alpha _{\rm B1950} = {19^{\rm h}21^{\rm m}22{^{\rm s}_{.}}2} , \delta _{\rm B1950} = {14^{\circ}25^{\prime }17{^{\prime\prime}_{.}}0}$|]. We used the chopper wheel method to convert the antenna temperature |$T_{\rm a}^*$| (Penzias & Burrus 1973; Ulich & Haas 1976; Kutner & Ulich 1981). The data in the antenna temperature (⁠|$T_{\rm a}^*$|⁠) scale was converted into main beam temperature (T_mb) by |$T_{\rm mb} = T_{\rm a}^* / \eta _{\rm mb}$|⁠, with a main beam efficiency (η_mb) of 0.43 for ¹²CO and 0.45 for ¹³CO and C¹⁸O (Minamidani et al. 2016; Umemoto et al. 2017). The typical rms noise levels after intensity calibration (T_mb scale) are ∼0.5 K, ∼0.2 K, and ∼0.2 K for ¹²CO, ¹³CO, and C¹⁸O J = 1–0, respectively. The intensity variation of ¹²CO, ¹³CO, and C¹⁸O J = 1–0 are 10%–20%, 10%, and 10%, respectively. The FUGIN project overview paper gives more detailed information (Umemoto et al. 2017). We summarize the observational parameters of the NANTEN2 and Nobeyama 45 m datasets in table 2.

2.3 Archive dataset

We use the ¹²CO J = 3–2 archival data of the CO High Resolution Survey (COHRS) obtained with JCMT (James Clark Maxwell Telescope: Dempsey et al. 2013). The spatial and velocity resolutions are 1.0 km s⁻¹ and |${16^{\prime \prime }}$|⁠, respectively. We smoothed the data to an angular resolution of |${30^{\prime \prime }}$|⁠. The data in the antenna temperature (⁠|$T_{\rm a}^*$|⁠) scale was converted into main beam temperature (T_mb) with the equation |$T_{\rm mb} = T_{\rm a}^* / \eta _{\rm mb}$|⁠, where the main beam efficiency (η_mb) of 0.61 was adopted by planet observations with an uncertainty of about 10%–15% (Dempsey et al. 2013; Buckle et al. 2009). The typical rms noise in the T_mb scale was ∼0.2 K for ¹²CO J = 3–2.

We use the following datasets to compare with the CO data: near- and mid-infrared data from the Spitzer space telescope (GLIMPSE in 3.6 μm and 8.0 μm: Benjamin et al. 2003; MIPSGAL in 24 μm: Carey et al. 2009); the 20 cm and 90 cm free–free radio continuum data from MAGPIS (A Multi-Array Galactic Plane Imaging Survey: Helfand et al. 2006) observed with the Very Large Array (VLA); and the H i 21 cm emission data from SGPS II (Southern Galactic Plane Survey: McClure-Griffiths et al. 2005) observed with the Australia Telescope Compact Array (ATCA) and the Parkes Radio Telescope. The angular resolution of the H i data is |$\sim {3{^{\prime}_{.}}3}$|⁠, and their spectrum resolution is 0.8 km s⁻¹.

3 Results

3.1 Distribution and properties of molecular gas

Figure 2a shows the longitude–velocity diagram of the NANTEN2 ¹²CO J = 1–0 data for a large area including W 33. There are four distinct velocity components toward W 33, as depicted by arrows. Goldsmith and Mao (1983) mentioned with CO observations that association of the 5–25 km s⁻¹ velocity component indicated with a black arrow with W 33 is not clear, as it is likely to be located in the Sagittarius arms (Reid et al. 2016). In this study, we therefore focus on the velocity components at 35 km s⁻¹, 45 km s⁻¹, and 58 km s⁻¹ as candidate clouds associated with W 33. We hereafter refer to these three clouds as “the 35 km s⁻¹ cloud,” “the 45 km s⁻¹ cloud,” and “the 58 km s⁻¹ cloud.” These three clouds appear to be connected with each other in the velocity space (figure 2a). The 35 km s⁻¹ cloud is continuously distributed along the Galactic longitude for the present region, while the CO emissions in the other two clouds at higher velocities are especially enhanced toward W 33 (dashed box in figure 2a). Figure 2b shows the longitude–velocity diagram of the three velocity clouds of W 33 using the FUGIN ¹²CO J = 1–0 data, which shows the detailed velocity distribution of the gas at high spatial resolution. The three velocity clouds can be separately identified at this spatial scale in figure 2b.

Figure 3 shows the spatial distributions of the three clouds in the four CO lines. The left, center, and right columns show the 30 km s⁻¹, 45 km s⁻¹, and 58 km s⁻¹ clouds, respectively, where the contours indicate the MAGPIS 90 cm radio continuum emission, and the positions of the dust clumps and OB stars identified by Messineo et al. (2015) are depicted by crosses and circles, respectively. We also present the velocity channel maps of the four CO transitions obtained with Nobeyama 45 m and JCMT in figures 14–17 in the Appendix as supplements. The CO emission in the 35 km s⁻¹ cloud is enhanced at the corresponding region of W 33 in all the transitions. W 33 Main shows the brightest CO emission in this region.

We found counterparts of W 33 A, W 33 Main, and W 33 B1 in the 35 km s⁻¹ cloud in the C¹⁸O emission (figure 3g). Figure 4 shows the C¹⁸O distributions of each dust clump (the C¹⁸O molecular clump properties are discussed in subsection 3.2.) Among these sources, W 33 Main and W 33 A are associated with CO molecular outflows, as discussed in subsection 3.6. The 35 km s⁻¹ cloud also shows morphological anti-correlations with the radio continuum emissions in H ii regions G012.745−00.153 and G012.820−00.238. The CO emission in G012.745−00.153 is enhanced at the eastern (l, b) = (12.78, −0.18) and southern rim (l, b) = (12.74, −0.18) of the H ii region, showing a steep intensity gradient in the ¹²CO emissions, while G012.820−00.238 is surrounded by molecular gas, especially in the ¹³CO and C¹⁸O emissions. W 33 Main is sandwiched by these two H ii regions. The 45 km s⁻¹ cloud has diffuse CO emission extended over the present region. The compact emissions at W 33 Main and W 33 A correspond to the wing features of the outflows (see subsection 3.6 for details), and are thus not related to the 45 km s⁻¹ cloud. Molecular gas in the 58 km s⁻¹ cloud is separated into the northern and southern components relative to W 33, and the central part corresponding to W 33 is weak in the CO emission. There are several clumpy structures embedded at the northern rim of the southern component, which are clearly seen in the ¹²CO emissions, and these clumps show spatial correlations with radio continuum emissions from the H ii regions G012.745−00.153 and G012.692−00.251 as well as W 33 B. In the C¹⁸O map in figure 3i and figure 4d, W 33 B is associated with the strong CO peak. There are several other clumpy molecular structures at the interspace between the northern and southern components of the 58 km s⁻¹ cloud, forming an arc-like molecular structure which looks to surround W 33. The size of the arc-like structure is roughly estimated to be ∼7 pc. On the other hand, clear associations of molecular clumps with W 33 A1 and W 33 Main1 are not recognized.

We derived the column densities and masses of the three velocity clouds using the ¹²CO integrated intensity maps shown in figures 3a–3c, where we defined the individual clouds by drawing contours at 5 σ noise levels in the integrated intensity of 8 K km s⁻¹ for the velocity interval of 10 km s⁻¹. By assuming an X(CO) factor of 2 × 10²⁰ (K km s⁻¹)⁻¹ cm⁻² (Strong et al. 1988), we estimated the mean column densities of the 35 km s⁻¹, 45 km s⁻¹, and 58 km s⁻¹ clouds as 1.7 × 10²² cm⁻², 1.7 × 10²² cm⁻², and 6.2 × 10²¹ cm⁻², respectively, with the total molecular masses derived as 1.1 × 10⁵ M_⊙, 1.0 × 10⁵ M_⊙, and 3.8 × 10⁴ M_⊙. The uncertainty of the mass estimation using the X factor is about ±30% (Bolatto et al. 2013). Lin et al. (2016) derived the mean column densities as 2.5 × 10²² cm⁻² using the infrared dust emission data obtained by Herschel, which is consistent with our estimate.

3.2 C¹⁸O molecular clump properties

We define C¹⁸O molecular clumps using the following procedures in order to investigate the physical properties of dense molecular gas belonging to the 35 km s⁻¹ and 58 km s⁻¹ clouds corresponding to the dust clumps.

1. Search for a peak integrated intensity toward the six dust clumps.

2. Define a clump boundary as the half level of its peak integrated intensity.

3. If the area enclosed by the boundary has multiple peaks, define the boundary as the contour of the 75% level of its peak integrated intensity.

We identified four molecular clumps corresponding to the dust clumps of W 33 Main, W 33 A, and W 33 B1 associated with the 35 km s⁻¹ cloud (figures 4a, 4b, and 4c), and W 33 B associated with the 58 km s⁻¹ cloud (figure 4d), while we could not define the boundary of the molecular clump toward W 33 Main1 and W 33 A1 because W 33 A1 does not have an intensity peak toward the dust clump and W 33 Main1 cannot be separated from the extended feature (figures 4e and 4f).

Then, we derived the physical parameters assuming local thermal equilibrium (LTE) using the following procedures (Wilson et al. 2009).

1. Derive the excitation temperature (T_ex) assuming that the ¹²CO J = 1–0 transition line is optically thick and T_ex is uniform throughout the molecular clump by following the equation from the ¹²CO peak intensity (⁠|$T_{\rm mb}(\rm ^{12}CO peak)$|⁠) at the peak position of the clump:

   \begin{eqnarray} T_{\rm ex} &=& 5.5 \bigg / \ln \left[1+ {5.5 \over T_{\rm mb}(\rm ^{12}CO peak) + 0.82 }\right]. \end{eqnarray}

   (1)
2. Estimate the optical depth of the C¹⁸O emission (τ₁₈) at each pixel and velocity channel from the C¹⁸O brightness temperature [T_mb(v)] at velocity v:

   \begin{eqnarray} \tau _{18} (v) &=& -\ln \left\lbrace 1-{T_{\rm mb}(v) \over 5.3} \left[ {1 \over \exp ({5.3 \over T_{\rm ex}})-1}-0.17 \right]^{-1} \right\rbrace . \nonumber\\ \end{eqnarray}

   (2)
3. Calculate the C¹⁸O column density [N(C¹⁸O)] at each pixel, summing the quantities of all v channels whose resolution is 1 km s⁻¹:

   \begin{eqnarray} N ({\rm C^{18}O}) &=& 2.4 \times 10^{14} \times \sum _v {T_{{\rm ex}} \tau _{18} (v) \Delta v \over 1-\exp \left(-{5.3 \over T_{\rm ex}} \right) } . \end{eqnarray}

   (3)
4. N(C¹⁸O)) is converted into H₂ column density (⁠|$N(\rm H_2)$|⁠) assuming the following conversion formula derived from the Ophiuchus molecular cloud (Frerking et al. 1982):

   \begin{eqnarray} N ({\rm H_2}) = \left[ {N (\rm C^{18}O) \over 1.7 \times 10^{14}} + 3.9 \right] \times 10^{21} . \end{eqnarray}

   (4)
5. Estimate the size (r) of each clump assuming that clumps are spherical shapes, using:

   \begin{eqnarray} r = \sqrt{{S \over \pi }}=\sqrt{D^2\Omega \over \pi }, \end{eqnarray}

   (5)

   where S is the area of the clump enclosed by the boundary contour, D is the distance to W 33, and Ω is the solid angle of the clump.
6. Calculate the mass of each clump:

   \begin{eqnarray} M = \mu m_{\rm H} D^2 \Omega \sum N(\rm H_2), \end{eqnarray}

   (6)

   where μ is the mean molecular weight 2.8, m_H is the proton mass, and the summation is performed over the each clump within the boundary.
7. The averaged number density of hydrogen molecules is calculated as

   \begin{eqnarray} n({\rm H_2}) ={ 3 M \over 4 \pi r^3 \mu m_{\rm H} } \end{eqnarray}

   (7)

   by assuming that the clumps have a uniform density.
8. The virial mass is estimated by

   \begin{eqnarray} M_{\rm vir} ={ 5 r \Delta v^2 \over 8 (\ln 2) G }=210 \left({r \over [\rm pc]}\right) \left({\Delta V_{\rm comp} \over [\rm km{\,\,}s^{-1}]}\right)^2 \end{eqnarray}

   (8)

   from the virial theorem using the clump size (r) and the full width at half maximum (FWHM) of the composite spectrum ΔV_comp in the clump estimated by a single Gaussian fitting.

Table 3 presents the physical properties of the molecular clumps. The typical values of the peak column density, mean column density, mass, size, hydrogen molecule number density, and virial mass are N(H₂)_peak ∼ 10²³ cm⁻², N(H₂)_mean ∼ 10²² cm⁻², M_clump ∼ 10³–10⁴ M_⊙, r ∼0.4–1 pc, |$n(\rm H_2) \sim 10^{4}$|–10⁵ cm⁻³, and M_vir ∼ 10²–10³ M_⊙. The peak and mean column densities are roughly consistent with those calculated from the dust continuum observations (Immer et al. 2014).

3.3 CO J = 3–2 / J = 1–0 intensity ratio

Figure 5 shows distributions of the ¹²CO J = 3–2/J = 1–0 intensity ratio (⁠|$R_{\rm 3\mbox{-}2/1\mbox{-}0}$|⁠) using |${50^{\prime \prime }}$| smoothing data (rms ∼0.26 K) in the 35 km s⁻¹, 45 km s⁻¹, and 58 km s⁻¹ clouds, respectively. We adopted the clipping level as 8 σ. Line intensity ratios between different J levels of CO provide useful diagnostics to investigate the physical association of molecular gas with high-mass stars, as these depend on the gas kinematic temperature (T_k) and number density of hydrogen [n(H₂)] following the Large Velocity Gradient (LVG) model (e.g., Goldreich & Kwan 1974). Figure 5a shows that the 35 km s⁻¹ cloud has a high |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of >1.0 around the central part of W 33, which includes W 33 Main, W 33 Main1, and W 33 B1. The 45 km s⁻¹ cloud shows lower |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of 0.4–0.5 except for the corresponding parts of W 33 Main and W 33 A, where |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| is locally elevated up to 0.8–1.2. As discussed in the next subsection, the high |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| in W 33 Main and W 33 A in the 45 km s⁻¹ cloud are due to the outflows emitted from the protostars embedded in the 35 km s⁻¹ cloud. The 58 km s⁻¹ cloud has high |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of 1.2 in the arc-like structure which seems to surround W 33, while the gas outside the arc-like structure shows low |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of ∼0.4. We also present the velocity channel map of |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| in figure 18 in the Appendix. The cause of the high |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of the gas in the 35 km s⁻¹ and 58 km s⁻¹ clouds can be understood as heating by the high-mass stars in W 33. The high |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of the molecular gas in the 35 km s⁻¹ and 58 km s⁻¹ clouds can be understood as high kinematic temperature heated by the high-mass stars in W 33. Therefore, the 35 km s⁻¹ and 58 km s⁻¹ clouds indicate physical associations of W 33, whereas association of the 45 km s⁻¹ is elusive.

3.4 Comparisons of the 35 km s⁻¹ and 58 km s⁻¹ clouds with infrared data

We focus here on the 35 km s⁻¹ and 58 km s⁻¹ clouds, as these are most likely associated with W 33. Figure 6 demonstrates comparisons between the ¹²CO J = 3–2 and infrared emissions, where the Spitzer 8 μm and 24 μm emissions are shown in grayscale in panels (a) and (b), respectively. The 8 μm emissions trace thermal emission from hot dust plus emission from polycyclic aromatic hydrocarbons (PAH), which are distributed in the photo-dissociation region (PDR; e.g., Churchwell et al. 2004). The 24 μm emission also traces warm dust grains heated by bright high-mass stars in the H ii region (e.g., Carey et al. 2009, Deharveng et al. 2010). The red and blue contours show the 35 km s⁻¹ and 58 km s⁻¹ clouds, respectively. The arc-like structure in the 58 km s⁻¹ cloud surrounds the strong emission part of the 35 km s⁻¹ cloud, showing complementary distributions between these two. The 8 μm emission, which is bright at the southeastern rim of the H ii region G012.745−00.153 (l, b) ∼ (12.75, −0.18), coincides with the steep intensity gradient of the CO emission in the 35 km s⁻¹ cloud. This lends more credence to association between the 35 km s⁻¹ cloud and G012.745−00.153.

The 24 μm emission shown in figure 6b is attributed to the warm dust grains and thus can be used to probe the region where the heating by the high-mass stars in W 33 is efficient. The 24 μm emission in W 33 is enhanced at the dust clumps and the H ii regions. In particular, G012.745−00.153 and G012.820−00.238, as well as W 33 Main, show bright 24 μm emissions. The arc-like structure in the 58 km s⁻¹ cloud overall traces the outline of the 24 μm distribution. This suggests that the high |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of the gas in the arc-like structure shown in figure 5c is due to interaction with the star-forming regions, which are spatially correlated with the arc-like structure.

The morphological correlations of the 58 km s⁻¹ cloud with infrared emissions, as well as the high |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of the gas in the 58 km s⁻¹ cloud, strongly suggest that not only the 35 km s⁻¹ cloud but also the 58 km s⁻¹ cloud are physically associated with W 33 despite the large velocity separation between the two clouds, ∼23 km s⁻¹. This is consistent with the previous studies by Gardner and Whiteoak (1972), Bieging, Pankonin, and Smith (1978), and Immer et al. (2013) based on observations of H₂CO absorption, radio recombination line, and water masers; Immer et al. (2013) discussed that the two velocity components are at the same annual parallax distance of 2.4 kpc.

On the other hand, association of the 45 km s⁻¹ cloud at the intermediate velocity range between the 35 km s⁻¹ and 58 km s⁻¹ clouds with W 33 is not evident in the present dataset.

3.5 Comparison with the H i 21 cm line

We analyzed the H i 21 cm line data in W 33 using the archival data obtained in the Southern Galactic Plane Survey (SGPS) project (McClure-Griffiths et al. 2005). Figure 7 shows the integrated intensity distributions of the three velocity clouds of the H i data. The contours show the MAGPIS 20 cm radio continuum emission, where the original image was spatially smoothed to have the same resolution as that of the SGPS H i data. In figure 7a the H i emission in the 35 km s⁻¹ cloud shows strong intensity depression at W 33 Main where the 20 cm continuum emission is strong. The depression can be confirmed in the H i 21 cm spectra shown in figure 8a. The H i profile has negative intensities at the velocity range of the 35 km s⁻¹ cloud, as indicated by the shading. This indicates absorption against the background continuum source, implying that W 33 Main is distributed at the inside or at the rear side of the 35 km s⁻¹ cloud. The H i maps of the 45 km s⁻¹ and 58 km s⁻¹ clouds in figures 7b and 7c show weak intensity depressions toward G012.745−00.153 and W 33 B, whereas it is not clear toward W 33 Main. The H i spectra toward G012.745−00.153 shown in figures 8b and 8c show absorption features at the velocity ranges where the ¹²CO emission appears, i.e., 40–50 km s⁻¹ in figure 8b and at 50–60 km s⁻¹ in figure 8c. We suggest that these absorptions are unrelated foreground components overlapping on the line of sight.

3.6 Molecular outflows

We identified bipolar outflows in the CO emissions toward W 33 Main and W 33 A. Figures 9 and 10 show the CO spectra and spatial distributions of the outflow lobes in W 33 Main and W 33 A, respectively. The outflows in both W 33 Main and W 33 A have total velocity widths as large as 40–50 km s⁻¹. The systemic velocities of these outflows are estimated to be ∼35 km s⁻¹ using the optically thin C¹⁸O J = 1–0 emission. W 33 A has been reported to be associated with outflows (Davies et al. 2010; de Wit et al. 2010; Galván-Madrid et al. 2010). As the blue and red lobes in figures 9 and 10 are not fully resolved in the present CO data due to the limitation of the spatial resolution, we estimate here the lengths of the outflow lobes to be 0.5 pc, defined by contours at half of the peak intensity level. If we assume an inclination angle of |${45^{\circ}}$|⁠, the dynamical timescales (t_dyn) of W 33 Main and W 33 A are calculated from the maximum velocity (V_max ∼ 23 km s⁻¹) and size (r) to be |$0.5/23 \times \sqrt{2} \sim 3 \times 10^{4}\:$|yr, which is consistent with previous JCMT ¹²CO and ¹³CO J = 3–2 observation results (Maud et al. 2015).

4 Discussion

The results of the present observations and analyses are summarized as follows:

1. We identified three molecular clouds toward W 33 at ∼35 km s⁻¹, ∼45 km s⁻¹, and ∼58 km s⁻¹ using the NANTEN2, FUGIN CO J = 1–0, and JCMT ¹²CO J = 3–2 data. The total molecular masses of each of the three clouds are derived as 1.1 × 10⁵ M_⊙, 1.0 × 10⁵ M_⊙, and 3.8 × 10⁴ M_⊙, respectively.

2. The 35 km s⁻¹ cloud is spatially coincident with W 33. Our CO data revealed a spatial correlation of the 35 km s⁻¹ cloud with the dust clumps of W 33 Main, W 33 A, W 33 B1, and the H ii regions G012.745−00.153 and G012.820−00.238, having |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of higher than 1.0. A strong absorption feature in the H i 21 cm line is seen in the velocity range of the 35 km s⁻¹ cloud. Therefore, W 33 Main is located inside or behind the 35 km s⁻¹ cloud.

3. The 45 km s⁻¹ cloud shows diffuse CO emission extended for the present region, and its association with W 33 is not clear in terms of spatial correlation.

4. In the 58 km s⁻¹ cloud the present CO dataset revealed an arc-like structure having a size of ∼7 pc. It shows complementary distributions with the 35 km s⁻¹ cloud and the infrared images along the line of sight, having |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of higher than 1.0. These observational properties suggest association of the arc-like structure with the dust clump W 33 B and the H ii regions G012.745−00.153 and G012.692−00.251 distributed at the southern part of W 33.

5. We identified two bipolar molecular outflows toward W 33 Main and W 33 A. The full velocity widths of the outflows are as large as 40–50 km s⁻¹. The dynamical timescales of the outflows can be estimated to be ∼3 × 10⁴ yr.

In this section, we discuss the origin of the observed properties of the multiple velocity clouds in W 33 over 23 km s⁻¹ and these relationships with the high-mass star formation in W 33.

4.1 Gravitational binding of the multiple velocity clouds

Based on the proper motion measurements, Immer et al. (2014) discussed whether W 33 A and W 33 B are gravitationally bound to W 33 Main, and the authors concluded that W 33 A and W 33 B are not gravitationally bound to W 33 Main, as the total speeds of W 33 A and W 33 B are larger than the derived escape velocity.

In addition to their calculations, we test here the dynamical binding of the 35 km s⁻¹ and 58 km s⁻¹ clouds, as these contain the high-mass star-forming regions other than W 33 Main, W 33 A, and W 33 B. If we assume that the two clouds are separated by 10 pc in space, the same order as the cloud size as a rough order-of-magnitude estimation, and by 23 km s⁻¹ in velocity, the total mass required to gravitationally bind these two clouds can be calculated as |$M = {rv^2 \over 2G}= 6 \times 10^{5}{\,\,}M_{\odot }$|⁠. This figure is a factor of four larger than the total molecular mass associated with W 33, 1.5 × 10⁵ M_⊙, calculated in subsection 3.1, indicating that the coexistence of the two velocity clouds in W 33 cannot be interpreted as a gravitationally bound system.

4.2 Expanding motion driven by feedback from high-mass stars

Another idea to interpret the multiple velocity clouds is expanding motion driven by feedback from high-mass stars. If we assume a spherical expansion of gas interacting with feedback, it would display a ring-like velocity distribution in both a spatial map and a position–velocity diagram (e.g., see figure 8 of Torii et al. 2015). The expanding velocity (v_exp) of an H ii region is limited by the sound speed of the ionized gas confined in the H ii region, which corresponds to v_exp ∼ 12 km s⁻¹ with an electron temperature of 10000 K (e.g., Ward-Thompson & Whitworth 2011). For the present case in W 33, the arc-like structure in the 58 km s⁻¹ cloud and velocity separation of 23 km s⁻¹ look consistent with this assumption. However, the longitude–velocity diagram presented in figure 2 shows that both the 35 km s⁻¹ and 58 km s⁻¹ clouds are separated, having uniform velocities for a large area including W 33. These velocity distributions in the position–velocity diagram are inconsistent with the prediction from the expansion assumption. Dale et al. (2013) discussed that the contribution of stellar wind to the expansion of a neutral medium is relatively minor compared with the expansion of an H ii region at gas density ∼10⁴ cm⁻¹. We thus conclude that the co-location of the multiple velocity clouds in W 33 over 23 km s⁻¹ was not formed by the feedback from the high-mass stars in W 33.

4.3 The cloud–cloud collision model

We postulate here a cloud–cloud collision scenario as an alternative idea to interpret the multiple velocity clouds associated with W 33. In a collision between two clouds with different sizes (Habe & Ohta 1992), a smaller cloud drives into a larger cloud, forming a cavity in the larger cloud, and a dense gas layer is formed at the interface of the collision, which corresponds to the bottom of the cavity, by strong compression of gas. High-mass stars are formed in this dense gas layer. In addition, a thin layer with turbulent gas is formed at the interface between the larger cloud and the dense layer, which has intermediate velocities between the larger cloud and the dense layer. This thin, turbulent layer is observed as “broad bridge features” in a position–velocity diagram, which is the emission between two velocity peaks with intermediate intensities (e.g., Haworth et al. 2015a, 2015b; Torii et al. 2017a). The turbulent gas is replenished as long as the collision continues.

Another observational signature of cloud–cloud collisions is “complementary distribution between two velocity clouds” (Fukui et al. 2017a; Torii et al. 2017a). For an observer viewing angle parallel to the collisional axis so that the two clouds are spatially coincident, the observer can see an anti-correlated or complementary distribution between the two clouds separated in velocity, as the larger cloud with a cavity displays a ring-like gas distribution on the sky. Fukui et al. (2017a) pointed out that, if the viewing angle has an inclination relative to the colliding axis, the complementary distribution has a spatial offset, which is proportional to the projected travel distance of the collision.

The gas distribution in the position–velocity diagram including the broad bridge feature depends on various parameters such as the cloud shapes, the density contrast of the clouds, the initial relative velocity between the clouds, etc. In figures 11 and 12, we show schematics of two extreme cases of a cloud–cloud collision and cartoons of the corresponding position–velocity diagrams, based on the discussions by Haworth et al. (2015a, 2015b), in which the authors post-processed the model data of cloud–cloud collisions calculated by Takahira, Tasker, and Habe (2014).

1. If the smaller cloud has a gas density much smaller than the larger cloud, which is the case shown in figure 11, the moving velocity of the dense layer is immediately decreased owing to the momentum conservation, and the collision halts in the middle of the larger cloud. In this case, the entirety of the smaller cloud undergoes collision quite quickly, streaming into the dense layer. The cartoon position–velocity diagram in the middle of the collision in figure 11 shows two velocity peaks separated by broad bridge features having intermediate intensities, whereas there is only one velocity peak after the collision is stopped.

2. On the other hand, if the smaller cloud has a density much higher than the larger cloud, deceleration of the dense layer is quite small, and the smaller cloud punches the larger cloud within a short time (figure 12). In this case, although the gas distribution in the position–velocity diagram during the collision is similar to that in case 1 in the zeroth-order approximation, the gas configuration after the collision finishes is quite different from case 1. In the position–velocity diagram at the final stage, as only some small fraction of the smaller cloud undergoes the collision, the two colliding clouds are separated in velocity without broad bridge features, and the larger cloud shows an intensity depression at the positions corresponding to the smaller cloud. In this case, spatially complementary distribution is an important diagnostic to investigate cloud–cloud collisions.

4.4 Cloud–cloud collision in W 33

We discuss here an application of the cloud–cloud collision model to the observed clouds in W 33. We propose a scenario that a collision between a smaller cloud having a size of ∼7 pc (which corresponds to the 35 km s⁻¹ cloud) and a larger cloud with a size of ∼17 pc (the 58 km s⁻¹ cloud) occurred in this region. In this scenario, it can be interpreted that the arc-like structure in the 58 km s⁻¹ cloud that shows complementary distribution with the 35 km s⁻¹ cloud was created through the collision. As shown in the position–velocity diagram in figure 2b, these two clouds are separated in velocity. In addition, since the 35 km s⁻¹ cloud has a similar velocity compared with its surroundings outside the sightline of W 33, as shown in figure 2b, the deceleration of the colliding velocity was as small as less than 4–5 km s⁻¹, which corresponds to the linewidth of the 35 km s⁻¹ cloud. These signatures suggest the case 2 collision discussed in the previous subsection and shown in figure 12. The mean density of the 35 km s⁻¹ and 58 km s⁻¹ clouds are roughly estimated to be 1.7 × 10²² cm⁻²/7 pc ∼ 2 × 10³ cm⁻³ and 6.2 × 10²¹ cm⁻²/17 pc ∼ 400 cm⁻³ from the ¹²CO J = 1–0 dataset, respectively. In addition, the density of the molecular clumps W 33 Main in the center of the cavity associated with the 35 km s⁻¹ cloud and W 33 B associated with the 58 km s⁻¹ cloud are ∼10⁵ cm⁻³ and ∼10⁴ cm⁻³ from the C¹⁸O J = 1–0 dataset, respectively. Therefore we assumed that the density of the smaller cloud is much higher than that in the larger cloud. Based on the schematic of the final stage of the collision in figure 12, the timescale of the collision can be estimated to be 17–(17 + 7) pc/23 km s⁻¹ ∼ 0.7–1.0 Myr, depending on the lengths of the dense layer and the smaller cloud in this stage.

4.5 High-mass star formation in W 33

In the previous subsection, we discussed that the observed gas properties in W 33 can be interpreted with the cloud–cloud collision model. We discuss here relationships between the collision and the high-mass star formation in W 33.

As already introduced, W 33 contains high-mass star-forming regions in various evolutionary stages. The six dust clumps in W 33, i.e., W 33 Main, W 33 A, W 33 B, W 33 Main1, W 33 A1, and W 33 B1, have ages younger than the estimated timescale of the cloud–cloud collision, 0.7–1.0 Myr, because the dynamical timescales of the outflows in W 33 Main and W 33 A are ∼10⁴ yr. We suggest that these sources were formed in the dense layer created though the collision between the 35 km s⁻¹ and 58 km s⁻¹ clouds. It is possible that continuous compression and accumulation of gas by the collision led to a mass accretion rate as high as 10⁻⁴–10⁻³ M_⊙ yr⁻¹ for these objects, which satisfies the condition of high-mass star formation (e.g., Wolfire & Cassinelli 1987; Hosokawa et al. 2010; McKee & Tan 2003).

In addition to these sources, there are several extended H ii regions in W 33. Messineo et al. (2015) identified four new O stars in W 33, as indicated in figure 1 by the orange circles. Among them, two O4–6I stars are distributed in G012.745−00.153, while one O4–6I star is in G012.820−00.238. The authors discussed that these three O stars have ages less than 6 Myr based on a stellar model by Ekström et al. (2012) or 2–4 Myr by comparing with the O4–6 supergiants in the Arches cluster in the Galactic center (Martins et al. 2008), which show similar K-band spectra to the present O stars in W 33. The estimate of 2–4 Myr is larger than our estimate of 0.7–1.0 Myr on the timescale of the cloud–cloud collision in W 33.

In order to investigate the formation timescales of these O stars in a different way, we estimate here the formation timescale of the H ii regions G012.745−00.153 and G012.820−00.238, that are excited by these O4–6I stars. In figure 13 evolutionary tracks of the H ii region radius |$r_{\rm H\,\small{II}}$| are plotted for different initial densities, based on the analytical model of the D-type expansion by Spitzer (1978). Here, Lyman continuum photons N_Ly of 10^49.55–10^49.93 s⁻¹, which correspond to an O6I–O4I star, are assumed (Panagia 1973). The horizontal dashed lines in figures 13a and 13b indicate the radii of G012.745−00.153 and G012.820−00.238 measured by the 90 cm data, respectively.

In figure 13a, in which the case for G012.745−00.153 is presented, the evolution timescales which correspond to the measured |$r_{\rm H{\small II}}$| are about 0.2 Myr, 0.7 Myr, and 2 Myr for initial densities of 10³ cm⁻³, 10⁴ cm⁻⁴, and 10⁵ cm⁻³, respectively, suggesting that quite a high density of ∼10⁵ cm⁻³ is required to be consistent with the age of 2–4 Myr estimated by Messineo et al. (2015). However, in the case of G012.820−00.238 shown in figure 13b, even with quite a high density of 10⁵ cm⁻³, the derived timescale is as short as 0.5 Myr, and it is inconsistent with the estimate by Messineo et al. (2015) based on the comparison with the Arches cluster in the Galactic center, although their other estimate of less than 6 Myr, which is based on the stellar model by Ekström et al. (2012), is consistent with the present results in figure 13.

If we tentatively assume that the ages of G012.745−00.153 and G012.820−00.238 are 0.7 Myr and 0.2 Myr, respectively, in which a reasonable initial gas density of 10⁴ cm⁻³ is assumed, these figures are within the estimated timescale of the cloud–cloud collision in W 33, 0.7–1 Myr, and it is therefore possible that formation of these O stars in G012.745−00.153 and G012.820−00.238 occurred in the dense layer formed through the cloud–cloud collision between the 35 km s⁻¹ and 58 km s⁻¹ clouds.

5 Conclusions

The conclusions of the present study are summarized as follows.

1. We carried out large-scale CO observations using the NANTEN2 and Nobeyama 45 m telescopes toward the galactic high-mass star-forming region W 33. Our dataset identified three velocity clouds at 35 km s⁻¹, 45 km s⁻¹, and 58 km s⁻¹ toward W 33. The 35 km s⁻¹ cloud spatially coincides with W 33, showing spatial correlations with the dust clumps W 33 Main, W 33 A, and W 33 B1 and the H ii regions G012.745−00.153 and G012.820−00.238, while the 58 km s⁻¹ cloud, which has an arc-like structure surrounding the 35 km s⁻¹ cloud and thus W 33, shows associations with the dust clump W 33 B and the H ii regions G012.745−00.153 and G012.692−00.251. The 45 km s⁻¹ cloud shows diffuse and extended CO emission throughout the corresponding region of W 33, although its association with W 33 is not clear. The total molecular masses of the three clouds, i.e., the 35 km s⁻¹, 45 km s⁻¹, and 58 km s⁻¹ clouds, are 1.1 × 10⁵ M_⊙, 1.0 × 10⁵ M_⊙, and 3.0 × 10⁴ M_⊙, respectively.

2. We analyzed the intensity ratio between ¹²CO J = 3–2 and ¹²CO J = 1–0 in these three velocity clouds, finding that the 35 km s⁻¹ and 58 km s⁻¹ clouds have |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of higher than 1.0 in the vicinity of W 33, whereas the 45 km s⁻¹ cloud has a low |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of ∼0.4 throughout the analyzed region. The cause of the high |$R_{\rm 3\mbox{-}2/1\mbox{-}0}$| of the gas in the 35 km s⁻¹ and 58 km s⁻¹ clouds can be understood as heating by the high-mass stars in W 33, suggesting physical associations of these two velocity clouds with W 33.

3. In the H i 21 cm line profiles, strong absorption features are detected at 35 km s⁻¹ toward W 33 Main. This result indicates that W 33 Main is located inside or behind the 35 km s⁻¹ cloud.

4. We identified two CO bipolar outflows toward W 33 Main and W 33 A. The full velocity widths of the outflows are as large as 40–50 km s⁻¹. The dynamical timescales of the outflows can be estimated to be ∼3 × 10⁴ yr.

5. In order to interpret the co-location of the two velocity components at W 33, we discussed three possibilities: a gravitationally bound system, expanding motion by feedback from high-mass stars in W 33, and a cloud–cloud collision model. We found that neither of the former two assumptions can explain the observed signatures of the molecular gas.

6. Finally, we proposed a cloud–cloud collision model to interpret the multiple velocity clouds in W 33. We assumed a collision of two molecular clouds with different sizes at a relative velocity of 23 km s⁻¹. In this model, the arc-like structure in the 58 km s⁻¹ cloud can be reasonably explained as a cavity created by the collision. The uniform velocity distributions of the 35 km s⁻¹ and 58 km s⁻¹ clouds in the position–velocity diagram suggest that the gas density of the smaller cloud, which corresponds to the 35 km s⁻¹ cloud, is much higher than that in the large cloud, the 58 km s⁻¹ cloud. This is consistent with the observations. The timescale of the collision was estimated to be 0.7–1.0 Myr. We also discussed that it is possible for the high-mass star-forming regions in W 33 to have been formed in the dense gas layer created through the collision.

Acknowledgements

We are grateful to Misaki Hanaoka for a useful discussion. NANTEN2 is an international collaboration of ten universities: Nagoya University, Osaka Prefecture University, University of Cologne, University of Bonn, Seoul National University, University of Chile, University of New South Wales, Macquarie University, University of Sydney, and Zurich Technical University. The Nobeyama 45 m radio telescope is operated by the Nobeyama Radio Observatory, a branch of the National Astronomical Observatory of Japan. Data analysis was carried out on the open use data analysis computer system at the Astronomy Data Center, ADC, of the National Astronomical Observatory of Japan.

The James Clerk Maxwell Telescope is operated by the East Asian Observatory on behalf of The National Astronomical Observatory of Japan, Academia Sinica Institute of Astronomy and Astrophysics, the Korea Astronomy and Space Science Institute, the National Astronomical Observatories of China and the Chinese Academy of Sciences (Grant No. XDB09000000), with additional funding support from the Science and Technology Facilities Council of the United Kingdom and participating universities in the United Kingdom and Canada.

This work is based [in part] on observations made with the Spitzer Space Telescope, which is operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA. Support for this work was provided by NASA.

The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc.

The Australia Telescope Compact Array is part of the Australia Telescope National Facility which is funded by the Australian Government for operation as a National Facility managed by CSIRO.

This work is financially supported by Grants-in-Aid for Scientific Research (KAKENHI Nos. 15K17607, 15H05694) from MEXT (the Ministry of Education, Culture, Sports, Science and Technology of Japan) and JSPS (Japan Society for the Promotion of Science).

Appendix. CO velocity channel maps of the Nobeyama and JCMT data sets

We show the velocity channel maps of the ¹²CO, ¹³CO, C¹⁸O J = 1–0, ¹²CO J = 3–2 emissions, and the ¹²CO J = 1–0/J = 3–2 ratio of W 33 in figures 14–18, respectively. The velocity range is between 27.5 and 63.5 km s⁻¹. The contours indicate the 90 cm radio continuum emission with VLA.

References