A Systematic Study of the Stochastic Gravitational-Wave Background due to Stellar Core Collapse

TL;DR

This study models the stochastic gravitational-wave background from stellar core collapse using empirical and Gaussian fits, assessing detectability with future gravitational-wave detectors like Advanced LIGO and Einstein Telescope.

Contribution

It introduces a systematic approach to model the gravitational-wave background from core collapse events and evaluates their detectability with upcoming detectors.

Findings

Models are 2-4 orders below Advanced LIGO sensitivity.

Models are 1-2 orders below Einstein Telescope sensitivity.

Systematic parameter space study of gravitational-wave spectra.

Abstract

Stellar core collapse events are expected to produce gravitational waves via several mechanisms, most of which are not yet fully understood due to the current limitations in the numerical simulations of these events. In this paper, we begin with an empirical functional form that fits the gravitational-wave spectra from existing simulations of stellar core collapse and integrate over all collapse events in the universe to estimate the resulting stochastic gravitational-wave background. We then use a Gaussian functional form to separately fit and model a low-frequency peak in the core-collapse strain spectra, which likely occurs due to prompt convection. We systematically study the parameter space of both models, as well as the combined case, and investigate their detectability by upcoming gravitational-wave detectors, such as Advanced LIGO and Einstein Telescope. Assuming realistic formation rates for progenitors of core-collapse supernovae, our results indicate that both models are 2--4 orders of magnitude below the expected sensitivity of Advanced LIGO, and 1--2 orders of magnitude below that of the Einstein Telescope.

Tables

Table 1. Table 1: Comparison of sensitivities for ξ 𝜉 \xi and A ′ superscript 𝐴 ′ A^{\prime} for the combined model and the individual low- and high-frequency models. The sensitivities for the combined model are better as compared to the individual models due to the overall increase in SGWB energy density for the combined model.

Model	$\xi$ [m2/s]		${𝐀}^{\prime }$ [m/kg1/2]
	aLIGO	ET	aLIGO	ET
Low frequency	—	—	$1.4\times {10}^{-14}$	$1.4\times {10}^{-15}$
High frequency	$1.8\times {10}^{11}$	$3.8\times {10}^6$	—	—
Combined	$2.9\times {10}^7$	$2.2\times {10}^5$	$1.1\times {10}^{-14}$	$1.1\times {10}^{-15}$

Equations

$${\Omega }_{\mathrm{GW}}(f)$$

$$R(z)={\lambda }_{\text{CC}}{R}_{\ast }(z),$$

$$E({\Omega }_{\text{m}},{\Omega }_{\Lambda },z)=\sqrt{{\Omega }_{\text{m}}(1+z{)}^3+{\Omega }_{\Lambda }},$$

$${\Omega }_{\mathrm{GW}}(f)$$

$$R_{\ast }(z)=\nu \frac{p{e}^{q(z-z_m)}}{p-q+q{e}^{p(z-z_m)}}$$

$$f_e|\tilde{h}(f_e)|=\frac{G}{\pi c^{4}D}{E}_{\nu }⟨q⟩{\left(1+\frac{f_e}{a}\right)}^3{e}^{-f_e/b},$$

$$\frac{dE}{d{f}_e}(f_e)$$

$$\begin{array}{rl}{\Omega }_{\mathrm{GW}}(f)& =\frac{8\pi Gf\xi }{3H_0^3{c}^2}\text{bigintssss}dz\frac{R_{\ast }(z)}{(1+z)E({\Omega }_{\text{m}},{\Omega }_{\Lambda },z)}\\ & {\left(1+\frac{f(1+z)}{a}\right)}^6{e}^{-2f(1+z)/b},\end{array}$$

$$\xi =\frac{G{\lambda }_{\mathrm{CC}}}{c^5}{E}_{\nu }^2⟨q{⟩}^2$$

$$\text{SNR}=\frac{3H_0^2}{10{\pi }^2}\sqrt{2T}{\left[{\int }_0^{\infty }df\frac{{\gamma }^2(f){\Omega }_{\mathrm{GW}}^2(f)}{f^6{P}_1(f)P_2(f)}\right]}^{1/2},$$

$$\displaystyle f_{e}|\tilde{h}(f_{e})|_{\text{low}}=\frac{A}{D}\exp{\bigg{(}-\frac{(f_{e}-\mu)^{2}}{2\sigma^{2}}\bigg{)}}$$

$$\displaystyle\frac{dE}{df_{e}}(f_{e})_{\text{low}}=\frac{\pi^{2}c^{3}A^{2}}{G}\exp{\bigg{(}-\frac{(f_{e}-\mu)^{2}}{2\sigma^{2}}\bigg{)}}$$

$${\Omega }_{\text{GW}}(f)=\frac{8{\pi }^{3}c}{3H_0^2}f{A}^{\prime 2}\int dz\frac{R_{\ast }(z)}{(1+z)E({\Omega }_{\text{m}},{\Omega }_{\Lambda },z)}$$

$$\exp \left(-\frac{(f(1+z)-\mu {)}^2}{{\sigma }^2}\right)$$

