11institutetext: Université Paris-Saclay, Université Paris Cité, CEA, CNRS,
AIM, 91191, Gif-sur-Yvette, France 22institutetext: Universität Heidelberg,
Zentrum für Astronomie, Institut für Theoretische Astrophysik, Albert-Ueberle-
Str 2, D-69120 Heidelberg, Germany 33institutetext: ENS Lyon, 15 parvis René
Descartes, BP 7000 69342 Lyon Cedex 07, France

# Large-scale turbulent driving regulates star formation in high-redshift gas-
rich galaxies II: Influence of the magnetic field and the turbulent
compressive fraction

Noé Brucy 1122 Patrick Hennebelle 11 Tine Colman 11 Simon Iteanu 1133

(Received 8 September 2022 / Accepted 28 May 2023)

###### Abstract

Context. The observed star formation rate (SFR) in galaxies is well below what
it should be if gravitational collapse alone were at play. There is still no
consensus about the main process that cause the regulation of the SFR.

Aims. It has recently been shown that one candidate that might regulate star
formation, the feedback from massive stars, is suitable only if the mean
column density at the kiloparsec scale is lower than $\approx
20\leavevmode\nobreak\ \text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$. On the
other hand, intense large-scale turbulent driving might slow down star
formation in high-density environments to values that are compatible with
observations. In this work, we explore the effect of the nature and strength
of the turbulent driving, as well as the effect of the magnetic field.

Methods. We performed a large series of feedback-regulated numerical
simulations of the interstellar medium in which bidimensional large-scale
turbulent driving was also applied. We determined the driving intensity needed
to reproduce the Schmidt-Kennicutt relation for several gas column densities,
magnetization, and driving compressibility.

Results. We confirm that in the absence of turbulent forcing and even with a
substantial magnetic field, the SFR is too high, particularly at a high column
density, compared to the Schmidt-Kennicutt relation. We find that the SFR
outcome strongly depends on the initial magnetic field and on the
compressibility of the turbulent driving. As a consequence, a higher magnetic
field in high column density environment may lower the energy necessary to
sustain a turbulence that is sufficiently intense to regulate star formation.

Conclusions. Stellar feedback does not seem to be sufficient to regulate star
formation in gas-rich galaxies where large-scale turbulent driving may be
needed. The sources of this large-scale turbulence as well as its
characteristics, such as its intensity, compressibility, and anisotropy, need
to be understood and quantified.

###### Key Words.:

Galaxies: kinematics and dynamics - ISM: general - ISM: magnetic fields -
stars: formation - Magnetohydrodynamics (MHD) - Radiative transfer

##  1 Introduction

Stars are complex objects that connect small- and large-scale physics. Key
characteristics of galaxies, such as their structure and shape, the amount of
gas, their chemistry, and the velocity field, depend on the efficiency and on
the rate of star formation. However, understanding and predicting the star
formation rate (SFR) remains both an observational and theoretical challenge.
On the theoretical side, the difficulty is the very broad span of scales
involved that cannot be captured in just one simulation (Dubois & Teyssier,
2008; Hopkins et al., 2011), and the high degeneracy between the multiple
processes at play.

A considerable amount of work has been done to understand the interstellar
medium (ISM) at the kiloparsec scale, which is an obligatory step on the path
towards determining how the various physical processes are combined and
influence the structure of the ISM. Since modelling a full galaxy with parsec
or sub-parsec resolution is very challenging, many studies have considered
computational boxes with a size of about 1 kpc. The advantage is that this
allows describing the disk stratification and the large molecular cloud
complexes while still resolving the clouds sufficiently well.

The complexity has been progressively increased in these simulations. The
first series of calculations (e.g. Rosen & Bregman, 1995; de Avillez &
Breitschwerdt, 2005; Joung & Mac Low, 2006; Hill et al., 2012) solved the
equations of Magnetohydrodynamics with an interstellar cooling curve in the
presence of an external gravitational potential representing the influence of
stars and dark matter. In these calculations, supernova (SN) explosions were
placed randomly at a rate aiming to reproduce the expected galactic SN rate.
These models produced appealing turbulent multi-phase ISM structures with a
velocity dispersion close to the observations. In a second series of
calculations (e.g. Kim et al., 2011, 2013; Hennebelle & Iffrig, 2014; Walch et
al., 2015; Iffrig & Hennebelle, 2017; Gatto et al., 2017), self-gravity was
introduced, together with sink particles that were to represent stellar
clusters. This important step allowed addressing the question of the SFR and
also the issue of the galactic thickness, which requires the correct balance
between gravity and turbulence. In these calculations, the stellar feedback
essentially due to SNe was spatially correlated to the sink particles, but no
time delay of a few to some tens of million years was introduced. That is to
say, the feedback was immediately delivered when enough mass was accreted to
form a massive star. Whereas these models could successively reproduce the
star formation rate typical of Milky Way-type galaxies, it was noted that the
scheme used for the SN explosions had a strong influence on the results. In
particular, if the SN exploded far away from its parent sink particle or
simply was forced to explode in the diffuse gas, it would become very
inefficient.

In the last generation of models (e.g. Padoan et al., 2016; Kim & Ostriker,
2017; Colling et al., 2018; Kannan et al., 2020; Ostriker & Kim, 2022; Rathjen
et al., 2023), the massive stars were tracked more self-consistently and the
SNe explode after a few millions to tens of million years, depending on the
mass of their progenitor. In combination with other processes such as stellar
wind, HII regions, or galactic shear, it has been concluded by various teams
that the SFR was compatible with Milky Way values.

A fundamental assumption made in these feedback-regulated models is that the
ISM regulation is relatively local. More precisely, the energy injected in the
medium at a given point is entirely due to stellar feedback, and it typically
comes from less than 1 kpc. However, several recent studies have shown that
injection of turbulence from large galactic scales has to be taken into
account in order to explain the observed velocity dispersion and SFR,
particularly in gas-rich galaxies (Bournaud et al., 2010; Renaud et al., 2012;
Goldbaum et al., 2015, 2016; Krumholz et al., 2018; Meidt et al., 2020; Nusser
& Silk, 2022). There is also evidence that large-scale driving is needed to
reproduce the statistics of the dense gas, in particular, the large-scale
structures in the Large Magellanic Cloud (e.g. Colman et al., 2022). Several
observational studies have indeed derived turbulence-injection scales of more
than 1 kpc (Dib et al., 2021; Chepurnov et al., 2015) and indicated large-
scale driving (Szotkowski et al., 2019; Besserglik & Goldman, 2021). Possible
sources of turbulence include the orbital energy (Wada et al., 2002; Agertz et
al., 2009) and the mass accretion onto the galaxies (Klessen & Hennebelle,
2010; Forbes et al., 2023). The former in particular requires a mechanism such
as an instability to convert this source of free energy into turbulent energy,
which might be the gravitational instability, for instance (Wada et al., 2002;
Agertz et al., 2009; Goldbaum et al., 2015; Fensch et al., 2023), or the
magneto-rotational instability (Piontek & Ostriker, 2007).

An important challenge of numerical simulations for the SFR in particular is
reproducing the Schmidt-Kennicutt relation (hereafter SK relation; Kennicutt,
1998; Kennicutt & Evans, 2012; Kennicutt & De Los Reyes, 2021) that links the
SFR to the column density of gas. In a recent paper, Brucy et al. (2020)
(hereafter Paper I, ) tested the effect of the injection of turbulence by
adding a large-scale turbulent driving similar to the one used by Schmidt et
al. (2009). Importantly, Brucy et al. (2020) have shown that whereas the SK
relation can be reproduced if the added turbulence is strong enough, the SFR
appears to be too high, particularly for gas-rich galaxies, when stellar
feedback alone is accounted for.

The goal of the present paper is to complement the study presented in Paper I
by studying the effects of the magnetic field strength and the compressive
fraction of the turbulent driving on the SFR. Both aspects are known to
significantly influence the gas distribution (e.g. Molina et al., 2012) and
likely influence the SFR values. We ran simulations of a local region of a
galactic disk within a cubic box of 1 kpc. We used a numerical setup very
similar to the one used by Colling et al. (2018) and Paper I. We simulated
regions of galaxies with a wide range of gas column densities that are
representative for Milky Way-like galaxies up to gas-rich galaxies at redshift
$z=$ 1–3 (Genzel et al., 2008, 2010; Daddi et al., 2010).

The paper is organised as follows: We start by discussing the orders of
magnitude of the turbulent energy that can be dissipated by various processes
in the ISM in section 2. Next, we present our numerical setup in section 3,
and we explore the effects of the turbulence and the magnetic field in
sections 4 and 5, respectively. We then discuss the main goal of this study,
that is, the reproduction of the SK relation in section 6. The caveats of this
study are discussed in section 7. Finally, we discuss the results and conclude
in section 8.

##  2 Energy dissipation and turbulent driving in the ISM, orders of magnitude

As discussed above, we further investigate whether large-scale turbulent
driving can help explain the SK relation. It is thus of primary importance to
discuss the possible sources of turbulent driving in galaxies. Below, we
discuss three mechanisms, namely gravitational instability (Wada et al., 2002;
Agertz et al., 2009), accretion-driven turbulence (Klessen & Hennebelle, 2010;
Forbes et al., 2023), and SN driving (e.g. Kim & Ostriker, 2017; Gatto et al.,
2017). We estimate the amount of energy $P$ that can be extracted per unit of
time through each of the mechanisms, as well as their efficiency $\epsilon$,
which is the fraction of that energy that is eventually converted into
turbulent kinetic energy.

We start by estimating the amount of turbulent energy that is dissipated per
unit of time in a region of size $L\simeq$ 1 kpc within a galactic disk, with
a column density $\Sigma=10-100$ M⊙ pc-2 and a velocity dispersion
$\sigma\simeq 10-100$ km s-1. Since the mass contained in a surface $L^{2}$ is
$\Sigma L^{2}$, whereas the specific kinetic energy is $\sigma^{2}/2$ and the
crossing time is $L/\sigma$, we obtain

$\displaystyle P_{\mathrm{diss}}$ $\displaystyle\simeq\Sigma
L^{2}\dfrac{\sigma^{3}}{2L},$ (1) $\displaystyle\simeq 3\times
10^{40}\mathrm{erg\,s^{-1}}\left(\dfrac{\Sigma}{100\,\mathrm{M_{\odot}pc^{-2}}}\right)\left(\dfrac{\sigma}{100\,\mathrm{km\,s^{-1}}}\right)^{3}\left(\dfrac{L}{\mathrm{1\,kpc}}\right).$

The largest amount of energy that can be extracted per unit of time from
galactic rotation is typically equal to the rotation energy
$V_{\mathrm{rot}}^{2}/2$, divided by a rotation time that is simply given
$2\pi R_{\mathrm{gal}}/V_{\mathrm{rot}}$ , where $R_{\mathrm{gal}}$ is the
galactic radius. This leads to

$\displaystyle P_{\mathrm{rot}}\simeq$ $\displaystyle\;\Sigma
L^{2}\dfrac{V_{\mathrm{rot}}^{3}}{4\pi R_{\mathrm{gal}}},$ (2)
$\displaystyle\simeq$ $\displaystyle\;2\times 10^{40}\mathrm{erg\,s^{-1}}$
$\displaystyle\left(\dfrac{\Sigma+\Sigma_{\star}}{100\,\mathrm{M_{\odot}pc^{-2}}}\right)\left(\dfrac{V_{\mathrm{rot}}}{200\,\mathrm{km\,s^{-1}}}\right)^{3}\left(\dfrac{L}{\mathrm{1\,kpc}}\right)^{2}\left(\dfrac{R_{\mathrm{gal}}}{\mathrm{2\,kpc}}\right)^{-1}$

where $\Sigma_{\star}$ is the stellar column density. The efficiency of
turbulent driving $\epsilon_{\mathrm{rot}}$ through gravitational instability
clearly depends on the disk stability, which is quantified by the Toomre
parameter (Toomre, 1964),

$Q=\dfrac{\sigma\kappa}{\pi G\Sigma},$ (3)

where $\kappa$ is the epicyclic frequency. Since a galactic disk is composed
of gas and stars, it is likely that the efficiency $\epsilon_{\mathrm{rot}}$
depends on the gas column density $\Sigma$ and the gas fraction
$\Sigma/(\Sigma+\Sigma_{\star})$ (Fensch & Bournaud, 2021). Galaxies with a
high gas fraction are more prone to gravitational instabilities, and we can
therefore expect the efficiency $\epsilon_{\mathrm{rot}}$ in these galaxies to
be high, probably close to 1. However, Milky Way-type galaxies can make
several rotations without shrinking significantly due to angular momentum
loss. For these galaxies, we expect the efficiency $\epsilon_{\mathrm{rot}}$
to be well below 1.

