The cosmic dark ages ended a few hundred million years after the Big Bang, when the first stars began to fill the universe with new light. It has generally been argued that these stars formed in isolation and were extremely massive -perhaps 100 times as massive as the Sun. In a recent study, Clark and collaborators showed that this picture requires revision. They demonstrated that the accretion disks that build up around Population III stars are strongly susceptible to fragmentation and that the first stars should therefore form in clusters rather than in isolation. We here use a series of high-resolution hydrodynamical simulations performed with the moving mesh code AREPO to follow up on this proposal and to study the influence of environmental parameters on the level of fragmentation. We model the collapse of five independent minihalos from cosmological initial conditions, through the runaway condensation of their central gas clouds, to the formation of the first protostar, and beyond for a further 1000 years. During this latter accretion phase, we represent the optically thick regions of protostars by sink particles. Gas accumulates rapidly in the circumstellar disk around the first protostar, fragmenting vigorously to produce a small group of protostars. After an initial burst, gravitational instability recurs periodically, forming additional protostars with masses ranging from ∼ 0.1 to 10 M ⊙ . Although the shape, multiplicity, and normalization of the protostellar mass function depend on the details of the sink-particle algorithm, fragmentation into protostars with diverse masses occurs in all cases, confirming earlier reports of Population III stars forming in clusters. Depending on the efficiency of later accretion and merging, Population III stars may enter the main sequence in clusters and with much more diverse masses than are commonly assumed.

We investigate the formation of the first stars at the end of the cosmic dark ages with a suite of three‐dimensional, moving‐mesh simulations that directly resolve the collapse of the gas beyond the formation of the first protostar at the centre of a dark matter minihalo. The simulations cover more than 25 orders of magnitude in density and have a maximum spatial resolution of 0.05 R⊙, which extends well below the radius of individual protostars and captures their interaction with the surrounding gas. In analogy to previous studies that employed sink particles, we find that the Keplerian disc around the primary protostar fragments into a number of secondary protostars, which is facilitated by H2 collisional dissociation cooling and collision‐induced emission. The further evolution of the protostellar system is characterized by strong gravitational torques that transfer angular momentum between the secondary protostars formed in the disc and the surrounding gas. This leads to the migration of about half of the secondary protostars to the centre of the cloud in a free‐fall time, where they merge with the primary protostar and enhance its growth to about five times the mass of the second most massive protostar. By the same token, a fraction of the protostars obtain angular momentum from other protostars via N‐body interactions and migrate to higher orbits. On average, only every third protostar survives until the end of the simulation. However, the number of protostars present at any given time increases monotonically, suggesting that the system will continue to grow beyond the limited period of time simulated here.

The very first stars to form in the universe heralded an end to the cosmic dark ages and introduced new physical processes that shaped early cosmic evolution. Until now, it was thought that these stars lived short, solitary lives, with only one extremely massive star, or possibly a very wide binary system, forming in each dark matter minihalo. Here we describe numerical simulations that show that these stars were, to the contrary, often members of tight multiple systems. Our results show that the disks that formed around the first young stars were unstable to gravitational fragmentation, possibly producing small binary and higher-order systems that had separations as small as the distance between the Earth and the Sun.

We use the moving mesh code arepo coupled to a time-dependent chemical network to investigate molecular gas in simulated spiral galaxies that is not traced by CO emission. We calculate H 2 and CO column densities, and estimate the CO emission and CO-H 2 conversion factor. We find that in conditions akin to those in the local interstellar medium, around 42% of the total molecular mass should be in CO-dark regions, in reasonable agreement with observational estimates. This fraction is almost insensitive to the CO integrated intensity threshold used to discriminate between CObright and CO-dark gas. The CO-dark molecular gas primarily resides in extremely long (> 100 pc) filaments that are stretched between spiral arms by galactic shear. Only the centres of these filaments are bright in CO, suggesting that filamentary molecular clouds observed in the Milky Way may only be small parts of much larger structures. The CO-dark molecular gas mainly exists in a partially molecular phase which accounts for a significant fraction of the total disc mass budget. The dark gas fraction is higher in simulations with higher ambient UV fields or lower surface densities, implying that external galaxies with these conditions might have a greater proportion of dark gas.

We use a suite of high resolution molecular cloud simulations carried out with the moving mesh code arepo to explore the nature of star-forming filaments. The simulated filaments are identified and categorised from column density maps in the same manner as for recent Herschel observations. When fit with a Plummer-like profile the filaments are in excellent agreement with observations, and have shallow power-law profiles of p ∼ 2.2 without the need for magnetic support. When data within 1 pc of the filament centre is fitted with a Gaussian function, the average FWHM is ∼ 0.3 pc, in agreement with predictions for accreting filaments. However, if the fit is constructed using only the inner regions, as in Herschel observations, the resulting FWHM is only ∼ 0.2 pc. This value is larger than that measured in IC 5146 and Taurus, but is similar to that found in the Planck Galactic Cold Cores and in Cygnus X. The simulated filaments have a range of widths rather than a constant value. When the column density maps are compared to the 3D gas densities, the filaments seen in column density do not belong to a single structure. Instead, they are made up of a network of short ribbon-like sub-filaments reminiscent of those seen in Taurus. The sub-filaments are pre-existing within the simulated clouds, have radii similar to their Jeans radius, and are not primarily formed through fragmentation of the larger filament seen in column density. Instead, small filamentary clumps are swept together into a single column density structure by the large-scale collapse of the cloud.