WHERE ARE THE LOW-MASS POPULATION III STARS?

The cosmological simulation consists of 2048³ dark matter particles in a comoving box of 8 $\;{h}^{-1}\;\mathrm{Mpc}$. The mass resolution is 5.13 × 10³$\;{h}^{-1}\;{M}_{\odot }$ and the gravitational softening length is 120 $\;{h}^{-1}\;\mathrm{pc}$, which enable us to handle minihalos with sufficient resolution. We generated the initial condition by a publicly available code, 2LPTic,⁵ using second-order Lagrangian perturbation theory (e.g., Crocce et al. 2006). We used the online version⁶ of CAMB (Lewis et al. 2000) to calculate the transfer function. The cosmological parameters adopted are consistent with an observation of the cosmic microwave background obtained by the Planck satellite (Planck Collaboration et al. 2014), namely, Ω₀ = 0.31, Ω_b = 0.048, λ₀ = 0.69, h = 0.68, n_s = 0.96, and σ₈ = 0.83. The initial and final redshifts are 127 and 0.

For the time integration, we used a massively parallel TreePM code, GreeM (Ishiyama et al. 2009a, 2012), with the Phantom-GRAPE software accelerator⁷ (Nitadori et al. 2006; Tanikawa et al. 2012, 2013), on the Aterui supercomputer at the Center for Computational Astrophysics, CfCA, of the National Astronomical Observatory of Japan. The snapshots were stored at the redshifts so that the logarithmic interval ${\rm{\Delta }}\mathrm{log}(1+z)$ is 0.01. We carefully checked how this interval affects the quality of merger trees. The value adopted here yields well-converged merger rates.

To identify halos, we used the Friends-of-Friends (FoF) algorithm (Davis et al. 1985) with a linking parameter of b = 0.2. The smallest halo consists of 32 particles, which set the minimum FoF halo mass to be 1.6 × 10⁵$\;{h}^{-1}\;{M}_{\odot }$. We extracted merger trees by the algorithm described in Ishiyama et al. (2015).

In our Pop III star formation model, Pop III stars are assumed to form in minihalos where virial temperature T_vir exceeds a threshold ${T}_{{\rm{vir}}}^{{\rm{crit}}}$. In such halos, H₂ molecules form efficiently and the gas collapses to form stars via H₂ line cooling. We also put an upper bound of T_vir so that we exclude the atomic cooling halos. We set the upper bound to be 2000 K in this work. We calculated T_vir of the minihalos using a function proposed by Kitayama et al. (2001),

$$\begin{eqnarray}{T}_{{\rm{vir}}} & = & 9.09\times {10}^{3}\left(\displaystyle \frac{\mu }{0.59}\right){\left[\displaystyle \frac{M({z}_{{\rm{c}}})}{{10}^{9}{h}^{-1}{M}_{\odot }}\right]}^{2/3}\\ & & \times {\left[\displaystyle \frac{{{\rm{\Delta }}}_{{\rm{c}}}({z}_{{\rm{c}}})}{18{\pi }^{2}}\right]}^{1/3}(1+{z}_{{\rm{c}}})\;{\rm{K}},\end{eqnarray} \tag{ 1 }$$

where μ is the mean molecular weight in units of the proton mass, $M({z}_{{\rm{c}}})$ is the minihalo mass at the collapse redshift z_c, and ${{\rm{\Delta }}}_{{\rm{c}}}({z}_{{\rm{c}}})$ is the mean overdensity of collapsed halos. We set μ = 1.22, and used the mean overdensity according to the spherical collapse model (Bryan & Norman 1998). Other values are directly derived from the simulation.

Once Pop III stars form, UV radiation in the Lyman–Werner (LW) bands from the stars causes photodissociations of H₂ molecules and suppresses Pop III star formation in relatively lower mass halos. We used the criterion of the virial temperature under the LW background proposed by Machacek et al. (2001),

