Numerical Methods for Simulating Star Formation

1. Introduction

Star formation is one of the main unsolved problems in astrophysics. Although our view of this fundamental process, at the nexus of galaxy formation, planetary science and stellar evolution, has considerably changed over the past decades, thanks to new observations and theoretical progress, many dark corners remain to be explored. One of the reasons why the true origin of stars still eludes us is the highly non-linear nature of the governing equations, describing self-gravitating, compressible, magnetized dust, and gas fluids interacting with radiation. In this review, we present these main governing equations, focusing on ideal magneto-hydrodynamics, radiation hydrodynamics, non-ideal effects and sink particle techniques. We describe the most popular discretisation schemes, focusing on possible sources of errors and their consequences on clouding our conclusions. We finally give a short description of the landscape in term of star formation simulations, with different set-up strategies addressing particular scales in this fundamentally multi-scale problem.

Star formation is believed to be the consequence of the collapse of mildly supersonic or transonic to subsonic molecular cores emerging out of supersonic turbulent flows in the interstellar medium (Mac Low and Klessen, 2004; McKee and Ostriker, 2007). The source of turbulence is probably to be found on large scales, as a consequence of large galactic shearing or colliding flows, but also on small scales, because of various sources of stellar feedback (McKee, 1989; Federrath et al., 2017). In this context, gravitational or thermal instabilities lead to the formation of dense gas clumps that undergo a more or less violent gravitational collapse, leading to the formation of a proto-star surrounded by a proto-stellar disk. Describing these processes using only simple analytical methods is almost impossible. Moreover, the traditional engineering methods in Computational Fluid Dynamics (CFD) are usually not robust enough to sustain highly supersonic and super-Alfvénic flows. Self-gravity, magnetic fields and radiation fields, taken together, define a very unique system of equations that has no equivalent in the industry. This has led astrophysicists to develop their own methods, largely inspired by traditional methods designed in the applied mathematics community, but significantly adapted to the specific problem at hand. In this context, Smoothed Particle Hydrodynamics (SPH) and Adaptive Mesh Refinement (AMR) techniques turned out to be particularly suitable to the problem of star formation. Both techniques have their own pros and cons, and comparing the two allows us to assess the robustness of our numerical results. New techniques have also been developed, that are specific to star formation, such as sink particles, a commonly adopted recipe to replace an otherwise collapsing fluid element by a collision-less particle, saving computational times and increasing the realism of the simulation.

In this review, we pay attention to the description of the equations, without necessarily discussing their physical foundations, such as the ideal MHD limit or the non-ideal diffusion processes. We describe the various numerical techniques, from low-order schemes to modern high-order methods, as well as from non-zero divergence schemes to exact divergence-free methods, etc. We refer to the corresponding literature, including all references that are relevant to the historical description of the discipline and that give a fair snapshot of the present state of the field. We apologize in advance for not having included all possible references on the topic.

2. Ideal MHD: Numerical Techniques

Developing new computational methods for solving the ideal MHD equations has generated a lot of interest within the applied mathematics and computational physics communities. Quite naturally, because of their success in solving the Euler equations, grid-based methods with flux upwinding, also known as Godunov's method, were applied to the MHD equations in the framework of finite-volume discretisation. In parallel, Smoothed Particle Hydrodynamics (SPH) generated a lot of interest in the astrophysics community because of its strict Galilean invariance. Both methods, however, quickly ran into difficulties trying to maintain the divergence-free property of the MHD equations.

2.1. The Ideal MHD Equations

Before we review the various improvements of MHD numerical solvers over the past decades, we briefly recall the ideal MHD equations, shown here in conservative form. We have first the mass conservation equation,

$\begin{array}{ll}\frac{\partial \rho }{\partial t}+\nabla ·(\rho \text{v})=0,& (1)\end{array}$

where ρ is the mass density and v is the fluid velocity vector. We also have the momentum conservation equation

$\begin{array}{ll}\frac{\partial }{\partial t}(\rho \text{v})+\nabla ·\left(\rho \text{v}\otimes \text{v}+P_{\text{tot}}𝕀-\frac{1}{4\pi }\text{B}\otimes \text{B}\right)=\rho \text{g},& (2)\end{array}$

where B is the magnetic field vector and P_tot is the total pressure, defined as the sum of the thermal pressure and the magnetic pressure

$\begin{array}{ll}{P}_{\text{tot}}=P+\frac{1}{8\pi }{B}^2.& (3)\end{array}$

Note that we work here in cgs units, hence the presence of the 4π term in these equations. We have also included the gravitational acceleration vector g as a source term on the right-hand side of the momentum conservation equation. Finally, we have the total energy conservation equation

$\begin{array}{ll}\frac{\partial E_{\text{tot}}}{\partial t}+\nabla ·\left((E_{\text{tot}}+P_{\text{tot}})\text{v}-\frac{1}{4\pi }\text{B}(\text{B}·\text{v})\right)=\rho \text{g}·\text{v},& (4)\end{array}$

where the total fluid energy is the sum of the kinetic energy, the internal energy and magnetic energy

$\begin{array}{ll}{E}_{\text{tot}}=\frac{1}{2}\rho v^2+e+\frac{1}{8\pi }{B}^2.& (5)\end{array}$

In order to close the system, we need the fluid equation of state, usually given be the ideal gas equation of state

$\begin{array}{ll}P=(\gamma -1)e,& (6)\end{array}$

and the induction equation for the evolution of the magnetic field

$\begin{array}{ll}\frac{\partial \text{B}}{\partial t}=-\nabla \times \text{E}\text{ with }\text{E}=-\text{v}\times \text{B},& (7)\end{array}$

where we have introduced the electric field E (also known as the electromotive force, emf). We need to add to these equations the most important property of the magnetic field, namely

$\begin{array}{ll}\nabla ·\text{B}=0,& (8)\end{array}$

also known as the solenoidal constraint or the divergence-free condition. Since we consider now in this section ideal MHD conditions, we have no dissipation terms in the induction equation, as well as in the fluid equations.

2.2. Some Important Mathematical Properties

The fluid equations appear to be in conservative form, which naturally calls for a finite volume representation, which will ensure exact conservation of mass, momentum and total energy, by construction, owing to the divergence theorem

$\begin{array}{ll}{\displaystyle {\int }_V}\nabla ·\text{F}\text{d}V={\displaystyle {\int }_S}\text{F}·\text{n}\text{d}S,& (9)\end{array}$

where the vector F can be the flux of mass, momentum or energy. Designing a numerical scheme boils down to computing the flux through the surface of the volume elements. On the other hand, the induction equation naturally calls for a finite surface representation, owing to the Stoke's theorem

$\begin{array}{ll}\frac{\partial }{\partial t}{\displaystyle {\int }_S\text{B}}·\text{n}\text{d}S={\displaystyle {\int }_S\nabla }\times \text{E}·\text{n}\text{d}S={\displaystyle {\int }_L\text{E}}·\text{d}\text{l}.& (10)\end{array}$

Similarly, the divergence-free condition, written in integral form as

$\begin{array}{ll}{\displaystyle {\int }_V}\nabla ·\text{B}\text{d}V={\displaystyle {\int }_S}\text{B}·\text{n}\text{d}S,& (11)\end{array}$

also calls for defining the magnetic field as a surface averaged quantity. This naturalness argument, together with the fact that the divergence-free condition can be maintained to machine precision accuracy, has led to the design of the constrained transport scheme (Evans and Hawley, 1988). This very popular method for grid-based techniques comes however with a price: the flow variables are not all co-located at cell centers, like volume-weighted fluid variables, but also at face centers, for the magnetic field components, and at edge centers for the electric field components.

A very important mathematical property of the magnetic field that emerges from the divergence-free condition is that the normal component of the field should be continuous across cell faces. The x-component of the field, for example, can be discontinuous in the y and z-directions, but has to vary smoothly in the x-direction. This property also naturally holds in the constrained transport method. It can be written as

$\begin{array}{ll}{\displaystyle {\int }_S}\left[\text{B}·\text{n}\right]\text{d}S=0,& (12)\end{array}$

where the operator [A] = A⁺−A⁻ denotes the jump of a quantity A across the surface element.

2.3. Preserving the Divergence-Free Condition

Historically, one of the first ideal MHD, general purpose codes, is ZEUS-2D developed specifically for astrophysics by Stone et al. (1992b). It is a finite-difference code using constrained transport and artificial viscosity to handle shocks, The continuity and divergence-free conditions are therefore satisfied by construction, but since artificial viscosity is used to handle shocks, instabilities can occur in fast rarefaction waves (Falle, 2002).

The other popular strategy for grid-based methods is to maintain all MHD variables as volume-averaged quantities, allowing for discontinuities across all cell faces. In the late 90s, a series of papers presented finite-volume MHD codes with proper upwinding of numerical fluxes using Riemann solvers (Dai and Woodward, 1994; Ryu et al., 1995; Tóth, 1996; Balsara, 1998; Lee and Deane, 2009). These methods are very powerful, because they are strictly conservative and because they satisfy the so-called entropy condition, meaning the entropy can only increase, a key property to maintain the stability of the numerical solution. These finite-volume codes all considered the magnetic field as a volume-averaged, piecewise-constant, cell-centered quantity, which violates the continuity condition of the normal component of the field. The resulting schemes are therefore not necessarily divergence-free anymore. In order to maintain the divergence as small as possible, an additional step is required that modifies the magnetic field components and decreases or nullifies the divergence: this operation is called divergence cleaning. In a seminal paper, Tóth (2000) has compared various schemes and showed that they all passed with mixed success a variety of MHD test problems. We will now review these grid-based methods that are using divergence cleaning.

The first method we discuss here is the projection scheme, introduced by Brackbill and Barnes (1980). The idea is to measure the spurious divergence after the main update of the MHD variables, and solve for a Poisson equation defined as

$\begin{array}{ll}\Delta \varphi =\nabla ·{\text{B}}^{*}.& (13)\end{array}$

where ϕ is a scalar potential. The new magnetic field is then corrected as

$\begin{array}{ll}{\text{B}}^{n+1}={\text{B}}^{*}-\nabla \varphi ,& (14)\end{array}$

so that by construction the divergence is now zero. This method works very well in many applications. It suffers however from two main issues: first, it is quite expensive as it requires to solve for a Poisson equation. Another consequence is that the correction is non-local: very localized divergence errors can be instantaneously transported across the grid to enforce the solenoidality condition. Second, the correction process modifies the magnetic field, without modifying the total energy. As a consequence, the resulting temperature can be modified, and the entropy condition might be violated. An easy fix is to remove the magnetic energy before the correction step and add the new magnetic energy after the correction, but the resulting scheme is not conservative anymore.

The second popular method is the so-called 8-waves formulation or Powell's method (Powell et al., 1999). The idea is to write more general ideal MHD equations allowing for the presence of magnetic monopoles. This results in the following non-conservative form

$\begin{array}{ll}\frac{\partial }{\partial t}(\rho \text{v})+\nabla ·\left(\rho \text{v}\otimes \text{v}+P_{\text{tot}}𝕀-\frac{1}{4\pi }\text{B}\otimes \text{B}\right)=-(\nabla ·\text{B})\text{B}+\rho \text{g},& (15)\end{array}$

$\begin{array}{ll}\frac{\partial E_{\text{tot}}}{\partial t}+\nabla ·\left((E_{\text{tot}}+P_{\text{tot}})\text{v}-\frac{1}{4\pi }\text{B}(\text{B}·\text{v})\right)=-(\nabla ·\text{B})\text{B}·\text{v}+\rho \text{g}·\text{v},& (16)\end{array}$

$\begin{array}{ll}\frac{\partial \text{B}}{\partial t}-\nabla \times (\text{v}\times \text{B})=-(\nabla ·\text{B})\text{v}.& (17)\end{array}$

This method proved very useful and robust for many applications. It is still widely used in astrophysics nowadays. The success of the method originates from the property that the spurious ∇ · B is advected away with the flow, using the so-called eighth MHD wave, so that divergence errors do not accumulate. There are however two problems: First, in stagnation points¹, the divergence errors will accumulate because the flow trajectories are converging. Second, the scheme is not strictly conservative. Shock waves lead to solutions that do not converge to the correct conservation laws anymore.

A very elegant solution to the first problem was proposed in Dedner et al. (2002) using the so-called hyperbolic-parabolic divergence cleaning technique, also known as a Generalized Lagrange Multiplier (GLM) formulation of the ideal MHD equations, in short, Dedner's scheme. The idea is to add a ninth wave to the previous Powell's modified MHD equations, introducing the cleaning field ψ that satisfies

$\begin{array}{ll}\frac{\partial \psi }{\partial t}+\text{v}·\nabla \psi +c_h^2\nabla ·\text{B}=-\frac{\psi }{\tau },& (18)\end{array}$

and is used as a source term in the induction equation.

$\begin{array}{ll}\frac{\partial \text{B}}{\partial t}-\nabla \times (\text{v}\times \text{B})=-\nabla \psi .& (19)\end{array}$

In the first equation, c_h is the divergence wave speed and τ is the divergence damping time scale. The latter is chosen equal to (or larger than) the fast magnetosonic wave speed, while the former is usually equal to (or larger than) the fast magnetosonic cell crossing time. At first sight, this new nine waves scheme can be seen as a combination of advection and damping of divergence errors, thus a clever combination of the projection method and the Powell's scheme. As shown recently by Derigs et al. (2017), it is in fact much more than that: this new divergence field ψ allows to restore the entropy condition, as this field can be interpreted as a divergence cleaning energy which stores temporarily the magnetic energy lost during the damping step. It is also conservative in a general sense, but still violates the Rankine-Hugoniot shock relations locally.

In parallel, accurate and stable MHD solvers have been the topics of many studies in the SPH framework. Probably because of its strict Galilean invariance, early SPH MHD solvers were quite oscillatory (see e.g., Dolag et al., 1999). Truncation errors associated to the non-zero divergence cannot be damped by numerical diffusion as they are advected away, as in the case of grid-based codes. Many regularization techniques have been proposed to provide stable SPH solvers for the ideal MHD equations: using the vector potential (e.g., Kotarba et al., 2009) or using artificial diffusivity (Price and Monaghan, 2005; Dolag and Stasyszyn, 2009). Interestingly, the work of Price and Monaghan (2005) was based on exploring various divergence cleaning methods used in grid-based codes for SPH. The authors concluded that this works in most cases, but in difficult cases, such as supersonic turbulent flows, this can lead to spurious effects, such as the violation of the entropy condition or wrong shock relations. More worrysome, these authors noticed that divergence cleaning can lead to an increase of the local magnetic energy, a price to pay to redistribute the truncation errors in the divergence. This results in spurious dynamo effects. Recently, however, Tricco and Price (2012) revisited the Dedner's scheme for SPH and found a formulation that guarantees that divergence cleaning leads to a decrease of the magnetic energy and an increase of the entropy, in the spirit of Derigs et al. (2017) for grid-based codes.

In light of these rather complex developments of divergence cleaning methods, the constrained transport (CT) approach we have introduced already seems quite appealing. Note however that this approach requires to have well defined cell-interfaces, which is not necessarily the case for SPH or the recently developed moving-mesh AREPO code (Springel, 2010). The CT scheme, introduced for astrophysical fluid flows by Evans and Hawley (1988) and used in the ZEUS-2D code (Stone et al., 1992b), features the nice property that the divergence of the magnetic field, defined in integral form as the net magnetic flux across the 6 faces of each cell, will always remain constant, So if it is initially zero, it will remain equal to zero within numerical round-off errors. This however requires to define the electric field on the cell edges, and this is the main difficulty of this approach (see Figure 1). Indeed, in the Godunov methodology, fluxes, defined at cell faces, are computed as the solution to the Riemann problem of the two neighboring cells meeting at their common cell interface. For the electric field, defined at cell edges, 4 neighboring cells meet and in order to properly upwind the electric field, we need to solve a two-dimensional Riemann problem, a rather daunting task (see Teyssier et al., 2006, for a discussion). Several solutions have been found to design 2D MHD Riemann solvers (Londrillo and del Zanna, 2004; Balsara et al., 2014), and this has been the main characteristic of several simulations codes with proper upwinding of both fluxes and electric fields, using the CT scheme within the Godunov methodology (Londrillo and del Zanna, 2004; Fromang et al., 2006; Stone et al., 2008; Li et al., 2012; Mocz et al., 2014a).

FIGURE 1

Figure 1. Schematic showing the geometry of a Cartesian cell in the Constrained Transport approach. The finite volume cell is shown as a gray cube. It is labeled i, j, k. The magnetic field components are defined perpendicular to the faces of the cube. The are shown in red, only in the rightermost face in all three directions. The electric field components are defined on the cell edges. They are shown in blue and only in one face for sake of simplicity.

2.4. Minimizing Numerical Errors: Higher Order vs. Mesh Refinement

When performing self-gravitating magnetized fluid simulations, it is of primary importance to quantify and understand numerical errors. These errors are often called truncation errors because they arise as leading order terms in the Taylor expansion of the adopted discrete numerical method. Usually, this Taylor expansion leads to the so-called modified equation, which encodes how the original ideal MHD equations have been transformed into a new set of equations with additional terms coming from the truncation errors. For first-order schemes, these terms are identical to a diffusion process, and are therefore called numerical viscosity or numerical diffusivity. In SPH, quite often, these terms are added explicitly to the equations as artificial viscosity and artificial diffusivity, while for grid-based Godunov solvers, these terms are only implicitly present, through this Taylor expansion of the discrete scheme. In both cases, however, these diffusion processes control how shock heating and magnetic reconnection proceed in the solution. They play a fundamental role in preserving the entropy condition, and in regulating the flow close to the grid scale. Unfortunately, many complex MHD processes, such as the small scale dynamo (Brandenburg and Subramanian, 2005) or the magneto-rotational instability (Balbus and Hawley, 1991) depend crucially on the so-called Prandtl number (Fromang et al., 2007; Federrath et al., 2014a), which is the ratio of the real magnetic diffusivity to the viscosity. In most case, this ratio is always close to 1 if one uses the numerical Prandtl number (Fromang et al., 2007; Federrath et al., 2011b), while in nature, it can vary widely. It is therefore crucial to adopt in some cases explicit viscosity and diffusivity coefficients and model the flow including these non-ideal processes (Fromang et al., 2007; Federrath, 2016).

As explained in the next section, in order to model these non-ideal effects, it is crucial to control and minimize the numerical diffusion as much as possible. There are two possible strategies to achieve this: refine the grid or increase the order of accuracy of the method. The first approach leads to the so-called Adaptive Mesh Refinement method, a very popular and successful technique in the context of star and galaxy formation. Popular AMR codes are available to the star formation community, such as RAMSES (Teyssier, 2002; Fromang et al., 2006), ORION (Klein, 1999; Krumholz et al., 2007b), ENZO (Bryan et al., 2014), or FLASH (Fryxell et al., 2000). Other codes, that used to be only unigrid, now propose an adapted grid or an adaptive grid extension, such as PLUTO (Mignone et al., 2007) and ATHENA (Stone et al., 2008). In all these codes, cells are adaptively refined according to various refinement criteria. In the context of star formation, the most popular approach is to always resolve the local Jeans length with 4 cells or more, the so-called Truelove criterion (Truelove et al., 1997)

$\begin{array}{ll}\text{If}{\lambda }_{\text{J}}=\frac{c_s}{\sqrt{4\pi G\rho }}<4\Delta x_{\ell }\text{ then refine to level }\ell +1.& (20)\end{array}$

This corresponds to a level-dependent density threshold that triggers new refinements. SPH methods are Lagrangian in nature, so they cannot refine as much as grid-based codes. Usually, much more SPH particles are needed in the initial conditions to reach a certain target Jeans mass, corresponding to the maximum resolution level of the corresponding AMR simulation. Particle splitting is an interesting alternative to classical SPH but is still under development (Kitsionas and Whitworth, 2002; Chiaki and Yoshida, 2015), the difficulty being to handle sharp transitions in particle mass and its interaction with the smoothing kernel. Note that if similar resolution requirements are met, AMR and SPH methods largely agree on the quantitative predictions on how the collapse proceeds (Commerçon et al., 2008). For magnetized flows, the Jeans length-based refinement strategy has to be more conservative, of the order of 30 cells per Jeans length, in order to capture properly the magnetic field amplification in collapsing cores (Sur et al., 2010; Federrath et al., 2011c; Turk et al., 2012).

A difficulty arises when one uses AMR for Constrained Transport. In this case, it is mandatory to be able to interpolate the magnetic field from the coarser level to the finer level and still satisfy the divergence-free condition. For this, divergence preserving interpolation schemes have been developed (Balsara, 2001; Tóth and Roe, 2002) and play an important role in the viability of the CT approach in the context of AMR.

To reduce truncation errors, the other option is to use higher order schemes. The solution inside each cell is not piecewise constant as in the traditional first-order Godunov method, but it is reconstructed using high-order polynomials of degree p as

$\begin{array}{ll}\rho (x)={\displaystyle \sum _{i=0,p}}{\alpha }_i{\psi }_i(x).& (21)\end{array}$

The ψ_i(x) are usually an orthonormal basis of polynomials of degree at most p. The coefficients α_i are computed using two different philosophies. The first option, the WENO approach computes the coefficients α_i at each time step using neighboring cells. The higher the polynomial degree, the more neighbors must be included in the stencil of the method (Jiang and Wu, 1999; Boscheri and Dumbser, 2014). The second option, the Discontinuous Galerkin (DG) approach, considers that the α_i are new variables, defined for each MHD variables, and updated using a projection of the ideal MHD equation onto the polynomial basis (Li and Xu, 2012; Mocz et al., 2014b; Klingenberg et al., 2017; Guillet et al., 2018). In the WENO case, one needs to access neighboring but possibly distant cells to compute the α_i, while in the DG case, one needs to store permanently the α_i and solve for their evolution at each time step. Note that other high-order methods are also being developed that do not strictly correspond to neither WENO nor DG (Felker and Stone, 2018).

In the context of high-order methods, the divergence free condition is particularly challenging. One can either implement high-order version of one of the above mentioned divergence cleaning techniques (Derigs et al., 2016), or one can try to preserve the divergence-free constraint within the method. The second approach could be seen as a generalization of the CT scheme to higher order. The key ingredients are: 1- the use of a divergence-free polynomial basis for the magnetic field (Li and Shu, 2005; Guillet et al., 2018), 2- enforcing the continuity of the normal field across cells boundaries (Li and Xu, 2012). Although some very promising solutions have been found recently, this is still a very active field of research. Interestingly, the traditional CT scheme with face-centered magnetic field variables can be re-interpreted as a DG scheme where each magnetic field component is piecewise-linear and continuous in the normal direction and piecewise constant in the transverse direction, so that two cell-centered coefficients are required (instead of one) for each field component.

3. Non-ideal MHD: Numerical Techniques

3.1. Equations and Basic Concepts

Star formation takes place in molecular clouds, which are made of a mixture of dust and gas, implying dust and gas collisions, and both constituents are far from being fully ionized. In the previous section, we presented the work done in the ideal MHD framework, which do not appear to be well suited to the ionization state in collapsing cores, in particular at the onset of disk formation. Recent works have emphasized the imperfect coupling of the dust and gas mixture with the magnetic fields at the transition between the envelop and the disk in collapsing cores (see the review by Wurster and Li, 2018 for the work done in the context of protostellar disk formation) and a lot of effort has been devoted over the past 10 years to include the so called non-ideal effects: the ambipolar diffusion, the Ohmic diffusion, and the Hall effect. The ambipolar diffusion is the common name to describe the interaction between neutrals and charged particles. It can be seen as a friction term, it enables the neutral field to respond to the magnetic forces, via collisions with charged particles. The Ohmic diffusion results from the collision of electrons with the neutrals. Last, the Hall effect is due to the drift between the positively and negatively charged species. As shown in Marchand et al. (2016) and Tsukamoto et al. (2015a), all these three terms can be dominant over the others at different scales within the envelop of collapsing dense cores. For a classical dust size distribution (Mathis et al., 1977), Marchand et al. (2016) find that amipolar diffusion and the Hall effect dominate at densities < 10¹² cm⁻³ and that Ohmic diffusion is the stronger resistive effect at higher densities.

The exact scale and density at which these resistive effects become dominant over the other dynamical processes depend on the chemistry, the ionization intensity and the dust grain size distribution (Zhao et al., 2016, 2018a; Zhao et al., 2018b; Dzyurkevich et al., 2017). Hennebelle et al. (2016) have shown that ambipolar diffusion regulates the flow dynamics over the other dynamical processes (induction, rotation and free-fall) at scale of a few 10s AU, which sets the initial size of protostellar disks. In addition, Masson et al. (2016) have shown that scales of a few 10s AU exhibit magnetic Reynolds numbers less than unity (see also Tomida et al., 2013; Tomida et al., 2015; Vaytet et al., 2018 for studies including the Ohmic diffusion). Recently, Koga et al. (2019) performed a similar analysis considering the Hall effect only and found coherent results as for ambipolar diffusion.

Before we describe the numerical implementation for the three aforementioned resistive effects, let us recall the main equations and define the necessary quantities. Ideally, the system should account for the different behavior of the neutrals, as well as negatively and positively charged particles. Among them the different constituents that participate to the ionization balance, we have neutral (molecule, atoms, dust grains), ions (molecular and atomic), electrons, and charged dust grains. The latter can be positively or negatively charged, and can hold multiple charges (Draine and Sutin, 1987). Most current works follow a one-fluid approximation to describe the evolution of this complex mixture, where the Ohm's law accounts for the non-ideal effects with some assumptions. In the following, we focus on the work done in the single and two-fluid approximations. We review briefly the work done toward a multi-fluid treatment of the charged and neutral particles in section 3.6.

In the low ionization limit, each kind of charged particle is scarce, so that the pressure terms, the gravity, and the inertia term can be neglected in their momentum evolution equation in comparison to the Lorentz force and the frictional force exerted by the neutrals. The generalized Ohm's law reads

$\begin{array}{ll}{\text{E}}_{\text{NIMHD}}=\text{v}\times \text{B}-{\eta }_{\Omega }\text{J}-\frac{{\eta }_{\text{H}}}{\left|\right|\text{B}\left|\right|}\text{J}\times \text{B}+\frac{{\eta }_{\text{AD}}}{\left|\right|\text{B}|{|}^2}\text{J}\times \text{B}\times \text{B},& (22)\end{array}$

where J = ∇ × B is the current, η_Ω, η_H, and η_AD are the Ohmic, Hall, and ambipolar resistivities, respectively (in units of cm² s⁻¹), and v is the velocity of the neutrals. The notation ||B|| represents the norm of the magnetic field vector B. The electric field is then replaced in Equation (7). It is worth noticing that all the resistive terms lead to parabolic partial differential equations to discretise. We note that only the Ohmic term can lead to reconnection here. It can be rewritten as a Laplacian operator η_ΩΔJ, whereas the two others cannot. The ambipolar diffusion and the Hall effect can be assimilated to drift velocities of the magnetic fields with respect to the neutral speed v. The total energy equation has to be modified accordingly, to account for the heating resulting from the friction terms, i.e., the ambipolar diffusion and the Ohmic diffusion

$\begin{array}{ll}\frac{\partial E_{\text{tot}}}{\partial t}+\nabla ·\left(\left(E_{\text{tot}}+P_{\text{tot}}\right)\text{v}-\frac{1}{4\pi }\text{B}(\text{B}·\text{v})\right)=\rho \text{g}·\text{v}+{\eta }_{\Omega }{\Vert \text{J}\Vert }^2+{\eta }_{\text{AD}}\frac{{\Vert \text{J}\times \text{B}\Vert }^2}{\Vert \text{B}\Vert }.& (23)\end{array}$

We refer readers to Shu et al. (1987), Nakano et al. (2002), and Balbus (2009) for more details on the derivation of the generalized Ohm's law in the one-fluid approximation. In the context of star formation, where a huge dynamical range has to be considered, the estimate of the resistivities is challenging. It requires to know the ionization state of the gas and dust mixture. The coupled chemistry of gas and dust has to be considered, which involves a huge variety of physical processes: gas phase and gas-dust chemistry, cosmic ray ionization, UV radiation, (in-)elastic collisions, grain growth. This domain remains the subject of intense research works that we do not detail in this review. We refer readers to section 4.0.1 in the review by Wurster and Li (2018) for a scan of the works that tackle the resistivity calculations.

Ambipolar diffusion and Ohmic diffusion were the first to be considered in star formation applications. Indeed, these two terms do not introduce new MHD waves in the system (they do not change the spectral properties of the hyberpoblic system), so that they do not require heavy modifications of the MHD solver. The induction equation integration is most of the time splitted in two steps: first the ideal MHD and second the resistive terms. If an energy equation is considered, the heating terms due to resistive effect in Equation (23) are integrated as source terms in most numerical implementations.

3.2. Ohmic Diffusion

The Ohmic diffusion is perhaps the simplest resistive term to introduce because of its diffusive nature, similar to the introduction of artificial viscosity and artificial resistivity to prevent any numerical artifact due to numerical diffusion. Assuming a constant resistivity, the Ohmic term can be rewritten as a Laplacian (which preserves the solenoidal constraint), leading to the following stability condition

$\begin{array}{ll}\Delta t_{\Omega }=\frac{\Delta x^2}{{\eta }_{\Omega }}.& (24)\end{array}$

It corresponds to the stability condition required in numerical schemes integrating parabolic equations such as the heat equation with an explicit scheme.