Accretion at a given mass rate $\dot{M}_{\mathrm{accr}}$ and a infall velocity
$v_{\mathrm{inf}}$ not only brings mass into the galaxy, but also delivers
energy at a rate $\dot{M}_{\mathrm{accr}}v_{\mathrm{inf}}^{2}/2$. To estimate
the mass accretion rate $\dot{M}_{\mathrm{accr}}$, we assume that the total
amount of gas is roughly constant, that is, the accretion of gas balances the
SFR,

$\dot{M}_{\mathrm{accr}}\simeq L^{2}\Sigma_{\mathrm{SFR}}.$ (4)

We estimate the surface density of the SFR by assuming that the SK relation
(Kennicutt & Evans, 2012) is verified,

$\Sigma_{\mathrm{SFR}}\simeq
0.1\mathrm{M_{\odot}\,yr^{-1}\,kpc^{-2}}\left(\dfrac{\Sigma}{100\mathrm{M_{\odot}\,pc^{-2}}}\right)^{1.4}.$
(5)

On the other-hand, we compute the infall velocity from the gravitational
potential,

$\dfrac{v_{\mathrm{inf}}^{2}}{2}\simeq
G\left(\Sigma+\Sigma_{\star}+\Sigma_{\mathrm{dm}}\right)R_{\mathrm{gal}},$ (6)

where $\Sigma_{\mathrm{dm}}$ is the column density of the dark matter. We thus
find that

$\displaystyle P_{\mathrm{acc}}\simeq$
$\displaystyle\;\Sigma_{\mathrm{SFR}}L^{2}G(\Sigma+\Sigma_{\star}+\Sigma_{\mathrm{dm}})R_{\mathrm{gal}}$
(7) $\displaystyle\simeq$ $\displaystyle\;5\times 10^{37}\mathrm{erg\,s^{-1}}$
$\displaystyle\left(\dfrac{\Sigma+\Sigma_{\star}+\Sigma_{\mathrm{dm}}}{100\,\mathrm{M_{\odot}pc^{-2}}}\right)\left(\dfrac{\Sigma}{100\,\mathrm{M_{\odot}pc^{-2}}}\right)^{1.4}\left(\dfrac{L}{\mathrm{1\,kpc}}\right)^{2}\left(\dfrac{R_{\mathrm{gal}}}{\mathrm{2\,kpc}}\right).$

It is complicated to determine the efficiency of turbulent driving by
accretion, $\epsilon_{\mathrm{acc}}$. Based on colliding-flow simulations,
Klessen & Hennebelle (2010) estimated it to be proportional to the density
contrast between the accreted material and the medium in which the turbulence
is driven. Forbes et al. (2023) have estimated it in the context of a galactic
disk and found that it varies with galactic radii, with a typical efficiency
of about $\epsilon_{\mathrm{acc}}\simeq$15-20$\%$. The stationary hypothesis
is rather conservative, and it is possible that the accretion rate is much
higher, in which case the resulting turbulent driving will be stronger as
well.

The amount of energy injected per unit of time by SN remnants is given by the
number of SNe produced per unit of time multiplied by the amount of energy a
single SN delivers to the ISM. This is roughly equal to
$(\Sigma_{\mathrm{SFR}}L^{2})/120$ M⊙, as given by the Salpeter initial mass
function, multiplied by $10^{51}$ erg.

By again assuming that the SFR is in line with the SK relation (Eq. (5)), we
can estimate

$\displaystyle P_{\mathrm{sn}}$
$\displaystyle\simeq\dfrac{\Sigma_{\mathrm{SFR}}L^{2}}{120\mathrm{M_{\odot}}}10^{51}\,\mathrm{erg}$
(8) $\displaystyle\simeq 2\times
10^{40}\mathrm{erg\,s^{-1}}\left(\dfrac{\Sigma}{100\,\mathrm{M_{\odot}pc^{-2}}}\right)^{1.4}\left(\dfrac{L}{\mathrm{1\,kpc}}\right)^{2}.$

The efficiency of turbulent driving by SN remnants, $\epsilon_{\mathrm{sn}}$,
has been estimated from kiloparsec turbulent box calculations (e.g. Iffrig &
Hennebelle, 2017) to be about a few percent.

To summarise, the three envisaged mechanisms could all provide a substantial
amount of energy to drive turbulence and compensate for turbulent dissipation.
They present different scaling with $\Sigma$, however, more precisely,
$P_{\mathrm{rot}}\propto\Sigma$, $P_{\mathrm{acc}}\propto\Sigma^{2.4}$ and
$P_{\mathrm{sn}}\propto\Sigma^{1.4}$. We stress that given the orders of
magnitude and the dependence on $\Sigma$, it is likely the case that for high
values of $\Sigma$, energy rotation of both gas and stars is likely the
dominant mechanism that triggers turbulence. This is also the source of energy
envisioned by Krumholz et al. (2018) and Nusser & Silk (2022).

##  3 Numerical setup

We used the RAMSES code (Teyssier, 2002) to solve the equations of
Magnetohydrodynamics with a Godunov solver (Fromang et al., 2006) on a uniform
cubic grid of $256^{3}$ cells with periodic boundaries in the $x$ and $y$
directions parallel to the disk and open boundaries in the vertical direction.
The box represents a cubic region of the galactic disk of size $L=1$ kpc, so
that the resolution is about 4 pc. The impact of the resolution on the results
is discussed in Appendix C.

###  3.1 Initial conditions

We used initial conditions similar to those described in previous works
(Colling et al., 2018; Brucy et al., 2020). The interstellar gas was initially
distributed as a Gaussian along the $z$-axis,

$n(z)=n_{0}\exp\left(-\frac{1}{2}\left(\frac{z}{z_{0}}\right)^{2}\right),$ (9)

with $n_{0}$ a free density parameter and $z_{0}=150\ \mathrm{pc}$ the typical
scale height. The column density of gas (hydrogen and helium), integrated
along the z-axis (perpendicular to the disk) is then

$\Sigma_{\text{gas},0}=\sqrt{2\pi}m_{p}n_{0}z_{0},$ (10)

where $m_{p}=1.4\times 1.66\cdot 10^{-24}$ g is the mean mass of a gas
particle. We set the initial temperature at $8000\ \mathrm{K}$ to match the
typical value of the temperature of the warm neutral medium phase of the ISM.
An initial synthetic turbulent velocity field with a 3D dispersion of
$5\sqrt{n_{0}/1.5}\leavevmode\nobreak\ \mathrm{km\cdot s}^{-1}$ was also
added. It presents a Kolmogorov power spectrum and was generated using a
random phase (Kolmogorov, 1941). The dependence on $n_{0}$ ensures that the
initial ratio of the kinetic and gravitational energy due to the gas stays the
same. Finally, we added a Gaussian magnetic field, oriented along the
$x$-axis,

$B_{x}(z)=B_{0}\exp\left(-\frac{1}{2}\left(\frac{z}{z_{0}}\right)^{2}\right),$
(11)

where $B_{0}$ is a parameter of our simulation with a value of a few
microgauss (see Table 1).

###  3.2 Numerical models

The gas is subject to an external gravitational potential, corresponding to
old stars and dark matter, of the form

$g_{\mathrm{ext}}(z)=-\frac{a_{1}z}{\sqrt{z^{2}+z_{0}^{2}}}-a_{2}z,$ (12)

with $a_{1}=1.42\times 10^{-3}$ kpc/Myr2, $a_{2}=5.49\times
10^{-4}\,\mathrm{Myr}^{-2}$ , and $z_{0}=0.18$ kpc (Kuijken & Gilmore, 1989;
Joung & Mac Low, 2006). This gravitational force adds up to the self-gravity
of the gas.

Sink particles (Bleuler & Teyssier, 2014) were used to follow the dense gas
and model star formation. Sink creation was triggered when the gas number
density passed a threshold of $10^{3}\ \mathrm{cm^{-3}}$. The simulation with
an initial column density higher than $50\leavevmode\nobreak\
\text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$ used a higher threshold of $5\times
10^{3}\ \mathrm{cm^{-3}}$. Since the SFR in these simulations is high, this
allowed us to reduce the number of sinks without changing the SFR much
(Colling et al., 2018). Gas with a density higher than the threshold within a
radius of four cells (16 pc) was then progressively accreted into the sink,
with no more than 10 % of the gas being accreted in one time step.

The gas cooling was as described in Audit & Hennebelle (2005), which takes
into account all the standard processes taking place in the atomic gas. It is
essentially identical to the cooling described in Koyama & Inutsuka (2000),
for instance.

###  3.3 Stellar feedback

The simulations included models for the formation and expansion of HII
regions, SN explosions, and the far-ultraviolet (FUV) feedback. They are the
same as in Paper I. The HII and SN feedbacks were attached to the sinks and
were described in length in Colling et al. (2018) and references therein. Each
time a sink had accreted a mass of 120$\leavevmode\nobreak\
\text{M}_{\astrosun}$, a massive star particle with a mass randomly determined
from the Salpeter initial mass function (Salpeter, 1955) between 8 and 120
$\leavevmode\nobreak\ \text{M}_{\astrosun}$ was created and attached to the
sink. The lifetime $\tau_{*}$ of this star was computed using the model

$\tau_{*}\left(M\right)=\tau_{0}\exp\left(-a\left(\log\left(\frac{M}{M_{0}}\right)\right)^{b}\right),$
(13)

with $\tau_{0}=3.265$ Myr, $M_{0}=148.16\leavevmode\nobreak\
\text{M}_{\astrosun}$, $a=0.238$, and $b=2.205$ (Woosley et al., 2002). When
this massive star reached the end of its lifetime $\tau_{*}$, it exploded in a
random location within a sphere with a radius $\tau_{*}\times$ 1 km $\cdot$
s-1. The gas that was located inside a sphere of 12 pc radius around the
location of the SN was heated up to inject a thermal energy of
$10^{51}\mathrm{erg}$. Since the cooling radius was almost never resolved, the
gas would cool down immediately and the SN would have no effect. To avoid
this, the explosion also injected $4\times 10^{43}$ g $\cdot$ cm $\cdot$ s -1
momentum into the same region. When SNe explode in a low-density environment,
a very high temperature and very high velocities can be generated, which
translate into very short time steps due to the Courant condition. A
limitation of the temperature to $10^{6}\leavevmode\nobreak\ \mathrm{K}$ and
of the velocity to $300\leavevmode\nobreak\ \mathrm{km}\cdot\mathrm{s}^{-1}$
generated by a SN explosion was therefore implemented. We discuss the impact
of this limitation on our results in Appendix B.

We also included self-consistent feedback from HII regions. Energy and
momentum were injected according to the flux of ionising photons emitted by
the massive star (Vacca et al., 1996). The evolution of the HII regions
themselves was computed via the so-called M1 radiative transfer method
(Rosdahl et al., 2013).

The FUV heating was uniform. However, it was not kept constant at the solar
neighbourhood value because young O-B stars contribute significantly to the
FUV emission. As a first approximation, we considered the UV heating to be
proportional to the SFR (Ostriker et al., 2010). The mean FUV density relative
to the solar neighbourhood value $G_{0}^{\prime}$ can then be written as

$G_{0}^{\prime}=\frac{\Sigma_{\mathrm{SFR}}}{\Sigma_{\mathrm{SFR,}\astrosun}}=\frac{\Sigma_{\mathrm{SFR}}}{2.5\times
10^{-9}\leavevmode\nobreak\
\mathrm{M}_{\astrosun}\cdot\mathrm{pc}^{-2}\cdot\mathrm{yr}^{-1}}.$ (14)

In our model, $G_{0}^{\prime}$ has a minimum value of $1$ (as a background
contribution) and follows equation 14 when the SFR increases.

###  3.4 Injection of turbulence

Bournaud et al. (2010), Krumholz & Burkhart (2016), and Krumholz et al. (2018)
showed that for galaxies with a high column density or high SFR, large-scale
gravitational instabilities are significant sources of turbulent energy and
might even dominate stellar feedback. The ISM in a kiloparsec region of the
galaxy is far from an isolated environment. It interacts with the rest of the
galaxy. Adding turbulent driving to simulations of restricted regions of a
galaxy enables us to take a part of these complex interactions into account
while keeping the computational cost of the simulation reasonable. The
turbulence added in the kiloparsec box is clearly not fully self-consistent.
We numerically investigated the effect of turbulent driving on star formation.
We used a model for turbulent driving that we adapted from the generalisation
of Ornstein-Uhlenbeck that was developed and used by several authors (Eswaran
& Pope, 1988; Schmidt et al., 2006, 2009; Federrath et al., 2010).