$$\begin{eqnarray}&&\left(\displaystyle \frac{{T}_{{\rm{vir}}}^{{\rm{crit}}}}{{\rm{1000}}\;{\rm{K}}}\right)=0.36{\left[{F}_{{\rm{LW}}}{({{\rm{\Omega }}}_{b}{h}^{2})}^{-1}{\left(\displaystyle \frac{1+z}{20}\right)}^{3/2}\right]}^{0.22},\end{eqnarray} \tag{ 2 }$$

where F_LW is the LW flux in a unit of 10⁻²¹ ergs⁻¹ cm⁻² Hz⁻¹. We assume a spatially uniform and time-dependent LW flux,

$$\begin{eqnarray}&&{F}_{{\rm{LW}}}=4\pi {J}_{{\rm{LW}}}=1.26\times {10}^{1.8(-1-\mathrm{tanh}(0.1(z-40)))},\end{eqnarray} \tag{ 3 }$$

which is a fitting function derived to reproduce cosmological reionization simulations (Ahn et al. 2012).

By tracking the merger tree in descending redshift, we listed halos with T^crit_vir < T_vir < 2000 K as minihalos where Pop III stars are assumed to be born with the formation redshift z_c. If any progenitors in the main branch of a halo are already listed, the halo is excluded from the list. The main branch of a halo is extracted by connecting the most massive progenitor of the most massive progenitor. This procedure ensures that gas in Pop III star-forming halos is metal free. We terminated the list at z = 10, when the reionization is assumed to start, since the photoheating accompanied by the reionization shut down the star formation in such minihalos. The upper limit of the virial temperature criterion is insensitive to the total number of minihalos. When we change the upper limit to 8000 K from 2000 K, the number of minihalos increases in only ∼1%.

The IMF of Pop III stars is still unknown, although there are some theoretical implications (e.g., Greif et al. 2011a; Hirano et al. 2014, 2015; Susa et al. 2014). As discussed in Section 1, the low-mass end is particularly uncertain even in theoretical models. Hence, we simply assume that the number ${n}_{{\rm{pop3}}}$ of low-mass Pop III stars with a main-sequence mass of 0.15–1.0 M_⊙ formed in a minihalo. Then we try to constrain the theoretical models in reverse.

We selected randomly ${n}_{{\rm{pop3}}}$ dark matter particles from each minihalo as tracers of the low-mass Pop III survivors. The spatial positions of the tracers at z = 0 are assumed to be those of Pop III survivors. In this study, we use ${n}_{{\rm{pop3}}}=1$ and 10. We adopt a Kroupa IMF (Kroupa 2001) for the low-mass Pop III stars and randomly set the mass of each Pop III star. Their magnitudes of various bands at z = 0 are calculated from their masses and ages using an isochrone model with Z = 0 for stars with M > 0.7 M_⊙ and Z = 0.0001 for stars with M < 0.7 M_⊙ (Girardi et al. 2000; Marigo et al. 2001).⁸ Since the lifetime of stars with masses larger than ∼0.8 M_⊙ is shorter than the cosmic time, these stars cannot survive and the number of Pop III survivors stars per minihalo is ∼0.9 ${n}_{{\rm{pop3}}}$.

The Pop III survivors in the Milky Way mix with other stellar components of the Milky Way. Thus, the field stars need to be excluded for the detection of the Pop III survivors. In order to evaluate their detectability in different fields, we extract the stellar distribution and metallicity distribution of the Milky Way using the online interface of the stellar population synthesis model for the Milky Way, Besançon⁹ (Robin et al. 2003). The effect of the ISM extinction in the disk was not taken into account.

We calculate the number of survivors and field stars and the sample size required to find one Pop III survivor in the following four different fields:

• 1.  
  high-latitude field: l = 0°–30°, b = 60°–90°;
• 2.  
  middle-latitude field: l = 0°–15°, b = 30°–60°;
• 3.  
  low-latitude field: l = 75–15°, b = 15°–30°; and
• 4.  
  central field: l = 0°–5°, b = 5°–75,