Numerous implementations for the Ohmic diffusion have been developed in the past 10 years, in all different kind of codes. Masson et al. (2012) present a fully explicit implementation in the AMR code RAMSES, where they update the emf to account for the resistive terms and then use the CT scheme as for ideal MHD. Because of the restrictive stability condition for the resistive term which can be smaller that the MHD one, some authors considered schemes that enable to relax the time step constraint. For SPH, Bonafede et al. (2011) implemented Ohmic diffusion with a constant resistivity in GADGET following the method used for artificial dissipation (e.g., Price and Monaghan, 2005; Dolag and Stasyszyn, 2009). Tomida et al. (2013) implemented an explicit scheme in a nested-grid code following the implementation of Matsumoto (2011) in the AMR code SFUMATO. They adopt the Super-Time-Stepping (STS) algorithm proposed by Alexiades et al. (1996) to accelerate the integration of the diffusion term and to relax the stringent explicit CFL condition. The STS consists in a discretisation as a classical explicit scheme, which does not require to satisfy the stability in every step but after a series of integration cycles which are made of N_STS sub-timesteps Δt_STS which are based on Chebyshev polynomials to guarantee the stability conditions over the super-timestep $\Delta t_{\text{tot}}=\sum _1^{N_{\text{STS}}}\Delta t_{\text{STS}}$. We also note the work of Mignone et al. (2007) and O'Sullivan and Downes (2006) where Ohmic diffusion is also integrated using the STS technique. Similarly, Tsukamoto et al. (2013) proposes a SPH discretisation of Ohmic dissipation using STS. Wurster et al. (2016) and Price et al. (2017) (in the PHANTOM code) use similar implementations on different SPH codes. Last, Marinacci et al. (2018) propose two implementations in AREPO: one using Powell divergence cleaning and another using CT. Both implementations can be integrated using an explicit or implicit scheme. Almost all the implementations mentioned so far for the Ohm diffusion are combined with ambipolar diffusion (see below).

Apart from resistive MHD, the STS technique is becoming increasingly popular in the computational astrophysics community to accelerate parabolic timestepping advancements in the context of radiation hydrodynamics with the FLD (Commerçon et al., 2011) and anisotropic diffusion (Meyer et al., 2012; Vaidya et al., 2017).

3.3. Ambipolar Diffusion

Ambipolar diffusion is not as straightforward to include since the associated term in the induction equation can not be rewritten as a Laplacian (Brandenburg and Zweibel, 1994). For the diffuse ISM, the drift between ions and neutrals is the effect of the ambipolar diffusion, that can be written as a friction force

$\begin{array}{ll}{\text{F}}_{\text{in}}={\rho }_{\text{i}}{\rho }_{\text{n}}{\gamma }_{\text{AD}}(\text{v}-{\text{v}}_{\text{i}}),& (25)\end{array}$

where ρ_i (resp. v_i) and ρ_n (v) are the mass density (resp. velocity) of ions and neutrals, and γ_AD is the collisional coupling constant, with ${\gamma }_{\text{AD}}~3\times 10^{13}$ cm³ s⁻¹ g⁻¹ in ISM mixture (Draine et al., 1983). In the one-fluid and strong coupling limit, the inertia of ions is neglected, as well as the pressure and gravitational forces (the ionization degree is very low). As a consequence, the Lorentz force equals the drag force so that the ion drift velocity is

$\begin{array}{ll}{\text{v}}_{\text{d}}\equiv {\text{v}}_{\text{i}}-\text{v}=\frac{1}{{\rho }_{\text{i}}{\rho }_{\text{n}}{\gamma }_{\text{AD}}}(\nabla \times \text{B})\times \text{B}.& (26)\end{array}$

The induction equation thus reads

$\begin{array}{ll}\frac{\partial \text{B}}{\partial t}=\nabla \times \left[{\text{v}}_{\text{i}}\times \text{B}\right]=\nabla \times \left[\text{v}\times \text{B}+\frac{1}{{\rho }_{\text{i}}{\rho }_{\text{n}}{\gamma }_{\text{AD}}}(\nabla \times \text{B})\times \text{B}\times \text{B}\right].& (27)\end{array}$

From Equation (22), we have ${\eta }_{\text{AD}}=\left|\right|B|{|}^2/{\rho }_i{\rho }_n{\gamma }_{\text{AD}}$. In the dense ISM, the heavy charged particles are ions and charged dust grains, and the expression of the resistivity is more complex (Marchand et al., 2016). In the two-fluid approximation, the ion inertia is not neglected and the system of equation accounts for both the ions and the neutrals momenta equation.

The first implementation of ambipolar diffusion in a grid code was presented in Black and Scott (1982) using a two-fluid approximation. They integrated the induction equation using the ion velocity and treated the collision between neutrals and ions as a friction force in the momenta equations. Mac Low et al. (1995) present a first 3D implementation of one-fluid ambipolar diffusion in the ZEUS code using a fully explicit scheme as well as CT. They propose to subcycle the magnetic fields evolution in the case when ambipolar diffusion timestep is much shorter than the dynamical one. The first SPH implementation can be found in Hosking and Whitworth (2004), where they used a two-fluid approximation, but their scheme did not preserve the solenoidal constraint. Duffin and Pudritz (2008) present a one-fluid implementation in the FLASH code using Powell's method (Powell et al., 1999). Masson et al. (2012) also implemented a one-fluid scheme for ambipolar diffusion in RAMSES with the same philosophy as for the Ohmic diffusion. A similar implementation using CT, operator splitting, and STS in the ATHENA code can be found in Bai and Stone (2011) and Chen and Ostriker (2014). Last, Wurster et al. (2014) propose a one-fluid implementation for SPH codes, which they implemented with success in PHANTOM.

Other one-fluid implementations in the strong coupling approximation can be found in Padoan et al. (2000) in a 3D explicit, finite difference grid code, Choi et al. (2009) in the MHD TVD code using a flux-interpolated CT scheme (Ryu et al., 1998; Balsara and Spicer, 1999), Christie et al. (2017) in the ENZO's Dedner MHD solver. Oishi and Mac Low (2006) and Li et al. (2006) present two simultaneous but independent implementations in the ZEUS code of a two-fluid solver in the heavy ion approximation.

For star formation applications, most of these methods are usually second-order accurate and at best second-order accurate in time (with a few exceptions though e.g., Meyer et al., 2014). Clearly, further work needs to be done to improve the accuracy of resistive MHD solvers. Alongside, multi-fluid approaches are also needed to account for more detailed physics (see section 3.6).

3.4. Hall Effect

The last resistive effect that has been implemented in star formation models is becoming increasingly popular. The Hall effect originates from the drift between ions and electrons, or, more generally the positively and negatively charged particles. The derivation of the generalized Ohm law (22) is not as straightforward as in the case of ambipolar diffusion only. We shall consider the momentum equations of ions and electrons, accounting for the Lorentz force as well as collisions between ions, electrons, and neutrals. If we restrict the resistive effects to the Hall effect only, the induction equation reads

$\begin{array}{ll}\frac{\partial \text{B}}{\partial t}=\nabla \times \left[\text{v}\times \text{B}-\frac{{\eta }_{\text{H}}}{\left|\right|\text{B}||}\text{J}\times \text{B}\right]=\nabla \times \left[(\text{v}+{\text{v}}_{\text{H}})\times \text{B}\right],& (28)\end{array}$

where we introduce the “Hall speed” ${\text{v}}_{\text{H}}=\frac{{\eta }_{\text{H}}}{\left|\right|\text{B}\left|\right|}\text{J}$. Contrary to the two previous resistive terms, the Hall effect introduces new waves in the MHD systems: the whistler waves. The dispersion relation for whistler waves reads (Balbus and Terquem, 2001)

$\begin{array}{ll}\omega =\pm \frac{{\eta }_{\text{H}}{k}^2}{2}+\sqrt{{\left(\frac{{\eta }_{\text{H}}{k}^2}{2}\right)}^2+k^2{c}_{\text{A}}^2}.& (29)\end{array}$

The whistler waves are dispersive, the phase speed c_w = ω/k depends on the wavenumber k and tends to infinity as the wavenumber tends to infinity. Since any discrete formulation cannot extend to infinity, a numerical solver cannot follow the whistler waves with very high frequencies. In practice, the whistler wave speed is chosen to be the one at the grid scale

$\begin{array}{ll}{c}_{\text{w}}=\left|\frac{{\eta }_{\text{H}}\pi }{2\Delta x}\right|+\sqrt{{\left(\frac{{\eta }_{\text{H}}\pi }{2\Delta x}\right)}^2+c_{\text{A}}^2}.& (30)\end{array}$

The time step is then constrained as

$\begin{array}{ll}\Delta t\le \frac{\Delta x}{c_{\text{w}}}\approx \frac{\Delta x^2}{{\eta }_{\text{H}}}.& (31)\end{array}$

The first way to implement the Hall effect is to use an operator split method, similarly to the Ohm and ambipolar diffusion. A pioneer implementation is presented in Sano and Stone (2002) in a second-order Godunov finite-difference code (Miyama et al., 1999), in the local shearing box framework. They use the CT method to update the induction equation with emf computed with the Hall term. The method has been implemented in 3D in ZEUS in Krasnopolsky et al. (2011). Falle (2003) presents an implicit scheme for multi-fluid MHD which alleviates the stringent time step conditions, where the time step goes to zero as the strength of the Hall effect increases. O'Sullivan and Downes (2006) designed a Hall Diffusion Scheme (HDS) in an explicit 3D code for multi-fluid MHD. They split the Hall diffusion in two parts, where the first part uses the STS time integration up to a critical value of the Hall resistivity (i.e., to satisfy the stability condition), and the second part, the HDS, diffuses the excess Hall resistivity with a different discretisation scheme. The balance between the ambipolar diffusion and the Hall resistivity determines the excess Hall resistivity. The HDS integrates the induction equation for the excess Hall using a dimensional splitting which is explicit for the first dimension, i.e., the x component $B_x^{n+1}$ of the magnetic fields is given by $B_y^n$ and $B_z^n$ at time n, explicit-implicit for the second dimension, $B_y^{n+1}$ is given by $B_x^{n+1}$ and $B_z^n$, and finally implicit for third dimension, $B_z^{n+1}$ is given by $B_x^{n+1}$ and $B_y^{n+1}$. Both STS and HDS can be subcycled to reach the time step limit given by the hyperbolic system. The HDS scheme has been also implemented in ATHENA by Bai (2014) using CT.

An alternative for finite volume methods is to include the Hall effect into the conservative integration scheme. Tóth et al. (2008) proposed a first implementation in the AMR code BATSRUS (Tóth et al., 2006). It uses block-adaptive grids with both explicit and implicit time discretisation, and various methods to control numerical errors in ∇ · B (8-waves scheme with diffusion, projection scheme). They achieve spatial second-order convergence by using symmetric slope limiters (like monotonized central) instead of asymmetric limiters (like minmod). In addition, they show that the first-order Lax-Friedrich Riemann solver is inconsistent for the Hall MHD equations. This comes from the fact that the whistler wave speed is proportional to the grid spacing, which makes a one order loss in space in the truncation errors of the scheme. The same applies for MUSCL schemes and asymmetric slope limiters which are only first-order accurate. Lesur et al. (2014) built upon Tóth et al. (2008) a Hall MHD solver in PLUTO for protoplanetary disc evolution studies. They implement a HLL solver to estimate the flux and then update the magnetic fields using CT. When written in a conservative form, the Hall induction equation reads

$\begin{array}{ll}\frac{\partial \text{B}}{\partial t}-\nabla ·\left[\text{v}\text{B}-\text{B}\text{v}-\frac{{\eta }_{\text{H}}}{\left|\right|\text{B}\left|\right|}(\text{J}\text{B}-\text{B}\text{J})\right]=0.& (32)\end{array}$

The current J = ∇ × B involves spatial derivatives of the magnetic fields, so that the flux depends on the conserved quantities, and on its derivative. This leads to a ill-defined Riemann problem. Lesur et al. (2014) solved this issue by assuming that J is an external parameter, which is computed at the cell interfaces using the predicted states. The left and right fluxes take this unique value as an input, with the predicted B. The Godunov fluxes are then computed using a whistler-modified HLL solver, which accounts for the whistler waves speed (30) in the characteristics speeds.

Last, Marchand et al. (2018) recently extended the Lesur et al. (2014) implementation to the CT scheme of RAMSES, which uses 2D Riemann solver to integrate the induction equation. They designed a 2D whistler-modified HLL solver to estimate the emf on cell border, assuming a uniform current in the 2D Riemann solver, again computed from the predicted states.

In SPH codes, the Hall effect has been implemented and used with success for star formation applications in Wurster et al. (2016), Tsukamoto et al. (2017), and Price et al. (2017). The standard method is based on an operator splitting and either a full explicit or a STS scheme based on the implementation for ambipolar diffusion (Wurster et al., 2014).

To date, only (Krasnopolsky et al., 2011; Wurster et al., 2016; Tsukamoto et al., 2017; Marchand et al., 2018) applied their methods for protostellar collapse. While the SPH methods seem to lead to a relatively accurate conservation of the angular momentum, the grid-based methods using CT are suffering from severe non-conservation of the angular momentum. The origin of this non-conservation is unclear. Krasnopolsky et al. (2011) invoked Alfvén waves going out of the simulation volume, while Marchand et al. (2018) did not find evidence of this. Instead, they propose that the non-conservation takes place after the accretion shock formation. Shocks indeed generate strong gradients and thus large Hall velocities. Further investigation is clearly needed to cope with this fundamental issue.

3.5. Full Resistive Implementation

Currently, a handful of 3D codes benefit of a full implementation of the three resistive terms: PLUTO (Lesur et al., 2014), ATHENA (Bai, 2014), RAMSES (Masson et al., 2012; Marchand et al., 2018), the SPH code by Tsukamoto et al. (2017), PHANTOM (Price et al., 2017), ZeusTW (Li et al., 2011), and GIZMO (Hopkins, 2017). A very recent full implementation for solar physics can be found in González-Morales et al. (2018) in the MANCHA3D grid code. It combines STS for ambipolar and ohmic diffusion, and HDS for the Hall effect.

Currently, the ATHENA, PLUTO, RAMSES, and ZeusTW implementations rely on CT algorithms for the resistive MHD equations integration. All other works use divergence cleaning algorithms. Given the variety of the methods developed in the literature, a quantitative comparison (in standard tests as well as in star formation applications) of all these implementations would be welcome in the coming years.

Last but not least, it is not clear which non-ideal process dominates over the others. The physical conditions in star-forming regions are so wide that there is room for regions where all three effects dominate. In addition, every effect needs to be tested and quantified before any conclusions on the relative importance of each of these effects (see the review by Wurster and Li, 2018).

3.6. Multifluid Approach

Multifluid approaches are designed to describe multiphase or multifluid flows. In the context of star formation, different fluids or phases are at stack: neutrals, ions, electrons, and dust grains. A few attempts have been made to describe all or a part of these components using multi-fluid approach. In this section, we consider only the approaches which do account for N momentum equations, where N is the number of phases or fluids. For instance, a popular two-fluid framework considers ions and neutrals and describes the friction between these two fluids, i.e., the ambipolar diffusion. Toth (1994) presented a first 2D implementation of two-fluid ions and neutrals dynamics, which is based on a flux-corrected transport finite difference method. Inoue and Inutsuka (2008) propose an unconditionally stable numerical method to solve the coupling between two fluids ions-neutrals (frictional forces/heatings, ionization, and recombination, Draine, 1986) for an ISM gas subject to the thermal instability. They split time integration as follows: 1/ the ideal hydrodynamical part for neutrals, 2/ the ideal MHD part for ions. The first part is integrated using a second-order Godunov scheme. The ideal MHD part is solved in two steps, similarly to what is done for one-fluid ideal MHD solvers, i.e., Stone et al. (1992b): the magnetic pressure terms are solved using a second-order Godunov method, and the magnetic tension term and the induction equation are solved using the method of characteristics for Alfvèn waves with the CT algorithm. Tilley and Balsara (2008) present a semi-implicit method for two-fluid ion-neutral ambipolar drift in the RIEMANN code (Balsara, 2001, 2004). Their scheme is second-order accurate in space and time and uses a mixed implicit-explicit formulation to deal with strong friction regimes. Pinto et al. (2008) consider a system composed of three fluids: positively-charged particles, negatively-charged particles, and neutrals. The charged particles include ions, electrons, as well as charged dust grains, depending on the main charge carrier. The self-consistent set of equations they designed remains to be implemented in 3D codes.

Last but not least, a complete multi-fluid approach should include dust grains dynamics on top of the neutrals, ions, and electrons. Falle (2003) presented a scheme for multi-fluid hydrodynamics in the limit of small mass densities of the charged particles, i.e., where the inertia of the charged particles are neglected. This is a similar to the heavy ion approximation but for neutrals, ions, electrons, and dust grains of various sizes. The HD equations are first integrated for the neutrals using a second-order Godunov scheme, and then the induction equation and the charged species velocities are updated using an implicit scheme. This scheme captures the three non-ideal effects and can deal with regime where the Hall effect dominates. We note that similar approaches were introduced in Ciolek and Roberge (2002) and O'Sullivan and Downes (2006).

The most complete work on multi-fluid designed for star formation was presented in Kunz and Mouschovias (2009). They considered a six-fluid set of equations (neutrals, electrons, molecular and atomic ions, positively charged, negatively charged, and neutral grains) and implemented it in a 2D version of the radiation-MHD ZEUS-MP code (Hayes et al., 2006). They modified the MHD solver to account for non-ideal effects (Ohmic dissipation and ambipolar diffusion). The abundances of the 6 species are calculated using a reduced chemical-equilibrium network. They applied their methods to protostellar collapse in Kunz and Mouschovias (2010). Unfortunately, the CPU cost of such a complete set of physics is prohibitive for 3D models, and no work has taken over since this first application.

3.7. What Is Next? Dust and Gas Mixture Dynamics

Currently, numerous studies report on the effect of non-ideal MHD in the star formation process. The wide majority of these works use the single-fluid approximation, assuming furthermore that the dust grains are perfectly coupled to the gas via collisions. Nevertheless, dust grains are observed to cover a wide range of sizes, from nanometer to micrometer in the ISM (Mathis et al., 1977) and even millimeter in protostellar disks environments (Chiang et al., 2012). Generally, it assumes that dust grains of different sizes coexist and follow a power-law dust distribution. Grains react to changes in the gas velocity via a drag force. For a given dust grain of mass m_g, the force is

$\begin{array}{ll}{\text{F}}_{\text{drag}}=-\frac{m_{\text{g}}}{t_{\text{s}}}({\text{v}}_{\text{g}}-\text{v}),& (33)\end{array}$

where v_g is the dust grain velocity, v the gas velocity, and t_s the stopping time, i.e., the characteristic decay timescale for the dust grain velocity relative to the gas. In star forming regions, the mean free path of the gas molecules is larger than the dust particles radius and the stopping time depends linearly on the dust grain size t_s ~ s (Epstein, 1924). The Stokes number St ≡ t_s/t_dyn characterizes the coupling of the dust grain relative to the gas, where well coupled dust grains follow St < 1. t_dyn represents a characteristic dynamical time, e.g., the eddy turbulent crossing time for molecular clouds or the orbital time for prostostellar discs. The Stokes number thus depends on the dust grain size so that part of the dust grains of a given dust size distribution will be coupled to the gas and the other part will be decoupled.

Several implementations of dust and gas mixtures have been developed in the literature to deal with the various Stokes regime. For St > 1, the bi-fluid formalism seems to be the most adapted (e.g., Bai and Stone, 2010; Laibe and Price, 2012; Meheut et al., 2012; Lorén-Aguilar and Bate, 2014; Booth et al., 2015). However, bi-fluid algorithms encounter sever limitations in the strong drag regime (small stopping time) (Laibe and Price, 2012), where prohibitive spatial and time discretisations may be required. In recent years, Laibe and Price (2014a) proposed a monofluid formalism for the dust and gas mixture which is well adapted to regimes with St ≤ 1. This monofluid formalism has been implemented with success in SPH and grid-based codes in the past 5 years (Laibe and Price, 2014b; Lebreuilly et al., 2019). Last, Price and Laibe (2015); Lin and Youdin (2017) present methods for simulating the dynamics of small grains with only one additional equation on the dust concentration on top of the Euler equations. These last methods have been recently implemented with success in PHANTOM (Price et al., 2017), as well as in PLUTO (Chen and Lin, 2018) and in RAMSES (Lebreuilly et al., 2019).

Recent works have emphasized the possible decoupling of micrometer dust grains in molecular clouds (Hopkins and Lee, 2016; Tricco et al., 2017) as well as of millimeter grains in collapsing dense cores (Bate and Lorén-Aguilar, 2017, see Figure 2). This decoupling is at the origin of a dynamical sorting which changes the dust size distribution. Other mechanisms may also affect the size distribution: coagulation, fragmentation, sublimation, etc. Further work should investigate the relative importance of each of these processes on the dust size distribution. In addition, the shape of the dust size distribution plays a critical role in the calculations of the resistivity as well as of the opacity. We can anticipate that coupled dust and gas dynamics will consider the evolution of the dust size distribution, coupled to non-ideal MHD and RHD algorithms. Last but not least, all current works on this field only account for the hydrodynamical (the pressure gradient) and gravitational forces. Adding the dynamics of charged dust grains is the natural next step.

FIGURE 2

Figure 2. SPH calculations of dust dynamics during the collapse of a 1 M_⊙ dense core. The left column shows the gas column density (edge-on view) and from left to right, the three other columns show the dust column density for dust grain sizes of 10 μm, 100 μm, and 1 mm. Large dust grains decouple from the gas and settle in the mid-plane faster into a large dusty disc with larger dust-to-gas ratio than in the envelop. Figure adapted from Bate and Lorén-Aguilar (2017) with permission from the authors.

4. Radiative Transfer: Numerical Techniques

4.1. A Poor's Man Approach: The Polytropic Equation of State

Radiative transfer plays a fundamental role in star formation and cannot be ignored, even for a review on magnetic fields. The importance of radiation fields is 2-fold. First, when the molecular cores collapse and reach a high enough density, the dust opacity becomes large enough to absorb the cooling radiation and the gas becomes optically thick. The traditional value for the dust opacity is

$\begin{array}{ll}{\kappa }_{\text{dust}}\simeq 0.1{\text{cm}}^2/\text{g at }T\simeq 10\text{K}.& (34)\end{array}$

Requiring that the Jeans length of the gas becomes optically thick to infrared radiation leads to

$\begin{array}{ll}\rho {\kappa }_{\text{dust}}{\lambda }_{\text{J}}>1\text{ or }\rho >{\rho }_{\text{crit}}=\frac{G}{\pi {\kappa }_{\text{dust}}^2{c}_{\text{s}}^2}\simeq 10^{-14}{\text{g}/\text{cm}}^3.& (35)\end{array}$

Many MHD simulations use this critical density to define a polytropic equation of state of the form

$\begin{array}{ll}T(\rho )=10\text{K}\left(1+{\left(\frac{\rho }{{\rho }_{\text{crit}}}\right)}^{2/3}\right),& (36)\end{array}$

to capture the transition from an optically thin, isothermal gas at low density to an optically thick, adiabatic gas at high density. Although this simplified approach is very useful for many already expensive simulations of collapse of turbulent clouds, this does not model properly the physics of radiative transfer (e.g., Commerçon et al., 2010). More importantly, once protostars have formed, higher energy radiation will emerge from the accretion shocks on the surface of the protostars, or from the main sequence stars after they ignited nuclear reactions. This high energy radiation in the optical or UV bands will then interact with the parent cloud, providing radiation feedback and influence the formation of the next generation of stars (Krumholz et al., 2007a; Myers et al., 2013). In conclusion, radiation plays a central role during the formation and the evolution of the first Larson core (see section 5), but later on acts as a self-regulating mechanism for star formation.

4.2. Ray Tracing and Long Characteristic Methods

The radiative transfer equation written in full generality using the radiation specific intensity I_ν(x, n, t)

$\begin{array}{ll}\frac{1}{c}\frac{\partial }{\partial t}{I}_{\nu }+\text{n}·\nabla I_{\nu }=j_{\nu }-{\alpha }_{\nu }{I}_{\nu }.& (37)\end{array}$

where j_ν is the emission coefficient and α_ν is the absorption coefficient. Note that the radiative transfer equation is written here in the laboratory frame. The absorption and emission coefficients are however defined in the comoving frame. Properly introducing the relativistic corrections to the radiative transfer equation and its various moments is a very important aspect of the problem we will only briefly touch upon in this review.

The most natural numerical technique to solve the radiative transfer equation is ray-tracing. The idea is to shoot a discrete set of light rays from a point source, solving the previous equation along the path of the light ray using the curvilinear coordinate s defined as ds = cdt. The equation of radiative transfer can be written as a Lagrangian time derivative following the trajectory of the photons as

$\begin{array}{ll}\frac{\text{d}}{\text{d}s}{I}_{\nu }=j_{\nu }-{\alpha }_{\nu }{I}_{\nu }.& (38)\end{array}$

The light ray is discretised along its path, matching the underlying gas distribution, whether it is a grid or a set of SPH particles. The difficulty is then to choose the appropriate number of rays, so that the angular distribution of the radiation field is properly sampled, but also so that enough light rays intersect the gas elements. In Abel and Wandelt (2002), an adaptive method was proposed, that splits or merges light rays adaptively, using the Healpix tessellation to sample the sphere. This method was first developed for the ENZO code (Wise and Abel, 2011) and used with success to describe photo-ionization regions at the Epoch of Reionization (Wise et al., 2012), as well as in turbulent molecular clouds (Shima et al., 2017). A similar scheme based on the Monte-Carlo approach has been developed for SPH in Altay et al. (2008) and Pawlik and Schaye (2008).

A similar approach is provided by the long characteristics method. This method can describe light rays emanating from a single point source like the ray-tracing method, but also diffuse radiation by adopting for each direction a set of parallel rays whose spacing matches the local grid resolution. This method has been implemented in the FLASH AMR code by Rijkhorst et al. (2006), with however an important addition we will discuss later. Diffuse radiation using long characteristics was included later by Buntemeyer et al. (2016) in the FLASH code, with an emphasis of describing both the optically thin and optically thick regimes. These methods exploit the particular property of the AMR implementations of FLASH and ENZO, namely large rectangular patches of various sizes and resolution. The long characteristic method has been developed for the graded octree AMR code RAMSES only recently by Frostholm et al. (2018).

4.3. Monte Carlo and Short Characteristic Methods

Although long characteristics offer the advantage of solving the exact attenuation along light rays, an especially accurate approach for point sources, they are tricky to parallelise because light rays propagate over multiple MPI domains. An alternative technique is the short characteristics method, for which light rays as short as the hydrodynamics cell size are considered. The radiation intensity from neighboring light rays is interpolated between neighboring cells, avoiding the use of long characteristics. Rijkhorst et al. (2006) use an hybrid long- and short-characteristics method in the FLASH code to overcome the difficulty to parallelise the long characteristics method. The short characteristic method also offers better scaling properties. The short characteristics method has also been used in the ATHENA code to solve the stationary radiative transfer equation by Davis et al. (2012) and in the C2-ray code for ionization problems (Mellema et al., 2006).

A third method which shares with the two previous ones the property of being accurate and able to capture complex radiation geometry is the Monte Carlo method. The Monte Carlo approach samples the radiation field with photon packets carrying information such as angle of propagation and frequency. This is as close as one can be from true photons propagating in the fluid and allows for complex scattering processes to be included in the models (Ercolano et al., 2003; Altay et al., 2008; Harries, 2011; Roth and Kasen, 2015). The method can be however quite expensive, as large particle numbers are required to avoid that the Poisson noise dominates the computed signal. In the diffusion limit, the mean free path becomes much smaller than the length scale of the system and the number of interactions becomes prohibitively large: It scales as the square of the optical depth (Nayakshin et al., 2009). It has also been mostly developed for stationary problems but many codes able to deal with non-stationary problems are now emerging (Nayakshin et al., 2009; Harries, 2011; Lomax and Whitworth, 2016).

4.4. Two Moments Methods

The common aspect in these three methods is that they can be quite expensive. This is probably not true for the short characteristics method, although the number of angular domains adopted plays a role in the final accuracy of the solution, and affects the cost. An alternative approach, called moment-based radiative transfer, can potentially solve this cost problem by integrating the radiative transfer equation over the angles. We obtain the radiation energy and the radiation momentum equation as

$\begin{array}{ll}\frac{\partial }{\partial t}{E}_{\nu }+\nabla ·{\text{F}}_{\nu }=4\pi j_{\nu }-{\alpha }_{\nu }c{E}_{\nu },& (39)\end{array}$

$\begin{array}{ll}\frac{1}{c}\frac{\partial }{\partial t}{\text{F}}_{\nu }+c\nabla ·{\mathbb{P}}_{\nu }=-{\alpha }_{\nu }{\text{F}}_{\nu },& (40)\end{array}$

where E_ν is the radiation energy, F_ν is the radiation flux and ℙ = 𝔻E_ν is the radiation pressure tensor and 𝔻 is the Eddington tensor. For slowly varying spectral absorption and emission coefficients, for which one has

$\begin{array}{ll}\nu \frac{\partial {\alpha }_{\nu }}{\partial \nu }\ll {\alpha }_{\nu }\text{ and }\nu \frac{\partial j_{\nu }}{\partial \nu }\ll j_{\nu },& (41)\end{array}$

one can write relativistic correction to first order in (v/c) as

$\begin{array}{ll}{\alpha }_{\nu }(\text{n})={\alpha }_{\nu }^0-(\text{n}·\text{v}/c){\alpha }_{\nu }^0,& (42)\end{array}$

$\begin{array}{ll}{j}_{\nu }(\text{n})=j_{\nu }^0+2(\text{n}·\text{v}/c)j_{\nu }^0,& (43)\end{array}$

where superscript 0 refers to quantities measured in the comoving frame (see Mihalas and Klein, 1982; Mihalas and Mihalas, 1984; Krumholz et al., 2007b, for a detailed discussion). Injecting this in the radiative transfer equation and taking the moments over the angles leads to new terms accounting for the relativistic corrections

$\begin{array}{ll}\frac{\partial }{\partial t}{E}_{\nu }+\nabla ·{\text{F}}_{\nu }=4\pi j_{\nu }^0-{\alpha }_{\nu }^{0}c{E}_{\nu }+{\alpha }_{\nu }^0{\text{F}}_{\nu }·\frac{\text{v}}{c},& (44)\end{array}$