The turbulent forcing was implemented by adding an external force density
$\bm{f}$ that accelerated the fluid on large scales. Our main hypothesis is
that this force is generated by large-scale processes at the galactic scale
within the galactic disk. As a consequence, the driving is bidimensional (the
force has no vertical component and no vertical mode), and it was applied only
at low altitude. The driving procedure is slightly different from the
procedure used in Paper I and is therefore described below at length.

The fluid equations that we solved were

$\displaystyle\dfrac{\partial\rho}{\partial
t}+\bm{\nabla}\cdot\left(\rho\bm{v}\right)$ $\displaystyle=0,$ (15)
$\displaystyle\rho\left(\dfrac{\partial\bm{v}}{\partial
t}+\left(\bm{v}\cdot\bm{\nabla}\right)\bm{v}\right)$
$\displaystyle=-\bm{\nabla}P+\dfrac{\left(\bm{\nabla}\times\bm{B}\right)\times\bm{B}}{4\pi}+\rho\bm{g}+\rho\bm{f},$
(16) $\displaystyle\rho\left(\dfrac{\partial e}{\partial
t}+\left(\bm{v}\cdot\bm{\nabla}\right)e\right)$
$\displaystyle=-P\left(\nabla\cdot\bm{v}\right)+\rho\bm{f}\cdot\bm{v}-\rho\mathcal{L},$
(17) $\displaystyle\dfrac{\partial\bm{B}}{\partial
t}-\bm{\nabla}\times\left(\bm{v}\times\bm{B}\right)$ $\displaystyle=0.$ (18)

In these equations, $\rho$ is the gas density, $\bm{v}$ is the gas speed, $P$
is the thermal pressure, $e$ is the internal energy, $\bm{B}$ is the magnetic
field, $\bm{g}$ is the gravitational potential, and $\mathcal{L}$ is the
energy loss function.

Table 1: List of simulations. Each group of simulations corresponds to an experiment discussed in the article. The simulations of the group noturb were discussed in Paper I. The parameter $n_{0}$ is the initial midplane density and sets the initial column density $\Sigma_{0,\mathrm{gas}}$ (see Eqs. (9) and (10)). The strength of the driving is given by parameter $f_{\mathrm{rms}}$ (see Eq. (23)). A value of $f_{\mathrm{rms}}=10^{4}$ in code units corresponds to an RMS acceleration due to the driving of $1.46\ \mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{Myr}^{-1}$. The parameter $\chi=1-\zeta$, introduced in Eqs. (21) and (22), is the compressibility fraction of the turbulent driving. Finally, $B_{0}$ is the initial value of the magnetic field, which was initialised as described in Eq. (11).  Group | $n_{0}$ | $\Sigma_{0,\mathrm{gas}}$ | $f_{\text{rms}}$ | $\chi$ | $B_{0}$
---|---|---|---|---|---
| $[\mathrm{cm}^{-3}]$ | $[\leavevmode\nobreak\ \text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}]$ | $[1.46\times 10^{-4}\leavevmode\nobreak\ \mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{Myr}^{-1}]$ |  | $[\mu\mathrm{G}]$
strength | 3 | 38.74 |  $10^{4}$ to $10^{5}$ | 0.25 | 3.8
6 | 77.4 |  $5\times 10^{4}$ to $8\times 10^{5}$ | 0.25 | 3.8
comp | 3 | 38.7 | $8\times 10^{4}$ | 0 to 1 | 7.6
6 | 77.4 | $2\times 10^{5}$ | 0 to 1 | 7.6
mag | 3 | 38.7 | $8\times 10^{4}$ | 0.25 | [0, 1.9, 3.8, 5.3, 7.6]
6 | 77.4 | $2\times 10^{5}$ | 0.25 | [0, 1.9, 3.8, 7.6, 15, 23]
noturb | 1 | 12.9 | 0 | N/A | 3.8
1.5 | 19.4 | 0 | N/A | 3.8
2 | 25.8 | 0 | N/A | 3.8
3 | 38.7 | 0 | N/A | 3.8
4 | 51.6 | 0 | N/A | 3.8
6 | 77.4 | 0 | N/A | 3.8
12 | 155 | 0 | N/A | 3.8
noturb_hB | 1 | 12.9 | 0 | N/A | 7.6
3 | 38.7 | 0 | N/A | 7.6
6 | 77.4 | 0 | N/A | 7.6
12 | 155 | 0 | N/A | 7.6
turb | 1.5 | 19.4 | $5\times 10^{3}$ | 0.25 | 3.8
3 | 38.7 | $6\times 10^{4}$ | 0.25 | 3.8
6 | 77.4 | $4\times 10^{5}$ | 0.25 | 3.8
12 | 155 | $3\times 10^{6}$ | 0.25 | 3.8
turb_hB | 1.5 | 19.4 | $1\times 10^{3}$ | 0.25 | 7.6
3 | 38.7 | $3\times 10^{4}$ | 0.25 | 7.6
6 | 77.4 | $10^{5}$ | 0.25 | 7.6
12 | 155 | $10^{6}$ | 0.25 | 7.6
turb_Bvar | 3 | 38.7 | $5\times 10^{4}$ | 0.25 | 5.3
6 | 77.4 | $2\times 10^{5}$ | 0.25 | 7.6
12 | 155 | $10^{6}$ | 0.25 | 10.7
turb3D | 3 | 38.74 | $6\times 10^{4}$ | 0.25 | 7.6
6 | 77.4 | $2\times 10^{5}$ | 0.25 | 7.6
| 12 | 155 | $8\times 10^{5}$ | 0.25 | 7.6

The evolution of the Fourier modes of the turbulence driving the acceleration
field $\bm{\hat{f}}(\bm{k},t)$ follows a stochastic differential equation
(Schmidt et al., 2006, 2009),

$\mathrm{d}\bm{\hat{f}}(\bm{k},t)=-\bm{\hat{f}}(\bm{k},t)\frac{\mathrm{d}t}{T}+F_{0}(\bm{k})\bm{P_{\zeta}}\left(\bm{k}\right)\mathrm{d}\bm{W}_{t}.$
(19)

In this equation, $\mathrm{d}t$ is the time step for integration, and $T$ is
the autocorrelation timescale. In our simulations, we chose $T=5$ Myr and
$\mathrm{d}t/T=1/100$. The value of $5$ Myr was chosen to be of the same order
of magnitude as the crossing time of the over-densities created by the
turbulent driving. $d\bm{W_{t}}$ is a small vector that was randomly chosen
following the Wiener process, as described in Schmidt et al. (2009). The main
modification to the Ornstein-Uhlenbeck process is the weighting function of
the driving modes $F_{0}$ and the projection operator $\bm{P_{\zeta}}$. We
selected only large 2D modes, and more weight was given to modes with a
wavelength of 500 pc (half the box size),

$F_{0}(\bm{k})=\begin{cases}1-\left(\dfrac{\bm{k}}{2\pi}-2\right)^{2}\text{ if
}1<\dfrac{|k|}{2\pi}<3\text{ and }k_{z}=0\\\ 0\text{ if not.}\end{cases}.$
(20)

In other words, the driving modes (divided by $2\pi$) were

$\left(\begin{array}[]{c}2\\\ 2\\\
0\end{array}\right),\left(\begin{array}[]{c}1\\\ 1\\\
0\end{array}\right),\left(\begin{array}[]{c}2\\\ 1\\\
0\end{array}\right),\left(\begin{array}[]{c}1\\\ 2\\\
0\end{array}\right),\left(\begin{array}[]{c}2\\\ 0\\\
0\end{array}\right),\left(\begin{array}[]{c}0\\\ 2\\\ 0\end{array}\right),$

and the respective weights were approximately

$0.31,\;0.65,\;0.65,\;0.94,\;0.94,\;1,\;1.$

Figure 1:  Column densities for the group strength (top) and comp (bottom) for
a medium column density. For each panel, the top subpanel is viewed edge-on
along the x-axis and the bottom subpanel is viewed face-on. White dots
represent sink particles, with a radius proportional to their mass. Top: RMS
acceleration of the turbulent driving is increased from left to right. Bottom:
Compressibilty of the turbulent driving is increased from left to right.
$\chi=0$ and $\chi=1$ refer to purely solenoidal and compressive driving,
respectively.

The projection operator $\bm{P_{\zeta}}$ is the weighted sum of the components
of the Helmholtz decomposition of compressive versus solenoidal modes,
projected in 2D,

$\bm{P_{\zeta}}(\bm{k})=\begin{pmatrix}\zeta+(1-2\zeta)\left(\dfrac{k_{x}^{2}}{k^{2}}\right)&(1-2\zeta)\left(\dfrac{k_{x}k_{y}}{k^{2}}\right)&0\\\
(1-2\zeta)\left(\dfrac{k_{x}k_{y}}{k^{2}}\right)&\zeta+(1-2\zeta)\left(\dfrac{k_{y}^{2}}{k^{2}}\right)&0\\\
0&0&0\\\ \end{pmatrix}.$ (21)

In other words,

$\bm{P_{\zeta}}(\bm{k})=\zeta\bm{P}^{\perp}(\bm{k})+(1-\zeta)\bm{P}^{\parallel}(\bm{k}),$
(22)

with $\bm{P}^{\perp}$ and $\bm{P}^{\parallel}$ the projection operators
perpendicular and parallel to $\bm{k},$ respectively (Federrath et al., 2010).
The weight $\zeta$ is the solenoidal fraction of the driving. In the
following, we also refer to the compressive fraction, $\chi=1-\zeta$.

Figure 2: Effect of the intensity of the turbulent driving on the SFR (group
strength). $f_{\mathrm{rms}}$ is a measure of the driving strength. More
precisely, $f_{\mathrm{rms}}$ is the RMS acceleration of the force density in
code units. Top: Mass accreted in sinks. The SFR computation started when 3%
were consumed or lost ($t_{\mathrm{3\%}}$, left triangle) and ended when 40%
of gas was consumed or lost ($t_{\mathrm{40\%}}$, right triangle). Bottom:
Surface density of the SFR as a function of the strength of the driving
$f_{\mathrm{rms}}$. The plotted value is the accretion rate between
$t_{\mathrm{3\%}}$ and $t_{\mathrm{40\%}}$ , and the error bars reflect the
standard deviation of all the accretion rates we obtained by choosing starting
and end points between these two values. The dashed line is the value of the
SFR given by the SK relation (Kennicutt & Evans, 2012). Figure 3: Relation
between the strength of the driving, the velocity dispersion, and the SFR. Top
panel: Mass-weighted velocity dispersion measured in the strength simulations
at 58 Myr at the kiloparsec scale. The velocity dispersion parallel to the
disk plane $\sigma_{\mathrm{disk}}$ is defined by
$\sigma_{\mathrm{disk}}=\sqrt{\sigma_{x}^{2}+\sigma_{y}^{2}}/\sqrt{2}$. Bottom
panel: SFR as a function of the velocity dispersion parallel to the disk. The
dashed line is the value of the SFR given by the SK relation.

The forcing field $\bm{f}(\bm{x},t)$ was then computed from the Fourier
transform,

$\bm{f}(\bm{x},t)=\mathrm{att}(z)g(\zeta)f_{\mathrm{rms}}\int\bm{\hat{f}}(\bm{k},t)e^{i\bm{k}\cdot
x}\mathrm{d}t^{3}\bm{k}.$ (23)

The parameter $f_{\mathrm{rms}}$ controls the total power that is injected by
the turbulent force into the simulation. The $g(\zeta)$ factor is an empirical
correction so that the resulting time-averaged root mean square (RMS) of the
power of the Fourier modes is equal to $f_{\mathrm{rms}}$, regardless of the
solenoidal fraction $\zeta$ (see the discussion in Appendix A). The
attenuation function $\mathrm{att}$ ensures that the driving occurs only on
the disk,

$\mathrm{att}(z)=\exp\left(\dfrac{-z^{2}}{2z_{t}^{2}}\right),$ (24)

where $z_{t}=75$ pc.

There are two main differences with the driving method used in Paper I. First,
the modes aligned with the $x$ \- and $y$ -axis were not used in Paper I,
leading to global motions than can be seen through the diagonal features in
the lower right panel of Fig. 1 of that paper. Second, the attenuation
function (24) was not used in Paper I. It is introduced here so that the
driving occurs only within the scale height of the disk.