$$\displaystyle\begin{split}\frac{dE}{df_{e}}(f_{e})_{\text{comb}}&=\frac{\xi}{\lambda_{\text{CC}}}\Bigg{[}\pi A^{\prime}\sqrt{\frac{c^{3}}{\xi G}}e^{-\frac{(f_{e}-\mu)^{2}}{2\sigma^{2}}}+\\ &\phantom{=}\left(1+\frac{f_{e}}{a}\right)^{3}e^{-f_{e}/b}\Bigg{]}^{2}.\end{split}$$

$${\Omega }_{\text{GW}}(f)$$

$$\mathcal{L}\propto \prod _i\exp \left[-\frac{(Y_i-{\Omega }_{\text{GW},i}(\xi ,a,b){)}^2}{2{\sigma }_i^2}\right],$$

Full text

I Introduction

The superposition of gravitational waves (GWs) generated by many uncorrelated and unresolved sources in the universe leads to a stochastic gravitational-wave background (SGWB). Such a background could be of cosmological or astrophysical origin. A cosmological SGWB may be produced by cosmic (super)strings Caldwell and Allen (1992); Damour and Vilenkin (2000, 2005); Siemens et al. (2007); Olmez et al. (2010, 2012), inflation Grishchuk (1975); Starobinskii (1979); Easther and Lim (2006); Barnaby et al. (2012), primordial black hole binaries Bird et al. (2016); Sasaki et al. (2016); Mandic et al. (2016), alternative cosmologies Gasperini and Veneziano (1993); Mandic and Buonanno (2006), and a variety of other phenomena. Astrophysical SGWB models include coalescences of compact objects in binary systems, like binary neutron stars or binary black holes Phinney (1991); Kosenko and Postnov (1998); Regimbau and de Freitas Pacheco (2006a); Zhu et al. (2011); Rosado (2011); Marassi et al. (2011a); Wu et al. (2012); Zhu et al. (2013); Evangelista and de Araujo (2013); Kowalska-Leszczynska et al. (2015); Evangelista and de Araujo (2015), rotating neutron stars Regimbau and de Freitas Pacheco (2001); Owen et al. (1998); Ferrari et al. (1999a); Chandrasekhar (1969); Houser et al. (1994); Lai and Shapiro (1995); Howell et al. (2004); Zhu et al. (2011); Rosado (2012); Lasky et al. (2013), magnetars Cutler (2002); Regimbau and de Freitas Pacheco (2006b); Regimbau and Mandic (2008); Marassi et al. (2011b); Howell et al. (2011); Wu et al. (2013); Lasky et al. (2013), the first stars Sandick et al. (2006), and white dwarf binaries Farmer and Phinney (2003).

Some of the proposed models have been constrained Siemens et al. (2007); Wu et al. (2012); Mandic et al. (2012); Crowder et al. (2013) by searches for the isotropic Abbott et al. (2005, 2007a, 2009a); Aasi et al. (2014, 2015) and anisotropic SGWB Abbott et al. (2007b, 2011). These searches have established upper limits on the energy density in the SGWB using data acquired by the first generation interferometric gravitational-wave detectors LIGO Abbott et al. (2004, 2009b) and Virgo Acernese et al. (2006). The second-generation Advanced LIGO (aLIGO) detectors have recently conducted their first observational run Abbott et al. (2016a); G. M. Harry (2010) (for the LIGO Scientific Collaboration), detecting the first gravitational-wave signals due to mergers of binary black hole systems Abbott et al. (2016b, c), and producing new estimates of the SGWB due to binary black hole mergers Abbott et al. (2016d). Other second-generation detectors, including Advanced Virgo aVi , GEO-HF Willke et al. (2006); Affeldt et al. (2014), and KAGRA Y. Aso and Y. Michimura and K. Somiya and M. Ando and O. Miyakawa and T. Sekiguchi and D. Tatsumi and and H. Yamamoto (2013); K. Somiya (2012) (for the KAGRA collaboration), are expected to become operational in the coming years. The third-generation gravitational-wave detectors, including the Einstein Telescope The Einstein Telescope science team (2011) and Cosmic Explorer Abbott et al. (2016e), are also currently being conceptualized.

Gravitational waves from stellar core collapse are expected to be produced via several mechanisms, including a quasi-periodic signal generated during the post-shock convection phase, hot-bubble convection and the standing accretion shock instability (SASI) Murphy et al. (2009); Marek et al. (2009); Kotake et al. (2009), anisotropic neutrino emission Marek et al. (2009); Kotake et al. (2009); Ott (2009), and the ringdown of the potentially newly formed black hole. Even though full three-dimensional simulations that include a complete set of relevant physical processes are not yet computationally feasible, predictions have been made about the gravitational-wave signals emitted during core collapse Ott et al. (2011, 2013); Müller et al. (2013); Yakunin et al. (2015). These predictions have been used to compute estimates of the corresponding SGWB from stellar core collapse due to both standard and early (Population III) stars Ferrari et al. (1999b); de Araujo et al. (2000, 2002); de Araujo et al. (2002); de Araujo et al. (2004a, b); Buonanno et al. (2005); Coward et al. (2002); Zhu et al. (2010); Marassi et al. (2009); Sandick et al. (2006); Pacucci et al. (2015). However, the dependency of the gravitational-wave signal on stellar progenitor properties, such as mass or spin, is not well known, and the rate of core collapse events is similarly uncertain. This means that the resulting SGWB estimates are necessarily only approximate.