$\begin{array}{ll}\frac{1}{c}\frac{\partial }{\partial t}{\text{F}}_{\nu }+c\nabla ·{\mathbb{P}}_{\nu }=-{\alpha }_{\nu }^0{\text{F}}_{\nu }+{\alpha }_{\nu }^0\mathbb{P}\text{v}+4\pi j_{\nu }^0\frac{\text{v}}{c}.& (45)\end{array}$

Since E_ν and F_ν are still expressed in the laboratory frame but ${\alpha }_{\nu }^0$ and $j_{\nu }^0$ are expressed in the comoving frame, these equations are often called mixed frame radiation moments equations. We can also Lorentz-transform the radiation moments from the laboratory frame to the comoving frame, which gives, to first order in (v/c)

$\begin{array}{ll}{E}_{\nu }^0=E_{\nu }-\frac{2}{c}{\text{F}}_{\nu }·\frac{\text{v}}{c},& (46)\end{array}$

$\begin{array}{ll}{\text{F}}_{\nu }^0={\text{F}}_{\nu }-E_{\nu }\text{v}-\mathbb{P}\text{v}.& (47)\end{array}$

Injecting these relations only in the right-hand side of the moments equations and neglecting all (v/c) terms gives the simpler form

$\begin{array}{ll}\frac{\partial }{\partial t}{E}_{\nu }+\nabla ·{\text{F}}_{\nu }=\left(4\pi j_{\nu }^0-{\alpha }_{\nu }^{0}c{E}_{\nu }^0\right)-{\alpha }_{\nu }^0{\text{F}}_{\nu }^0·\frac{\text{v}}{c},& (48)\end{array}$

$\begin{array}{ll}\frac{1}{c}\frac{\partial }{\partial t}{\text{F}}_{\nu }+c\nabla ·{\mathbb{P}}_{\nu }=-{\alpha }_{\nu }^0{\text{F}}_{\nu }^0+\left(4\pi j_{\nu }^0-{\alpha }_{\nu }^{0}c{E}_{\nu }^0\right)\frac{\text{v}}{c},& (49)\end{array}$

where all right-hand side terms are now expressed in the comoving frame. In practice, after the radiation transport step, encoded in the left-hand side, one then converts the laboratory frame radiation variables into the comoving frame and solve the thermo-chemistry encoded in the right-hand side in the comoving frame. Local thermodynamical equilibrium conditions can then be reached in the comoving frame, eventually reaching the diffusion limit (see below). At the end of the thermo-chemistry step, one finally converts back the radiation variables into the laboratory frame before entering the next transport step (see Rosdahl and Teyssier, 2015, for a detailed implementation).

These equations are quite handy because of their simplicity: We have no angular dependence anymore. There is a catch however, since everything depends on the Eddington tensor, for which we have no equation at this order of the moment's hierarchy. At this point, two different strategies have been explored in the star formation literature: (1) compute an accurate Eddington tensor by solving exactly the stationary radiative transfer equation using one of the methods presented above (short or long characteristics, ray-tracing or Monte Carlo), (2) compute an approximate form of the Eddington tensor based on a particular closure of the moment's hierarchy.

The first approach is often referred to a the Variable Eddington Tensor (VET) method. In Gnedin and Abel (2001), the authors developed a moment-based solver coupled to a cosmological hydrodynamics code, called Optically Thin Variable Eddington Tensor or OTVET. A similar technique has also been implemented for the SPH code GADGET (Petkova and Springel, 2011). The idea is to solve the stationary radiative transfer equation in the full domain assuming an optically thin medium, so that the problem boils down to a collection of 1/r² sources combined together. This can be achieved efficiently using any gravity solver. Once the corresponding Eddington tensor is computed, the moments equations are solved using a finite difference scheme. In Jiang et al. (2012), on the other hand, the authors developed a moment-based solver for the ATHENA MHD code, using VET together with the short characteristic method of Davis et al. (2012) to compute the Eddington tensor. Obviously, this leads to significantly more accurate results in optically thicker environments. An intermediate solution has been implemented in the TreeCol code (Clark et al., 2012), for which a gravity tree code is used to collect column densities between particles.

The second approach relies on simple but approximate models for the Eddington tensor. The simplest one, called M0, assumes the Eddington tensor is isotropic, with

$\begin{array}{ll}\text{M}0\text{ closure }:𝔻=\frac{1}{3}𝕀,& (50)\end{array}$

as it is the case for example in the diffusion limit (see below). A slightly more elaborate model is called the M1 closure (Dubroca and Feugeas, 1999; Ripoll et al., 2001). It assumes that the radiation intensity can be fitted by a Lorentz-boosted Planckian distribution, in other words an ellipse in angular space. The parameters of this distribution are determined by matching the radiation energy and the radiation flux. This leads to

$\begin{array}{ll}\text{M}1\text{ closure }:𝔻=\frac{1-\chi }{2}𝕀+\frac{3\chi -1}{2}\text{n}\otimes \text{n}.& (51)\end{array}$

This closure captures the optically thick regime, but also optically thin conditions, in case there is one dominant source of radiation. Indeed, in this case, one gets the free-streaming regime for radiation with 𝔻 ≃ n ⊗ n. The unit vector n and the parameter χ are determined using the radiation energy and the radiation flux as

$\begin{array}{ll}\text{n}=\frac{{\text{F}}_{\nu }}{\left|{\text{F}}_{\nu }\right|},\chi (f)=\frac{3+4f^2}{5+2\sqrt{4-3f^2}},f=\frac{\left|{\text{F}}_{\nu }\right|}{c{E}_{\nu }}& (52)\end{array}$

A particularly interesting property of the M1 model is that it leads to a system of hyperbolic conservation laws, that can be numerically integrated using Godunov schemes. The M1 model has also some major caveats, such as radiation shocks, in case multiple strong sources are present, or incorrect radiation geometry, in case of sharp shadows. This model was introduced for the first time in astrophysics by González et al. (2007) in the HERACLES code. It was then adapted to cosmic reionization by Aubert and Teyssier (2008) in the ATON code, and ported to GPU into the RAMSES-CUDATON code (Aubert and Teyssier, 2010), used now routinely for simulations of the Epoch of Reionization (EoR) (e.g., Ocvirk et al., 2016). Later on, the M1 method has been ported to the AMR code RAMSES, leading to the design of the RAMSES-RT solver (Rosdahl et al., 2013; Rosdahl and Teyssier, 2015). The M1 method has also been ported to ATHENA by Skinner and Ostriker (2013) and recently to AREPO by Kannan et al. (2018), the latter based on a new higher-order implementation of the M1 closure.

4.5. Flux-Limited Diffusion Methods

A definitive advantage of moment-based radiative transfer methods is the possibility to model both optically thin and optically thick regimes. For optically thick conditions, the radiation field follows the diffusion limit of the radiative transfer equation. In that limit, the radiation field in the comoving frame is quasi-isotropic, so that $𝔻\simeq \frac{1}{3}𝕀$ and the radiation flux can be approximated by

$\begin{array}{ll}{\text{F}}_{\nu }^0\simeq -\frac{c}{3{\alpha }_{\nu }^0}\nabla E_{\nu }^0.& (53)\end{array}$

Since now the flux in the comoving frame is a direct function of the energy density in the comoving frame, we do not need the radiation flux equation anymore and we are left with only the radiation energy equation. The Lorentz transform can be written to leading order as

$\begin{array}{ll}{E}_{\nu }\simeq E_{\nu }^0& (54)\end{array}$

$\begin{array}{ll}{\text{F}}_{\nu }\simeq -\frac{c}{3{\alpha }_{\nu }^0}\nabla E_{\nu }^0+\frac{4}{3}{E}_{\nu }^0\text{v}.& (55)\end{array}$

Injecting these into the radiation energy equation leads to a fully covariant formulation in the diffusion limit

$\begin{array}{l}\frac{\partial }{\partial t}{E}_{\nu }^0+\nabla ·\left(\frac{4}{3}{E}_{\nu }^0\text{v}\right)-\nabla ·\left(\frac{c}{3{\alpha }_{\nu }^0}\nabla E_{\nu }^0\right)\\ =\left(4\pi j_{\nu }^0-{\alpha }_{\nu }^{0}c{E}_{\nu }^0\right)+\frac{1}{3}\nabla E_{\nu }^0·\text{v},& (56)\end{array}$

which can be written in the familiar form

$\begin{array}{ll}\frac{\partial }{\partial t}{E}_{\nu }^0+\nabla ·\left(E_{\nu }^0\text{v}\right)+\frac{1}{3}{E}_{\nu }^0\nabla ·\text{v}=\nabla ·\left(\frac{c}{3{\alpha }_{\nu }^0}\nabla E_{\nu }^0\right)+\left(4\pi j_{\nu }^0-{\alpha }_{\nu }^{0}c{E}_{\nu }^0\right).& (57)\end{array}$

This is the radiation energy equation in the diffusion limit using only comoving radiation variables. If one wants to extend this equation outside its validity range (high optical depth), then one can use Flux Limited Diffusion (FLD), for which the flux function is modified as

$\begin{array}{ll}{\text{F}}_{\nu }^0\simeq -\frac{c\lambda (R)}{{\alpha }_{\nu }^0}\nabla E_{\nu }^0\text{ and }R=\frac{|\nabla E_{\nu }^0|}{{\alpha }_{\nu }^0{E}_{\nu }^0},& (58)\end{array}$

where λ(R) is the flux limiter, a function that has to connect the diffusion limit with λ ≃ 1/3 for R ≃ 0 to the free-streaming regime where λ ≃ 1/R for R → + ∞. A possible form for the flux limiter has been proposed by Minerbo (1978) and reads

$\begin{array}{ll}\lambda =\frac{2}{3+\sqrt{9+12R^2}}\text{ for }0\le R\le \frac{3}{2},& (59)\end{array}$

$\begin{array}{ll}\lambda =\frac{1}{1+R+\sqrt{1+2R^2}}\text{ for }\frac{3}{2}\le R\le +\infty .& (60)\end{array}$

while the Levermore (1984) flux limiter reads

$\begin{array}{ll}\lambda =\frac{1}{R}\left(cothR-\frac{1}{R}\right)& (61)\end{array}$

Historically, FLD has been the first method implemented in radiation hydrodynamics codes. Here again, the ZEUS 2D code was a precursor (Stone et al., 1992a). For AMR, the ORION code was developed specifically for FLD in collapsing star forming core by Krumholz et al. (2007b), using the mixed frame formulation outlined above or in Mihalas and Klein (1982). A version of the RAMSES code with FLD was developed later by Commerçon et al. (2011), with a particular emphasis on adaptive time stepping together with implicit time integration described in Commercon et al. (2014). Because of its simplicity, FLD has been used for many detailed studies of gravitational collapse with the effect of radiation included. The current frontier for FLD is probably the multigroup treatment of radiation, allowing for more accurate models of the full spectral energy distribution (Shestakov and Offner, 2008; Zhang et al., 2013; González et al., 2015). In the context of multifrequency radiative transfer, hybrid solutions have also been proposed, with ray-tracing techniques dedicated to the high-energy radiation field, and FLD for the infrared, dust-absorbed radiation field (Kuiper et al., 2010b; Klassen et al., 2014).

FLD has also been developed for SPH by Whitehouse and Bate (2004), Whitehouse et al. (2005), and Mayer et al. (2007), where most of the difference with the many grid-based versions lies in the radiation energy gradient that uses the SPH kernel. A simplified version based on a local estimate of the column density has also been proposed by Stamatellos et al. (2007) which avoids entirely the need for an implicit radiation solver. This method can be considered as an intermediate solution between the polytropic equation of state that we introduced at the beginning of this section and FLD. The main caveat here is the complete lack of any proper transport mechanism for radiation. In a similar spirit, Dale et al. (2007) also developed an approximate method for ionizing radiation, with Stromgren sphere iteratively grown around each star to locate photo-ionized, photo-heated gas.

5. Second Collapse

From first principles, the ultimate goal of star formation is the formation of the protostar itself. The major difficulty sits in the huge dynamical range in physical and temporal scales that have to be described: protostars with radii of about 1 R_⊙ form within dense cores of sizes ≃ 0.1 pc. So there are two ways to deal with protostars in star formation. The first one is to try to follow the collapse down to the stellar scales with the best numerical model accounting for the complex physical processes at stack combined with a very high numerical resolution at a cost of very short horizon of predictability. The second is the opposite: the physics and numerical resolutions are degraded, but the models are integrated over longer dynamical timescales to study the impact of protostellar feedback on the ISM dynamics. In the next three sections, we review the work that has been done using these two approaches, as well as attempts to account for protostellar evolution in large-scale models.

5.1. Historical Work and 1D Studies

Contemporary numerical star formation studies started in the late 60s with the outstanding pioneering work of Larson (1969) who first computed numerically the collapse of a dense core down to the formation of the protostar. He used a 1D spherically symmetric model, accounting for coupled gas dynamics and radiative transfer in a modified Eulerian scheme. Larson identified two distinct stages during the protostellar collapse, called the first and the second collapse. Each of these two stages are followed by the formation of a hydrostatic object, commonly referred to as the first and second Larson cores in post (Larson, 1969) studies. Larson's work established an empirical evolutionary sequence as follows. Dense cores first collapse isothermally because dust thermal emission is very efficient at radiating away the gas compression energy (dust and gas are thermally coupled within dense cores, see for instance Galli et al., 2002). The first collapse stops at densities of the order of ≃ 10⁻¹³ g cm⁻³ at which dust grains become opaque to their own radiation so that the gas begins to heat up quasi-adiabatically. This is the formation of the first Larson core. The first collapse lasts roughly a free-fall time. The first core accretes matter and its temperature increases up about 2000K where H₂ dissociation starts. This endothermic reaction enables the gas to behave as a fluid with an effective polytropic index smaller than the critical value 4/3 for collapse. The second collapse thus starts at typical densities of ≃ 10⁻⁸ g cm⁻³ until stellar densities are reached within ≃ 100 yr to form the second core, i.e., the protostar. The typical properties of first Larson cores at the start of the second collapse are a size ≃ 10 AU, a mass ≃ 0.01 M_⊙, and a lifetime ≃ 1000 yr. We note that first Larson cores are predicted by theoretical and numerical studies, but there is no observational confirmation of such objects even though a few sources are considered as good candidates (Tsitali et al., 2013; Gerin et al., 2015; Maureira et al., 2017).

Since Larson's work numerous 1D calculations using spherical symmetry have been conducted. We note the work in the 80s by Stahler et al. (1980) and Winkler and Newman (1980) in which Larson predictions were confirmed to establish the current empirical evolution sequence of collapsing low-mass dense cores. This evolutionary picture is still currently the commonly admitted scenario for low-mass star formation. For massive star formation, recent work indicates that first cores may not have time to form because of the large accretion rates (Bhandare et al., 2018). Note that in the reminder of the section dedicated to second collapse, we mention only the work which is consistent with the Larson evolutionary picture, i.e., where the gas thermal and chemical budgets are taken into account, and we do not mention the work that has been done using an isothermal approximation.

Modern studies began with the calculations performed in Masunaga and Inutsuka (2000) who incorporated a more accurate radiation transport scheme as well as a realistic gas equation of state which accounts for H₂ dissociation. Masunaga and Inutsuka (2000) integrated their models throughout the Class 0 and Class I phases up to an age of 1.3 × 10⁵ yr. Nowadays, 1D models are still used either as a first step toward describing more accurately the physics of the protostellar collapse, e.g., with multigroup radiative transfer (Vaytet et al., 2013, 2014) or with a view to provide quantitative predictions for the first and second Larson core properties (Vaytet and Haugbølle, 2017; Bhandare et al., 2018). These models are not pushed in time as far as the ones by Masunaga and Inutsuka (2000).

5.2. From 1D to 3D

Pioneered studies acknowledged up front that 1D spherical models were limited because the effect of rotation, turbulence, and magnetic fields cannot be taken into account. Three decades after Larson work, the first 3D numerical simulation able to describe the formation of the protostar was performed by Bate (1998) using SPH. The model accounted for initial rotation, but magnetic fields were neglected and radiative transfer was crudely mimicked by a piecewise polytropic equation of state. It is worth to notice at this point that grid-based codes took some time to reproduce (Bate, 1998) results for almost one decade. SPH is indeed naturally well-suited to collapse problems thanks to its Lagrangian nature while innovative techniques, such as AMR and nested grids, were not mature enough to capture the entire dynamical range and to cover eight orders of magnitude in physical length given the computer capabilities at that time. As mentioned in the section 2, ideal MHD was introduced early in protostellar collapse models using grid-based codes (e.g., Dorfi, 1982; Phillips and Monaghan, 1985). The first calculations of second collapse incorporating magnetic fields were thus done on a grid-based code by Tomisaka (2002) using a nested grid with cylindrical coordinates assuming axisymmetry (Tomisaka and Bregman, 1993). Tomisaka (2002) code solved the ideal MHD equations using a second-order “monotonic scheme” (van Leer, 1977; Winkler and Norman, 1986), a constrained transport scheme for the induction equation (Evans and Hawley, 1988), as well as an angular momentum preserving scheme. To capture the huge range in spatial scales, Tomisaka (2002) used up to 16 levels of refinement. He found that two different types of outflows were launched at the first and second cores scales. He suggested that the second core outflow corresponds to the optical jets observed in young stellar objects (YSOs), while the first core outflow corresponds to the molecular bipolar outflow. The first full 3D MHD models were presented in Machida et al. (2006) using nested grid with Cartesian coordinates and 21 levels of refinement. They compared the results obtained with resistive (Ohmic diffusion) and ideal MHD. They confirmed the results of Tomisaka (2002) even with non-ideal MHD included. After Machida et al. (2006) work, the numerical challenge of forming a protostar in 3D with MHD was faced, but the horizon of predictability remained very short since they had to stop their calculations 20 days after the protostar formation because of extremely small timesteps resulting from the high velocity of the jet propagating from the finest level. Banerjee and Pudritz (2006) also performed 3D collapse calculations down to stellar scales, using ideal MHD as well as molecular cooling, and they confirmed the results of Tomisaka (2002) and Machida et al. (2006).

All the multidimensional studies we mentioned so far were done assuming a piecewise polytropic equation of state to mimic the thermal behavior of the gas throughout the evolutionary sequence. Such barotropic equation of states are parameterized using the results of 1D spherical symmetry simulations previously mentioned, so that they cannot account for multidimensional effects such as the optical depth drop in the vertical direction that allows efficient cooling of the nascent protostellar disks and thus leads to incorrect results for the fragmentation (e.g., Boss et al., 2000; Commerçon et al., 2010; Tomida et al., 2010). The next level of complexity is to add a more accurate model for the thermal balance of the gas and dust mixture, i.e., to model radiative transfer (see section 4). The addition of a more accurate model for radiation transport can be crudely resumed to an additional heating/cooling term in the gas internal energy equation due to matter/radiation interactions. For the models which are restricted to the first collapse and first core formation, i.e., for temperatures less than ≃ 2000 K, the addition is conceptually straightforward for the gas evolution equations. An usual ideal gas equation of state can be used, with the only modification on the adiabatic index γ which varies with the temperature. At low temperature (T ≲ 150 K), the H₂ molecule behaves like a monoatomic gas (γ = 5/3) because only the vibrational degrees of freedom are excited. At higher temperature, the rotational degrees start to be excited and the adiabatic index decreases to γ = 7/5. Accounting for this non-constant adiabatic index is important regarding the fragmentation properties (Commerçon et al., 2010) and the first core properties (Tomida et al., 2010; Lee and Hennebelle, 2018b). The driver of the second collapse is the endothermic dissociation of the H₂ molecule, which modifies considerably the gas equation of state. Currently, there are two ways to integrate non-ideal EOS in second collapse calculations. The first one is to compute on-the-fly the thermodynamic quantities from equilibrium chemical abundances and to use an ideal gas equation of state, including the effects of dissociation and ionization, to compute the gas pressure and the internal energy (Black and Bodenheimer, 1975). Most of the second-collapse including RHD are done with this approximation (Masunaga and Inutsuka, 2000; Whitehouse and Bate, 2006; Stamatellos et al., 2007; Forgan et al., 2009; Tomida et al., 2013; Bhandare et al., 2018). The second one consists in using tabulated EOS table coming from detailed studies to account for non-ideal effects (e.g., ionization by pressure, interaction between particles) which are important at high density/pressure (Saumon et al., 1995). Saumon et al. (1995) EOS is used in the work by Vaytet et al. (2013), Vaytet et al. (2014), Vaytet and Haugbølle (2017), and Vaytet et al. (2018), and also in the SPH public code PHANTOM (Price et al., 2017). We note that these developments imply that the dependency of the specific heat capacity as a function of the physical conditions has to be taken into account. Caution should thus be taken within the radiation solvers which are often integrated using implicit schemes. The most common and simple way to deal with it is to compute a heat capacity-like factor, $\widetilde{C_{\text{V}}}\equiv e/T$ where e is the gas internal energy density, which is kept constant in the radiation solve. Readers can refer to Tomida et al. (2013) for a detailed implementation. Unfortunately, details on the numerical implementations are usually not reported in the literature. The first full 3D RHD calculations of the second collapse of an initial 1 M_⊙ dense core were performed by Whitehouse and Bate (2006) using the FLD approximation in a SPH code, followed by Stamatellos et al. (2007) who used a local radiative cooling approximation and SPH. Bate (2010) extended the work of Whitehouse and Bate (2006) beyond the formation of the stellar core for about 50 yr, using 3 × 10⁶ equal-mass SPH particles to satisfy the Bate and Burkert (1997) mass resolution criterion throughout the entire collapse. Interestingly, there is to date no 3D RHD calculations of the second collapse, i.e., without magnetic field, performed with a grid-based code.

Last but not least, it has been demonstrated in Commerçon et al. (2010) that the interplay between magnetic fields and radiative transfer is of primary importance for the thermal budget of the collapsing gas and of the first and second core accretion shocks (Commerçon et al., 2011a; Vaytet et al., 2012; Vaytet et al., 2013, 2018). The infall velocity is indeed greatly modified in the presence of magnetic fields which tends to focus the collapsing gas on a smaller area compared to the case without magnetic fields (Commerçon et al., 2011b; Myers et al., 2013). The magnetic braking mechanism transports angular momentum from the inner parts of the collapsing cloud to the envelop. As a consequence, the infall velocity increases along with the incident kinetic energy, which then modifies the thermal budget through the accretion shocks. The next step in the increasing complexity is to perform full 3D RMHD calculations. The first SPH RMHD calculations of protostellar core formation were performed in Bate et al. (2014), extending their previous work to ideal MHD. Tomida et al. (2013) and Tomida et al. (2015) performed the first 3D full RMHD calculations using a grid-based code, extending the work of Machida et al. (2006) to account for radiative transfer with the FLD approximation, as well as the Ohmic dissipation and the ambipolar diffusion. They used a finite-volume nested-grid code with standard MUSCL predictor-corrector scheme to achieve second-order accuracy.More details on Tomida et al. (2015) implementation for non-ideal MHD are given in section 3. For the implicit radiation update, they used a linear system solver with a combination of the BiCGStab solver and the incomplete LU decomposition preconditionner. Tomida et al. (2013) used a nested-grid constructed to ensure that the Jeans length is always resolved by at least 16 cells, with typical spatial resolution at the maximum level of refinement of Δx ≃ 6.6 × 10⁻⁵ AU (23 levels). When this resolution condition cannot be satisfied, they stopped the calculations. This is particularly the case when discs form and grow with time to spread over different level of refinements. As a consequence, they cannot follow the second collapse in their model with ambipolar diffusion since a large disc is formed. More flexible refinement techniques, such as AMR or SPH are thus required to continue the calculations with a reasonable numerical resolution. Currently, the state-of-the-art includes 3D RMHD calculations with resistive effects, done with both SPH (Tsukamoto et al., 2015b; Wurster et al., 2018) and AMR (Vaytet et al., 2018). The typical resolution in SPH calculations currently achievable with limited CPU resources is of 3 × 10⁶ equal-mass particles per M_⊙. In AMR, Vaytet et al. (2018) used a coarse grid of 64³ with 21 additional levels of refinement, and a refinement criterion of 32 points per Jeans length. While the SPH calculations are able to integrate the three non-ideal effects (Wurster et al., 2018) down to stellar scales, there is no AMR work describing the Hall effect in second collapse calculations (Vaytet et al., 2018 accounts for the ambipolar diffusion and the Ohmic diffusion). The calculations by Wurster et al. (2018) were integrated 17 yr after the stellar core formation, while Vaytet et al. (2018) could only integrate for 24 days. Beside the major achievement which represents the formation of a protostar on a computer, we must admit that these kind of studies are not yet capable to give long term predictability. In addition, the CPU time required to integrate full 3D RMHD models is relatively big. Clearly, second collapse study cannot be done on a frame covering a large parameter space. In addition, even at the first core scale, the integration timestep becomes so short in the adiabatic fragments that long time evolution calculations becomes prohibitive. In that view, the sink particles are a good alternative to reach a long horizon of predictability and to cover a wide parameter space.

6. Sink Particles

Sink particles are Lagrangian particles, and can accrete matter as well as angular momentum in order to ensure mass and angular momentum conservation. Sink particles should not be confounded with sink cells which are usually introduced in grid-based codes which use a spherical grid to deal with the singularity in the center (e.g., Boss and Black, 1982; Kuiper et al., 2010a; Li et al., 2011). The differences between sink particles and sink cells are discussed in section 6.4. Sink particles are often erroneously referred as protostars or stars in the literature whereas the resolution at which they are introduced is generally much larger than the stellar scales, typically 1−10 AU for collapse calculations, and a few 100 AU for cluster formation and ISM evolution calculations. Multiple systems (small clusters, binaries), as well as systems made of a disk around a star can thus be encompassed within a single sink particle depending on the resolution.

6.1. Standard Implementations

Sink particles were introduced in star formation calculations by Bate et al. (1995) to enable long time integration to study the fragmentation of collapsing clouds. Sink particles were introduced to mimic the second collapse and the protostars at scales that could not be reached in numerical simulations. They are treated as accreting non-collisional point masses and enable one to follow the evolution through the main accretion phase, with a view to compare with observations of clusters in which most of the gas has been accreted. Their original implementation was done in a SPH code, with three main parts in the algorithm: (i) sink particle creation, (ii) sink particle accretion, (iii) boundary conditions for sink particles. A SPH gaseous particle is tagged for sink particle creation if its density exceeds a specific density threshold (defined by the user) and satisfied a number of tests. First, the smoothing length of the SPH particle that is a candidate for sink particle creation has to be less than half the accretion radius r_acc of the sink particle. This ensures that the sink particle is formed from at least N_neigh. r_acc is defined before the calculations start and remains fixed. Its value is chosen by the user depending on the required level of resolution. It sets the smallest scale that can be resolved in the calculations. The knowledge of the flow within r_acc is lost and the gas is assumed to collapse beyond this scale to form protostars. Then, a series of tests is performed on the system composed of the particle and its neighbors to decide if it should create a sink particle. The tests ensure that the gas particles would continue to collapse if the sink particle was not created. First, the ratios α and β of the thermal energy and rotational energy to the gravitational energy must satisfy

$\begin{array}{ll}\alpha \le 1/2;\alpha +\beta \le 1.& (62)\end{array}$

Additionally, the total energy must be negative, and finally, the system must be contracting ($\nabla ·\overrightarrow{v}<0$). If all the tests are passed, a sink particle is formed from the system at the center of mass. Initially, the sink particle has thus a mass of ~N_neigh times the standard SPH particle mass. Afterwards, sink particles interact only with the gas particles via gravity. The sink particle is allowed to accrete mass at each integration timestep if a gas particle enters r_acc and passes several criteria: the particle must be the most bound to the sink particle (and not to another more distant but more massive sink particle), and the particle must have a moderate specific angular momentum. When a gas particle is accreted, its mass and linear momentum are added to the sink particle, as well as the angular momentum which modifies the spin of the sink particle. The sink particle is moved to the center of mass of the sink and the accreted gas particle. Sink particles also contribute to the computation of the total gravitational potential in the same way as standard gaseous particles. Bate et al. (1995) pointed out that the introduction of sink particles in SPH calculations may affect the gas outside the accretion radius, in particular because there is a discontinuity in the number of particles across the accretion radius. They proposed different types of boundary conditions to account for missing neighbors. Note that the use of boundary conditions in modern sink particles implementation is not used anymore thanks to the improvements that have been made, in particular in the case of optically thick flows (e.g., Hubber et al., 2013).

Sink particles were introduced in Eulerian calculations for star formation purposes by Krumholz et al. (2004) in the AMR code ORION. They introduce a Jeans length criterion for sink creation, based on the resolution study of Truelove et al. (1997). A sink is created within a cell at the maximum level of refinement if the density ρ exceeds the maximum Jeans density ρ_J that can be resolved

$\begin{array}{ll}{\rho }_{\text{J}}<J^2\frac{\pi c_{\text{s}}^2}{G\Delta x_{\text{min}}^2},& (63)\end{array}$