###  3.5 Computation of the SFR and the column density

The total mass of gas accreted by the sinks was considered to form stars and
is denoted $\mathrm{M}_{\star}(t)$. As described in section 3.2, only a
fraction of the dense gas is accreted onto the stars so that the efficiency of
the conversion of dense gas into stars is strictly lower than 1. The gas in
the box is continually consumed by sinks or expelled from the box through the
vertical boundaries. As a consequence, the gas available for star formation
decreases with time. Since we used the simulations to test the star formation
relation that links the SFR with the column density of gas, and to avoid too
huge variations of the gas reservoir, we restricted the analysis to a short
enough time span. In detail, we started the analysis at $t_{\mathrm{3\%}}$
when 3% of the gas was consumed or lost, and we stopped it when this depletion
fraction reached 40% ($t_{\mathrm{40\%}}$). The surface density of SFR was
then computed via

$\Sigma_{\mathrm{SFR}}(t_{3\%},t_{40\%})=\dfrac{\mathrm{M}_{\star}(t_{40\%})-\mathrm{M}_{\star}(t_{3\%})}{L^{2}\left(t_{40\%}-t_{3\%}\right)}.$
(25)

The uncertainties on this value were computed by taking the standard deviation
of all values of the SFR between $t_{3\%}$ and $t_{40\%}$, that is,
$\Sigma_{\mathrm{SFR}}(t_{i},t_{j})$ for $t_{3\%}\leq t_{i}<t_{j}\leq
t_{40\%}$.

The mean column density in the considered time span varies between
$0.97\Sigma_{0}$ and $0.6\Sigma_{0}$, where $\Sigma_{0}$ is the initial column
density. For each simulation, the associated column density is then
$\Sigma=0.8\Sigma_{0}\pm 0.2\Sigma_{0}$. The determination of
$\Sigma_{\mathrm{SFR}}$ and of $\Sigma$ is different from the determinations
used in Paper I, where the SFR was averaged over a fixed period of time and
only the initial column density was considered. This means that the depletion
of the gas during the time span covered by the simulation was not taken into
account.

###  3.6 List of simulations

We ran several simulations that we list in Table 1. We separated them into
groups, each of which corresponds to an experiment. There are two families of
experiments.

The first family of experiments (strength, comp, and mag) aimed to quantify
the effect of one given parameter on the SFR, and the properties of the ISM in
general. For these experiments, we set all parameters except for the tested
parameter to fiducial values. The tested parameters were the strength of the
turbulent driving, its compressibility, and the initial magnetic field.

The goal of the second family (noturb, noturb_hB, turb, turb_hB, turb_Bvar,
and turb3D) was to test whether it is possible to derive an SFR that is
compatible with the SK relation in simulation boxes of 1 kpc3 even in a dense
environment under various conditions. Each of the groups contained three to
seven simulations with increasing column density. We computed the SFR for each
of the simulations and compared it to what is expected from the SK relation.
The main parameters we varied are the strength of the turbulent driving and
the intensity of the initial magnetic field. For the latter, several values
were considered for each column density. For the simulation with turbulent
driving, we chose the driving strength by trial and error to find an
acceptable SFR. The goals here were to show that it is possible to obtain an
SFR that was similar to the observations when turbulent driving was added, but
also to analyse the resulting velocity dispersion and gas distribution.

##  4 Effects of turbulence

In the following experiments, we consider dense regions of galaxies with an
initial column density of $38$ and $77\leavevmode\nobreak\
\text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$. They are referred to as medium and
high column density in what follows. We kept all parameters constant except
for one parameter.

We started by investigating the properties of turbulence, in particular, the
strength of the driving and the fraction of compressible modes. Whereas we
expect that higher compressibility leads to a higher SFR, the impact of higher
Mach numbers is less straightforward. One of the difficulties that arises when
trying to understand the effect of the turbulent driving on the SFR is that
turbulence can both trigger and quench star formation (Mac Low & Klessen,
2004). In classical models for star formation (Padoan & Nordlund, 2011;
Hennebelle & Chabrier, 2011; Federrath & Klessen, 2012), the SFR is
principally set by the width of the density probabilty distribution function
(PDF), the density threshold above which the gas collapses under its own
gravity, and the time needed for the collapse. When the Mach number is
increased, the width of the PDF as well as the density threshold also
increase.

###  4.1 Influence of the strength of the turbulent driving

Figure 4: Effect of the compressibility of the turbulent driving on the SFR
(group comp). Surface density of the SFR as a function of the compressibility
$\chi$ of the driving. The dashed line is the value of the SFR given by the SK
relation. Figure 5: Mass-weighted velocity dispersion measured in the comp
simulations at 58 Myr at the kiloparsec scale.  Figure 6: Effect of the
intensity of the initial magnetic field on the SFR (group mag). $B_{0}$ is the
initial midplane intensity of the magnetic field. The legend is the same as in
Fig. 2. The dashed line is the value of the SFR given by the SK relation.

Figure 7: Effect of the magnetic field on the structure of the disk and the
outflows (group mag)´ for a high column density. _Top left panel:_ Averaged
density profile (slightly smoothed for readability) at $t\approx
58\leavevmode\nobreak\ \mathrm{Myr}$. For runs in which the midplane density
was close to the maximum density, the position at which the density is half of
the midplane density is marked with a plus. The dashed black line corresponds
to the initial profile. The colours of the lines refer to the same colour bar
as in the plot of the top right panel. _Top right panel_ : Mass loss via
outflows through the open up and down boundaries. _Lower panel:_ Edge-on (top)
and face-on (bottom) column densities.  Figure 8: Velocity dispersion measured
in the MAG simulations at 58 Myr.

We used our numerical setup to probe the effect of the strength of turbulence
(simulation set strength in Table 1) by varying the RMS acceleration of the
turbulence between $1\times 10^{4}$ to $1\times 10^{5}$ in code units111the
code units for the RMS acceleration of the turbulence convert into 1.46
$\times 10^{-4}\leavevmode\nobreak\
\mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{Myr}^{-1}$. for medium column
densities and from $5\times 10^{4}$ to $8\times 10^{5}$ for a high column
density.

These values were selected because they lead to SFR values that are not too
far to the value inferred from the SK relation. Figure 1 portrays column
density images for the medium column density simulation with several values of
the turbulent forcing strength, projected along the x-axis (top panel) and
z-axis (bottom panel). When the turbulent forcing is weak, the gas and the
recently formed stars are organised in a few dense regions. That is to say,
the gas is organised as few giant molecular complexes. As turbulent driving is
increased, the density field is less strongly concentrated and the stars tend
to be distributed more uniformly. We also note that their number is clearly
reduced. Figure 2 shows that with stronger turbulent driving, accretion onto
the sinks is indeed slower, and the SFR decreases.

Fig. 1 shows that for the highest forcing values, the anisotropy that results
from the stratification is even more pronounced. This is because not only the
forcing alone applies in the xy-plane, but also because the star formation and
therefore the stellar feedback are reduced.

Increasing the strength of the driving naturally leads to an increase in the
velocity dispersion (Fig. 3). Because the driving is 2D, projected on the
plane of the disk, only the velocity dispersion along the disk plane, defined
as $\sigma_{disk}=\sqrt{\sigma_{x}^{2}+\sigma_{y}^{2}}/\sqrt{2}$, is
increased. We would expect the vertical speed dispersion, which is assumed to
be mainly due to the stellar feedback, to be reduced with stronger driving,
but no such trend is clearly observed. The velocity dispersion scales linearly
with the driving strength, but saturates for very strong driving.
Interestingly, the velocity field is generally significantly non-isotropic,
with a ratio of the velocity dispersion parallel to the disk $\sigma_{disk}$
and the vertical velocity dispersion $\sigma_{z}$ that can reach a factor of
5.

The outcome of this increased velocity dispersion is a reduction in the SFR
(Figs. 2 and 3). An increase in the velocity dispersion by a factor of a few
leads to a decrease in the SFR by nearly a factor 10, which is a very
significant drop. This illustrates the fact that the externally driven
turbulence can substantially modify the ISM evolution in this regime.

###  4.2 Influence of the compressive fraction on the turbulence

The dual effect of turbulence, which both triggers and quenches star
formation, is particularly clear when the compressive fraction of the driving
is varied. The velocity field can be divided through the Helmotz decomposition
into compressible modes $\bm{v_{c}}$ and solenoidal modes $\bm{v_{s}}$.
Compressible modes verify $\nabla\times\bm{v_{c}}=0$ and trigger compression
and dilation of the density field. By contrast, solenoidal modes verify
$\nabla\cdot\bm{v_{s}}=0$ and thus do not affect the density field. When the
compressive fraction $\chi$ of the turbulent driving is changed (see Eq.
(22)), the amount of energy in the different velocities modes is changed
accordingly. However, it is important to note that the correspondence between
the compressive fraction of the turbulent driving $\chi$ and the resulting
ratio of the compressive and solenoidal modes in the velocity field is not
perfect. Rotation motions induced by solenoidal driving eventually generate
compression as well (Vazquez-Semadeni et al., 1996), especially as the
turbulence cascades to lower scales (Federrath, 2013). With solenoidal
driving, a smaller fraction of the kinetic energy is globally in compressive
modes that can produce density fluctuations. This is well illustrated in the
study of the PDF in purely hydrodynamic simulations (Federrath et al., 2008).
The density PDF is broader and explores higher densities when the compressive
fraction $\chi$ is closer to one (compressive driving) than when it is closer
to zero (solenoidal driving). To quantify the effect of the compressibility of
the turbulent driving on the SFR, we ran the comp experiment. We again
selected galaxies with an initial column density of $38$ and
$77\leavevmode\nobreak\ \text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$. With every
other parameter fixed, we modified the compressive fraction $\chi$ from $0$ to
$1$.

The column density maps of snapshots at 58 Myr of these simulations are shown
in the lower panel of Fig. 1. From left to right, the compressibility is
increased from 0 to 1. In the edge-on view (top), the purely solenoidal case
produces horizontal filaments, while more clumpy structures are produced when
the compressibility increases. The face-on view shows that dense structures
tend to be smaller when the compressive fraction is smaller. These larger and
denser structures are actively forming stars. As a result, and this is
depicted in Fig. 4, increasing the compressibility increases the SFR by a
factor of about 10. However, even a purely compressive run produces fewer
stars than a run without turbulent driving (compared with Fig. 9). Therefore,
the effect of turbulent driving at this strength is always to quench star
formation, regardless of the compressibility.

Interestingly, even though the driving was calibrated to produce an identical
velocity dispersion regardless of the compressibility (see Appendix A), a
dependence of the velocity dispersion on the compressibility (see Fig. 5) is
nevertheless observed. The variations in the velocity dispersion remain well
below the variations induced by the variations of the driving strength,
however.

Overall, these results clearly indicate that the nature and thus the sources
of the turbulence are important to identify, not only regarding the driving
intensity of the turbulence, but also regarding the compressibility of the ISM
turbulence. Both play a significant role in setting the SFR.

Figure 9: Stellar mass for the group noturb. The red line is the stellar mass
creation when the SFR in the kiloparsec box is equal to $2.5\times
10^{-3}(\Sigma/10\leavevmode\nobreak\
\text{M}_{\astrosun}\cdot\mathrm{pc}^{-2})^{1.4}\leavevmode\nobreak\
\text{M}_{\astrosun}\cdot\mathrm{yr}^{-1}$ (Kennicutt & Evans, 2012). The red
zone corresponds to an SFR that is three times higher and lower and aims to
represent the spread of the SK relation. The orange histogram is the mass
(multiplied by 100 for readability) of stars that blow up into SNe in a time
span of 20 Myr. The SFR computation for Fig. 10 starts when 3% of the gas is
consumed or lost ($t_{\mathrm{3\%}}$, left triangle) and ends when 40% of the
gas is consumed or lost ($t_{\mathrm{40\%}}$, right triangle).

##  5 Influence of the magnetic field

The magnetic field plays an important role on the gas dynamics and in the
structure of the ISM (Hennebelle & Inutsuka, 2019). In the following, we study
the impact of stronger initial magnetic fields on the SFR in our simulations
when all other parameters remain the same. We considered the simulations from
the group mag (cf Table 1). We varied the initial midplane magnetic field
$B_{0}$ from 0 $\mu\mathrm{G}$ to 15 $\mu\mathrm{G}$. As a reminder, the
initial magnetic field is uniform in the galactic plane and follows a Gaussian
distribution in the vertical direction (see Eq. (11)). The initial column
density for these runs was also $33.7$ (medium) or $77.4\leavevmode\nobreak\
\text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$ (high).

