A NEW GRAVITATIONAL-WAVE SIGNATURE FROM STANDING ACCRETION
SHOCK INSTABILITY IN SUPERNOVAE
Takami Kuroda
1
, Kei Kotake
2,3
, and Tomoya Takiwaki
3
1
Department of Physics, University of Basel, Klingelbergstrasse 82, 4056 Basel, Switzerland
2
Department of Applied Physics, Fukuoka University, 8-19-1, Jonan, Nanakuma, Fukuoka 814-0180, Japan
3
Division of Theoretical Astronomy, National Astronomical Observatory of Japan, 2-21-1, Osawa, Mitaka, Tokyo 181-8588, Japan
Received 2016 May 30; revised 2016 July 29; accepted 2016 August 24; published 2016 September 20
ABSTRACT
We present results from fully relativistic three-dimensional core-collapse supernova simulations of a non-rotating
M15
star using three different nuclear equations of state (EoSs). From our simulations covering up to ∼350 ms
after bounce, we show that the development of the standing accretion shock instability (SASI) differs significantly
depending on the stiffness of nuclear EoS. Generally, the SASI activity occurs more vigorously in models with
softer EoS. By evaluating the gravitational-wave (GW) emission, we find a new GW signature on top of the
previously identified one, in which the typical GW frequency increases with time due to an accumulating accretion
to the proto-neutron star (PNS). The newly observed quasi-periodic signal appears in the frequency range from
∼100 to 200 Hz and persists for ∼150 ms before neutrino-driven convection dominates over the SASI. By
analyzing the cycle frequency of the SASI sloshing and spiral modes as well as the mass accretion rate to the
emission region, we show that the SASI frequency is correlated with the GW frequency. This is because the SASI-
induced temporary perturbed mass accretion strikes the PNS surface, leading to the quasi-periodic GW emission.
Our results show that the GW signal, which could be a smoking-gun signature of the SASI, is within the detection
limits of LIGO, advanced Virgo, and KAGRA for Galactic events.
Key words: gravitational waves – hydrodynamics – supernovae: general
1. INTRODUCTION
Clarifying a correspondence between core-collapse super-
nova (CCSN) dynamics and the gravitational-wave (GW)
signals is a time-honored attempt since the 1980s
(Müller
1982). Very recently the observational horizon of
GW astronomy extends far enough to allow the first detection
coined by LIGO for the black hole merger event (Abbott
et al. 2016). Extensive research over the decades has
strengthened our confidence that CCSNe, next to compact
binary mergers, could also be one of the most promising
astrophysical sources of GWs (see Ott
2009; Kotake 2013 for
reviews).
Traditionally, most of the theoretical predictions have
focused on the GW signals from rotational core-collapse and
bounce (see, e.g., Dimmelmeier et al.
2002; Scheidegger et al.
2010; Ott et al. 2012; Kuroda et al. 2014; Yokozawa et al.
2015). In the postbounce phase, a variety of GW emission
processes have been proposed, including convection inside the
proto-neutron star (PNS)
and in the postshock region (Burrows
& Hayes
1996), the standing accretion shock instability (SASI;
Kotake et al.
2007, 2009; Marek & Janka 2009; Murphy et al.
2009) and nonaxisymmetric instabilities (Ott et al. 2005;
Scheidegger et al.
2010; Kuroda et al. 2014).
In the non-rotating core, Murphy et al. (
2009) first showed in
their two-dimensional (2D) models that the evolution of
convective activities in the PNS surface regions can be imprinted
in the GW spectrogram. The characteristic GW frequency is
considered as a result of the g-mode oscillation excited by the
downflows to the PNS (Marek et al.
2009) and by the
deceleration of convection plumes hitting the surface (Murphy
et al.
2009). These features have also been identified in more
recent 2D models with the best available neutrino transport
scheme (Yakunin et al. 2010;Mülleretal.2013;Yakunin
et al.
2015). Furthermore, Müller et al. (2013) showed in their
self-consistent 2D models that the SASI motions become
generally more violent for more massive progenitors, which
tends to make the GW amplitudes and frequencies higher.
Not to mention the explosion dynamics (e.g., Hanke et al.
2012; Couch 2013; Takiwaki et al. 2014; Janka et al. 2016), the
GW signatures are very sensitive to the spatial dimension
employed in the numerical modeling (e.g., Kotake et al.
2009;
Müller et al.
2012). Due to the high numerical cost, however,
only a few full three-dimensional (3D
) models have been
reported so far to study the postbounce GW features (without
any symmetry constraints and excision of the PNS; e.g.,
Scheidegger et al.
2010; Ott et al. 2012; Kuroda et al. 2014).
Using a prescribed boundary condition of the PNS contraction,
Hanke et al. (
2013) showed in their 3D models that a rapid
shrinking of the PNS fosters the development of the SASI.
General relativity (GR) should play a crucial role because the
SASI is favored by smaller shock radii due to the short SASI’s
growth rate (Foglizzo et al. 2006). To have a final word on
recent hot debates about the impacts of neutrino-driven
convection versus the SASI on the supernova mechanism
(e.g., Burrows
2013), full 3D-GR models are needed, which is
also the case for clarifying the GW emission processes.
In this Letter, we study the GW emission from a non-rotating
M15
star by performing 3D-GR hydrodynamic simulations
with an approximate neutrino transport. Using three modern
nuclear equations of states (EoSs), we investigate its impacts on
both the postbounce dynamics and the GW emission. Our
results reveal a new GW signature where the SASI activity is
imprinted. We discuss how the detectability of the signals, if
detected, could provide the live broadcast that shows how the
supernova shock is dancing in the core.
The Astrophysical Journal Letters, 829:L14 (6pp), 2016 September 20
doi:10.3847/2041-8205/829/1/L14
© 2016. The American Astronomical Society. All rights reserved.
1
2. NUMERICAL METHODS
In our full GR radiation-hydrodynamics simulations, we
solve the evolution equations of metric, hydrodynamics, and
neutrino radiation. Each of them is solved in an operator-
splitting manner, but the system evolves self-consistently as a
whole satisfying the Hamiltonian and momentum constraints
(Kuroda et al.
2012, 2014).
Regarding the metric evolution, we evolve the standard BSSN
variables
˜
g
ij
, f,
˜
A
ij
, K,and
˜
G
i
(Shibata & Nakamura 1995;
Baumgarte & Shapiro
1999). The gauge is specified by the “1
+log” lapse and by the Gamma-driver-shift condition.
In the radiation-hydrodynamic part, the total stress-energy
tensor
()
ab
T
total
is expressed as
()
() ()
¯
()
å
=+
ab ab
nnnn
n
ab
Î
TT T,1
total fluid
,,
eex
where
()
ab
T
fluid
and
()n
a
b
T
are the stress-energy tensor of fluid and
the neutrino radiation field, respectively. All radiation and
hydrodynamical variables are evolved in conservative ways.
We consider all three flavors of neutrinos (
¯
n
nn,,
eex
) with
n
x
representing heavy-lepton neutrinos (i.e.,
n
n
m
t
,
and their anti-
particles). To follow the 3D hydrodynamics up to
400
ms
postbounce, we shall omit the energy dependence of the
radiation in this work (see, however, Kuroda et al.
2016).
We use three EoSs based on the relativistic-mean-field
theory with different nuclear interaction treatments, which are
DD2 and TM1 of Hempel & Schaffner-Bielich (
2010) and
SFHx of Steiner et al. (
2013). For SFHx, DD2, and TM1
4
, the
maximum gravitational mass
M
max
and the radius of cold NS R
in the vertical part of the mass–radius relationship are
=M 2.13
max
, 2.42, and 2.21
M
and
~
R
1
2
, 13, and, 14.5
km, respectively (Fischer et al.
2014). SFHx is thus softest
followed in order by DD2 and TM1. Among these three, while
DD2 is consistent with nuclear experiments, such as for its
symmetry energy (Lattimer & Lim
2013), SFHx is the best-fit
model with the observational mass – radius relationship. All
EoSs are compatible with NS mass measurement ∼2.04
M
(Demorest et al. 2010). Our 3D-GR models are named DD2,
TM1, and SFHx, which simply reflects the EoS used.
We study a frequently used 15 M
e
star of Woosley &
Weaver (
1995). The 3D computational domain is a cubic box
with 15,000 km width, and nested boxes with eight refinement
levels are embedded. Each box contains 128
3
cells, and the
minimum grid size near the origin is
D
=x 458
m. In the
vicinity of the stalled shock front
~
R
100
km, our resolution
achieves
D
~x 1.9
km, i.e., the effective angular resolution
becomes
~
1
.
Extraction of GWs from our simulations is done by the
conventional quadrupole formula in which the transverse and
the trace-free gravitational field h
ij
is expressed by (Misner
et al.
1973)
()
() ()
()qf
qf qf
=
+
++´´
h
AeAe
D
,
,,
.2
ij
In Equation (2),
()qf
+´
A
,
represents the amplitude of
orthogonally polarized wave components with emission angle
(
)qf,
dependence (Scheidegger et al. 2010; Kuroda
et al.
2014),
+´
e
denotes unit polarization tensors, and D is
the source distance where we set D=10 kpc in this Letter.
3. RESULTS
We start by describing the hydrodynamics at bounce. The
central rest mass density
r
c
reaches
r
= 3.69,
c
3.75 and 4.50
×10
14
gcm
−3
for TM1, DD2, and SFHx, which is higher, as
expected, for the softer EOS (e.g., Fischer et al. 2014).
Figure 1. In each set of panels, we plot (top) the gravitational-wave amplitude of plus mode
+
A
[cm] and (bottom) the characteristic wave strain in the frequency-time
domain
˜
h
in a logarithmic scale that is overplotted by the expected peak frequency
F
peak
(black line denoted by “A”). “B” indicates the low-frequency component. The
component “A” is originated from the PNS g-mode oscillation (Marek & Janka
2009; Müller et al. 2013). The component “B” is considered to be associated with the
SASI activities (see Section
3). Left and right panels are for TM1 and SFHx, respectively. We note that SFHx (left) and TM1 (right) are the softer and stiffer EoS
models, respectively.
4
The symmetry energy S at nuclear saturation density is S=28.67, 31.67,
and 36.95 MeV, respectively (e.g., Fischer et al. 2014).
2
The Astrophysical Journal Letters, 829:L14 (6pp), 2016 September 20 Kuroda, Kotake, & Takiwaki
After bounce, the non-spherical matter motion develops and
starts GW emission. In Figure
1, we plot time evolution of the
angle-dependent GW amplitude (only plus mode
()qf
+
A
,
,
black line) in the top panels and the characteristic wave strain
in the frequency-time domain
˜
()qfhF,,
(see Equation (44) in
Kuroda et al.
2014) in the bottom ones. Here F denotes the GW
frequency. We extract GWs along the north pole
(
)(
)
qf=,0,0
. The postbounce hydrodynamics evolutions in
DD2 are rather similar to TM1 and we mainly focus on the
comparison between SFHx and TM1 in the following.
The GW amplitude (
+
A
, top panels) shows a consistent
behavior as reported in Müller et al. (
2013), Ott et al. (2013),
and Yakunin et al. (
2015). It shows an initial low frequency and
slightly larger amplitude until
~T 60
pb
ms, which is followed
by a quiescent phase with a higher frequency until
~T 150
pb
ms. Afterward, the amplitude and frequency become
larger with time.
From the spectrograms (bottom panels), we see a narrow-
band spectrum (labeled “A” in both models) that shows an
increasing trend in its peak frequency. Müller et al. (2013) and
Murphy et al. (2009) showed that this peak shift can be
explained by properties of PNS, such as its compactness and
surface temperature. By following Equation (17) in Müller
et al. (
2013), we overplot
F
peak
in the bottom panels (black
line). In both models,
F
peak
indeed tracks spectral peak quite
well, although there are some exceptions in the late phase of
SFHx (
T 200
pb
ms) when the other strong component
appears at
F100 200 Hz
(labeled “B”). The component
“A” is thus actually originated from the g-mode oscillation of
the PNS surface.
Before going into detail to explain the origin of the low-
frequency component “B,” we briefly focus on several key
differences in the hydrodynamic evolution between SHFx and
TM1. In Figure
2, SFHx experiences violent sloshing (top left)
and spiral motions of the SASI (top right) before neutrino-
driven convection dominates over the SASI (bottom left),
whereas the SASI activities are less developed in TM1. For
SFHx, the clear SASI motions are observed after the prompt
convection phase ceases at
~T 50
pb
ms.
In Figure
3, we plot time evolutions of maximum, average,
an minimum shock radii
R
shoc
k
(top, solid) and normalized
mode amplitudes
∣
∣∣∣∣
∣
ºAcc
lm lm 00
(see Burrows et al. 2012
for c
lm
) of spherical polar expansion of the shock surface
(
)
qf
R
,
shock
. For A
lm
, we plot models SFHx (middle) and TM1
(bottom) with focusing a period of
T120 300
pb
ms that
corresponds to the appearance of component “B.” We also plot
Figure 2. Snapshots of the entropy distribution (
k
B
baryon
−1
) for models SFHx and TM1 (top left,
=T 150
pb
ms of SFHx; top right,
=T 237
pb
ms of SFHx; bottom
left,
=T 358
pb
ms of SFHx; bottom right,
=T 358
pb
ms of TM1). The contours on the cross sections in the x=0 (back right), y=0 (back left), and z=0 (bottom)
planes are, respectively, projected on the sidewalls of the graphs. The 90° wedge on the near side is excised to see the internal structure. Note that to see the entropy
structure clearly in each dynamical phase, we change the maximum entropy in the color bar as
=s 16
max
, 20, and 22
k
B
baryon
−1
for
=T 150
pb
, 237, and 358 ms,
respectively.
3
The Astrophysical Journal Letters, 829:L14 (6pp), 2016 September 20 Kuroda, Kotake, & Takiwaki
the spherically averaged gain radius
R
gain
(dashed) in the top
panel.
The characteristic SASI motions seen in Figure 2 are
reflected in the evolution of
∣
∣A
lm
. For SFHx, the most dominant
mode during the first phase of the SASI (50 ms
T 150
pb
ms) is the sloshing mode, i.e.,
(
)( )=lm,1,0
,
which is in accord with the clear one-sided shock-heated region
(top left panel of Figure
2). Regarding the EoS dependence,
although we do not see any qualitative differences between the
stiffest EoS model TM1 and the softest EoS model SFHx, TM1
shows less SASI development, i.e., smaller values of
∣
∣A
lm
,
during the SASI development phase. DD2 also shows less
SASI development compared to SFHx. Such a quantitative
difference can be explained by the shock radius. In the top
panel of Figure
3, TM1 shows more extended shock radii until
~T 150
pb
ms. This is because, depending on the stiffness of
nuclear EoS, the bounce shock can be formed at larger radius
that can sometimes amount to
~ M0.1
difference in mass
coordinates (Suwa et al.
2013; Fischer et al. 2014). Conse-
quently, the prompt shock has to plunge into more material and
stalls at smaller radius in our softest EoS model SFHx. The
smaller shock radius is a favorable condition for the SASI
development due to the shorter advective–acoustic cycle
(Foglizzo
2002; Scheck et al. 2008). Initial SASI activities
reach their maxima when the shock expansion occurs due to
sudden drop of the mass accretion rate at
~T 150
pb
ms.
Afterward, the spiral mode becomes dominant as seen in
A
11
(see also the top right panel of Figure 2), which lasts another
~150 200
ms (SFHx/TM1).
In the final phase, the core experiences neutrino-driven
convection until the end of our calculation time
~T 350
pb
ms.
During this phase, matters in the gain region are exposed
intensively to neutrino radiations and form high entropy
(
)~s 20
k
B
smaller-scale convection plumes (middle and
bottom panels of Figure
2). Following Foglizzo et al. (2006),
we check the parameter χ. Although
c
3
is expected to be
satisfied for convection to develop, we find that χ stays ∼0.5
until
T 350
pb
ms in both models despite the appearance of
convection plumes. As already pointed out in Ott et al. (2013)
and Hanke et al. (
2013), this is because the initial perturbations
in the gain region are already not small when the neutrino
convection phase initiates. The gain radius (
R
gain
in Figure 3)
appears more inward in SFHx, which leads to higher entropic
convection plumes compared to those in TM1 (compare the
bottom two panels of Figure
2).
Now, we discuss how these hydrodynamical evolutions
affect the GW emission “B” in Figure
1. By spatially
decomposing the quadrupole moment of matters into several
spherical shells, we roughly localize this emission at
R10 20
km (Figure 4).
Before going into further discussion, we present a back-of-
the-envelope estimation of the GW amplitude as
∣∣
˙
()
~~~Dh MR T M R R M M222,3
2
dyn
2
22
2
///
where M, R, and
T
dy
n
represent the mass, size, and dynamical
timescale of the system, respectively, in geometrized unit.
Figure 3. Top: time evolution of maximum, average, and minimum shock radii
(solid) and spherically averaged gain radius (dashed) for models SFHx (red)
and TM1 (black). Two vertical dotted lines represent the period when the low-
frequency component “B” appears (Figure
1). Time evolution of normalized
mode amplitudes
∣
∣A
lm
for several representative modes (l, m) of SFHx (middle)
and TM1 (bottom). We show the period bounded by two vertical dotted lines in
the top panel.
Figure 4. Rough measurement of contributions from each spherical shell to (a)
the GW amplitude and (b1–4) their spectrogram
˜
h
in a logarithmic scale. We
show the contributions from four spherical shells with intervals of [0, 10], [10,
20], [20, 30], and [30, 100] km. Black contours overplotted on spectrograms
for
˜
h
represent the half-maximum of spectrograms for the mass accretion rate
measured at R=17 (b2),23(b3), and 48 (b4) km.
4
The Astrophysical Journal Letters, 829:L14 (6pp), 2016 September 20 Kuroda, Kotake, & Takiwaki
Here, we have used the following reasonable assumptions:
˙
()~TMM 4
dyn
or
()~~TRVRM,5
dyn
3
//
with
~
V
MR
being the velocity derived by the energy
conservation. From the last relation in Equation (
3), we expect
that significant time variation in the mass accumulation onto
the PNS can potentially lead to the GW emission. In Figure
4,
we superimpose the spectrogram of the mass accretion rate
˙
()MR
(the black contour at half-maximum) measured at
R=17, 23, and 48 km on top of the GW spectrogram. While
˙
()=MR 48 km
starts quasi-periodic oscillation at
F∼100–200 Hz around
~T 120
pb
ms, we find a time delay
of ∼60 ms for their appearance at a deeper region (R = 17 and
23 km). Since the density averaged mean radial velocity
between the lepton-driven (
R10 20
km) and the entropy-
driven (
R
40
km) convection layers is
~´510c
m
7
s
−1
, the
time delay is consistent with the advection timescale over the
stable layer (
R
2
040
km). Furthermore, coincidence of
time modulation in
˙
()MR
and the GW component “B” is
obvious from panel (b2).
Finally, to connect the SASI activities with the GW
component B, we plot spectrograms of the normalized mode
amplitude of the sloshing-SASI mode
∣
˜
∣
A
10
, the mass accretion
rate
∣
˙
˜
∣
M
measured at R=17 km, normalized quadrupole
deformation of the isodensity surface
˜
l
for l=2, and a rough
measurement of the GW energy spectrum in Figure
5.
˜
l
denotes
a Fourier component of normalized mode amplitude
l
defined by
() ()
å
º
=-
RR,6
l
mll
lm
,
,
14
2
0,0
14
/
where
R
lm,
14
is evaluated by the spherical polar expansion of the
isodensity surface R
14
extracted at
r
= 10
14
gcm
−1
as the same
way as for the shock surface. Although several other modes are
excited at the surface, only the leading contribution (l = 2
mode) to the GW emission is shown in the panel. As a
reference, the isodensity surface R
14
locates ∼13.5 km during
T150 300
pb
ms in SFHx. From the last relation in
Equation (3), we plot
∣∣
˙
~+hMlog log const.
10 10
2
in panels
(d) of Figure
5 with assuming
=MM0.5
, a mass contained in
R10 20
km, and
=
R
13.5
14
km stays nearly constant.
During
T140 180
pb
ms in SFHx, we see a strong
sloshing motion that has its peak frequency at
F100 200 Hz
(a1). With some time delay (∼50 ms) from
the appearance of it, the mass accretion rate
˙
M
starts showing a
quasi-periodic oscillation at the same frequency range
F100 200 Hz
(b1) and it excites oscillation on the
isodensity surface (c1). A combination of large
˙
M
and
2
causes us to expect GW emissions to appear in panel (d1), and
it can well explain Figure
1. During
T
2
00 300
pb
ms,
2
stays at ∼3×10
−4
in SFHx. A rough measurement of the GW
amplitude due to this deformation,
~
-
A
MR2
2
21
, deduces
~
A
2c
m
, which is consistent with the actual amplitude
(Figure
4).
4. SUMMARY AND DISCUSSION
We have presented relativistic 3D SN simulations with three
different nuclear EoSs. The overall pictures of SN dynamics are
qualitatively the same among all three models, although the
development of the SASI differs quantitatively. The softer the
EoS is, the more the SASI develops, since the prompt shock
stalls at smaller radii. The evolution shows the first prompt
convection phase, the sloshing-SASI phase, which shifts to the
spiral mode and finally to the neutrino-driven convection
phase.
Regarding the GWs, we have also confirmed previously
reported emissions originated from the PNS surface g-mode
oscillation (Murphy et al.
2009; Müller et al. 2013).
Additionally, in the softest EoS model SFHx, in which the
most vigorous SASI motion was observed, we have found
another low-frequency (
F100 200
Hz) quasi-periodic
emission. This emission was spatially localized at
Figure 5. Spectrograms of (a) Fourier decomposed normalized mode amplitude
∣
˜
∣A
10
of the shock surface for the sloshing-SASI mode; (b) the mass accretion rate
˙
˜
M
(with a dimension of
M
), through surface of a sphere with radius of R=20 km; (c) deformation of the isodensity surface
˜
l
for l=2 mode; and (d) a rough
measurement of the GW energy spectrum that is proportional to
˙
~
-
RMM
2
2
1
(see the text). Top and bottom rows are for SFHx and TM1, respectively.
5
The Astrophysical Journal Letters, 829:L14 (6pp), 2016 September 20 Kuroda, Kotake, & Takiwaki