where J < 0.25 is the Jeans length resolution criterion, c_s the gas sound speed, and Δx_min the minimum cell size. The initial mass of the sink is $m_{\text{sink}}=(\rho -{\rho }_{\text{J}})\Delta x_{\text{min}}^3$. In Krumholz et al. (2004), no additional checks are performed to validate the sink creation. Their algorithm is coupled to a sink merging scheme based on a friends-of-friends (FOF) algorithm to handle situations where a block of contiguous cells create sink particles in a single time step. The FOF is performed after each time step to group all the sink particles (old and new), with a linking length equal to r_acc. All the groups of sink particles found by the FOF are then merged and replaced by a single particle at the center of mass of the group. All the merged sink particles quantities (mass, momentum, angular momentum) are also added conservatively to the new sink particle. Krumholz et al. (2004) sink accretion scheme was designed to handle situations where the flow onto the sink particle is subsonic (which is analogous to the case where boundary conditions have to be introduced in standard SPH implementation). They set the accretion rate using an approximate formula from Bondi-Hoyle accretion which handles also the regime where the flow is supersonic. The accretion rate is estimated using average properties in the region within r_acc. An important point in Krumholz et al. (2004) implementation is the description of the accretion zone. The mass accreted by the sink particle comes from all the cells in this accretion zone. Based on their experiments, Krumholz et al. (2004) use a value of r_acc = 4Δx. In the case where the Bondi-Hoyle radius is smaller than the accretion zone, it is incorrect to estimate the accretion rate from all the cells in the entire accretion zone. Krumholz et al. (2004) thus use an accretion kernel where each cell within r_acc is assigned a weight from a Gaussian-like function. The kernel is used to compute the average quantities to determine the mass accretion rate. The mass to be accreted is redistributed within the accretion zone onto virtual cloud particles (each cell in the accretion zone is divided into 8³ cloud particles). Additional checks can be then performed to avoid to accrete gas that has for instance too much angular momentum to be gravitationally bound. When a parcel of gas is accreted, it also transfers its linear and angular momentum to the sink. An absolute mass cap is also used to limit to 25% the amount of mass accreted from a unique cell in a single time step. After accretion, the position of the sink particle is changed in two steps according to its momentum after gas accretion and through gravity. First the momentum is updated to account for the gravitational interaction with the gas. The total force on a sink particle is computed by a direct summation over all the cells in the computational domain, except the ones in the accretion zone. This method remains practicable as long as the number of sink particles created remains limited. For the accretion zone, the interaction is computed between the sink and the cloud particles created during the accretion step. A Plummer law is used to soften the gravitational interaction at short distance, i.e., $F_{\text{grav}}~1/(r^2+{ϵ}^2)$ where ϵ = 2Δx_min is the softening length. Second, the acceleration due to particle-particle interaction is added to the sink particle momentum using a Bulirsch-Stoer algorithm with adaptive time step and error control. Note that the stability of the temporal discretisation scheme required that the sink particle velocity satisfies a CFL-like condition (v_sinkΔt < CΔx, with C < 1).

Sink particles algorithms have been complemented in the past 10 years to improve their versatility, in particular in the transition between the sub- and super-sonic regimes. All newer implementations have been built upon these two pioneering works.

6.2. Contemporary Implementations

In all current implementations, the sink creation and accretion are the most sensitive algorithms. Of course, the sink particles can be introduced in the initial conditions. Otherwise, they are created after the creation checks.

Each of the commonly used codes in star formation has its own implementation of sink particles. For the grid-based codes, we note the work of Wang et al. (2010) in the ENZO code, which is very similar to Krumholz et al. (2004) but with a different merging scheme which involves a maximum mass sink mass for merging. Federrath et al. (2010) presented an implementation in the AMR code FLASH, where the creation checks were similar to the ones by Bate et al. (1995) for SPH with two additional checks: the cell tagged for sink creation should be at the maximum level of refinement, and it should sit in a local gravitational potential minimum. They also account for magnetic energy in the energy budget done in the check algorithm. The checks are performed on a spherical region of radius r_acc around the cell candidate for sink creation. Their scheme for accretion is also different from the previous grid-based works. If a cell i within r_acc exceeds the density threshold ρ_res for sink creation, the mass increment $\Delta M=({\rho }_i-{\rho }_{\text{res}})\Delta x_i^3$ is considered for accretion onto the sink particle. The mass increment has to be gravitationally bound to the sink, and in the case of overlapping sink particles, it is accreted to the most gravitationally bound sink particle, the latter being moved to the center of mass of the particle-gas configuration before the accretion step. They also improve the treatment of the gas-sink gravitational interaction. Similarly to Krumholz et al. (2004), the sink-sink interaction is done by direct summation over all the sink particles. The gravitational acceleration of the sink due to the gas is estimated from the gravitational potential of the gas which is handled by the Poisson solver of FLASH to handle the gas-gas interaction. The acceleration is interpolated using a first-order cloud-in-cell method for each sink particle. Last, the acceleration of the gas due to the sink is done by a direct summation. Each sink particle contributes to the acceleration of each cell in the computational domain. The direct summations for the sink-sink and sink-gas interactions requires a gravitational softening. Federrath et al. (2010) use a cubic spline softening similar to what is used in SPH, which is less aggressive within the softening length than a simple Plummer profile. Federrath et al. (2011a) updated their first implementation to ensure exact momentum conservation using a direct summation between all cells and sink particles for the gas-sink interaction as well.

In addition, Federrath et al. (2010) performed a quantitative comparison between their implementation and the one by Bate et al. (1995) by comparing the results of two star formation calculations performed with FLASH and with a SPH code. They find a good quantitative agreement which was quite encouraging given the fundamental differences between SPH and AMR sink implementations.

Hubber et al. (2013) presented an improved algorithm for the sink particles in SPH, which can handle situations where sinks are introduced after the gas has become adiabatic. They propose a different combination of creation criteria which are based on the studies we just mentioned: (i) density threshold, (ii) no overlapping with another sink accretion volume, (iii) the gaseous particle flagged for sink creation sits in a minimum of gravitational potential (on a volume defined by its neighbors), (iv) the gaseous particle is outside the Hill sphere of fragments harboring another already formed sink.

Last, Bleuler and Teyssier (2014) presented an innovative sink creation scheme based on a clump finder algorithm used on-the-fly in the RAMSES code. The clump finder is performed in a first step to identify peaks and their associated regions. Second, the peaks are considered for sink creation. A virial theorem type analysis is performed on the clumps, accounting for surface pressure as well as tidal forces. The virial check is passed if the second derivative of the moment of inertia in the center of mass frame is negative. Then a collapse check is performed, which corresponds to the contraction check of Krumholz et al. (2004) and Federrath et al. (2010), but adapted according to the virial analysis. Last, the final check is a proximity check. We also note that Bleuler and Teyssier (2014) offer also the possibility to use the alternative sink creation and accretion schemes reported in the literature. They provided an excellent review of the different creation checks, accretion, sink merging and trajectories schemes. In particular, they tested different schemes to perform sink accretion: Bondi-Hoyle accretion, flux accretion, and density threshold accretion. They show that all methods give good results for spherical Bondi accretion (i.e., when the flow is supersonic). In the subsonic regime, only the Bondi formula gives a correct accretion rate. They recommend to automatically adapt the accretion scheme depending on the surrounding flow properties.

A number of other sink particle implementations have been reported in the literature that we cannot detail in this review. We give in the following a non-exhaustive list, restricted to the field of star formation. In SPH, sink particles have been implemented in the codes GADGET (Jappsen et al., 2005), DRAGON (Goodwin et al., 2004), and SEREN (Hubber et al., 2011). For grid-based codes, sink particle are used in the uniform grid code ATHENA by Gong and Ostriker (2013), in RAMSES by Padoan and Nordlund (2011), and in the ORION2 code by Lee et al. (2014). Last, Greif et al. (2011) also proposed a sink particle algorithm within the AREPO moving-mesh code. Their implementation is very similar to the one used in SPH codes. We also note that sink particle merging is now widely used in most of the implementations.

Twenty years after their introduction, sink particles remain a sensible subject in the community. It is commonly accepted that an accreting sink particle should have a minimal impact on the collapsing gas around it. The easiest way to deal with it is to create the sink particle during the (quasi-) isothermal phases of the collapse (either the first or the second collapse), when the gas velocity is supersonic. The perturbation created by the sink particles on the gas dynamics cannot propagate since the CFL timestep is limited by the supersonic motions. An example of bad behaviors observed in protostellar collapse SPH calculations using sink particles introduced in the vicinity of the first cores can be found in appendix D of Commerçon et al. (2008). When a barotropic law is used to mimic the gas thermal budget, a simple method used to introduce sink particles in the vicinity of hydrostatic object can be found in Price and Bate (2008). They artificially change the gas equation of state to isothermal at an arbitrary density (e.g., 10⁻¹¹ g cm⁻³) and set the density threshold for sink creation to a value of more than two orders of magnitudes larger.

6.3. Sink Particles and Magnetic Fields

In the previous section, we did not discuss the sink particle algorithms in the presence of magnetic fields. We try to describe here some attempts done to provide a good description of MHD flows for sink creation and accretion. Sink particles can of course be introduced in MHD calculations, but the implementation is not as straightforward as in the hydrodynamical case. While it is relatively easy to accrete gas and momentum, magnetic flux accretion cannot bypass the solenoidal constraint $\nabla ·\overrightarrow{B}=0$.

Price and Bate (2008) and Wang et al. (2010) performed star formation calculations in turbulent and magnetized clouds, but do not account for any contribution of the magnetic fields in their sink algorithms. Federrath et al. (2010) have introduced the magnetic energy in their negative total energy check for sink creation, but not in the accretion scheme. Lee et al. (2014) have extended the work of Krumholz et al. (2004) to deal with magnetized flows. For sink creation, they follow the work of Myers et al. (2013) who derived a magnetic Truelove criterion. They derive simple steady state accretion rates as a function of the plasma beta β and the Mach number ${M}$.

In all implementations, gas is accreted by the sink particles but not magnetic flux. It stays the gas in the surrounding of the sink particles, and magnetic flux accumulates as collapse goes on. While at scales of a few 10-100s AU this flux accumulation is not problematic (Federrath and Klessen, 2012; Li et al., 2018), it eventually leads to a magnetic explosion and magnetic flux redistribution in the surrounding of the sink (Zhao et al., 2011) due to the development of an interchange instability (Spruit and Taam, 1990; Spruit et al., 1995; Li and McKee, 1996), which is also observed in high resolution models without sink particles and ideal MHD (Li et al., 2014; Masson et al., 2016). In addition, since density stays roughly constant while magnetic intensity increases in the sink accretion zone, the Alfvén speed increases significantly, leading to non-physical accretion rates, as well as short integration time steps. Lee et al. (2014) proposes to cap the accretion rate at a maximum value corresponding to a constant Alfvén speed after mass accretion. This trick tempers the accumulation of magnetic flux as well.

In the recent years, resistive MHD calculations have shown that all non-ideal effects leads to a decrease of the magnetic flux at the first core scale. Magnetic fields decouple from the gas and dust mixture dynamics. The magnetic flux does not accumulate at the first core border, but rather at a larger distance which corresponds to the radius at which non-ideal effects start to dominate (Hennebelle et al., 2016; Masson et al., 2016). Matter is accreted but leaves its associated magnetic flux outside the decoupling region. Adding sink particles in restive MHD calculations with resolution capable of resolving the decoupling region has not yet been investigated in detail (Machida et al., 2009; Matsumoto et al., 2017; Tomida et al., 2017), but one can foresee that introducing a sink within the decoupling region should prevent magnetic flux accumulation. It is expected that non-ideal MHD can alleviate the conditions of interchange instability, but is then necessary to resolve the decoupling scale, i.e., the disk. To date, the most detailed study of sink particles introduction in MHD calculations can be found in Machida et al. (2014).

6.4. Sink Particles vs. Sink Cells

A fundamental difference between sink cells and sink particles is that sink cells are fixed in position and have physical boundaries. There are pros and cons for the two methods, which we briefly summarize in the following.

First, the sink cell can be considered as a boundary in codes using a spherical grid. In collapse calculations, a sink cell is most often introduced as a passive boundary with zero gradients, or also called outflow (Yorke et al., 2002; Kuiper et al., 2010a; Li et al., 2011). The mass flux across the boundary is added to the central point mass of the sink cell and used to compute the gravitational acceleration. The sink cell enables one to make more detailed models for the magnetic field configuration at the boundary, as for instance in star-disk interaction studies (e.g., Mellon and Li, 2008; Takasao et al., 2018). In principle, this is a quite powerful numerical set-up to overstep the issue of time integration after the second core formation. The current state-of-the-art is not yet mature to design such numerical studies, but this should certainly be the purpose of future works.

On the other hand, a major inconvenience of sink cells resides in their fixed position in time, which forces the center of mass to be equal to the geometrical center of the disk. As a consequence, the development of non-axisymmetric azimuthal low wave-number modes may be artificially damped, which alters the development of eccentric gravitational instabilities (Adams et al., 1989; Shu et al., 1990). This effect can be particularly problematic in the case of young stellar objects with disks of mass comparable to the stellar mass, M_disk/M_star > 1/3 (Krumholz et al., 2007a; Kuiper et al., 2011; Sigalotti et al., 2018).

Last, the use of “outflow” boundary conditions for the mass accretion on the sink cell is very different from the algorithms we just presented for the sink particles. Machida et al. (2014) have shown that applying similar criteria for the sink accretion than the ones employed in studies with a sink cell (Li et al., 2011) may lead to very different results, i.e., no disk formation, compared to the same calculations done with a classical sink particle algorithm, i.e., formation of a disk.

7. Sub-Grid Models for Sink Particles

Beyond the practical advantage of saving computational time, the sink particle has been shown to be very useful in works where (proto-)stellar feedback is considered. Similar to cosmological and galactic evolution studies, the sub-grid models can be associated to the sink particle evolution. We distinguish two types of sub-grid physics, depending on the feedback origin: radiative (luminosity, accretion shock, ionization), and dynamical (jets/outflow, winds). In any case, the source terms introduced by sub-grid physics are treated using operator splitting schemes.

7.1. Radiative Feedback

The first class of sub-grid models were designed to account for the radiative feedback of the nascent protostars. The radiative feedback term includes the accretion, the internal, and the ionization luminosities. When integrating over discrete timesteps, the luminosity is an energy or radiative flux input. Radiative feedback was introduced early in 2D models (e.g., Bodenheimer et al., 1990). They used a 2D axisymmetric hydrodynamical code (Rozyczka, 1985), which was primarily designed for stellar winds. Bodenheimer et al. (1990) used the gray FLD approximation for the radiation transport. They combined the evolution of their central zone, or sink cell, with a pre-main sequence evolution model based on Maclaurin spheroids in hydrostatic equilibrium. They computed the structure of the protostar according to the central (ρ, T) path tracks reported by Winkler and Newman (1980) from 1D spherical calculations. The protostar mass M_⋆ and equatorial radius R_⋆ were then used to compute the accretion luminosity L_acc = 0.5GM_⋆Ṁ/R_⋆, where Ṁ is the mass flux crossing the inner boundary. Bodenheimer et al. (1990) further assumed that the radiative energy input LΔt was only due to the accretion luminosity and that all the infall kinetic energy was converted into radiation at the stellar core accretion shock. They found that most of the energy that went for heating in the envelope comes from this protostellar luminosity. Yorke and Bodenheimer (1999) and Yorke et al. (2002) used a more sophisticated protostellar evolution model which accounts for the intrinsic luminosity of the protostar L_int. They realized that the central luminosity evolution sets the thermal budget within the envelop, but could not resolve the innermost regions of the star-disk system. They designed a protostellar accretion luminosity model which accounts for the luminosity coming from Keplerian motion of the disk being converted into heat and radiated away at the disk-core boundary layer following (Adams and Shu, 1986)

$\begin{array}{ll}{L}_{\text{tot}}=L_{\text{int}}+\frac{3G{M}_{\star }Ṁ}{4R_{\star }},& (64)\end{array}$

where they assume that a fourth of the total potential energy of the accreted material is released into heat within the disk, the remaining being radiated away within the unresolved inner disk region and at the accretion shock. In addition, Yorke and Bodenheimer (1999) and Yorke et al. (2002) followed radiation transport using a frequency dependent ray-tracing algorithm (with 64 frequency bins). In particular, Yorke et al. (2002) demonstrated the importance to handle the protostellar irradiation using a frequency dependent irradiation scheme in the context of massive star formation. The flashlight effect, i.e., radiation escapes in the vertical direction and the radiative acceleration is reduced within the disk to allow accretion, is enhanced compared to a simple gray model. In addition, Krumholz et al. (2005) showed how protostellar outflows allow radiation to escape in the outflow cavities to allow a continued accretion onto the central massive stars. These results have been further confirmed by 3D dynamical simulations of Krumholz et al. (2007a), Kuiper et al. (2011, 2016), and Rosen et al. (2016).

Currently, most 3D numerical RHD calculations include stellar radiative feedback sub-grid models attached to sink particles/cells evolution. Most of the implementations account for both the internal and the accretion luminosities, with some modulation factors on the amount of potential energy radiated away and on the efficiency of the mass transfer for the sink accretion radius to the protostar. The simpler sub-grid models use a gray FLD model for the radiation transport where the protostellar luminosity is simply a source term in the radiative energy equation (Krumholz et al., 2007a; Offner et al., 2009; Stamatellos et al., 2011; Fontani et al., 2018; Jones and Bate, 2018). Frequency dependent protostellar irradiation modules have been developed using ray-tracing to compute the radiation fields, combined with moment models for radiation hydrodynamics (Kuiper et al., 2010a; Wise and Abel, 2011; Klassen et al., 2014; Ramsey and Dullemond, 2015; Buntemeyer et al., 2016; Rosen et al., 2017). These frequency dependent modules have been used primarily in the context of massive star formation (Klassen et al., 2016; Rosen et al., 2016).

Last, we note the developments made to capture ionizing radiation feedback from massive stars: Dale et al. (2007) used a Strömgren volume method in SPH, Peters et al. (2010) ray tracing in the FLASH code, Kuiper and Hosokawa (2018) an hybrid ray-tracing and FLD irradiation module in the PLUTO code, Geen et al. (2015) the M1 moment method in the RAMSES code, and Harries (2015) Monte Carlo radiative transfer in the AMR code TORUS. These schemes are applied first in isolated collapse calculations (Kuiper and Hosokawa, 2018) and, mostly, in cluster formation studies (Dale and Bonnell, 2011; Peters et al., 2011; Dale et al., 2012; Geen et al., 2015, 2018; Ali et al., 2018).

7.2. Dynamical Feedback

The most advanced MHD collapse calculations have shown that outflows and jets naturally develop at different scales during the collapse (Tomisaka, 2002; Banerjee and Pudritz, 2006; Machida et al., 2006; Hennebelle and Fromang, 2008; Ciardi and Hennebelle, 2010; Commerçon et al., 2010; Tomida et al., 2010; Price et al., 2012). In addition, outflows and jets are commonly observed in YSO's (it is for instance a selection criterion for Class 0 sources, Andre et al., 1993) and they are considered as possible sources of turbulence driving in star forming clouds and participate in the regulation of the star formation rate (Matzner, 2007; Nakamura and Li, 2007; Krumholz, 2014; Federrath, 2015). Outflows and jets launching scales are nevertheless not captured in most studies because it requires a very high numerical resolution to reach sub-AU scales. Besides, studies neglecting magnetic fields cannot generate MHD winds driven centrifugally or by toroidal magnetic pressure. Consequently, very few numerical calculations have been able to launch self-consistently MHD outflows (e.g., Hennebelle et al., 2011). Sub-grid models have thus been developed to account for outflows/jets in numerical works with unresolved launching scales and missing physics, similarly to what has just been presented for the radiative feedback. The numerical implementation are based on analytical works dedicated to MHD wind launching (Blandford and Payne, 1982; Pelletier and Pudritz, 1992; Shu et al., 1994; Ferreira, 1997; Matzner and McKee, 1999).

The first works on sub-grid models have been reported in Li and Nakamura (2006) and Nakamura and Li (2007) who employed ideal MHD and self-gravity. They attached a protostellar outflow model on sink particles following (Matzner and McKee, 2000) because they were unable to describe the self-consistent launch of outflows with the crude resolution they use (128³ for a ≈ 2 pc box). Each star injects in the ambient medium a momentum that is proportional to the stellar mass M_⋆ with a two-component outflow model. The volume of injection is given by the sink accretion volume (which they refer to as “supercell”). A fraction η of the total outflow momentum is put in a conical collimated jet component to facilitate the transport of energy and momentum at large distances. The rest of the momentum is put into a spherical component around the protostar. The direction of the jet is given by the orientation of the magnetic field in the central cell, and they assume an opening angle of 30° about this direction.

Most of the implementation of outflows and jets have a similar construction, with flavors depending on the sink particle algorithm. Wang et al. (2010) implemented in the MHD ENZO code a simplified version for protostellar outflows where they only account for the collimated jet component. They adopt a continuous momentum injection ΔP = P_⋆ΔM where P_⋆ is a proportionality constant which depends on the stellar mass $~M_{\star }^{1/2}$). If the momentum injection occurs in a cell with a too low density, they take 10% of ΔM out of the sink particle and redistribute it evenly between the injection cells in order to avoid too large outflow speeds (and too small integration timestep). Along the same line, Cunningham et al. (2011); Hansen et al. (2012) extended the work of Krumholz et al. (2005) and Krumholz et al. (2007a) and presented for the first time a combined sub-grid model which includes protostellar outflows (Matzner and McKee, 1999) and protostellar radiative heating (Offner et al., 2009), ignoring magnetic fields though. Their sub-grid model accounts for pre-main sequence evolution as in Offner et al. (2009) and the calculated protostellar radius is used in the outflow and radiative heating models. They parameterized their outflow model with dimensionless parameters which, for instance, sets the fraction of mass accreted by the protostar or launched in the outflow, and the ejection velocity as a function of the Kepler speed at the surface. Contrary to the case of radiative heating where energy is put within the accretion volume, they inject outflows at a distance comprised between 4Δx and 8Δx from the sink particle, and used angular dependency from Matzner and McKee (1999). We note that they not only inject momentum, but also mass and thermal energy. They set the wind temperature to 10⁴ K, which is appropriate for an ionized wind. Using the same tool, Myers et al. (2014) added magnetic field in their simulations with ORION, but because of the low resolution they used (> 20 AU), they where not able to launch self-consistent outflows and used the sub-grid model we just mentioned. The outflows they produce may form strong shocks at very high temperature, higher than the dust sublimation temperature so that the dust opacity drops. Their radiative transfer scheme, primarily designed for dust thermal emission does not allow the gas to cool efficiently. They thus needed to add a line cooling function for the thermal budget of the gas.

Federrath et al. (2014b) implemented a sub-grid-scale outflow model upon the sink particles algorithm they implemented in FLASH (Federrath et al., 2010) for MHD collapse. Their model combines different features of the Li and Nakamura (2006), Nakamura and Li (2007), and Cunningham et al. (2011) implementations. They used a normalized velocity profile for momentum injection, which reproduces the two components outflow/jet with opening angles of 30° and 5°, respectively. They also improved the algorithm for the outflow orientation by recording the angular momentum transfer from the accreted gas to the sink particle. Their algorithm also conserves mass and momentum exactly in the sink accretion step. Murray et al. (2018) proposed a similar implementation in RAMSES. They follow (Offner et al., 2009; Cunningham et al., 2011; Federrath et al., 2014b; Myers et al., 2014) for the outflow injection as well as the protostellar evolution. Outflows are there injected within a conical volume about the spin axis of the sink particle in a radial extent comprised between 4 and 8 cells (at the highest refinement level) away from the sink particle. Federrath et al. (2014b) performed a parameter study on the number of cells for the radial outflow direction and found that convergence on the maximum outflow velocity requires 16 cells per outflow radius. Convergence on the mass, linear and angular momentum of the outflow is achieved though with 8 cells.

All these protostellar outflow sub-grid models are currently widely used in the community for astrophysical applications on the star formation rate/efficiency (e.g., Federrath, 2015; Kuiper et al., 2016; Offner and Chaban, 2017; Murray et al., 2018), stellar initial mass function (Krumholz et al., 2012), origin of turbulence driving (Nakamura and Li, 2011; Offner and Liu, 2018), as well as cluster and massive star formation (Wang et al., 2010; Cunningham et al., 2011; Kuiper et al., 2015; Li et al., 2018).

In the context of massive star formation, another class of outflows models are based on the luminosity and ionization. Sub-grids models of outflows generated by massive star luminosity have been developed (e.g., Krumholz et al., 2005; Peters et al., 2014), but the physical mechanisms driving this type of outflows are not yet well understood. No work in the context of massive star formation has accounted for the combined effects of magnetic fields and radiative force with a resolution sufficient to resolve the magnetic outflow launch. There is certainly a combination of magneto-centrifugal process and from radiative pressure, depending on the mass of the forming protostar.

7.3. Second Collapse and Pre-main Sequences

Ideally, the sub-grid models should be designed to reproduce the results of high resolution studies, which resolve the protostar scales. The short time integration after stellar core formation, dominated by the adiabatic evolution, does not allow to design sub-grid models for protostellar feedback from full 3D RMHD calculations with a high level of confidence. For these reason, all recent works mostly rely on the results of 1D spherical symmetry results to set their protostellar core evolution models. For instance, Kuiper and Yorke (2013) use the pre-main sequence evolution code STELLAR (Bodenheimer et al., 2007; Beuther et al., 2008) to determine the protostellar luminosity of the non-resolved protostars. Pelupessy et al. (2013) developed the open source Astrophysical Multi-purpose Software Environment (AMUSE), a component library of simulations involving different physical domains and scales. For instance, it enables one to couple hydrodynamical codes (GADGET, ATHENA, AMRVAC) with pre-main sequence and binary evolution codes such as MESA (Paxton et al., 2011). Recently, Wall et al. (2019) used AMUSE to perform a star cluster formation model that includes individual star formation from self-gravitating, magnetized gas, coupled to collisional stellar dynamics. It couples different tools: the AMR code FLASH (Fryxell et al., 2000), the N-body code ph4 (Capuzzo-Dolcetta et al., 2012), and the stellar evolution code SeBa (Portegies Zwart and Verbunt, 1996). We anticipate to see more works using this strategy of coupling tools and scales in the near future.

Recent work from Vaytet et al. (2018) have managed to resolve the second core formation with non-ideal MHD (Ohm diffusion and ambiploar diffusion) from collapse 1 M_⊙ dense core (see Figure 3). Interestingly, they show that the magnetic field lines topology is very different at the first and second core scales. Ambipolar diffusion is dominating at the first core scale. The magnetic fields direction is essentially vertical and the magnetic fields and the gas evolution are decoupled. At the second core scale, the ionization rises due to high temperature and the system is back to ideal MHD, but with a high plasma beta β = P_therm/P_mag. The weak magnetic field is efficiently wrapped up by the gas rotation and the resulting fields direction is essentially toroidal. On the opposite, when ideal MHD is preserved all the way from the dense core to the stellar core, the magnetic fields lines are strongly pinched and the field remains poloidal at the second core scale. The energy balance at the first and second core accretion shocks is also different. While the accretion shock at the first core has been classified as a radiative supercritical shock, with all the incident kinetic energy radiated away (or often referred to as cold accretion), in 1D spherical core collapse studies (Commerçon et al., 2011a; Vaytet et al., 2012), recent 3D works have shown that the first core accretion experiences both cold and hot accretion at the same time at its surface, depending on the accretion flows morphology. On the surface of the stellar core, the accretion shock has been found to be subcritical (Vaytet et al., 2012; Vaytet et al., 2018; Tomida et al., 2013) with a significant part of the incident kinetic energy transferred to the protostars. Although these results have been investigated only for the very first stages of the protostar evolution (only 1 month in Vaytet et al., 2018), they indicate that it is not easy to derive general properties of small scales phenomena to design robust sub-grid models. In particular, Wurster et al. (2018) show that the magnetic fields evolution of the protostars they formed in their SPH models is largely affected by numerical diffusion. Last, the radiative efficiency of the stellar accretion shock has been shown to be a key process for the pre-main sequence evolution and sets the radius of the young protostars (Baraffe et al., 2012). We recall that the sub-grids models designed for the radiative and dynamical feedback take the protostellar radius as an input. The launching speed of the protostellar jets as well as the accretion luminosity depends on $R_{\star }^{-1}$.

FIGURE 3

Figure 3. Volume rendering from collapsing 1 M_⊙ dense cores at the protostellar core scale. The left panel shows the results of an ideal MHD 3D simulations, and the right panel the comparison with non-ideal MHD including ambipolar diffusion and Ohmic diffusion. The color lines indicate the 3D magnetic field lines, and the color coding the magnetic field amplitude. Figure reproduced from Vaytet et al. (2018) with permission of A&A.

8. Discussion

We have described almost all the numerical methods that are required to model the star formation process within the turbulent interstellar medium. These discretised equations must now be solved for a given set of initial and boundary conditions, over some prescribed time. This defines the numerical set up and the overall simulation strategy. Obviously, given the wide range of spatial and temporal scales involved, one cannot realistically model an entire galaxy all the way down to the formation of brown dwarves and proto-planetary disks. Unfortunately, given the complexity of the star formation process, and the fact that large and small scales are strongly coupled, this is required by the physics of the problem. Indeed, supersonic turbulence is seeded on large scales, where kinetic energy is injected through shearing or colliding flows induced by galactic rotation, spiral waves and extra-galactic accretion flows. Moreover, turbulence and gas ejection can be triggered by stellar feedback processes occurring on very small scales, such as radiation driven flows, jets and ultimately supernovae explosions.

Modeling star formation is therefore fundamentally a multi-scale, multi-physics problem, and a very difficult one. Past or present day simulations always rely on some compromises. We can decompose them into several categories, corresponding roughly to the adopted box size and initial conditions. First, full galaxy simulations, where spiral waves and spiral shocks are used as the seed for turbulence. The box size must be of several tens of kiloparsecs to contain the entire rotating disk and possibly the galactic corona. The gas is compressed through the spiral shocks or various colliding flows. Thermal and gravitational instabilities fragment the gas into small molecular clouds. The spatial resolution never really drops below 0.1 or 1 pc in these simulations, so that stars cannot be modeled individually, at least for massive galaxies like the Milky Way (see Figure 4 for illustration). Simulators rely on simplified, sub-grid star formation recipes based on stochastically spawning star cluster particles of the same mass, or on the sink particle technique where in this case sink particles represent a star cluster, rather than a single star. Bournaud et al. (2010) studied the properties the ISM substructure and turbulence in galaxy simulations with resolutions up to 0.8 pc and 5 × 10³ M_⊙ with RAMSES. Renaud et al. (2013) and Kraljic et al. (2014) were able to capture the transition from turbulence-supported to self-gravitating gas with resolution up to 0.05 pc in simulations of Milky Way-like galaxies using RAMSES. Hopkins et al. (2012) presented SPH simulations dwarfs galaxies and Milky Way (MW) analogs to massive star-forming galactic discs with pc-scale resolution. Later, Hopkins et al. (2014) wrapped their implementation of feedback in the FIRE (Feedback In Realistic Environments) simulations suite, aiming, among others objectives, at resolving the formation of giant molecular clouds and the multi-phase interstellar medium (ISM). More recently, Hopkins et al. (2018) updated the FIRE implementation in the GIZMO code and performed sub-parsec resolution simulations of galactic discs. We also refer readers to the work of Dobbs et al. (2018) using SPH, Tasker (2011) using ENZO, and Smith et al. (2018) using AREPO. In addition, we note that for dwarf galaxies, it becomes possible to model much smaller scales so that at least massive stars can be modeled individually: Hu et al. (2016); Hu (2019) using SPH and a resolution of 1 − 5 M_⊙ by gas particle, Emerick et al. (2018) and Emerick et al. (2019) using ENZO with a maximum resolution of 1.8 pc.

