{'fraction_non_alphanumeric': 0.03567972530045263, 'fraction_numerical': 0.02025909161854222, 'mean_word_length': 4.441820961440462, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The major mechanism to drive turbulence in star-forming galaxies

Two competing models, gravitational instability-driven transport and stellar feedback, have been proposed to interpret the high velocity dispersions observed in high-redshift galaxies. We study the major mechanisms to drive the turbulence in star-forming galaxies using a sample of galaxies from the xCOLD GASS survey, selected based on their star-formation rate (SFR) and gas fraction to be in the regime that can best distinguish between the proposed models. We perform Wide Field Spectrograph (WiFeS) integral field spectroscopic (IFS) observations to measure the intrinsic gas velocity dispersions, circular velocities and orbital periods in these galaxies. Comparing the relation between the SFR, velocity dispersion, and gas fraction with predictions of these two theoretical models, we find that our results are most consistent with a model that includes both transport and feedback as drivers of turbulence in the interstellar medium. By contrast, a model where stellar feedback alone drives turbulence under-predicts the observed velocity dispersion in our galaxies, and does not reproduce the observed trend with gas fraction. These observations therefore support the idea that gravitational instability makes a substantial contribution to turbulence in high redshift and high SFR galaxies.


INTRODUCTION
High redshift star-forming galaxies have higher star formation rates (SFRs) than low redshift star-forming galaxies. The SFR density peaked around the redshift range of = 1 − 3 (e.g., Madau et al. 1996;Lilly et al. 1996;Madau & Dickinson 2014). High redshift star-forming galaxies also differ from local counterparts in their small size, high gas fraction and clumpy, thick star-forming discs (e.g. Genzel et al. 2011). One of the most intriguing properties of high-redshift disc galaxies is that their gas velocity dispersions are higher by factors of 2-5 compared to local star-forming galaxies (e.g. Förster Schreiber et al. 2006Schreiber et al. , 2009Cresci et al. 2009;Wisnioski et al. 2011Wisnioski et al. , 2015Übler et al. 2019). These high velocity dispersions likely indicate a highly turbulent ionised interstellar medium (ISM) (e.g. Law et al. 2007;Förster Schreiber et al. 2009;Green et al. 2014;Krumholz & Burkhart 2016;Simons et al. 2017). However, high velocity dispersions do not correlate solely with redshift. Instead, samples of large numbers of galaxies show that the gas velocity dispersion is well correlated with galactic SFR (e.g. Lehnert et al. as other secondary parameters, including gas fraction and circular velocity (Krumholz & Burkhart 2016;Krumholz et al. 2018). Combining the observations with theoretical models, Krumholz & Burkhart (2016) argue that stellar feedback alone cannot explain the high velocity dispersions. They propose that for galaxies with high SFR and high velocity dispersion, the gravity-driven model agrees better with observations. Krumholz et al. (2018) extend this work by proposing a unified model that can link mass transport, star formation fuelling and gravitational instability with the energy and momentum balance of star formation feedback. By comparing their models with collected observation data from the literature, they suggest that transport+feedback can explain the observed velocity dispersions for both local and high galaxies. The model predicts that transport dominates the turbulence of high star forming galaxies and feedback dominates the turbulence of low SFR star forming galaxies.
While these results are suggestive, the data available in the Krumholz & Burkhart (2016) and Krumholz et al. (2018) studies were limited and heterogeneous. It is therefore of interest to carry out higher precision tests using more targeted and homogeneous data sets. There are two main ways to approach this task. One is using only H data, the approach followed by most authors to date (e.g., Johnson et al. 2018;Übler et al. 2019;Varidel et al. 2020). The advantage of this approach is the large sample size one can obtain if only H data are required. The disadvantage is there is no information about molecular and atomic gas, and thus one cannot test any model predictions that depend on gas fraction or other gas properties. The alternative approach that we pursue here is to accept smaller sample sizes, but make use of data on the neutral interstellar gas. The GASS and xCOLD GASS survey provide well-measured molecular gas and atomic gas masses based on CO and H 21 cm observations (Saintonge et al. 2011(Saintonge et al. , 2012Catinella et al. 2012). In this work, we select a representative sample of star-forming galaxies based on their SFRs and molecular and atomic gas masses from the xCOLD GASS survey (Section 2.1). We carry out IFS observations of these star-forming galaxies. Such observations provide robust measurements of velocity dispersion, circular velocity, and orbital period in these galaxies, which can be used to better test these two theoretical models. In Section 2, we describe the details of sample selection, observations, data reduction and data analysis. Section 3 shows discusses the results and their implications. We summarise our main conclusions in Section 4. We use the cosmological parameters 0 = 70 km s −1 Mpc −1 , Ω = 0.3, Ω Λ = 0.7 throughout this paper.

Sample Selection
We draw our sample of galaxies from the xCOLD GASS survey, selecting targets with SFRs and gas fractions in the regime that provides the greatest sensitivity to whether turbulence is driven primarily by gravity or by feedback. These galaxies have typical gas fraction of g 0.3 (where g ≡ g /( g + * ), for gas mass g and stellar mass * ) and SFR in the range 4 M yr −1 to 13 M yr −1 . We focus on this range because Krumholz et al. (2018) show that for this range of gas fraction and SFR, the feedback-only and feedback plus gravity models make very different predictions: the feedback plus gravity model predicts both higher velocity dispersion and a strong relationship between gas fraction and velocity dispersion that is absent in the feedback-only model. Thus both the absolute value of the velocity dispersion and the scaling between velocity dispersion and gas fraction provide means of distinguishing between the models in this regime. We exclude the galaxies with close companions to minimize the impact on the gas turbulence from the merger process. We only include galaxies with declination (Decl.) less than 20 degrees to be able to observe them from Siding Spring Observatory (SSO). These cuts yield a sample of 14 galaxies for which we conduct integrate field observations, of which 7 reach a signalto-noise ratio sufficient for our further analysis (see See Section 2.2 and Section 2.3 for more details) In Figure 1, we show the relationship between gas velocity dispersion ( ) and SFR, with galaxies drawn from the literature (both local and high redshift) shown in the background. The fuchsia foreground data points show the galaxies used in this work. Figure 1 shows that there is a positive trend between and SFR although with a large scatter. A potential interpretation of the large scatter in the versus SFR plane, which we intend to test here, is that it is partly due to the diversity of gas fraction and rotation curve speed among the galaxies plotted.
Our sample galaxies have well measured SFR, molecular gas mass ( H2 ), atomic gas mass ( HI ), and stellar mass ( * ), all of which we take from the xCOLD GASS catalog ; see details in Table 1). The SFRs are measured from WISE + GALEX for galaxies detected in both data sets, and are based on SED modeling for galaxies that lack GALEX detections (Janowiecki et al. 2017). The total molecular gas mass ( H2 ) is derived from the total CO (1-0) line luminosity measured from the IRAM 30 meter telescope (Saintonge et al. 2012. The stellar mass is from the SDSS DR7 MPA/JHU catalog 1 , and is measured through photometry (Kauffmann et al. 2003;Salim et al. 2007). Finally, Catinella et al. (2012) provides measurements of atomic gas mass ( HI ) from observations of the H 21 cm line.

Observations and Data Reduction
From October of 2016 to March of 2017, we observed the 14 galaxies over 6 nights using the WiFeS instrument. WiFeS is an integral field, double-beam, concentric, image-slicing spectrograph mounted on the 2.3-m telescope at SSO. It provides 25 slitlets, each 38 arcsec long and 1 arcsec wide, yielding a 25 × 38 arcsec Field of View (FoV) (Dopita et al. 2007(Dopita et al. , 2010. The galaxies were observed using the B3000 grating on the blue arm and the R7000 grating on the red arm. The typical seeing conditions during the observations were 1.5 to 2 . The data are reduced by the PyWiFeS data reduction pipeline (Childress et al. 2014). The reduced data are reconstructed to three dimensional data cubes. In this study, we only focus our analysis on the red arm, which covers the H emission line. We subtract the sky using the sky spectra extracted from the spaxels without galaxy emission. We also use the sky spectra to characterise the instrument resolution, which is crucial for the velocity dispersion estimation. We measure the velocity dispersion of four skylines that are close to H emission line for each galaxy. The mean value of the measured velocity dispersion is adopted as the instrument resolution for each galaxy, which is about instr ∼ 17 km s −1 , consistent with spectral resolution = /Δ in the WiFeS manual. We summarise the instrument resolution ( instr ) measurements for each galaxy in Table 2. ,and The HI of nearby galaxies are from THINGS Leroy et al. 2008;Ianjamasimanana et al. 2012) and Stilp et al. (2013). The data on nearby ULIRGs are from Downes & Solomon (1998); Sanders et al. (2003); Veilleux et al. (2009);Scoville et al. (2015Scoville et al. ( , 2017. The data on high-redshift galaxies are from Law et al. (2009)

Spectral Analysis
Among the 14 galaxies selected in Section 2.1, we discard one that is classified as an AGN in the xCOLD GASS catalogue , and analyse the WiFeS data cubes for the remaining 13 galaxies.
We study galaxy kinematics using the H emission line. First, the spectral continuum is subtracted using spectra at the rest-frame wavelength range between 6500 Å and 6530 Å. Then we fit the H , [N ] 6549 and [N ] 6583 emission lines simultaneously using three Gaussian profiles in each spaxel. In this process, we obtain the H flux, the line of sight velocity and observed velocity dispersion ( obs ) for each spaxel ( Figure 2).
We compute the intrinsic gas velocity dispersion for each spaxel ( pix ) from the observed ( obs ) velocity dispersion by correcting the instrument broadening as follows: pix = ( 2 obs − 2 instr ) 1/2 . To have a reliable (S/N > 5) intrinsic gas velocity dispersion measurement for each spaxel, we require that the S/N of the H be greater than 20 -see the detailed discussion in Zhou et al. (2017). Seven galaxies in our sample have S/N ratios in the H line that are high enough for further analysis.

Gas velocity dispersions
We measure the typical gas velocity dispersion for each galaxy by minimizing the "beam smearing" effect. The "beam smearing" effect blends the information from one spaxel with its neighbouring spaxels in the IFU observations due to the limited spatial resolution. Following the method proposed by Wisnioski et al. (2015), we measure velocity dispersion along the major kinematic axis from the outer regions, where the "beam smearing" effect is negligible. We first extract the observed velocity dispersion ( obs ) for each galaxy along the major kinematic axis from the observed velocity dispersion map (see Section 2.3.2 for details on how we derive the major kinematic axis). We then estimate the typical obs of our galaxies by averaging the velocity dispersion in the region of the flat rotation curve (red dots in the fifth column of Figure 2). Last, we remove the instrument resolution ( instr ) from the observed H gas velocity dispersion ( obs ) and obtain the intrinsic ionised gas velocity dispersion ( ) through = ( 2 obs − 2 instr ) 1/2 (Table 2). Before proceeding, we must make one final correction. The measured velocity dispersions in this work are from H emission lines, which are dominated by the gas from H regions. This gas has a thermal velocity dispersion th = √︁ / H , where is the gas temperature and is the mean mass of free particles, normalised to the hydrogen mass H . This thermal velocity will be added in quadrature to the bulk velocity dispersion of the neutral ISM ( 0 ) that is of interest for testing theoretical models (Krumholz & Burkhart 2016;Krumholz et al. 2018). The characteristic temperature of H regions is ≈ 10 4 K (Andrews & Martini 2013), and for fully ionised gas that is 73% H and 25% He by mass, = 0.61, so this corresponds to a typical thermal broadening th ≈ 12 km s −1 (Zhou et al. 2017) 2 . We therefore estimate the bulk ISM gas velocity dispersion 0 as 2 0 = 2 − 2 th . We use 0 for all the analysis we present below; however, we note that, given the relatively large values of reported in Table 2, the difference between and 0 is in all cases 20%.

Circular Velocity and Orbital Period
Besides the SFR, the velocity dispersion is also predicted to depend on the circular velocity and orbital period (Krumholz et al. 2018). Therefore, we also calculate these two quantities from the WiFeS IFU observations. To do so, we measure the rotation curve along the major axis using the H velocity map. First, we estimate the kinematic position angle (PA) of gas and determine the major kinematic axis of each galaxy using the 'FIT KINEMATIC PA' routine (Krajnović et al. 2006). The kinematic PA of gas is defined as the counter-clockwise angle between north and a line that bisects the velocity field of the gas, measured on the receding side. The first and second columns of Figure 2 show the resulting fits for direction of the major kinematic axis. Then we measure the line of sight velocity along the major kinematic axis and correct the velocity for the inclination of each galaxy, using the inclination angle given in the xCOLD GASS catalog ). The fourth column of Figure 2 shows the resulting rotation curves along the major kinematic axis. We estimate the circular velocity ( c ) from the flat part of the rotation curves, which we indicate by the red points in the fourth column of Figure 2. We adopt the largest value from the red dots as the circular velocity for each galaxy. We summarise the resulting values for c in Table 2. We find that the mean circular velocity of our sample of galaxies is 191 km s −1 .
We also measure the orbital period using the rotation curves shown in Figure 2. We first define a typical radius, half-light radius (r Re ), so that to estimate the orbital period ( orb ) for each galaxy. We used SDSS r-band image 3 to run GALFIT (Peng et al. 2002(Peng et al. , 2010 and obtain the half-light radius for each galaxy. One Sérsic model is enough to fit the surface brightness profile for each galaxy. The r Re are listed in Table 2. From the forth column of Figure 2, we find that the circular velocity corresponding to half-light radius larger Kinematic maps and axis profiles for the seven disc galaxies with S/N > 20 that constitute our final sample. From left to right we show for each galaxy: observed H velocity dispersion map, H velocity field (corrected the inclination angle), H emission map, velocity and velocity dispersion profile along the major kinematic axis, which is shown by the white solid lines over plotted on the velocity dispersion and velocity maps. The white cross over-plotted on the velocity dispersion and velocity maps are the data points from which we extract the velocity and velocity dispersion profiles. The blue dashed lines in the right two columns indicate the radii beyond which we deem the rotation curve to be flat; points outside the radius are shown in red, and we use these points to measure the circular velocity ( ) and the observed velocity dispersion ( ). The red dashed lines indicate the half-light radius(r Re ) we measured in section 2.3.2. than the rotation curve turn-flatten radius. So when we estimate the orbital period ( orb ), we used r Re and c to estimate the orb . We show our estimated values of orb in Table 2. The mean value of orb for our galaxies is about 211 Myr. Krumholz et al. (2018) propose a galactic disc model based on the requirement that the disc be in both vertical hydrostatic and energy equilibrium. In this model, the sources of energy input include both star formation feedback and gravitational potential energy released   by inward flow of gas through the disc. They show that these requirements lead to a relationship between star formation rate and gas velocity dispersion

Theoretical models
where is rotation curve index, is the offset between resolved and unresolved star formation law, is the stability parameter (Romeo & Wiegert 2011;Romeo & Falstad 2013), mp is the ratio of total pressure to turbulent pressure at mid-plane, , and , are the fractional contributions of gas self-gravity to mid-plane pressure and gravitational stability, respectively, ff is the star formation efficiency per free-fall time, sf is the fraction of the ISM in the star-forming (molecular) phase, and orb are the circular velocity and orbital period at the outer edge of the star-forming disc, and sf,max = 2 Gyr is the maximum star formation timescale. Eq. 1 is the prediction for a model in which both mass transport and star formation feedback contribute to energy equilibrium, and it applies for galaxies where the gas velocity dispersion is larger than sf , the maximum velocity dispersion that can be sustained by star formation feedback alone. On the other hand, one can also use the same framework for a model in which one assumes that star formation feedback alone contributes to turbulence, and thus sf = 0 by assumption; in order for this to occur, one must allow ff to vary freely. Doing so yields the pure feedback model (in the terminology of Krumholz et al. 2018 Here is a factor of order unity that measures the fraction of turbulent energy dissipated per crossing time, nt is the fraction of velocity dispersion that is non-thermal (close to unity for all our galaxies), is the one plus ratio of gas to stellar , and * / * is the momentum injected per unit mass of stars formed.
More details about these models are in the work of Krumholz et al. (2018). We emphasise that, although we use the precise formulations of these models derived by Krumholz et al., the predictions for the scaling between 0 , SFR, , and gas fraction derived there are generic to all models in the same class; for example, the scalings in Eq. 2 are identical to those obtained by, e.g., Faucher-Giguère et al. (2013) or Orr et al. (2020), who also present models where feedback is the sole source of turbuelence. We can therefore use the scalings in Eq. 1 and 2 as generic tests of all models of that class.

RESULTS AND DISCUSSION
In this section, we compare our observational results with the theoretical models described in Section 2.4 and discuss the results.

Measurement of the model input parameters
Our IFU observations provide straightforward measurements of some of the parameters that enter the two theoretical models presented in Section 2.4: the gas velocity dispersion ( 0 ), orbital period  The relationship between 2 and the ratio of SFR to the square of velocity dispersion (SFR/ 2 ). The red solid lines represent the pure feedback model. Both left and right panels, the shade regions represent the circular velocity range from = 153 − 248 km s −1 , which is the typical range for our galaxies. The we use in this work has been subtracted the thermal broadening in quadrature of HII region, about 12 km s −1 .
( orb ), and circular velocity ( ). Similarly, the xCOLD GASS survey provides the SFR and gas content in these galaxies, which allow us to estimate the gas fraction in the star-forming molecular phase as sf = ( H2 )/( H2 + HI ). We report our measured values of sf in Table 1; the mean value of sf is about 0.35.
Other parameters appearing in the theoretical models require somewhat more indirect estimation procedures. Following Krumholz et al. (2018), we estimate , ≈ , , and we compute , , the effective gas fraction in the disc for the purpose of computing gravitational stability, from its definition where Σ is gas surface density, Σ * is stellar surface density, is the intrinsic gas velocity dispersion, and * is the stellar velocity dispersion. We lack direct measurements of * , and we therefore adopt the ratio * / 0 = 2.3 measured in the Milky Way (Kalberla & Kerp 2009;McKee et al. 2015;Krumholz et al. 2018). This ratio may well be somewhat smaller in our more gas-rich galaxies, which presumably have dynamically "hotter" gas discs. However, we have verified that even using an extreme value * / 0 = 1 (i.e., gas and stars have the same velocity dispersion) does not make a substantial difference to our results. For this reason we simply adopt the Milky Way value for the remainder of the paper. Similarly, we do not have independent measures of the stellar and gas areas, and we therefore adopt Σ g /Σ * = ( H2 + HI )/ * . Table 1 shows our measured values of , . The mean value is about 0.60, slightly higher than for the Milky Way ( , ≈ 0.5 - Krumholz et al. 2018), which is not surprising given that we have intentionally selected more gas-rich galaxies.
Similarly, we estimate the stability parameters for our discs from where = 0 / Σ . Here, = √︁ 2( + 1)Ω is the epicyclic frequency and Ω = 2 / orb is the galaxy angular velocity (Krumholz et al. 2018). We can measure Ω directly from the orbital period, and all our galaxies have relatively flat rotation curves outside of their central regions, so we adopt = 0, which provides an estimate of . The gas surface density, Σ , requires somewhat more care because xCOLD GASS does not resolve the target galaxies. For the H component, we rely on the well-established observational and theoretical result galactic H surface densities never exceed ≈ 10 M pc −2 , because higher gas column densities lead to the formation of self-shielded regions that convert to H 2 (Krumholz et al. 2008(Krumholz et al. , 2009Kalberla & Kerp 2009;McKee & Krumholz 2010;Wong et al. 2013). Since our galaxies are gas-rich and more rapidly starforming than the Milky Way, we expect them to lie near this limit, and adopt Σ HI = 10 M pc −2 as the surface density of the atomic component. For the molecular gas, we assume that the H emission traces star-forming molecular gas, consistent with the observation that H and CO extremely well-correlated when averaged on kpc scales (e.g., Leroy et al. 2008Leroy et al. , 2013. We therefore use the resolved H emission as a proxy for the radius of the unresolved CO. We compute a third-order polynomial to fit the spatial H flux as a function of radius, and use this fit to obtain the half-light radius of ionized gas emission ( eff,H ), which we list in the Table 2, and derive the surface density of H 2 from H2 / 2 eff,H . Inserting these values into Eq. 3 yields estimates of for each galaxy, which we report Table 2. The mean value of Q is 0.77. The final parameter we estimate from observations is where * = * / Σ * . We can estimate * , and thence , using the same assumed ratio * / 0 adopted above. Re-arranging the definitions of , * , and , , we have and plugging our estimates of (Krumholz et al. 2018;Romeo & Wiegert 2011;Romeo & Falstad 2013), = , , and * / 0 into this expression yields * , and then subsituting into Eq. 5 gives . We report our estimate for this parameter in Table 2. The mean value is is 2.03. For all the remaining parameters that enter Eq. 1 and Eq. 2, we adopt the fiducial values of Krumholz et al. (2018). We list these for convenience in Table 3.

Comparison of theoretical models to observations
We are now in a position to compare our observational results with the theoretical models discussed above: the trans-port+feedback model (Eq. 1) that includes both gravitationallydriven and feedback-driven contributions, and the feedback-only model (Eq. 2). Figure 3 shows the relationship between SFR and the intrinsic gas velocity dispersion. The solid blue line and solid red lines show the transport+feedback model and pure feedback model, respectively. For the purposes of generating the theoretical lines, we adopt the mean values of measured parameters, which are sf = 0.35, = 191 km s −1 , orb = 211 Myr, = 0.77, and , = , = 0.60 (see Table 1, Table 2, and Section 2.3.2). The blue and red shaded regions represent the range from = 153 − 248 km s −1 , which is the range of measured circular velocity in our galaxies. We find that the transport+feedback model prediction is in better agreement with the observations for our galaxies. The stellar feedback only model significantly underestimates the velocity dispersion for a given SFR.
Because we have individual measurements for the gas fractions, we can use also these values instead of the means, thereby performing a test of the "secondary" dependence of the SFR on gas fraction (as opposed to the "primary" dependence on velocity dispersion). To this end, we re-arrange Eq. 1 so as to isolate the dependence on gas fraction. We therefore define the quantity which captures the dependence of the model prediction on gas fraction (and Toomre ). Examining Eq. 1, we can see that, in the transport+feedback model, we expect SFR/ ∝ 1 . Performing an analogous operation for the pure feedback model (Eq. 2), we find that the gas fraction dependence is captured by the quantity Again, comparison to Eq. 2 shows that, in the pure feedback model, we expect SFR/ 2 ∝ 2 . We plot these two predictions against the data in the left and right panels of Figure 4, respectively. First examining the transport+feedback model, we find that, given the observational uncertainties, the model prediction agrees well with the data, both in normalisation and in trend. We emphasise that this is a non-trivial result: the model predicts that, at fixed star formation rate, more gas-rich galaxies (higher 1 ) tend to have smaller velocity dispersion, and the data bear out this prediction. By contrast the feedback-only model in the right panel does not match the data in normalisation, which is not surprising, since we have already seen this effect in Figure 3. While the trend in 2 with SFR/ 2 is qualitatively consistent with the prediction within the errors, the observed galaxies are much more gas-rich (higher 2 ) at a given value of SFR/ 2 than would be expected for the feedback-only model.

Sample selection
As discussed in Section 2.1, our sample selection is largely driven by the need to explore parts of parameter space where feedbackonly and feedback-plus-gravity models make measurably different predictions, as well as where we also have access to data on the gas content of galaxies. These two goals are in tension, because the predictions of the two models differ most strongly for high-SFR, highly-turbulent galaxies of the type found primarily at high-redshift, but for such galaxies we have little information on their molecular content and essentially none on their H content. Conversely, our knowledge of gas properties is the richest for nearby galaxies, but the vast majority of these fall into the low-SFR, low-turbulence regime where the feedback-only and feedback-plusgravity models make nearly identical predictions. As a result of these competing imperatives, our sample winds up covering a fairly narrow range of SFR and velocity dispersion. While our exploration of this regime favours the feedback-plus-gravity model, we cannot from this sample alone rule out the possibility that a feedback-only model might provide a better fit to data in other parts of parameter space, for example at SFRs 30 M yr −1 or velocity dispersions 40 km s −1 , larger than those we have been able to explore.

Bias in H velocity dispersions
A second caveat to our work concerns our use of H , which traces ionised gas, as our tool for measuring the velocity dispersion. In deriving the velocity dispersion of the dominant neutral ISM component, we have corrected for thermal broadening, but this correction is relatively small (≈ 20%) for the range of velocity dispersions found in our sample galaxies. Our approach is consistent with the observational results of Übler et al. (2018) and Girard et al. (2019), who find that the ionised gas velocity dispersion differs little from the molecular gas velocity dispersion in galaxies where the total velocity dispersion 20 km s −1 . More recently, however, Girard et al. (2021) obtained much larger differences between molecular and ionised gas velocity dispersions, a factor of ≈ 2.5, in a sample of 9 local galaxies, several with velocity dispersions similar to the range we are exploring, 20 km s −1 . These results are inconsistent with the previous findings in the literature, and the origin of the discrepancy is unknown. Further complicating this story, recent simulations by Kretschmer et al. (2021) suggest that galaxies with large velocity dispersions can show a significant systematic offset in velocity dispersion between molecular and atomic phases of the ISM, such that the ionised gas velocity disperison is actually closer than the molecular gas value to the velocity dispersion of the neutral ISM as a whole. Given the divergent and contradictory results in the literature, we have elected to adopt a minimal approach of not making a larger correction to the H velocity dispersion. However, it is important to consider how such a larger correction would affect our results. As an extreme, adopting a simple factor of 2.5 reduction to all our velocity dispersions, consistent with Girard et al. (2021), would have the effect of bringing our measurements much more closely into alignment with the feedback-only model, and pushing them further from the feedback-plus-gravity one. However, such a change would also have other effects; for example, such a correction would reduce our estimated values of Toomre to ≈ 0.3, which would likely invalidate the feedback-only model as well, since even in this model it is assumed that galaxies self-regulate to ≈ 1.

CONCLUSIONS
We summarise our main results are as follows: (i) We select a sample of star-forming galaxies in the xCOLD GASS survey, chosen to lie in a range of star formation rate and the gas fraction that offers maximum power to distinguish between alternative models for the origin of turbulence in the ISM: on where stellar feedback alone provides the required energy, and where both stellar and gravitational instability-driven transport of mass down the potential well supply energy to the turbulence.
(ii) We carry out IFU observations of these galaxies with the WiFeS instrument, obtaining robust measurements of the gas velocity dispersion, circular velocity, and orbital period in our galaxies. These parameters, together the gas fractions taken from the xCOLD GASS survey, provide the inputs required by the competing theoretical models.
(iii) We compare the relation between SFR, velocity dispersion, and gas fraction with the two models. We find that our results are consistent with the transport+feedback model. In particular, not only does this model match the mean relationship between velocity dispersion and star formation rate seen in the data, it reproduces the secondary trend with gas fraction. By contrast, the feedback-only model both underpredicts the velocity dispersion and fails to match the dependence on gas fraction.
Overall, we find that the transport+feedback model is a better match to the observations. This suggests that stellar feedback only cannot provide the energy to support the high velocity dispersion found in our galaxies. Our results suggest that the gravitational instability dominates the high turbulence of high and high SFR star-forming galaxies.