In this paper, we try to take a more systematic look at the stellar core collapse SGWB model. We begin with an empirical functional form that has been used in the past to describe the gravitational-wave spectrum emitted by a single stellar core collapse. We fit this functional form to the available GW spectra obtained in numerical simulations of the core collapse process to determine the plausible range of the free parameters. In the process of studying these simulations, we observe that many simulations predict a low-frequency peak that our functional form is unable to capture. As a result, we model this feature using a second functional form. For both of these models, we compute the SGWB due to all stellar core collapse events in the universe, assuming that the redshift distribution of core collapse events follows the star formation rate. We use the star formation rate model from Vangioni et al. (2015), which is shown to agree with observations of metallicity, optical depth, and the reionization redshift, as obtained from CMB measurements. We then scan the parameter space and determine the regions that will be accessible to Advanced LIGO (at design sensitivity) and Einstein Telescope detectors as the representative second and third generation detectors, respectively. We note that this study is a follow up on the study presented in Crocker et al. (2015), which focused only on gravitational waves produced by the ringdown of the black hole formed at the end of the core collapse.

In Section 2, we present a model of the SGWB due to stellar core collapse, identify the plausible ranges of the model parameter space, and discuss the low-frequency feature in the spectrum. In Section 3, we perform a systematic scan of the parameter space of this SGWB model and discuss its accessibility to the future detectors. We summarize our results in Section 4.

II SGWB Model due to Stellar Core Collapse

Following Regimbau (2011); Regimbau and Mandic (2008); Wu et al. (2012); Abbott et al. (2016d), we define the GW energy density normalized by the critical energy density:

[TABLE]

where $\rho_{\rm GW}(f)$ is the energy density in gravitational waves in the frequency band $(f,f+df)$. The critical energy density is the energy density necessary to close the universe, $\rho_{c}=\frac{3H_{0}^{2}c^{2}}{8\pi G}$, where $H_{0}$ is the current value of the Hubble constant, taken to be 67.7 km/s/Mpc, $G$ is the gravitational constant, and $c$ is the speed of light. Then, following Abbott et al. (2016d), we have:

[TABLE]

where the energy spectrum emitted by a single astrophysical source is denoted $dE/df_{e}(f_{e})$ and it is evaluated at the source (emitted) frequency $f_{e}$ which is related to the observed frequency as $f_{e}=f(1+z)$. The rate of stellar core collapse events is denoted $R(z)$ and it can be expressed in terms of the star formation rate (SFR) $R_{*}(z)$:

[TABLE]

where $\lambda_{\text{CC}}$ denotes the mass fraction of stars which undergo the core-collapse process. In our analysis, we will treat $\lambda_{\text{CC}}$ as a free parameter, but we note that it can be estimated from the initial mass function (IMF). Assuming the Salpeter IMF $\phi(m)=Nm^{-2.35}$ normalized such that $\int_{0.1M_{\odot}}^{\infty}\phi(m)\,m\,dm=1$, we can estimate $\lambda_{\text{CC}}$ by assuming that stars with masses larger than $8M_{\odot}$ undergo core collapse: $\lambda_{\text{CC}}=\int_{8M_{\odot}}^{\infty}\phi(m)dm\approx 0.007M_{\odot}^{-1}$.

The factor of $1+z$ in Eq. 2 corrects for the cosmic expansion, and converts the time in the source frame to the detector frame. Finally, $E(\Omega_{\text{m}},\Omega_{\Lambda},z)$ in Eq. 2 captures the redshift dependence of the comoving volume:

[TABLE]

where $\Omega_{\text{m}}=0.309$ and $\Omega_{\Lambda}=0.691$ denote, respectively, the energy density in matter and in dark energy Ade et al. (2015). Combining the above, we have

[TABLE]

Finally, we comment on our choice of the star formation rate model. The SFR and its redshift dependence have been examined by multiple authors, who have proposed several functional forms Lilly et al. (1996); Hopkins and Beacom (2006); Fardal et al. (2007); Wilkins et al. (2008); Nagamine et al. (2006); Hernquist and Springel (2003); Vangioni et al. (2015). We use the functional form proposed by Vangioni et al. (2015) that is based on measurements of the luminosity function in high-redshift galaxies Behroozi et al. (2013). As noted in Vangioni et al. (2015), this SFR model is consistent with metallicity observations, as well as with the optical depth and reionization redshift inferred from cosmic microwave background observations Hinshaw et al. (2013); Ade et al. (2015). This SFR model is given by the Springel & Hernquist functional form Hernquist and Springel (2003):

[TABLE]