FIGURE 4

Figure 4. Simulation of the Milky Way: column density of the gas disc in a sub-parsec resolution simulation. The color table only applies to the main panel: the table has been changed in each zoom-in view to enhance contrast. Figure adapted from Renaud et al. (2013) with permission from the authors.

The second approach relies on simulating only a small portion of vertically stratified galactic disks. The box size is usually around 1 kilo-parsec, with periodic boundary conditions in the direction of the disk plane and outflow conditions in the direction perpendicular to the disk plane. The geometry is clearly heavily constrained by this elongated, vertical and stratified layer of gas but it captures most of the phenomenon at work locally in the disk. Using the shearing box technique, one can also add more realistic shearing conditions, so that turbulence can be maintained both by stellar feedback and a large scale galactic shearing flow (Kim et al., 2002). In this case, the resolution can drop significantly below 0.1 parsec, may be 0.01 parsec or even less for the highest resolution simulations. Unfortunately, it is still impossible to model individual stars at this resolution so that here again most papers are based on the sink particle techniques for representing star clusters. It is however almost possible to resolve massive molecular cores (and their associated sink particles) individually, so that individual supernovae and HII regions can be resolved. These simulations are probably the most realistic models of the interstellar medium of the Galaxy, although still far from resolving the entire stellar population. In standard setups, the models account for the thermal instability and feedback (mechanical by SN and radiative by massive stars). First grid-based models were presented in Korpi et al. (1999) and included magnetic fields, SN heating and radiative cooling. 3D MHD AMR models of a kpc-scale ISM driven by SNe were presented in de Avillez and Breitschwerdt (2005) and Joung and Mac Low (2006). Currently, numerous studies and projects are based on a similar setup, with increasing physics put on over several years (SNe, radiative feedback, cosmic-rays, shear, etc.): the series of papers by Hennebelle and Iffrig (2014), Iffrig and Hennebelle (2017), and Colling et al. (2018) as well as the work of Martizzi et al. (2016) and Butler et al. (2015, 2017) using RAMSES, by Kim et al. (2011), Kim et al. (2013), Kim and Ostriker (2015), Kim and Ostriker (2017), and Kim and Ostriker (2018) using ATHENA, the SILCC project papers using FLASH (Walch et al., 2015; Girichidis et al., 2016, 2018; Gatto et al., 2017; Peters et al., 2017).

The third approach aims at simulating individual molecular cores. The box size ranges from 100 parsec down to 0.1 parsec for the highest resolution simulations. The mass typically ranges from 50 to a few 10³ M_⊙. Boundary conditions are either periodic or isolated. In the latter case, the cloud as a whole can collapse, while the former set of simulations represents a more uniform background with collapsing regions only at very small scales. The smallest box sizes represents internal regions of larger clouds and are the only simulations for which the entire star population, down to the brown dwarf limit, can be modeled. These simulations are directly tackling the origin of the stellar IMF, and the formation of proto-planetary disks. We mention in the following a non-exhaustive list of the work done in this intense field of research. The series of papers by Bate & collaborators have investigated the collapse of molecular clouds of mass ranging from 50 to 500 M_⊙ with increasing physics and numerical resolution since their pioneer work (Bate et al., 2003) using SPH. In particular, they have investigated the effect of radiative transfer (Bate, 2009, 2012) and magnetic fields (Price and Bate, 2008; Price and Bate, 2009). The work of Dale et al. (2007), Dale and Bonnell (2011), and Dale et al. (2012) have also used SPH simulations to study the effect of ionizing radiation from the forming massive stars. AMR codes have also been widely used in this context, accounting for a lot of physics and initial and boundary conditions. We note the work of Girichidis et al. (2011), Girichidis et al. (2012b), and Girichidis et al. (2012a) using FLASH who studied the importance of isolated initial conditions in isothermal cloud core collapse (without stellar feedback). Federrath et al. (2014b) also used FLASH to perform simulations of isolated cores with protostellar feedback (jets). Krumholz et al. (2007a), Krumholz et al. (2012), Hansen et al. (2012), Myers et al. (2013), Myers et al. (2014), Cunningham et al. (2011), Li et al. (2018), and Offner and Chaban (2017) presented 3D collapse models using ORION with radiative transfer and/or ideal MHD, as well as protostellar feedback (luminosity, outflows). RAMSES has also been extensively used in this approach: Hennebelle et al. (2011) and Commerçon et al. (2011b) for 100 M_⊙ isolated core collapse models with (R)MHD, Lee and Hennebelle (2018a), Lee and Hennebelle (2018b), and Lee and Hennebelle (2018c) focused on the peak of the IMF using (M)HD 1000 M_⊙ isolated core collapse models, Geen et al. (2015), Geen et al. (2016), and Gavagnin et al. (2017) studied the effect of ionizing radiation using the M1 moment method for radiative transfer. Last, we note the work of Wang et al. (2010) using ENZO with isolated boundary conditions, ideal MHD and protostellar jets. This non-exhaustive list of works done using this approach demonstrates the importance and the utility of such models. These simulations, however, rarely follow the formation of the second Larson core. They still rely on the sink particle technique, and on sub-grid models as well, when it comes to modeling proto-stellar jet launching, energy budget at accretion shocks and their associated radiation.

The last category of simulations deals with isolated collapsing low-mass dense cores, with the goal of following the formation of protostellar disks and ultimately the entire second collapse without relying on any sub-grid model. This ambitious strategy has already been discussed at length in the previous sections, so we will not describe it here again. We also refer the reader to the review by Wurster and Li (2018) on the formation of protostellar disks.

An interesting intermediate strategy has been developed recently with the zoom-in technique. This method is used now routinely in galaxy formation simulations, for which only a single dark matter halo can be re-simulated at much higher resolution, keeping the entire cosmological environment alive but at a coarser resolution. Here, the idea is to model the turbulent molecular cloud on large scales, say with a box of size 100 parsecs, to identify a single collapsing core, and to zoom-in on the core with much higher resolution and follow the first and possibly the second collapse and the formation of the proto-planetary disk. For instance, Kuffmeier et al. (2017) and Kuffmeier et al. (2018) presented RAMSES AMR zoom-in calculations from an outer scale of 40 pc down to cell sizes of 2 au to study the effect of the environment on the formation and the evolution of protoplanetary disks. Zhang et al. (2018) performed high resolution calculations of the formation and evolution of a star-forming core, obtained by running larger scale calculations of molecular cloud formation (Zamora-Avilés et al., 2018). They used FLASH with a maximum resolution of 25 AU and were able to cover scales from 256 pc to 25 AU (see Figure 5). Mocz et al. (2017) performed moving mesh AREPO calculations from a box size of about 5 pc down to a maximum resolution of ≃ 4 AU (≃ 5 × 10⁻⁵ pc) to study the magnetic fields morphology from large to dense core scales. At larger scale, Padoan et al. (2017) performed a simulation of supernova-driven turbulence using RAMSES, with a box size of 250 pc and a maximum resolution of 0.0076 pc. At kpc scales, Hennebelle (2018) (FRIGG project with RAMSES, simulation of the formation of self-gravitating cores) and Seifried et al. (2017) (SLICC-Zoom project with FLASH, simulations of the formation of molecular clouds) presented zoom calculations of stratified Galactic disks down to resolution of 10⁻²−10⁻³ pc. Thanks to the steady increase in CPU power, we expect to see more and more work using the zoom-in technique in the coming years.

FIGURE 5

Figure 5. Simulation of protostellar disk formation within molecular clouds. The column density maps show the entire molecular (upper left panel) and successive zoom within a star-forming dense core. Figure adapted from Zhang et al. (2018) with permission from the authors.

These different strategies are all in different ways very ambitious and address different problems in the theory of star formation. In addition, a compromise needs to be found between physical realism and computational efficiency. On one hand, we need to include magnetic fields, radiation fields and complex chemistry, but on the other hand, we need many small resolution elements and many time steps. Moreover, in the context of high performance computing, modern supercomputers are very hard to use at full efficiency when complex grid geometries is used in conjunction with expensive and demanding algorithms. There is a clear tendency for simulation projects deployed on the largest supercomputers in the world to use simple grid geometry, like Cartesian meshes and periodic boxes, with a highly simplified physical model. These large scale simulations are however very interesting to explore statistical aspects and large inertial range for turbulent flows. For example, Federrath (2013) performed 4096³ hydrodynamic isothermal turbulence simulations using FLASH with periodic boundary condition for more than 40,000 time-steps on 32,768 CPU cores. In the framework of the magnetorotational turbulence, Fromang (2010); Ryan et al. (2017) performed high resolution simulations of isothermal shearing boxes using, respectively, ZEUS and a GPU version of RAMSES.

In conclusion, all simulations of star formation published in the literature, many of which have been discussed here, explore various corners of the numerical parameter space: resolution vs. box size, statistics vs. internal structure, physical realism vs. computational speed. It is interesting that at the smallest scales of interest here, a “first principle” approach is in principle possible and pursued by several groups. At larger scales, however, star formation simulations share the same kind of limitations as galaxy formation simulations, or climate models and weather forecasting simulations, namely a strong dependence of the results on small, unresolved scales. Although bigger computers with more efficient, higher order codes and more realistic models will certainly help shed light of the mysteries of star formation, a robust methodology to implement sub-grid models and couple them properly to the fluid equations still needs to be invented in our field. Various attempts have been proposed in the context of unresolved turbulence and star formation from analytical works (Krumholz and McKee, 2005; Hennebelle and Chabrier, 2011; Padoan and Nordlund, 2011; Federrath and Klessen, 2012) but the methodology, in the context of the full spectrum of required physical processes, is still at its infancy.

Author Contributions

RT has written the introduction, the ideal MHD part and the radiative transfer part, as well as the discussion. BC has written the non-ideal MHD part, the second collapse part, the sink particle part and the subgrid models part.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

We would like to warmly thank our two referees, Richard Klein and Christoph Federrath, for their help in reviewing the manuscript. We thank Matthew Bate, Florent Renaud, and Shangjia Zhang for providing adapted figures from their work.

Footnotes

1. ^Stagnation points are regions where the flow is brought to rest in the frame of the system under study. More generally, even if the velocity does not vanish, stagnation points are regions in the flow where the streamlines are converging and the flow becomes compressive.

References

Abel, T., and Wandelt, B. D. (2002). Adaptive ray tracing for radiative transfer around point sources. Month. Notices R. Astron. Soc. 330, L53–L56. doi: 10.1046/j.1365-8711.2002.05206.x

CrossRef Full Text | Google Scholar

Adams, F. C., Ruden, S. P., and Shu, F. H. (1989). Eccentric gravitational instabilities in nearly keplerian disks. Astrophys. J. 347:959. doi: 10.1086/168187

CrossRef Full Text | Google Scholar

Adams, F. C., and Shu, F. H. (1986). Infrared spectra of rotating protostars. Astrophys. J. 308:836. doi: 10.1086/164555

CrossRef Full Text | Google Scholar

Alexiades, V., Amiez, G., and Gremaud, P. A. (1996). Super-time-stepping acceleration of explicit schemes for parabolic problems. Commun. Numer. Methods Eng. 12, 31–42.

Google Scholar

Ali, A., Harries, T. J., and Douglas, T. A. (2018). Modelling massive star feedback with Monte Carlo radiation hydrodynamics: photoionization and radiation pressure in a turbulent cloud. Month. Notices R. Astron. Soc. 477, 5422–5436. doi: 10.1093/mnras/sty1001

CrossRef Full Text | Google Scholar

Altay, G., Croft, R. A. C., and Pelupessy, I. (2008). SPHRAY: a smoothed particle hydrodynamics ray tracer for radiative transfer. Month. Notices R. Astron. Soc. 386, 1931–1946. doi: 10.1111/j.1365-2966.2008.13212.x

CrossRef Full Text | Google Scholar

Andre, P., Ward-Thompson, D., and Barsony, M. (1993). Submillimeter continuum observations of rho ophiuchi A: the candidate protostar VLA 1623 and prestellar clumps. Astrophys. J. 406:122. doi: 10.1086/172425

CrossRef Full Text | Google Scholar

Aubert, D., and Teyssier, R. (2008). A radiative transfer scheme for cosmological reionization based on a local Eddington tensor. Month. Notices R. Astron. Soc. 387, 295–307. doi: 10.1111/j.1365-2966.2008.13223.x

CrossRef Full Text | Google Scholar

Aubert, D., and Teyssier, R. (2010). Reionization simulations powered by graphics processing units. I. On the structure of the ultraviolet radiation field. Astrophys. J. 724, 244–266. doi: 10.1088/0004-637X/724/1/244

CrossRef Full Text | Google Scholar

Bai, X.-N. (2014). Hall-effect-Controlled Gas Dynamics in Protoplanetary Disks. I. Wind Solutions at the Inner Disk. Astrophys. J. 791:137. doi: 10.1088/0004-637X/791/2/137

CrossRef Full Text | Google Scholar

Bai, X.-N., and Stone, J. M. (2010). Particle-gas dynamics with athena: method and convergence. Astrophys. J. Suppl. Ser. 190, 297–310. doi: 10.1088/0067-0049/190/2/297

CrossRef Full Text | Google Scholar

Bai, X.-N., and Stone, J. M. (2011). Effect of ambipolar diffusion on the nonlinear evolution of magnetorotational instability in weakly ionized disks. Astrophys. J. 736:144. doi: 10.1088/0004-637X/736/2/144

CrossRef Full Text | Google Scholar

Balbus, S. A. (2009). Magnetohydrodynamics of protostellar disks. arXiv:0906.0854.

Google Scholar

Balbus, S. A., and Hawley, J. F. (1991). A powerful local shear instability in weakly magnetized disks. I - Linear analysis. II - Nonlinear evolution. Astrophys. J. 376, 214–233. doi: 10.1086/170270

CrossRef Full Text | Google Scholar

Balbus, S. A., and Terquem, C. (2001). Linear analysis of the hall effect in protostellar disks. Astrophys. J. 552, 235–247. doi: 10.1086/320452

CrossRef Full Text | Google Scholar

Balsara, D. S. (1998). Total variation diminishing scheme for adiabatic and isothermal magnetohydrodynamics. Astrophys. J. Suppl. Ser. 116, 133–153. doi: 10.1086/313093

CrossRef Full Text | Google Scholar

Balsara, D. S. (2001). Divergence-free adaptive mesh refinement for magnetohydrodynamics. J. Comput. Phys. 174, 614–648. doi: 10.1006/jcph.2001.6917

CrossRef Full Text | Google Scholar

Balsara, D. S. (2004). Second-order-accurate schemes for magnetohydrodynamics with divergence- free reconstruction. Astrophys. J. Suppl. Ser. 151, 149–184. doi: 10.1086/381377

CrossRef Full Text | Google Scholar

Balsara, D. S., Dumbser, M., and Abgrall, R. (2014). Multidimensional HLLC Riemann solver for unstructured meshes - With application to Euler and MHD flows. J. Comput. Phys. 261, 172–208. doi: 10.1016/j.jcp.2013.12.029

CrossRef Full Text | Google Scholar

Balsara, D. S., and Spicer, D. S. (1999). A staggered mesh algorithm using high order godunov fluxes to ensure solenoidal magnetic fields in magnetohydrodynamic simulations. J. Comput. Phys. 149, 270–292. doi: 10.1006/jcph.1998.6153

CrossRef Full Text | Google Scholar

Banerjee, R., and Pudritz, R. E. (2006). Outflows and jets from collapsing magnetized cloud cores. Astrophys. J. 641, 949–960. doi: 10.1086/500496

CrossRef Full Text | Google Scholar

Baraffe, I., Vorobyov, E., and Chabrier, G. (2012). Observed luminosity spread in young clusters and FU Ori stars: a unified picture. Astrophys. J. 756:118. doi: 10.1088/0004-637X/756/2/118

CrossRef Full Text | Google Scholar

Bate, M. R. (1998). Collapse of a molecular cloud core to stellar densities: the first three-dimensional calculations. Astrophys. J. 508, L95–L98. doi: 10.1086/311719

CrossRef Full Text | Google Scholar

Bate, M. R. (2009). The importance of radiative feedback for the stellar initial mass function. Month. Notices R. Astron. Soc. 392, 1363–1380. doi: 10.1111/j.1365-2966.2008.14165.x

CrossRef Full Text | Google Scholar

Bate, M. R. (2010). Collapse of a molecular cloud core to stellar densities: the radiative impact of stellar core formation on the circumstellar disc. Month. Notices R. Astron. Soc. 404, 79–83. doi: 10.1111/j.1745-3933.2010.00839.x

CrossRef Full Text | Google Scholar

Bate, M. R. (2012). Stellar, brown dwarf and multiple star properties from a radiation hydrodynamical simulation of star cluster formation. Month. Notices R. Astron. Soc. 419, 3115–3146. doi: 10.1111/j.1365-2966.2011.19955.x

CrossRef Full Text | Google Scholar

Bate, M. R., Bonnell, I. A., and Bromm, V. (2003). The formation of a star cluster: predicting the properties of stars and brown dwarfs. Month. Notices R. Astron. Soc. 339, 577–599. doi: 10.1046/j.1365-8711.2003.06210.x

CrossRef Full Text | Google Scholar

Bate, M. R., Bonnell, I. A., and Price, N. M. (1995). Modelling accretion in protobinary systems. Month. Notices R. Astron. Soc. 277, 362–376. doi: 10.1093/mnras/277.2.362

CrossRef Full Text | Google Scholar

Bate, M. R., and Burkert, A. (1997). Resolution requirements for smoothed particle hydrodynamics calculations with self-gravity. Month. Notices R. Astron. Soc. 288, 1060–1072. doi: 10.1093/mnras/288.4.1060

CrossRef Full Text | Google Scholar

Bate, M. R., and Lorén-Aguilar, P. (2017). On the dynamics of dust during protostellar collapse. Month. Notices R. Astron. Soc. 465, 1089–1094. doi: 10.1093/mnras/stw2853

CrossRef Full Text | Google Scholar

Bate, M. R., Tricco, T. S., and Price, D. J. (2014). Collapse of a molecular cloud core to stellar densities: stellar-core and outflow formation in radiation magnetohydrodynamic simulations. Month. Notices R. Astron. Soc. 437, 77–95. doi: 10.1093/mnras/stt1865

CrossRef Full Text | Google Scholar

Beuther, H., Linz, H., and Henning, T. (2008). “Massive star formation: observations confront theory,” in ASP Conference Series, Vol. 387 (San Francisco, CA: Astronomical Society of the Pacific), 189.

Google Scholar

Bhandare, A., Kuiper, R., Henning, T., Fendt, C., Marleau, G.-D., and Kölligan, A. (2018). First core properties: from low- to high-mass star formation. arXiv:1807.06597. doi: 10.1051/0004-6361/201832635

CrossRef Full Text | Google Scholar

Black, D. C., and Bodenheimer, P. (1975). Evolution of rotating interstellar clouds. I - Numerical techniques. Astrophys. J. 199, 619–632. doi: 10.1086/153729

CrossRef Full Text | Google Scholar

Black, D. C., and Scott, E. H. (1982). A numerical study of the effects of ambipolar diffusion on the collapse of magnetic gas clouds. Astrophys. J. 263, 696–715. doi: 10.1086/160541

CrossRef Full Text | Google Scholar

Blandford, R. D., and Payne, D. G. (1982). Hydromagnetic flows from accretion discs and the production of radio jets. Month. Notices R. Astron. Soc. 199, 883–903. doi: 10.1093/mnras/199.4.883

CrossRef Full Text | Google Scholar

Bleuler, A., and Teyssier, R. (2014). Towards a more realistic sink particle algorithm for the RAMSES CODE. Month. Notices R. Astron. Soc. 445, 4015–4036. doi: 10.1093/mnras/stu2005

CrossRef Full Text | Google Scholar

Bodenheimer, P., Laughlin, G. P., Rózyczka, M., and Yorke, H. W. (eds). (2007). Astronomy and Astrophysics. Numerical Methods in Astrophysics: An Introduction. Boca Raton, FL: Taylor & Francis Group.

Google Scholar

Bodenheimer, P., Yorke, H. W., Rozyczka, M., and Tohline, J. E. (1990). The formation phase of the solar nebula. Astrophys. J. 355:651. doi: 10.1086/168798

CrossRef Full Text | Google Scholar

Bonafede, A., Dolag, K., Stasyszyn, F., Murante, G., and Borgani, S. (2011). A non-ideal magnetohydrodynamic GADGET: simulating massive galaxy clusters. Month. Notices R. Astron. Soc. 418, 2234–2250. doi: 10.1111/j.1365-2966.2011.19523.x

CrossRef Full Text | Google Scholar

Booth, R. A., Sijacki, D., and Clarke, C. J. (2015). Smoothed particle hydrodynamics simulations of gas and dust mixtures. Month. Notices R. Astron. Soc. 452, 3932–3947. doi: 10.1093/mnras/stv1486

CrossRef Full Text | Google Scholar

Boscheri, W., and Dumbser, M. (2014). A direct Arbitrary-Lagrangian-Eulerian ADER-WENO finite volume scheme on unstructured tetrahedral meshes for conservative and non-conservative hyperbolic systems in 3D. J. Comput. Phys. 275, 484–523. doi: 10.1016/j.jcp.2014.06.059

CrossRef Full Text | Google Scholar

Boss, A. P., and Black, D. C. (1982). Collapse of accreting, rotating, isothermal, interstellar clouds. Astrophys. J. 258, 270–279. doi: 10.1086/160077

CrossRef Full Text | Google Scholar

Boss, A. P., Fisher, R. T., Klein, R. I., and McKee, C. F. (2000). The jeans condition and collapsing molecular cloud cores: filaments or binaries? Astrophys. J. 528, 325–335. doi: 10.1086/308160

CrossRef Full Text | Google Scholar

Bournaud, F., Elmegreen, B. G., Teyssier, R., Block, D. L., and Puerari, I. (2010). ISM properties in hydrodynamic galaxy simulations: turbulence cascades, cloud formation, role of gravity and feedback. Month. Notices R. Astron. Soc. 409, 1088–1099. doi: 10.1111/j.1365-2966.2010.17370.x

CrossRef Full Text | Google Scholar

Brackbill, J. U., and Barnes, D. C. (1980). The effect of nonzero product of magnetic gradient and B on the numerical solution of the magnetohydrodynamic equations. J. Comput. Phys. 35, 426–430. doi: 10.1016/0021-9991(80)90079-0

CrossRef Full Text

Brandenburg, A., and Subramanian, K. (2005). Astrophysical magnetic fields and nonlinear dynamo theory. Phys. Rep. 417, 1–209. doi: 10.1016/j.physrep.2005.06.005

CrossRef Full Text | Google Scholar

Brandenburg, A., and Zweibel, E. G. (1994). The formation of sharp structures by ambipolar diffusion. Astrophys. J. 427:L91. doi: 10.1086/187372

CrossRef Full Text | Google Scholar

Bryan, G. L., Norman, M. L., O'Shea, B. W., Abel, T., Wise, J. H., Turk, M. J., et al. (2014). ENZO: an adaptive mesh refinement code for astrophysics. Astrophys. J. Suppl. 211:19. doi: 10.1088/0067-0049/211/2/19

CrossRef Full Text | Google Scholar

Buntemeyer, L., Banerjee, R., Peters, T., Klassen, M., and Pudritz, R. E. (2016). Radiation hydrodynamics using characteristics on adaptive decomposed domains for massively parallel star formation simulations. New Astron. 43, 49–69. doi: 10.1016/j.newast.2015.07.002

CrossRef Full Text | Google Scholar

Butler, M. J., Tan, J. C., Teyssier, R., Rosdahl, J., Van Loo, S., and Nickerson, S. (2017). Kiloparsec-scale simulations of star formation in disk galaxies. IV. regulation of galactic star formation rates by stellar feedback. Astrophys. J. 841:82. doi: 10.3847/1538-4357/aa7054

CrossRef Full Text | Google Scholar

Butler, M. J., Tan, J. C., and Van Loo, S. (2015). Kiloparsec-scale simulations of star formation in disk galaxies III. Structure and dynamics of filaments and clumps in giant molecular clouds. Astrophys. J. 805:1. doi: 10.1088/0004-637X/805/1/1

CrossRef Full Text | Google Scholar

Capuzzo-Dolcetta, R., Limongi, M., and Tornambè, A. (2012). “Advances in computational astrophysics: methods, tools, and outcome,” in ASP Conference Proceedings, Vol. 453 (San Francisco, CA: Astronomical Society of the Pacific), 129.

Google Scholar

Chen, C.-Y., and Ostriker, E. C. (2014). Formation of magnetized prestellar cores with ambipolar diffusion and turbulence. Astrophys. J. 785:69. doi: 10.1088/0004-637X/785/1/69

CrossRef Full Text | Google Scholar

Chen, J.-W., and Lin, M.-K. (2018). Dusty disc-planet interaction with dust-free simulations. Month. Notices R. Astron. Soc. 478, 2737–2752. doi: 10.1093/mnras/sty1166

CrossRef Full Text | Google Scholar

Chiaki, G., and Yoshida, N. (2015). Particle splitting in smoothed particle hydrodynamics based on Voronoi diagram. Month. Notices R. Astron. Soc. 451, 3955–3963. doi: 10.1093/mnras/stv1227

CrossRef Full Text | Google Scholar

Chiang, H.-F., Looney, L. W., and Tobin, J. J. (2012). The envelope and embedded disk around the class 0 protostar L1157-mm: dual-wavelength interferometric observations and modeling. Astrophys. J. 756:168. doi: 10.1088/0004-637X/756/2/168

CrossRef Full Text | Google Scholar

Choi, E., Kim, J., and Wiita, P. J. (2009). An explicit scheme for incorporating ambipolar diffusion in a magnetohydrodynamics code. Astrophys. J. 181, 413–420. doi: 10.1088/0067-0049/181/2/413

CrossRef Full Text | Google Scholar

Christie, D., Wu, B., and Tan, J. C. (2017). GMC collisions as triggers of star formation. IV. The role of ambipolar diffusion. Astrophys. J. 848:50. doi: 10.3847/1538-4357/aa8a99

CrossRef Full Text | Google Scholar

Ciardi, A., and Hennebelle, P. (2010). Outflows and mass accretion in collapsing dense cores with misaligned rotation axis and magnetic field. Month. Notices R. Astron. Soc. 409, L39–L43. doi: 10.1111/j.1745-3933.2010.00942.x

CrossRef Full Text | Google Scholar

Ciolek, G. E. L. E. C., and Roberge, W. A. G. R. (2002). Time-dependent, multifluid, magnetohydrodynamic shock waves with grain dynamics I. Formulation and numerical tests. Astrophys. J. 567, 947–961. doi: 10.1086/338591

CrossRef Full Text | Google Scholar

Clark, P. C., Glover, S. C. O., and Klessen, R. S. (2012). TreeCol: a novel approach to estimating column densities in astrophysical simulations. Month. Notices R. Astron. Soc. 420, 745–756. doi: 10.1111/j.1365-2966.2011.20087.x

CrossRef Full Text | Google Scholar

Colling, C., Hennebelle, P., Geen, S., Iffrig, O., and Bournaud, F. (2018). Impact of galactic shear and stellar feedback on star formation. Astron. Astrophys. 620:A21. doi: 10.1051/0004-6361/201833161

CrossRef Full Text | Google Scholar

Commerçon, B., Audit, E., Chabrier, G., and Chièze, J.-P. (2011a). Physical and radiative properties of the first-core accretion shock. Astron. Astrophys. 530:A13. doi: 10.1051/0004-6361/201016213

CrossRef Full Text | Google Scholar

Commercon, B., Debout, V., and Teyssier, R. (2014). A fast, robust, and simple implicit method for adaptive time-stepping on adaptive mesh-refinement grids. Astron. Astrophys. 563:A11. doi: 10.1051/0004-6361/201322858

CrossRef Full Text | Google Scholar

Commerçon, B., Hennebelle, P., Audit, E., Chabrier, G., and Teyssier, R. (2008). Protostellar collapse: a comparison between smoothed particle hydrodynamics and adaptative mesh refinement calculations. Astron. Astrophys. 482:371. doi: 10.1051/0004-6361:20078591

CrossRef Full Text | Google Scholar

Commerçon, B., Hennebelle, P., Audit, E., Chabrier, G., and Teyssier, R. (2010). Protostellar collapse: radiative and magnetic feedbacks on small-scale fragmentation. Astron. Astrophys. 510:L3. doi: 10.1051/0004-6361/200913597

CrossRef Full Text | Google Scholar

Commerçon, B., Hennebelle, P., and Henning, T. (2011b). Collapse of massive magnetized dense cores using radiation magnetohydrodynamics: early fragmentation inhibition. Astrophys. J. Lett. 742:L9. doi: 10.1088/2041-8205/742/1/L9

CrossRef Full Text | Google Scholar

Commerçon, B., Teyssier, R., Audit, E., Hennebelle, P., and Chabrier, G. (2011). Radiation hydrodynamics with adaptive mesh refinement and application to prestellar core collapse. I. Methods. Astron. Astrophys. 529:A35. doi: 10.1051/0004-6361/201015880

