A 33 GHz Survey of Local Major Mergers: Estimating the Sizes of the Energetically Dominant Regions from High-resolution Measurements of the Radio Continuum

We used the Karl G. Jansky VLA to observe radio continuum emission from the most luminous nearby LIRGs and ULIRGs. Our sample (see Table 1) consists of 22 sources from the IRAS Revised Bright Galaxy Sample (Sanders et al. 2003). These galaxies have IR luminosities ${L}_{\mathrm{IR}}[8-1000\,\mu {\rm{m}}]={10}^{11.6}-{10}^{12.6}\,{L}_{\odot }$ and were selected to be northern enough to be observed by the VLA, i.e., $\delta \gt -15^\circ$. These systems are also a subset of the Great Observatories All-sky LIRG Survey (GOALS; Armus et al. 2009), for which multiwavelength data are available.

Table 1.  Characteristics of the Sample Galaxies

Galaxy Name	Alternate Name	R.A. (J2000)	Decl. (J2000)	${\mathrm{log}}_{10}({L}_{\mathrm{IR}})$	DL	Scale (kpc/arcsec)	ID Number
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)
CGCG 436−030	MCG +02-04-025	${01}^{{\rm{h}}}{20}^{{\rm{m}}}02\buildrel{\rm{s}}\over{.} 722$	$+14^\circ 21^{\prime} 42\buildrel{\prime\prime}\over{.} 94$	11.64	127	0.601	1
IRAS F01364−1042	2MASX J01385289−1027113	${01}^{{\rm{h}}}{38}^{{\rm{m}}}52\buildrel{\rm{s}}\over{.} 921$	$-10^\circ 27^{\prime} 11\buildrel{\prime\prime}\over{.} 42$	11.81	201	0.942	2
III Zw 035	⋯	${01}^{{\rm{h}}}{44}^{{\rm{m}}}30\buildrel{\rm{s}}\over{.} 500$	$+17^\circ 06^{\prime} 05\buildrel{\prime\prime}\over{.} 00$	11.58	111	0.526	3
VII Zw 031	⋯	${05}^{{\rm{h}}}{16}^{{\rm{m}}}46\buildrel{\rm{s}}\over{.} 096$	$+79^\circ 40^{\prime} 13\buildrel{\prime\prime}\over{.} 28$	11.95	229	1.066	4
IRAS 08572+3915	⋯	${09}^{{\rm{h}}}{00}^{{\rm{m}}}25\buildrel{\rm{s}}\over{.} 390$	$+39^\circ 03^{\prime} 54\buildrel{\prime\prime}\over{.} 40$	12.13	254	1.176	5
UGC 04881	Arp 55	${09}^{{\rm{h}}}{15}^{{\rm{m}}}55\buildrel{\rm{s}}\over{.} 100$	$+44^\circ 19^{\prime} 55\buildrel{\prime\prime}\over{.} 00$	11.70	169	0.796	6
UGC 05101	⋯	${09}^{{\rm{h}}}{35}^{{\rm{m}}}51\buildrel{\rm{s}}\over{.} 595$	$+61^\circ 21^{\prime} 11\buildrel{\prime\prime}\over{.} 45$	11.97	168	0.792	7
MCG +07-23-019	Arp 148	${11}^{{\rm{h}}}{03}^{{\rm{m}}}53\buildrel{\rm{s}}\over{.} 200$	$+40^\circ 50^{\prime} 57\buildrel{\prime\prime}\over{.} 00$	11.61	149	0.704	8
NGC 3690	Arp 299	${11}^{{\rm{h}}}{28}^{{\rm{m}}}32\buildrel{\rm{s}}\over{.} 300$	$+58^\circ 33^{\prime} 43\buildrel{\prime\prime}\over{.} 00$	11.82	45.2	0.217	9
UGC 08058	Mrk 231	${12}^{{\rm{h}}}{56}^{{\rm{m}}}14\buildrel{\rm{s}}\over{.} 234$	$+56^\circ 52^{\prime} 25\buildrel{\prime\prime}\over{.} 24$	12.52	181	0.849	10
VV 250	UGC 08335 NED02	${13}^{{\rm{h}}}{15}^{{\rm{m}}}34\buildrel{\rm{s}}\over{.} 980$	$+62^\circ 07^{\prime} 28\buildrel{\prime\prime}\over{.} 66$	11.77	132	0.621	11
UGC 08387	Arp 193, IC 883	${13}^{{\rm{h}}}{20}^{{\rm{m}}}35\buildrel{\rm{s}}\over{.} 300$	$+34^\circ 08^{\prime} 21\buildrel{\prime\prime}\over{.} 00$	11.65	101	0.479	12
UGC 08696	Mrk 273	${13}^{{\rm{h}}}{44}^{{\rm{m}}}42\buildrel{\rm{s}}\over{.} 111$	$+55^\circ 53^{\prime} 12\buildrel{\prime\prime}\over{.} 65$	12.15	162	0.761	13
VV 340a	UGC 09618 NED02	${14}^{{\rm{h}}}{57}^{{\rm{m}}}00\buildrel{\rm{s}}\over{.} 826$	$+24^\circ 37^{\prime} 04\buildrel{\prime\prime}\over{.} 12$	11.67	144	0.665	14
VV 705	I Zw 107	${15}^{{\rm{h}}}{18}^{{\rm{m}}}06\buildrel{\rm{s}}\over{.} 344$	$+42^\circ 44^{\prime} 36\buildrel{\prime\prime}\over{.} 69$	11.87	172	0.807	15
IRAS 15250+3608	⋯	${15}^{{\rm{h}}}{26}^{{\rm{m}}}59\buildrel{\rm{s}}\over{.} 404$	$+35^\circ 58^{\prime} 37\buildrel{\prime\prime}\over{.} 53$	12.02	238	1.105	16
UGC 09913	Arp 220	${15}^{{\rm{h}}}{34}^{{\rm{m}}}57\buildrel{\rm{s}}\over{.} 116$	$+23^\circ 30^{\prime} 11\buildrel{\prime\prime}\over{.} 47$	12.16	77.2	0.369	17
IRAS 17132+5313	⋯	${17}^{{\rm{h}}}{14}^{{\rm{m}}}20\buildrel{\rm{s}}\over{.} 000$	$+53^\circ 10^{\prime} 30\buildrel{\prime\prime}\over{.} 00$	11.90	217	1.012	18
IRAS 19542+1110	⋯	${19}^{{\rm{h}}}{56}^{{\rm{m}}}35\buildrel{\rm{s}}\over{.} 440$	$+11^\circ 19^{\prime} 02\buildrel{\prime\prime}\over{.} 60$	12.07	277	1.277	19
CGCG 448−020	II Zw 096	${20}^{{\rm{h}}}{57}^{{\rm{m}}}23\buildrel{\rm{s}}\over{.} 900$	$+17^\circ 07^{\prime} 39\buildrel{\prime\prime}\over{.} 00$	11.79	148	0.698	20
IRAS 21101+5810	2MASX J21112926+5823074	${21}^{{\rm{h}}}{11}^{{\rm{m}}}30\buildrel{\rm{s}}\over{.} 400$	$+58^\circ 23^{\prime} 03\buildrel{\prime\prime}\over{.} 20$	11.74	162	0.764	21
IRAS F23365+3604	2MASX J23390127+3621087	${23}^{{\rm{h}}}{39}^{{\rm{m}}}01\buildrel{\rm{s}}\over{.} 273$	$+36^\circ 21^{\prime} 08\buildrel{\rm{s}}\over{.} 31$	12.16	273	1.262	22

Note. Column (1): name of the galaxy. Column (2): alternate name. Column (3): source right ascension (J2000) from NED. Column (4): source declination (J2000) from NED. Column (5): total IR luminosity from 8–1000 μm in log₁₀ solar units computed from the IRAS flux densities from Sanders et al. (2003) and following the equation from Sanders & Mirabel (1996). Column (6): luminosity distance from NED. Column (7): scale at Hubble flow distances from NED, corrected for the CMB, used to convert from arcseconds to kiloparsecs. Column (8): number used to identify each system in some of the figures presented in this paper.

Download table as:  ASCIITypeset image

As part of the resident shared risk project AL746, we observed the radio continuum emission from each source at C band (4–8 GHz) and Ka band (26.5–40 GHz). For each observation we used dual-polarization mode with two 1 GHz-wide basebands. Each band was split into eight 128 MHz spectral windows with 64 channels each. We centered the 1 GHz basebands at ∼4.7 and 7.2 GHz in C band and ∼29 and 36 GHz in Ka band.

In order to recover emission across a wide range of angular scales, we observed our sample in each frequency range in separate sessions using each of the four VLA configurations (A, B, C, and D, from highest to lowest angular resolution). Observations spanned the period from 2010 August 2 to 2011 August 16. In the D and C configurations, we observed each source for 5 minutes. In the B configuration, we observed each source for 10 minutes split between two 5-minute scans. In the A configuration we observed most sources for 20 minutes, split into four 5-minute scans. Due to scheduling constraints, eight sources were not observed in the A configuration at Ka band; these are identified with an asterisk in Table 3. Thus, the total time on source for most targets was ∼40 minutes per band.

At the beginning of each session, we observed either 3C 48 or 3C 286, which was used to set the flux density scale and calibrate the bandpass. Through the rest of the session we alternated between observations of science targets and a secondary calibrator within a few degrees of each science target. We used observations of these secondary calibrators to measure phase and gain variations due to atmospheric/ionospheric and instrumental fluctuations. Table 2 summarizes the calibrators used for each science target.

These data have also appeared in Leroy et al. (2011) and Barcos-Muñoz et al. (2015). Leroy et al. (2011) presented first results from our observations at both C and Ka band but used only observations from the two most compact VLA configurations. Barcos-Muñoz et al. (2015) presented C- and Ka-band observations using all four configurations for the specific case of Arp 220. In this paper, we report on the full survey, emphasizing the Ka-band observations and the combination of all four array configurations. These represent the highest-resolution, highest-sensitivity radio observations for these galaxies published to date. The C-band observations combining all four array configurations will be reported in an upcoming paper focused on the resolved spectral energy distribution (SED), i.e., across the disks of the systems in our sample (L. Barcos-Muñoz et al. 2017, in preparation).

We used the Common Astronomy Software Application (CASA; McMullin et al. 2007) to calibrate, inspect, and analyze the data. We followed a standard VLA reduction procedure, including calibrating the bandpass, phase, and amplitude response of each antenna. We set the overall flux density scale using "Perley-Butler 2010" models for the primary calibrators and assuming that the Ka-band emission shares the same structure as the VLA-provided Q-band model.