with the following parameters Vangioni et al. (2015): $\nu=0.178$ M${}_{\odot}/{\rm yr}/{\rm Mpc}^{3}$, $z_{m}=2.00$, $p=2.37$, and $q=1.80$. We note that at low redshifts, this SFR model agrees well with a model based on GRB observations at high redshift Vangioni et al. (2015).At high redshifts, the difference is more substantial, but as noted in Crocker et al. (2015), the contribution to the SGWB from high-redshift stellar core collapse events is not dominant, and therefore the choice between these two SFR models does not make a significant difference in the resulting SGWB. As a result, we continue our analysis assuming the luminosity function-based SFR model, since it matches CMB observations more closely in terms of optical depth and redshift of reionization.

II.1 High Frequency Model

As noted above, the physics of the stellar core collapse process is not yet fully understood, and three-dimensional numerical simulations which include the complete physical description of the process are yet to be conducted. Past simulations have primarily been in two dimensions Müller et al. (2013); although some are three-dimensional, these simulations typically use a coarser physical description Ott et al. (2013). In Buonanno et al. (2005); Sandick et al. (2006), it was shown that the following functional form could describe the GW amplitude spectra $\tilde{h}(f_{e})$ emitted during the core collapse process (in the local frame of the star):

[TABLE]

where $a$ and $b$ are free parameters of the model, $D$ is the distance to the star (assumed small enough that redshifting effects can be ignored), $E_{\nu}$ is the energy carried away by neutrinos during the core collapse and $\langle q\rangle$ is the luminosity-weighted averaged neutrino anisotropy Buonanno et al. (2005).

The GW energy spectrum from a single core collapse event can then be computed as:

[TABLE]

Inserting this spectrum into Eq. LABEL:eqomega, we obtain the explicit form of the SGWB spectrum:

[TABLE]

where $\xi$ is a combination of unknown scaling factors, defined as

[TABLE]

In anticipation of the low-frequency peak which is discussed and modeled below, we will refer to this as the high-frequency model. Figure 1 shows examples of spectra generated with this model. See the caption for discussion of the effect of the choice of model parameters on the morphology of the spectrum.

Figure 1 also shows the SGWB sensitivities of the Advanced LIGO aLI and Einstein Telescope The Einstein Telescope science team (2011) detectors, assuming the detection statistic defined by Allen and Romano (1999) for which the signal to noise ratio is defined as:

[TABLE]

where $T$ is the observation time (set to 1 year in our case), $\gamma(f)$ is the overlap reduction function for the chosen detector pair (set to 1, assuming colocated detector pairs) arising from the different locations and orientations of the detectors Allen and Romano (1999), and $P_{1}(f)$ and $P_{2}(f)$ are the strain power spectral densities of the two detectors. In Figure 1, we plot the $\text{SNR}=2$ curves.

In Figure 2, we show the effect of $a$ and $b$ on the spectrum by plotting $\Omega_{\text{GW}}/\xi$ at 100 Hz as a function of these parameters. For a fixed value of $b$, increasing $a$ leads to decreased $\Omega_{\text{GW}}$; this is apparent from the fact that $\Omega_{\text{GW}}$ goes approximately as $a^{-6}$ (see Eq. 9). The converse is true for $a$: fixing $a$ and allowing $b$ to grow leads to higher $\Omega_{\text{GW}}$, since increasing $b$ pushes the exponential term in Eq. 9 closer to 1.

Similar results were found in Sandick et al. (2006) and Dvorkin et al. (2016). In Sandick et al. (2006), two modes of star formation were considered: a normal mode of star formation as considered here, and several possibilities for an additional population of massive stars (Population III) which were deemed necessary as initial reports of the optical depth were quite high. Using $a=200$ and $b=300$, they found a peak value of $\Omega_{\rm GW}h^{2}=3-4\times 10^{-10}$ at $f\approx 300$ Hz. In general, the massive modes made a relatively small contribution (no more than a factor of 30) at lower frequencies ($f<100$ Hz). Above 300 Hz, all models were indistinguishable. More recently, using the same SFR considered here, Ref. Dvorkin et al. (2016) compared different models of black hole formation from single star collapse and mergers. Single star collapse also showed peak values of $\Omega_{\rm GW}h^{2}$ between $10^{-10}$ and $10^{-9}$, though the peak frequency was model dependent.

In order to perform a scan of the parameter space of this model, it is necessary to determine the plausible ranges of the parameters $a$ and $b$. Since these are empirical parameters, we determine the ranges of their values by fitting the functional form from Eq. 7 to the spectra generated in various simulations of the core collapse. We performed fits systematically to all waveforms provided by Ott et al. (2013); Abdikamalov et al. (2014); Yakunin et al. (2015); Müller et al. (2013) and show below several examples to illustrate the quality of the fits and the range of the fitted parameters. We limit the fits to the frequencies below $f_{e}\approx 400$ Hz: since the SGWB will be dominated by progenitors at redshifts $z\lesssim 3$, imposing the 400 Hz cut-off in the emitted frequency will impact the spectrum only at observed frequencies $\gtrsim 100$ Hz, where the sensitivity of terrestrial GW detectors is significantly reduced (c.f. Figure 1).