CrossRef Full Text | Google Scholar

Cunningham, A. J., Klein, R. I., Krumholz, M. R., and McKee, C. F. (2011). Radiation-hydrodynamic simulations of massive star formation with protostellar outflows. Astrophys. J. 740:107. doi: 10.1088/0004-637X/740/2/107

CrossRef Full Text | Google Scholar

Dai, W., and Woodward, P. R. (1994). Extension of the piecewise parabolic method to multidimensional ideal magnetohydrodynamics. J. Comput. Phys. 115, 485–514. doi: 10.1006/jcph.1994.1212

CrossRef Full Text | Google Scholar

Dale, J. E., and Bonnell, I. (2011). Ionizing feedback from massive stars in massive clusters: fake bubbles and untriggered star formation. Month. Notices R. Astron. Soc. 414, 321–328. doi: 10.1111/j.1365-2966.2011.18392.x

CrossRef Full Text | Google Scholar

Dale, J. E., Ercolano, B., and Bonnell, I. A. (2012). Ionizing feedback from massive stars in massive clusters - II. Disruption of bound clusters by photoionization. Month. Notices R. Astron. Soc. 424, 377–392. doi: 10.1111/j.1365-2966.2012.21205.x

CrossRef Full Text | Google Scholar

Dale, J. E., Ercolano, B., and Clarke, C. J. (2007). A new algorithm for modelling photoionizing radiation in smoothed particle hydrodynamics. Month. Notices R. Astron. Soc. 382, 1759–1767. doi: 10.1111/j.1365-2966.2007.12486.x

CrossRef Full Text | Google Scholar

Davis, S. W., Stone, J. M., and Jiang, Y.-F. (2012). A radiation transfer solver for athena using short characteristics. Astrophys. J. Suppl. 199:9. doi: 10.1088/0067-0049/199/1/9

CrossRef Full Text | Google Scholar

de Avillez, M. A., and Breitschwerdt, D. (2005). Global dynamical evolution of the ISM in star forming galaxies. I. High resolution 3D simulations: effect of the magnetic field. Astron. Astrophys. 436, 585–600. doi: 10.1051/0004-6361:20042146

CrossRef Full Text | Google Scholar

Dedner, A., Kemm, F., Kröner, D., Munz, C. D., Schnitzer, T., and Wesenberg, M. (2002). Hyperbolic divergence cleaning for the MHD equations. J. Comput. Phys. 175, 645–673. doi: 10.1006/jcph.2001.6961

CrossRef Full Text | Google Scholar

Derigs, D., Winters, A. R., Gassner, G. J., and Walch, S. (2016). A novel high-order, entropy stable, 3D AMR MHD solver with guaranteed positive pressure. J. Comput. Phys. 317, 223–256. doi: 10.1016/j.jcp.2016.04.048

CrossRef Full Text | Google Scholar

Derigs, D., Winters, A. R., Gassner, G. J., Walch, S., and Bohm, M. (2017). Ideal GLM-MHD: about the entropy consistent nine-wave magnetic field divergence diminishing ideal magnetohydrodynamics equations. arXiv, 420–467. doi: 10.1016/j.jcp.2018.03.002

CrossRef Full Text | Google Scholar

Dobbs, C. L., Pettitt, A. R., Corbelli, E., and Pringle, J. E. (2018). Simulations of the flocculent spiral M33: what drives the spiral structure? Month. Notices R. Astron. Soc. 478, 3793–3808. doi: 10.1093/mnras/sty1231

CrossRef Full Text | Google Scholar

Dolag, K., Bartelmann, M., and Lesch, H. (1999). SPH simulations of magnetic fields in galaxy clusters. Astron. Astrophys. 348, 351–363.

Google Scholar

Dolag, K., and Stasyszyn, F. (2009). An MHD GADGET for cosmological simulations. Month. Notices R. Astron. Soc. 398, 1678–1697. doi: 10.1111/j.1365-2966.2009.15181.x

CrossRef Full Text | Google Scholar

Dorfi, E. (1982). 3D models for self-gravitating, rotating magnetic interstellar clouds. Astron. Astrophys. 114, 151–164.

Google Scholar

Draine, B. T. (1986). Multicomponent, reacting MHD flows. Month. Notices R. Astron. Soc. 220, 133–148. doi: 10.1093/mnras/220.1.133

CrossRef Full Text | Google Scholar

Draine, B. T., Roberge, W. G., and Dalgarno, A. (1983). Magnetohydrodynamic shock waves in molecular clouds. Astrophys. J. 264, 485–507. doi: 10.1086/160617

CrossRef Full Text | Google Scholar

Draine, B. T., and Sutin, B. (1987). Collisional charging of interstellar grains. Astrophys. J. 320, 803–817. doi: 10.1086/165596

CrossRef Full Text | Google Scholar

Dubroca, B., and Feugeas, J. (1999). Etude théorique et numérique d'une hiérarchie de modèles aux moments pour le transfert radiatif. Academie des Sciences Paris Comptes Rendus Serie Sciences Mathematiques 329, 915–920. doi: 10.1016/S0764-4442(00)87499-6

CrossRef Full Text | Google Scholar

Duffin, D. F., and Pudritz, R. E. (2008). Simulating hydromagnetic processes in star formation: introducing ambipolar diffusion into an adaptive mesh refinement code. Month. Notices R. Astron. Soc. 391, 1659–1673. doi: 10.1111/j.1365-2966.2008.14026.x

CrossRef Full Text | Google Scholar

Dzyurkevich, N., Commerçon, B., Lesaffre, P., and Semenov, D. (2017). Magnetic diffusivities in 3D radiative chemo-hydrodynamic simulations of protostellar collapse. Astron. Astrophys. 603:A105. doi: 10.1051/0004-6361/201628995

CrossRef Full Text | Google Scholar

Emerick, A., Bryan, G. L., and Mac Low, M.-M. (2018). Stellar radiation is critical for regulating star formation and driving outflows in low-mass dwarf galaxies. Astrophys. J. 865:L22. doi: 10.3847/2041-8213/aae315

CrossRef Full Text | Google Scholar

Emerick, A., Bryan, G. L., and Mac Low, M.-M. (2019). Simulating an isolated dwarf galaxy with multichannel feedback and chemical yields from individual stars. Month. Notices R. Astron. Soc. 482, 1304–1329. doi: 10.1093/mnras/sty2689

CrossRef Full Text | Google Scholar

Epstein, P. S. (1924). On the resistance experienced by spheres in their motion through gases. Phys. Rev. 23, 710–733. doi: 10.1103/PhysRev.23.710

CrossRef Full Text | Google Scholar

Ercolano, B., Barlow, M. J., Storey, P. J., and Liu, X. W. (2003). MOCASSIN: a fully three-dimensional Monte Carlo photoionization code. Month. Notice R. Astron. Soc. 340, 1136–1152. doi: 10.1046/j.1365-8711.2003.06371.x

CrossRef Full Text | Google Scholar

Evans, C. R., and Hawley, J. F. (1988). Simulation of magnetohydrodynamic flows - A constrained transport method. Astrophys. J. 332, 659–677. doi: 10.1086/166684

CrossRef Full Text | Google Scholar

Falle, S. (2003). A numerical scheme for multifluid magnetohydrodynamics. Month. Notices R. Astron. Soc. 344, 1210–1218. doi: 10.1046/j.1365-8711.2003.06908.x

CrossRef Full Text | Google Scholar

Falle, S. A. E. G. (2002). Rarefaction shocks, shock errors, and low order of accuracy in ZEUS. Astrophys. J. 577, L123–L126. doi: 10.1086/344336

CrossRef Full Text | Google Scholar

Federrath, C. (2013). On the universality of supersonic turbulence. Month. Notices R. Astron. Soc. 436, 1245–1257. doi: 10.1093/mnras/stt1644

CrossRef Full Text | Google Scholar

Federrath, C. (2015). Inefficient star formation through turbulence, magnetic fields and feedback. Month. Notices R. Astron. Soc. 450, 4035–4042. doi: 10.1093/mnras/stv941

CrossRef Full Text | Google Scholar

Federrath, C. (2016). Magnetic field amplification in turbulent astrophysical plasmas. J. Plasma Phys. 82:535820601. doi: 10.1017/S0022377816001069

CrossRef Full Text | Google Scholar

Federrath, C., Banerjee, R., Clark, P., and Klessen, R. (2010). Modeling collapse and accretion in turbulent gas clouds: implementation and comparison of sink particles in AMR and SPH. Astrophys. J. 713, 269–290. doi: 10.1088/0004-637X/713/1/269

CrossRef Full Text | Google Scholar

Federrath, C., Banerjee, R., Seifried, D., Clark, P. C., and Klessen, R. S. (2011a). “Implementing and comparing sink particles in AMR and SPH,” in Computational Star Formation, Proceedings of the International Astronomical Union, IAU Symposium, vol. 270, 425–428. doi: 10.1017/S1743921311000755

CrossRef Full Text | Google Scholar

Federrath, C., Chabrier, G., Schober, J., Banerjee, R., Klessen, R. S., and Schleicher, D. R. (2011b). Mach number dependence of turbulent magnetic field amplification: solenoidal versus compressive flows. Phys. Rev. Lett. 107:114504. doi: 10.1103/PhysRevLett.107.114504

PubMed Abstract | CrossRef Full Text | Google Scholar

Federrath, C., and Klessen, R. S. (2012). The star formation rate of turbulent magnetized clouds: comparing theory, simulations, and observations. Astrophys. J. 761:156. doi: 10.1088/0004-637X/761/2/156

CrossRef Full Text | Google Scholar

Federrath, C., Rathborne, J. M., Longmore, S. N., Kruijssen, J. M. D., Bally, J., Contreras, Y., et al. (2017). The link between solenoidal turbulence and slow star formation in G0.253+0.016. Multi-Messenger Astrophys. Galactic Centre 322, 123–128. doi: 10.1017/S1743921316012357

CrossRef Full Text | Google Scholar

Federrath, C., Schober, J., Bovino, S., and Schleicher, D. R. G. (2014a). The turbulent dynamo in highly compressible supersonic plasmas. Astrophys. J. Lett. 797:L19. doi: 10.1088/2041-8205/797/2/L19

CrossRef Full Text | Google Scholar

Federrath, C., Schrön, M., Banerjee, R., and Klessen, R. S. (2014b). Modeling jet and outflow feedback during star cluster formation. Astrophys. J. 790:128. doi: 10.1088/0004-637X/790/2/128

CrossRef Full Text | Google Scholar

Federrath, C., Sur, S., Schleicher, D. R. G., Banerjee, R., and Klessen, R. S. (2011c). A new jeans resolution criterion for (M)HD simulations of self-gravitating gas: application to magnetic field amplification by gravity-driven turbulence. Astrophys. J. 731:62. doi: 10.1088/0004-637X/731/1/62

CrossRef Full Text | Google Scholar

Felker, K. G., and Stone, J. M. (2018). A fourth-order accurate finite volume method for ideal MHD via upwind constrained transport. J. Comput. Phys. 375, 1365–1400. doi: 10.1016/j.jcp.2018.08.025

CrossRef Full Text | Google Scholar

Ferreira, J. (1997). Magnetically-driven jets from Keplerian accretion discs. Astron. Astrophys. 319, 340–359.

Google Scholar

Fontani, F., Commerçon, B., Giannetti, A., Beltrán, M. T., Sánchez-Monge, Á., Testi, L., et al. (2018). Fragmentation properties of massive protocluster gas clumps: an ALMA study. Astron. Astrophys. 615:A94. doi: 10.1051/0004-6361/201832672

CrossRef Full Text | Google Scholar

Forgan, D., Rice, K., Stamatellos, D., and Whitworth, A. (2009). Introducing a hybrid radiative transfer method for smoothed particle hydrodynamics. Month. Notices R. Astron. Soc. 394, 882–891. doi: 10.1111/j.1365-2966.2008.14373.x

CrossRef Full Text | Google Scholar

Fromang, S. (2010). MHD simulations of the magnetorotational instability in a shearing box with zero net flux: the case Pm = 4. Astron. Astrophys. 514:L5. doi: 10.1051/0004-6361/201014284

CrossRef Full Text | Google Scholar

Fromang, S., Hennebelle, P., and Teyssier, R. (2006). A high order Godunov scheme with constrained transport and adaptive mesh refinement for astrophysical magnetohydrodynamics. Astron. Astrophys. 457, 371–384. doi: 10.1051/0004-6361:20065371

CrossRef Full Text | Google Scholar

Fromang, S., Papaloizou, J., Lesur, G., and Heinemann, T. (2007). MHD simulations of the magnetorotational instability in a shearing box with zero net flux. Astron. Astrophys. 476, 1123–1132. doi: 10.1051/0004-6361:20077943

CrossRef Full Text | Google Scholar

Frostholm, T., Haugbølle, T., and Grassi, T. (2018). Lampray: multi-group long characteristics ray tracing for adaptive mesh radiation hydrodynamics. arXiv:1809.05541.

Google Scholar

Fryxell, B., Olson, K., Ricker, P., Timmes, F. X., Zingale, M., Lamb, D. Q., et al. (2000). FLASH: an adaptive mesh hydrodynamics code for modeling astrophysical thermonuclear flashes. Astrophys. J. Suppl. Ser. 131, 273–334. doi: 10.1086/317361

CrossRef Full Text | Google Scholar

Galli, D., Walmsley, M., and Gonçalves, J. (2002). The structure and stability of molecular cloud cores in external radiation fields. Astron. Astrophys. 394, 275–284. doi: 10.1051/0004-6361:20021125

CrossRef Full Text | Google Scholar

Gatto, A., Walch, S., Naab, T., Girichidis, P., Wünsch, R., Glover, S. C. O., et al. (2017). The SILCC project - III. Regulation of star formation and outflows by stellar winds and supernovae. Month. Notices R. Astron. Soc. 466, 1903–1924. doi: 10.1093/mnras/stw3209

CrossRef Full Text | Google Scholar

Gavagnin, E., Bleuler, A., Rosdahl, J., and Teyssier, R. (2017). Star cluster formation in a turbulent molecular cloud self-regulated by photoionization feedback. Month. Notices R. Astron. Soc. 472, 4155–4172. doi: 10.1093/mnras/stx2222

CrossRef Full Text | Google Scholar

Geen, S., Hennebelle, P., Tremblin, P., and Rosdahl, J. (2015). Photoionization feedback in a self-gravitating, magnetized, turbulent cloud. Month. Notices R. Astron. Soc. 454, 4484–4502. doi: 10.1093/mnras/stv2272

CrossRef Full Text | Google Scholar

Geen, S., Hennebelle, P., Tremblin, P., and Rosdahl, J. (2016). Feedback in Clouds II: UV photoionization and the first supernova in a massive cloud. Month. Notices R. Astron. Soc. 463, 3129–3142. doi: 10.1093/mnras/stw2235

CrossRef Full Text | Google Scholar

Geen, S., Watson, S. K., Rosdahl, J., Bieri, R., Klessen, R. S., and Hennebelle, P. (2018). On the indeterministic nature of star formation on the cloud scale. Month. Notices R. Astron. Soc. 481, 2548–2569. doi: 10.1093/mnras/sty2439

CrossRef Full Text | Google Scholar

Gerin, M., Pety, J., Fuente, A., Cernicharo, J., Commerçon, B., and Marcelino, N. (2015). Nascent bipolar outflows associated with the first hydrostatic core candidates Barnard 1b-N and 1b-S. Astron. Astrophys. 577:L2. doi: 10.1051/0004-6361/201525777

CrossRef Full Text | Google Scholar

Girichidis, P., Federrath, C., Allison, R., Banerjee, R., and Klessen, R. S. (2012a). Importance of the initial conditions for star formation - III. Statistical properties of embedded protostellar clusters. Month. Notices R. Astron. Soc. 420, 3264–3280. doi: 10.1111/j.1365-2966.2011.20250.x

CrossRef Full Text | Google Scholar

Girichidis, P., Federrath, C., Banerjee, R., and Klessen, R. S. (2011). Importance of the initial conditions for star formation - I. Cloud evolution and morphology. Month. Notices R. Astron. Soc. 413, 2741–2759. doi: 10.1111/j.1365-2966.2011.18348.x

CrossRef Full Text | Google Scholar

Girichidis, P., Federrath, C., Banerjee, R., and Klessen, R. S. (2012b). Importance of the initial conditions for star formation - II. Fragmentation-induced starvation and accretion shielding. Month. Notices R. Astron. Soc. 420, 613–626. doi: 10.1111/j.1365-2966.2011.20073.x

CrossRef Full Text | Google Scholar

Girichidis, P., Seifried, D., Naab, T., Peters, T., Walch, S., Wünsch, R., et al. (2018). The SILCC project - V. The impact of magnetic fields on the chemistry and the formation of molecular clouds. Month. Notices R. Astron. Soc. 480, 3511–3540. doi: 10.1093/mnras/sty2016

CrossRef Full Text | Google Scholar

Girichidis, P., Walch, S., Naab, T., Gatto, A., Wünsch, R., Glover, S. C. O., et al. (2016). The SILCC (SImulating the lifeCycle of molecular Clouds) project - II. Dynamical evolution of the supernova-driven ISM and the launching of outflows. Month. Notices R. Astron. Soc. 456, 3432–3455. doi: 10.1093/mnras/stv2742

CrossRef Full Text | Google Scholar

Gnedin, N. Y., and Abel, T. (2001). Multi-dimensional cosmological radiative transfer with a Variable Eddington Tensor formalism. New Astron. 6, 437–455. doi: 10.1016/S1384-1076(01)00068-9

CrossRef Full Text | Google Scholar

Gong, H., and Ostriker, E. C. (2013). Implementation of sink particles in the athena code. Astrophys. J. Suppl. Ser. 204:8. doi: 10.1088/0067-0049/204/1/8

CrossRef Full Text | Google Scholar

González, M., Audit, E., and Huynh, P. (2007). HERACLES: a three-dimensional radiation hydrodynamics code. Astron. Astrophys. 464, 429–435. doi: 10.1051/0004-6361:20065486

CrossRef Full Text | Google Scholar

González, M., Vaytet, N., Commercon, B., and Masson, J. (2015). Multigroup radiation hydrodynamics with flux-limited diffusion and adaptive mesh refinement. Astron. Astrophys. 578:A12. doi: 10.1051/0004-6361/201525971

CrossRef Full Text | Google Scholar

González-Morales, P., Khomenko, E., Downes, T., and De Vicente, A. (2018). MHDSTS: a new explicit numerical scheme for simulations of partially ionised solar plasma. Astron. Astrophys. 615:67. doi: 10.1051/0004-6361/201731916

CrossRef Full Text | Google Scholar

Goodwin, S. P., Whitworth, A. P., and Ward-Thompson, D. (2004). Astronomy & Astrophysics Simulating star formation in molecular cloud cores I. The influence of low levels of turbulence on fragmentation and multiplicity. Astron. Astrophys. 414, 633–650. doi: 10.1051/0004-6361:20031594

CrossRef Full Text | Google Scholar

Greif, T. H., Springel, V., White, S. D., Glover, S. C., Clark, P. C., Smith, R. J., et al. (2011). Simulations on a moving mesh: the clustered formation of population III protostars. Astrophys. J. 737:75. doi: 10.1088/0004-637X/737/2/75

CrossRef Full Text | Google Scholar

Guillet, T., Pakmor, R., Springel, V., Chandrashekar, P., and Klingenberg, C. (2018). High-order magnetohydrodynamics for astrophysics with an adaptive mesh refinement discontinuous galerkin scheme. arXiv. doi: 10.1093/mnras/stz314

CrossRef Full Text | Google Scholar

Hansen, C. E., Klein, R. I., McKee, C. F., and Fisher, R. T. (2012). Feedback effects on low-mass star formation. Astrophys. J. 747:22. doi: 10.1088/0004-637X/747/1/22

CrossRef Full Text | Google Scholar

Harries, T. J. (2011). An algorithm for Monte Carlo time-dependent radiation transfer. Month. Notices R. Astron. Soc. 416, 1500–1508. doi: 10.1111/j.1365-2966.2011.19147.x

CrossRef Full Text | Google Scholar

Harries, T. J. (2015). Radiation-hydrodynamical simulations of massive star formation using Monte Carlo radiative transfer - I. Algorithms and numerical methods. Month. Notices R. Astron. Soc. 448, 3156–3166. doi: 10.1093/mnras/stv158

CrossRef Full Text | Google Scholar

Hayes, J. C., Norman, M. L., Fiedler, R. A., Bordner, J. O., Li, P. S., Clark, S. E., et al. (2006). Simulating radiating and magnetized flows in multiple dimensions with ZEUS-MP. Astrophys. J. Suppl. Ser. 165, 188–228. doi: 10.1086/504594

CrossRef Full Text | Google Scholar

Hennebelle, P. (2018). The FRIGG project: from intermediate galactic scales to self-gravitating cores. Astron. Astrophys. 611:A24. doi: 10.1051/0004-6361/201731071

CrossRef Full Text | Google Scholar

Hennebelle, P., and Chabrier, G. (2011). Analytical star formation rate from gravoturbulent fragmentation. Astrophys. J. 743:L29. doi: 10.1088/2041-8205/743/2/L29

CrossRef Full Text | Google Scholar

Hennebelle, P., Commerçon, B., Chabrier, G., and Marchand, P. (2016). Magnetically self-regulated formation of early protoplanetary disks. Astrophys. J. Lett. 830:L8. doi: 10.3847/2041-8205/830/1/L8

CrossRef Full Text | Google Scholar

Hennebelle, P., Commerçon, B., Joos, M., Klessen, R. S., Krumholz, M., Tan, J. C., et al. (2011). Collapse, outflows and fragmentation of massive, turbulent and magnetized prestellar barotropic cores. Astron. Astrophys. 528:A72. doi: 10.1051/0004-6361/201016052

CrossRef Full Text | Google Scholar

Hennebelle, P., and Fromang, S. (2008). Magnetic processes in a collapsing dense core. I. Accretion and ejection. Astron. Astrophys. 477, 9–24. doi: 10.1051/0004-6361:20078309

CrossRef Full Text | Google Scholar

Hennebelle, P., and Iffrig, O. (2014). Simulations of magnetized multiphase galactic disc regulated by supernovae explosions. Astron. Astrophys. 570:A81. doi: 10.1051/0004-6361/201423392

CrossRef Full Text | Google Scholar

Hopkins, P. F. (2017). Anisotropic diffusion in mesh-free numerical magnetohydrodynamics. Month. Notices R. Astron. Soc. 466, 3387–3405. doi: 10.1093/mnras/stw3306

CrossRef Full Text | Google Scholar

Hopkins, P. F., Kereš, D., Oñorbe, J., Faucher-Giguère, C.-A., Quataert, E., Murray, N., et al. (2014). Galaxies on FIRE (feedback in realistic environments): stellar feedback explains cosmologically inefficient star formation. Month. Notices R. Astron. Soc. 445, 581–603. doi: 10.1093/mnras/stu1738

CrossRef Full Text | Google Scholar

Hopkins, P. F., and Lee, H. (2016). The fundamentally different dynamics of dust and gas in molecular clouds. Month. Notices R. Astron. Soc. 456, 4174–4190. doi: 10.1093/mnras/stv2745

CrossRef Full Text | Google Scholar

Hopkins, P. F., Quataert, E., and Murray, N. (2012). The structure of the interstellar medium of star-forming galaxies. Month. Notices R. Astron. Soc. 421, 3488–3521. doi: 10.1111/j.1365-2966.2012.20578.x

CrossRef Full Text | Google Scholar

Hopkins, P. F., Wetzel, A., Kereš, D., Faucher-Giguère, C.-A., Quataert, E., Boylan-Kolchin, M., et al. (2018). FIRE-2 simulations: physics versus numerics in galaxy formation. Month. Notices R. Astron. Soc. 480, 800–863. doi: 10.1093/mnras/sty1690

CrossRef Full Text | Google Scholar

Hosking, J. G., and Whitworth, a. P. (2004). Modelling ambipolar diffusion with two-fluid smoothed particle hydrodynamics. October 1000, 994–1000. doi: 10.1111/j.1365-2966.2004.07273.x

CrossRef Full Text | Google Scholar

Hu, C.-Y. (2019). Supernova-driven winds in simulated dwarf galaxies. Month. Notices R. Astron. Soc. 483, 3363–3381. doi: 10.1093/mnras/sty3252

CrossRef Full Text | Google Scholar

Hu, C.-Y., Naab, T., Walch, S., Glover, S. C. O., and Clark, P. C. (2016). Star formation and molecular hydrogen in dwarf galaxies: a non-equilibrium view. Month. Notices R. Astron. Soc. 458, 3528–3553. doi: 10.1093/mnras/stw544

CrossRef Full Text | Google Scholar

Hubber, D. A., Batty, C. P., McLeod, A., and Whitworth, A. P. (2011). SEREN - a new SPH code for star and planet formation simulations. Algorithms and tests. Astron. Astrophys. 529:A27. doi: 10.1051/0004-6361/201014949

CrossRef Full Text | Google Scholar

Hubber, D. A., Walch, S., and Whitworth, A. P. (2013). An improved sink particle algorithm for SPH simulations. Month. Notices R. Astron. Soc. 430, 3261–3275. doi: 10.1093/mnras/stt128

CrossRef Full Text | Google Scholar

Iffrig, O., and Hennebelle, P. (2017). Structure distribution and turbulence in self-consistently supernova-driven ISM of multiphase magnetized galactic discs. arXiv 70, 1–22. doi: 10.1051/0004-6361/201630290

CrossRef Full Text | Google Scholar

Inoue, T., and Inutsuka, S. (2008). Two-fluid magnetohydrodynamic simulations of converging H I flows in the interstellar medium. I. Methodology and basic results. Astrophys. J. 687, 303–310. doi: 10.1086/590528

CrossRef Full Text | Google Scholar

Jappsen, A. K., Klessen, R. S., Larson, R. B., Li, Y., and Mac Low, M. M. (2005). The stellar mass spectrum from non-isothermal gravoturbulent fragmentation. Astron. Astrophys. 435, 611–623. doi: 10.1051/0004-6361:20042178

CrossRef Full Text | Google Scholar

Jiang, G.-S., and Wu, C.-c. (1999). A high-order WENO finite difference scheme for the equations of ideal magnetohydrodynamics. J. Comput. Phys. 150, 561–594. doi: 10.1006/jcph.1999.6207

CrossRef Full Text | Google Scholar

Jiang, Y.-F., Stone, J. M., and Davis, S. W. (2012). A godunov method for multidimensional radiation magnetohydrodynamics based on a variable eddington tensor. Astrophys. J. Suppl. 199:14. doi: 10.1088/0067-0049/199/1/14

CrossRef Full Text | Google Scholar

Jones, M. O., and Bate, M. R. (2018). Sink particle radiative feedback in smoothed particle hydrodynamics models of star formation. Month. Notices R. Astron. Soc. 480, 2562–2577. doi: 10.1093/mnras/sty1969

CrossRef Full Text | Google Scholar

Joung, M. K. R., and Mac Low, M.-M. (2006). Turbulent structure of a stratified supernova-driven interstellar medium. Astrophys. J. 653, 1266–1279. doi: 10.1086/508795

CrossRef Full Text | Google Scholar

Kannan, R., Vogelsberger, M., Marinacci, F., McKinnon, R., Pakmor, R., and Springel, V. (2018). AREPO-RT: radiation hydrodynamics on a moving mesh. arXiv:1804.01987. doi: 10.1093/mnras/stz287

CrossRef Full Text | Google Scholar

Kim, C.-G., Kim, W.-T., and Ostriker, E. C. (2011). Regulation of star formation rates in multiphase galactic disks: numerical tests of the thermal/dynamical equilibrium model. Astrophys. J. 743:25. doi: 10.1088/0004-637X/743/1/25

CrossRef Full Text | Google Scholar

Kim, C.-G., and Ostriker, E. C. (2015). Vertical equilibrium, energetics, and star formation rates in magnetized galactic disks regulated by momentum feedback from supernovae. Astrophys. J. 815:67. doi: 10.1088/0004-637X/815/1/67

CrossRef Full Text | Google Scholar

Kim, C.-G., and Ostriker, E. C. (2017). Three-phase interstellar medium in galaxies resolving evolution with star formation and supernova feedback (TIGRESS): algorithms, fiducial model, and convergence. Astrophys. J. 846:133. doi: 10.3847/1538-4357/aa8599

CrossRef Full Text | Google Scholar

Kim, C.-G., and Ostriker, E. C. (2018). Numerical simulations of multiphase winds and fountains from star-forming galactic disks. I. Solar neighborhood TIGRESS model. Astrophys. J. 853:173. doi: 10.3847/1538-4357/aaa5ff

CrossRef Full Text | Google Scholar

Kim, C.-G., Ostriker, E. C., and Kim, W.-T. (2013). Three-dimensional hydrodynamic simulations of multiphase galactic disks with star formation Feedback. I. Regulation of star formation rates. Astrophys. J. 776:1. doi: 10.1088/0004-637X/776/1/1

CrossRef Full Text | Google Scholar

Kim, W.-T., Ostriker, E. C., and Stone, J. M. (2002). Three-dimensional simulations of parker, magneto-jeans, and swing instabilities in shearing galactic gas disks. Astrophys. J. 581, 1080–1100. doi: 10.1086/344367

CrossRef Full Text | Google Scholar

Kitsionas, S., and Whitworth, A. P. (2002). Smoothed particle hydrodynamics with particle splitting, applied to self-gravitating collapse. Month. Notices R. Astron. Soc. 330, 129–136. doi: 10.1046/j.1365-8711.2002.05115.x

CrossRef Full Text | Google Scholar

Klassen, M., Kuiper, R., Pudritz, R. E., Peters, T., Banerjee, R., and Buntemeyer, L. (2014). A general hybrid radiation transport scheme for star formation simulations on an adaptive grid. Astrophys. J. 797:4. doi: 10.1088/0004-637X/797/1/4

CrossRef Full Text | Google Scholar