which are displayed in Figure 4. These four fields represent the galactic halo, disk, bulge, and center of the Milky Way. The required sample size is defined as the ratio of the number of field stars to that of survivors. Since most survivors have a V-band magnitude, m_V > 17, we do not count field stars brighter than this value.

Zoom In Zoom Out Reset image size

Figure 5 shows the required sample size in four Milky-Way-sized halos as a function of V-band limiting magnitude for four different target fields. The result of the ${n}_{{\rm{pop3}}}=10$ model is shown (solid curves). Dashed curves are the results of this model multiplied by a factor of 10, which mimics the ${n}_{{\rm{pop3}}}=1$ model. The dependence of the number of survivors on ${n}_{{\rm{pop3}}}$ shown in Equation (4) justifies this simple scaling.

Zoom In Zoom Out Reset image size

In spite of the larger number of survivors toward the galactic center (Figure 3), higher latitude fields require lower sample sizes. This is simply because the field stars are more concentrated at the low-latitude field than the survivors. The higher number density of field stars at lower latitude makes the detection of Pop III survivors less efficient. In the low-latitude and central fields, the sample size required is by a factor of 3 ∼ 5 and 30 ∼ 200 larger than those in middle- and high-latitude fields, respectively. For all halos, the high-latitude field with the direction to the galactic plane (l = 0°–30°, b = 60°–90°) are most efficient.

The required sample size is slightly reduced with a deep limiting magnitude. This is because the deep observation reaches the stars in the galactic halo. The fraction of survivors in the galactic halo is higher than the galactic disk since the disk star shows more centrally concentrated distribution than that of dark matter. Furthermore, the number of survivors per square degree is increasing as we observe deeper, as indicated in Figure 6, which gives the number density of survivors brighter than the given V-band apparent magnitude per square degree in H5 for four different target fields. However, the weak dependence on the limiting magnitude indicates that a wide and shallow survey detects more survivors than a deep and narrow survey with a given survey power, i.e., telescope and instrument, as long as the number of targets is sufficient.

Zoom In Zoom Out Reset image size

When we take into account the ISM extinction, the expected number of both survivors and field stars becomes smaller, in particular, for lower latitude fields. We also extract the stellar distribution of the Milky Way with the extinction using the Besançon (Robin et al. 2003).¹⁰ Figure 7 shows the effect of the ISM extinction on the number of field stars for the central and the low-latitude fields, respectively; there is little effect on the low-latitude field (l = 75–15°, b = 15°–30°). On the other hand, the number of stars at the bright end is largely reduced for the central field (l = 0°–5°, b = 5°–75°). Assuming that the extinction acts on survivors in the same manner as the field stars, the number density of bright survivors (m_V < 19) is reduced by a factor ∼2 in only the central field.

Zoom In Zoom Out Reset image size

The observer position does not change the total number of survivors in Milky-Way-sized halos but can influence those in these fields. To investigate the effect, we calculate the number of survivors seen from randomly placed 1000 observers on a spherical surface lying 8.5 kpc from the halo center, for high- and low-latitudes fields of H5, l = 0°–30°, b = 60°–90°, and l = 75–15°, b = 15°–30°. In Table 2, we summarize the median, 5%, 25%, 75%, and 95% of the sample size for these fields. Regardless of the V-band limiting magnitude m_V, the sample size of the high-latitude field (l = 0°–30°, b = 60°–90°) spread by a factor of nearly two and four from 25% to 75% and from 5% to 95%. Those of the low-latitude field (l = 75–15°, b = 15°–30°) spread by a factor of nearly 1.4 and 2 from 25% to 75% and from 5% to 95%. The 95% values in the high-latitude field are always larger than the 5% values in the low-latitude field.

Table 2.  Uncertainty of the Sample Size by Various Observational Positions