The first set of simulated spectra comes from a study performed by Ott et al. Ott et al. (2013). They perform multiple simulations in which a $27{\rm\;M_{\odot}}$ star is made to undergo core collapse and examine the post-core bounce phase. These simulations are carried out in three dimensions, use general relativity, and apply a neutrino leakage scheme that is consistent with radiation-hydrodynamics simulations. The neutrino heating rate and simulation duration are varied in each simulation Ott et al. (2013) (see Table 1 in Ott et al. (2013) for details regarding each simulation). Figure 3 shows $f|\tilde{h}(f)|$ from a cross-polarized gravitational-wave signal, along with a fit of our functional form. This signal was generated from a 184 ms simulation which used a neutrino heating parameter of 1.00 (see Ott et al. (2013) for more detail).

We also use simulations performed by Abdikamalov et al. Abdikamalov et al. (2014), which examine the angular momentum dependence of GW signals coming from stellar core collapse. In particular, the collapse of a $12{\rm\;M_{\odot}}$ progenitor model described in Woosley and Heger (2007) is carried out in two dimensions under the assumption of axisymmetry and with the use of the general relativistic hydrodynamics code CoCoNuT Müller et al. (2010). The simulation also uses an approximate electron capture model, the Lattimer-Swesty equation of state (EOS) Lattimer and Swesty (1991) with bulk modulus $K=220$ MeV, and a neutrino leakage scheme for postbounce evolution. The simulation is carried out many times, with the angular momentum in the core varied in each iteration (see Tables I and II in Abdikamalov et al. (2014) for details). Figure 4 shows the strain spectrum from a plus-polarized gravitational wave signal, rescaled by observer distance $D$, for simulation A4O06.5 (see Table II in Abdikamalov et al. (2014)), along with our fit.

Next, we use spectra from simulations performed by Yakunin et al. Yakunin et al. (2015). These simulations use the CHIMERA code Bruenn et al. (2016) to study the collapse of four Woosley-Heger non-rotating stellar progenitors between 12 M*⊙* and 25 M*⊙* Woosley and Heger (2007). These simulations are in two dimensions, use multifrequency neutrino transport in ray-by-ray approximation with relativistic corrections, and a Lattimer-Swesty EOS Lattimer and Swesty (1991) with $K=220$ MeV Yakunin et al. (2015). Figure 5 shows the strain spectrum from this simulation carried out for the 20 M*⊙* progenitor, as well as our fit to these data.

Finally, we use spectra from simulations performed by Müller et al. Müller et al. (2013). These simulations were carried out using the Vertex-CoCoNuT Müller et al. (2010) or Vertex-Prometheus Buras et al. (2006); Rampp and H.-T.Janka (2002) codes. These are similar neutrino hydrodynamics codes, but Vertex-CoCoNuT uses general relativity while Vertex-Prometheus is pseudo-Newtonian. The simulations model core-collapse supernovae in two dimensions, and include neutrino transport using the ray-by-ray plus approximation Buras et al. (2006) and a variable Eddington factor technique via Vertex, the neutrino transport module Rampp and H.-T.Janka (2002). Six zero-age main sequence progenitor masses were assumed, from $8.1{\rm\;M_{\odot}}$ to $27{\rm\;M_{\odot}}$. All models use the Lattimer-Swesty EOS Lattimer and Swesty (1991), three with $K=220$ MeV and three with $K=180$ MeV Müller et al. (2013). Figure 6 shows the strain spectrum from the $8.1{\rm\;M_{\odot}}$ progenitor evolved using GR and an EOS with $K=180$ MeV, as well as our fit to this data.

Studies of these simulated GW spectra suggest that the following parameter ranges encapsulate most of the variation in the GW strain spectra:

• •

$5{\rm\;Hz}<a<150{\rm\;Hz}$

• •

$10{\rm\;Hz}<b<400{\rm\;Hz}$.

The fits are qualitative, but generally fall within a factor of two of the simulation data. The scaling parameter $\xi$ is not restricted by these fits since it is related to the overall amplitude of the SGWB energy density spectrum rather than to the individual source strain spectrum morphology.

II.2 Low-Frequency Model

The spectra shown in Figures 3-6 were fitted well with the functional form of Eq. 7. However, we observed that in some simulations the emitted GW spectra feature two peaks, one at high frequencies, which is well-modeled by Eq. 7, and another at $\approx$ 60–120 Hz. Of the 64 numerical core collapse simulations we considered, over 80% exhibited signs of the low-frequency peak, some more prominently than others.

The precise physical process that produces this low-frequency peak in the GW spectrum is unclear, but several papers have made note of it and attempted to identify its source. Müller et al. Müller et al. (2013) and Yakunin et al. Yakunin et al. (2015) consider the deflection of infalling matter through the shock as a major contributor to this peak, but ultimately rule it out. Müller et al. Müller et al. (2013) go on to indicate that a quasi-periodic signal around $100$ Hz arises due to prompt post-shock convection during the first tens of milliseconds following core bounce. However, they rule out convection as the direct source of this signal because the entropy and lepton gradients which drive the convection quickly disperse, while the GW emission lasts for tens of milliseconds; instead, they speculate that the convection leads to the development of acoustic waves which produce a GW signal due to the resulting mass motions. Finally, they note that this signal always peaks near $100$ Hz in their simulations; there appears to be no significant dependence on the properties of the progenitor.