Once the data were calibrated, we imaged each science target. To do this, we used the task CLEAN in mode mfs (multifrequency synthesis; Sault & Wieringa 1994), with Briggs weighting setting robust = 0.5. For each array configuration, we imaged each baseband independently. Whenever possible, we iterated this imaging with phase and amplitude self-calibration. The number of self-calibration iterations varied from zero to eight based on the signal-to-noise ratio of the data, with four iterations typical. After several iterations of phase-only self-calibration, when possible, we also performed amplitude self-calibration. We always solved for only relative variations in the amplitude gains among the antennas (solnorm = True in CASA's gaincal) and hence avoided forcing the flux of the source to some value.

After self-calibrating the two basebands independently, we combined both into a single image using clean's multifrequency synthesis mode and nterms = 2. The latter allows us to model the frequency dependence of the sky emission with two Taylor coefficients. After the described combination, we ended up with four images per source (one per array configuration). Finally, we jointly imaged all self-calibrated data, combining all eight measurement sets (four configurations and two basebands). This combined image represents our best data product, using all of our observations with sensitivity to a wide range of angular scales. In the highest signal-to-noise cases, for example, UGC 08058 (Mrk 231) and UGC 09913 (Arp 220), we performed further self-calibration during this final imaging step.

These final images have a nominal frequency $\nu =32.5$ GHz¹⁶ and a typical rms noise $26\,\mu \mathrm{Jy}\,{\mathrm{beam}}^{-1}$. Table 3 reports the exact beam size and rms noise for the combined image for each target.

We combine our survey with previous observations of our sample at 1.49 GHz (beam FWHM ∼ 15'') from Condon et al. (1990). We also use the 5.95 GHz flux densities (beam FWHM ∼ 04) from Leroy et al. (2011) and CO (1−0) flux densities, obtained using the ARO 12 m antenna (FWHM = 1'), the latter of which will be reported in G. Privon et al. (2017, in preparation). We present a compilation of the flux densities at these different frequencies, along with 32.5 GHz flux densities measured from our new images, in Table 3. The uncertainties in the 1.49 GHz flux density values are assumed to be dominated by flux density calibration errors (∼5%; see Section 3 in Condon et al. 1990).

Five of our sources lack flux density measurements at 1.49 GHz. For three of these—VII Zw 031, CGCG 448−020, and IRAS F23365+3608—we use the 1.4 GHz flux density from the NVSS catalog (Condon et al. 1998). For IRAS 19542+1110 and IRAS 21101+5810, we take the values at 1.425 GHz measured by Condon et al. (1996). We use these flux densities interchangeably with the 1.49 GHz fluxes, but assign them a larger (10%) uncertainty in these cases to reflect some uncertainty in the spectral index.

We measure integrated flux densities for each source from the lowest angular resolution observations, which were obtained using the VLA in its D configuration. The maximum recoverable scale for the D configuration, $\approx 22^{\prime\prime}$, corresponds to ∼16 kpc at the 165 Mpc median distance of our sample. This would cover most of the star-forming activity in a local disk galaxy (e.g., Schruba et al. 2011). U/LIRGs are observed to be much more compact, with sizes less than a few kiloparsecs based on previous radio (Liu et al. 2015), near-IR (Haan et al. 2011), mid-IR (Díaz-Santos et al. 2010), and far-IR (FIR) observations (Lutz et al. 2016). Therefore, we expect negligible missing flux in the D-configuration-based flux densities.

Confirming this expectation, most of our targets appear unresolved in the D-configuration-only images, which have beam sizes $\approx 2\buildrel{\prime\prime}\over{.} 7$. We obtained the flux densities reported in Table 3 using CASA task imfit to fit 2D Gaussians to these mostly unresolved point sources. A few targets, including NGC 3690, CGCG 448−020, IRAS 17132+5313, VII Zw 031, VV 250, VV 340, and VV 705, showed some extent or multiple components in the D configuration maps. In most of these cases, we tapered the D configuration data to a lower resolution until the morphology became a single point-like source. Then we fit a 2D Gaussian to this degraded image. NGC 3690 and VV 250 show well-separated components that can only be fit using two Gaussians, even in the tapered images. We report the sum of these two components as the integrated flux density.

The uncertainties that we report sum (in quadrature) the statistical error calculated by imfit with uncertainty in the calibration of the flux scale, which we estimated to be ∼12% in Barcos-Muñoz et al. (2015). For the two faintest galaxies in our sample, UGC 04881 and IRAS 08572+3915, the signal-to-noise ratio of the D configuration data only was not high enough to recover integrated flux densities. For these two systems, we instead report results from the combined data using all configurations, which we tapered until we recovered point-like structures that could be fit using Gaussians.

In addition to our new 32.5 GHz flux densities, Table 3 reports literature flux densities for our sources at $\nu =1.49$ and 5.95 GHz. We combine these with the $\nu =32.5\,\,\mathrm{GHz}$ measurements to calculate the galaxy-integrated spectral index between 1.49 and 5.95 GHz (${\alpha }_{1.5-6}$) and between 5.95 and 32.5 GHz (${\alpha }_{6-33}$). Here, we define the spectral index, α, by ${S}_{\nu }\propto {\nu }^{\alpha }$. Note that because we use the flux density integrated over the whole galaxy, we do not expect the different angular resolutions at different frequencies to affect these calculations.

In Figure 2, we show the derived spectral indices. We plot ${\alpha }_{1.5-6}$ as a function of ${\alpha }_{6-33}$. Here the solid line shows equal spectral indices for both pairs of bands, which we would expect if a single spectral index holds across the entire radio regime (from 1.5 to 33 GHz). Dashed lines show $\alpha =-0.8$, a typical spectral index for synchrotron emission without any opacity effects (e.g., Condon 1992).

Zoom In Zoom Out Reset image size

A main goal of our study is to measure the extent of the radio continuum emission in our targets with the purpose of constraining the size of the energetically dominant region.

To do this, we analyzed the final images combining data from all the array configurations. These high-resolution images are sensitive to the brightest compact cores, but they have lower surface brightness sensitivity than the D configuration data that we used to determine the total flux density. Therefore, they may miss extended, low surface brightness emission. To take this into account, we measure the size of the energetically dominant region from the half-light area (A₅₀). This is the area enclosed by the highest-intensity isophote that includes half of the total integrated flux density of the system, which we measured from the lower-resolution data above and expect to be complete. Note that this approach measures the observed A₅₀, which reflects the true size of the source convolved with the synthesized beam of the array.

We require the intensity of the isophote enclosing the half-light area, or C₅₀, to be at least 5 times the rms noise in the image. If C₅₀ would be less than 5σ in the combined image, we interpret this to indicate an important component of extended, low surface brightness emission. In order to recover this emission, we measure A₅₀ for these systems from lower-resolution versions of the data that have better surface brightness sensitivity. In these cases, we first tried using natural weighting instead of Briggs (see Section 2). If we still could not recover half of the light within an S/N > 5 contour, then we produced progressively lower resolution images by applying larger and larger u–v-tapers to the data. We stepped the size of the taper by $0\buildrel{\prime\prime}\over{.} 2$ and used Briggs weighting schemes with robust parameter 0.5 at each step. In this way, we measure A₅₀ from the highest-resolution image where C₅₀ can be reliably measured, i.e., where C₅₀ $\geqslant$ 5σ.

The following systems showed extended, low surface brightness emission and required u–v-tapering: CGCG 436−030, CGCG 448−020, IRAS 21101+5810, IRAS 17132+5313, VV 340a, and VV 705. For NGC 3690, the natural weighting approach was sufficient.

Once we identified a reliable half-light contour, C₅₀, we calculated A₅₀ by multiplying the number of pixels within the C₅₀ contour by the pixel area. Figure 1 shows the images that were used to measure A₅₀ and the C₅₀ contour (in red) for each source.

Many of our sources show sizes close to that of the synthesized beam. We show this in Figure 3. There, we plot the ratio of the observed A₅₀ to the beam area, A_beam, as a function of the beam area in units of arcsec² (top left panel) and kpc² (top right panel).

Zoom In Zoom Out Reset image size

The quantity of physical interest is the true size of the 33 GHz emission with the beam deconvolved, A_50d. In the top and bottom left panels of Figure 3, a dashed line indicates a value of $2\times \ {{\rm{A}}}_{\mathrm{beam}}$, which we consider a practical threshold for the emission to be viewed as resolved. Here ${A}_{\mathrm{beam}}\,=\,\pi {\theta }_{\mathrm{maj}}{\theta }_{\min }/4$, with ${\theta }_{\mathrm{maj}}\,\mathrm{and}\,{\theta }_{\min }$ the FWHM of the synthesized beam along its major and minor axis. In this definition, A_beam refers to the area expected to enclose half the total power in the beam. This definition is consistent with our measured area A₅₀ and appropriate for deconvolution.

We treat the sources that show extent larger than the beam but size smaller than $2\times {A}_{\mathrm{beam}}$ as marginally resolved (region between the solid and dashed lines in Figure 3). In these cases, we assume that the intrinsic shape (deconvolved size) of the source follows a Gaussian distribution. We then estimate the deconvolved size of the source by ${A}_{50}(\mathrm{deconvolved})\,={A}_{50}(\mathrm{observed})-{A}_{\mathrm{beam}}\equiv {A}_{50{\rm{d}}}$, equivalent to deconvolving the FWHM in quadrature.

In the top panels of Figure 3, two sources lie below the solid line, indicating an observed size smaller than the beam. These are IRAS F08572+3915 and UGC 04881NE. Although statistical fluctuations could produce this situation, the signal-to-noise ratio of the data appears to be too high for this explanation to hold. The most likely culprit is a calibration issue when combining observations using the different array configurations. We adopt a conservative upper limit of ${A}_{50{\rm{d}}}\lt {A}_{\mathrm{beam}}$ for these two systems.¹⁷

In order to determine the best estimate of ${A}_{50{\rm{d}}}$ for "resolved" sources with ${A}_{50}\gt 2\times {A}_{\mathrm{beam}}$, we inspected the shape of the C₅₀ contour (red in Figure 1) to determine whether the source exhibits a Gaussian shape. If it did, then we apply the same approach used for the marginally resolved sources to each component and summed the results to find the total ${A}_{50{\rm{d}}}$. This tended to be the case when more than one component is present, such as VV 705 and CGCG 448−020.

If C₅₀ showed a more complex morphology, our simple Gaussian treatment becomes invalid. In these cases we instead assume that the measured (not deconvolved) A₅₀ represents an upper limit to the true size. This is true for the following galaxies: IRAS 19542+1110, IRAS F23365+3604, UGC 08387, VII Zw 031, VV 250a, and VV 340a.

For two sources, UGC 04881 and VV 250, a second, faint component could be recovered only in the low-resolution map used to assess the integrated flux density. In both cases, the individual components are unresolved in this integrated map. Here, we had to lower our conservative limit of 5σ in order to recover the half-light area. In these two systems, we measure C₅₀ from a contour with ${\rm{S}}/{\rm{N}}\approx 3$ and treat the size estimate as an upper limit (see Table 4).

For NGC 3690 and IRAS 17132+5313, one component of the C₅₀ contour shows a Gaussian distribution, while others show more complex morphology. In both cases, we performed the deconvolution on the Gaussian components. Then we have a partially deconvolved estimate, ${A}_{50{\rm{d}}}\lt {A}_{50}(\mathrm{observed})$, which is still an upper limit because of the un-deconvolved, more complex structure. We report values for ${A}_{50{\rm{d}}}$ and C₅₀ in Table 4, along with an equivalent ${R}_{50{\rm{d}}}$ value where ${A}_{50{\rm{d}}}=\pi {R}_{50{\rm{d}}}^{2}$. We caution, however, that R_50d is only a representative number reflective of the upper limit to the area in these cases.

In Table 4, we also report the degree of Gaussianity, defined as the ratio between the flux density level of the C₅₀ contour and the peak flux density. For a 2D Gaussian, this value is 0.5.

In addition to the integrated flux density and a characteristic size, we measured the contribution of compact sources to the overall flux density of each target using the maps of Figure 1. For our purposes, compact sources are those that clearly belong to the system and show a Gaussian morphology.

For each target, we identify these sources by eye and fit them using imfit, providing estimates of the rms noise and reasonable starting guesses for the sizes and peak intensity and position. The locations of the fit compact sources appear as white crosses in Figure 1. Their sizes, which are often comparable to the size of the beam, are shown in the bottom left panel of Figure 3. We also calculated the flux density that is originating from all the compact sources in a system and compared it to the integrated flux density (see top panels in Figure 4). We note that such comparisons may be affected by the different physical resolutions achieved from the observations; however, we find no trend relating the fraction of flux in compact sources to beam physical area. In the bottom panel of Figure 4 we show instead the contribution of each point source—especially important when more than one is present—to the integrated flux density at 33 GHz.

Zoom In Zoom Out Reset image size

We identified compact sources in each of our targets except the northeast component in IRAS F17132+5313, which shows mostly extended emission. For the cases of the faint components in the systems UGC 04881 and VV 250, the Gaussian fit was performed on the low-resolution image that was used to obtain the integrated flux density of the system.

A subset of our sources show most of their emission concentrated into a very small area, consistent with a point source producing much of the flux density even at our highest angular resolution. To make the strongest possible measurement of the compactness of these targets, we used our highest-resolution images. This is usually the A configuration image ($\sim 0\buildrel{\prime\prime}\over{.} 1$), except in those cases with B as the longest baseline array configuration observed ($\sim 0\buildrel{\prime\prime}\over{.} 2;$ see Table 3).

From this highest-resolution image, we measured the flux density detected at S/N ≥ 5, which corresponds well to the total flux density in the compact core of the image. We compared this flux density in the bright core at the highest resolution to the integrated flux density of the system, ${f}_{A(\mathrm{or}B)}$. Most of the U/LIRGs in this sample show single bright point sources in the highest-resolution image, although a few, including NGC 3690, UGC 08387, Arp 220, and VV 705, show more than one compact core.

We also measured the size of the 33 GHz emission showing significant detection, as set by the $5{\sigma }_{A(\mathrm{or}B)}$ ¹⁸ contour, at this highest angular resolution image. We report the beam size of the A, or B, array configuration images along with the sizes of the $5{\sigma }_{A(\mathrm{or}B)}$ contour and ${f}_{A(\mathrm{or}B)}$ of each system in Table 5. We highlight those sources with most of their emission being contributed by a single bright compact source, being good potential AGN candidates. These include IRAS F01364−1042, III Zw 035, and IRAS 15250+3609. Arp 220 should also be on this list as it shows ${f}_{A(\mathrm{or}B)}\gt 50 \%$; however, we refer the reader to a more exhaustive discussion on the morphology of its 33 GHz emission presented in Barcos-Muñoz et al. (2015). There are six other sources with ${f}_{A(\mathrm{or}B)}\gt 50 \%$, but unfortunately the highest resolution achieved was only $\sim 0\buildrel{\prime\prime}\over{.} 2$ and the constraint on their compactness is then weaker. However, note that Mrk 231, a known AGN (e.g., Ulvestad et al. 1999; Lonsdale et al. 2003), belongs to that group.

Table 5.  Analysis of A (or B)^a Array Configuration-only Images

Galaxy Name	${\theta }_{M}\times {\theta }_{m}$	${\sigma }_{A(\mathrm{or}B)}$ (μJy beam−1)	${f}_{A(\mathrm{or}B)}$ (%)	Log10(${A}_{5\sigma }$) (arcsec2)	${R}_{5\sigma }$ (kpc)
(1)	(2)	(3)	(4)	(5)	(6)
CGCG 436−030	$0\buildrel{\prime\prime}\over{.} 072\times 0\buildrel{\prime\prime}\over{.} 061$	36.5	16.4	−1.873	0.039
IRAS F01364−1042	$0\buildrel{\prime\prime}\over{.} 101\times 0\buildrel{\prime\prime}\over{.} 060$	44.0	59.3	−1.530	0.091
III Zw 035	$0\buildrel{\prime\prime}\over{.} 073\times 0\buildrel{\prime\prime}\over{.} 062$	45.6	61.3	−1.476	0.054
VII Zw 031	$0\buildrel{\prime\prime}\over{.} 119\times 0\buildrel{\prime\prime}\over{.} 062$	24.8	0.3	−3.222	0.015
IRAS 0857+3915	$0\buildrel{\prime\prime}\over{.} 241\times 0\buildrel{\prime\prime}\over{.} 180$	28.8	97.5	−0.831	0.255
UGC 04881	$0\buildrel{\prime\prime}\over{.} 247\times 0\buildrel{\prime\prime}\over{.} 184$	27.7	84.3	−0.739	0.192
UGC 05101	$0\buildrel{\prime\prime}\over{.} 259\times 0\buildrel{\prime\prime}\over{.} 240$	27.2	74.4	−0.112	0.393
MCG +07-23-019	$0\buildrel{\prime\prime}\over{.} 216\times 0\buildrel{\prime\prime}\over{.} 189$	31.5	67.8	−0.615	0.196
NGC 3690	$0\buildrel{\prime\prime}\over{.} 239\times 0\buildrel{\prime\prime}\over{.} 218$	26.6	42.3	0.312	0.175
UGC 08058	$0\buildrel{\prime\prime}\over{.} 227\times 0\buildrel{\prime\prime}\over{.} 202$	27.6	80.4	−0.299	0.340
VV 250	$0\buildrel{\prime\prime}\over{.} 219\times 0\buildrel{\prime\prime}\over{.} 202$	26.5	38.2	−0.305	0.249
UGC 08387	$0\buildrel{\prime\prime}\over{.} 073\times 0\buildrel{\prime\prime}\over{.} 051$	24.4	39.8	−1.301	0.060
UGC 08696	$0\buildrel{\prime\prime}\over{.} 210\times 0\buildrel{\prime\prime}\over{.} 204$	28.4	65.5	−0.233	0.328
VV 340a	$0\buildrel{\prime\prime}\over{.} 085\times 0\buildrel{\prime\prime}\over{.} 065$	17.8	0.9	−2.959	0.013
VV 705	$0\buildrel{\prime\prime}\over{.} 059\times 0\buildrel{\prime\prime}\over{.} 051$	20.9	17.2	−1.712	0.063
IRAS 15250+3609	$0\buildrel{\prime\prime}\over{.} 058\times 0\buildrel{\prime\prime}\over{.} 051$	25.9	73.2	−1.767	0.082
Arp 220	$0\buildrel{\prime\prime}\over{.} 066\times 0\buildrel{\prime\prime}\over{.} 052$	24.8	66.5	−0.589	0.106
IRAS 17132+5313	$0\buildrel{\prime\prime}\over{.} 060\times 0\buildrel{\prime\prime}\over{.} 053$	24.5	10.5	−2.137	0.049
IRAS 19542+1110	$0\buildrel{\prime\prime}\over{.} 072\times 0\buildrel{\prime\prime}\over{.} 063$	25.7	49.4	−1.472	0.132
CGCG 448−020	$0\buildrel{\prime\prime}\over{.} 073\times 0\buildrel{\prime\prime}\over{.} 063$	27.8	23.1	−1.752	0.052
IRAS 21101+5810	$0\buildrel{\prime\prime}\over{.} 075\times 0\buildrel{\prime\prime}\over{.} 052$	27.0	26.9	−1.785	0.055
IRAS F23365+3604	$0\buildrel{\prime\prime}\over{.} 069\times 0\buildrel{\prime\prime}\over{.} 062$	37.1	33.9	−1.818	0.088

Note. Column (1): galaxy name, with those galaxies that represent good AGN candidates in bold face. Those with a weaker argument to be potential AGNs are underlined (see Section 3.4 for more details). Column (2): beam size of the A (or B) array configuration-only images (major × minor axis). Column (3): rms of the A (or B) array configuration image. Column (4): percentage of the total flux density recovered at 33 GHz from the A (or B) array configuration-only image. This A (or B) array configuration-only flux density was obtained adding pixels with emission above $5{\sigma }_{A(\mathrm{or}B)}$. Column (5): observed area of the 5${\sigma }_{A}$ contour in arcseconds². Column (6): equivalent radii for column (5), assuming ${A}_{5\sigma }=\pi \,{R}_{5\sigma }^{2}$ in units of kpc obtained using the scale conversion from Table 1.

^aEight systems were not observed with the A configuration of the VLA. See Table 3 for reference.

Download table as:  ASCIITypeset image

In Table 5 we also note two systems, VII Zw 031 and VV 340a, with ${f}_{A(\mathrm{or}B)}\lt 1 \%$ indicating that most of their emission at $0\buildrel{\prime\prime}\over{.} 1$ resolution is filtered out and then is mostly extended in nature.

For a resolved or nearly resolved source, where beam filling is a minor consideration, the brightness temperature, T_b, offers the prospect to constrain the emission mechanism and opacity of the source (e.g., Condon et al. 1991). At radio frequencies, the brightness temperature, ${T}_{{\rm{b}}}$, follows the Rayleigh–Jeans approximation where

$$\begin{eqnarray}&&{T}_{{\rm{b}}}=\left(\displaystyle \frac{{S}_{\nu }}{{{\rm{\Omega }}}_{\mathrm{source}}}\right)\displaystyle \frac{{c}^{2}}{2{k}_{B}{\nu }^{2}}\,,\end{eqnarray} \tag{ 1 }$$

with ${S}_{\nu }$ the flux density at frequency ν and ${{\rm{\Omega }}}_{\mathrm{source}}$ the area of the source.

Most of our targets are resolved. Thus, an "averaged nuclear ${T}_{{\rm{b}}}$" at 32.5 GHz can be derived using ${{\rm{\Omega }}}_{\mathrm{source}}={A}_{50{\rm{d}}}$ and ${S}_{\nu }=0.5\times {S}_{32.5}$ (see above for the explanation of the aperture correction). We also calculate ${T}_{{\rm{b}}}$ from the point of highest intensity in the highest-resolution image for each target, peak ${T}_{{\rm{b}}}$, where ${{\rm{\Omega }}}_{\mathrm{source}}={{\rm{\Omega }}}_{\mathrm{beam}}$ in that case. Figure 5 shows histograms of these peak and averaged nuclear ${T}_{{\rm{b}}}$ at $\nu =32.5\,\,\mathrm{GHz}$.

Zoom In Zoom Out Reset image size

The averaged nuclear ${T}_{{\rm{b}}}$ for our targets is typically a few tens of Kelvin to a few times 100 K, reaching up to a few thousand Kelvin in the brightest targets.

For only free–free emission filling the beam, we would expect ${T}_{{\rm{b}}}$ for optically thick emission to approach ${T}_{{\rm{e}}}$ for the H ii regions. For physical conditions like those present in our sample, the expected electron temperature ${T}_{{\rm{e}}}\sim 5000-{10}^{4}$ K (Hummer & Storey 1987; Condon 1992). In metal-rich environments, such as the central regions of ULIRGs (Veilleux et al. 2009), the cooling is more efficient and ${T}_{{\rm{e}}}$ may tend toward the low end of this range, ∼5000 K (e.g., Puxley et al. 1989), though note that Anantharamaiah et al. (2000) found ${T}_{{\rm{e}}}$ of 7500 K for Arp 220 from integrated measurements of radio recombination lines.

In Figure 5 we observe that ${T}_{{\rm{b}}}$ does not exceed either 10⁴ K or 5000 K for any galaxy. In theory the unresolved, or marginally resolved, sources could be optically thick and highly clumped at scales much smaller than the beam size. However, both the observed spectral index (which would be positive with $\alpha \sim 2$ for the free–free emission if optically thick) and the relative smoothness of the images argue against such a scenario. Instead, low opacity at 32.5 GHz appears to be the natural explanation for the ${T}_{{\rm{b}}}$ that we observe.

In Figure 5 we observe that the peak brightness temperatures do not exceed the likely ${T}_{{\rm{e}}}$. However, ${T}_{{\rm{b}}}$ (peak) should be treated as a lower limit for the unresolved and marginally resolved sources. Are these sources likely to be optically thick? Excluding the case of Mrk 231 since it hosts an AGN, the lower limits for the peak ${T}_{{\rm{b}}}$ go from 20 K up to 690 K, with the unresolved case, UGC 04881, having a temperature of 22 K. In the marginally resolved cases, we can gain insight into the likely size of the source by contrasting the peak and average ${T}_{{\rm{b}}}$. Figure 3 shows that for most of these marginally resolved with T_b(peak) < T_b(average), we would expect to be able to resolve them with a beam area that is 2 times smaller at most. This would imply a true ${T}_{{\rm{b}}}$ peak ∼2 times larger than what we measure, still not enough for these sources to reach the optically thick regime. In these marginally resolved cases, in particular, the substructure of the emission remains unclear. Our data offer limited insight into whether the data may be structured into smaller optically thick regions beneath the beam.

For the unresolved source, the situation is less clear. With the size unconstrained, the source could be optically thick at 33 GHz and heavily beam diluted. However, we note again that the spectral index that we observe does not appear consistent with optically thick free–free emission. We proceed assuming that we observe optically thin 33 GHz emission for this source.

The flux densities of many of our targets have been measured at 1.4 GHz (Table 3), but even in its most extended configuration, the VLA reaches only $\approx 1^{\prime\prime}$ resolution at this frequency. Using the measured 1.4 GHz flux densities, we calculate the averaged nuclear ${T}_{{\rm{b}}}$ at 1.4 GHz assuming that the 32.5 GHz sizes also describe the true extent of the 1.4 GHz emission. These span from 10³ up to ${10}^{6.5}$ K.

These are high values. Values of ${T}_{{\rm{b}}}\gt 5\times {10}^{3}$ or 10⁴ K imply that the emission at 1.49 GHz is mostly synchrotron in nature, because the source function of the free–free emitting ionized gas is a blackbody at $T\sim 5\times {10}^{3}\mbox{--}{10}^{4}\,{\rm{K}}$, as explained above.

Dominant synchrotron emission may be expected at 1.4 GHz, but the values that we find may in fact be too high for the standard mixture of free–free and synchrotron emission seen in starburst galaxies. Considering such a mixture, Condon et al. (1991) suggested a maximum ${T}_{{\rm{b}}}$ for a starburst of ∼10${}^{4.6}$ K at 1.49 GHz (their Equation (9), using ${T}_{{\rm{e}}}\sim 5000\,{\rm{K}}$). At least 12 sources in the sample show ${T}_{{\rm{b}}1.4\mathrm{GHz}}\gt {10}^{4.6}$ K when we combine the 33 GHz sizes and the 1.49 GHz flux densities (see Table 6). This could imply that the 1.49 GHz emission from these sources includes a significant AGN contribution. One of those sources, Mrk 231, is well known to be dominated by an AGN, which explains why it has the highest predicted averaged nuclear ${T}_{{\rm{b}}}$ at 1.49 GHz.

Based on this line of argument, for these high-brightness 1.49 GHz sources, we would expect much of the flux density to be confined to an unresolved core in VLA 1.49 GHz imaging. In Mrk 231, most of the emission is unresolved at 1.49 GHz; however, other sources show resolved emission at 1.49 GHz. In these cases, the 33 GHz sizes, which are small compared to the $\sim 1\buildrel{\prime\prime}\over{.} 5$ VLA beam, may not be representative of the true 1.49 GHz emission. Indeed, we might worry that the 33 GHz size will underestimate the size at 1.49 GHz if the system is optically thick at these lower frequencies. In such a case, the emission will emerge from a photosphere larger than the emitting (optically thin) region at 33 GHz and the true brightness temperature at 1.49 GHz will be lower than our estimate.

Another alternative is that an extended synchrotron component may contribute to the integrated flux density. This component would have to have a spectrum steep enough that it does not contribute much to the flux density at 33 GHz, implying substantial variations in the resolved spectral index.

On the other hand, several sources have high T_b and remain barely resolved even at 8.44 GHz (see maps in Condon et al. 1991): IRAS 08572+3915, IRAS 17132+5313, IRAS 15250+3609, and III Zw 035. These are our best AGN candidates based on T_b arguments. Here extra information is needed to determine whether they are powered by an AGN and/or starbursts. In an upcoming paper, we investigate this possibility by combining the current observations with the lower-frequency ($\nu =4-8$ GHz) part of our survey (L. Barcos-Muñoz et al. 2017, in preparation).

IR luminosity, ${L}_{\mathrm{IR}}[8-1000\,\mu {\rm{m}}]$, and radio emission both trace recent star formation in starburst galaxies. IR luminosity reflects reprocessed light from young stars, while the 33 GHz continuum predominately captures a mix of synchrotron and thermal emission, both of which originate indirectly from young stars.

Considering a mix of synchrotron and thermal emission, Murphy et al. (2012) relate the recent SFR to the 33 GHz luminosity, ${L}_{33\mathrm{GHz}}$, via

$$\begin{eqnarray}\left(\displaystyle \frac{{\mathrm{SFR}}_{\nu }}{{M}_{\odot }\,{\mathrm{yr}}^{-1}}\right) & = & {10}^{-27}\left[2.18{\left(\displaystyle \frac{{T}_{{\rm{e}}}}{{10}^{4}{\rm{K}}}\right)}^{0.45}{\left(\displaystyle \frac{\nu }{\mathrm{GHz}}\right)}^{-0.1}\right.\\ & & +{\left.15.1{\left(\displaystyle \frac{\nu }{\mathrm{GHz}}\right)}^{{\alpha }^{\mathrm{NT}}}\right]}^{-1}\left(\displaystyle \frac{{L}_{\nu }}{\mathrm{erg}\,{{\rm{s}}}^{-1}\,{\mathrm{Hz}}^{-1}}\right),\end{eqnarray} \tag{ 2 }$$

where T_e is the electron temperature and ${\alpha }_{\mathrm{NT}}$ is the nonthermal spectral index. Murphy et al. (2012) relate the IR luminosity to the recent SFR via

$$\begin{eqnarray}&&\left(\displaystyle \frac{{\mathrm{SFR}}_{\mathrm{IR}}}{{M}_{\odot }\,{\mathrm{yr}}^{-1}}\right)=3.15\times {10}^{-44}\left(\displaystyle \frac{{L}_{\mathrm{IR}}[8-1000\,\mu {\rm{m}}]}{\mathrm{ergs}\,{{\rm{s}}}^{-1}}\right).\end{eqnarray} \tag{ 3 }$$

In the left panel of Figure 6 we compare IR-based and 33 GHz based SFRs estimated for each U/LIRG in our sample. Following Murphy et al. (2012), we adopt ${T}_{{\rm{e}}}={10}^{4}$ K and ${\alpha }_{\mathrm{NT}}=-0.8$ at $\nu =32.5\,\,\mathrm{GHz}$, but note both as a source of uncertainty. If we use ${T}_{{\rm{e}}}=5000$ K, the SFR based on 33 GHz increases by ∼37%.

Zoom In Zoom Out Reset image size

The left panel of Figure 6 shows that these two simple SFR estimates agree in our sample. The strong outlier, source #10, is Mrk 231. This system is known to be dominated by an AGN that appears to contribute substantially to the 33 GHz emission. The other sources are consistent with a simple radio–IR correlation that has a normalization in agreement with the Murphy et al. (2012) relations.

If the assumption is made that the 33 GHz size, ${A}_{50,{\rm{d}}}$, reflects the distribution of star formation, we can derive a SFR surface density, ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}$. As above, we take ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}\,=0.5\times {\mathrm{SFR}}_{\mathrm{IR}}/{A}_{50,{\rm{d}}}$.¹⁹