Area	mv	Median	5%	25%	75%	95%
l = 0°–30°, b = 60°–90°	22	4.9 × 103	2.3 × 103	3.3 × 103	6.5 × 103	8.8 × 103
l = 0°–30°, b = 60°–90°	20	7.8 × 103	4.4 × 103	5.7 × 103	1.1 × 104	1.6 × 104
l = 75–15°, b = 15°–30°	22	3.3 × 104	2.2 × 104	2.7 × 104	3.9 × 104	4.8 × 104
l = 75–15°, b = 15°–30°	20	6.0 × 104	4.0 × 104	5.1 × 104	7.4 × 104	9.4 × 104

Note. The median, 5%, 25%, 75%, and 95% of the sample size from 1000 randomly placed observers are shown. m_v is the V-band magnitude cutoff.

Download table as:  ASCIITypeset image

The spectroscopic metal-poor surveys are mainly conducted at high-latitude fields so far. The results of cosmological simulations can be directly compared with outcomes of the past and ongoing/planning surveys. Figure 5 also shows the sample size and limiting magnitude of the following surveys with crosses; the Hamburg/ESO survey (Christlieb et al. 2008), SEGUE (Yanny et al. 2009), and future prospects by PFS (Takada et al. 2014) and LAMOST (Li et al. 2015). If we assume that the Hamburg/ESO survey and SEGUE do not detect any survivors and take into account the incompleteness of follow-up moderate dispersion spectroscopy for candidates found by the Hamburg/ESO survey (Christlieb 2006), the formation model of Pop III stars with ${n}_{{\rm{pop3}}}=10$ is already excluded by these observations because the main targets of SEGUE are in relatively higher latitude fields (Figure 1 in Yanny et al. 2009).

There is a study that suggests that the surface pollution can enhance the surface metal abundance of Pop III survivors up to [Fe/H] ∼ −5 (Komiya et al. 2015), below which the Hamburg/ESO survey had found two stars (Christlieb et al. 2002; Frebel et al. 2005). Even if their claim is real, the outcomes of the past surveys still rule out the model with ${n}_{{\rm{pop3}}}=10$ and favor ${n}_{{\rm{pop3}}}=1$.

Recently, Skymapper performs a photometric metal-poor survey with a narrow-band filter with a bandpass corresponding to Ca II HK lines and broad-band filters (Keller et al. 2007). Skymapper photometrically excludes the metal-rich field stars and successfully found the most iron-deficient star (Keller et al. 2014) and an extremely metal-poor star in the bulge (Howes et al. 2015). Although the efficiency and completeness of the photometric classification have not been presented, the stars with [Fe/H] > −1.5 are largely excluded (Figure 2 in Howes et al. 2014). Assuming that the photometric classification could exclude 95% of stars with [Fe/H] > −1.5 and adopting the metallicity distribution of the Milky Way model, the numbers of field stars in the high-, middle-, low-latitude, and central fields are reduced by factors of 2.8, 3.8, 7.8, and 17, respectively. The photometric classification enhances the efficiency of follow-up spectroscopy by these factors. As a result, the required sample sizes are comparable in the high- and middle-latitude fields, while the required sample sizes in the low-latitude and central fields are still 2 times and 20 times larger than that of the high-latitude field, respectively. Practically, the most efficient field and depth depend on the number of fibers and the field of view of instruments for spectroscopy. If the photometric classification works well, while the high-latitude field is most efficient for instruments with low fiber density, e.g., HERMES (∼120 fibers deg⁻², Sheinis et al. 2014) for the GALAH survey (e.g., De Silva et al. 2015) and LAMOST (∼800 fibers deg⁻², Li et al. 2015), the middle-latitude field is most efficient for instruments with high fiber density, e.g., PFS (∼1800 fibers deg⁻², Takada et al. 2014). By performing spectroscopic observation of a million of stars after the photometric classification in the future, we can constrain the low-mass Pop III star IMF with unprecedented accuracy.