Figure 6 shows that an increased magnetic field has a quasi-linear effect on
the SFR. For the two cases (medium and high column density), as $B_{0}$
increased from 0 to 7 $\mu$G, the SFR has dropped by a factor of nearly 10.
Since high magnetic intensities lead to a substantial reduction in SFR, this
implies that a strong magnetic field could reduce the turbulent energy that
would otherwise be required to sustain a low SFR. This point is further
investigated in section 6.3.

Figure 10: Star formation relation for the groups noturb and turb. The SFR and
the associated uncertainties are computed as explained in section 3.5. We
added data from observations of normal and irregular galaxies, and the blue
line represents the SK relation (Kennicutt & Evans, 2012). We also compare our
results with the results from the R models of Kim et al. (2020).

The effect of the magnetic field on the SFR is a consequence of magnetic
support, which locally resists the impact of self-gravity, but also of the
global magnetic support that leads to a large disk scale height (Iffrig &
Hennebelle, 2017) and therefore to reduced densities. To further verify this
point, we investigate the effect of the magnetic intensity on the density
profile in Fig. 7. For low magnetic fields, the SFR is so high that the gas in
the midplane of the disk is accreted in stars or is swept up by the SNe, and
this depletes the gas near the equatorial plane, as we show in the column
density map and in the upper left panel, which features the averaged density
profile.

The scale height, denoted with plus marks, increases when we increase the
initial magnetic from $3$ to $15\mu$G. However, the density profile is not
smooth, and the effect of the magnetic field is better understood by comparing
the profile at $t\approx 58\leavevmode\nobreak\ \mathrm{Myr}$ (in color) with
the initial profile (dashed black line). The gas is lifted, and a significant
amount of gas is lost through outflows (upper right panel). A stronger
magnetic field triggers stronger outflows. This effect is also seen in
simulations reported by Girichidis et al. (2018). Magnetically driven outflows
can be distinguished from the SN outflows that are observed in the low
magnetic field cases because they occur earlier in the simulation. It is
important to note that increased outflows also contribute to the diminution of
the SFR. Less gas remains available for the formation of gas.

The bottom panel of Fig. 7 shows that the value of the initial magnetic field
does not seem to affect the distribution of the gas in the horizontal
direction, where turbulent motion generated by the driving dominates the
dynamics of the gas. For the same reason, it has very little effect on the
velocity dispersion, as we show in Fig. 8. The regulation effect of the
magnetic field is thus mainly due to how it changes the vertical structure.

Because the driving is not computed self-consistently, hypothetical
stabilising effects of the magnetic field analogous to what was reported by
Heitsch et al. (2007) and Zamora-Avilés et al. (2018) at smaller scales would
not be captured by our model. However, the turbulent driving aims to model
turbulence that is generated by instabilities at the galactic scale. Even if
the magnetic field intensity is locally high, the averaged magnetic field at
the galactic scale is probably not strong enough to have a significant
influence on the development of the instability.

##  6 Reproducing the Schmidt-Kennicutt relation

###  6.1 Self-regulation by stellar feedback?

Figure 11: Star formation relation for the simulations with a high magnetic
field.

In Paper I we found that stellar feedback alone cannot regulate the SFR to
match the observed values, especially when the column density is high. In this
section, we focus on the simulation of the noturb group, without injection of
turbulence from the large scale.

Self-regulation works well for galaxies with column densities similar to the
mean value of the Milky Way ($\approx 10\leavevmode\nobreak\
\text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$). In our simulation with a starting
column density of approximately $12\leavevmode\nobreak\
\text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$ ($n_{0}=1$), star formation is
moderate and well in line with the observations. However, for higher column
densities, Fig. 9 reveals that the SFR, that is, the slope of the blue solid
line, is far above the observationally determined values represented by the
slope of the thick red line. Moreover, the star formation versus column
density relation that can be derived is much steeper (see Fig. 10), with
typically ${\dot{\Sigma}}\propto\Sigma^{2.5}$. Interestingly, a similar result
was found by Kim et al. (2020). However, this comparision should be taken
carefully because the simulations from Kim et al. (2020) have a different box
size and shape, a different integration time, a different feedback
implementation, and a higher external potential for high column density.

Figure 9 shows that high SFR generates a peak in SN mass. After this peak,
star formation is reduced for a while. However, the accumulation of gas at the
large scale is not destroyed, and the very dense gas reservoir is replenished
over a period of a few tens of million years. As a consequence, star formation
can continue at a fast rate. This adds up to the fact that was already
highlighted in previous studies (Brucy et al., 2020; Nusser & Silk, 2022),
that the energy injection via stellar feedback is lower than what would be
needed to maintain the SFR at a rate compatible with the SK relation, as well
as the energy injected via gravo-turbulence to maintain the Toomre parameter
$Q\sim 1$. Therefore, this constitutes strong evidence that stellar feedback
cannot regulate star formation sufficiently in high column density
environments.

Figure 12: Markers represent the time-averaged power injected by the turbulent
driving force in order to reproduce the SK relation. The power is computed in
the whole simulation box and averaged between 4 and 60 Myr. In these
simulations, the driving strength was calibrated so that the SFR was close to
the SK relation (see section 6.2 and 6.3, as well as Figs. 10 and 11). The
mean value is the total energy injected during the period divided by the time
span, while the error bar is the standard deviation from the power computed at
each time step. The solid lines refer to the estimates of the possible source
of turbulent driving, and the notation is the same as in section 2. The
efficiency coefficients ($\epsilon_{\mathrm{sn}}\sim
0.05,\epsilon_{\mathrm{acc}}\sim 0.2,\epsilon_{\mathrm{rot}}\in[0.1,1]$) are
reasonable estimates, as discussed in section 2 as well. To draw the blue
line, we assumed that the velocity dispersion was proportional to the column
density.

###  6.2 Impact of turbulent driving on the Schmidt-Kennicutt relation

Paper I has shown that it was possible to obtain a star formation relation
that was compatible with the SK relation when the turbulent driving was scaled
with the initial column density. In this work, we reproduced this result with
a slightly different method for turbulent injection that we presented in
section 3.4 and with the improved SFR computation presented in section 3.5.
For the group of simulations turb, we kept the same magnetic field for all the
simulations, and we studied the increasing column densities. The compressive
fraction was set at 0.25. This choice can be discussed: it is lower than the
natural mix, which would be 0.5 (Federrath et al., 2008), but the velocity
field in the galactic scale simulation of Jin et al. (2017) indicates that the
solenoidal mode may dominate at large scales.

For each column density, we ran several simulations with different driving
strengths (as in section 4.1) and selected the simulation in which the SFR was
more similar to the SK relation. Figure 10 shows the result, where the
differences between the SFR from simulations with and without driving are
rather clear.

This proves that the turbulent driving from the large scale is a valid
candidate for closing the gap between the observations and the simulations
with stellar feedback alone. The question clearly is whether the added driving
is realistic in terms of injected energy and resulting speed dispersion, and
also regarding the consequence on the structure of the ISM.

One possible diagnostic is to compute the energy injected within the
simulation, as has been done in Paper I. Another diagnostic is to consider the
velocity dispersion. Both are discussed below.

###  6.3 Impact of the magnetic field on the Schmidt-Kennicutt relation

Figure 13: Mass-weighted velocity dispersion in the full computational volume
of 1 kpc3 computed at $t=58$ Myr. Left panel: Velocity parallel to the
galactic plane. Right panel: Vertical velocity.

The actual value of the magnetic field in galaxies (Beck, 2015; Han, 2017)
remains controversial because it is very difficult to measure. However, based
on energy equipartition, it is reasonable to think that gas-rich galaxies are
subject to stronger magnetic fields. As shown in section 5, a stronger
magnetic field results in a lower SFR.

It is thus interesting to investigate how the magnetic field, in conjunction
with turbulence, may contribute to shaping the SK relation. To investigate
this, we ran three series of simulations called noturb_hB, turb_hB, and
turb_Bvar. The description of all simulations can be found in Table 1. We
recall that noturb_hB and turb_hB are the same as noturb and turb,
respectively, but with doubled magnetic intensities. For the group turb_Bvar,
the initial magnetic field, parametrised by its midplane value $B_{0}$, scales
as the square root of the column density. In this way, the ratio of the
kinetic and magnetic energy at the beginning of the simulation remains the
same.

Figure 11 displays the SFR for the runs . In spite of the initial field, which
is twice as intense, the trend remains similar to the SFR of the noturb runs.
In particular, the SFR dependence on the column density is clearly too steep.

To investigate the impact that magnetic fields have on the SK relation in
conjunction with turbulence, we proceeded as with the turb group. For each
value of the initial column density and magnetisation, we made an educated
guess of what would be a good strength of the turbulent RMS acceleration to
obtain an SFR similar to the SK (equal to within a factor of 2). We then
launched the simulation, computed the SFR, and corrected the RMS acceleration
if necessary. The simulations we finally selected are described in Table 1 .
They produce a mean SFR that is more similar to the SK value. The goal is to
determine for an initial given magnetic field intensity and column density how
much turbulent energy is needed to reproduce SK.

First, Fig. 11 shows that with our choices of RMS acceleration, we reproduce
the SK relation for all groups, as expected. The SFR is globally slightly
lower for the turb_Bvar group, meaning that we slightly overestimated the
turbulent energy required for this group.

Second, Fig. 12 shows how much power is injected into the simulation by the
driven turbulence for all groups as a function of the column density. The
injected power was computed at each step by computing the difference between
the kinetic energy before and after applying the turbulent driving and
dividing by the time step. This reveals that the required power is
significantly reduced in the presence of substantial magnetic fields. Compared
to the series of runs turb, the power with which it is necessary to drive the
turbulence is more than a factor of ten lower when the magnetic field is
doubled, and the slope is roughly the same (group turb_hB). When the magnetic
field scales as the square root of the column density (group turb_Bvar), the
power-law slope of the power needed to reproduce the SK relation versus the
column density is lowered from almost 5 to 4. This means that magnetic fields
may significantly contribute to the origin of the SK relation. The values of
the magnetic field we considered are likely reasonable, but it is entirely
possible that stronger intensities have to be considered, which would further
reduce the required turbulent driving.

Figure 12 also shows the various estimates for the energy injection made in
section 2. In gas-rich galaxies, rotation could in principle provide enough
energy up to $\Sigma\simeq 100$ M⊙ pc-2. For higher values, however, this does
not seem to be the case, meaning that a source of energy that could explain
the SK relation for column densities higher than $100$ M⊙ pc-2 is currently
lacking.

###  6.4 Velocity dispersion

The velocity dispersion is also a clear and simple diagnostic to assess the
realism of simulations. Figure 13 displays the mass-weighted velocity
dispersion computed at the scale of the full kiloparsec box as a function of
column density for the various runs we performed.

Several trends are worth noticing. First, in the absence of turbulent driving
(runs noturb and noturb_hB), the velocity dispersion remains broadly isotropic
while in the presence of driving, particularly at high column densities, the
velocity dispersion in the disk plane is several times higher than the
velocity along the z-axis. Interestingly, in the presence of driving, the
velocity dispersion in the disk plane is higher than without driving, whereas
the reverse is true for the z-velocity. This is because in the latter case,
self-gravity is strong and triggers isotropic collapse in various places. In
the presence of strong magnetic fields (runs turb_hB and turb_Bvar), the
velocity dispersion necessary to reproduce the SK relations is lower by a
factor of about 2. Since the energy dissipation rate is expected to be
proportional to $\sigma^{3}$, this agrees well with the factor of about ten
induced for the driving power presented in Fig. 12.

It is important to stress that the estimated injected power of Fig. 12 and the
velocity dispersion of Fig. 13 are not predictions, but an estimate of the
strength of the turbulence needed to quench star formation down to the SK
relation. The way in which the turbulence is injected and the parameters of
the model, such as the strength of the magnetic field, may change these
requirement estimates (see sections 6.5 and 6.3, respectively).