Ott et al. Ott et al. (2013) found that early GW emission due to prompt convection sets in at $\approx 10$ ms after core bounce and lasts for about $30$ ms. They also note that the GW signal produced by prompt convection is particularly sensitive to the perturbations which drive the convection; this is compatible with the statement from Müller et al. Müller et al. (2013) that the signal from prompt convection may depend on the width of the layer in which the convection occurs.

Kuroda et al. Kuroda et al. (2014) found that the speed of core rotation during the collapse process may impact the frequency of the signal from prompt convection. Simulations in which the core was rapidly rotating led to the emission of spiral waves; the acoustic waves which produce the low-frequency GW signal were Doppler shifted by emission on top of the spiral waves, leading to a low-frequency peak closer to 200 Hz. They also note that accurate neutrino transport is a necessary component of these simulations in order to have a reliable prediction of the GW signal; this agrees with similar findings from Müller et al. Müller et al. (2013) which found significantly reduced GW emission at low frequencies when simplified neutrino rates were used.

Finally, recent 3D simulations by Andresen et al. Andresen et al. (2016) indicate the presence of an additional low-frequency component of the GW signal, which appears to be related to SASI. They speculate that the convection-related GW signal in the 100–200 Hz range may be actually driven by the coupling of SASI dynamics to the stellar interior. They also note that the signals from their 3D simulations are considerably different from those of 2D simulations, and amplitudes are significantly lower, especially at higher frequencies.

To summarize, there is significant evidence that the low-frequency peak is a realistic feature of the GW signal produced during core-collapse; it is present in many different simulations and there is a general agreement that the signal is related to prompt postbounce convection or SASI. Although it seems to vary for the different simulations, the peak in the GW frequency spectrum tends to occur between $\approx$60–120 Hz, or even as high as 200 Hz for the rapidly-rotating models of Kuroda et al. (2014). We also note that not all of the simulations included in our study exhibited this peak; however, it seems that the more modern simulations, which feature general relativistic hydrodynamics and sophisticated neutrino transport, generally include this feature Prestegard (2016).

Figure 7 shows that the functional form of Eq. 7 is not sufficient to capture the low-frequency behavior of the simulated spectrum. This issue is particularly important because the terrestrial GW detectors have optimal sensitivities below 100 Hz, as shown in Figure 1. We model the low-frequency behavior with a Gaussian function in the strain spectrum governed by three parameters. The mean, denoted $\mu$, defines the center of the Gaussian peak, the standard deviation, denoted $\sigma$, defines the width of the peak, and $A$ denotes the amplitude of this peak:

[TABLE]

where $A$ is scaled by $D$, the distance from the star at which the signal is observed (assumed to be sufficiently small for redshifting effects to be neglected). The corresponding GW energy spectrum is given by:

[TABLE]

which gives GW energy density:

[TABLE]

This equation combines the two scaling parameters $A$ and $\lambda_{\text{CC}}$ into $A^{\prime}=A\lambda_{\text{CC}}^{1/2}$, with units of m/kg1/2 Prestegard (2016). Figure 8 shows this model for several choices of parameters.

II.3 Combined High and Low-Frequency Model

We also consider a combined case, which includes both the low-frequency peak and the higher-frequency functional form for the core-collapse spectrum. Following Prestegard (2016), this results in an energy spectrum of

[TABLE]

[TABLE]

III Results

III.1 High Frequency Model

As discussed above, the free parameters of the high-frequency SGWB model are $\xi$, $a$, and $b$. The $\xi$ parameter is simply a scaling factor related to the total energy emitted in GWs, and the $a$ and $b$ parameters effectively determine the shape of the spectrum: its rise, peak position, and drop-off at high frequencies. We scan the parameter space $(a,b,\xi)$ of this model, restricting $5{\rm\;Hz}<a<150{\rm\;Hz}$ and $10{\rm\;Hz}<b<400{\rm\;Hz}$ as discussed in the previous section. The range of $\xi$ is tuned so as to optimally probe the accessibility of the model to each detector pair. For Advanced LIGO, we use a range of $\xi=10^{6}\text{--}10^{13}{\rm\;m^{2}/s}$, and for ET we use $\xi=10^{2}\text{--}10^{9}{\rm\;m^{2}/s}$ Prestegard (2016). For each point in this parameter space, we compute the spectrum $\Omega_{\rm GW}(f)$ and compare it to the sensitivities of the Advanced LIGO and ET detectors, computing the likelihood function:

[TABLE]

where the index $i$ runs over frequency bins, $Y_{i}$ is the expected measurement of the GW energy density in the bin $i$, $\sigma_{i}$ is the corresponding measurement error, and $\Omega_{{\rm GW},i}(\xi,a,b)$ is the modelled energy density in the bin $i$ for the given free parameters $\xi$, $a$, and $b$. For projecting the future experimental sensitivities we set $Y_{i}=0$. In order to determine the accessibility of this three-dimensional parameter space to future detectors, we marginalize (integrate) the likelihood function over one of the parameters, and then compute the expected sensitivity contours at 95% confidence in the plane of the remaining two parameters.