Only bright survivors are reachable in dwarf galaxies but the distribution of survivors is concentrated (Figure 3). Thus, the detection of survivors in dwarf galaxies could be more efficient than that in the Milky Way. In this subsection, we investigate the possibility of the detection of Pop III survivors in dwarf galaxies. However, the known missing satellite problem (e.g., Klypin et al. 1999; Moore et al. 1999; Ishiyama et al. 2009b) deters the direct prediction from the cosmological simulations and thus we evaluate the average number of Pop III survivors in the model with ${n}_{{\rm{pop3}}}=10$, brighter than given absolute magnitudes in a subhalo (left panel of Figure 8). The thresholds correspond to absolute magnitudes of the tip of red giant stars (M_V = 0), the red giant stars (M_V = 2), and the turn-off stars (M_V = 4). We identify subhalos in each halo with the ROCKSTAR phase space halo/subhalo finder (Behroozi et al. 2013).

Zoom In Zoom Out Reset image size

While the average number is larger for fainter thresholds, dependence of the average number on the distance from the galactic center is similar. The average number decreases with the distance because the distant subhalos were recently formed and do not satisfy the criteria of Pop III star formation (Section 2.2). Furthermore, the average number per subhalo is larger for more massive subhalos (right panel of Figure 8) because the more massive subhalos include more minihalos with Pop III star formation. The inverse of the average number per subhalo is the required number of dwarf galaxies to find one Pop III survivor, if subhalos hosting dwarf galaxies do not preferentially contain a larger number of survivors than all other subhalos. The required number of dwarf galaxies is estimated in less than 10 at <100 kpc for the tip of red giant stars (corresponding to m_V ∼ 20), if the dwarf galaxies have a common mass of ∼10⁷ M_⊙ (e.g., Strigari et al. 2008; Okamoto & Frenk 2009), although there is a debate (Hayashi & Chiba 2012).

The number of dwarf galaxies discovered so far is nearly 30 (e.g., Bechtol et al. 2015; Koposov et al. 2015). Although the extremely metal-poor stars with [Fe/H] < −3 have been discovered especially in ultra-faint dwarf galaxies (e.g., Frebel et al. 2014), no Pop III survivor has been detected. This constraint disfavors n_pop3 = 10 but is consistent with n_pop3 = 1.

Recently, the wide-field photometric surveys newly found dwarf galaxies including ultra-faint dwarf galaxies at <100 kpc (e.g., Belokurov et al. 2007; Laevens et al. 2015). In addition to the large average number of Pop III survivors in nearby dwarf galaxies, the faint threshold can be realized for nearby dwarf galaxies. Thus, the all-sky survey of nearby dwarf galaxies is highly demanded for the detection of Pop III survivors. The increasing number of nearby dwarf galaxies refines the constraint on n_pop3.

In contrast to the Pop III survivors in the Milky Way, the survivors in dwarf galaxies are as faint as m_V ∼ 20 and the numbers of the detectable Pop III survivors and dwarf galaxies strongly depend on the limiting magnitude. Therefore, the light collecting power is essential on the detection of survivors in dwarf galaxies. Although the number of stars that are bright enough for the high-dispersion spectroscopy with 8 m class telescopes is small, candidates of the survivors are reached with the narrow-band imaging and/or low/medium-dispersion spectroscopy with 8 m class telescopes, and will be good targets for the follow-up high-dispersion spectroscopy with 30 m class telescopes.

Hartwig et al. (2015) provided the expected number of Pop III survivors in the Milky Way. They reported that an unbiased survey of 4 × 10⁶ halo stars could impose a constraint on the low-mass end of the Pop III IMF. This value is more than two orders of magnitude larger than our ${n}_{{\rm{pop3}}}=10$ model (Figure 5). Even in the ${n}_{{\rm{pop3}}}=1$ model, more than an order-of-magnitude difference exists. They used a simple model to distinguish halo and bulge stars because they adopted a semi-analytic model based on merger trees extracted from extended Press–Schechter theory, which could not give the spatial distribution of halos.

In the model of Hartwig et al. (2015), there are physical processes we do not include in our model, dynamical heating of gas during mass accretion and mergers. However, they reported that the main mechanisms to suppress Pop III star formation are metal enrichment at low redshifts and the threshold mass of the minihalo at high redshifts, both of which are included in our model by different manners.