Comparisons with observations are difficult because postprocessing steps
beyond the scope of the work are needed to ensure that the same quantities are
compared in terms of scale and tracked gas. Data from the ShiZELS galaxies
(Swinbank et al., 2012), where properties also measured at scales close to 1
kpc in the H$\alpha$ band, feature regions with a surface density of SFR
$\Sigma_{\mathrm{sfr}}$ around 0.1 $\leavevmode\nobreak\
\text{M}_{\astrosun}\cdot\mathrm{yr}^{-1}\mathrm{kpc}^{-2}$ and measured
velocity dispersion that can reach up to 200
$\mathrm{km}\cdot\mathrm{s}^{-1}$. Similar results were found by Law et al.
(2009) and in the WiggleZ survey (Wisnioski et al., 2011). The values from the
turb group, with a initial magnetic field at $3.8\leavevmode\nobreak\
\mu\mathrm{G}$, are quite extreme, indicating that for these column densities,
the magnetic field is probably stronger.

Figure 14: SFR surface density for the turb3D group compared to the turb_hB.
The legend is the same as in Fig. 10.

###  6.5 Three-dimensional driving

The high anisotropy of the velocity field found in the driven simulations is a
consequence of the 2D turbulent driving. The driving force is parallel to the
disk. Furthermore, it does not depend on the vertical coordinate because the
wave vectors of the Fourier decomposition of the force are also aligned with
the disk. Both of these assumptions are quite strong and can be relaxed.

Even though large-scale motions are expected to be aligned with the galactic
plane, turbulent eddies can spread in the vertical direction if their size is
below the scale height of the disk. In this last part, we try to modify our
driving model to take this into account. The driving force follows the same
equations as in section 3.4 with the following modifications.

We changed the power spectrum $F_{0}$ to add vertical modes with a wavelength
of about 32 pc (the disk scale height is initially 150 pc),

$F_{0}(\bm{k})=\begin{cases}\left(\dfrac{k}{2\pi}\right)^{-5/3}\text{ if
}2\leq\dfrac{k_{x}}{2\pi},\dfrac{k_{y}}{2\pi}\leq 3\text{ and
}30\leq\dfrac{k_{z}}{2\pi}\leq 32,\\\ \\\ 0\text{ if not.}\end{cases}$ (26)

The projection operator $\bm{P_{\zeta}}$ was the normal 3D projection operator
defined in Schmidt et al. (2009). However, we applied a factor of $0.3$ to the
vertical motion of the resulting force to avoid motions that were too strong
in the vertical direction that may destroy the disk.

We were able to find the value of the turbulent driving for the 2D driving
that leads to SFR compatible with the SK relation (Fig 14, left panel).
Compared with the most similar set of simulation, turb_hB, the energy needed
to achieve SK is slightly increased, as we show in Fig. 12. However, Fig. 13
shows that the velocity dispersion in the disk is lower, as is the anisotropy.
The vertical velocity dispersion does not increase significantly, probably
because of the strong attenuation of the vertical component of the force. With
this experiment, a less anisotropic velocity dispersion can be produced by
releasing the assumption on 2D turbulence. When more accurate constraints on
the turbulent driving are available, the presence and strength of a vertical
component are an important item to monitor.

##  7 Caveats

###  7.1 Shear and Coriolis effect

The fiducial model includes neither shear nor the Coriolis effect because the
driven turbulence module and the shearing box used in aprevious study (
Colling et al. (2018)) require incompatible boundary conditions. Our 2D and
mainly solenoidal turbulent driving reproduces some features of the shear and
the Coriolis effect, but in a less self-consistent way. In its current state,
it does not allow us to distinguish between the effect of large-scale
unordered turbulence and ordered rotation.

The effect of shear and the Coriolis effect were studied by Colling et al.
(2018) and Kim et al. (2020). Both studies showed that a shear value
comparable to it did not reduce the SFR enough to match the SK relation and
just slightly improved the discrepancy observed in group turb in Fig. 10.
Colling et al. (2018) showed that for an initial column density of
$19.4\leavevmode\nobreak\ \leavevmode\nobreak\
\text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$, the shear associated with an
angular speed of $28\leavevmode\nobreak\
\mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{kpc}^{-1}$ , which is typical of
the solar neighbourhood, only reduced the SFR by a factor $1.2$. This is
almost negligible given the high uncertainties of the type of computations. On
the other hand, a higher angular speed of $56\leavevmode\nobreak\
\mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{kpc}^{-1}$, which may be found in
denser environments at a galactocentric radius around $4\leavevmode\nobreak\
\mathrm{kpc}$, leads to a reduction of the SFR of a factor $3.6$. The
turbulent driving we chose for the same initial column density of
$19.4\leavevmode\nobreak\ \leavevmode\nobreak\
\text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$ reduced the SFR by a factor $4.5$
and led to an SFR below the SK relation. This means that a region with this
moderate column density undergoing a strong shear would not need additional
driving. The question that arises is whether the shear may reduce the much
higher energy needed to reproduce the SK law (see Fig. 12) for the cases with
a high column density. While we cannot answer this question for the shear with
our current setup, it is possible to answer this for the Coriolis effect
alone, as we did in Appendix D. The conclusion is that for a high column
density of $77\leavevmode\nobreak\ \text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$,
the Coriolis effect has a limited impact on the SFR and cannot replace or
complete the large-scale turbulent driving to match the SK relation.

###  7.2 Nature of the turbulence

Following this line of thought, an important point to keep in mind is that the
turbulence injected in groups turb, turb_Bvar, turb_hB, and turb3D is
calibrated so that the resulting SFR stays close to the observed value at a
given column density and is not generated self-consistently. We showed in
section 4 that the exact nature of the turbulence matters greatly for the SFR.
The choice of the Ornstein-Uhlenbeck process makes it hard to model a very
coherent driving. These drivings may produce velocity fields that are even
more solenoidal that the field yielded by a purely solenoidal Ornstein-
Uhlenbeck driving (e.g. with the magneto-rotational instability in Gong et
al., 2020, Fig. 2). To proceed with more precise studies of gas-rich galaxies,
a better modelling of the driving source is required.

###  7.3 Sources of feedback

Several sources of stellar feedback are not included in the model, such as
stellar winds and cosmic rays. Both have been studied in other works (Gatto et
al. (2017) and Rathjen et al. (2021), respectively). Stellar feedback effects
do not add up linearly for the SFR: including at least one kind of early
feedback (photo-ionisation, stellar winds, or infrared irradiation) may change
the result dramatically, but combining them does not change the SFR by more
than a factor of three.

###  7.4 Integration time

Another possible caveat is that unlike other similar studies (e.g. Kim &
Ostriker (2017)), the time span over which the simulations were run was only
150 - 200 Myr. This ensured that the total mass of gas within the box did not
vary by more than 40% because of outflows and star formation. Previous studies
have often solved the issue of outflowing gas by using very elongated boxes
that kept a large fraction of the outflowing gas within the computational
domain, maybe artificially.

###  7.5 Numerical parameters

Finally, as always with numerical simulations, the influence of the choice of
the numerical parameters may be a source of concern. This includes the
influence of the sink creation criterion, which was discussed in Colling et
al. (2018). They showed in their figure 11b that the SFR is unchanged even
when they modified the threshold for sink accretion by a factor 16, but the
addition of a test to verify the boundedness of the gas was able to increase
the SFR by almost a factor of two. We discuss the details of the SN recipe in
Appendix B. The limitation we applied for numerical reasons to the SN subgrid
model altered the properties of the low-density gas, but had a limited impact
on the dense gas and the SFR. We present the convergence study with resolution
in Appendix C and show that an increase in the resolution only modifies the
SFR by about 10%.

##  8 Discussion and conclusions

We extended the study of Paper I of what sets the SFR at the kiloparsec scale.
Using a very similar numerical setup of a kiloparsec cube region of a galaxy,
we explored the parameter space, and more particularly, the effect of the
compressibility of large-scale turbulence and the strength of the magnetic
field. In Paper I we showed that stellar feedback alone was not sufficient to
quench star formation to match the observed SFRs, and that large-scale
turbulent driving may be the main effect that is missing, explaining the gap.
To refine the study, both the mechanism of turbulence injection (section 3.4)
and the computation of the SFR (section 3.5) were improved to allow a better
comparison with the observed SFR. Compared to the reworked method for
turbulence injection with the method used in Paper I, quenching via large-
scale turbulence is more efficient in low column density regions and slightly
less efficient in high column density regions.

The main results of this study are the following:

  1. 1.

In high column density environments, stellar feedback is unable to destroy a
large accumulation of dense gas and quench star formation efficiently enough
to yield a rate that is compatible with observations (section 6.1).

  2. 2.

If it is strong enough, 2D large-scale turbulent driving can destroy these
giant clouds and can efficiently quench star formation down to levels that are
compatible with observations (section 6.2).

  3. 3.

Increasing the 2D turbulent driving force linearly increases the velocity
dispersion parallel to the disk up to a maximum limit at which it saturates.
Below this saturation limit, stronger turbulent driving translates into a
lower SFR (section 4.1). In the case of a low magnetic field, the velocity
dispersion associated with the turbulent-driving strength that is required to
match the SK relation is high compared with observations (section 6.4).

  4. 4.

The compressibility of the turbulent driving force matters as more compressive
driving schemes are ten times less efficient than solenoidal ones in quenching
star formation (section 4.2).

  5. 5.

An increased magnetic field can also have a dramatic effect on the SFR, where
an increase of roughly $10\leavevmode\nobreak\ \mu G$ reduces the SFR by a
factor of ten (section 5). It cannot explain the observed SFR with stellar
feedback alone, but it reduces the amount of energy needed from large-scale
driving (section 6.3) and the associated velocity dispersion (section 6.4).

With this work, we provided an overview of how the large-scale driving and the
mean value of the magnetic field can influence the SFR. The next step towards
comprehending what regulates the SFR at the kiloparsec scale in general and
the role of the turbulent driving in particular is to obtain more precise
constraints on the energy that can be tapped from the large scale, as well as
the nature of the driving.

Several analytical works have investigated the possible source of this
turbulence, gravitational instabilities (Nusser & Silk, 2022), mass transfer
within the galactic disk (Krumholz et al., 2018; Meidt et al., 2020), or
accretion onto the disk from the circumgalactic medium. Observations from
local (e.g. Sun et al., 2022) or high-redshift galaxies, as well as isolated
and cosmological simulations, are other ways to better know how energy is
transferred from large to small scales.

## Acknowledgements

The authors thank the referee, Enrique Vazquez-Semadeni, for his useful
comments which helped to significantly improve the article. The authors
acknowledge Interstellar Institute’s program ”With Two Eyes” and the Paris-
Saclay University’s Institut Pascal for hosting discussions that nourished the
development of the ideas behind this work. NB thanks Jérémy Fensch, Chang Goo
Kim and Munan Gong for interesting and fruitful discussions. This work was
granted access to HPC resources of CINES and CCRT under the allocation
x2020047023 made by GENCI (Grand Equipement National de Calcul Intensif) and
the special allocation 2021SA10spe00010. NB, PH and TC acknowledges financial
support from the European Research Council (ERC) via the ERC Synergy Grant
”ECOGAL: Understanding our Galactic ecosystem – From the disk of the Milky Way
to the formation sites of stars and planets” (grant 855130).

#### Software

We made use of the following software and analysis tools: GNU/Linux, Ramses
(Teyssier, 2002), Python, Matplotlib (Hunter, 2007), Numpy (van der Walt et
al., 2011), Pymses (Guillet et al., 2013), Astrophysix, Astropy (Astropy
Collaboration et al., 2013, 2018). Many thanks to the authors for making them
publicly available. Special thanks to all the contributors of ramses, and in
particular to Andrew McLeod for the introduction of the turbulence driving
module that has been heavily used in this work.

#### Data availability

The data underlying this article are available in the Galactica Database at
http://www.galactica-simulations.eu, and can be accessed with the unique
identifier ISMFEED. Additional data and the source code used to run
simulations and perform analysis will be shared on reasonable request to the
corresponding author.