The right panel in Figure 6 shows our calculated ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}$. These span from ${10}^{0.6}$ up to ${10}^{4.1}\,{M}_{\odot }\,{\mathrm{yr}}^{-1}\,{\mathrm{kpc}}^{-2}$ (right panel, bottom axis, in Figure 6). The high end of this range represents the highest ${{\rm{\Sigma }}}_{\mathrm{SFR}}$ found for any galaxy in the local universe. The wide range indicates diverse conditions. Even though we have observed the brightest and closest U/LIRGs, these span about four orders of magnitude in ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}$.

The IR surface brightness is also of interest. In local U/LIRGs, most of the bolometric luminosity is emitted in the 8−1000 μm range. By assuming that half of ${L}_{\mathrm{IR}}$ is concentrated within ${A}_{50,{\rm{d}}}$, we estimate ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ for this inner region. For our approach from Equation (3), ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ is identical to ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}$ within a constant factor. Therefore, we show the ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ axis along the top of the right panel of Figure 6.

The U/LIRGs in this sample have ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ ranging from ${10}^{10.5}\,\mathrm{to}\,{10}^{14.1}\,{L}_{\odot }\,{\mathrm{kpc}}^{-2}$. The high end of this range is of particular interest. The dashed vertical line in Figure 6 indicates ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}={10}^{13}\,{L}_{\odot }\,{\mathrm{kpc}}^{-2}$. This value of ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ has been argued to correspond to the characteristic Eddington limit set by radiation pressure on dust in self-regulated, optically thick disks (Thompson et al. 2005). Some sources in our sample show ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}\geqslant {10}^{13}\,{L}_{\odot }\,{\mathrm{kpc}}^{-2}$, indicating that they may be Eddington-limited starbursts (see Section 5.5 for further discussion).