The main reason for such large difference is the adoption of the IMF of the Pop III stars. Whereas Hartwig et al. (2015) used a logarithmically flat IMF from 0.01 to 100 M_⊙, we adopted the Kroupa IMF from 0.15 to 1.0 M_⊙ for low-mass Pop III stars. Consequently, the number of low-mass Pop III stars per minihalo is also largely different. In their model, the average number was about 0.1 for an IMF in the mass range 0.65–100 M_⊙ (T. Hartwig 2015, private communication), which is nearly two orders of magnitude smaller than our fixed ${n}_{{\rm{pop3}}}=10$ model. This large difference is sufficient to explain large gaps in the estimation of the sample sizes needed to constrain the low-mass end of the IMF. From Figure 5, it is easy to infer that the large difference in the estimation of sample sizes would disappear when we used the ${n}_{{\rm{pop3}}}=0.1$ model, consistent with Hartwig et al. (2015). This means that our and their results qualitatively agree with each other although there is quantitative disagreement.

Why is the average number of Pop III stars per minihalo in both models so different? The decisive difference is that the model of Hartwig et al. (2015) used the lowest mass of the IMF as an arbitrary parameter and tried to constrain it whereas the number of Pop III stars is used in our model. In fact, the lower bound of the IMF 0.65 M_⊙ in their model was obtained to match the sample size of Hamburg/ESO survey.

The radial number density of survivors (Figure 2) distributes more concentrated than the dark matter and qualitatively agrees with Gao et al. (2010) (see also White & Springel 2000; Diemand et al. 2005; Tumlinson 2010), although the concentration of our results are less than that of Gao et al. (2010). This quantitative disagreement can be explained by the difference of cosmological parameters. The parameter σ₈ of Gao et al. (2010) is 0.9, while 0.83 is used in our simulation. The Press–Schechter theory predicts the number of halos of M ∼ 10⁶ M_⊙ at z = 25 in σ₈ = 0.9 by a factor ∼3 larger than σ₈ = 0.83. Because earlier formed progenitor halos tend to concentrate on the center of halos, larger σ₈ should result in a higher concentration of survivors.

One may imagine that the difference can be explained by the way to select tracers of the low-mass Pop III stars. Whereas we use randomly selected dark matter particles in minihalos, the most bound particles were used in Gao et al. (2010). Since Pop III stars should be born in the central dense regions of minihalos, they might stay there with a high probability. However, in such situations, Pop III stars could not be long-lived low-mass stars because of subsequent mass accretion. On the other hand, some of stars could be kicked away from the central regions of minihalos to regions with shallower potential, via the gravitational many-body interaction (e.g., Clark et al. 2011a, 2011b; Greif et al. 2011a, 2012; Smith et al. 2011; Umemura et al. 2012; Machida & Doi 2013; Susa 2013), As a consequence, they could be long-lived low-mass stars because of poor mass accretion, and easily stripped from minihalos by tides of larger halos. Our method to select Pop III tracers mimics this latter process.

To investigate how the distribution of survivors is sensitive to the choice of tracers, we selected the most bound particle from each minihalo as tracers and calculated the radial number density distributions of survivors. Compared with the randomly selected model (${n}_{{\rm{pop3}}}=1$), we confirm that both distributions agree well. This suggests that the criteria of Pop III forming minihalos and adopted cosmological parameters are more important for the distributions of survivors, and cosmological radiation hydrodynamical simulations with large volume are necessary to calibrate models.

In this study, we neglect some physical processes that can affect the number of survivors. We will take them into account for our model in future.