## References

  * Agertz et al. (2009) Agertz, O., Lake, G., Teyssier, R., et al. 2009, MNRAS, 392, 294
  * Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, ApJ, 156, 123
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33
  * Audit & Hennebelle (2005) Audit, E. & Hennebelle, P. 2005, A&A, 433, 1
  * Beck (2015) Beck, R. 2015, A&AR, 24, 4
  * Besserglik & Goldman (2021) Besserglik, D. & Goldman, I. 2021, ApJ, 915, 117
  * Bleuler & Teyssier (2014) Bleuler, A. & Teyssier, R. 2014, MNRAS, 445, 4015
  * Bournaud et al. (2010) Bournaud, F., Elmegreen, B. G., Teyssier, R., Block, D. L., & Puerari, I. 2010, MNRAS, 409, 1088
  * Brucy et al. (2020) Brucy, N., Hennebelle, P., Bournaud, F., & Colling, C. 2020, ApJL, 896, L34
  * Chepurnov et al. (2015) Chepurnov, A., Burkhart, B., Lazarian, A., & Stanimirovic, S. 2015, ApJ, 810, 33
  * Colling et al. (2018) Colling, C., Hennebelle, P., Geen, S., Iffrig, O., & Bournaud, F. 2018, A&A, 620, A21
  * Colman et al. (2022) Colman, T., Robitaille, J.-F., Hennebelle, P., et al. 2022, MNRAS, 514, 3670
  * Daddi et al. (2010) Daddi, E., Bournaud, F., Walter, F., et al. 2010, ApJ, 713, 686
  * de Avillez & Breitschwerdt (2005) de Avillez, M. A. & Breitschwerdt, D. 2005, A&A, 436, 585
  * Dib et al. (2021) Dib, S., Braine, J., Gopinathan, M., et al. 2021, A&A, 655, A101
  * Dubois & Teyssier (2008) Dubois, Y. & Teyssier, R. 2008, A&A, 477, 79
  * Eswaran & Pope (1988) Eswaran, V. & Pope, S. B. 1988, Computers and Fluids, 16, 257
  * Federrath (2013) Federrath, C. 2013, MNRAS, 436, 1245
  * Federrath & Klessen (2012) Federrath, C. & Klessen, R. S. 2012, ApJ, 761, 156
  * Federrath et al. (2008) Federrath, C., Klessen, R. S., & Schmidt, W. 2008, ApJ, 688, L79
  * Federrath et al. (2010) Federrath, C., Roman-Duval, J., Klessen, R. S., Schmidt, W., & Mac Low, M. M. 2010, A&A, 512, A81
  * Fensch & Bournaud (2021) Fensch, J. & Bournaud, F. 2021, MNRAS, 505, 3579
  * Fensch et al. (2023) Fensch, J., Bournaud, F., Brucy, N., et al. 2023, A&A, 672, A193
  * Forbes et al. (2023) Forbes, J. C., Emami, R., Somerville, R. S., et al. 2023, ApJ, 948, 107
  * Fromang et al. (2006) Fromang, S., Hennebelle, P., & Teyssier, R. 2006, A&A, 457, 371
  * Gatto et al. (2017) Gatto, A., Walch, S., Naab, T., et al. 2017, MNRAS, 466, 1903
  * Genzel et al. (2008) Genzel, R., Burkert, A., Bouché, N., et al. 2008, ApJ, 687, 59
  * Genzel et al. (2010) Genzel, R., Tacconi, L. J., Gracia-Carpio, J., et al. 2010, MNRAS, 407, 2091
  * Girichidis et al. (2018) Girichidis, P., Seifried, D., Naab, T., et al. 2018, MNRAS, 480, 3511
  * Goldbaum et al. (2015) Goldbaum, N. J., Krumholz, M. R., & Forbes, J. C. 2015, ApJ, 814, 131
  * Goldbaum et al. (2016) Goldbaum, N. J., Krumholz, M. R., & Forbes, J. C. 2016, ApJ, 827, 28
  * Gong et al. (2020) Gong, M., Ivlev, A. V., Zhao, B., & Caselli, P. 2020, ApJ, 891, 172
  * Guillet et al. (2013) Guillet, T., Chapon, D., & Labadens, M. 2013, Astrophysics Source Code Library, ascl:1310.002
  * Han (2017) Han, J. L. 2017, ARA&A, 55, 111
  * Heitsch et al. (2007) Heitsch, F., Slyz, A. D., Devriendt, J. E. G., Hartmann, L. W., & Burkert, A. 2007, ApJ, 665, 445
  * Hennebelle & Chabrier (2011) Hennebelle, P. & Chabrier, G. 2011, ApJ, 743, L29
  * Hennebelle & Iffrig (2014) Hennebelle, P. & Iffrig, O. 2014, A&A, 570, A81
  * Hennebelle & Inutsuka (2019) Hennebelle, P. & Inutsuka, S.-i. 2019, Frontiers in Astronomy and Space Sciences, 6, 5
  * Hill et al. (2012) Hill, A. S., Joung, M. R., Mac Low, M.-M., et al. 2012, ApJ, 750, 104
  * Hopkins et al. (2011) Hopkins, P. F., Quataert, E., & Murray, N. 2011, MNRAS, 417, 950
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science and Engineering, 9, 90
  * Iffrig & Hennebelle (2015) Iffrig, O. & Hennebelle, P. 2015, A&A, 576, A95
  * Iffrig & Hennebelle (2017) Iffrig, O. & Hennebelle, P. 2017, A&A, 604, A70
  * Jin et al. (2017) Jin, K., Salim, D. M., Federrath, C., et al. 2017, MNRAS, 469, 383
  * Joung & Mac Low (2006) Joung, M. K. R. & Mac Low, M.-M. 2006, ApJ, 653, 1266
  * Kannan et al. (2020) Kannan, R., Marinacci, F., Simpson, C. M., Glover, S. C. O., & Hernquist, L. 2020, MNRAS, 491, 2088
  * Kennicutt & Evans (2012) Kennicutt, R. C. & Evans, N. J. 2012, ARA&A, 50, 531
  * Kennicutt (1998) Kennicutt, Jr., R. C. 1998, ApJ, 498, 541
  * Kennicutt & De Los Reyes (2021) Kennicutt, Jr., R. C. & De Los Reyes, M. A. C. 2021, ApJ, 908, 61
  * Kim et al. (2011) Kim, C.-G., Kim, W.-T., & Ostriker, E. C. 2011, ApJ, 743, 25
  * Kim & Ostriker (2017) Kim, C.-G. & Ostriker, E. C. 2017, ApJ, 846, 133
  * Kim et al. (2020) Kim, C.-G., Ostriker, E. C., Fielding, D. B., et al. 2020, ApJ, 903, L34
  * Kim et al. (2013) Kim, C.-G., Ostriker, E. C., & Kim, W.-T. 2013, ApJ, 776, 1
  * Klessen & Hennebelle (2010) Klessen, R. S. & Hennebelle, P. 2010, A&A, 520, A17
  * Kolmogorov (1941) Kolmogorov, A. 1941, Akademiia Nauk SSSR Doklady, 30, 301
  * Koyama & Inutsuka (2000) Koyama, H. & Inutsuka, S.-I. 2000, ApJ, 532, 980
  * Krumholz & Burkhart (2016) Krumholz, M. R. & Burkhart, B. 2016, MNRAS, 458, 1671
  * Krumholz et al. (2018) Krumholz, M. R., Burkhart, B., Forbes, J. C., & Crocker, R. M. 2018, MNRAS, 477, 2716
  * Kuijken & Gilmore (1989) Kuijken, K. & Gilmore, G. 1989, MNRAS, 239, 605
  * Law et al. (2009) Law, D. R., Steidel, C. C., Erb, D. K., et al. 2009, ApJ, 697, 2057
  * Mac Low & Klessen (2004) Mac Low, M.-M. & Klessen, R. S. 2004, Reviews of Modern Physics, 76, 125
  * Meidt et al. (2020) Meidt, S. E., Glover, S. C. O., Kruijssen, J. M. D., et al. 2020, ApJ, 892, 73
  * Molina et al. (2012) Molina, F. Z., Glover, S. C. O., Federrath, C., & Klessen, R. S. 2012, MNRAS, 423, 2680
  * Nusser & Silk (2022) Nusser, A. & Silk, J. 2022, MNRAS, 509, 2979
  * Ostriker & Kim (2022) Ostriker, E. C. & Kim, C.-G. 2022, ApJ, 936, 137
  * Ostriker et al. (2010) Ostriker, E. C., McKee, C. F., & Leroy, A. K. 2010, ApJ, 721, 975
  * Padoan & Nordlund (2011) Padoan, P. & Nordlund, A. 2011, ApJ, 730, 40
  * Padoan et al. (2016) Padoan, P., Pan, L., Haugbolle, T., & Nordlund, A. 2016, ApJ, 822, 11
  * Piontek & Ostriker (2007) Piontek, R. A. & Ostriker, E. C. 2007, ApJ, 663, 183
  * Rathjen et al. (2021) Rathjen, T.-E., Naab, T., Girichidis, P., et al. 2021, MNRAS, 504, 1039
  * Rathjen et al. (2023) Rathjen, T.-E., Naab, T., Walch, S., et al. 2023, MNRAS, 522, 1843
  * Renaud et al. (2012) Renaud, F., Kraljic, K., & Bournaud, F. 2012, ApJL, 760, L16
  * Rosdahl et al. (2013) Rosdahl, J., Blaizot, J., Aubert, D., Stranex, T., & Teyssier, R. 2013, MNRAS, 436, 2188
  * Rosen & Bregman (1995) Rosen, A. & Bregman, J. N. 1995, ApJ, 440, 634
  * Salpeter (1955) Salpeter, E. E. 1955, ApJ, 121, 161
  * Schmidt et al. (2009) Schmidt, W., Federrath, C., Hupp, M., Kern, S., & Niemeyer, J. C. 2009, A&A, 494, 127
  * Schmidt et al. (2006) Schmidt, W., Hillebrandt, W., & Niemeyer, J. C. 2006, Computers & Fluids, 35, 353
  * Sun et al. (2022) Sun, J., Leroy, A. K., Rosolowsky, E., et al. 2022, ApJ, 164, 43
  * Swinbank et al. (2012) Swinbank, A. M., Sobral, D., Smail, I., et al. 2012, MNRAS, 426, 935
  * Szotkowski et al. (2019) Szotkowski, S., Yoder, D., Stanimirović, S., et al. 2019, ApJ, 887, 111
  * Teyssier (2002) Teyssier, R. 2002, A&A, 385, 337
  * Toomre (1964) Toomre, A. 1964, ApJ, 139, 1217
  * Vacca et al. (1996) Vacca, W. D., Garmany, C. D., & Shull, J. M. 1996, ApJ, 460, 914
  * van der Walt et al. (2011) van der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Computing in Science and Engineering, 13, 22
  * Vazquez-Semadeni et al. (1996) Vazquez-Semadeni, E., Passot, T., & Pouquet, A. 1996, ApJ, 473, 881
  * Wada et al. (2002) Wada, K., Meurer, G., & Norman, C. A. 2002, ApJ, 577, 197
  * Walch et al. (2015) Walch, S., Girichidis, P., Naab, T., et al. 2015, MNRAS, 454, 238
  * Wisnioski et al. (2011) Wisnioski, E., Glazebrook, K., Blake, C., et al. 2011, MNRAS, 417, 2601
  * Woosley et al. (2002) Woosley, S. E., Heger, A., & Weaver, T. A. 2002, Reviews of Modern Physics, 74, 1015
  * Zamora-Avilés et al. (2018) Zamora-Avilés, M., Vázquez-Semadeni, E., Körtgen, B., Banerjee, R., & Hartmann, L. 2018, MNRAS, 474, 4824

##  Appendix A Calibration of the turbulent driving

In section 4.2 we studied the effect of changing the solenoidal fraction
$\zeta$ of the driving on the SFR. In order to have a fair comparison, it is
important to ensure that the strength of the driving is the same, regardless
of $\zeta$. A way to check this is to consider the RMS power of the Fourier
modes,

$\mathrm{RMS(t)}=\sqrt{\int\left|\bm{\hat{f}}(\bm{k},t)\right|^{2}\mathrm{d}^{3}\bm{k}}.$
(27)

The time-averaged RMS depends on the solenoidal fraction (see Fig. 15), and
the role of the normalization factor $g_{\zeta}$ of Eq. (23) is to correct for
this. Because we used 2D driving instead of 3D driving, we cannot use the
normalisation factor given by Schmidt et al. (2009) in their Equation (9).
Instead, we empirically computed the relevant factor by running the turbulent
generation code without any normalisation and compared the measured RMS in the
Fourier space with the input $f_{\mathrm{rms}}$ (see Fig. 15). A fit of the
obtained curve then gives the $g_{\zeta}$ factor we used in the code,

$g_{\zeta}=0.93(1-\zeta)^{2}-0.87(1-\zeta)+0.76.$ (28)

Figure 15: Ratio of the time-averaged RMS (as in Eq. (27)) over the
$f_{\mathrm{rms}}$ parameter for unnormalised simulations. The fit of the
curve gives the normalisation factor $g_{\zeta}$.

##  Appendix B Limitation of the supernova feedback

Figure 16: Averaged ratio of the injected momentum over the expected reference
momentum of $4\times 10^{43}$ g $\cdot$ cm $\cdot$ s -1 for all the SNe that
blew up during the simulation. For simulations with higher initial column
densities, the injected momentum approaches the reference value.