Note that for systems that are optically thick in the IR, the ${\tau }_{\mathrm{IR}}\sim 1$ photosphere may be larger than the 33 GHz size. In this case, the ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ that we calculate would never be observed, even if very high resolution FIR observations were available. This does not mean that this quantity lacks physical meaning, however. These systems are incredibly opaque to UV and optical light, which we expect to be generated in the region of active star formation traced by our 33 GHz data. This will then be quickly reprocessed into IR light, which then scatters out to the photosphere before leaving the system.

In this case, that inner region captured by the 33 GHz emission and ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ and ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}$ are the quantities directly related to the region of most intense feedback and the immediate sites of star formation.

Are our sources optically thick in the IR? IR observations are limited to relatively coarse resolution, so direct size measurements in this range provide only modest constraints. Díaz-Santos et al. (2010) and Lutz et al. (2016) measured IR sizes²⁰ for the systems in our sample at 13 μm (Spitzer) and 70 μm (Herschel).²¹ For most systems, the sizes at 13 μm are larger than those at 70 μm, and the latter are larger than those we measured at 33 GHz.

A more powerful constraint comes from comparing our measured size to that implied by the measured dust temperature and luminosity. To do this, we consider the emission emitted in the IR, specifically between 8 and 1000 μm, and the dust temperature found comparing 63 and 158 μm emission (for more details, see Diaz-Santos et al. submitted). For an optically thick blackbody of temperature T_dust,

$$\begin{eqnarray}&&0.5\times {L}_{\mathrm{IR}}[80-1000\,\mu {\rm{m}}]=4\pi {R}_{50\mathrm{IR}}^{2}\sigma {T}_{\mathrm{dust}}^{4}\,,\end{eqnarray} \tag{ 4 }$$

where L_IR is shown in Table 1. ${L}_{\mathrm{IR}}$ and T_dust are measured, and this approach allows for the size expected for a photosphere with ${T}_{\mathrm{dust}}$ to produce ${L}_{\mathrm{IR}}$.

In Figure 7, we compare the sizes measured at 13, 70 μm, and 33 GHz and those calculated for a blackbody (assuming ${A}_{50}^{\mathrm{IR}}=\pi {R}_{50\mathrm{IR}}^{2}$). We see that at least half of the sources in our sample are optically thick at IR wavelengths, with our measured 33 GHz size smaller than the blackbody size. Thus, IR opacity appears significant in our sample, which might be expected considering that we are studying the most obscured systems in the local universe.

Zoom In Zoom Out Reset image size

Our sample consists of gas-rich mergers. In these systems, large masses of gas are funneled to the center, where they become mostly molecular (e.g., Larson et al. 2016). The surface and volume densities of this gas relate closely to its self-gravity and ability to form stars. Again, we assume that the 33 GHz size is characteristic of the system, and by combining this with half of the integrated CO (1−0) measurements, we estimate these quantities for the sample.

Both the assumption of the 33 GHz characteristic size and the conversion between CO luminosity and mass ("conversion factor") introduce uncertainties into the calculation. Our calculation assumes that the gas shares a characteristic size with the star formation traced by the radio. If our targets harbor large amounts of non-star-forming gas or the internal relationship between gas and star formation is strongly nonlinear, e.g., with stars forming much faster in a subset of very dense gas, the calculation will yield biased results. We do expect the approximation to hold, at least to first order. On larger scales star formation traced by IR and CO emission do track one another approximately one-to-one in major mergers (Daddi et al. 2010). Moreover, interferometric CO measurements find that nearly half of the total CO mass is enclosed in the central few kiloparsecs in local U/LIRGs (e.g., Downes & Solomon 1998; Wilson et al. 2008).

For a starburst ${\alpha }_{\mathrm{CO}}=0.8\,{M}_{\odot }\,{({\rm{K}}\mathrm{km}{{\rm{s}}}^{-1}{\mathrm{pc}}^{2})}^{-1}$ (including helium) and coexisting gas and radio emission, we infer values for the molecular gas surface density, ${{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}}$, from ${10}^{2.5}\,\mathrm{to}\,{10}^{5.3}\,{M}_{\odot }\,{\mathrm{pc}}^{-2}$. Even the low end of this range corresponds to source-averaged surface densities in excess of many Local Group molecular clouds (e.g., Bolatto et al. 2008; Fukui & Kawamura 2010). The high end is far in excess of ∼ 1 g cm⁻², which is commonly invoked as an immediate precondition for star formation considering dense substructure inside molecular clouds. Here this gas column density is the average value across the whole energetically dominant area of a galaxy.