Surfaces of Pop III survivors could be polluted with metals (e.g., Yoshii 1981; Shigeyama et al. 2003; Frebel et al. 2009; Komiya et al. 2010, 2015). Frebel et al. (2009) studied the accretion of the ISM on a star in the galactic disk and concluded that it is generally negligible. On the other hand, Komiya et al. (2015) investigated the effect of the ISM accretion in minihalos by a semi-analytic model based on the hierarchical clustering scenario. They demonstrated that the pollution is effective in the early stage of the hierarchical formation of halos and the surface iron abundance of survivors could be enhanced up to [Fe/H] ∼ −5. The Fe abundance of hyper metal-poor stars (HMP) can be explained by this pollution scenario. Even if we take into account this pollution effect, the models with ${n}_{{\rm{pop3}}}=10$ are still ruled out because only two stars of [Fe/H] < −5 have been discovered so far. Thus, the main constraint of the Pop III star formation model in our study is not much affected by this effect as discussed in Section 4.1.

The abundance of minihalos could be suppressed by streaming velocity, that is, the supersonic relative velocity of baryon and dark matter that arises at the time of recombination (Tseliakhovich & Hirata 2010; Tseliakhovich et al. 2011). Recent numerical simulations including streaming velocity indicate that typical minihalo mass increases by a factor of three and Pop III star formation is delayed by Δz ∼ 4 (Greif et al. 2011b). On the other hand, Stacy et al. (2011) suggest that there is little effect on the gas evolution by the typical streaming velocity. Currently the effect of the streaming velocity is highly uncertain, and thus, it is needed to evaluate it accurately through large Pop III formation simulations.

In the present paper, we have discussed the observational possibility of finding Pop III survivors. The effort to search these stars will help to constrain the low-mass end of the Pop III IMF at M ≲ 1 M_⊙ severely. On the other hand, the theoretical IMF of Pop III stars extends to $\sim 1000\;{M}_{\odot }$ (Susa et al. 2014; Hirano et al. 2014, 2015). Thus, it is worth mentioning the observations to be compared with the high-mass part of the Pop III IMF.

Metal-poor stars could have been born in the remnants of Pop III stars. Hence the comparison of the abundance ratios in the atmospheres of these stars with the theoretical predictions of the nucleosynthesis in the Pop III stars could have great significance in constraining the high-mass part of the IMF. In fact, the observed abundance ratios provide little evidence of pair instability supernovae (PISNe), which should be found if some Pop III stars form in the range of $140\;{M}_{\odot }\lesssim M\lesssim 260\;{M}_{\odot }$. It is rather consistent with the assumption that most Pop III stars are less massive than 100 M_⊙ to supply the metals by core collapse supernovae (Susa et al. 2014). We have to keep in mind that we cannot directly conclude the less massive (≲100 M_⊙) IMF is favored, because the lack of PISNe pattern is only evident in very metal-poor stars with [Fe/H] ≦ −3, and has not been proved for [Fe/H] > −3. In fact, at [Fe/H] = −2.5 Aoki et al. (2014) found a possible candidate of a second-generation star formed in a remnant of a massive star with ≳100 M_⊙. Further observations of the stars with higher metallicity will give us more information on Pop III IMF. Present high-resolution simulations will be coupled with semi-analytical models of low-metallicity star formation to be compared with these observations in the near future.

Some theoretical calculations on the Pop III star formation naturally predict formation of multiple stellar systems including massive binaries (Stacy et al. 2012; Susa 2013; Stacy & Bromm 2014; Susa et al. 2014). Such systems will evolve into black hole binaries, which would merge to form a more massive black hole by emitting gravitational waves. Kinugawa et al. (2016) estimated the detection rate of such events to be pretty high (∼180 events yr⁻¹). If the gravitational wave from such objects are detected at a predicted rate, it will be circumstantial evidence that we are witnessing the merging of Pop III black holes. In fact, the recent discovery of the gravitational wave from 29 M_⊙ to 36 M_⊙ black hole binary coalescence (Abbott et al. 2016) suggests a high rate of such events. The forthcoming data release will provide a better estimate of the frequency. Meanwhile, the high-resolution cosmological simulations as presented in this paper will be coupled with Pop III binary formation/evolution theory to give a better prediction of the event rate from the theoretical side.