Klassen, M., Pudritz, R., Kuiper, R., Peters, T., and Banerjee, R. (2016). Simulating the formation of massive protostars. I. Radiative feedback and accretion disks. Astrophys. J. 823:28. doi: 10.3847/0004-637X/823/1/28

CrossRef Full Text | Google Scholar

Klein, R. I. (1999). Star formation with 3-D adaptive mesh refinement: the collapse and fragmentation of molecular clouds. J. Comput. Appl. Math. 109, 123–152. doi: 10.1016/S0377-0427(99)00156-9

CrossRef Full Text | Google Scholar

Klingenberg, C., Pörner, F., and Xia, Y. (2017). An efficient implementation of the divergence free constraint in a discontinuous galerkin method for magnetohydrodynamics on unstructured meshes. Commun. Comput. Phys. 21, 423–442. doi: 10.4208/cicp.180515.230616a

CrossRef Full Text | Google Scholar

Koga, S., Tsukamoto, Y., Okuzumi, S., and Machida, M. N. (2019). Dependence of Hall coefficient on grain size and cosmic ray rate and implication for circumstellar disc formation. Month. Notices R. Astron. Soc. 484, 2119–2136. doi: 10.1093/mnras/sty3524

CrossRef Full Text | Google Scholar

Korpi, M. J., Brandenburg, A., Shukurov, A., Tuominen, I., and Nordlund, Å. (1999). A supernova-regulated interstellar medium: simulations of the turbulent multiphase medium. Astrophys. J. 514, L99–L102. doi: 10.1086/311954

CrossRef Full Text | Google Scholar

Kotarba, H., Lesch, H., Dolag, K., Naab, T., Johansson, P. H., and Stasyszyn, F. A. (2009). Magnetic field structure due to the global velocity field in spiral galaxies. Month. Notices R. Astron. Soc. 397, 733–747. doi: 10.1111/j.1365-2966.2009.15030.x

CrossRef Full Text | Google Scholar

Kraljic, K., Renaud, F., Bournaud, F., Combes, F., Elmegreen, B., Emsellem, E., et al. (2014). The role of turbulence in star formation laws and thresholds. Astrophys. J. 784:112. doi: 10.1088/0004-637X/784/2/112

CrossRef Full Text | Google Scholar

Krasnopolsky, R., Li, Z. Y., and Shang, H. (2011). Disk formation in magnetized clouds enabled by the hall effect. Astrophys. J. 733, 1–6. doi: 10.1088/0004-637X/733/1/54

CrossRef Full Text | Google Scholar

Krumholz, M. R. (2014). The big problems in star formation: the star formation rate, stellar clustering, and the initial mass function. Phys. Rep. 539, 49–134. doi: 10.1016/j.physrep.2014.02.001

CrossRef Full Text | Google Scholar

Krumholz, M. R., Klein, R. I., and McKee, C. F. (2007a). Radiation-hydrodynamic simulations of collapse and fragmentation in massive protostellar cores. Astrophys. J. 656, 959–979. doi: 10.1086/510664

CrossRef Full Text | Google Scholar

Krumholz, M. R., Klein, R. I., and McKee, C. F. (2012). Radiation-hydrodynamic simulations of the formation of orion-like star clusters. II. The initial mass function from winds, turbulence, and radiation. Astrophys. J. 754:71. doi: 10.1088/0004-637X/754/1/71

CrossRef Full Text | Google Scholar

Krumholz, M. R., Klein, R. I., McKee, C. F., and Bolstad, J. (2007b). Equations and algorithms for mixed-frame flux-limited diffusion radiation hydrodynamics. Astrophys. J. 667, 626–643. doi: 10.1086/520791

CrossRef Full Text | Google Scholar

Krumholz, M. R., and McKee, C. F. (2005). A general theory of turbulence-regulated star formation, from spirals to ultraluminous infrared galaxies. Astrophys. J. 630, 250–268. doi: 10.1086/431734

CrossRef Full Text | Google Scholar

Krumholz, M. R., McKee, C. F., and Klein, R. I. (2004). Embedding lagrangian sink particles in eulerian grids. Astrophys. J. 611, 399–412. doi: 10.1086/421935

CrossRef Full Text | Google Scholar

Krumholz, M. R., McKee, C. F., and Klein, R. I. (2005). How protostellar outflows help massive stars form. Astrophys. J. 618, L33–L36. doi: 10.1086/427555

CrossRef Full Text | Google Scholar

Kuffmeier, M., Frimann, S., Jensen, S. S., and Haugbølle, T. (2018). Episodic accretion: the interplay of infall and disc instabilities. Month. Notices R. Astron. Soc. 475, 2642–2658. doi: 10.1093/mnras/sty024

CrossRef Full Text | Google Scholar

Kuffmeier, M., Haugbølle, T., and Nordlund, Å. (2017). Zoom-in simulations of protoplanetary disks starting from GMC scales. Astrophys. J. 846:7. doi: 10.3847/1538-4357/aa7c64

CrossRef Full Text | Google Scholar

Kuiper, R., and Hosokawa, T. (2018). First hydrodynamics simulations of radiation forces and photoionization feedback in massive star formation. Astron. Astrophys. 616:A101. doi: 10.1051/0004-6361/201832638

CrossRef Full Text | Google Scholar

Kuiper, R., Klahr, H., Beuther, H., and Henning, T. (2010a). Circumventing the radiation pressure barrier in the formation of massive stars via disk accretion. Astrophys. J. 722, 1556–1576. doi: 10.1088/0004-637X/722/2/1556

CrossRef Full Text | Google Scholar

Kuiper, R., Klahr, H., Beuther, H., and Henning, T. (2011). Three-dimensional simulation of massive star formation in the disk accretion scenario. Astrophys. J. 732:20. doi: 10.1088/0004-637X/732/1/20

CrossRef Full Text | Google Scholar

Kuiper, R., Kuiper, R., Klahr, H., Klahr, H., Dullemond, C., Dullemond, C., et al. (2010b). Fast and accurate frequency-dependent radiation transport for hydrodynamics simulations in massive star formation. Astron. Astrophys. 511:81. doi: 10.1051/0004-6361/200912355

CrossRef Full Text | Google Scholar

Kuiper, R., Turner, N. J., and Yorke, H. W. (2016). Protostellar outflows and radiative feedback from massive stars. II. Feedback, star-formation efficiency, and outflow broadening. Astrophys. J. 832:40. doi: 10.3847/0004-637X/832/1/40

CrossRef Full Text | Google Scholar

Kuiper, R., and Yorke, H. (2013). On the simultaneous evolution of massive protostars and their host cores. Astrophys. J. 772:61. doi: 10.1088/0004-637X/772/1/61

CrossRef Full Text | Google Scholar

Kuiper, R., Yorke, H. W., and Turner, N. J. (2015). Protostellar outflows and radiative feedback from massive stars. Astrophys. J. 800:86. doi: 10.1088/0004-637X/800/2/86

CrossRef Full Text | Google Scholar

Kunz, M. W., and Mouschovias, T. C. (2009). The nonisothermal stage of magnetic star formation. I. Formulation of the problem and method of solution. Astrophys. J. 693, 1895–1911. doi: 10.1088/0004-637X/693/2/1895

CrossRef Full Text | Google Scholar

Kunz, M. W., and Mouschovias, T. C. (2010). The non-isothermal stage of magnetic star formation - II. Results. Month. Notices R. Astron. Soc. 408, 322–341. doi: 10.1111/j.1365-2966.2010.17110.x

CrossRef Full Text | Google Scholar

Laibe, G., and Price, D. J. (2012). Dusty gas with smoothed particle hydrodynamics - I. Algorithm and test suite. Month. Notices R. Astron. Soc. 420, 2345–2364. doi: 10.1111/j.1365-2966.2011.20202.x

CrossRef Full Text | Google Scholar

Laibe, G., and Price, D. J. (2014a). Dusty gas with one fluid. Month. Notices R. Astron. Soc. 440, 2136–2146. doi: 10.1093/mnras/stu355

CrossRef Full Text | Google Scholar

Laibe, G., and Price, D. J. (2014b). Dusty gas with one fluid in smoothed particle hydrodynamics. Month. Notices R. Astron. Soc. 440, 2147–2163. doi: 10.1093/mnras/stu359

CrossRef Full Text | Google Scholar

Larson, R. B. (1969). Numerical calculations of the dynamics of collapsing proto-star. Month. Notices R. Astron. Soc. 145:271. doi: 10.1093/mnras/145.3.271

CrossRef Full Text | Google Scholar

Lebreuilly, U., Commerçon, B., and Laibe, G. (2019). Small dust grain dynamics on adaptive mesh refinement grids I. Methods. Astron. Astrophys. 626:96. doi: 10.1051/0004-6361/201834147

CrossRef Full Text | Google Scholar

Lee, A. T., Cunningham, A. J., McKee, C. F., and Klein, R. I. (2014). Bondi-Hoyle accretion in an isothermal magnetized plasma. Astrophys. J. 783:50. doi: 10.1088/0004-637X/783/1/50

CrossRef Full Text | Google Scholar

Lee, D., and Deane, A. E. (2009). An unsplit staggered mesh scheme for multidimensional magnetohydrodynamics. J. Comput. Phys. 228, 952–975. doi: 10.1016/j.jcp.2008.08.026

CrossRef Full Text | Google Scholar

Lee, Y.-N., and Hennebelle, P. (2018a). Stellar mass spectrum within massive collapsing clumps. I. Influence of the initial conditions. Astron. Astrophys. 611:A88. doi: 10.1051/0004-6361/201731522

CrossRef Full Text | Google Scholar

Lee, Y.-N., and Hennebelle, P. (2018b). Stellar mass spectrum within massive collapsing clumps. II. Thermodynamics and tidal forces of the first Larson core. A robust mechanism for the peak of the IMF. Astron. Astrophys. 611:A89. doi: 10.1051/0004-6361/201731523

CrossRef Full Text | Google Scholar

Lee, Y.-N., and Hennebelle, P. (2018c). Stellar mass spectrum within massive collapsing clumps III. Effects of temperature and magnetic field. arXiv:1812.05508. doi: 10.1051/0004-6361/201834428

CrossRef Full Text | Google Scholar

Lesur, G., Kunz, M. W., and Fromang, S. (2014). Thanatology in protoplanetary discs. The combined influence of Ohmic, Hall, and ambipolar diffusion on dead zones. Astron. Astrophys. 566:A56. doi: 10.1051/0004-6361/201423660

CrossRef Full Text | Google Scholar

Levermore, C. D. (1984). Relating Eddington factors to flux limiters. J. Quant. Spectrosc. Radiat. Transfer 31, 149–160. doi: 10.1016/0022-4073(84)90112-2

CrossRef Full Text | Google Scholar

Li, F., and Shu, C.-W. (2005). Locally divergence-free discontinuous galerkin methods for MHD equations. J. Sci. Comput. 22–23, 413–442. doi: 10.1007/s10915-004-4146-4

CrossRef Full Text | Google Scholar

Li, F., and Xu, L. (2012). Arbitrary order exactly divergence-free central discontinuous Galerkin methods for ideal MHD equations. J. Comput. Phys. 231, 2655–2675. doi: 10.1016/j.jcp.2011.12.016

CrossRef Full Text | Google Scholar

Li, P. S., Klein, R. I., and McKee, C. F. (2018). Formation of stellar clusters in magnetized, filamentary infrared dark clouds. Month. Notices R. Astron. Soc. 473, 4220–4241. doi: 10.1093/mnras/stx2611

CrossRef Full Text | Google Scholar

Li, P. S., Martin, D. F., Klein, R. I., and McKee, C. F. (2012). A stable, accurate methodology for high mach number, strong magnetic field MHD turbulence with adaptive mesh refinement: resolution and refinement studies. Astrophys. J. 745:139. doi: 10.1088/0004-637X/745/2/139

CrossRef Full Text | Google Scholar

Li, P. S., McKee, C. F., and Klein, R. I. (2006). The heavy-ion approximation for ambipolar diffusion calculations for weakly ionized plasmas. Astrophys. J. 653, 1280–1291. doi: 10.1086/508977

CrossRef Full Text | Google Scholar

Li, Z.-Y., Krasnopolsky, R., and Shang, H. (2011). Non-ideal Mhd effects and magnetic braking catastrophe in protostellar disk formation. Astrophys. J. 738:180. doi: 10.1088/0004-637X/738/2/180

CrossRef Full Text | Google Scholar

Li, Z.-Y., Krasnopolsky, R., Shang, H., and Zhao, B. (2014). On the role of pseudodisk warping and reconnection in protostellar disk formation in turbulent magnetized cores. Astrophys. J. 793:130. doi: 10.1088/0004-637X/793/2/130

CrossRef Full Text | Google Scholar

Li, Z.-Y., and McKee, C. F. (1996). Hydromagnetic accretion shocks around low-mass protostars. Astrophys. J. 464:373. doi: 10.1086/177329

CrossRef Full Text | Google Scholar

Li, Z.-Y., and Nakamura, F. (2006). Cluster formation in protostellar outflow-driven turbulence. Astrophys. J. 640, L187–L190. doi: 10.1086/503419

CrossRef Full Text | Google Scholar

Lin, M.-K., and Youdin, A. N. (2017). A thermodynamic view of dusty protoplanetary disks. Astrophys. J. 849:129. doi: 10.3847/1538-4357/aa92cd

CrossRef Full Text | Google Scholar

Lomax, O., and Whitworth, A. P. (2016). SPAMCART: a code for smoothed particle Monte Carlo radiative transfer. Month. Notices R. Astron. Soc. 461, 3542–3551. doi: 10.1093/mnras/stw1568

CrossRef Full Text | Google Scholar

Londrillo, P., and del Zanna, L. (2004). On the divergence-free condition in Godunov-type schemes for ideal magnetohydrodynamics: the upwind constrained transport method. J. Comput. Phys. 195, 17–48. doi: 10.1016/j.jcp.2003.09.016

CrossRef Full Text | Google Scholar

Lorén-Aguilar, P., and Bate, M. R. (2014). Two-fluid dust and gas mixtures in smoothed particle hydrodynamics: a semi-implicit approach. Month. Notices R. Astron. Soc. 443, 927–945. doi: 10.1093/mnras/stu1173

CrossRef Full Text | Google Scholar

Mac Low, M.-M., and Klessen, R. S. (2004). Control of star formation by supersonic turbulence. Rev. Mod. Phys. 76, 125–194. doi: 10.1103/RevModPhys.76.125

CrossRef Full Text | Google Scholar

Mac Low, M.-M., Norman, M. L., Konigl, A., and Wardle, M. (1995). Incorporation of ambipolar diffusion into the ZEUS magnetohydrodynamics code. Astrophys. J. 442:726. doi: 10.1086/175477

CrossRef Full Text | Google Scholar

Machida, M. N., Inutsuka, S.-i., and Matsumoto, T. (2006). Second core formation and high-speed jets: resistive magnetohydrodynamic nested grid simulations. Astrophys. J. Lett. 2, 151–154. doi: 10.1086/507179

CrossRef Full Text | Google Scholar

Machida, M. N., Inutsuka, S.-I., and Matsumoto, T. (2009). First direct simulation of brown dwarf formation in a compact cloud core. Astrophys. J. 699, L157–L160. doi: 10.1088/0004-637X/699/2/L157

CrossRef Full Text | Google Scholar

Machida, M. N., Inutsuka, S.-i., and Matsumoto, T. (2014). Conditions for circumstellar disc formation: effects of initial cloud configuration and sink treatment. Month. Notices R. Astron. Soc. 438, 2278–2306. doi: 10.1093/mnras/stt2343

CrossRef Full Text | Google Scholar

Marchand, P., Commerçon, B., and Chabrier, G. (2018). Impact of the Hall effect in star formation and the issue of angular momentum conservation. arXiv:1808.08731, 37. doi: 10.1051/0004-6361/201832907

CrossRef Full Text | Google Scholar

Marchand, P., Masson, J., Chabrier, G., Hennebelle, P., Commerçon, B., and Vaytet, N. (2016). Chemical solver to compute molecule and grain abundances and non-ideal MHD resistivities in prestellar core-collapse calculations. Astron. Astrophys. 592:A18. doi: 10.1051/0004-6361/201526780

CrossRef Full Text | Google Scholar

Marinacci, F., Vogelsberger, M., Kannan, R., Mocz, P., Pakmor, R., and Springel, V. (2018). Non-ideal magnetohydrodynamics on a moving mesh. Month. Notices R. Astron. Soc. 476, 2476–2492. doi: 10.1093/mnras/sty397

CrossRef Full Text | Google Scholar

Martizzi, D., Fielding, D., Faucher-Giguère, C.-A., and Quataert, E. (2016). Supernova feedback in a local vertically stratified medium: interstellar turbulence and galactic winds. Month. Notices R. Astron. Soc. 459, 2311–2326. doi: 10.1093/mnras/stw745

CrossRef Full Text | Google Scholar

Masson, J., Chabrier, G., Hennebelle, P., Vaytet, N., and Commerçon, B. (2016). Ambipolar diffusion in low-mass star formation. I. General comparison with the ideal magnetohydrodynamic case. Astron. Astrophys. 587:A32. doi: 10.1051/0004-6361/201526371

CrossRef Full Text | Google Scholar

Masson, J., Teyssier, R., Mulet-Marquis, C., Hennebelle, P., and Chabrier, G. (2012). Incorporating ambipolar and ohmic diffusion in the AMR MHD code RAMSES. Astrophys. J. Suppl. 201:24. doi: 10.1088/0067-0049/201/2/24

CrossRef Full Text | Google Scholar

Masunaga, H., and Inutsuka, S.-i. (2000). A radiation hydrodynamic model for protostellar collapse. II. The second collapse and the birth of a protostar. Astrophys. J. 531, 350–365. doi: 10.1086/308439

CrossRef Full Text | Google Scholar

Mathis, J. S., Rumpl, W., and Nordsieck, K. H. (1977). The size distribution of interstellar grains. Astrophys. J. 217, 425–433. doi: 10.1086/155591

CrossRef Full Text | Google Scholar

Matsumoto, T. (2011). An implicit scheme for ohmic dissipation with adaptive mesh refinement. Publ. Astron. Soc. Jpn. 63, 317–323. doi: 10.1093/pasj/63.2.317

CrossRef Full Text | Google Scholar

Matsumoto, T., Machida, M. N., and ichiro Inutsuka, S. (2017). Circumstellar disks and outflows in turbulent molecular cloud cores: possible formation mechanism for misaligned systems. Astrophys. J. 839:69. doi: 10.3847/1538-4357/aa6a1c

CrossRef Full Text | Google Scholar

Matzner, C. D. (2007). Protostellar outflow-driven turbulence. Astrophys. J. 659, 1394–1403. doi: 10.1086/512361

CrossRef Full Text | Google Scholar

Matzner, C. D., and McKee, C. F. (1999). Bipolar molecular outflows driven by hydromagnetic protostellar winds. Astrophys. J. 526, L109–L112. doi: 10.1086/312376

PubMed Abstract | CrossRef Full Text | Google Scholar

Matzner, C. D., and McKee, C. F. (2000). Efficiencies of low-mass star and star cluster formation. Astrophys. J. 545, 364–378. doi: 10.1086/317785

CrossRef Full Text | Google Scholar

Maureira, M. J., Arce, H. G., Dunham, M. M., Pineda, J. E., Fernández-López, M., Chen, X., et al. (2017). Kinematics of a young low-mass star-forming core: understanding the evolutionary state of the first-core candidate L1451-mm. Astrophys. J. 838:60. doi: 10.3847/1538-4357/838/1/60

CrossRef Full Text | Google Scholar

Mayer, L., Lufkin, G., Quinn, T., and Wadsley, J. (2007). Fragmentation of gravitationally unstable gaseous protoplanetary disks with radiative transfer. Astrophys. J. 661, L77–L80. doi: 10.1086/518433

CrossRef Full Text | Google Scholar

McKee, C. F. (1989). Photoionization-regulated star formation and the structure of molecular clouds. Astrophys. J. 345, 782–801. doi: 10.1086/167950

CrossRef Full Text | Google Scholar

McKee, C. F., and Ostriker, E. C. (2007). Theory of star formation. Annu. Rev. Astro. Astrophys. 45, 565–687. doi: 10.1146/annurev.astro.45.051806.110602

CrossRef Full Text | Google Scholar

Meheut, H., Meliani, Z., Varniere, P., and Benz, W. (2012). Dust-trapping Rossby vortices in protoplanetary disks. Astron. Astrophys. 545:A134. doi: 10.1051/0004-6361/201219794

CrossRef Full Text | Google Scholar

Mellema, G., Iliev, I. T., Alvarez, M. A., and Shapiro, P. R. (2006). C 2-ray: a new method for photon-conserving transport of ionizing radiation. New Astron. 11, 374–395. doi: 10.1016/j.newast.2005.09.004

CrossRef Full Text | Google Scholar

Mellon, R. R., and Li, Z. (2008). Magnetic braking and protostellar disk formation: the ideal MHD limit. Astrophys. J. 681, 1356–1376. doi: 10.1086/587542

CrossRef Full Text | Google Scholar

Meyer, C. D., Balsara, D. S., and Aslam, T. D. (2012). A second-order accurate Super TimeStepping formulation for anisotropic thermal conduction. Month. Notices R. Astron. Soc. 422, 2102–2115. doi: 10.1111/j.1365-2966.2012.20744.x

CrossRef Full Text | Google Scholar

Meyer, C. D., Balsara, D. S., and Aslam, T. D. (2014). A stabilized Runge-Kutta-Legendre method for explicit super-time-stepping of parabolic and mixed equations. J. Comput. Phys. 257, 594–626. doi: 10.1016/j.jcp.2013.08.021

CrossRef Full Text | Google Scholar

Mignone, A., Bodo, G., Massaglia, S., Matsakos, T., Tesileanu, O., Zanni, C., and Ferrari, A. (2007). PLUTO: a numerical code for computational astrophysics. Astrophys. J. Suppl. Ser. 170, 228–242. doi: 10.1086/513316

CrossRef Full Text | Google Scholar

Mihalas, D., and Klein, R. I. (1982). On the solution of the time-dependent inertial-frame equation of radiative transfer in moving media to O(vc). J. Comput. Phys. 46, 97–137. doi: 10.1016/0021-9991(82)90007-9

CrossRef Full Text | Google Scholar

Mihalas, D., and Mihalas, B. W. (1984). Foundations of Radiation Hydrodynamics. New York, NY: Oxford University Press.

Google Scholar

Minerbo, G. N. (1978). Maximum entropy Eddington factors. J. Quant. Spectr. Radiat. Transf. 20, 541–545. doi: 10.1016/0022-4073(78)90024-9

CrossRef Full Text | Google Scholar

Miyama, S. M., Tomisaka, K., and Hanawa, T. (1999). “Numerical astrophysics,” in Proceedings of the International Conference on Numerical Astrophysics 1998 (NAP98) (Boston, MA: Kluwer Academic), 383.

Google Scholar

Mocz, P., Burkhart, B., Hernquist, L., McKee, C. F., and Springel, V. (2017). Moving-mesh Simulations of Star-forming Cores in Magneto-gravo-turbulence. Astrophys. J. 838:40. doi: 10.3847/1538-4357/aa6475

CrossRef Full Text | Google Scholar

Mocz, P., Vogelsberger, M., and Hernquist, L. (2014a). A constrained transport scheme for MHD on unstructured static and moving meshes. Month. Notices R. Astron. Soc. 442, 43–55. doi: 10.1093/mnras/stu865

CrossRef Full Text | Google Scholar

Mocz, P., Vogelsberger, M., Sijacki, D., Pakmor, R., and Hernquist, L. (2014b). A discontinuous Galerkin method for solving the fluid and magnetohydrodynamic equations in astrophysical simulations. Month. Notices R. Astron. Soc. 437, 397–414. doi: 10.1093/mnras/stt1890

CrossRef Full Text | Google Scholar

Murray, D., Goyal, S., and Chang, P. (2018). The effects of protostellar jet feedback on turbulent collapse. Month. Notices R. Astron. Soc. 475, 1023–1035. doi: 10.1093/mnras/stx3153

CrossRef Full Text | Google Scholar

Myers, A. T., Klein, R. I., Krumholz, M. R., and McKee, C. F. (2014). Star cluster formation in turbulent, magnetized dense clumps with radiative and outflow feedback. Month. Notices R. Astron. Soc. 439, 3420–3438. doi: 10.1093/mnras/stu190

CrossRef Full Text | Google Scholar

Myers, A. T., McKee, C. F., Cunningham, A. J., Klein, R. I., and Krumholz, M. R. (2013). The fragmentation of magnetized, massive star-forming cores with radiative feedback. Astrophys. J. 766:97. doi: 10.1088/0004-637X/766/2/97

CrossRef Full Text | Google Scholar

Nakamura, F., and Li, Z.-Y. (2007). Protostellar turbulence driven by collimated outflows. Astrophys. J. 662, 395–412. doi: 10.1086/517515

CrossRef Full Text | Google Scholar

Nakamura, F., and Li, Z.-Y. (2011). Clustered star formation in magnetic clouds: properties of dense cores formed in outflow-driven turbulence. Astrophys. J. 740:36. doi: 10.1088/0004-637X/740/1/36

CrossRef Full Text | Google Scholar

Nakano, T., Nishi, R., and Umebayashi, T. (2002). Mechanism of magnetic flux loss in molecular clouds. Astrophys. J. 573, 199–214. doi: 10.1086/340587

CrossRef Full Text | Google Scholar

Nayakshin, S., Cha, S.-H., and Hobbs, A. (2009). Dynamic Monte Carlo radiation transfer in SPH: radiation pressure force implementation. Month. Notices R. Astron. Soc. 397, 1314–1325. doi: 10.1111/j.1365-2966.2009.15091.x

CrossRef Full Text | Google Scholar

Ocvirk, P., Gillet, N., Shapiro, P. R., Aubert, D., Iliev, I. T., Teyssier, R., et al. (2016). Cosmic Dawn (CoDa): the first radiation-hydrodynamics simulation of reionization and galaxy formation in the local universe. Month. Notices R. Astron. Soc. 463, 1462–1485. doi: 10.1093/mnras/stw2036

CrossRef Full Text | Google Scholar

Offner, S. S. R., and Chaban, J. (2017). Impact of protostellar outflows on turbulence and star formation efficiency in magnetized dense cores. Astrophys. J. 847:104. doi: 10.3847/1538-4357/aa8996

CrossRef Full Text | Google Scholar

Offner, S. S. R., Klein, R. I., McKee, C. F., and Krumholz, M. R. (2009). The effects of radiative transfer on low-mass star formation. Astrophys. J. 703, 131–149. doi: 10.1088/0004-637X/703/1/131

CrossRef Full Text | Google Scholar

Offner, S. S. R., and Liu, Y. (2018). Turbulent action at a distance due to stellar feedback in magnetized clouds. Nat. Astron. 2, 896. doi: 10.1038/s41550-018-0566-1

CrossRef Full Text | Google Scholar

Oishi, J. S., and Mac Low, M.-M. (2006). The inability of ambipolar diffusion to set a characteristic mass scale in molecular clouds. Astrophys. J. 638, 281–285. doi: 10.1086/498818

CrossRef Full Text | Google Scholar

O'Sullivan, S., and Downes, T. P. (2006). An explicit scheme for multifluid magnetohydrodynamics. Month. Notices R. Astron. Soc. 366, 1329–1336. doi: 10.1111/j.1365-2966.2005.09898.x

CrossRef Full Text | Google Scholar

Padoan, P., Haugbølle, T., Nordlund, Å., and Frimann, S. (2017). Supernova driving. IV. The star formation rate of molecular clouds. Astrophys. J. 840:48. doi: 10.3847/1538-4357/aa6afa

CrossRef Full Text | Google Scholar

Padoan, P., and Nordlund, Å. (2011). The star formation rate of supersonic magnetohydrodynamic turbulence. Astrophys. J. 730:40. doi: 10.1088/0004-637X/730/1/40

CrossRef Full Text | Google Scholar

Padoan, P., Zweibel, E., and Nordlund, Å. (2000). Ambipolar drift heating in turbulent molecular clouds. Astrophys. J. 540, 332–341. doi: 10.1086/309299

CrossRef Full Text | Google Scholar

Pawlik, A. H., and Schaye, J. (2008). TRAPHIC - radiative transfer for smoothed particle hydrodynamics simulations. Month. Notices R. Astron. Soc. 389, 651–677. doi: 10.1111/j.1365-2966.2008.13601.x

CrossRef Full Text | Google Scholar

Paxton, B., Bildsten, L., Dotter, A., Herwig, F., Lesaffre, P., and Timmes, F. (2011). Modules for experiments in stellar astrophysics (MESA). Astrophys. J. Suppl. Ser. 192:3. doi: 10.1088/0067-0049/192/1/3

CrossRef Full Text | Google Scholar

Pelletier, G., and Pudritz, R. E. (1992). Hydromagnetic disk winds in young stellar objects and active galactic nuclei. Astrophys. J. 394, 117–138. doi: 10.1086/171565

CrossRef Full Text | Google Scholar

Pelupessy, F. I., van Elteren, A., de Vries, N., McMillan, S. L. W., Drost, N., and Portegies Zwart, S. F. (2013). The astrophysical multipurpose software environment. Astron. Astrophys. 557:A84. doi: 10.1051/0004-6361/201321252

CrossRef Full Text | Google Scholar

Peters, T., Banerjee, R., Klessen, R. S., and Mac Low, M.-M. (2011). The interplay of magnetic fields, fragmentation, and ionization feedback in high-mass star formation. Astrophys. J. 729:72. doi: 10.1088/0004-637X/729/1/72

CrossRef Full Text | Google Scholar

Peters, T., Banerjee, R., Klessen, R. S., Mac Low, M.-M., Galván-Madrid, R., and Keto, E. R. (2010). H II regions: witnesses to massive star formation. Astrophys. J. 711, 1017–1028. doi: 10.1088/0004-637X/711/2/1017

CrossRef Full Text | Google Scholar