These values obviously depend on the mass-to-light ratio adopted to convert CO luminosity to mass. The appropriate conversion factor for starburst galaxies has been a matter of debate, with suggestions ranging from approximately Galactic (e.g., Papadopoulos et al. 2012; Scoville et al. 2014) to low (e.g., Downes & Solomon 1998) and highly environment-dependent (Shetty et al. 2011) values. To see the effect of a higher, Milky Way ${\alpha }_{\mathrm{CO}}$, one should multiply our nominal surface and volume densities by $\approx 5.4$.

Also note that this assumption of matched ${{\rm{\Sigma }}}_{\mathrm{gas}}$ and ${{\rm{\Sigma }}}_{\mathrm{SFR}}$ distributions does not hold for some ULIRGs. For example, for IRAS 13120−5453 the measured starburst size derived from submillimeter continuum is found to be more compact than the emission from dense (Privon et al. 2017) and more diffuse (Sliwa et al. 2017) molecular gas tracers. Moreover, recent high-resolution observations of the CO emission in Arp 220 (Scoville et al. 2017) suggest that the gas is distributed in a larger area compared to the star formation area traced by the 33 GHz emission (see Figure 1 and Barcos-Muñoz et al. 2015). Only $\sim 20 \%$ of the total CO emission is coming from the nuclei. At the moment, Arp 220 is uniquely well studied. These results argue that high-resolution interferometric observations of the gas to match our SFR-tracing continuum will yield important information on how the SFR per unit gas varies across the system. Lacking such information, we proceed assuming matched gas and SFR. If these ULIRGs represent the general case, the reader may think of our ${{\rm{\Sigma }}}_{\mathrm{gas}}$ as an upper limit, with tens of percent of the material in an extended, comparatively non-star-forming disk. This will imply even higher SFR per unit gas mass in the nuclear regions than we calculate below.

One class of models considers the total mass surface density a main driver of the conversion factor, largely via its effect on the line width (e.g., Shetty et al. 2011; Narayanan et al. 2012). Our measured sizes give us the opportunity to illustrate the effect of such a dependence on derived surface densities. To do this, we use the prescription in Bolatto et al. (2013, their Equation (31)), which follows Shetty et al. (2011). Neglecting any metallicity dependence and considering only the regime where ${\rm{\Sigma }}\gt 100$ ${M}_{\odot }$ pc⁻², their prescription is

$$\begin{eqnarray}&&\left(\displaystyle \frac{{\alpha }_{\mathrm{CO}}}{\tfrac{{M}_{\odot }}{({\rm{K}}\,\mathrm{km}\,{{\rm{s}}}^{-1}\,{\mathrm{pc}}^{2})}}\right)\approx 4.35{\left(\displaystyle \frac{{{\rm{\Sigma }}}_{\mathrm{total}}}{100{M}_{\odot }{\mathrm{pc}}^{-2}}\right)}^{-0.5},\end{eqnarray} \tag{ 5 }$$

where ${{\rm{\Sigma }}}_{\mathrm{total}}$ is the total mass surface density driving the potential well. We will assume the systems studied here to be gas dominated in the main CO-emitting region and take ${{\rm{\Sigma }}}_{\mathrm{total}}\sim {{\rm{\Sigma }}}_{\mathrm{gas}}$. The overall gas mass fraction in local U/LIRGs is closer to $\sim 30 \%$ (Larson et al. 2016). However, we expect the gas to be concentrated relative to the stars, so that we can assume the systems to be locally gas dominated in the emitting region. We assume that in the dense, well-shielded central regions of U/LIRGs, the H i content is negligible, and we consider ${{\rm{\Sigma }}}_{\mathrm{gas}}\sim {{\rm{\Sigma }}}_{\mathrm{mol}}$.

We calculate the conversion factor from Equation (5) iteratively, because ${{\rm{\Sigma }}}_{\mathrm{mol}}$ changes as ${\alpha }_{\mathrm{CO}}$ changes. Numerically iterating, we reach a value of ${\alpha }_{\mathrm{CO}}$ that converges to within 0.1%. These values go from 0.2 up to 1.65, with a median value of 0.43 for our sample. We report the gas properties derived using this surface-density-dependent ${\alpha }_{\mathrm{CO}}$ in brackets, along with ${\alpha }_{\mathrm{CO}}$ for each source, in Table 6. The effect of applying this correction is to narrow the range of derived gas surface densities, as the high surface density systems have low ${\alpha }_{\mathrm{CO}}$.

The gas surface density values derived here translate to average hydrogen column densities that range from ${10}^{22.6}\,\mathrm{to}\,{10}^{25.4}\,{\mathrm{cm}}^{-2}$ when using ${\alpha }_{\mathrm{CO}}=0.8$, and ${10}^{22.9}\,\mathrm{to}\,{10}^{24.8}\,{\mathrm{cm}}^{-2}$ when using the surface-density-dependent conversion factor. Assuming a Galactic dust-to-gas ratio (Bohlin et al. 1978), which may be roughly appropriate (Rupke et al. 2008; Iono et al. 2009), these column densities imply line-of-sight extinctions of ${A}_{V}\sim 22$ to 12,000 mag, for a starburst conversion factor, and ${A}_{V}\sim$ 48 to 3200 mag, for a surface-density-dependent conversion factor.

The gas volume density and the corresponding freefall time are central quantities for many theories of star formation (e.g., Krumholz et al. 2012). We estimate the gas volume density from the measured sizes and the integrated CO luminosities. This requires additional geometric assumptions. We consider the most basic approach and assume that our sources are 3D Gaussians. In this case, ∼30% of the mass exists inside the FWHM of the Gaussian,²² ${R}_{50,{\rm{d}}}$.

Adopting this geometry, we find ${n}_{{{\rm{H}}}_{2}}$ from ${10}^{0.5}\,\mathrm{to}\,{10}^{4.9}\,{\mathrm{cm}}^{-3}$ for a fixed starburst ${\alpha }_{\mathrm{CO}}$. Using the variable, surface-density-dependent ${\alpha }_{\mathrm{CO}}$, we find a narrower range of ${10}^{0.9}\,\mathrm{to}\,{10}^{4.3}\,{\mathrm{cm}}^{-3}$. The freefall collapse times associated with these densities range from 6 to 100 (2–190) Myr with the fixed (variable) ${\alpha }_{\mathrm{CO}}$.

In models like those of Condon (1992) and Murphy et al. (2012), the radio SED reflects a mixture of thermal and nonthermal emission. What powers the emission that we observe from U/LIRGs at 33 GHz? In the Condon (1992) model for a starburst galaxy like M82, about 50% of the total 33 GHz continuum is produced by free–free ("thermal") emission; for comparison, $\lt 12 \%$ of the emission is expected to be produced by free–free emission at 1.5 GHz.

We have several constraints on the nature of the emission mechanism in our targets: the SED shape, the brightness temperature, and the comparison with the SFR implied by the IR. Together, these indicate some tens of percent contribution of thermal emission to the 33 GHz flux density, with the balance being synchrotron. However, a detailed understanding of the emission mechanism will need to wait for better coverage of the radio SED in these targets (L. Barcos-Muñoz et al. 2017, in preparation).

Brightness Temperature and Optical Depth: The brightness temperature of optically thick free–free emission is expected to be ∼5 × 10³–10⁴ K. If the 33 GHz T_b exceeded this value, this would provide evidence that synchrotron dominates the emission. Figure 5 shows that the averaged nuclear ${T}_{{\rm{b}}}$ does not exceed this limit. Either the emission is patchy within the beam, or the emission at 33 GHz is optically thin. Thus, the brightness temperature in the sources allows for a normal mix of emission mechanisms and is consistent with optically thin free–free emission making up a large part (or all) of the observed 33 GHz flux density.

If we neglect filling factor effects and assume that $\approx 50 \%$ of the total ${T}_{{\rm{b}}}$ is due to thermal emission, then we can estimate the optical depth of the free–free emission. We derive ${\tau }_{\mathrm{thermal}}\sim {T}_{{\rm{b}}}/{T}_{{\rm{e}}}\leqslant 0.2$ for all our sample. This number is still less than 1, and therefore optically thin, even if we assume that 100% of the 33 GHz flux density is due to thermal emission.

Spectral Index: For a mixture of synchrotron ("nonthermal") emission and optically thin free–free emission, Condon & Yin (1990) give the following approximation to the fraction of emission that is thermal:

$$\begin{eqnarray}&&\displaystyle \frac{S}{{S}_{{\rm{T}}}}\sim 1+10{\left(\displaystyle \frac{\nu }{\mathrm{GHz}}\right)}^{0.1+{\alpha }_{\mathrm{NT}}}\,.\end{eqnarray} \tag{ 6 }$$

Here S is the total flux density, ${S}_{{\rm{T}}}$ is the flux density from thermal emission, and ${\alpha }_{\mathrm{NT}}$ is the typical nonthermal spectral index $\sim -0.8$. The formula assumes a power-law SED for the nonthermal emission.

We combine Equation (6) with the ${S}_{5.95}$ from Table 3 to calculate $S/{S}_{{\rm{T}}}$ at 5.95 GHz. Then, knowing that ${S}_{{\rm{T}}}\propto {\nu }^{-0.1}$, we predict the spectral index between ${\alpha }_{6-33}$. Based on this, we expect an average ${\alpha }_{6-33}=-0.53$. We expect ${\alpha }_{6-33}$ to approach ${\alpha }_{\mathrm{NT}}=-0.8$ as the thermal fraction decreases to zero, while if the thermal fraction is higher than this estimate, ${\alpha }_{6-33}$ will be $\gt -0.53$.

Figure 2 shows that 17 out of the 22 systems in our sample have ${\alpha }_{6-33}\lt -0.53$, implying that in most of our sample, nonthermal emission is stronger relative to thermal emission than predicted by Equation (6). Barcos-Muñoz et al. (2015) found a similar result comparing 6 and 33 GHz emission in Arp 220. We caution that our assumed ${\alpha }_{\mathrm{NT}}$ affects this result and that we cannot, at present, distinguish between variations in the thermal fraction and ${\alpha }_{\mathrm{NT}}$ from only two frequencies. Indeed, multifrequency observations, particularly at high frequency, suggest curvature in the radio SED (e.g., see Clemens et al. 2008, 2010; Leroy et al. 2011; Marvil et al. 2015), so that the power-law assumption for the nonthermal emission model in Equation (6) represents a simplification. Observations that cover a wide band will allow for a more complex treatment for a better disentanglement of the contribution of the two components at these frequencies (S. Linden et al. 2017, in preparation).

Spectral Index and Implied Opacity at Lower Frequencies: Following the same approach, we use Equation (6) and the flux at 5.95 GHz to predict an integrated ${\alpha }_{1.5-6}$ of −0.71. However, less than half of the sample shows spectral indices that agree with this predicted value. Most of our targets show shallower spectral indices. This is most likely due to opacity affecting the low-frequency emission, especially the observations at 1.49 GHz where free–free absorption is known to play a major role in compact starbursts (see, e.g., Condon et al. 1991; Murphy et al. 2013). In fact, in Figure 2 we also observe a change in slope as frequency increases for several sources, from shallower to steeper in most cases. Mrk 231 even shows a change from positive ${\alpha }_{1.5-6}$ to a negative ${\alpha }_{6-33}$. For a compact starburst this would indicate that ${\tau }_{\mathrm{thermal}}$ becomes one at some frequency between 1.5 and 33 GHz.²³ However, we know that Mrk 231 has a very compact core (e.g., Lonsdale et al. 2003; Helmboldt et al. 2007), which suggests instead that the change in slope is most likely due to synchrotron self-absorption at low frequencies. In addition, it is also possible that the flattening in the observed ${\alpha }_{1.5-6}$ could be caused by ionization and bremsstrahlung losses (Thompson et al. 2006; Lacki et al. 2010), which become important at low frequencies in high-density environments such as those found in our sample (see Section 4.3).