Figure 9 shows the resulting 95% expected sensivitity contours for both Advanced LIGO and Einstein Telescope in the $\xi$-$b$ and $\xi$-$a$ planes. The curves in the $\xi$-$b$ plane are decreasing with $b$, which is a consequence of the fact that increasing $b$ cuts off the spectrum at higher frequencies, hence making the model more detectable by the GW detectors. Similarly, lowering the value of $a$ allows the spectrum to rise at lower frequencies (c.f. Figure 1), making the model more accessible to GW detectors and causing the increasing trend (with $a$) of the curves in the $\xi$-$a$ plane. It is also evident that the Einstein Telescope will provide a substantially better probe of this model than the second-generation detectors.

Marginalizing over the $a$ and $b$ parameters yields the posterior distribution for the amplitude parameter $\xi$. The 95% expected sensitivity on $\xi$ is $1.8\times 10^{11}$ m2/s and $3.8\times 10^{6}$ m2/s for Advanced LIGO and Einstein Telescope, respectively. For the representative value of the mass fraction parameter $\lambda_{\text{CC}}\approx 0.01$ M${}_{\odot}^{-1}$, these sensitivities can be translated into sensitivity curves in the $E_{\nu}$-$\langle q\rangle$ plane, also shown in Figure 9 (bottom). For a neutrino asymmetry of $\langle q\rangle=0.0045$ (following Sandick et al. (2006); Buonanno et al. (2005)), we see that Advanced LIGO would require the total neutrino energy to be approximately $2.5\times 10^{57}$ ergs for a detectable SGWB, which is rather high. Einstein Telescope would require $E_{\nu}\approx 1.2\times 10^{55}$ ergs, which is still about two orders of magnitude above expected neutrino energies from the literature ($\approx 3\times 10^{53}$ ergs Sandick et al. (2006); Buonanno et al. (2005)). It is important to note, however, that these sensitivity numbers are obtained after integrating over a large parameter space. For example, it is evident from the top two panels of Figure 9 that there is a substantial part of the parameter space with $\xi<3.8\times 10^{6}$ m2/s that would still be accessible to Einstein Telescope, hence potentially including the required neutrino energy level.

III.2 Low Frequency Model

The ranges for the low frequency model parameters are based on the previously described work with the individual simulations. For most of the simulations, the low frequency peak occurs between 60–120 Hz; however, there are some cases described in the literature where this peak may occur as high as 200 Hz due to rapid rotation during the collapse process Kuroda et al. (2014). Thus, we have used a liberally defined range on $\mu$, 30–200 Hz. Most of the waveforms that we studied from Ott et al. (2013) had relatively narrow low-frequency peaks, with typical widths between 40–60 Hz. The low-frequency peaks in the waveforms from Abdikamalov et al. (2014) tended to be wider, with some as broad as $140$ Hz. We assume the range of $\sigma$ to be 10–80 Hz. Because $A^{\prime}$ is a scaling parameter like $\xi$, its scanned range was tuned to each detector pair: $A^{\prime}=10^{-15}$–$10^{-12}$ m/kg1/2 for Advanced LIGO and $A^{\prime}=10^{-17}$–$10^{-13}$ m/kg1/2 for Einstein Telescope.

We have performed a scan of the parameter space using an identical method to that described in the previous section. Figure 10 shows 95% sensitivities in two-dimensional parameter spaces where we have marginalized over the third model parameter. Regions of parameter space above each of the curves should be detectable by Advanced LIGO (black) or the Einstein Telescope (gray). In the $A^{\prime}$-$\mu$ plane, the Advanced LIGO curve reaches a minimum near $\mu=80$ Hz, since this is where the stochastic search achieves its best sensitivity with this detector (see Figure 1). However, the Einstein Telescope is most sensitive below 10 Hz, which is not included in our parameter space for $\mu$. As a result, this contour decreases with decreasing $\mu$. For both detectors, the contours in the $A^{\prime}$-$\sigma$ parameter space decrease as $\sigma$ increases. This is because $\sigma$ acts similarly to a scaling parameter for $\Omega_{\text{GW}}(f)$, as previously discussed.

We also marginalize over the $\mu$ and $\sigma$ parameters to estimate an expected sensitivity on $A^{\prime}$ at 95% confidence. These results indicate that Advanced LIGO should be sensitive to cases with $A^{\prime}>1.4\times 10^{-14}{\rm\;m/kg^{1/2}}$ and Einstein Telescope should be sensitive to cases with $A^{\prime}>1.4\times 10^{-15}{\rm\;m/kg^{1/2}}$. In the bottom panel of Figure 10, we have used these expected sensitivities on $A^{\prime}$ to calculate contours in the $A$-$\lambda_{\text{CC}}$ parameter space. Taking $\lambda_{\text{CC}}$ to be 0.01 M${}_{\odot}^{-1}$ gives sensitivities of $A>$ 197 m for Advanced LIGO and $A>$ 19.5 m for the Einstein Telescope.

