MNRAS
507,
3672–3689 (2021)
https://doi.org/10.1093/mnras/stab2388
Advance Access publication 2021 August 19
Infrared dust echoes from neutron star mergers
Wenbin Lu
,
1
,
2
‹
Christopher F. McKee
3
and Kunal P. Mooley
4
,
5
1
Department of Astrophysical Sciences, Princeton University, Princeton, NJ 08544, USA
2
TAPIR, Walter Burke Institute for Theoretical Physics, Mail Code 350-17, Caltech, Pasadena, CA 91125, USA
3
Physics Department and Astronomy Department, University of California at Berkeley, Berkeley, CA 94720, USA
4
National Radio Astronomy Observatory, P.O. Box O, Socorro, NM 87801, USA
5
Cahill Center for Astronomy, MC 249-17, Caltech, Pasadena, CA 91125, USA
Accepted 2021 August 5. Received 2021 July 28; in original form 2021 April 17
ABSTRACT
A significant fraction of binary neutron star mergers occur in star-forming galaxies where the UV-optical and soft X-ray emission
from the relativistic jet may be absorbed by dust and re-emitted at longer wavelengths. We show that, for mergers occurring in
gas-rich environment (
n
H
0
.
5cm
−
3
at a few to tens of pc) and when the viewing angle is less than about 30
◦
, the emission
from heated dust should be detectable by
James Webb Space Telescope
(
JWST
), with a detection rate of
∼
1yr
−
1
. The spatial
separation between the dust emission and the merger site is a few to 10 milli-arcsecs (for a source distance of 150 Mpc),
which may be astrometrically resolved by
JWST
for sufficiently high signal-noise-ratio detections. Measuring the superluminal
apparent speed of the flux centroid directly gives the orbital inclination of the merger, which can be combined with gravitational
wave data to measure the Hubble constant. For a line of sight within the jet opening angle, the dust echoes are much brighter
and may contaminate the search for kilonova candidates from short gamma-ray bursts, such as the case of GRB 130603B.
Key words:
gravitational waves – dust, extinction – infrared: general – gamma-ray bursts – neutron star mergers.
1 INTRODUCTION
Observations of GW170817 provided unambiguous support that short gamma-ray bursts (GRBs) are associated with binary neutron star
(bNS) mergers (Abbott et al.
2017
; Goldstein et al.
2017
; Kasliwal et al.
2017
; Fong et al.
2019
; Mooley et al.
2018a
; Margutti et al.
2018
),
in agreement with theoretical expectations (Eichler et al.
1989
; Narayan, Paczynski & Piran
1992
; Narayan, Piran & Kumar
2001
;Lee&
Ramirez-Ruiz
2007
; Nakar
2007
;Berger
2014
; Kumar & Zhang
2015
). Most (
50 per cent) of the short GRBs observed so far are from
late-type galaxies whereas
∼
30 per cent are from early-type galaxies (Fong et al.
2013
), which indicates that the probability of short GRB
occurrence is not only related to the stellar mass but also strongly influenced by star formation activity (Berger
2014
). Another independent
evidence for short delay times between star formation and bNS merger is that
∼
30 per cent of known merging Galactic bNSs have merger time
less than 100 Myr (e.g. Tauris, Langer & Podsiadlowski
2015
; Farrow, Zhu & Thrane
2019
), and the cosmological bNS merger rate inferred
from the Galactic bNS population (e.g. Kim, Perera & McLaughlin
2015
) is comparable to that measured by LIGO-Virgo gravitational wave
(GW) observations (Abbott et al.
2020a
). Furthermore, modeling of the Galactic bNS population shows that there is a statistically significant
excess of systems with short delay times compared to the canonical power-law of
t
−
1
(Beniamini & Piran
2019
), and such a steep delay-time
distribution is also consistent with the redshift distribution of short GRBs (Wanderman & Piran
2015
). Therefore, we expect a significant
fraction of bNS mergers/short GRBs to occur in gas-rich environment, where the UV-optical and soft X-ray emission from the relativistic jet
will be reprocessed by gas and dust along the beaming cone into longer wavelengths in the form of an echo. Indeed, the broad-band afterglows
of many short GRBs show evidence of dust extinction
A
V
∼
0.3–1
.
5 mag (Fong et al.
2015
) along the jet axis, including e.g. GRB 130603B
(Berger, Fong & Chornock
2013
;Tanviretal.
2013
; de Ugarte Postigo et al.
2014
; Fong et al.
2014
).
On the other hand, it is well-known that long GRBs of core-collapse origin are located in intensively star-forming environment. A large
fraction (10–50 per cent) of them occur in dusty environment such that their optical afterglows are attenuated by
A
V
>
0
.
5 mag (Cenko et al.
2009
;Melandrietal.
2012
; Perley et al.
2013
). Thus, their afterglow emission are more susceptible to gas and dust reprocessing, which has
been studied in the literature (Esin & Blandford
2000
; Madau, Blandford & Rees
2000
; Waxman & Draine
2000
; Fruchter, Krolik & Rhoads
2001
; Draine & Hao
2002
; Perna & Lazzati
2002
; Heng, Lazzati & Perna
2007
;Evansetal.
2014
). However, for the practical cases where
the observer’s line of sight (LOS) is close to the jet axis, the GRB-associated supernova (Woosley & Bloom
2006
) typically outshines the dust
echo in the near-infrared band.
E-mail:
wenbinlu@astro.princeton.edu
C
2021 The Author(s)
Published by Oxford University Press on behalf of Royal Astronomical Society
Downloaded from https://academic.oup.com/mnras/article/507/3/3672/6354800 by California Institute of Technology user on 03 November 2021
Infrared dust echo
3673
Figure 1.
Schematic picture of the model. A relativistic jet, launched from a bNS merger, generates bright UV-optical and soft X-ray emission when it is
decelerated by the circum-stellar medium at a distance of the deceleration radius
r
dec
∼
0
.
1 pc. The jet emission is beamed in a cone of half opening angle of
θ
j
,
and the dust grains at distances of a few to tens of parsecs within this cone are heated by the radiation to high temperatures. All dust grains are deplete
d below
the sublimation radius
r
sub
∼
3 pc. The energy absorbed by dust at
r
>
r
sub
is re-radiated in the infrared according to the local dust temperature, which roughly
follows a power-law scaling with radius
T
∝
r
−
1/3
. For a LOS at an angle
θ
obs
from the jet axis, light-travel delay causes longer wavelength dust emission to
arrive at the observer at later time.
In this paper, we propose that dust echoes from a fraction of bNS mergers whose jet axis are generally misaligned with the observer’s
LOS should be detectable by
James Webb Space Telescope
(
JWST
). This provides a new observable electromagnetic (EM) counterpart
1
to bNS
mergers, in addition to the kilonova/macronova (Li & Paczy
́
nski
1998
; Kulkarni
2005
; Metzger et al.
2010
; Barnes & Kasen
2013
; Tanaka &
Hotokezaka
2013
;Yu,Zhang&Gao
2013
; Metzger
2017
), short GRB, and synchrotron afterglow from the jet/ejecta (Nakar & Piran
2011
;
Metzger & Berger
2012
; Piran, Nakar & Rosswog
2013
; Hotokezaka et al.
2018
), as well as non-thermal emission from the wind nebula of
a possible long-lived NS remnant (Zhang
2013
; Metzger & Piro
2014
). We further demonstrate that, for sufficiently bright dust echoes, the
separation between the echo flux centroid and the position of the merger may be astrometrically resolved by
JWST
, and it can be used to infer
the inclination angle of the bNS orbital axis and measure the Hubble constant (Hotokezaka et al.
2019
).
The broad-brush picture considered in this work is shown in Fig.
1
. Before presenting the model, we list a number of potential physical
sources of dust extinction in the vicinity of bNS mergers in late-type galaxies.
(1) Diffuse interstellar medium. For typical Galactic gas-to-dust ratio, the V-band extinction is related to the gas column density by
A
V
5
.
3
×
10
−
22
N
H
mag cm
2
0
.
5mag(
H
/
10 M
pc
−
2
), where
N
H
(or
H
) is the number (or mass) column density, and typical non-star-
burst spiral galaxies have
H
∼
10
M
pc
−
2
(Kennicutt & Evans
2012
). The vertical scale-height of atomic gas in the disc is a few hundred
parsecs, whereas the molecular gas is more vertically confined. Thus, for those short GRBs with significant V-band extinction
A
V
=
0
.
5mag
and LOS path-length
H
=
300 pc through the gas-rich region, the mean number density of the intervening gas is
̄
n
H
N
H
/
H
1cm
−
3
.
However, due to the inhomogeneity of the interstellar medium, the gas density in the vicinity of a bNS merger may be smaller than
̄
n
H
(see
later discussions). Generally, the extinction is stronger at shorter wavelengths, with the rough scaling of
A
λ
∝
λ
−
1
for an average Galactic LOS
(Cardelli, Clayton & Mathis
1989
), so UV photons
λ
∼
0.1
μ
m may be significantly attenuated within a distance
H
from the source (and
ionizing photons above 13.6 eV may be largely absorbed by neutral gas). UV extinction is dominated by very small grains of size
λ
/2
π
∼
0.02
μ
m. These smaller grains have lower sublimation temperature. The fact that they cool less efficiently by infrared (IR) emission makes them
easier to heat up and sublimate. These two effects lead to the depletion of smaller grains in a larger volume than for bigger grains, which will
later be modeled in detail.
(2) Cradle giant molecular clouds (GMCs). The average mass of molecular clouds hosting OB stars, potentially multiple generations of
them, is about 3
×
10
5
M
. These massive GMCs are in turn destroyed mainly by photoionization on a time-scale of between 10 and 30 Myr
(Williams & McKee
1997
; Matzner
2002
). It is likely that some bNS mergers occur with delay times shorter than the cloud destruction time,
and these mergers may still be embedded in their cradle GMC. However, the fraction of bNS mergers with delay time less than 30 Myr cannot
be estimated with great confidence, given limited observational clues on the existence of such rapid mergers (they may be identified by future
1
An earlier study on the Compton scattering echo of short GRBs found the signal to be too faint except for very nearby sources where the observer’s LOS is
only slightly away from the jet beaming cone (Beniamini et al.
2018
).
MNRAS
507,
3672–3689 (2021)
Downloaded from https://academic.oup.com/mnras/article/507/3/3672/6354800 by California Institute of Technology user on 03 November 2021
3674
W. Lu, C. F. McKee and K. P. Mooley
LISA-like GW missions; Lau et al.
2020
). Theoretically, there is a ‘fast channel’
2
where the envelope expansion of a low-mass (2–4 M
)He
star leads to Roche–Lobe overflow towards a NS companion, and this leads to further shrinkage of the binary orbit before an ultra-stripped
supernova and potentially very short GW inspiral time
t
GW
a few Myrs (Belczynski et al.
2006
;Taurisetal.
2015
). The other consideration
is the motion of the bNS system. Roughly 20 Myr after the binary formation (for a primary mass of 12 M
for instance), the natal kick on the
first-born NS,
v
k1
, will induce a binary centre-of-mass speed of
v
k1
/10, because the other member is still a massive main-sequence star that is
10 times more massive than the neutron star. We see that
v
k1
100 km s
−
1
is required for the binary to remain bound to the GMC, which has
escape speed of
∼
10 km s
−
1
(similar to that of the ultra-faint dwarf galaxy Rec II, which has been enriched by r-process elements likely by a
bNS merger; Ji et al.
2016
). Such a small kick is possible because the primary star has lost its hydrogen envelope due to Roche–Lobe overflow
in a close binary before exploding as a Type Ib supernova, which may generate a substantially weaker NS kick than that of isolated pulsars
(e.g. Tauris et al.
2017
). Then, the second-born NS causes another a centre-of-mass kick of
v
k2
/2, that could be as small as 10 km s
−
1
in the
case (ultra-stripped) electron capture supernova (Beniamini & Piran
2016
;Taurisetal.
2017
). Observationally, the short GRB 101219A (
z
=
0.718) occurred in a host galaxy with estimated star-formation rate of
∼
16 M
yr
−
1
and stellar population age of
∼
50 Myr (Fong et al.
2013
).
The averaged star-formation rate surface density of this galaxy can be roughly converted to a molecular column density of
H
2
∼
10
2
M
pc
−
2
(Kennicutt & Evans
2012
). We note that the majority of long GRBs from core-collapse of massive stars should be in this category.
(3) Unassociated GMCs. The demographics of GMCs in the Milky Way have been studied by wide field surveys of molecular emission lines.
Their mass and size distributions are found to be d
N
/d
M
∝
M
−
1.8
and d
N/
d
r
c
∝
r
−
3
.
2
c
(e.g. Williams & McKee
1997
; Heyer, Carpenter & Snell
2001
), meaning that most molecular mass [
∫
d
MM
(d
N/
d
M
)
∝
M
0
.
2
] is in the most massive GMCs with
M
max
∼
5
×
10
6
M
and that most
geometric cross-section for LOS-intersection [
∫
d
r
c
r
2
c
(d
N/
d
r
c
)
∝
r
−
0
.
2
c
] is provided by the least massive GMCs with mass
M
min
∼
10
3
M
and size
r
c
,
min
∼
3 pc. These small GMCs contains a fraction
f
M
∼
(
M
min
/M
max
)
0
.
2
∼
20 per cent of the total molecular mass. Suppose a
bNS merger occurs in a molecular region with column density
H
2
=
10
2
M
pc
−
2
and vertical scale height of
z
H
2
=
50 pc (Heyer & Dame
2015
), then the mean number density of GMCs can be estimated by
̄
n
GMC
=
f
M
H
2
/
(
M
min
2
z
H
2
)
∼
2
×
10
−
4
pc
−
3
. Then, the probability
that a narrow radiation beam from a GRB jet intersects a molecular cloud within a distance
r
∼
10 pc is given by
P
∼
̄
n
GMC
πr
2
c
,
min
r
6 per cent (
r/
10 pc)(
H
2
/
10
2
M
pc
−
2
). When such an intersection occurs, nearly all UV photons (including those above 13.6 eV) will be
absorbed/scattered by dust grains, because at high gas density
n
H
∼
10
3
cm
−
3
, dust attenuation occurs within a small thickness where all the
gas is ionized by a very small fraction of the ionizing photons.
The above considerations suggest that a significant fraction of bNS mergers may be surrounded by dense gas on lenthscales of a few to
tens of parsecs. This motivates us to model the interaction between the jet radiation beam and the surrounding gas/dust, as well as to explore
what can be learned from the detection of such dust echoes.
This paper is organized as follows. In Section 2, we provide an order-of-magnitude estimate of the luminosity and duration of the dust
echo for arbitrary viewing angles, under the assumption that the jet generates a short bright pulse of optical-UV emission within its beaming
angle. Then, we show the light curves and flux centroid motions for a number of physically motivated cases and provide an estimate of the
detection rate of dust echoes in Section 3. Detailed modeling of the dust sublimation and emission processes is presented in Section 4, which
is the basis of Section 3. There we focus on dust heating by UV-optical emission; the effects of X-ray photons are discussed in Section 5. We
summarize in Section 6. The UV-optical emission from the reverse shock formed when the jet runs into the circum-stellar medium is calculated
Appendix A. We use the convention
Q
=
10
x
Q
x
in CGS units.
2 ORDER-OF-MAGNITUDE ESTIMATE
In this section, we estimate the luminosity and duration of the dust echo in a number of simple scenarios.
The first case is that the observer’s viewing angle
θ
obs
wrt. the jet axis is much smaller than the opening angle of the jet
θ
j
–theLOS
is effectively coincident with the jet axis. As a result of light travel delay, a spherical shell of dust at radius
r
heated from an isotropic short
burst of radiation produces a square-function echo light curve of duration 2
r
/
c
(this is called the ‘response function’). Here, the jet emission
only illuminates a small patch of solid angle
πθ
2
j
on the sphere, so the duration of the echo is
rθ
2
j
/
(2
c
) but the luminosity is the same as in the
case of a full sphere. Then we consider a large number of spherical shells, each of which has absorption optical depth d
τ
abs
(the fraction of
the incident UV-optical energy
E
UV
absorbed by this shell), and the total bolometric echo luminosity contributed by dust grains in all shells is
given by
L
on
-
axis
d
∼
∫
r
sub
d
τ
abs
E
UV
e
−
τ
ext
2
r/c
,
(1)
where
τ
ext
2
τ
abs
is the extinction (absorption plus scattering) optical depth below radius
r
,and
r
sub
is the sublimation radius below which
all dust grains are heated to above the sublimation temperature
T
sub
∼
2200 K (see Section 4) and have evaporated. Here for simplicity,
we ignore the attenuation due to gas ionization, because the isotropic equivalent energy it takes to ionize all the gas up to radius
r
is
∼
10
47
erg (
n
H
/
1cm
−
3
)(
r/
3pc)
3
E
UV
∼
10
50
erg. As shown in Appendix A, when the jet is decelerated by the gas in the circum-stellar
2
Another channel is that if the second-born NS receives a large natal kick comparable to the pre-explosion binary orbital speed (a few hundred km s
−
1
), there is
a finite chance that the resulting bNS system may be in a highly eccentric orbit. In the high-eccentricity limit, the GW merger time distribution as a res
ult of such
a strong kick is given by d
P/
dlog
t
GW
∝
t
2
/
7
GW
(Lu, Beniamini & Bonnerot
2021
), and it is possible to achieve extremely short GW inspiral time
t
GW
Myr.
MNRAS
507,
3672–3689 (2021)
Downloaded from https://academic.oup.com/mnras/article/507/3/3672/6354800 by California Institute of Technology user on 03 November 2021
Infrared dust echo
3675
medium (CSM) at a distance
∼
0
.
1 pc (the deceleration radius), the electrons accelerated by the reverse shock generate bright synchrotron
emission in the UV-optical frequency range with isotropic equivalent energy
E
UV
∼
10
50
erg, for an isotropic equivalent jet kinetic energy
E
j
∼
a
few
×
10
52
erg (as inferred from the luminosity function of short GRBs; e.g. Beniamini et al.
2019
), and CSM density
0
.
1cm
−
3
. For isotropic
equivalent UV luminosity
L
UV
=
E
UV
/
t
UV
and typical grain size
a
=
0
.
1
μ
m, the sublimation radius is estimated to be
r
sub
3pc
L
1
/
2
UV
,
48
a
−
1
/
2
0
.
1
μ
m
(see Section 4). At much larger distances
r
r
sub
, the grain temperature drops as
T
∝
r
−
1/3
(see Section 4) and the peak frequency of the dust
emission goes to longer and longer wavelengths.
In the optically thin limit
τ
ext
1, even though each radial shell with equal intervals in log(
r
) gives the same amount of bolometric echo
luminosity, the contribution to the near-IR (
λ
obs
∼
afew
μ
m) flux is dominated by the region near the sublimation radius. Thus, the specific
luminosity of the dust echo near the peak wavelength
λ
obs
∼
hc/
3
kT
sub
∼
2
μ
m(
T
sub
/
2200 K)
−
1
is given by
L
on
-
axis
d
,ν
∼
h
3
kT
sub
min
[
τ
abs
(
r
sub
)
,
1
2
]
E
UV
2
r
sub
/c
∼
(6
×
10
26
erg s
−
1
Hz
−
1
)min
[
2
τ
abs
(
r
sub
)
,
1
]
L
1
/
2
UV
,
48
t
UV
,
2
a
1
/
2
0
.
1
μ
m
,
(2)
where we have used fiducial UV-optical luminosity
L
UV
=
10
48
erg s
−
1
and duration
t
UV
=
100 s (cf. Appendix A). Note that equation (2)
applies to both optically thin (
τ
abs
1) and optically thick (
τ
abs
1) dust columns. The duration of the dust echo is given by
t
on
-
axis
d
,
obs
∼
θ
2
j
2
r
sub
c
∼
(9 d) (
θ
j
/
4
o
)
2
L
1
/
2
UV
,
48
a
−
1
/
2
0
.
1
μ
m
.
(3)
The above estimate shows that dust echo may be responsible for the famous kilonova candidate from GRB 130603B (Berger et al.
2013
;
Tanvir et al.
2013
), which was based on one H-band (rest-frame
λ
obs
=
1.2
μ
m) detection with
L
ν
8
×
10
26
erg s
−
1
Hz
−
1
at
t
obs
=
6
.
9din
the host-galaxy rest frame (redshift
z
=
0.356). This will be further discussed in Section 3.1.
The second case is that the observer’s viewing angle is much larger than the jet opening angle,
θ
obs
θ
j
. In this case, the dust echo arrives
at the observer after a delay
t
off
-
axis
d
,
obs
∼
r
sub
c
(1
−
cos
θ
obs
)
1
.
3yr(
θ
obs
/
30
o
)
2
L
1
/
2
UV
,
48
a
−
1
/
2
0
.
1
μ
m
,
(4)
where the second expression is valid for
θ
obs
1 rad (assumed to be true hereafter). The peak flux in the off-axis case is much smaller than
the on-axis case by a factor that is of the order (
θ
j
/
θ
obs
)
2
, for the following two reasons. First, only a fraction of the entire azimuth around
the LOS is occupied by the jet-illuminated patch,
θ
j
/
[2
π
sin(
θ
obs
)]
0
.
02 (
θ
j
/
4
o
)(
θ
obs
/
30
o
)
−
1
. Secondly, since the echo from an infinitesimal
radial shell at radius
r
sub
lasts for a duration of
t
d, obs
(
r
sub
/
c
)
θ
j
sin
θ
obs
, the total flux only comes from the contribution from a smaller
radial thickness of
r/r
sub
θ
j
sin
θ
obs
/
(1
−
cos
θ
obs
)
1
/
4(
θ
j
/
4
o
)(
θ
obs
/
30
o
)
−
1
. Combining this two factors, we obtain the peak flux for the
off-axis case
L
off
-
axis
d
,ν
∼
θ
j
2
π
sin
θ
obs
θ
j
sin
θ
obs
1
−
cos
θ
obs
L
on
-
axis
d
,ν
h
3
kT
sub
min
[
τ
(
r
sub
)
,
1
2
]
E
UV
θ
2
j
4
π
(1
−
cos
θ
obs
)
r
sub
/c
∼
(3
×
10
24
erg s
−
1
Hz
−
1
)(
θ
j
/
4
o
)
2
(
θ
obs
/
30
o
)
−
2
min
[
2
τ
abs
(
r
sub
)
,
1
]
L
1
/
2
UV
,
48
t
UV
,
2
a
1
/
2
0
.
1
μ
m
,
(5)
Another way of understanding the
θ
−
2
obs
scaling of the peak luminosity is that the total energy from dust emission is independent of the viewing
angle and that the flux is spread out over a longer time for larger viewing angles with the dependence
L
off
-
axis
d
,ν
∝
(1
−
cos
θ
obs
)
−
1
≈
θ
−
2
obs
/
2.
After reaching a peak at
t
off
-
axis
d
,
obs
(equation 4), the flux will stay near the peak for a duration comparable to the peak time, and then the flux will
decline gradually, as the contribution from colder dust at larger radii
r
r
sub
becomes dominant. Generally, at longer wavelengths, the flux
reaches the peak at later times and the post-peak decline is slower.
So far we have thrown away the dust-scattered light, which may be observable in the optical band by ground-based observatories (e.g. future
Extremely Large Telescopes
). Since the scattering cross-sections are similar to absorption, we expect the specific luminosity or photon number
flux of the scattered light in the optical band (at wavelength
λ
opt
) to be smaller than that from dust emission by a factor of (
λ
opt
/
λ
IR
)(
E
opt
/
E
UV
)
∼
0.1, where we have taken
λ
opt
∼
0
.
6
μ
m (V-band),
λ
IR
∼
2
μ
m (peak wavelength for dust near the sublimation temperature), and a fraction
(
E
opt
/
E
UV
)
∼
1/3 of the jet emission in the 0.02–1
μ
m range to be in the optical band near
λ
opt
. At small viewing angles
θ
obs
30
◦
,this
flux-density reduction is compensated (by a factor of
∼
2) by the fact that forward-scattering is slightly favoured (Draine
2011
). Additionally,
the kilonova emission (especially in the case of a long-live NS remnant; Yu et al.
2013
; Metzger & Piro
2014
) may also contribute to the
flux of the scattering echo, whose spectrum may have line features from atomic transitions. However, there may be additional dust extinction
along the LOS which further reduces the flux of the scattered light. In this work, we focus on the infrared emission from the dust heated by
UV-optical emission from the jet and do not consider scattered light.
Finally, the flux centroid of the dust emission is located at a projected distance of
r
sub
sin
θ
obs
from the centre of explosion. For a source
(angular-diameter) distance of
D
, this gives an angular separation of
̃
θ
proj
∼
r
sub
sin
θ
obs
D
3mas(
θ
obs
/
30
o
)
L
1
/
2
UV
,
48
a
−
1
/
2
0
.
1
μ
m
(
D/
100 Mpc)
−
1
.
(6)
In our numerical model, we find the projected angular separation to be slightly larger than the estimate above because of the flux contribution
from radii slightly larger than
r
sub
. The spatial extent of the dust echo is given by the transverse size of the illuminated patch
θ
j
r
sub
,which
gives an angular size of less than 1 mas. This means that the dust echo will be almost like a point source. Since the physical position of the
MNRAS
507,
3672–3689 (2021)
Downloaded from https://academic.oup.com/mnras/article/507/3/3672/6354800 by California Institute of Technology user on 03 November 2021
3676
W. Lu, C. F. McKee and K. P. Mooley
dust emission region moves at speed of light, the observed flux centroid of the dust echo will have apparent proper motion (in units of speed
of light)
β
app
=
sin
θ
obs
1
−
cos
θ
obs
3
.
7(
θ
obs
/
30
o
)
−
1
,
(7)
which should in principle be measurable given two epochs of astrometric dust echo data, or one epoch of dust echo data plus an earlier epoch
of kilonova emission (as the reference point).
3 DUST ECHO LIGHT CURVES, FLUX CENTROID POSITIONS, AND DETECTION RATE
We defer the detailed description of our numerical model to Section 4. In this section, we show the dust echo light curves and flux centroid
positions for a number of representative cases (the data can be downloaded from this URL).
3
Then, we provide an estimate of the detection
rate of dust echoes by
JWST
.
We consider a physical situation where the UV-optical emission from the jet is a square pulse of isotropic equivalent luminosity
L
UV
=
3
×
10
47
erg s
−
1
and duration
t
UV
=
300 s. The spectrum is taken to be a power law from synchrotron emission with
L
ν
∝
ν
(1
−
p
)/2
and
p
=
2.2, and the normalization is
L
ν
max
=
(3
−
p
)
L
UV
/
(2
ν
max
) at maximum frequency
hν
max
=
50 eV. As shown in Appendix A, these
parameters are expected from a bNS merger-like GW170817 with isotropic equivalent jet energy of
E
j
∼
afew
×
10
52
erg (e.g. Gill & Granot
2018
;Lazzatietal.
2018
; Mooley et al.
2018a
; Beniamini et al.
2019
; Hajela et al.
2019
) but embedded in a denser CSM of density
10
−
1
cm
−
3
in a late-type galaxy [see O’Connor, Beniamini & Kouveliotou (
2020
) for constraints on CSM densities for a sample of short GRBs]. Note
that the jet optical emission depends on the CSM density at a distance of
∼
0
.
1 pc, which could be unrelated to the density of the gas at a few
to tens of parsecs that is responsible for the dust echo. Photons above 50 eV may contribute comparable amount of heating (despite reduced
absorption cross-sections, see Fig.
B1
), but we do not consider soft X-ray photons in this work because rapid cooling of the radiating electrons
accelerated by the reverse shock may suppress the emission above 50 eV. The effects of X-ray photons will be discussed in Section 5. We
assume the half opening angle of the jet emission to be
θ
j
=
4
◦
, which is consistent with the jet in GW170817 (Mooley et al.
2018a
,
b
).
We calculate the sublimation radius for dust grains of different sizes and the extinction in an iterative way, and then the volumetric
emissivity is calculated from the dust temperatures for all grain sizes, and finally, we account for light-travel delays and obtain the light curve
of dust echo from a given observer’s viewing angle
θ
obs
and wavelength
λ
obs
. Below, we show the results for the off-axis (
θ
obs
θ
j
)casein
Section 3.1 and then the on-axis (
θ
obs
<θ
j
) case in Section 3.2.
3.1 Off-axis case
The results of the off-axis case are shown in Figs
2
,
3
,and
4
for three different viewing angles of
θ
obs
=
20
◦
,10
◦
, and 30
◦
, respectively, and
in each figure, we consider three different hydrogen number densities
n
H
=
0.3, 1, 3 cm
−
3
for the gas-dust mixture and different observing
wavelengths
λ
obs
=
2
,
4
μ
m. We find that the dust echo brightens on a time-scale of months to a few years after the merger, strongly depending
on the viewing angle. For a source distance of 150 Mpc, gas density
n
H
0
.
5cm
−
3
and viewing angle
30
◦
, the peak flux is photometrically
detectable at 3
σ
level by
JWST
with 20 ks exposure.
We also show that the late-time jet afterglow emission is subdominant compared to the dust echo, regardless of the density of the CSM
that the jet runs into. This is because, at a higher CSM density, the jet decelerates more rapidly, the afterglow flux peaks earlier and at a
larger value, but then the flux drops rapidly as
F
ν
∝
t
−
p
obs
after the peak time, where
p
is the power-law index of the electron Lorentz factor
distribution. Let us consider a jet with a fixed total energy interacting with CSM of different densities. For an off-axis LOS
θ
obs
θ
j
,the
time of the afterglow flux peak has the analytical scaling
t
obs
,
p
∝
n
−
1
/
3
θ
2
obs
, the peak flux scales as
F
ν,
p
∝
n
(
p
+
1)
/
4
θ
−
2
p
obs
(Nakar, Piran & Granot
2002
; Gottlieb, Nakar & Piran
2019
). The afterglow flux after the peak is given by
F
ν
F
ν,
p
(
t
obs
/t
obs
,
p
)
−
p
∝
n
(3
−
p
)
/
12
t
−
p
obs
, which depends on
the CSM density
n
extremely weakly (for
p
=
2.2, the dependence is
n
1/15
) and is independent of the viewing angle
θ
obs
. Therefore, for a jet
with similar energy and structure as GW170817 interacting with higher-density CSM, we expect the afterglow light curve after the flux peak to
be close to the extrapolation of the post-peak segment of the GW170817 afterglow. Such an extrapolation is shown as thin dash-double-dotted
lines in Figs
2
,
3
,and
4
, where the flux normalization at two different wavelengths
λ
obs
=
2
μ
m (red) and 4
μ
m (blue) are obtained from the
precisely measured power-law spectrum of GW170817 afterglow
F
ν
∝
ν
−
0.58
(e.g. Margutti et al.
2018
;Makhathinietal.
2020
; Troja et al.
2020
). This means that, for any CSM density, the late-time jet afterglow emission should be near or to the left of the extrapolated lines.
On the right
-
hand panels of the figures, we show the time evolution of the position of the flux centroid projected on the sky. For a source
distance of 150 Mpc, the angular separation between the echo flux centroid and the centre of explosion is a few to 10 milli-arcseconds, which
may be resolvable by
JWST
when the signal-noise-ratio (SNR) is above 10. The time-dependent intensity maps of the dust emission for a
representative case are shown in Fig.
5
. The astrometric error of
JWST
is estimated to be 1
/
SNR times the half width half-maximum (HWHM)
of the point spread function and the HWHM near 2
μ
m is about 30 mas (which is also the pixel scale of NIRCam), so the 1
σ
astrometric
precision is given by 3 mas (SNR
/
10)
−
1
(Mooley et al., in preparation). The echo has superluminal apparent motion at a speed given by
equation (7), which is independent of the dust grain model and the details of the jet optical emission. Thus, measuring this apparent speed
3
https://github.com/wenbinlu/dustecho.git
MNRAS
507,
3672–3689 (2021)
Downloaded from https://academic.oup.com/mnras/article/507/3/3672/6354800 by California Institute of Technology user on 03 November 2021
Infrared dust echo
3677
Figure 2.
Left-hand panel
: Light curves of the dust echo emission at different wavelengths and for three different scenarios of hydrogen number densities
n
H
=
0
.
3cm
−
3
(dash–dotted lines), 1 cm
−
3
(solid lines), and 3 cm
−
3
(dashed lines). The grey dotted (dashed) line denotes 10
σ
(3
σ
) sensitivity of
JWST
for 20
ks exposure time at a luminosity distance of 150 Mpc. Different line colours corresponding to wavelengths of
λ
obs
=
2
μ
m (red) and 4
μ
m (blue). For all cases,
we have fixed the observer’s viewing angle
θ
obs
=
20
◦
, jet opening angle
θ
j
=
4
◦
, and dust-to-gas ratio parameter
n
0
/
n
H
=
1.45
×
10
−
15
(for a typical LOS
in the Milky Way). The emission from the jet is taken to be a square pulse with isotropic equivalent UV-optical luminosity
L
UV
=
3
×
10
47
erg s
−
1
, spectrum
L
ν
∝
ν
(1
−
p
)/2
(
p
=
2.2), and duration
t
UV
=
300 s. The red diamonds show the afterglow light curve of GW170817, which are obtained by scaling the R-band
flux measurements to wavelength of 2
μ
m according to the known broad-band power-law spectrum (Fong et al.
2019
). The red and blue dashed-double-dotted
lines denote
F
ν
∝
t
−
p
obs
extrapolations of the GW170817 afterglow at
λ
obs
=
2and4
μ
m, respectively. For any CSM density, the late-time jet afterglow emission
from a GW170817-like bNS merger should be near or to the left of these two lines.
Right-hand panel
: Time evolution of the angular separation between the dust
echo flux centroid and the centre of explosion, for an assumed angular diameter distance of 150 Mpc. The results in all cases agree with the analytical re
sult of
apparent speed
β
app
=
sin
θ
obs
/(1
−
cos
θ
obs
). The astrometric precision of
JWST
is estimated to be about 3 mas for SNR
=
10, which is shown by a horizontal
grey dotted line.
Figure 3.
Same as Fig.
2
but for a viewing angle of
θ
obs
=
10
◦
and the LOS is not far from the jet beaming cone of opening angle
θ
j
=
4
◦
. The flux reaches
the peak at an earlier time than the
θ
obs
=
20
◦
case and the peak flux is higher.
provides a model-independent way of obtaining the LOS inclination angle that can be used to measure the Hubble constant (Hotokezaka et al.
2019
).
It should also be noted that, since we are mainly interested in the cases with
θ
obs
30
◦
, the dust echo from the counter jet (which
propagates along the jet axis away from the observer) is expected to be temporally separate from that from the forward jet (propagating towards
the observer) by 2
r
sub
/c
∼
20 yr or longer. Furthermore, since the dust echo from the counter jet is spread out over a much longer duration of
∼
r
sub
/
c
, the flux is much smaller than the forward jet echo by at least an order of magnitude (e.g. Heng et al.
2007
). Therefore, astrometric
measurement of the flux centroid motion within 20 years after the merger is practically unaffected by the counter jet.
In the following, we provide a rough estimate of the rate of GW-selected bNS mergers that have detectable dust echoes. The two main
factors affecting the detectability of dust echoes are the gas density near the jet axis
n
H
at distances of a few to tens of parsecs and the viewing
angle
θ
obs
.
The interstellar medium is known to be highly inhomogeneous, with
∼
half the volume occupied by very tenuous hot ionized medium
(
T
HIM
∼
10
6
K; McKee & Ostriker
1977
) and the other
∼
half is the warm neutral medium (WNM;
T
WNM
8000 K, Wolfire et al.
2003
). It is
MNRAS
507,
3672–3689 (2021)
Downloaded from https://academic.oup.com/mnras/article/507/3/3672/6354800 by California Institute of Technology user on 03 November 2021
3678
W. Lu, C. F. McKee and K. P. Mooley
Figure 4.
Same as Fig.
2
but for a viewing angle of
θ
obs
=
30
◦
. The flux reaches the peak at a later time than the
θ
obs
=
20
◦
case and the peak flux is lower.
Figure 5.
Multi-epoch intensity maps of the dust echo, showing the flux centroid (marked by an orange cross in each panel) moving across the sky at a
superluminal speed. For this case, we consider hydrogen density
n
H
=
1cm
−
3
(unimportant, as long as the dust column is optically thin), viewing angle
θ
obs
=
20
◦
, observing wavelength
λ
obs
=
2
μ
m, and source angular diameter distance of 150 Mpc. The position of the centre of explosion (marked by a red circle) is at
(
̃
x
=
0
,
̃
y
=
0), which can be determined to a precision much better than 1 mas, provided that the bright kilonova emission is detected by
JWST
within the first
month after the merger. The astrometric precision of
JWST
is estimated to be about 3 mas (SNR
/
10)
−
1
, so the proper motion of the dust echo may be resolved
at
t
obs
2 yr after the merger.
likely that a fraction of bNS mergers, those with very short delay times
30 Myr, may be preferentially located in denser regions (the cold
neutral medium), but here we conservatively assume that bNS mergers are uniformly distributed within the gaseous disc. For the
∼
half of them
embedded in WNM, the density of the surrounding gas can be estimated following the model by Ostriker, McKee & Leroy (
2010
)basedon
thermal equilibrium as well as dynamical equilibrium in the direction perpendicular to the galactic disc. If the gravity in the disc is dominated
by stars and dark matter, the thermal pressure of WNM is given by (Ostriker et al.
2010
)
P
th
H
(
2
π
Gζ
d
c
2
w
f
w
ρ
sd
α
)
1
/
2
,
(8)
where
H
is the column density of diffuse H
I
gas,
ζ
d
0.33 is a constant that depends weakly on the vertical distribution of the gas,
c
w
8kms
−
1
is the thermal sound speed (since the
T
WNM
is set by thermal balance and is insensitive to local galactic conditions),
f
w
∼
0.5 is
the fraction of H
I
in WNM phase,
ρ
sd
is the density of stars and dark matter in the disc mid-plane, and
α
∼
5 is the ratio of total pressure (turbulent
plus magnetic) to thermal pressure. Adopting the relation between V-band extinction and gas column density
A
V
0
.
5mag(
H
/
10 M
pc
−
2
)
MNRAS
507,
3672–3689 (2021)
Downloaded from https://academic.oup.com/mnras/article/507/3/3672/6354800 by California Institute of Technology user on 03 November 2021
Infrared dust echo
3679
and
ρ
sd
∼
0
.
1M
pc
−
3
(as appropriate for the Galactic inner disc),
4
we obtain the typical density of the WNM
n
WNM
0
.
8cm
−
3
A
V
mag
(
ρ
sd
0
.
1M
pc
−
3
)
1
/
2
.
(9)
We also note that, if the self-gravity of the gas in the disc is significant, then the density of the WNM will be larger.
Therefore, we see that bNS mergers with extinction
5
A
V
0
.
5 mag, LOS viewing angle
θ
obs
30
◦
, and distance
D
150 Mpc should
yield dust echoes detectable by
JWST
. These requirements are described by: (1)
P
1
∼
20–50 per cent of bNS mergers occur in star-forming
galaxies with significant dust extinction (
A
V
0
.
5 mag) along the jet axis (Fong et al.
2013
), (2)
P
2
∼
0.5 of the gaseous disc volume is
occupied by WNM, (3)
P
3
40 per cent of GW-selected mergers have viewing angle less than 30
◦
(Schutz
2011
). Thus, we infer the detection
rate of dust echoes within 150 Mpc to be in the range 0.5–1
.
4yr
−
1
(
R
bNS
/
10
3
Gpc
−
3
yr
−
1
), where
R
bNS
is the (highly uncertain) bNS merger
rate.
6
There is also an interplay between viewing angle and source distance. Since the peak flux of dust echo scales as
F
peak
∝
θ
−
2
obs
∝
D
−
2
,
we see that the detection horizon for dust echoes scales as
D
max
∝
θ
−
1
obs
. For small viewing angles
θ
obs
30
◦
, the detection rate scales as
D
3
max
θ
1
.
7
obs
∝
D
1
.
3
max
, where the factor of
θ
1
.
7
obs
comes from GW selection (Schutz
2011
). Thus, events at larger distances up to the GW detection
horizon of
∼
400 Mpc for the upcoming O4 observing run (Chen et al.
2021
) increases the detection rate of dust echoes by a factor of about
3.6 – the total detection rate of dust echoes by
JWST
is in the range 1.7–5 yr
−
1
(
R
bNS
/
10
3
Gpc
−
3
yr
−
1
). Furthermore, the kilonova emission
(especially the blue component), as well as the jet afterglow, from events with smaller inclination angles are generally brighter (Gottlieb
et al.
2019
; Darbha & Kasen
2020
; Kawaguchi, Shibata & Tanaka
2020
;Zhuetal.
2020
). Thus, the inclinations of the events with precise
localization by their EM counterparts are biased towards lower angles, and hence the fraction of them with detectable dust echoes is higher
than the purely GW-selected sample. However, we note that, for source distances larger than 150 Mpc, it may be very challenging to measure
the angular separation between the merger position and the flux centroid of the dust echo.
We conclude that it should be possible for
JWST
to detect the dust echoes from bNS mergers and that the detection rate is of the order once
per year. We suggest the following strategy: (1) after the GW alert of each bNS merger, conduct rapid ground-based wide-field optical and/or
infrared surveys (e.g. Andreoni et al.
2020
) to identify the host galaxy and focus on the cases of star-forming galaxies (non-star-forming-host
cases may be interesting for other purposes); (2) use
JWST
to take a deep exposure of the kilonova emission within the first month of the
merger (Kasliwal et al.
2019
), with an exposure time such that photometric SNR is high enough (
10) to accurately establish the astrometric
position of the merger, which is registered in the GAIA frame; (3) when the jet afterglow has faded away, take one more deep (
∼
20 ks)
JWST
exposure to search for the dust echo; (4) in the case of a successful detection of the dust echo, take additional exposures to increase the SNR
to greater than 10 and then look for astrometric shift between the echo and the position of the kilonova.
3.2 On-axis case
When the observer’s LOS is close to the jet axis (in the case where a short GRB is observed), the dust echoes may be sufficiently bright in the
near-IR to be confused with kilonovae. Here we show that the H-band excess from GRB 130603B at host-galaxy rest-frame time
t
obs
=
6
.
9d
(Berger et al.
2013
;Tanviretal.
2013
) may indeed be due to emission from heated dust, which is an alternative explanation to the kilonova
model. The specific luminosity of the excessive emission is 8
×
10
26
erg s
−
1
Hz
−
1
at rest-frame wavelength
λ
obs
1.2
μ
m, which is slightly
shorter than the peak of the emission spectrum for dust temperature
T
sub
2200 K. From equation (2), we see that the dust echo may be
consistent with the excess, as long as (1) the isotropic equivalent UV energy from the jet is
E
UV
∼
afew
×
10
50
erg and (2) an optically thick
dust layer exists at a distance of
r
sub
∼
3 pc. The first requirement can be satisfied if the jet energy is
E
j
∼
3
×
10
52
erg and the CSM density at
a distance of 0
.
1 pc is higher than 0
.
3cm
−
3
(cf. Appendix A). An optically thick (
A
V
∼
1 mag) dust layer of thickness
∼
3 pc means that the
hydrogen number density is
n
H
∼
300 cm
−
3
. This is possible if the source is embedded in or near a molecular cloud. In fact, the steep spectrum
of the optical afterglow from this event indeed showed LOS dust extinction of
A
V
1 mag in the star-forming host galaxy (de Ugarte Postigo
et al.
2014
; Fong et al.
2014
).
As a concrete example, we consider the UV-optical emission from the jet to be a square pulse with isotropic equivalent UV-optical
luminosity
L
UV
=
10
48
erg s
−
1
and duration
t
UV
=
500 s. The dust echo light curves for an on-axis observer are shown in Fig.
6
, for three
different hydrogen number densities of
n
H
=
100, 300, 1000 cm
−
3
. We find that the H-band excess in GRB 130603B is consistent with the dust
echo model. We note that the H-band excess can also be reasonably explained by a kilonova from lanthanide-rich ejecta of mass 0.05–0
.
1M
4
McKee, Parravano & Hollenbach (
2015
) provided an estimate of
ρ
sd
0
.
05 M
pc
−
3
in the solar neighbourhood. Studies of Galactic supernova remnants
show that most supernovae occur near galactocentric radius of
∼
5 kpc (Green
2015
), where the density of stars and dark matter is higher than that of the solar
neighbourhood by a factor of
∼
2.
5
Note that, technically speaking,
A
V
is the extinction for a LOS that goes through the entire disc in the vertical direction. Suppose a bNS merger occurs, on
average, near the disc mid-plane, then the LOS only goes through half of the disc, but inclinations away from the vertical direction makes the path-len
gth
∼
twice
longer than the half thickness.
6
We note that Abbott et al. (
2020b
) inferred
R
bNS
∼
300 Gpc
−
3
yr
−
1
based on a uniform prior of [1, 2.5] M
on the masses of the merging NSs. If most NS
masses are closer to 1.35 M
(as inferred from the Galactic bNS population; Farrow et al.
2019
) which correspond to a much smaller detectable volume than
the 2.5 M
ones, then the rate from Abbott et al. (
2020b
) is likely an underestimate. The bNS merger rate will be accurately measured by the upcoming O4 run
(see Abbott et al.
2018
).
MNRAS
507,
3672–3689 (2021)
Downloaded from https://academic.oup.com/mnras/article/507/3/3672/6354800 by California Institute of Technology user on 03 November 2021
3680
W. Lu, C. F. McKee and K. P. Mooley
Figure 6.
Dust echo light curves viewed from an on-axis LOS (
θ
obs
=
0) at host-galaxy rest-frame wavelength
λ
obs
=
1.2
μ
m (red) and 2.2
μ
m (blue lines),
and for different hydrogen number densities
n
H
=
100 cm
−
3
(dash–dotted lines), 300 cm
−
3
(solid lines), and 1000 cm
−
3
(dotted lines). The red circle and arrow
denotes the H-band (rest-frame
λ
obs
1.2
μ
m) flux excess at
t
obs
=
6
.
9 d and H-band upper limit at
t
obs
=
21
.
8 d from GRB 130603B (Berger et al.
2013
;
Tanvir et al.
2013
), respectively. For all cases, we have fixed the jet opening angle
θ
j
=
4
◦
, and dust-to-gas ratio parameter
n
0
/
n
H
=
1.45
×
10
−
15
(which
is calibrated by the extinction law of a typical Galactic LOS). The emission from the jet is taken to be a square pulse with isotropic equivalent UV-opti
cal
luminosity
L
UV
=
10
48
erg s
−
1
and duration
t
UV
=
500 s. For comparison, we also show the J-band (
λ
obs
1.2
μ
m, red diamonds) and H-band (
λ
obs
2.2
μ
m, blue squares) light curves of the kilonova associated with GW170817 (Villar et al.
2017
).
and expansion speed
∼
0.2
c
, as modeled by the original authors (Berger et al.
2013
;Tanviretal.
2013
). A near-IR spectrum can be used to
distinguish between these two scenarios: a kilonova spectrum has features due to atomic transitions (e.g. Watson et al.
2019
) whereas the
spectrum from dust emission is expected to be featureless. The high grain temperature means that a large number of vibrational modes are
excited and the energy levels occupied are those with very large quantum numbers above the ground state (Draine & Li
2001
), so the emission
spectrum is continuous. The silicate feature near 10
μ
m, due to O-Si-O strentching mode which has a large dipole matrix element, may stand
out above the continuum when the grain is hotter than the excitation temperature of
∼
200 K. We conclude that future kilonova searches in
short GRB afterglow photometric data should be aware of the possible contamination from dust echo, especially when the optical-IR afterglow
spectrum shows evidence for strong dust extinction (
A
V
0
.
5 mag) or when the X-ray spectrum shows large neutral-hydrogen column density
(
N
H
10
21
cm
−
3
).
4 DUST ECHO MODEL
In this section, we construct a detailed model to calculate the light curves of dust echoes as well as the intensity distribution as viewed from
an arbitrary LOS.
To predict the emission from the UV-heated dust, we need a model for the distribution of grain sizes and how they interact with (i.e. absorb,
scatter, and emit) radiation from the infrared to the UV bands. In the last few decades, significant progress has been made in understanding the
size distribution, composition, and various frequency-dependent features in the dielectric functions, based on measurements of the extinction
curve, spectra of reflection nebulae, polarization of starlight, near-infrared emission features, thermal emission from cold dust at longer
wavelengths, as well as
in situ
measurement of interstellar dust in the Solar System (see Draine
2003a
,
2011
, and references therein). However,
since the interstellar medium is dynamic, multiphase, and highly complex, even the ‘best-fitting’ models are only crude approximations to
the real physical situations. In this work, we adopt the following simplified model that captures the main physical processes relevant to our
discussion.
We consider the Mathis, Rumpl & Nordsieck (
1977
, hereafter
MRN
) distribution of grain sizes,
d
n
g
d
a
μ
m
=
n
0
a
−
3
.
5
μ
m
,a
min
<a
μ
m
<a
max
,
(10)
where
n
0
is the normalization in units of cm
−
3
,
a
is the radius of an equal-volume sphere for each grain in units of
μ
m, and we take minimum
size
a
min
=
0
.
01
μ
m and maximum size
a
max
=
0
.
3
μ
m. Much smaller grains
a<
0
.
01
μ
m (e.g. polycyclic aromatic hydrocarbons or PAHs)
are optically thin, evaporate quickly and hence do not contribute significantly to dust echo. Much larger grains
a>
0
.
3
μ
m (if they exist)
make a minor contribution to the extinction of optical and UV radiation from the source. In the long wavelength limit
λ
2
π
a
, the extinction
is dominated by absorption, and the wavelength dependence of the absorption cross-section is given by
C
abs,
λ
∝
a
3
λ
−
2
under the electric
dipole approximation (equation 22.27 in Draine
2011
). At shorter wavelengths,
7
both absorption and scattering may be important, and the
7
In the special case of a spherical grain, the absorption/scattering cross-sections can be computed by the Mie theory, which has been widely used to mod
el
interstellar dust (e.g. Draine & Lee
1984
). However, physical dust particles are non-spherical and are much more difficult to model accurately.
MNRAS
507,
3672–3689 (2021)
Downloaded from https://academic.oup.com/mnras/article/507/3/3672/6354800 by California Institute of Technology user on 03 November 2021
Infrared dust echo
3681
total extinction cross-section
C
ext,
λ
asymptotically approaches twice the geometric area of
π
a
2
in the limit
λ
2
π
a
, as long as the grain is
optically thick. This motivates us to consider the following model for the efficiency factors for absorption and extinction
Q
abs
,λ
≡
C
abs
,λ
π
a
2
1
1
+
(
λ/λ
0
)
2
a
−
1
μ
m
,Q
ext
,λ
≡
C
ext
,λ
π
a
2
1
1
/
2
+
(
λ/λ
0
)
2
a
−
1
μ
m
,
(11)
where
λ
0
is a critical wavelength (typically a few to 10
μ
m) that depends on grain composition and internal structures. This simple model
captures the broad-brush picture of dust absorption/emission and makes the (otherwise complicated) integrals over the grain size distribution
analytically tractable. We adopt
λ
0
2
μ
m which gives Planck-averaged absorption efficiency (see equation 16) in between graphite and
silicate models by Laor & Draine (
1993
) in the temperature range (500–2500 K) of interest for dust echo.
The optical depth due to dust extinction can be integrated analytically
τ
d
,λ
N
H
=
πμ
m
2
n
H
∫
a
max
a
min
d
a
μ
m
d
n
g
d
a
μ
m
a
2
μ
m
Q
ext
,λ
=
2
√
2
πμ
m
2
n
0
n
H
λ
0
λ
(
atan
√
0
.
5
x
max
−
atan
√
0
.
5
x
min
)
,
(12)
where
x
min/max
=
a
min/max
(
λ
0
/
λ
)
2
,
n
H
is the H number density (including atomic, ionic and molecular forms), and
N
H
is the H column
density. For
λ
0
of a few
μ
m, this roughly reproduces the observed extinction curve of
A
λ
/mag
=
1.086
τ
d,
λ
∝
λ
−
1
in the optical and UV
for a typical Galactic LOS with
R
V
=
3.1 (Cardelli et al.
1989
; Fitzpatrick
1999
). We normalize the extinction curve by V-band extinction
A
V
/N
H
5
.
3
×
10
−
22
mag cm
2
at
λ
0.55
μ
m. This gives
n
0
/
n
H
1.45
×
10
−
15
and hence
τ
d,λ
0
.
49
A
V
/
mag
λ
μ
m
(
atan
√
0
.
5
x
max
−
atan
√
0
.
5
x
min
)
.
(13)
In this way, our simplified dust model (equations 10 and 11) gives qualitatively similar extinction curve as the more sophisticated model by
e.g. Weingartner & Draine (
2001
). A comparison is shown in Fig.
B1
in the Appendix. We note that for
A
V
∼
1 mag, the gas-dust column is
highly optically thick for UV photons, which are strongly attenuated by a small fraction of the column density.
At a given time, the specific luminosity
L
ν
is given by afterglow emission from the jet-driven external shocks, and the heating rate for a
grain of size
a
at a distance
r
from the source is given by
̇
q
h
=
∫
ν
max
d
ν
L
ν
e
−
τ
λ
4
π
r
2
π
a
2
Q
abs
,λ
,
(14)
where we consider maximum photon energy of
hν
max
=
50 eV (or
λ
min
0.025
μ
m) because higher energy photons have smaller absorption
cross-sections (see Fig.
B1
). The cooling rate due to thermal emission at grain temperature
T
is given by (using Kirkhoff’s theorem)
̇
q
c
=
4
π
a
2
∫
Q
abs
,λ
π
B
ν
d
ν
=
4
π
a
2
Q
abs
P
σ
SB
T
4
,
(15)
where
σ
SB
is the Stefan-Boltzmann constant,
B
ν
(
T
) is the Planck function, and the Planck-averaged absorption efficiency is defined as
Q
abs
P
≡
π
σ
SB
T
4
∫
Q
abs
,λ
B
ν
d
ν.
(16)
For absorption efficiency in equation (11), we find that the Planck-averaged
Q
abs
P
can be approximated by
Q
abs
P
1
1
+
3
.
5(
T/
1000 K)
−
2
a
−
1
μ
m
,
(17)
which is in between the graphite and silicate models of Laor & Draine (
1993
) and is also similar to the approximation adopted by Waxman &
Draine (
2000
, their equation 13). A comparison is shown in Fig.
B2
in the Appendix. Therefore, the dust temperature
T
(
a
,
r
) at a given time
8
t
is given by the balance between heating and cooling,
̇
q
h
=
̇
q
c
⇒
T
(
a, r, t
)
.
(18)
Using the approximated Planck-averaged
Q
abs
P
, the energy balance can be written into a simple cubic equation whose solution is given by
T
1000 K
=
[
(
2
y
2
3
ξ
)
1
3
+
(
ξy
18
)
1
3
]
1
/
2
,y
=
̇
q
h
[cgs]
7
.
13
a
2
μ
m
,ξ
=
[
(
31
.
5
a
μ
m
)
2
−
12
y
]
1
/
2
+
31
.
5
a
μ
m
.
(19)
The above expression breaks down at extremely high heating rate when
y
>
(31.5/
a
μ
m
)
2
/12, and in that situation, the dust temperature is higher
than 3240
a
−
1
/
2
μ
m
K and the physical consequence is that the grain evaporates almost immediately. In the limit of
y
ξ
2
/12, the second term
in the temperature expression dominates, and we have
ξ
62/
a
μ
m
and
y
3
.
5
L
48
/r
2
19
(ignoring attenuation), so a rough estimate of the dust
temperature is given by
T
2
.
2
×
10
3
K
L
1
/
6
48
r
−
1
/
3
19
a
−
1
/
6
0
.
1
μ
m
.
The grain sublimation rate at a given temperature is determined by the balance between atoms evaporating from the surface and those
returning from the gas phase, which are poorly known (given our limited knowledge about grain composition, surface properties, and shape).
8
This is the retarded time, i.e. the time since the arrival of the first causal signal from the centre of explosion.
MNRAS
507,
3672–3689 (2021)
Downloaded from https://academic.oup.com/mnras/article/507/3/3672/6354800 by California Institute of Technology user on 03 November 2021