Several systems show the opposite trend, exhibiting steep ${\alpha }_{1.5-6}$ and a shallower ${\alpha }_{6-33}$. The simplest explanation for these measurements is that these systems have a higher thermal fraction than the other targets. Alternatively, some other source may contribute to the 33 GHz emission, e.g., anomalous dust emission (Draine & Lazarian 1998; Ali-Haïmoud et al. 2009; Murphy et al. 2010). More detailed SED coverage could confirm this interpretation. Another possible explanation includes contribution from thermal dust, which is normally important only at much higher frequencies, ≳100 GHz. Again, better frequency coverage will play a key role.

In Figure 2, we find a tentative correlation between ${\alpha }_{1.5-6\mathrm{GHz}}$ and ${A}_{50,{\rm{d}}}$, showing a shallower spectral index for more compact sources. This trend makes sense if more compact sources are also more opaque. In this case, 1.5 GHz emission in more opaque systems will be suppressed owing to a higher opacity at 1.5 GHz relative to 6 GHz.

Integrated spectral indices only give us a partial view of the processes that are powering star formation in our sample. We require more detailed spectral index maps to dissect the distribution of the radio emission. We will report resolved spectral index maps between 6 and 33 GHz in a future paper (L. Barcos-Munõz et al. 2017, in preparation). These results will be greatly complemented by spectral indices maps between 1.49 and 8.44 GHz reported in Vardoulaki et al. (2015) using the Condon et al. (1990) and Condon et al. (1991) observations.

Expectations from IR-Based SFRs: The contrast of the 33 GHz flux density with the total IR emission also sheds some light on the emission mechanism. Inasmuch as the IR tells us about the SFR, it also makes a prediction for the expected thermal emission, along with some simplifying assumptions.

We derive the expected free–free emission, ${S}_{{\rm{T}}}$, and then thermal fraction, ${S}_{{\rm{T}}}$/S, at 33 GHz by assuming that all the IR luminosity is due to star formation and that none of the ionizing photons (that will potentially produce free–free emission) are absorbed by dust. Note that if an AGN is present and contributes significantly to the IR luminosity, then the SFR derived by this method will be overestimated (see Armus et al. 2007; Petric et al. 2011, for an estimation of the AGN contribution to L_IR in local U/LIRGs). We use Equation (3) and the thermal SFR from Table 8 in Murphy et al. (2012), which relates SFR and the thermal luminosity, L^T, by the following equation:

$$\begin{eqnarray}\left(\displaystyle \frac{{\mathrm{SFR}}_{\nu }^{{\rm{T}}}}{{M}_{\odot }\,{\mathrm{yr}}^{-1}}\right) & = & 4.6\times {10}^{-28}{\left(\displaystyle \frac{{T}_{{\rm{e}}}}{{10}^{4}{\rm{K}}}\right)}^{-0.45}\\ & & \times {\left(\displaystyle \frac{\nu }{\mathrm{GHz}}\right)}^{0.1}\left(\displaystyle \frac{{L}_{\nu }^{{\rm{T}}}}{\mathrm{erg}\,{{\rm{s}}}^{-1}\,{\mathrm{Hz}}^{-1}}\right),\end{eqnarray} \tag{ 7 }$$

where we assume ${T}_{{\rm{e}}}\sim {10}^{4}$ K (see Sections 4.1 and 4.2 for further discussion on this assumption). In this way, we predict the thermal radio emission expected given the IR luminosity. Comparing it to ${L}_{\mathrm{IR}}$, we derive the thermal fractions shown in the top right panel of Figure 2. We see no clear trend; however, note that ${T}_{{\rm{e}}}$ is uncertain, and the derived thermal fractions depend on it. Lower values of ${T}_{{\rm{e}}}$, or higher thermal optical depths, imply lower thermal fractions. We also observe that Mrk 231 shows the lowest predicted thermal fraction in our sample. This is expected since it does not follow the radio–IR correlation (see Figure 6), with SFR${}_{33\mathrm{GHz}}$ being ∼4 times higher than SFR${}_{\mathrm{IR}}$. By comparing the thermal fractions shown in Figure 2 with the radio–IR correlation shown in Figure 6, we see that all 11 sources with low thermal fraction (i.e., thermal fractions $\lt 60 \%$) are below the equality line in Figure 6. This is consistent with Equation (2) underestimating the SFR due to a more dominant nonthermal component (i.e., a plausible shallower ${\alpha }_{\mathrm{NT}}$) than what is assumed for the equation (−0.8).

From our analysis of the spectral index, we expect thermal fractions ≤50%. Figure 2 shows that, based on the prediction from the IR, most of the sources have thermal fractions ≈50%–100%. We expect that this is the combination of three effects. First, even if the IR is all powered by star formation, some of the ionizing photons produced by young stars that could otherwise produce free–free emission will be absorbed dust and thus not produce free–free emission. These should not be counted in our prediction for the thermal emission, and the true thermal fraction would be accordingly smaller. We highlighted a similar situation in Arp 220, where the predicted thermal fraction is ∼50% but SED analysis shows that it should be closer to 35% (see Barcos-Muñoz et al. 2015). Second, as noted above, the SED-based estimates remain hampered by the lack of sensitive, wide-band coverage of the SED. As long as the adopted nonthermal spectral index (or SED) remains uncertain, so will the thermal fractions estimated in this way. Third, if an AGN contributes a substantial amount to the IR emission, then the thermal fraction would be overestimated because the AGN will not contribute to the free–free emission in the same way as star formation.

Two sources, UGC 04881NE and IRAS 08572+3915, show thermal fractions >100%, meaning that they have very high ratios of IR to radio emission (see Figure 6). This IR excess has been reported before for IRAS 08572+3915 (see discussion in Yun et al. 2004), and this system was already noted as an interesting source in discussion of first results from this survey (Leroy et al. 2011). See the Appendix for further discussion of these two sources.

Radio–FIR Correlation at 33 GHz: As a more observational restatement of the previous result, we derive q₃₃, the ratio of FIR flux (between ∼42 and ∼122 μm) to radio flux density at 33 GHz:

$$\begin{eqnarray}&&{q}_{\nu }={\mathrm{Log}}_{10}(({S}_{\mathrm{FIR}}/3.75\times {10}^{12}\,\mathrm{Hz})/{S}_{\nu })\,.\end{eqnarray} \tag{ 8 }$$

Here ${S}_{\nu }$ is the flux density at frequency ν in units of ${\rm{W}}\,{{\rm{m}}}^{-2}\,{\mathrm{Hz}}^{-1}$, and ${S}_{\mathrm{FIR}}[42-122\,\mu {\rm{m}}]=1.26\times {10}^{-14}\,(2.58\,{S}_{60\mu {\rm{m}}}+{S}_{100\mu {\rm{m}}})$, in units of ${\rm{W}}\,{{\rm{m}}}^{-2}$, is the FIR flux, with the flux density at 60 and 100 μm measured in Jy.

We show a histogram of ${q}_{32.5\mathrm{GHz}}$ in Figure 8. We find a median ${q}_{33}\approx 3.32$ and a dispersion of 0.19 dex. q₃₃ is similar to that found by Rabidoux et al. (2014) studying regions in local star-forming galaxies, but we find a tighter correlation. Their measured dispersion is 0.1 dex larger than ours. In fact, the 0.19 dex in dispersion we observed for q₃₃ is similar to that found in Condon et al. (1991) at 1.49 GHz. The tighter dispersion found for our global measurements compared to the local ones of Rabidoux et al. (2014) appears to corroborate the global nature of the IR–radio correlation. Note as well that ${q}_{32.5\mathrm{GHz}}$ does not appear to correlate with ${{\rm{\Sigma }}}_{\mathrm{SFR}}$.

Zoom In Zoom Out Reset image size

Our size estimates imply that a large part of the star-forming activity, and hence presumably also the gas that fuels it, is concentrated in areas with half-light radii from 30 pc up to 1.7 kpc.²⁴ Applying these sizes to global quantities using the proper aperture corrections, we estimate ${{\rm{\Sigma }}}_{\mathrm{SFR}}$, ${{\rm{\Sigma }}}_{\mathrm{IR}}$, ${N}_{{\rm{H}}}$, and n_mol.

The resulting values span a wide range, typically 4 dex. The high end of the range for each property is among the highest average gas, SFR, or luminosity surface density measured for any galaxy. The low end of the range is still high compared to values found in "normal" disk galaxies: the lowest-density systems have ${{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}}\sim {10}^{2}-{10}^{3}\,{M}_{\odot }$ pc⁻² and ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}\sim {10}^{0}$–10¹ ${M}_{\odot }$ yr⁻¹ kpc⁻². These already resemble the highest kiloparsec-resolution values (which come from active galaxy centers) found in Leroy et al. (2013) (see bottom panel in Figure 10). Moreover, the gas surface densities in our sample, even the lowest values, resemble those found for individual molecular clouds, but here they extend over the whole energetically dominant region of a galaxy. This implies average interstellar gas pressures that match or exceed those found inside individual clouds. Because of this high pressure, a Milky Way GMC dropped into any of the targets would not remain an isolated, self-gravitating object. Self-gravitating, overpressured clouds in these targets must be more extreme and denser than clouds in normal galaxies, a conjecture born out by observations of nearby starburst galaxies (e.g., Keto et al. 2005; Wei et al. 2012; Johnson et al. 2015; Leroy et al. 2015).

About half (13) of the 22 targets studied here show galaxy-averaged ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}\geqslant {10}^{2.7}\,{M}_{\odot }\,{\mathrm{yr}}^{-1}\,{\mathrm{kpc}}^{2}$. This corresponds to ≥2 times higher than the ${{\rm{\Sigma }}}_{\mathrm{SFR}}$ that would be inferred based on the IR emission from the Orion core (Soifer et al. 2000). Several (7) sources show ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}\gt {10}^{3}\,{M}_{\odot }\,{\mathrm{yr}}^{-1}\,{\mathrm{kpc}}^{-2}$, corresponding to ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}\gt {10}^{13}\,{L}_{\odot }\,{\mathrm{kpc}}^{-2}$. This value has been put forward as the characteristic Eddington limit for ${{\rm{\Sigma }}}_{\mathrm{SFR}}$ in a radiation-pressure-supported, optically thick disk (Scoville 2003; Thompson et al. 2005; see Section 5.5 for further discussion).

The high column densities obscure the energetically dominant regions at nonradio wavelengths. Assuming a "starburst" conversion factor, 13 U/LIRGs show hydrogen column densities consistent with being Compton thick, ${N}_{{\rm{H}}}\gt 1.5\times {10}^{24}\,{\mathrm{cm}}^{-2}$ (e.g., Comastri 2004), which would directly affect the ability of X-ray diagnostics to detect the presence of AGNs in these systems. As mentioned above, the implied optical extinctions are extreme, 22−12,000 mag for our sample assuming a Galactic dust-to-gas ratio. Even IR wavelengths, at which a normal star-forming galaxy is usually optically thin, will show significant opacity for these dust columns. At 100 μm, for a mass absorption coefficient of ${\kappa }_{100}=31.3\,{\mathrm{cm}}^{2}\,{{\rm{g}}}^{-1}$ (Li & Draine 2001), the dust opacity of these targets is τ₁₀₀ ∼ 0.02−12, with those same 13 (except for one) Compton-thick sources also being optically thick at 100 μm, i.e., ${\tau }_{100}\gt 1$.

Several studies have reported a "deficit" in the [C ii] 158 μm-to-FIR luminosity (from 40 to 120 μm) ratio, ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{\mathrm{FIR}}$, in U/LIRGs relative to lower-luminosity star-forming galaxies (e.g., Malhotra et al. 2001; Díaz-Santos et al. 2013; Lutz et al. 2016). The ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{\mathrm{FIR}}$ decreases with increasing dust temperature, mid-IR opacity, star formation efficiency (${L}_{\mathrm{IR}}/{M}_{{{\rm{H}}}_{2}}$), and IR surface density (where Spitzer and Herschel data are utilized to measure sizes). The deficit arises because the collisional energy required to produce [C ii] is suppressed in the compact, dense starburst environments of U/LIRGs and/or because the IR luminosity is increased.