For comparison, we examine the s27f${}_{\text{heat}}$1.00 simulation of Ott et al. (2013) and the A1O05.5 simulation of Abdikamalov et al. (2014). For each simulation, we compute the strain spectrum $f|\tilde{h}(f)|$ and fit the functional form from Eq. 7 to the simulation (not including the frequencies where the low-frequency structure is prevalent). Next, we subtract this fit from the data and then fit the Gaussian functional form given in Eq. 12 to the low-frequency residual. The best fits for this functional form correspond to values of $A=0.18$ m for the s27f${}_{\text{heat}}$1.00 simulation and $A=0.85$ m for the A1O05.5 simulation. These values are approximately 3 and 2 orders of magnitude below the expected Advanced LIGO sensitivity to $A$, respectively, and 2 and 1 orders of magnitude below the expected Einstein Telescope sensitivity, implying that these experiments are unlikely to detect this signal.

III.3 Combined Case

We have also performed a parameter estimation study for the combined model, which is the sum of the low-frequency peak and the higher frequency model, as shown in Eq. 15. In this case, the model has six parameters, and a brute force exploration of the model parameter space is computationally infeasible. In order to study this model, we employ the MultiNest algorithm Feroz and Hobson (2008); Feroz et al. (2009a, b). The same parameter ranges given in the previous two sections are used in this combined study.

Figure 11 shows 95% sensitivity contours in the $b$-$a$ and $A^{\prime}$-$\xi$ planes for both Advanced LIGO and the Einstein Telescope, obtained by marginalizing over the other 4 parameters. In the $b$-$a$ plane, the region to the right of the curves is enclosed by the contours; thus, these detectors should be sensitive to parts of the parameter space to the left of their respective contours. This behavior is apparent because $a$ acts as an inverse scaling parameter in our model and $b$ acts similarly to a scaling parameter.

In the $A^{\prime}$-$\xi$ plane, the region below the curves is enclosed by the contours, so these detectors should be sensitive to areas of the parameter space above the contours. The shape of the contours is immediately apparent since both of the parameters illustrated here are scaling parameters.

The 2D posterior in the $\sigma$-$\mu$ plane spans most of the parameter space and does not yield significant insights into these parameters.

In Table 1, we compare the expected sensitivities on the $\xi$ and $A^{\prime}$ parameters for the combined model to those from the individual low- and high-frequency models. The sensitivities for the combined model tend to be better because the combined model includes both low- and high-frequency contributions, increasing the overall GW energy density and improving the detectability of such an SGWB.

IV Conclusions

In this paper, we have studied the stochastic gravitational wave background generated by stellar core collapse events occurring throughout the universe. Since the energy spectrum of gravitational waves emitted in a single core collapse event is not well understood, we have modeled this spectrum with an empirical functional form, using model parameter ranges determined from fits to a number of numerical simulations of core collapse events. In addition, we have noted that some of the simulations predict a low-frequency peak in the emitted spectrum, possibly due to acoustic waves generated by prompt convection occurring immediately after core bounce. The peak often appears in the vicinity of 100 Hz, and is therefore in or near the most sensitive frequency band for stochastic searches with terrestrial gravitational-wave detectors.

We have performed systematic scans of the model parameter space, including both the empirical broadband spectrum and the low-frequency peak, and compared the resulting stochastic background spectra to the expected sensitivities of the second- and third-generation detectors. While this background is unlikely to be detected by the second-generation, advanced detectors, Einstein Telescope may be able to detect the background in optimistic scenarios. However, in the majority of the parameter space we examined, even the Einstein Telescope would not be sensitive enough to detect this background.

It should be noted that the recent detections of gravitational waves from the binary black hole (BBH) mergers Abbott et al. (2016b, c) have led to new estimates of the stochastic background due to the BBH mergersAbbott et al. (2016d). The BBH background is estimated to be significantly stronger than the core collapse background for most of the parameter space examined here, and is potentially detectable by advanced detectors. The two backgrounds are predicted to have different spectral shapes (the BBH stochastic background spectrum is proportional to $f^{2/3}$), which offers a potential handle for distinguishing between them. However, distinguishing between the two types of stochastic background will require high signal-to-noise measurement of the stochastic background across the sensitive frequency band 1–100 Hz, as well as a better understanding of the core collapse process and of the gravitational wave spectrum emitted in this process, so that both statistical and systematic uncertainties are smaller than the amplitude of the core collapse background. Further development of numerical simulations of the core collapse is therefore critical for such studies.

Acknowledgments

The work K. A. O. was supported in part by DOE grant DE-SC0011842 at the University of Minnesota. The work of K. C., T. P., and V. M. was supported in part by NSF grant PHY-1204944 at the University of Minnesota. The work of E. V. has been carried out at the ILP LABEX (under reference ANR-10-LABX-63) supported by French state funds managed by the ANR within the Investissements d’Avenir programme under reference ANR-11-IDEX-0004-02. LIGO-P1600324.