Peters, T., Klaassen, P. D., Mac Low, M.-M., Schrön, M., Federrath, C., Smith, M. D., et al. (2014). Collective outflow from a small multiple stellar system. Astrophys. J. 788:14. doi: 10.1088/0004-637X/788/1/14

CrossRef Full Text | Google Scholar

Peters, T., Naab, T., Walch, S., Glover, S. C. O., Girichidis, P., Pellegrini, E., et al. (2017). The SILCC project - IV. Impact of dissociating and ionizing radiation on the interstellar medium and Hα emission as a tracer of the star formation rate. Month. Notices R. Astron. Soc. 466, 3293–3308. doi: 10.1093/mnras/stw3216

CrossRef Full Text | Google Scholar

Petkova, M., and Springel, V. (2011). A novel approach for accurate radiative transfer in cosmological hydrodynamic simulations. Month. Notices R. Astron. Soc. 415, 3731–3749. doi: 10.1111/j.1365-2966.2011.18986.x

CrossRef Full Text | Google Scholar

Phillips, G. J., and Monaghan, J. J. (1985). A numerical method for three-dimensional simulations of collapsing, isothermal, magnetic gas clouds. Month. Notices R. Astron. Soc. 216, 883–895. doi: 10.1093/mnras/216.4.883

CrossRef Full Text | Google Scholar

Pinto, C., Galli, D., and Bacciotti, F. (2008). Three-fluid plasmas in star formation. I. Magneto-hydrodynamic equations. Astron. Astrophys. 484, 1–15. doi: 10.1051/0004-6361:20078818

CrossRef Full Text | Google Scholar

Portegies Zwart, S. F., and Verbunt, F. (1996). Population synthesis of high-mass binaries. Astron. Astrophys. 309, 179–196.

Google Scholar

Powell, K. G., Roe, P. L., Linde, T. J., Gombosi, T. I., and De Zeeuw, D. L. (1999). A solution-adaptive upwind scheme for ideal magnetohydrodynamics. J. Comput. Phys. 154, 284–309. doi: 10.1006/jcph.1999.6299

CrossRef Full Text | Google Scholar

Price, D. J., and Bate, M. R. (2008). The effect of magnetic fields on star cluster formation. Month. Notices R. Astron. Soc. 385, 1820–1834. doi: 10.1111/j.1365-2966.2008.12976.x

CrossRef Full Text | Google Scholar

Price, D. J., and Bate, M. R. (2009). Inefficient star formation: the combined effects of magnetic fields and radiative feedback. Month. Notices R. Astron. Soc. 398, 33–46. doi: 10.1111/j.1365-2966.2009.14969.x

CrossRef Full Text | Google Scholar

Price, D. J., and Laibe, G. (2015). A fast and explicit algorithm for simulating the dynamics of small dust grains with smoothed particle hydrodynamics. Month. Notices R. Astron. Soc. 451, 813–826. doi: 10.1093/mnras/stv996

CrossRef Full Text | Google Scholar

Price, D. J., and Monaghan, J. J. (2005). Smoothed particle magnetohydrodynamics — III. multidimensional tests and the div B= 0 constraint. Month. Notices R. Astron. Soc. 364, 384–406. doi: 10.1111/j.1365-2966.2005.09576.x

CrossRef Full Text | Google Scholar

Price, D. J., Tricco, T. S., and Bate, M. R. (2012). Collimated jets from the first core. Month. Notices R. Astron. Soc. 423, L45–L49. doi: 10.1111/j.1745-3933.2012.01254.x

CrossRef Full Text | Google Scholar

Price, D. J., Wurster, J., Tricco, T. S., Nixon, C., Toupin, S., Pettitt, A., et al. (2017). Phantom: a smoothed particle hydrodynamics and magnetohydrodynamics code for astrophysics. arXiv:1702.03930. doi: 10.1017/pasa.2018.25

CrossRef Full Text | Google Scholar

Ramsey, J. P., and Dullemond, C. P. (2015). Radiation hydrodynamics including irradiation and adaptive mesh refinement with AZEuS. I. Methods. Astron. Astrophys. 574:A81. doi: 10.1051/0004-6361/201424954

CrossRef Full Text | Google Scholar

Renaud, F., Bournaud, F., Emsellem, E., Elmegreen, B., Teyssier, R., Alves, J., et al. (2013). A sub-parsec resolution simulation of the Milky Way: global structure of the interstellar medium and properties of molecular clouds. Month. Notices R. Astron. Soc. 436, 1836–1851. doi: 10.1093/mnras/stt1698

CrossRef Full Text | Google Scholar

Rijkhorst, E. J., Plewa, T., Dubey, A., and Mellema, G. (2006). Hybrid characteristics: 3D radiative transfer for parallel adaptive mesh refinement hydrodynamics. Astron. Astrophys. 452, 907–920. doi: 10.1051/0004-6361:20053401

CrossRef Full Text | Google Scholar

Ripoll, J. F., Dubroca, B., and Duffa, G. (2001). Modelling radiative mean absorption coefficients. Combust. Theory Model. 5, 261–274. doi: 10.1088/1364-7830/5/3/301

CrossRef Full Text | Google Scholar

Rosdahl, J., Blaizot, J., Aubert, D., Stranex, T., and Teyssier, R. (2013). RAMSES-RT: radiation hydrodynamics in the cosmological context. Month. Notices R. Astron. Soc. 436, 2188–2231. doi: 10.1093/mnras/stt1722

CrossRef Full Text | Google Scholar

Rosdahl, J., and Teyssier, R. (2015). A scheme for radiation pressure and photon diffusion with the M1 closure in RAMSES-RT. Month. Notices R. Astron. Soc. 449, 4380–4403. doi: 10.1093/mnras/stv567

CrossRef Full Text | Google Scholar

Rosen, A. L., Krumholz, M. R., McKee, C. F., and Klein, R. I. (2016). An unstable truth: how massive stars get their mass. Month. Notices R. Astron. Soc. 463, 2553–2573. doi: 10.1093/mnras/stw2153

CrossRef Full Text | Google Scholar

Rosen, A. L., Krumholz, M. R., Oishi, J. S., Lee, A. T., and Klein, R. I. (2017). Hybrid adaptive ray-moment method (HARM²): a highly parallel method for radiation hydrodynamics on adaptive grids. J. Comput. Phys. 330, 924–942. doi: 10.1016/j.jcp.2016.10.048

CrossRef Full Text | Google Scholar

Roth, N., and Kasen, D. (2015). Monte Carlo radiation-hydrodynamics with implicit methods. Astrophys. J. Suppl. 217:9. doi: 10.1088/0067-0049/217/1/9

CrossRef Full Text | Google Scholar

Rozyczka, M. (1985). Two-dimensional models of stellar wind bubbles. I - Numerical methods and their application to the investigation of outer shell instabilities. Astron. Astrophys. 143, 59–71.

Google Scholar

Ryan, B. R., Gammie, C. F., Fromang, S., and Kestener, P. (2017). Resolution dependence of magnetorotational turbulence in the isothermal stratified shearing box. Astrophys. J. 840:6. doi: 10.3847/1538-4357/aa6a52

CrossRef Full Text | Google Scholar

Ryu, D., Jones, T. W., and Frank, A. (1995). Numerical magnetohydrodynamics in astrophysics: algorithm and tests for multidimensional flow. Astrophys. J. 452:785. doi: 10.1086/176347

CrossRef Full Text | Google Scholar

Ryu, D., Miniati, F., Jones, T. W., and Frank, A. (1998). A divergence-free upwind code for multidimensional magnetohydrodynamic flows. Astrophys. J. 509, 244–255. doi: 10.1086/306481

CrossRef Full Text | Google Scholar

Sano, T., and Stone, J. M. (2002). The effect of the hall term on the nonlinear evolution of the magnetorotational instability. I. Local axisymmetric simulations. Astrophys. J. 570, 314–328. doi: 10.1086/339504

CrossRef Full Text | Google Scholar

Saumon, D., Chabrier, G., and van Horn, H. M. (1995). An equation of state for low-mass stars and giant planets. Astrophys. J. Suppl. Ser. 99:713. doi: 10.1086/192204

CrossRef Full Text | Google Scholar

Seifried, D., Walch, S., Girichidis, P., Naab, T., Wünsch, R., Klessen, R. S., et al. (2017). SILCC-Zoom: the dynamic and chemical evolution of molecular clouds. Month. Notices R. Astron. Soc. 472, 4797–4818. doi: 10.1093/mnras/stx2343

CrossRef Full Text | Google Scholar

Shestakov, A. I., and Offner, S. S. R. (2008). A multigroup diffusion solver using pseudo transient continuation for a radiation-hydrodynamic code with patch-based AMR. J. Comput. Phys. 227, 2154–2186. doi: 10.1016/j.jcp.2007.09.019

CrossRef Full Text | Google Scholar

Shima, K., Tasker, E. J., and Habe, A. (2017). Does feedback help or hinder star formation? The effect of photoionization on star formation in giant molecular clouds. Month. Notices R. Astron. Soc. 467, 512–523. doi: 10.1093/mnras/stw3279

CrossRef Full Text | Google Scholar

Shu, F., Najita, J., Ostriker, E., Wilkin, F., Ruden, S., and Lizano, S. (1994). Magnetocentrifugally driven flows from young stars and disks. 1: a generalized model. Astrophys. J. 429, 781–796. doi: 10.1086/174363

CrossRef Full Text | Google Scholar

Shu, F. H., Adams, F. C., and Lizano, S. (1987). Star formation in molecular clouds: observation and theory. Annu. Rev. Astron. Astrophys. 25, 23–81. doi: 10.1146/annurev.aa.25.090187.000323

CrossRef Full Text | Google Scholar

Shu, F. H., Tremaine, S., Adams, F. C., and Ruden, S. P. (1990). SLING amplification and eccentric gravitational instabilities in gaseous disks. Astrophys. J. 358:495. doi: 10.1086/169003

CrossRef Full Text | Google Scholar

Sigalotti, L. D. G., Cruz, F., Gabbasov, R., Klapp, J., and Ramírez-Velasquez, J. (2018). From large-scale to protostellar disk fragmentation into close binary stars. Astrophys. J. 857:40. doi: 10.3847/1538-4357/aab619

CrossRef Full Text | Google Scholar

Skinner, M. A., and Ostriker, E. C. (2013). A Two-moment radiation hydrodynamics module in athena using a time-explicit godunov method. Astrophys. J. Suppl. 206:21. doi: 10.1088/0067-0049/206/2/21

CrossRef Full Text | Google Scholar

Smith, M. C., Sijacki, D., and Shen, S. (2018). Supernova feedback in numerical simulations of galaxy formation: separating physics from numerics. Month. Notices R. Astron. Soc. 478, 302–331. doi: 10.1093/mnras/sty994

CrossRef Full Text | Google Scholar

Springel, V. (2010). E pur si muove: galilean-invariant cosmological hydrodynamical simulations on a moving mesh. Month. Notices R. Astron. Soc. 401, 791–851. doi: 10.1111/j.1365-2966.2009.15715.x

CrossRef Full Text | Google Scholar

Spruit, H. C., Stehle, R., and Papaloizou, J. C. B. (1995). Interchange instability in and accretion disc with a poloidal magnetic field. Month. Notices R. Astron. Soc. 275, 1223–1231. doi: 10.1093/mnras/275.4.1223

CrossRef Full Text | Google Scholar

Spruit, H. C., and Taam, R. E. (1990). Mass transport in a neutron star magnetosphere. Astron. Astrophys. 229, 475–493.

Google Scholar

Stahler, S. W., Shu, F. H., and Taam, R. E. (1980). The evolution of protostars. I - Global formulation and results. Astrophys. J. 241, 637–654. doi: 10.1086/158377

CrossRef Full Text | Google Scholar

Stamatellos, D., Whitworth, A., Bisbas, T., and Goodwin, S. (2007). Radiative transfer and the energy equation in SPH simulations of star formation. Astron. Astrophys. 475, 37–49. doi: 10.1051/0004-6361:20077373

CrossRef Full Text | Google Scholar

Stamatellos, D., Whitworth, A. P., and Hubber, D. A. (2011). The importance of episodic accretion for low-mass star formation. Astrophys. J. 730:32. doi: 10.1088/0004-637X/730/1/32

CrossRef Full Text | Google Scholar

Stone, J. M., Gardiner, T. A., Teuben, P., Hawley, J. F., and Simon, J. B. (2008). Athena: a new code for astrophysical MHD. Astrophys. J. Suppl. Ser. 178, 137–177. doi: 10.1086/588755

CrossRef Full Text | Google Scholar

Stone, J. M., Mihalas, D., and Norman, M. L. (1992a). ZEUS-2D: a radiation magnetohydrodynamics code for astrophysical flows in two space dimensions. III - The radiation hydrodynamic algorithms and tests. Astrophys. J. Suppl. Ser. 80, 819–845. doi: 10.1086/191682

CrossRef Full Text | Google Scholar

Stone, J. M., Stone, J. M., and Norman, M. L. (1992b). ZEUS-2D: a radiation magnetohydrodynamics code for astrophysical flows in two space dimensions. II. The magnetohydrodynamic algorithms and tests. Astrophys. J. Suppl. 80:791. doi: 10.1086/191681

CrossRef Full Text | Google Scholar

Sur, S., Schleicher, D. R. G., Banerjee, R., Federrath, C., and Klessen, R. S. (2010). The generation of strong magnetic fields during the formation of the first stars. Astrophys. J. Lett. 721, L134–L138. doi: 10.1088/2041-8205/721/2/L134

CrossRef Full Text | Google Scholar

Takasao, S., Tomida, K., Iwasaki, K., and Suzuki, T. K. (2018). A three-dimensional simulation of a magnetized accretion disk: fast funnel accretion onto a weakly magnetized star. Astrophys. J. 857:4. doi: 10.3847/1538-4357/aab5b3

CrossRef Full Text | Google Scholar

Tasker, E. J. (2011). Star formation in disk galaxies. II. The effect Of star formation and photoelectric heating on the formation and evolution of giant molecular clouds. Astrophys. J. 730:11. doi: 10.1088/0004-637X/730/1/11

CrossRef Full Text | Google Scholar

Teyssier, R. (2002). Cosmological hydrodynamics with adaptive mesh refinement. A new high resolution code called RAMSES. Astron. Astrophys. 385, 337–364. doi: 10.1051/0004-6361:20011817

CrossRef Full Text | Google Scholar

Teyssier, R., Fromang, S., and Dormy, E. (2006). Kinematic dynamos using constrained transport with high order Godunov schemes and adaptive mesh refinement. J. Comput. Phys. 218, 44–67. doi: 10.1016/j.jcp.2006.01.042

CrossRef Full Text | Google Scholar

Tilley, D. A., and Balsara, D. S. (2008). A two-fluid method for ambipolar diffusion. Month. Notices R. Astron. Soc. 389, 1058–1073. doi: 10.1111/j.1365-2966.2008.13636.x

CrossRef Full Text | Google Scholar

Tomida, K., Machida, M. N., Hosokawa, T., Sakurai, Y., and Lin, C. H. (2017). Grand-design spiral arms in a young forming circumstellar disk. Astrophys. J. 835:L11. doi: 10.3847/2041-8213/835/1/L11

CrossRef Full Text | Google Scholar

Tomida, K., Okuzumi, S., and Machida, M. N. (2015). Radiation magnetohydrodynamic simulations of protostellar collapse: nonideal magnetohydrodynamic effects and early formation of circumstellar disks. Astrophys. J. 801:117. doi: 10.1088/0004-637X/801/2/117

CrossRef Full Text | Google Scholar

Tomida, K., Tomisaka, K., Matsumoto, T., Hori, Y., Okuzumi, S., Machida, M. N., et al. (2013). Radiation magnetohydrodynamic simulations of protostellar collapse: protostellar core formation. Astrophys. J. 763:6. doi: 10.1088/0004-637X/763/1/6

CrossRef Full Text | Google Scholar

Tomida, K., Tomisaka, K., Matsumoto, T., Ohsuga, K., MacHida, M. N., and Saigo, K. (2010). Radiation magnetohydrodynamics simulation of proto-stellar collapse: two-component molecular outflow. Astrophys. J. Lett. 714(1 Pt 2):58–63. doi: 10.1088/2041-8205/714/1/L58

CrossRef Full Text | Google Scholar

Tomisaka, K. (2002). Collapse of rotating magnetized molecular cloud cores and mass outflows. Astrophys. J. 575, 306–326. doi: 10.1086/341133

CrossRef Full Text | Google Scholar

Tomisaka, K., and Bregman, J. N. (1993). Extended hot-gas halos around starburst galaxies. Publ. Astron. Soc. Jpn. 45, 513–528.

Google Scholar

Toth, G. (1994). Numerical study of two-fluid C-type shock waves. Astrophys. J. 425, 171–194. doi: 10.1086/173973

CrossRef Full Text | Google Scholar

Tóth, G. (1996). A general code for modeling MHD flows on parallel computers: versatile advection code. Astrophys. Lett. Commun. 34:245. doi: 10.1007/978-94-009-0315-9_101

CrossRef Full Text | Google Scholar

Tóth, G. (2000). The div B=0 constraint in shock-capturing magnetohydrodynamics codes. J. Comput. Phys. 161, 605–652. doi: 10.1006/jcph.2000.6519

CrossRef Full Text | Google Scholar

Tóth, G., De Zeeuw, D. L., Gombosi, T. I., and Powell, K. G. (2006). A parallel explicit/implicit time stepping scheme on block-adaptive grids. J. Comput. Phys. 217, 722–758. doi: 10.1016/j.jcp.2006.01.029

CrossRef Full Text | Google Scholar

Tóth, G., Ma, Y., and Gombosi, T. I. (2008). Hall magnetohydrodynamics on block-adaptive grids. J. Comput. Phys. 227, 6967–6984. doi: 10.1016/j.jcp.2008.04.010

CrossRef Full Text | Google Scholar

Tóth, G., and Roe, P. L. (2002). Divergence- and curl-preserving prolongation and restriction formulas. J. Comput. Phys. 180, 736–750. doi: 10.1006/jcph.2002.7120

CrossRef Full Text | Google Scholar

Tricco, T. S., and Price, D. J. (2012). Constrained hyperbolic divergence cleaning for smoothed particle magnetohydrodynamics. J. Comput. Phys. 231, 7214–7236. doi: 10.1016/j.jcp.2012.06.039

CrossRef Full Text | Google Scholar

Tricco, T. S., Price, D. J., and Laibe, G. (2017). Is the dust-to-gas ratio constant in molecular clouds? Month. Notices R. Astron. Soc. 471, L52–L56. doi: 10.1093/mnrasl/slx096

CrossRef Full Text | Google Scholar

Truelove, J. K., Klein, R. I., McKee, C. F., Holliman, J. H. I., Howell, L. H., and Greenough, J. A. (1997). The jeans condition: a new constraint on spatial resolution in simulations of isothermal self-gravitational hydrodynamics. Astrophys. J. Lett. 489:L179. doi: 10.1086/310975

CrossRef Full Text | Google Scholar

Tsitali, A. E., Belloche, A., Commerçon, B., and Menten, K. M. (2013). The dynamical state of the first hydrostatic core candidate Chamaeleon-MMS1. Astron. Astrophys. 557:A98. doi: 10.1051/0004-6361/201321204

CrossRef Full Text | Google Scholar

Tsukamoto, Y., Iwasaki, K., and Inutsuka, S.-i. (2013). An explicit scheme for ohmic dissipation with smoothed particle magnetohydrodynamics. Month. Notices R. Astron. Soc. 434, 2593–2599. doi: 10.1093/mnras/stt1205

CrossRef Full Text | Google Scholar

Tsukamoto, Y., Iwasaki, K., Okuzumi, S., Machida, M. N., and Inutsuka, S. (2015a). Bimodality OF circumstellar disk evolution induced by the hall current. Astrophys. J. Lett. 810:L26. doi: 10.1088/2041-8205/810/2/L26

CrossRef Full Text | Google Scholar

Tsukamoto, Y., Iwasaki, K., Okuzumi, S., Machida, M. N., and Inutsuka, S. (2015b). Effects of Ohmic and ambipolar diffusion on formation and evolution of first cores, protostars, and circumstellar discs. Month. Notices R. Astron. Soc. 452, 278–288. doi: 10.1093/mnras/stv1290

CrossRef Full Text | Google Scholar

Tsukamoto, Y., Okuzumi, S., Iwasaki, K., Machida, M. N., and Inutsuka, S.-I. (2017). The impact of the Hall effect during cloud core collapse: implications for circumstellar disk evolution. Publ. Astron. Soc. Jpn. 69:95. doi: 10.1093/pasj/psx113

CrossRef Full Text | Google Scholar

Turk, M. J., Oishi, J. S., Abel, T., and Bryan, G. L. (2012). Magnetic fields in population III star formation. Astrophys. J. 745:154. doi: 10.1088/0004-637X/745/2/154

CrossRef Full Text | Google Scholar

Vaidya, B., Prasad, D., Mignone, A., Sharma, P., and Rickler, L. (2017). Scalable explicit implementation of anisotropic diffusion with Runge- Kutta-Legendre super-time stepping. Month. Notices R. Astron. Soc. 472, 3147–3160. doi: 10.1093/mnras/stx2176

CrossRef Full Text | Google Scholar

van Leer, B. (1977). Towards the ultimate conservative difference scheme. IV. A new approach to numerical convection. J. Comput. Phys. 23:276. doi: 10.1016/0021-9991(77)90095-X

CrossRef Full Text | Google Scholar

Vaytet, N., Audit, E., Chabrier, G., Commerçon, B., and Masson, J. (2012). Simulations of protostellar collapse using multigroup radiation hydrodynamics: I. the first collapse. Astron. Astrophys. 543:60. doi: 10.1051/0004-6361/201219427

CrossRef Full Text | Google Scholar

Vaytet, N., Chabrier, G., Audit, E., Commerçon, B., Masson, J., Ferguson, J., et al. (2013). Simulations of protostellar collapse using multigroup radiation hydrodynamics. II. The second collapse. Astron. Astrophys. 557:A90. doi: 10.1051/0004-6361/201321423

CrossRef Full Text | Google Scholar

Vaytet, N., Commerçon, B., Masson, J., González, M., and Chabrier, G. (2018). Protostellar birth with ambipolar and ohmic diffusion. Astron. Astrophys. 615:A5. doi: 10.1051/0004-6361/201732075

CrossRef Full Text | Google Scholar

Vaytet, N., and Haugbølle, T. (2017). A grid of one-dimensional low-mass star formation collapse models. Astron. Astrophys. 598:A116. doi: 10.1051/0004-6361/201628194

CrossRef Full Text | Google Scholar

Vaytet, N., Tomida, K., and Chabrier, G. (2014). On the role of the H₂ ortho: para ratio in gravitational collapse during star formation. Astron. Astrophys. 563:A85. doi: 10.1051/0004-6361/201322855

CrossRef Full Text | Google Scholar

Walch, S., Girichidis, P., Naab, T., Gatto, A., Glover, S. C. O., Wünsch, R., et al. (2015). The SILCC (SImulating the LifeCycle of molecular Clouds) project - I. Chemical evolution of the supernova-driven ISM. Month. Notices R. Astron. Soc. 454, 238–268. doi: 10.1093/mnras/stv1975

CrossRef Full Text | Google Scholar

Wall, J. E., McMillan, S. L. W., Mac Low, M.-M., Klessen, R. S., and Portegies Zwart, S. (2019). Collisional N-body dynamics coupled to self-gravitating magnetohydrodynamics reveals dynamical binary formation. arXiv:1901.01132.

Google Scholar

Wang, P., Li, Z.-Y., Abel, T., and Nakamura, F. (2010). Outflow feedback regulated massive star formation in parsec-scale cluster-forming clumps. Astrophys. J. 709, 27–41. doi: 10.1088/0004-637X/709/1/27

CrossRef Full Text | Google Scholar

Whitehouse, S. C., and Bate, M. R. (2004). Smoothed particle hydrodynamics with radiative transfer in the flux-limited diffusion approximation. Month. Notices R. Astron. Soc. 353, 1078–1094. doi: 10.1111/j.1365-2966.2004.08131.x

CrossRef Full Text | Google Scholar

Whitehouse, S. C., and Bate, M. R. (2006). The thermodynamics of collapsing molecular cloud cores using smoothed particle hydrodynamics with radiative transfer. Month. Notices R. Astron. Soc. 367, 32–38. doi: 10.1111/j.1365-2966.2005.09950.x

CrossRef Full Text | Google Scholar

Whitehouse, S. C., Bate, M. R., and Monaghan, J. J. (2005). A faster algorithm for smoothed particle hydrodynamics with radiative transfer in the flux-limited diffusion approximation. Month. Notices R. Astron. Soc. 364, 1367–1377. doi: 10.1111/j.1365-2966.2005.09683.x

CrossRef Full Text | Google Scholar

Winkler, K. A., and Norman, M. L. (1986). “Astrophysical radiation hydrodynamics,” in NATO Advanced Research Workshop on Astrophysical Radiation Hydrodynamics (Dordrecht; Holland; Boston, MA: D. Reidel Publishing Co.), 187.

Winkler, K. H. A., and Newman, M. J. (1980). Formation of solar-type stars in spherical symmetry. I - The key role of the accretion shock. Astrophys. J. 236, 201–211. doi: 10.1086/157734

CrossRef Full Text | Google Scholar

Wise, J. H., and Abel, T. (2011). ENZO+MORAY: radiation hydrodynamics adaptive mesh refinement simulations with adaptive ray tracing. Month. Notices R. Astron. Soc. 414, 3458–3491. doi: 10.1111/j.1365-2966.2011.18646.x

CrossRef Full Text | Google Scholar

Wise, J. H., Turk, M. J., Norman, M. L., and Abel, T. (2012). The birth of a galaxy: primordial metal enrichment and stellar populations. Astrophys. J. 745:50. doi: 10.1088/0004-637X/745/1/50

CrossRef Full Text | Google Scholar

Wurster, J., Bate, M. R., and Price, D. J. (2018). On the origin of magnetic fields in stars. Month. Notices R. Astron. Soc. 481, 2450–2457. doi: 10.1093/mnras/sty2438

CrossRef Full Text | Google Scholar

Wurster, J., Bate, M. R., and Price, D. J. (2018). The collapse of a molecular cloud core to stellar densities using radiation non-ideal magnetohydrodynamics. Month. Notices R. Astron. Soc. 475, 1859–1880. doi: 10.1093/mnras/stx3339

CrossRef Full Text | Google Scholar

Wurster, J., and Li, Z.-Y. (2018). The role of magnetic fields in the formation of protostellar discs. Front. Astron. Space Sci. 5:39. doi: 10.3389/fspas.2018.00039

CrossRef Full Text | Google Scholar

Wurster, J., Price, D., and Ayliffe, B. (2014). Ambipolar diffusion in smoothed particle magnetohydrodynamics. Month. Notices R. Astron. Soc. 444, 1104–1112. doi: 10.1093/mnras/stu1524

CrossRef Full Text | Google Scholar

Wurster, J., Price, D. J., and Bate, M. R. (2016). Can non-ideal magnetohydrodynamics solve the magnetic braking catastrophe? Month. Notices R. Astron. Soc. 457, 1037–1061. doi: 10.1093/mnras/stw013

CrossRef Full Text | Google Scholar

Yorke, H. W., and Bodenheimer, P. (1999). The formation of protostellar disks. III. The influence of gravitationally induced angular momentum transport on disk structure and appearance. Astrophys. J. 525, 330–342. doi: 10.1086/307867

CrossRef Full Text | Google Scholar

Yorke, H. W., Yorke, H. W., Sonnhalter, C., and Sonnhalter, C. (2002). On the formation of massive stars. Astrophys. J. 569, 846–862. doi: 10.1086/339264

CrossRef Full Text | Google Scholar

Zamora-Avilés, M., Vázquez-Semadeni, E., Körtgen, B., Banerjee, R., and Hartmann, L. (2018). Magnetic suppression of turbulence and the star formation activity of molecular clouds. Month. Notices R. Astron. Soc. 474, 4824–4836. doi: 10.1093/mnras/stx3080

CrossRef Full Text | Google Scholar

Zhang, S., Hartmann, L., Zamora-Avilés, M., and Kuznetsova, A. (2018). On estimating angular momenta of infalling protostellar cores from observations. Month. Notices R. Astron. Soc. 480, 5495–5503. doi: 10.1093/mnras/sty2244

CrossRef Full Text | Google Scholar

Zhang, W., Howell, L., Almgren, A., Burrows, A., Dolence, J., and Bell, J. (2013). Castro: a new compressible astrophysical solver. III. Multigroup radiation hydrodynamics. Astrophys. J. Suppl. Ser. 204:7. doi: 10.1088/0067-0049/204/1/7

CrossRef Full Text | Google Scholar

Zhao, B., Caselli, P., and Li, Z.-Y. (2018a). Effect of grain size on differential desorption of volatile species and on non-ideal MHD diffusivity. Month. Notices R. Astron. Soc. 478, 2723–2736. doi: 10.1093/mnras/sty1165

CrossRef Full Text | Google Scholar

Zhao, B., Caselli, P., Li, Z.-Y., Krasnopolsky, R., Shang, H., and Nakamura, F. (2016). Protostellar disc formation enabled by removal of small dust grains. Month. Notices R. Astron. Soc. 460, 2050–2076. doi: 10.1093/mnras/stw1124

CrossRef Full Text | Google Scholar

Zhao, B., Caselli, P., Li, Z. Y., and Krasnopolsky, R. (2018b). Decoupling of magnetic fields in collapsing protostellar envelopes and disc formation and fragmentation. Month. Notices R. Astron. Soc. 473, 4868–4889. doi: 10.1093/mnras/stx2617

CrossRef Full Text | Google Scholar

Zhao, B., Li, Z.-Y., Nakamura, F., Krasnopolsky, R., and Shang, H. (2011). Magnetic flux expulsion in star formation. Astrophys. J. 742:10. doi: 10.1088/0004-637X/742/1/10

CrossRef Full Text | Google Scholar