The sizes used to gauge the IR surface brightness in Díaz-Santos et al. (2013) come from IR space telescopes, which have much coarser angular resolution than our maps. In Figure 10 (top left panel), we plot ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{\mathrm{FIR}}$ from Díaz-Santos et al. (2013) as a function of the SFR surface density inferred using our sizes, ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}$. The plot shows clear, strong anticorrelation between ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{\mathrm{FIR}}$ and ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}$. The top right panel of Figure 9 shows ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{\mathrm{FIR}}$ as a function of A_50d. Both plots show that more compact systems with more locally intense star formation show stronger ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{\mathrm{FIR}}$ deficits (lower ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{\mathrm{FIR}}$). This is strong corroboration, using the best size measurements to date, of the correlation found by Díaz-Santos et al. (2013) of higher deficit for systems with higher luminosity densities.

Zoom In Zoom Out Reset image size

The spectral index between 1.5 and 6 GHz may give some indication of the opacity at low frequencies. In the bottom panel of Figure 9, we plot ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{\mathrm{FIR}}$ as a function of this spectral index, ${\alpha }_{1.5-6\mathrm{GHz}}$. ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{\mathrm{FIR}}$ is lower, and thus the [C ii] deficit is larger, for systems with flatter (more nearly 0) spectral indices. This flattening is believed to be due to increasing opacity (e.g., see Murphy et al. 2013), so the bottom panel of Figure 9 shows that the ${L}_{[{\rm{C}}{\rm{II}}]}/{L}_{\mathrm{FIR}}$ ratio is lowest for U/LIRGs that are most obscured at radio, as well as IR, wavelengths.

With the exception of IRAS F08572+3915, the five U/LIRGs (Mrk 231, IRAS 15250+3608, III Zw 035, IRAS F01364−1042, and Arp 220) with the flattest ${\alpha }_{1.5-6\mathrm{GHz}}$, and among the largest [C ii] deficit, also have the lowest estimated thermal fraction at 33 GHz in our sample. These results are broadly consistent with our detailed study of Arp 220 (Barcos-Muñoz et al. 2015), where we presented evidence of suppressed 33 GHz thermal emission and speculated that the suppression is due to the absorption of ionizing UV photons by dust concentrated within the H ii regions (see also Luhman et al. 2003; Fischer et al. 2014). Such a scenario would also imply a lack of heating of photodissociated regions and thus a suppression of the amount of collisional energy available to produce [C ii].

The observed scaling between SFR surface density, ${{\rm{\Sigma }}}_{\mathrm{SFR}}$, and gas surface density, ${{\rm{\Sigma }}}_{\mathrm{gas}}$, is often used as a main diagnostic of the physics of star formation in galaxies (e.g., Kennicutt 1998). Kennicutt (1998) fit a scaling between galaxy-averaged ${{\rm{\Sigma }}}_{\mathrm{SFR}}$ and ${{\rm{\Sigma }}}_{\mathrm{gas}}$ that describes both normal disk galaxies and starbursts. The starbursts in Kennicutt (1998) have high ${{\rm{\Sigma }}}_{\mathrm{SFR}}$ and ${{\rm{\Sigma }}}_{\mathrm{gas}}$ and include U/LIRGs like those studied here.

The contrast between the normal disks (low ${{\rm{\Sigma }}}_{\mathrm{SFR}}$, ${{\rm{\Sigma }}}_{\mathrm{gas}}$) and the starburst galaxies (high ${{\rm{\Sigma }}}_{\mathrm{SFR}}$, ${{\rm{\Sigma }}}_{\mathrm{gas}}$) played a main role in driving the best overall fit of Kennicutt (1998), ${{\rm{\Sigma }}}_{\mathrm{SFR}}\sim {{\rm{\Sigma }}}_{\mathrm{gas}}^{1.4}$. This contrast depends on the sizes adopted for the starburst galaxies. Changing the size affects both surface densities by the same factor, but because the overall relationship between ${{\rm{\Sigma }}}_{\mathrm{SFR}}$ and ${{\rm{\Sigma }}}_{\mathrm{gas}}$ is nonlinear, the adopted size affects the slope.

In Figure 10 we place each of our targets in the ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}$–${{\rm{\Sigma }}}_{\mathrm{gas}}^{33\,\mathrm{GHz}}$ (or ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}$–${{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}}$) plane (see Section 4.2 and 4.3 for details on the derivation of ${{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}$ and ${{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}}$). In the top left panel, we show only the U/LIRGs from our sample and adopt a fixed ${\alpha }_{\mathrm{CO}}$ = 0.8 ${M}_{\odot }$ pc⁻² (K km s⁻¹). These U/LIRGs show high surface densities and an approximately linear relationship. A nonlinear least-squares fit²⁵ yields

$$\begin{eqnarray}{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}) & = & (1.02\pm 0.10)\,{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}})\\ & & -(1.33\pm 0.47).\end{eqnarray} \tag{ 9 }$$

This slope is in good agreement with the results found by Liu et al. (2015) for disk galaxies and for U/LIRGs. Genzel et al. (2010) also noted that the internal relationship for starburst galaxies was more nearly linear than the relationship using both types of galaxies, giving rise to the idea of "two sequences" of star formation. A similar conclusion of "two sequences" of star formation is also derived by Daddi et al. (2010), although they obtained a steeper slope (∼1.4) for each type of galaxy that approaches unity within the uncertainty of their measurements. With a slope close to unity, another way to express Equation (9) is that for a "starburst" conversion factor, we find a typical gas depletion time, ${\tau }_{\mathrm{dep}}\equiv {M}_{\mathrm{mol}}/\mathrm{SFR}$, of ${\tau }_{\mathrm{dep}}\sim 20\,\mathrm{Myr}$ for the targets studied here. Note that this short timescale would potentially lead to a relatively flat-spectrum radio source inconsistent with the observed FIR/radio correlation (see Section 5.1); however, the uncertainty in the calculated ${\tau }_{\mathrm{dep}}$ is at least a factor of a few.

Zoom In Zoom Out Reset image size

In addition to the size, the adopted conversion factor can have a large effect on the results. Because we find an approximately linear relationship within our sample, shifting from one constant ${\alpha }_{\mathrm{CO}}$ to another will not affect the slope. For example, if we use a Galactic ${\alpha }_{\mathrm{CO}}$ = 4.35 ${M}_{\odot }$ pc⁻² (K km s⁻¹) instead, the coefficient would shift to −2.08 ± 0.55, raising the depletion time to ${\tau }_{\mathrm{dep}}=125\,\mathrm{Myr}$. For comparison, Leroy et al. (2013) find a significantly longer ${\tau }_{\mathrm{dep}}$, ∼1.6 Gyr, in the disks of nearby normal galaxies.

Several suggestions posit a continuous variation in ${\alpha }_{\mathrm{CO}}$ that depends on surface density (see Equation (5)). Adopting such a prescription affects the slope of the derived relation. If we adopt the surface-density-dependent slope discussed in Section 4.3, the best fit shifts to

$$\begin{eqnarray}{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}) & = & (1.52\pm 0.16)\,{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}})\\ & & -(3.09\pm 0.66).\end{eqnarray} \tag{ 10 }$$

The top right panel in Figure 10 shows our data for two cases: a fixed "starburst" conversion factor and the mass surface-density-dependent value. Internal to the starburst sample, the linearity or nonlinearity of the slope depends entirely on the treatment of the conversion factor and the assumption of the cospatiality between CO and radio emission; the apparent relationship between ${{\rm{\Sigma }}}_{\mathrm{SFR}}$ and CO luminosity surface brightness is approximately linear.

As mentioned above, the contrast between normal disk galaxies and starbursts played a large role in determining the Kennicutt (1998) fit. The bottom panel of Figure 10 explores this contrast. There, we compare our results to those found for kiloparsec-size regions drawn from 30 nearby disk galaxies by Leroy et al. (2013). Individual regions appear as green squares, and the median and scatter in ${{\rm{\Sigma }}}_{\mathrm{SFR}}$, in bins of fixed ${{\rm{\Sigma }}}_{\mathrm{mol}}$, appear as red points with error bars. Note that, in contrast to Kennicutt (1998), we consider only the molecular gas component of the ISM, and, in the normal galaxies, we consider individual kiloparsec-sized regions. Kennicutt (1998) include atomic gas and consider whole-disk averages. We chose our approach to focus on star-forming (molecular) gas in comparably sized regions in order to contrast the ability of gas to form stars in the two types of systems.

Figure 10 shows a significant contrast between disks and our starburst sample, even for matched ${\alpha }_{\mathrm{CO}}$ (a similar contrast was seen when comparing ${\tau }_{\mathrm{dep}}$). In that case, ${\alpha }_{\mathrm{CO}}\,=4.35$ ${M}_{\odot }$ pc⁻² (K km s⁻¹) for both samples, a fit to our sample and the Leroy et al. (2013) bins yields

$$\begin{eqnarray}{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}) & = & (1.35\pm 0.04)\,{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}})\\ & & -(3.85\pm 0.13).\end{eqnarray} \tag{ 11 }$$

Meanwhile, adopting the starburst ${\alpha }_{\mathrm{CO}}=0.8$ ${M}_{\odot }$ pc⁻² (K km s⁻¹) for our sample only yields

$$\begin{eqnarray}{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}) & = & (1.63\pm 0.07)\,{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}})\\ & & -(4.18\pm 0.22).\end{eqnarray} \tag{ 12 }$$

In both cases, the data appear to support the "two sequences" idea, at least to some degree, with internal relationships in the two subsamples that are more nearly linear, and a steep slope when contrasting both populations (but see below). This is particularly the case when we use a starburst conversion factor for our sample.

Adopting ${\alpha }_{\mathrm{CO}}\propto {{\rm{\Sigma }}}_{\mathrm{mol}}^{-0.5}$ (see Equation (5)), we find instead

$$\begin{eqnarray}{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{{\mathrm{SFR}}_{\mathrm{IR}}}^{33\,\mathrm{GHz}}) & = & (1.87\pm 0.06)\,{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}})\\ & & -(4.63\pm 0.19).\end{eqnarray} \tag{ 13 }$$

In this case we find an even steeper slope when fitting the combined data, from the U/LIRGs studied here and the normal spirals from Leroy et al. (2013), than when we use a starburst conversion factor for our sample only, and even more so when we fit either sample alone. To some degree, this reinforces the "two sequences" view, but with a strong caveat. Our results are consistent with the idea that the depletion time is multivalued at a fixed gas surface density, but they do not offer any strong evidence regarding a true bimodality. The data that we use are not complete in any meaningful sense. Therefore, the absence of intermediate ${\tau }_{\mathrm{dep}}$ points near where the two samples would overlap can easily be a selection effect. That is, there may be plenty of parts of galaxies that fill in apparently empty space in Figure 10; our samples are simply not constructed to reveal this. Indeed, Saintonge et al. (2011), Huang & Kauffmann (2015), Genzel et al. (2015), and others have convincingly shown that a continuous range of gas depletion times appear to exist within the population (see also Scoville et al. 2016, for further discussion on continuous and bimodal star formation scaling relations).

Our results do strongly reinforce the idea that the disk–starburst contrast is essential to probe the nonlinear nature of star formation scaling relations. We also show, following a number of others (e.g., see Bouché et al. 2007; Ostriker & Shetty 2011), that the adopted conversion factor, in addition to the starburst sizes, plays a large role in the results. We summarize all the different fits to the gas star formation law using the different conversion factors in Table 7.

Table 7.  Summary of Best Nonlinear Least-squares Fit to Equation ${\mathrm{log}}_{10}({{\rm{\Sigma }}}_{\mathrm{SFR}}^{33\,\mathrm{GHz}})=A\,{\mathrm{log}}_{10}({{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}})+B$