Figure 17: Comparison of the SFR with (blue line) and without (orange line)
limiting the velocity of the SN ejecta for $n_{0}=1.5$ (top) and $n_{0}=6$
(bottom).

Figure 18: Comparison of the edge-on (top) and face-on (bottom) column density
map with and without (left and right of each panel, respectively) limiting the
velocity of the SN ejecta for $n_{0}=1.5$ (left panel) and $n_{0}=6$ (right
panel). Figure 19: Mass-weighted probability distribution function of the gas
density for the same snapshots as in Fig. 18.

When SNe explode in a low-density environment, they can generate very high
velocities and temperatures. This leads to a numerical issue because it
implies very short time-steps, to a point that the simulation can hardly
progress. In order to avoid this in our simulation, the energy from SNe was
mainly injected as a momentum, and both the gas temperature and the gas
velocity cannot exceed $\mathrm{T}_{\mathrm{sat}}=10^{6}\leavevmode\nobreak\
\mathrm{K}$ and $\mathrm{V}_{\mathrm{sat}}=300\leavevmode\nobreak\
\mathrm{km}\cdot\mathrm{s}^{-1}$. As a consequence, we did not model the very
hot phase of the gas, which, despite representing only a small fraction of the
total mass of the ISM, fills the majority of the volume. However, the hot gas
cools down very slowly and thus cannot collapse to form stars at the timescale
we studied. Furthermore, with our 4 pc resolution, the cooling radius is not
resolved for most of the SNe, especially for the most massive simulations.

Limiting the velocity means that overall, the momentum injected is lower than
prescribed by the SN model, as we show in Fig. 16. The figure depicts the
ratio of the total momentum effectively injected by SNe over the total
expected momentum, which is the number of SNe multiplied by the reference
value of $4\times 10^{43}$ g $\cdot$ cm $\cdot$ $\mathrm{s}^{-1}$. On average,
$30$ to $60\%$ of the reference momentum is injected, with a higher ratio in a
dense environment. This may seem very low, but it is worth noting that SNe
exploding in a low-density environment have a quasi-negligible impact on dense
gas, as was demonstrated by Iffrig & Hennebelle (2015).

In order to verify that the limitation of the velocity of the gas does not
change the outcome of the SFR significantly, we ran two simulations analogous
to noturb with $n_{0}=1.5$ and $6$ without the limitation. When the SNe
ignited, the time step of the simulation without the limited velocity was ten
times lower than the fiducial time step. Therefore, the simulations without a
limited velocity evolved for a shorter time. Since our main focus is the
computation SFR, we compared the amount of gas accreted by the sink particles
in both cases in Fig. 17. It shows that for the period we considered, the
limitation of the velocity due to SN feedback leads to an overestimation of
the SFR by a factor $1.6$. This is comparable to the uncertainty due to other
factors, such as specific choices of the sink creation recipe, implementation
of other stellar feedback, or resolution. We were unable to run the simulation
longer, but similar simulations run by Kim & Ostriker (2017) and Ostriker &
Kim (2022) with an elongated box in the vertical direction indicates that the
SFR is bound to oscillate, with bursts of similar or slightly lower amplitude
as the first event of star formation depicted in Fig. 17.

The column density map (Fig. 18) reveals that our fiducial simulation fails to
capture an important feature of the gas structure: superbubbles of hot gas.
However, the dense gas does not seem to be strongly affected by the change in
the SN feedback recipe: it is displaced, but its internal structure and its
quantity are roughly the same. The latter point can be verified through the
mass-weighted probability distribution function of the gas density shown in
Fig. 19. This figure shows that the amount of very low density gas increases
when the SN speed is not limited, but the amount of dense gas stay the same.

##  Appendix C Convergence with the resolution

All the simulations we discussed were run with a uniform grid resolution of
about $4$ pc, with $256^{3}$ cells. This allowed the wide parameter study
presented here, with 63 simulations discussed out of a total of more than 200.
It is important, however, to determine the influence of the resolution and
whether convergence is reached. To this end, we selected two simulations
(group turb with $n_{0}=6$ and group noturb with $n_{0}=6$) and reran them
with a doubled resolution ($512^{3}$ cells of size $\sim$ 2 pc) and half
resolution ($128^{3}$ cells of size $\sim$ 8 pc). The sink creation threshold
$n_{\mathrm{sink}}$ was adapted so that
$n_{\mathrm{sink}}\propto\dfrac{1}{|\Delta x|^{2}}$ , where $\Delta x$ is the
size of a cell (Kim & Ostriker 2017).

Figure 20 shows the mass accreted by sinks in these simulations. In the case
without turbulent driving (top), increasing the resolution results in a
increased SFR. The case with driving is less clear, especially because star
formation is more stochastic in this configuration. In both cases, the
variations in the amount of mass accreted by the sinks and the resulting SFR
(Fig. 21) due to the resolution are about 10%. This is small compared with the
variations measured when changing the physical parameters (sections 4 and 5).
It is within the same range as the other source of uncertainties, such as the
choice of other numerical parameters (sink creation threshold and seed for the
turbulence) or temporal variations. This convergence result is in line with
resolution studies that were carried out for a similar setup at a lower column
density ((Colling et al. 2018, Fig. 9) and (Kim & Ostriker 2017, Fig. 14)).

Figure 22 shows that similarly, the variation in the mass-weighted velocity
dispersion are small between our fiducial resolution of $4$ pc and a doubled
resolution of $2$ pc.

Figure 20: Comparison of the total mass accreted by sinks as a function of
time for different resolutions for a simulation with turbulence driving (top)
and a simulation without it (bottom).

Figure 21: Comparison of the computed SFR for different resolutions for a
simulation with turbulence driving (top) and a simulation without it (bottom).
The dashed line is the value of the SFR given by the SK relation.  Figure 22:
Comparison of the mass-weighted velocity dispersion for different resolutions
for a simulation with turbulence driving at $t=72$ Myr.

##  Appendix D Impact of the Coriolis effect on the SFR

As stated in section 7, our model does not include shear, mainly for numerical
reasons. However, it is rather easy to implement the Coriolis effect in our
simulations and assess how much it impacts the SFR. To make this experiment,
we added a source term $\rho\bm{f_{\mathrm{cor}}}$ in the momentum equation
(16). We have

$\bm{f_{\mathrm{cor}}}=-2\Omega\bm{e_{z}}\times\bm{v},$ (29)

where $\bm{e_{z}}$ is a unit vector in the vertical direction, and $\Omega$ is
the angular speed, a constant and uniform parameter in our setup. In practice,
we updated the velocity field after each hydro and gravity step using an
implicit Crank Nicholson scheme,

$\displaystyle v_{\mathrm{x,\mathrm{new}}}$
$\displaystyle=\dfrac{v_{\mathrm{x,\mathrm{old}}}\left(1-\left(\Omega
dt\right)^{2}\right)+2v_{\mathrm{y,\mathrm{old}}}\Omega dt}{\left(\Omega
dt\right)^{2}},$ (30) $\displaystyle v_{\mathrm{y,\mathrm{new}}}$
$\displaystyle=\dfrac{v_{\mathrm{y,\mathrm{old}}}\left(1-\left(\Omega
dt\right)^{2}\right)-2v_{\mathrm{x,\mathrm{old}}}\Omega dt}{\left(\Omega
dt\right)^{2}},$ (31) $\displaystyle v_{\mathrm{z,\mathrm{new}}}$
$\displaystyle=v_{\mathrm{z,\mathrm{old}}}.$ (32)

The implementation was similar to that used in Colling et al. (2018) and was
dutifully tested.

We chose an initial column density of $\Sigma_{0}=77.4\leavevmode\nobreak\
\leavevmode\nobreak\ \text{M}_{\astrosun}\cdot\mathrm{pc}^{-2}$ and an initial
magnetic intensity $B_{0}=3.8\leavevmode\nobreak\ \mu\mathrm{G}$. We
considered values of the angular speed of $27.5$ and $110\leavevmode\nobreak\
\mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{kpc}^{-1}$. When we assume a
constant rotation velocity of $220\leavevmode\nobreak\
\mathrm{km}\cdot\mathrm{s}^{-1}$, these values would correspond to the angular
speed at $8$ and $2\leavevmode\nobreak\ \mathrm{kpc}$, respectively. As a
first step, we ran the simulations without additional large-scale turbulent
driving222Self-consistent SN driving is still present. to determine whether
the Coriolis effect constituted an alternative mechanism to the large-scale
driving.

Fig. 23 (top) shows that the Coriolis effect slightly stabilises the disk,
with a reduction of the SFR of nearly a factor 1.5 with
$\Omega=110\leavevmode\nobreak\
\mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{kpc}^{-1}$ compared to the
simulation without the Coriolis effect. However, the resulting SFR is still
much higher than the observed values. This means that even with a very strong
Coriolis effect, an additional source of driving is still needed.

This result is slightly surprising. Because the angular speed $\Omega$ is
constant and uniform in our simulation domain, the implied Toomre parameter is
given by

$Q=\dfrac{2\Omega\sigma}{\pi G\Sigma}.$ (33)

If it is computed at the box scale ($1$ kpc), the velocity dispersion
$\sigma=\sqrt{\sigma_{x}^{2}+\sigma_{y}^{2}+\sigma_{z}^{2}}/\sqrt{3}$ is
about$10\leavevmode\nobreak\ \mathrm{km}\cdot\mathrm{s}^{-1}$ at
$20\leavevmode\nobreak\ \mathrm{Myr}$ for the two runs. For the chosen values
of the angular speed $\Omega=27.5$ and $110\leavevmode\nobreak\
\mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{kpc}^{-1}$, we obtain $Q=0.5$ and
$2.1$, respectively. We would expect the disk to be stabilised against
gravitational instabilities for the simulation with
$\Omega=110\leavevmode\nobreak\
\mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{kpc}^{-1}$, which has a value of
the Toomre parameter $Q$ above one, the stability limit. However, because of
the density fluctuations and because the velocity dispersion is lower when it
is computed at a smaller scale, it is possible to reach lower values of Q
locally. Fig. 24 shows the value of Q when it is computed in squares of
$100\leavevmode\nobreak\ \mathrm{pc}\times 100\leavevmode\nobreak\
\mathrm{pc}$ at $t=20\mathrm{Myr}$. For all the tested values of the angular
speed $\Omega$, at least one region has a Toomre parameter well below $1$. It
is also important to have in mind that the value of the velocity dispersion
$\sigma$ used in Eq. (33) also encompasses compressive motions, which tend to
enhance star formation.

Because of the slight stabilisation effect observed, we may also wonder
whether the additional turbulent energy required to reproduce the SK relation
discussed in section 6 is lowered when the Coriolis force is added. To test
this, we ran the same simulations, but with a turbulent driving
$f_{\mathrm{rms}}$ of $10^{5}$ (code units) with a compressibility
$\chi=0.25$. We compared them with the run with the same parameters in the
group strength, which has no Coriolis effect and yields an SFR that is ten
times higher that the SFR predicted by the SK relation. Fig. 23 (bottom) shows
that the Coriolis effect has almost no impact on the SFR when it is combined
with the large-scale turbulent driving. Despite a slightly different star
formation history, the slope of the curve is indeed roughly similar,
regardless of the angular speed, without any systematic trend. As a
consequence, we find that the additional turbulent energy required to
reproduce the SK relation is not affected by the presence or absence of the
Coriolis effect. A more complete study of the exact role of the Coriolis force
and why in this case it does not seem to play a major role in the stability of
the disk is beyond the scope of this article, but can be the topic of further
research work.

Figure 23: Comparison of the accreted mass in sink particles for different
angular speeds $\Omega$ (in
$\mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{kpc}^{-1}$) for the simulations
without turbulent driving (top) and with turbulent driving $f_{\mathrm{rms}}$
of $10^{5}$ in code units (bottom). The time is rescaled so that 0 Myr
corresponds to the time at which $10^{5}\leavevmode\nobreak\
\text{M}_{\astrosun}$ were accreted into sinks.

Figure 24: Value of the Toomre parameter $Q$ (Eq. (33)) computed in a square
of $100\leavevmode\nobreak\ \mathrm{pc}\times 100\leavevmode\nobreak\
\mathrm{pc}$ for the simulations without turbulent driving (top) and with
turbulent driving (bottom) at $t=20$ Myr and for different angular speeds
$\Omega$ (in $\mathrm{km}\cdot\mathrm{s}^{-1}\cdot\mathrm{kpc}^{-1}$).