${\alpha }_{\mathrm{CO}}$	Sample Included in Fit	A	B
${M}_{\odot }\,{\mathrm{pc}}^{-2}({\rm{K}}\,\mathrm{km}\,{{\rm{s}}}^{-1})$
0.8	U/LIRGs only (this paper)	1.02 ± 0.10	−1.33 ± 0.47
4.35	U/LIRGs only (this paper)	1.02 ± 0.11	−2.08 ± 0.55
∝${{\rm{\Sigma }}}_{\mathrm{gas}}^{-0.5}$	U/LIRGs only (this paper)	1.52 ± 0.16	−3.09 ± 0.66
4.35	U/LIRGs (this paper) + kiloparsec-size regions (Leroy et al. 2013)	1.35 ± 0.04	−3.85 ± 0.13
0.8a	U/LIRGs (this paper) + kiloparsec-size regions (Leroy et al. 2013)	1.63 ± 0.07	−4.18 ± 0.22
∝${{\rm{\Sigma }}}_{\mathrm{gas}}^{-0.5}$ a	U/LIRGs only (this paper) + kiloparsec-size regions (Leroy et al. 2013)	1.87 ± 0.06	−4.63 ± 0.19
4.35	kiloparsec-size regions (Leroy et al. 2013)	1.06 ± 0.02	−3.48 ± 0.03

Note.

^aOnly applied to the U/LIRGs.

Download table as:  ASCIITypeset image

Efficiency per Freefall Time: A popular class of models posits an approximately fixed fraction of gas converts to stars per gravitational freefall time, ${\tau }_{\mathrm{ff}}^{\mathrm{mol}}=\sqrt{3\pi /(32{\rm{G}}{\rho }_{\mathrm{mol}})}$ (e.g., Krumholz & McKee 2005; Krumholz et al. 2012). If we adopt a simple, spherical (with radius ${R}_{50,{\rm{d}}}$) view of the geometry of the systems, we can estimate ${\tau }_{\mathrm{ff}}^{\mathrm{mol}}$. For a 3D Gaussian, this implies an aperture correction of $\sim 1/3.4$ for the total gas mass (or SFR) within that volume.

Comparing ${\tau }_{\mathrm{ff}}^{\mathrm{mol}}$ to the depletion time of the molecular gas mass, ${\tau }_{\mathrm{dep}}^{\mathrm{mol}}$, we estimate the efficiency of the conversion of the gas mass into stars per freefall time, or ${\epsilon }_{\mathrm{ff}}$ = ${\tau }_{\mathrm{ff}}^{\mathrm{mol}}/{\tau }_{\mathrm{dep}}^{\mathrm{mol}}$. We find median values for ${\tau }_{\mathrm{ff}}^{\mathrm{mol}}$ of 1.1, 1.5, and 0.5 Myr for "starburst," surface-density-dependent, and Galactic conversion factors. These numbers imply median ${\epsilon }_{\mathrm{ff}}$ of 8%, 15%, and 0.6%. The first two numbers appear high compared to the universal $\sim 1 \%$ ${\epsilon }_{\mathrm{ff}}$ assumed in the Krumholz et al. (2012) model, and in more agreement with a nonuniversal star formation efficiency (Semenov et al. 2016), but we emphasize the uncertainty in the adopted geometry.

The high density of star formation and luminosity in the inner parts of our targets undoubtedly creates strong feedback on the gas. This can suppress or even halt ongoing star formation, and in equilibrium we might expect this feedback to counterbalance the force of gravity, leading to some degree of self-regulation. Radiation pressure on dust has been proposed as the main feedback mechanism for compact, optically thick starbursts (Scoville 2003; Murray et al. 2005; Thompson et al. 2005; Andrews & Thompson 2011). Momentum injection by supernova (SN) explosions (Thompson et al. 2005; Kim & Ostriker 2015) and cosmic-ray pressure (e.g., Socrates et al. 2008; Faucher-Giguère et al. 2013) also likely play a key role.

The high ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ values derived for our targets and their very dusty nature make them excellent candidates to be "Eddington-limited" starbursts. In such a system, the star formation surface density will increase until it yields a radiation pressure on dust that balances the force of gravitational collapse. Because we expect that the force from radiation pressure must be present, if a source shows a luminosity surface density above this equilibrium value, then some other assumption in the calculation must break down. This could be the assumption of equilibrium, as the pressure exerted by radiation might temporarily or permanently suppress star formation and/or expel gas from the system in a galactic wind. Alternatively, the source of the luminosity could be something other than star formation. One common inference when this "maximal starburst" case is exceeded is that an appreciable part of the luminosity in the system may arise from an AGN. Alternatively, the assumptions about disk structure used to calculate the force of gravity may be wrong. For example, in the models of Thompson et al. (2005) the gas fraction and velocity dispersion play a key role.

We have already seen some evidence that this case may apply to our systems. Thompson et al. (2005) noted an IR luminosity surface density of ${{\rm{\Sigma }}}_{\mathrm{IR}}\sim {10}^{13}\,{L}_{\odot }\,{\mathrm{kpc}}^{-2}$ as characteristic for dense, optically thick Eddington-limited starbursts. We showed above that a subset of our targets exhibit ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ near, or even above, this limit.

In detail, the exact limiting ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ depends on the detailed structure of the starburst disk, including its size, stellar velocity dispersion (σ), gas mass fraction (${f}_{{\rm{g}}}$), dust-to-gas ratio, and the Rosseland mean opacity (κ) of the system. Thus, the Eddington limit varies from source to source. Taking this into account, we compare our inferred ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}$ (or ${F}_{\mathrm{obs}}$) for each target to the predicted Eddington flux. For hydrostatic equilibrium in a disk, the Eddington flux ${F}_{\mathrm{edd}}$ is

$$\begin{eqnarray}&&{F}_{\mathrm{edd}}=\displaystyle \frac{4\pi {Gc}{\rm{\Sigma }}}{\kappa }\,,\end{eqnarray} \tag{ 14 }$$

where Σ is the surface density of the mass that dominates the gravitational potential involved in the star-forming region and κ is the effective opacity.

The effective opacity depends on the characteristics of the system under study. Following Thompson et al. (2005) and Andrews & Thompson (2011), for systems that are optically thick to the UV radiation, but optically thin to the re-radiated FIR emission, $\kappa (\mathrm{thin})\sim {{\rm{\Sigma }}}_{\mathrm{gas}}^{-1}$. For systems that are optically thick to the reprocessed FIR emission, i.e., when ${{\rm{\Sigma }}}_{\mathrm{gas}}\,\gtrsim \,1\,{\rm{g}}\,{\mathrm{cm}}^{-2}$, $\kappa (\mathrm{thick})\approx {\kappa }_{o}{{\rm{T}}}^{2}$, where T is the temperature of the central star-forming disk and ${\kappa }_{o}\,\approx 2.4\times {10}^{-4}\,{\mathrm{cm}}^{2}\,{\rm{g}}\,{{\rm{K}}}^{-2}$ (Semenov et al. 2003). The transition between regimes is expected to occur when ${{\rm{\Sigma }}}_{\mathrm{gas}}\,\sim \,1\,{\rm{g}}\,{\mathrm{cm}}^{-2}$. Note that in systems without large IR optical depths, the momentum and turbulence from SNe are expected to dominate support of the disk, rather than radiation pressure.

For a Milky Way gas-to-dust ratio and ${\rm{\Sigma }}={{\rm{\Sigma }}}_{\mathrm{mol}}/{f}_{{\rm{g}}}$ a version of Equation (14) that captures all three possible regimes is

$$\begin{eqnarray}&&{F}_{\mathrm{edd}}=\displaystyle \frac{\pi {Gc}{{\rm{\Sigma }}}_{\mathrm{mol}}^{2}}{{{\rm{f}}}_{{\rm{g}}}(1+{\tau }_{\mathrm{IR}}+10{{\rm{n}}}_{\mathrm{mol}}^{-1/7}-\exp [-{\tau }_{\mathrm{UV}}])}\,.\end{eqnarray} \tag{ 15 }$$

Here ${{\rm{\Sigma }}}_{\mathrm{mol}}={{\rm{\Sigma }}}_{\mathrm{mol}}^{33\,\mathrm{GHz}}$ is the gas surface density (see Table 6), and f_g is the gas mass fraction in the core of the galaxy. ${\tau }_{\mathrm{IR}}={\kappa }_{\mathrm{IR}}({\rm{T}}){{\rm{\Sigma }}}_{\mathrm{mol}}/2$ is the IR optical depth and ${\tau }_{\mathrm{UV}}\sim 500\,-1000({\mathrm{cm}}^{2}\,{{\rm{g}}}^{-1})\times {{\rm{\Sigma }}}_{\mathrm{mol}}/2$ the ultraviolet optical depth. We approximate the contribution to support by SNe as $10{n}_{\mathrm{mol}}^{-1/7}$ (see the Appendix of Faucher-Giguère et al. 2013; Kim & Ostriker 2015); the numerical prefactor can vary by a factor of several, up to $\sim 30{n}_{\mathrm{mol}}^{-1/7}$, where ${n}_{\mathrm{mol}}\equiv {n}_{\mathrm{mol}}^{33\,\mathrm{GHz}}$ is the number density of the gas (see Table 6).

In order to derive ${\kappa }_{\mathrm{IR}}({\rm{T}})$, we assume that Equation (40) from Thompson et al. (2005) describes the relation between T, ${T}_{\mathrm{eff}}$, and the vertical IR optical depth. We then solve the implicit equation for T assuming that

$$\begin{eqnarray}\displaystyle \frac{{\kappa }_{\mathrm{IR}}(T)}{({\mathrm{cm}}^{2}\,{{\rm{g}}}^{-1})}=\left\{\begin{array}{ll}2.4\times {10}^{-4}\,{T}^{2}, & {\rm{if}}\,{\rm{}}T{\rm{}}\,{\rm{\lt }}\,{\rm{180}}\,{\rm{K}}\\ 2.4\times {10}^{-4}\,{180}^{2}\approx 7.8, & {\rm{if}}\,{\rm{}}T\geqslant 180\,{\rm{K}}\end{array}\right.,\end{eqnarray} \tag{ 16 }$$

where ${{\rm{\Sigma }}}_{\mathrm{IR}}^{33\,\mathrm{GHz}}={\sigma }_{\mathrm{SB}}{T}_{\mathrm{eff}}^{4}$ and ${\sigma }_{\mathrm{SB}}$ is the Stefan–Boltzmann constant.

In Figure 11 we show the resulting ${F}_{\mathrm{obs}}/({f}_{{\rm{g}}}\ast {F}_{\mathrm{edd}}$ estimated for each galaxy. We draw lines at the Eddington limit where ${F}_{\mathrm{obs}}={F}_{\mathrm{edd}}$ for different gas fractions, and above which the systems are considered super-Eddington. We include the cases for which we consider SN feedback (filled black circles) and where we do not (open red circles).

Zoom In Zoom Out Reset image size

If we assume that the gas fraction in the center of the sources is closer to 1 and neglect SN feedback, we observe five systems that are super-Eddington, including Mrk 231 (UGC 08058), Arp 220, and CGCG 448−020. As noted earlier, Mrk 231 is known to host a strong AGN, and if this drives the IR luminosity, then this Eddington calculation for a starburst disk does not apply. Arp 220 has mid-IR evidence of an energetic AGN—based on low PAH equivalent widths (=0.03−0.17) and/or high 30–15 μm flux density ratios (=10−20) indicative of very warm, Seyfert-like mid-IR dust emission (Stierwalt et al. 2013). Arp 220 is a special case, where the mid-IR diagnostics potentially break owing to the high dust opacity of the system. This is also applicable to CGCG 448−020, for which the source dominating the emission at IR and radio wavelengths (northeast component; see Figure 1 and the Appendix for more information) is highly obscured.

The sources in our sample are extreme starbursts for which SN feedback is most likely important, especially in systems that are more extended and warm (T < 180 K). If we include SN feedback in the calculation of F_edd (filled black circles in Figure 11), we observe that for a gas fraction of 1, 11 out of 22 systems in our sample show super-Eddington values, including the systems mentioned above. Assuming a more conservative gas fraction of 0.3, which is about the system-averaged gas mass fraction based on Larson et al. (2016), we find that five systems are super-Eddington. Note that CGCG 448−020 is a special case since it shows super-Eddington values independent of the gas fraction, indicating the potential presence of an AGN.

We note that our results highly depend on the adopted gas fraction, and while we might expect some funneling of gas to the center to raise the gas fraction to higher values locally, the best way to further improve our constraints is by resolved observations of the disk dynamics, which can yield the total (dynamical) mass, velocity dispersion, and gas mass.