Late Engine Activity in Neutron Star Mergers and Its Cocoon: An Alternative Scenario for the Blue Kilonova

Follow-up observations of short gamma-ray bursts (sGRBs) have continuously unveiled late extended/plateau emissions, attributed to jet launch due to late engine activity, the nature of which remains enigmatic. Observations of GW170817 have confirmed that sGRBs are linked to neutron star (NS) mergers, and discovered a kilonova (KN) transient. Nevertheless, the origin of the early blue KN in GW170817 remains unclear. Here, we investigate the propagation of late jets in the merger ejecta. By analytically modeling jet dynamics, we determine the properties of the jet-heated cocoon, and estimate its cooling emission. Our results reveal that late jets generate significantly brighter cocoons compared to prompt jets, primarily due to reduced energy loss by adiabatic cooling. Notably, with typical late jets, emission from the cocoon trapped inside the ejecta can reproduce the blue KN emission. We estimate that the forthcoming Einstein Probe mission will detect the early cocoon emission at a rate of ∼2.1−1.6+3.2 yr−1, and that optical/UV follow-ups in the LIGO-Virgo-KAGRA O5 run will be able to detect ∼1.0−0.7+1.5 cocoon emission events. As an electromagnetic counterpart, this emission provides an independent tool to probe NS mergers in the Universe, complementing insights from sGRBs and gravitational waves.


INTRODUCTION
Short gamma-ray bursts (sGRBs) have traditionally been classified based on the short duration (< 2 s) of their hard gamma-ray emission (Kouveliotou et al. 1993), also referred to as "prompt" emission.sGRBs have theoretically been associated with binary neutron star (BNS; also black hole neutron star: BH-NS) mergers (Paczynski 1986;Goodman 1986;Eichler et al. 1989).In this scenario, the merger produces a compact object surrounded by an accretion disk.This system is short lived, and is considered as the central engine that powers relativistic jets and explain the short nature of sGRBs' prompt emission.Hence, sGRBs have been re-Corresponding author: Hamid Hamidani hhamidani@astr.tohoku.ac.jp garded as a tool to study merging BNS systems in the Universe.
Additionally, BNS mergers have been proposed as a site for the rapid neutron-capture process (r-process) nucleosynthesis thanks to the neutron rich content of NSs (Eichler et al. 1989).Simulations showed that a substantial amount of mass (∼ 10 −3 − 10 −2 M ⊙ ) is ejected during the merger and ideal conditions for r-process are achieved (Shibata 1999;Shibata & Uryū 2000;Hotokezaka et al. 2013;Bauswein et al. 2013;Radice et al. 2016;etc.).Later on, it has been suggested that radioactive decay of freshly synthesised heavy elements can power a unique ultraviolet/optical/infrared transient called macronova/kilonova (KN hereafter) [Li & Paczyński 1998;Kulkarni 2005;Metzger et al. 2010].The KN was expected to take a red color, due to the high opacity of heavy elements synthesized (i.e., lanthanides; Barnes & Kasen 2013;Tanaka & Hotokezaka 2013;Kasen et al. 2013).
In 2017, a BNS merger event (GW170817) was observed for the first time with two LIGO GW detectors assisted by the VIRGO detector (Abbott et al. 2017a).About 1.7 s after the merger, a sGRB (GRB 170817A) was detected (Abbott et al. 2017b).Follow-up observations in X-ray (Margutti et al. 2017), optical (Lyman et al. 2018), and radio (Mooley et al. 2018;Troja et al. 2019) confirmed the presence of a relativistic jet.This further confirmed that sGRBs are produced by BNS mergers (Paczynski 1986;Goodman 1986;Eichler et al. 1989).
In addition to the merger ejecta produced by the merger, there is the sGRB-jet launched by the merger remnant.The sGRB's prompt emission is produced by emission from this jet in the line-of-sight (Ruderman 1975).In NS mergers, mass is ejected dynamically before the merger (tidal interaction) and immediately after the merger (≲ 10 ms; Hotokezaka et al. 2013;Shibata & Hotokezaka 2019).Hence, as shown in Figure 1 (panel B), the jet is initially surrounded by the dense ejecta, and it had to propagate through this dense ejecta to produce the observed sGRB prompt emission (panel C).This is known to produce a jet shock that heats up a part of the ejecta forming a "cocoon" component (in AGNs: Begelman & Cioffi 1989; in collapsars Martí et al. 1997;MacFadyen & Woosley 1999;Matzner 2003;and in sGRBs Nagakura et al. 2014;Murguia-Berthier et al. 2014;Duffell et al. 2015;Nathanail et al. 2020Nathanail et al. , 2021;;Hamidani & Ioka 2021;Gottlieb & Nakar 2022).Hence, the jet, with its large energy budget, is another element to take into account in NS mergers (Nakar & Piran 2017;Kasliwal et al. 2017;Gottlieb et al. 2018).
This prompt jet heated cocoon has been proposed to explain the blue KN in GW170817/AT2017gfo (Piro & Kollmeier 2018;Ioka & Nakamura 2018).Piro & Kollmeier (2018), claimed that AT2017gfo shows evidence for the prompt jet cocoon.However, this contradicts simulation results, showing that for prompt sGRBjets most of the cocoon is trapped inside the ejecta (see Figure 1) and only a small fraction of its mass (< 10%) manage to escape (∼ 10 −4 − 10 −6 M ⊙ ; see Table 2 in Hamidani & Ioka 2023a).This small mass implies that the timescale of emission from the escaped cocoon is too short (∼ 1 h) to explain the blue KN (∼ 1 d), and by the time the escaped cocoon is transparent, the trapped cocoon has adiabatically cooled down to the point of being too dim to be relevant (Hamidani & Ioka 2023a,b).Therefore, parameters of the prompt jet's cocoon assumed in Piro & Kollmeier (2018) are questionable (in particular mass as M ∼ 0.01M ⊙ ).
After the prompt phase, interaction of the prompt-jet with the circumstellar medium produces a softer emission, from X-ray to radio, called "afterglow" (Sari et al. 1998).This emission is expected to decay over time following a simple power-law function (e.g., see Figure 1 in Zhang et al. 2006).However, follow up observations in X-ray show an excess luminosity after the prompt phase in most sGRBs (Barthelmy et al. 2005;Norris & Bonnell 2006;etc.).As shown in Figure 2, these phases (referred to as "extended" and "plateau") are unusually flat over time, in terms of observed luminosity, before decaying (also see Figure 1 in Gompertz et al. 2013; Figure 1 in Kisaka et al. 2017; Figures 1 and 2 in Kagawa et al. 2019; although time evolution can be different for some events: Kaneko et al. 2015;Lien et al. 2016;Zhu et al. 2022).As a result, explaining these phases with afterglow modeling is challenging, and their origin has generally been attributed to late engine activity (Ioka et al. 2005).Kisaka & Ioka (2015) showed that late engine activity can be explained with late time mass accretion (see their Figure 2).In a recent study, Gottlieb et al. (2023) showed using self-consistent end-to-end simulations that late engine activity (extended phase) naturally arises as the disk transitions into a magnetically arrested disk (MAD) state after the prompt phase (see their Figure 1).These late phases have been identified in many sGRBs (∼ 70% of Swift/BAT sGRBs according to Kisaka et al. 2017).Also, recent observations of GRB 211211A and GRB 230307A (long GRBs but with potential kilonova candidate) show evidence of late emission in the form of a soft tail (for GRB 211211A: Rastinejad et al. 2022;Troja et al. 2022;Yang et al. 2022; and for GRB 230307A: Levan et al. 2023;Sun et al. 2023;Gillanders et al. 2023;Yang et al. 2023).Hence, in BNS mergers, after the prompt phase, late engine activity is expected to continue in most cases.
The observed luminosity and duration of both the prompt phase and these subsequent late phases are shown in Figure 2 (right panel) for a sample of sGRBs (taken from: Kisaka et al. 2017 andZou et al. 2018).The observed luminosities and duration show some scatter, however, in general, compared to the prompt emission that shines at ∼ 10 51 −10 52 erg s −1 for ≲ 2 s, the observed extended emission is on the order of ∼ 10 48 −10 49 erg s −1 and lasts for ∼ 10 2.5 s; while the plateau emission is about ∼ 10 46 erg s −1 for ∼ 10 4 − 10 6 s.
Details about this late engine activity are not clear.The flat (or shallow) temporal behaviour of the observed luminosity during these two late phases (extended and plateau) requires extra energy.Energy injection at the decelerating forward shock has been suggested to explain this phase (Zhang et al. 2006;Nousek et al. 2006;Panaitescu et al. 2006).However, this cannot explain cases where this phase is followed by a steep decay (Troja et al. 2007;Liang et al. 2007).Also, this would require extreme radiation efficiencies in the prompt phase (e.g., see Granot et al. 2006;Ioka et al. 2006).Additionally, these phases (the extended emission in particular) show highly variable and flaring features (e.g., see Figure 3 in Barthelmy et al. 2005; Figure 1 in Kisaka et al. 2017).In Matsumoto et al. (2020), using the compactness argument, it has been found that the late emission (extended phase in particular) must have a Lorentz factor Γ ≳ 10 (see their Table 1).For the plateau emission, the situation is more enigmatic.The external forward shock model has been invoked to explain the plateau emission (Zhang et al. 2006; interaction between the prompt and the extended outflow: Matsumoto et al. 2020; also see Lamb et al. 2019 for GRB 160821B; and other scenarios: Oganesyan et al. 2020;Beniamini et al. 2020), but in many cases it is interpreted as late central engine (BH or magnetar in some cases) activity (e.g., Gompertz et al. 2013;Kisaka & Ioka 2015;Becerra et al. 2019;Gottlieb et al. 2023;etc.).Therefore, based on these arguments, late emission phases are typically interpreted as indicators of late engine activity launching jet outflow (Ioka et al. 2005;Kisaka & Ioka 2015;Gottlieb et al. 2023).
Despite observations hinting toward late engine activity in many sGRBs, with a few exceptions, the impact of such late engine activity on the KN has not been subject to thorough investigation (Kisaka et al. 2015(Kisaka et al. , 2016;;Ioka & Nakamura 2018;Matsumoto et al. 2018).In particular, late jet launch is expected to shock and heat-up the ejecta, producing a cocoon, at a much larger radius (than in the prompt jet case).Adiabatic cooling is much less efficient in the case of later activity, and the cocoon is expected to be much brighter at later times.Such emission powered by late engine activity could give hints about the central engine of sGRBs, which is still poorly understood.
In this paper, we investigate late jet launch and propagation through BNS merger ejecta at later times (see Figure 1).We analytically estimate the properties of late jet-cocoons (using Hamidani & Ioka 2021, 2023a), and calculate their cooling emission analytically (using Hamidani & Ioka 2023b).We found that late jets' cocoon can produce a bright blue KN component, especially for wide jets, and for jets that fail to penetrate the ejecta (i.e., choked).We also confirm that such a cocoon is a bright source of soft X-rays immediately after the breakout (Hamidani & Ioka 2023b) [see Figure 1].Our estimates suggest that the cocoon emission will be detected on a yearly basis, and in the near future, emission from the cocoon can be a new set of eyes to probe BNS mergers.
This paper is organized as follows.Analytic modeling for jet propagation in the merger ejecta is presented in §2.Analytic modeling for the cocoon cooling emission is presented in §3.Results are presented in §4, and discussed in §5.Finally, a conclusion is presented in Section 6.Details related to our analytic model can be found in appendices §A, B, C, D, and E.

Overview
As shown in Figure 1, the picture of a BNS after the merger time (t m ) is as follows.A compact merger remnant with an accretion disk is formed; we refer to it as the central engine.This central engine is surrounded by the expanding merger ejecta.The prompt sGRBjet is launched soon after the merger, delayed for a certain time (∼ 1 s for GRB 170817A; Abbott et al. 2017b).This jet is active for a certain duration (∼ 2 s for sGRBs), and is characterized by its opening angle and power.This jet propagates through the merger ejecta.This propagation produces a cocoon.Once the jet reaches the outer edge of the ejecta, both the jet and the cocoon break out.

Gravitational Waves
Figure 1.Schematic illustration of a late jet heated cocoon in BNS mergers, and its cooling emission.After the merger (A), mass is ejected and a sGRB's prompt jet is launched within ∼ 1 s (B).At much later times, ∼ 10 − 10 4 s, jets associated with extended/plateau emission, commonly observed in sGRBs (see Figure 2), are launched (C).These jets (white) shock the ejecta (red), heating it up producing a hot cocoon.There are two cocoon components: i) a component that is fast enough to escape the ejecta, named "escaped cocoon" (purple); and ii) a component that is slower and is trapped inside the ejecta, named "trapped cocoon" (light blue).Cooling emission from the escaped cocoon is expected to produce soft X-ray to UV photons.For the trapped cocoon heated by late jets, cooling emission could be similar to a blue KN component.The r-process heated ejecta and its red KN emission (red) is also shown.
At later times (∼ 1 − 10 s), more mass is ejected (i.e., post-merger ejecta; Fujibayashi et al. 2018).Also, after the prompt phase, late jets are launched in two phases: the extended phase and the plateau phase (see Figure 2).Jet parameters are different in each phase.These late jets are expected to propagate through the ejecta in a manner similar to the prompt jets; forming a cocoon and breaking out to produce the extended/plateau emission (full arguments behind this reasoning are presented in §A.2).A small fraction of the cocoon is fast enough to break out of the ejecta, and the rest is trapped inside the ejecta (Hamidani & Ioka 2023b) [see panel C in Figure 1].
In GW170817/GRB 170817A, late emission has not been observed.It should be noted that this could be due to the large off-axis viewing angle (∼ 20 • − 30 • ≫ Γ −1 jet ).Therefore, a late jet in GW170817/GRB 170817A can not be excluded from the lack of extended/plateau emission.Also, the presence of a late jet, with its smaller energy budget relative to the prompt jet, would not necessarily contradict afterglow (and VLBI; Mooley et al. 2018)  Figure 2. On left: a schematic light curve of an observed sGRB with late engine activity.Prompt (light grey), extended (green), and plateau (blue) phases are shown, in terms of their observed isotropic equivalent luminosities, and time (since the prompt burst trigger time).Dashed lines illustrate the transition between the different phases, which is often interpreted as the end of energy injection (Ioka et al. 2005).On the right: A sample of observed (squares) sGRB's prompt (grey), extended (green), and plateau (blue) events as a function of their observed luminosity and duration (with redshift measurement; taken from Kisaka et al. 2017 andZou et al. 2018).Typical luminosities (and duration), indicated with dashed lines and star symbols, are adopted throughout this paper (assuming ηγ ∼ 0.1; see §3 and Table 1).
calization in GW170817/GRB 170817A is another issue that could explain the lack of EM counterparts related to e.g., late jet breakout.While, the two recent nearby GRBs (GRB 211211A and GRB 230307A) associated with a KN, show evidence of late emission in the form of a soft tail (Rastinejad et al. 2022;Levan et al. 2023; also see Sun et al. 2023).The high energy emission in GRB 211211A, delayed by ∼ 10 3 s, further suggest late engine activity and late jet launch (Mei et al. 2022).
Parameters of the ejecta and the jet are set in the following ( §2.2 and §2.3; respectively).Given these parameters, the dynamics of the jet-cocoon are determined analytically ( §2.4 and §2.5) following the analytic model (and approximations) in Hamidani & Ioka (2021, 2023a,b).In particular, we derive the breakout time [see equation ( 5)], and use it to find the cocoon's total mass and energy (in §2.4), and the mass and energy of the cocoon that escapes the ejecta (in §2.5).

Initial conditions: the ejecta
The time since the merger is taken as t − t m = t (t m = 0 for simplicity).The ejecta is assumed to be expanding homologously: The ejecta is assumed to have a total mass M e .In terms of mass-density, we assume a simple spherical wind-like function ρ(β, t) ∝ r −n with n = 2, and M e = βm β0 4π(ct) 3 ρ(β, t)β 2 dβ; where the maximum velocity is β m = √ 3/5 (so that ⟨β⟩ ∼ 0.2; see Hotokezaka et al. 2013), and the minimum velocity is taken as β 0 = β m /100 ≪ β m .For the prompt jet we take the ejecta mass as ∼ 0.01 M ⊙ .This is comparable to the expected dynamical ejecta immediately after the merger (Hotokezaka et al. 2013).At later times, post-merger ejecta is expected to dominate (∼ 1 − 10 s after the merger; Fujibayashi et al. 2018).Considering measurements of the ejecta mass from the late time KN emission in AT2017gfo, indicating a mass of ∼ 0.05 M ⊙ (also similar to the KNe associated with GRB 211211A and GRB 130603B; Berger et al. 2013;Tanvir et al. 2013;Rastinejad et al. 2022), we adopt the same ejecta mass at the time of late jet launch (see Table 1)1 .

Initial conditions: the jet
The sGRB-engine-powered polar jet can be defined by its opening angle θ 0 ; and its power L j (2L j for two polar jets), or more practically its isotropic equivalent luminosity L iso,0 = 2L j /(1 − cos θ 0 ).We consider three independent engine phases: prompt, extended, and plateau; each producing a distinct relativistic jet.The jet launch time is t 0 ; and the jet launch lasts until t e for a duration t e − t 0 .As a convention, the subscript " 0 " is used to indicate parameters that are related to the central engine.
For simplicity, we assume an initial top-hat structure for the jet as launched by the central engine.This is reasonable considering that any initial complex jet structure (e.g., Gaussian) is expected to converge to a top-hat structure as the jet will strongly be collimated during its propagation across the dense medium (e.g., see Urrutia et al. 2023).It should be noted that this is not inconsistent with observations of GW170817/GRB 170817A, as the more complex "structured jet" is expected to emerge after the jet breakout (also see §2.3.2).Also, during the engine activity (t 0 < t < t e ) at a given phase (prompt, extended, or plateau), we assume that the jet is continuously launched with a constant jet power and constant jet opening angle.
In the following ( §2.3.1;§2.3.2;§2.3.3),these parameters are inferred from typical values observed for these three phases.

Jet launch time and duration
Considering observations of GW170817 (Abbott et al. 2017a) and GRB 170817A (Abbott et al. 2017b), the delay between the BNS merger and the prompt emission should be < 1.7 s (more accurately ≲ 1.3 s according to Hamidani et al. 2020).Considering that jet propagation through the ejecta contributes to this delay, it is reasonable to take the delay between the merger and the jet launch as t 0 − t m = 1 s.Also, based on the bimodal distribution of GRBs, the duration of sGRB's prompt emission (t e − t b ≲ t e − t 0 ) is ≲ 2 s, and with the time it takes the jet to breakout (t b − t 0 ) typically ∼ 0.5 − 1 s, we roughly assume that the prompt jet duration as t e − t 0 = 2 s.
For the later extended and plateau phases, the jet launch and quenching times are not well understood.However, one can roughly use the time of the start of these phases from observations to infer the jet breakout.Observations show that the extended emission starts ∼ 10 s after the prompt phase, and the plateau emission starts ∼ 10 3 s after the extended emission phase (not observed after GW170817, but this is most likely due to the off-axis line-of-sight; see §2.1).Hence, based on this data and considering the time it takes the jet to breakout, we take the starting time of the extended emission's jet as a few seconds after the merger t 0 − t m = 10 0.5 s, and for the plateau phase we take t 0 − t m = 10 2.5 s.In other words, it has been assumed that the engine activity is mostly continuous throughout these three hypothetical jet phases (although the luminosity drops by orders of magnitudes; see Table 1).
Finally, the jet duration can simply be inferred from the observed duration of these phases as it can be seen in right panel of Figure 2. We estimated the average observed duration values and adopted them as "typical" values (as shown with the dashed lines and the star symbols).One exception is that for the prompt phase we adopt a typical duration of 2 s, which is regarded as the typical duration of sGRB jets (Kouveliotou et al. 1993).This is because the breakout time can be comparable to this timescale, making the observed duration of sGRBs much shorter than the actual timescale of the engine activity (Bromberg et al. 2013).

Jet opening angle
First, we assume that the post-breakout jet opening angle θ j eventually converges to the jet's initial opening angle θ 0 .In reality, after the breakout, there is a complex jet structure, composed of a jet-core and the cocoon in the wings (Mizuta & Ioka 2013;Lazzati et al. 2017;Lamb & Kobayashi 2017, 2018;Lazzati et al. 2018;Gottlieb et al. 2021).However, it has been shown that the jet-core can be approximated to have a top-hat structure, and its opening angle is comparable to θ 0 at sufficiently late times after the breakout (Gottlieb et al. 2021;Urrutia et al. 2023 for prompt jets).Hence, here, with the jet opening angle θ j ∼ θ 0 we refer to the jetcore structure.
With measurements of "jet breaks" in the afterglow light curve, the initial jet opening angle can be inferred.In Rouco Escorial et al. (2023), afterglow analysis for a sample of sGRBs shows that while the typical value of the jet opening angle is θ j ∼ 6 • , several sGRBs are associated with narrow jets θ j ∼ 3 • , and the jet opening angle can be as large as θ j ≳ 15 • (e.g., GRB 211106A; see Laskar et al. 2022;Rouco Escorial et al. 2023).With the intention of representing these observed features, we adopt three jet models with different opening angles: narrow (θ 0 = 3 • ), typical (θ 0 = 6 • ) and wide (θ 0 = 15 • ) [see Table 1].It should be noted that, as prompt jets dominate in terms of jet energy, these opening angles (as inferred from observations) are expected to reflect prompt jets.Therefore, the actual opening angles of late jets are not clear.Here, for simplicity, similarity to the prompt jet have been assumed.

Jet luminosity
Not all the jet energy is released in the form of emission.Hence, in order to infer the intrinsic jet power from observations, we introduce radiation efficiency η γ defined as: where E γ,iso is the isotropic equivalent radiated energy, and E K,iso is the (isotropic equivalent) kinetic energy of the jet after photon radiation.
For the extended and plateau phases, the radiation efficiency is not well studied.However, Matsumoto et al. (2020) presented an evaluation of the radiation efficiency of late engine jets as ∼ 0.03 on average, although the standard deviation is large.
Here, for simplicity, we assume the radiation efficiency for prompt, extended, and plateau jets to be the same, and we set it as: (3) Hence, from the above radiation efficiency [equation (2)] the initial luminosity of the jet can be found as a function of the observed luminosity L obs ≈ E γ,iso /(t e − t 0 ) as (using the assumption θ j ∼ θ 0 ; see §2.3.2): Using equation (3), and adopting typical values for the observed luminosity and duration for each of the three phases (see right panel in Figure 2), the typical initial parameters of the jet can be found for the three phases as presented in Table 1.Note that the afterglow modeling requires several free parameters that are difficult to constrain (i.e., n, ϵ e and ϵ B ; jet geometry as well Lamb et al. 2022) and the results are often sensitive to the assumed values.In particular, the estimated values of η γ are highly sensitive to the assumption that all the electrons are non-thermal, i.e., f e ∼ 1 (as We follow the analytic model in Hamidani et al. (2020) and Hamidani & Ioka (2021) to solve the jet propagation, in particular the jet collimation by the ejecta and the cocoon self-consistently.This model has been demonstrated to give analytic results consistent with numerical simulations (see Figures 3, 4, and 5 in Hamidani & Ioka 2021; also see Figure 4 in Hamidani et al. 2020).Eventually, the jet breakout time t b and cocoon properties at t b are determined.
The breakout time since the merger time (t m = 0) in an expanding medium can be written as [for more details see equation ( 44) in Hamidani & Ioka (2021)]: where c is the speed of light, β m is the velocity of the outer ejecta, n is the power-law index of the ejecta's density profile, t 0 is the delay between the merger and the jet launch.A 1 is a constant that can be written as a function of the parameters of the jet and the ejecta, as N s is a calibration coefficient determined through comparison with numerical simulations (N s = 0.46 here)2 ; ⟨η ′ ⟩ = (1 − ⟨β h ⟩)η with η = 0.5 here (equipartition of energy), and is the average jet head velocity3 ; ⟨χ⟩ is another free-parameter to account for the effect of the lateral expansion of ejecta 4 ; L j is the jet power (one side) and θ 0 is its initial opening angle; and M e is the total mass of the ejecta (see §2.2; §2.3; and Table 1 for the typical values).The breakout radius can be found [using equation ( 1)] as: r b /2 is the semi-major axis of the ellipsoidally shaped cocoon.The semi-minor axis r c can also be found (see §A).Hence the cocoon volume at the breakout time can be found as V c = 4π 3 r 2 c r b (for both hemispheres).For more details see §A.
For late time jets, where t b ≫ t 0 , the breakout time can approximately be found as With the typical jet opening angle in the extended emission phase as a reference (model EX-T in Table 1), the parameter dependency of the breakout time can be found as: It should be noted that this is a simplified expression that does not include all dependencies (e.g., η ′ ∝ 1 − ⟨β h ⟩, and χ) and is expected to result in a factor of ∼ 2 error; compared to the accurate expression in equation ( 5) that is used in the following.

The pre-breakout cocoon
The cocoon is distinguished from the relativistic jet as being slower, having Γ < 10, and from the ejecta as being heated by the jet shock (Hamidani & Ioka 2023a).As explained, before the breakout, the cocoon is assumed to take an ellipsoidal shape.Also, the pressure inside the cocoon P c is assumed to be dominated by radiation pressure (adiabatic index of 4/3), and roughly constant throughout the cocoon's volume.
With these assumptions, the cocoon mass (M c ), energy (E c ), and internal energy (E c,i ) can be found until the breakout following the analytic model in Hamidani & Ioka (2023a).Details are given in §B.

Analytic model: post-breakout
We follow the same model in Hamidani & Ioka (2023a) to determine the fraction of the escaped cocoon mass (M es c /M c ) and internal energy (E es c,i /E c,i ) analytically.The main argument is that, after the jet breakout, the cocoon fluid accelerates (due to its highly thermalized energy content), until it reaches its terminal velocity determined by Bernoulli equation as β = 1 − Γ 2 inf , where Γ inf ≈ hΓ (h is enthalpy).This terminal velocity is assumed to be reached roughly at t 1 ∼ 2t b (∼ 1.5t b or shorter, for the plateau jet case where t b is significantly long).Then, the cocoon fluid is designated as escaped or trapped inside the ejecta based on its velocity, so that: As the cocoon expansion converges to homologous (reached at t 1 ), the escaped cocoon's mass and energy can be found using the above condition (see Hamidani & Ioka 2023a).Analytic results are given in §B (also, for more details see Hamidani & Ioka 2023a).
The escaped cocoon mass (M es c and energy) distribution is approximated as conical (see Figure 4 in Hamidani & Ioka 2023a), with an opening angle θ es c [see equation (B10)], and roughly with an axial symmetric distribution.
The escaped cocoon is separated into a relativistic part and a non-relativistic part, and each is modeled (in terms of mass-density and internal energy density) with a single power-law function (see §B; and Hamidani & Ioka 2023b).
For the trapped cocoon (M tr c ∼ M c ≫ M es c ), the mass-density is approximately considered as the same as that of the ejecta.

Main approximations
The cocoon is assumed to radiate as a blackbody in thermal equilibrium.We employ the diffusion approximation, where, for a given shell, once the diffusion time becomes comparable to the dynamical time, the internal energy is approximated as instantly transforming into photons that are emitted (as in Nakar & Sari 2012;and similarly in Kisaka et al. 2015).Radiation behind this shell is considered as unable to diffuse out.
We use the grey opacity to estimate the optical depth.For a typical range of density and temperature of the considered cocoon fluid (in the comoving frame), we solved the ionization degrees of the fluid under the assumption of local thermodynamic equilibrium.We confirmed that the ionization edge is always located at an energy higher than hν peak by at least a factor of 3; where h here is Planck's constant and ν peak is the frequency at which the spectrum peaks following Wien's displacement law.Hence, the strong attenuation by the boundfree process (e.g, see Figure 8 in Metzger 2019) is not expected to affect the bulk of the emission discussed in this paper.Thus, we adopt the same (grey) opacity both in UV and soft X-ray range, as κ = 1 cm 2 g −1 (as for lanthanide-free material because the cocoon outflow mostly originates from the polar region; see Shibata et al. 2017;Fujibayashi et al. 2018), which is derived from the detailed calculations of bound-bound opacities (Tanaka et al. 2020;Banerjee et al. 2020Banerjee et al. , 2023)).
For late jets, the mass-density at the jet shock is significantly lower compared to the case of prompt jets.Therefore, we estimated the time required to realize the photon blackbody spectrum from the thermalized cocoon fluid.Considering free-free process for electrons, and a composition of fully ionized iron-peak elements (Z ∼ 30), this time can be estimated as the ratio of the internal energy density of the fluid aT 4 over the free-free cooling rate Rybicki & Lightman 1979).This timescale is found, at most, as comparable to the breakout time.Therefore, our assumption of blackbody spectrum is reasonable; especially considering that bound-free (and boundbound) process is even more effective at thermalizing the shocked cocoon fluid.
Finally, as suggested in Granot et al. (2023), radiation-mediated shocks can affect the chemical composition (via fusion or fission) and the opacity in the jet-shocked cocoon.However, the significance of this effect is not fully understood.

Cooling emission from the escaped cocoon
We use an updated version of the model presented in Hamidani & Ioka (2023b), which was intended for the Table 1.Typical parameters of prompt emission's, extended emission's, and plateau emission's jets, based on observations (see Figure 2).In total, nine models are considered, with three jet models for each of the three phases (for more information see §2.3).Breakout times and escaped cocoon mass fractions, calculated using our analytic model, are also shown (see §A; and B).case of prompt jets.In the case of late jets, the breakout time (and radius) can be large and non-negligible (e.g., compared to the photospheric radius of the emitting shell).Therefore, in the following, quantities are expressed using alternative coordinates for laboratory time, radius, and observed time, with the breakout time (t b ), breakout radius (r b ), and observed breakout time [t b,obs = (1 − β m )t b ] as their reference, so that: ∆t = t − t b , ∆r = r − r b , and ∆t obs = t obs − t b,obs (respectively) [see §C].

Jet opening
With the cocoon mass determined (see §2.5 and §B), the optical depth and the velocity of the diffusion shell as a function of the observed time can be calculated [see equation (C16)].Also, with the total internal energy determined (see §B), the internal energy density can be calculated.
Then, cooling emission can be found as the product of i) the leaking internal energy over time (i.e., power due to the change of the diffusion shell's velocity), ii) the effect of the observer's time, and iii) the effect of relativistic beaming as: is the solid angle of the emission coming from the shell moving with Γ d β d .
With the photospheric radius found as: the observed temperature can be found using Stefan-Boltzmann law (in the relativistic limit) as where σ is the Stefan-Boltzmann constant.
It should be noted that, here, by cooling emission we refer exclusively to jet heating as the source of internal energy.More details are given in §C (for details on other types of emission, see Gottlieb & Loeb 2020;Hamidani & Ioka 2023b).

Cooling emission from the trapped cocoon
As discussed in Hamidani & Ioka (2023b), for the prompt jet case, the trapped cocoon [defined as having β ⩽ β m ; see equation ( 8)] is not a relevant component for the cocoon cooling emission.This is because i) internal energy is deposited at an early times (∼ t b ∼ 1 s after the merger), and ii) as it is concealed by both the escaped cocoon and the ejecta, by the time photons can escape out of the trapped cocoon (usually at ∼ 10 3 s; as the escaped cocoon becomes optically thin), its internal energy would have been significantly reduced due to the adiabatic expansion; roughly by a factor of ∼ t b /t ∼ 10 −3 .However, the situation is fundamentally different for the late jet case; as the jet launch and breakout is expected at much later time (t b ≫ 1 s), and the reduction of internal energy due to adiabatic cooling is much less severe.
The same method used for the cooling emission of the escaped cocoon is followed here.Notable differences are: i) the power-law index of the mass-density profile is approximated here to the initial mass-density profile of the ejecta (n = 2); ii) as the velocity of the trapped cocoon (and the ejecta) is significantly smaller (β < β m ), the approximation β ≪ 1 is used, allowing us to express the observed time and the laboratory time much simply as t obs ∼ t; and iii) the trapped cocoon's structure is assumed as conical for simplicity.Assuming a simplified conical shape for the trapped cocoon, we found that its overall opening angle [defined by the mass fraction; see equation (D23)] is larger than that of the escaped cocoon (due to the large mass fraction), particularly for wide jets, making it easier to detect.The detailed analytic model is presented in §D.

R-process powered KN emission
Here we refer to all emission powered by r-process elements' radioactive heating as "KN emission" (regardless of where it is taking place; e.g., even if it is within the cocoon).
In short, the ejecta, expanding homologously, is defined with its mass M e , power-law index of its density profile n, its maximum velocity β m , and its minimum velocity β 0 .The optical depth can then be calculated from the shell moving with β m and inward (assuming grey opacity and neglecting the escaped cocoon's contribution).Then, as the r-process heating can be modeled with a power-law function, the energy deposition is found, allowing us to determine the diffused internal energy (i.e., luminosity) in the same manner as in §3.2 and §3.3 (also see Kisaka et al. 2015;Hamidani & Ioka 2023b).The full analytic model for the KN emission is presented in §E.

Dynamics
In Table 1 we show the jet breakout time for our sample of models.We found that, typically, the jet breakout time is of the order of ∼ 1 s (prompt jets), ∼ 10 s (extended emission's jets), and ∼ 10 3 s (plateau emission's jets).These breakout time are roughly consistent with the observed timing of late emission, considering that the actual "observed" delay is shorter by a factor of 1 − β m (∼ 0.65 here).This delay can easily be evaluated from observation of sGRBs based on the start of the late emission phase (or based on the start of the cocoon emission; see §4.2) and could be useful to evaluate the parameters of the jet and the ejecta [see equation (7)].
We find that plateau jets (if they exist) are inherently difficult to break out.This is different from the results of Matsumoto & Kimura (2018), where jet propagation did not take into account the effect of the expanding ejecta on the cocoon and on the jet collimation (see Hamidani & Ioka 2021).In one case (model PL-W), the jet failed to breakout (i.e., choked).Therefore, it is difficult to explain the observed fainter plateau emissions (≪ 10 46 erg s −1 ; see right panel in Figure 2) with this scenario, unless the ejecta mass is somehow significantly smaller than what has been suggested (∼ 0.05 M ⊙ ) for AT2017gfo (e.g., in cases of equal mass BNS where the disk/ejecta mass is small; see Kiuchi et al. 2019); or the prompt/extended jet is wide and strong enough to permanently disrupt the ejecta, leaving a cavity for the plateau jet to propagate through (although this is not trivial due to the trapped cocoon; see §A.2).
Table 1 also shows the fraction of the escaped cocoon mass for each model is shown.These mass fractions have been calculated using equation (B7) [also see equations ( 5) and ( B6)].The larger α (= Energy deposited by the jet/kinetic energy of the ejecta) is, the larger the fraction of the escaped cocoon should be.
Two trends are apparent.First, the wider the jet opening angle, the larger the mass fraction of the escaped cocoon.This is due to the larger jet energy (∝ θ 2 0 ), longer breakout times (also ∝ θ 2 0 ), resulting in a larger jet energy-density deposited into the cocoon, boosting its velocity (more details in Hamidani & Ioka 2023a).Second, the later the jet phase is, the much smaller the mass fraction of the escaped cocoon is.This trend has not been discussed before, and is mainly explained by the smaller energy budget [L iso,0 (t e − t 0 )] of late jets (1/2, and 1/20 that of the prompt jet, for the extended and plateau jet, respectively), while the kinetic mass of the ejecta (∝ M e ) is larger at later times due to post-merger mass ejection.Hence, for late jets, most of the cocoon mass (and energy) is concentrated in the trapped part, making its contribution to the cocoon emission more significant.

Components of the cocoon's cooling emission
In Figure 3, cocoon cooling emission calculated using our analytic model for a late jet model (EX-W in Table 1) is presented for reference.There are three cocoon components to consider: relativistic (escaped), non-relativistic (escaped), and trapped.
Bolometric luminosities have been calculated for each of the three cocoon components separately (see §C and §D, respectively, for more details), as shown in Figure 3.Then, contributing components are connected to produce the overall light curve.
Transition time from the relativistic cocoon to the non-relativistic cocoon is set with the optical depth.Transition time from the non-relativistic to the trapped cocoon is set with the luminosity; i.e., once the trapped cocoon's luminosity dominates over the non-relativistic cocoon's luminosity, we assume that transition takes place.This is a crude simplification.One alternative is to rather use the optical depth to determine the time of this transition.This alternative procedure would produce a bump feature in luminosity at the transition time (sudden rise in luminosity due to the trapped cocoon's contribution; e.g., at t obs ∼ 3000 s in Figure 3).However, this choice overlooks the pressure gradient at β m that may result in internal energy transferred from the outer trapped cocoon to the inner non-relativistic cocoon.Future numerical simulations are needed to clarify this point.
As pointed out in §4.1, the start of the cocoon cooling emission from a late jet is delayed from the merger time and from the prompt phase.This is one major difference compared to cocoon emission from prompt jets (also see Figures 4, 5 and 6).
The very early emission comes from the escaped cocoon (relativistic or non-relativistic).It is the brightest part of the emission, although quite short; ∼ 10 3 s.As explained in Hamidani & Ioka (2023b) the luminosity, then, drops roughly as ∝ t −2 obs for prompt jets.But considering the late start of the emission here, time dependency of this early emission can be generalized rather as ∝ ∆t −2 obs .Internal energy available for cooling emission can be found as There are two important effects to highlight: i) for late jets the average jet head velocity ⟨β h ⟩ ≈ t b t b −t0 β m is much smaller (∼ β m ) compared to late prompt jets (∼ 2β m ); Ioka & Nakamura 2018;Hamidani et al. 2020) [see Table 1] making late jets more efficient at heating the ejecta; and ii) the energy (and mass; as explained in §4.1) fraction of the trapped cocoon is significantly larger for late jets compared to prompt jets [by a factor of ∼ 3 for typical jet models; see equation (B9)].The combination of i) and ii) boosts the internal energy of late jet's cocoon (e.g., by ∼ 2 × 3 for the extended phase).Hence, the internal energy of the trapped cocoon at the breakout time for the prompt jet case is roughly comparable to that of the late jets, despite the large difference in jet luminosity.
Not all this energy is converted to photons.One crucial point to highlight is that, in the cocoon, there are two competing processes for internal energy: i) adiabatic cooling and ii) cooling emission.Due to the late breakout time, adiabatic cooling ∝ t b t [see equation (C17)] is much less efficient (than in the case of prompt jet's cocoon) making the overall cocoon emission more luminous for later jets.
As explained in Hamidani & Ioka (2023b), the trapped cocoon is not relevant for prompt jets; as by the time the escaped cocoon is optically thin (at ∼ 10 3 s) only a small fraction (∼ t b t ∼ 1 1000 ) is left to power the cooling emission.This is much different for late jets as the adiabatic loss (∼ t b ∆t ) is orders of magnitude less severe [see equation (C17)].
Also, the cocoon emission from late jets (trapped in particular) is released at times when the r-process powered KN emission is dimmer (compared to prompt jets).The bright early afterglow emission (from the prompt jet) is also dimmer at late times.Hence, the late-jet's cocoon suffers less from contamination and is easy to identify (compared to prompt jets).
Cooling emission from the trapped cocoon dominates in terms of bolometric luminosity at later times.As shown in Figure 3, for EX-W, the trapped cocoon shines at decent luminosities (∼ 10 43 erg s −1 at ≳ 10 3 s), much higher than the early KN.Considering that the trapped cocoon emission is blue (see §4.3), this model offers an alternative explanation to the origin of the early blue KN in AT2017gfo (Drout et al. 2017;Villar et al. 2017;Waxman et al. 2018).

Observational properties
Figure 4 presents our analytic results for the jet models in Table 1.Cocoon emission from prompt jets, the r-process powered KN emission (see §E), and observations of the early KN in AT2017gfo (Villar et al. 2017;Waxman et al. 2018) are shown for comparison (in grey solid lines, grey dotted lines, and grey symbols, respectively).

Bolometric luminosity
Bolometric luminosities in Figure 4 (top) show two general trends in regards to the cocoon cooling emission: i) the later the jet (launch) phase is, the brighter the emission; and ii) the wider the jet opening angle is, the brighter the emission [assuming constant observed isotropic equivalent luminosities; see equation (4) and Table 1].To our best knowledge, this is the first time that the first trend (i) has been found.The second point (ii) has already been discussed [see equations ( 59) and (60) in Hamidani & Ioka 2023b].
The delay between the merger time and the start of the cocoon cooling emission is t b,obs = (1 − β m )t b ; it is a simple function of the outer ejecta velocity β m , and the jet breakout t b ; also a function of the ejecta and the jet parameters [see equation ( 7)].Since it is easy to evaluate such a timescale from observations, it can prove useful to evaluate the time it took the jet to break out (t b ) and the ejecta velocity (β m ).Hence, a combination of the parameters of the jet and the ejecta can also be constrained with this timescale [e.g., using equation ( 7)].
The early luminosity decreases rapidly.This is due to its dependency on the time since the breakout as ∝ (t obs − t b,obs ) −2 = ∆t −2 obs ; this is because E c,i ∝ 1/∆t due to adiabatic cooling and L ∝ E c,i /∆t obs [see §C and equation ( C21)].This dependency is consistent with Hamidani & Ioka (2023b) [for prompt jets, where t obs − t b,obs ∼ t obs ].
For late jets, the cocoon luminosity shows a late time excess, powered by the trapped cocoon.Contribution of this component is much less significant in prompt jets (due to the reasons mentioned in §4.2).This component is not as bright as the early cocoon emission, coming from the escaped cocoon, but it lasts much longer thanks to the larger trapped cocoon mass.This feature would make it easier to identify (see Figure 5) especially in wide jet cases, where it is more luminous than the KN; this is particularly helpful in future GW170817like events where follow-up observations are expected to be delayed due to the time it takes to cover the region of the sky determined by the GW signal.For plateau jets, this component can be extremely bright, even overshadowing the r-process powered KN.Therefore, this component can offer a natural explanation to the blue KN in AT2017gfo (Drout et al. 2017;Villar et al. 2017;Waxman et al. 2018).It can also offer an explanation to the KN candidates, associated with certain sGRBs, that somehow seem to have a luminous blue KN component (compared to AT2017gfo) and a relatively typical red KN component [e.g., GRB 050724 (Gao et al. 2017), GRB 060614 (Jin et al. 2015;Yang et al. 2015;Rossi et al. 2020), GRB 070714B (Gao et al. 2017), GRB 070809 (Jin et al. 2020), and GRB 150101B (Troja et al. 2018); also see (Gompertz et al. 2018), and Table 2 in (Rossi et al. 2020) for a summary].
Emission from the trapped cocoon, with its delayed start and long timescale, contributes to making the cocoon emission from late jets easier to detect and to exploit, as a tool to probe BNS and sGRBs (also its color; see §4.3.3).Also, the trapped cocoon is wider than the escaped cocoon, and can be as wide as ∼ 60 • − 70 • (for wide jets), making it even more accessible (see §5.3.1 and §D).Considering the model PL-T (see Table 1), the parameter dependency of this trapped cocoon emission can be found as: where t obs ∝ (κM c ) 1 2 [see equation (D25)], t obs ∼ 2t b is roughly considered as the starting time, and p is the power-law index of the time dependency of this luminosity, with p ∼ 1.1 − 1.3 at the early phase of this emission; roughly between t obs ∼ 2t b and 20 × t b .This shallow power-law index p (∼ 2 for the escaped cocoon; see Hamidani & Ioka 2023b) is due to the less severe adiabatic cooling at early times (as ∆t1 ∆t ∼ 1 for t ≳ t 1 ∼ 2t b ).For more details see §D and equation (D29) in particular.Dependency of this luminosity can be used to estimate the combination L iso,0 θ 2 0 ∝ L j (i.e., jet power), if this emission is observed (i.e., both L bol and t b ).

Photospheric four-velocity
Photospheric four-velocity is shown in Figure 4 (middle panels).The late time evolution is similar for all jet models, as the photosphere moves inward, across the escaped cocoon, and to the trapped cocoon.As explained in Hamidani & Ioka (2023b), the photospheric velocity can be inferred from observations, and is potentially a powerful tool to reveal the relativistic nature of the cocoon.Also, it can be used to, in theory, measure the product (κM e ) 1 7 .Although it is challenging to break the apparent model degeneracy by only relaying on the observed photospheric velocity at later times, as it can be seen from Figure 4; distinguishing between prompt, extended, and plateau jets is feasible provided that early observations (≲ 10 2 s after the jet breakout) of the photospheric velocity are available5 .
At early times, and immediately after the breakout, the photospheric velocity drops from Γ ph β ph ∼ 10 to ∼ 1 in ∼ 10 2 − 10 3 s.This unique evolution of Γ ph β ph is a very useful quantity to monitor in order to identify the cocoon and extract information about its environment (e.g., parameters of the engine/jet and the ejecta).High cadence observations are ideal to monitor this fast evolution.

Observed temperature
With our blackbody approximation, the color of the cocoon emission is determined by the observed temperature.In Figure 4 (bottom panel), T obs is plotted as a function of the observed time since the merger t obs − t m,obs .The evolution of T obs is determined by the time since the breakout t obs − t b,obs .The difference in t obs −t b,obs causes an apparent difference between T obs for prompt and for late jet models in Figure 4. Apart from this, the overall color of the cocoon emission from prompt and late jets is intrinsically similar.
In all cases, immediately after the breakout, for a few tens of seconds, the observed temperature from the escaped cocoon is so high (≳ 10 6 K), that blackbody emission is expected to peak in soft X-ray.The temperature dependency, of the peak of the emission, can be expressed using Wien's displacement law as hν peak ≈ 0.24 keV T obs 10 6 K .
As explained in §3.1, strong attenuation by the boundfree process was found as located at an energy higher than hν peak (by at least a factor of 3) hence its irrelevance.The idea that early cocoon emission can power a bright soft X-ray transients has already been presented in Hamidani & Ioka (2023b).This X-ray emission part is discussed in more details later in §5.2 (also see Figure 6).At later times, emission from the escaped cocoon is expected to peak at UV bands (∼ 10 4 − 10 5 K) [see Figure 5].
In the case of prompt jets, after the escaped cocoon becomes optically thin at late times (≳ 10 3 − 10 4 s), the trapped cocoon's emission is irrelevant due to its extreme faintness (e.g., relative to the r-process powered KN).However, for late jets, the trapped cocoon is brighter and lasts longer (due to the arguments presented in §4.2).Temperature of the trapped cocoon part is on the order of ∼ a few 10 4 K, making it bright in optical bands too.This is comparable to the temperature range of the blue KN (e.g., in AT2017gfo), and suggests that the trapped cocoon component, heated by late engine activity, can be one possible explanation of the blue KNe component (in AT2017gfo; Kasliwal et al. 2017;Drout et al. 2017;Villar et al. 2017).Also, being bright in optical makes this cocoon emission very accessible from an observational point of view (compared to UV observations) as optical facilities are more abundant, and optical emission is less affected by extinction.
In all cases, the fast temperature evolution shows the multi-messenger and unique nature of the cocoon emission in sGRBs; as it shifts from soft X-ray to UV, and then to optical bands in a matter of hours.This makes it easy to identify, although high cadence observations are required.In late jets, this emission come at later times (e.g., after the merger) making it even more accessible.

UV/Optical signal
In Figure 5 we show the magnitudes (absolute, and apparent for a distance of 40 Mpc as for GW170817/AT2017gfo) of the cocoon emission, in comparison to that of the r-process powered KN and observations of AT2017gfo.
The pattern of late engine-powered cocoon emission, compared to that of the prompt cocoon, show three main differences.First, the emission is delayed due to the late jet breakout.Second, the early part (from the mildly relativistic escaped cocoon) is slightly dimmer.Third, the cocoon emission is still luminous at late times (∼ 10 3 −10 5 s; ∼ 1 h -1 d), thanks to the contribution of the trapped cocoon (due to the same arguments presented in §4.2 and §4.3.1).
For the plateau jet case (also for extended wide jet), in terms of magnitude, the late time peak (trapped cocoon) coincides with the peak of the r-process powered KN component (grey dotted line in Figure 5).Again, this shows that the late-jet's trapped cocoon is capable of explaining the observed blue KN component by itself.With the wide plateau jet, this component is far brighter, mainly due to the fact that the jet failed to breakout (t b > t e ).This case is already excluded for AT2017gfo, but could be a possible explanation of KNe candidates claimed by Rossi et al. (2020) that have typical red KN components but their blue KN components are exceptionally bright (e.g., sGRBs 050724, 060614, and 070714B; also see Jin et al. 2015;Gao et al. 2017;Troja et al. 2018;Gompertz et al. 2018;Jin et al. 2020 and the references within).(Brown et al. 2010; for UVM2 band) and (Schlegel et al. 1998; for g band).A jet with parameters similar to those of the typical plateau jet model (PL-T) appears to reproduce the blue kilonova (KN) of AT2017gfo effectively.There is a slight difference between the analytic models in the UV band and the g band, e.g., in comparison to AT2017gfo.This discrepancy, although the emission process is a blackbody in both cases, arises from our use of a simplified, non-wavelength-dependent mean grey opacity (κ = 1 cm 2 g −1 ; see Banerjee et al. 2023).
Comparison with the r-process powered KN (in Figure 5), and considering limiting magnitudes of different instruments (in UV: UVOT and the upcoming ULTRA-SAT; in optical: ZTF and the upcoming LSST), suggest that cocoon emission from wider late jets dominates the KN for longer timescales, and can easily be detected.For typical jet opening angles, the emission is dimmer, especially for the extended emission's jet case, as the trapped cocoon is not powerful enough to compete with the r-process powered KN.For narrow jets, the cocoon emission seems extremely difficult to detect.
Nakar & Piran (2017), followed by Gottlieb et al. (2018) and Hamidani & Ioka (2023b), explained that the cocoon emission is mainly a UV transient.This is confirmed here, in more specific details; while the escaped part of the cocoon shines in UV, in the form of a fast transient, the trapped cocoon part (only for late jets) is visible in optical at later times for longer.
These findings stress the importance of a wide field UV survey, such as the ULTRASAT mission (∼ 200 deg 2 ; Sagiv et al. 2014); which could detect the cocoon emission without the need of a GW signal.Optical wide field surveys (e.g., LSST) are also important.However, high cadence (≲ 10 3 s) is necessary, otherwise, the detection would come at later times when the uniqueness of the cocoon emission (e.g., in terms of the fast photospheric evolution) could be missed, making it difficult to identify due to contamination by other fast optical transients (e.g., peculiar ultra-stripped supernovae Ho et al. 2023;Agudo et al. 2023; also, galactic dwarf novae, and variable quasars Oshikiri et al. 2023).

X-ray cocoon signal
In Figure 6 we present the cocoon emission in soft X-ray.This X-ray luminosity has been calculated by integrating Plank law in the range 0.1 − 3 keV.Note that the possibility of strong attenuation by the boundfree process has been considered, but we found that the effect is weak in the escaped cocoon at the considered soft X-ray range (as discussed in §3.1).
Depending on the engine phase, the onset of the Xray emission is different relative to the merger time (top) [or to the prompt emission time].Emission since the jet breakout time (bottom) show that the X-ray emission is similar for the different models; with two exceptions: the wider the jet opening angle the longer the emission (due to the large escaped cocoon mass), and the model PL-W is an outlier due to the different fate of its jet (failed; i.e., no relativistic cocoon) and the non-relativistic cocoon being the emitting component, making its time evolution slower.
This shows that the cocoon breakout can trigger a soft X-ray transient (in addition to emission in a multitude of other bands; see Figure 5).This X-ray emission is different from other discussed emissions that can be observed in X-ray, such as shock breakout, and cocoon afterglow (Nakar & Sari 2010;Nakar & Piran 2017).Also, this emission is similar in nature to that expected in the case of long GRBs (De Colle et al. 2018), although here it is fainter (and shorter) by about one order of magnitude (each) due to the weaker jet and the less massive ejecta.
The upcoming Einstein Probe (EP) mission (expected in 2024; Yuan et al. 2022a) in particular, with its very wide field of view (∼ 3600 deg 2 ) would be ideal to detect such transient, independently of GW/GRB signals (also the HiZ-GUNDAM mission planned for the next decade, Yonetoku et al. 2020; SVOM telescope with a FoV of ∼1 degree (Atteia et al. 2022); and KOYOH, small satellites with wide-field X-ray launched recently6 ).With EP's deep sensitivity, detection could be possible for distances much larger than 40 Mpc.Spectroscopic analysis will allow us to identify it as cocoon emission; given that the spectrum is blackbody with its peak shifting to lower frequencies with time (see Figure 4; for more details see the dependency of T obs in Hamidani & Ioka 2023b).
Another relevant soft X-ray detector for this emission is the all sky monitor MAXI (Matsuoka et al. 2009).In fact, as argued in Hamidani & Ioka (2023b), such transients might have already been detected by MAXI in the past but remained unidentified (MAXI Unidentified Short Soft Transient;e.g., D'Ai et al. 2016;Morokuma et al. 2015;Honda et al. 2015;Uchida et al. 2014;etc.)detected by MAXI (in 2 − 10 keV), but most likely missed by Swift/BAT (only sensitive to harder Xray, 15 − 50 keV; same for Fermi/GBM) and also missed by late Swift/XRT follow-ups due to their fast decay.This is supported by data indicating that the rate of MAXI's X-ray transients is three times higher than the GRB rate, which might suggest that contribution from off-axis GRB events have been captured (although contamination from galactic transients is another possibility; Serino et al. 2014).As a note, with MAXI's X-ray band and sensitivity, the duration of this emission is expected to be shorter than in Figure 6.Also, this X-ray emission could be present in the sample of fast X-ray transient candidates observed by Chandra in 0.5 − 7.0 keV (Bauer et al. 2017;Quirola-Vásquez et al. 2022, 2023) although the duration of most of the candidates is much longer (∼ 10 3 −10 4 s; with a few exceptions, e.g., XRT 191127).Also, it is worth noting that this X-ray emission can also explain some X-Ray Flashes (XRFs; Yamazaki et al. 2002;Sakamoto et al. 2005Sakamoto et al. , 2008)).
Our modeling of this X-ray emission is simplistic, assuming LTE, grey opacity, etc.These are crude but fairly reasonable considerations (as discussed in §3.1).Also, extinction in the host galaxy has not been taken into account and is assumed to have a limited effect on this emission (mostly Galactic as BNS mergers are expected to occur at the outskirts of their host galaxies).Lastly, one challenge related to detecting this X-ray emission is contamination by early jet afterglow (nonthermal) emission.Therefore, ideal events are sGRBs seen with a LOS slightly outside the jet core.

Future detection
Here we estimate the number of cocoon-powered events that could be detected in the near future.For the rate of BNS mergers, we take ∼ 320 +490 −240 Gpc −3 yr −1 as our fiducial value (Abbott et al. 2021; based on the first half of O3; also see Mandel & Broekgaarden 2022).We roughly assume that every BNS merger will launch a jet (as in GW170817/GRB 170817A; Abbott et al. 2017b).We consider the cocoon from the wide jet, as our model  1).Sensitivity of MAXI (∼ 10 −9 erg cm −2 s −1 per orbit in 2 − 10 keV) and Einstein probe's WXT (∼ 1.2 × 10 −10 erg cm −2 s −1 for 100 s in 0.5 − 4 keV) instruments are shown (Matsuoka et al. 2009;Yuan et al. 2022b).The starting time of this cocoon emission is set as the start of the free-expansion phase t1 (see §C).
(with θ 0 ∼ θ j ∼ 15 • ); and we take the probability of such wide jet among sGRBs as ∼ 0.2 (see Figure 4 in Rouco Escorial et al. 2023; also see Fong et al. 2015).We consider the case of late jet launch after the prompt jet; we assume that there is a probability of ∼ 0.7 for this to happen (following Kisaka et al. 2017 results suggesting that ∼ 70% of Swift-sGRBs show signs of late emission).

Detection in soft X-ray
Let us estimate the detection rate of the cocoon emission by Einstein Probe (EP) mission.The opening angle of the escaped cocoon is on the order of θ es c ∼ 30 • for the wide jet model [see equation (B10); also see simula-tion results in Hamidani & Ioka 2023a].This is also the expected opening angle for the soft X-ray emission (as θ c > 1/Γ).This reduces the number of detectable events by a factor [1 − cos(θ es c )] ∼ 0.13.Then, considering the sensitivity of EP (see Figure 6), and considering late jet models (for D = 40 Mpc), the cocoon emission can be detected up to ∼ 25 times further away (up to ∼ 1 Gpc away in an optimistic scenario) 7 .EP's field of view as a fraction of all sky reduces the number of events by the factor ∼ 3600 41253 ∼ 0.087.The yearly number of events to be expected can be found as the product of, BNS merger rate (∼ 320 +490 −240 Gpc −3 yr −1 ), volume (4π/3 Gpc 3 ), the probability of sGRB-jet launch in a BNS merger (taken as 1), the probability of a wide jet launch (∼ 0.2), the probability of a late jet launch (∼ 0.7), the probability of the cocoon emission in the line-of-sight (∼ 0.13), and the probability of the emission being in the field of view of EP's WXT instrument (∼ 0.087).Hence, the rate of yearly detection of this soft X-ray emission can be found as: This detection rate is equivalent to ∼ 1/600 BNS mergers within 1 Gpc being detected.Considering the more recent, but much less constrained, estimate of BNS event rate (as 10 − 1700 Gpc −3 yr −1 ; Abbott et al. 2023), the yearly detection rate will take a wider range of values as 0.07 − 11 yr −1 .The BNS event rate has also been estimated using knowledge of the jet structure and GRB emission as 2 − 680 Gpc −3 yr −1 (Hayes et al. 2023), but more observations of BNS mergers (and sGRBs) are needed to further constrain it.Hence, we expect that EP will be able to independently detect cocoon cooling emission in soft X-ray (also HiZGUNDAM and other soft X-ray instruments).With the properties of this emission known (see §5.2), identification will be possible.Such detection will have several important uses, such as constraining the rate of BNS merger and the fraction of sGRB-jets that are choked (failed).Such detection could also be used to trigger follow-up campaigns and identify KNe independently.Therefore, this detection process is important as it offers us a new channel to probe BNS mergers in the Universe, independently of sGRBs and GWs.
As the opening angle of the escaped cocoon (θ es c ) [see equation (B10)] is ∼ 2 times larger than the jet opening angle (θ j ), simple estimation indicates that the expected number of these X-ray transients is about ∼ 4 times larger than the number of sGRBs in the local Universe (< 1 Gpc).

Detection as an EM counterpart to GW
Here, we roughly evaluate the detection rate of the cocoon emission as a counterpart of GW events, in optical or UV bands, and for BNS merger events observed with at least 2 GW detectors.We focus on the upcoming O5 run of LIGO-VIRGO-KAGRA.With O5 run, LIGO is expected to detect BNS mergers up to ∼200 Mpc (angleaveraged sensitivity; Abbott et al. 2020).We assume that, at such distance, early EM follow-up and localiza-tion will be successful (e.g., with LSST's wide field of view).The duration of O5 is set for 2 years.However, as each detector is expected to be active for ∼ 70 − 80% of the time, effectively the duration will be ∼ 1 year (Abbott et al. 2020;Davis et al. 2021).
We take into account contribution from the trapped cocoon, as its emission has the widest opening angle and can be as large as θ c ∼ 70 • (see §D) [considering that it is concealed by the ejecta on the sides making it dimmer].Hence, considering that we are using an angle-averaged sensitivity, the probability of detecting BNS with a viewing angle within θ c ∼ 70 • of the polar axis8 can be found as 1 − cos θ c ≈ 0.66.
The number of events to be expected in the O5 run can be found as the product of, BNS merger rate (∼ 320 +490 −240 Gpc −3 yr −1 ), observable volume (∼ 0.03 Mpc 3 ), the probability of sGRB-jet launch (taken as 1), the probability of a wide jet launch (∼ 0.2), the probability of a late jet launch (∼ 0.7), the probability of BNS merging so that an on-axis line-of-sight of the trapped cocoon, within θ c , is possible (∼ 0.66), and the time of O5 multiplied by the duty cycle (∼ 2 yr×0.5).This gives: This is equivalent to ∼ 1/20 of every BNS mergers within 200 Mpc being detected.This rate of detection is slightly lower than that found in equation ( 14) [due to the limitation in distance D] but overall similar remarks can be made.Considering recent analysis by Abbott et al. (2023) suggesting a BNS merger rate of 10 − 1700 Gpc −3 yr −1 , the above number of events (for 2 years) is found within the range ∼ 0.03 − 5.3 Gpc −3 .Therefore, although these estimations are based on preliminary parameters (i.e., of O5 that will most likely be subject to changes), we estimate that the odds are favorable for the cocoon emission to be detected and well monitored in O5 (as an EM counterpart to GWs) for the first time.Hence, O5 will play an important role at answering questions related to GW170817/AT2017gfo, sGRBs, and late engine activity.
For the LIGO-VIRGO-KAGRA O4B run, at the time of writing, VIRGO's sensitivity is planned to reach 80 Mpc for the last 3 quarters of 2024 (although this is constantly subject to change).If so, the expected detection rate in O4B would be much lower than that found for O5 [see equation ( 15)], by a factor of ∼ (80/200) 3 (0.75/2) ∼ 1/40; i.e., ∼ 0.024 ≪ 1 event in O4B.
It should be noted that combining these detection channels with targets of opportunity from GRB detectors, is expected to increase even further the number of events.Also, as the detection rate is ∝ D 3 , future 3rd generation of GW detectors (Cosmic Explorer Reitze et al. 2019;and Einstein Telescope Punturo et al. 2010a,b;Maggiore et al. 2020) expected in the coming decade and with future wide field detectors (e.g., LSST; Ivezić et al. 2019) will dramatically increase the detection rate by ∼ 3 orders of magnitude (e.g., Zhu et al. 2023).

SUMMARY & CONCLUSION
Considering the observational evidence of late extended/plateau emission phases after the prompt emission in most sGRBs (e.g., see Kisaka et al. 2017), and that these late emission phases are most likely due to late engine activity, we explored the scenario of propagation of late jet outflow in BNS merger ejecta (see Figure 1).Following observational analysis of sGRBs (Fong et al. 2015;Kisaka et al. 2017;Rouco Escorial et al. 2023), we estimated the typical parameters of prompt jets and late jets in consistency with the observed prompt, extended, and plateau emission phases (see §2.3; and Figure 2).For each of these three phases, we considered three types of jets according to their width: narrow (3 • ), typical (6 • ), and wide (15 • ) [see Table 1].
Based on the analytic model for prompt jet propagation in expanding medium Hamidani et al. (2020); Hamidani & Ioka (2021), we analytically solved late jets' propagation and breakout [see equation (7); §A].Following the same arguments as in Hamidani & Ioka (2023a), we estimated the properties of the jet-shock heated cocoon before and after the breakout (see §B).Then, in §3, we presented an analytic model to estimate the cocoon cooling emission.Our model was based on Hamidani & Ioka (2023b) with additional improvements as required for late jets (see §C).To our best knowledge, this is the first time that the cocoon of late jets have been considered, and its emission have been calculated.
Our results indicate that typical plateau jets can barely breakout (Table 1).This indicates that, although it is possible to explain plateau emission of ∼ 10 46 erg s −1 with jet luminosities ∼ 10 47 erg s −1 and radiation efficiency η γ ∼ 0.1, it would be difficult to explain fainter plateau emissions that have been observed (see Figure 2), unless the ejecta mass is much smaller than what has been observed in AT2017gfo (∼ 0.05 M ⊙ or prompt/extended jets are strong enough to permanently leave a cavity (see §A.2).Such low mass is plausible in systems where the disk/ejecta mass is small, which is favored by mergers of BNS of equal masses (see Kiuchi et al. 2019).Another possibility is plateau jets with η γ ≳ 0.1.
We explained that the cocoon emission initiates upon the jet breakout.Consequently, the delay between the merger and the onset of the cocoon emission can be easily measured, rendering it a useful quantity as a proxy for the jet breakout time [see equation ( 7)].
We found that for later jets, the cocoon emission is more luminous and lasts longer, compared to prompt jets (see top panels in Figure 4).We showed that late jet head velocity ( Hamidani et al. 2020) making late jets more efficient at depositing energy into the cocoon.Furthermore, because breakout times are longer in late jets, adiabatic cooling (∝ t b /t) is less efficient, and more internal energy is available for the cooling emission.In prompt jets, compared to late jets, shock heating takes place much earlier (∼ 1 s; compared ∼ 10 − 10 4 s), and adiabatic cooling consumes much of the internal energy before cooling emission start to diffuse out.Consequently, for certain late jets, the trapped cocoon (irrelevant in prompt jets) can be brighter than the expected r-process powered KN at late times (∼ 10 4 − 10 5 s) [see Figures 3  and 4 (top panels); also see equation ( 12)].
We also showed that the brightness of the cocoon emission is strongly dependent on the jet opening angle, in agreement with Hamidani & Ioka (2023b) [see their equations ( 59) and ( 60)].
The cocoon emission was found to peak in the UV during most of its emission, ∼ 10 2 − 10 5 s after the merger (see bottom panels in Figure 4 and top panels in Figure 5), making it an excellent target for the upcoming UL-TRAST mission.The trapped cocoon (of late jets) is found to peak with ∼ 10 4 − 10 3 K, at ∼ 10 3 − 10 5 s after the merger, making it bright in optical bands as well.This makes it easier to identify due to the less extinction in optical, and the large number of facilities/surveys in optical bands with more depth in their observations (e.g., LSST).
At late times (∼ 10 4 − 10 5 s), we found that late jets (plateau in particular) with their cocoons (e.g., model PL-T) can produce a blue KN-like emission, similar to the one in AT2017gfo (see Figures 4 and 5).For wider jets (e.g., model PL-W), an extremely bright blue KN can be produced (see Figures 4 and 5).This might be one possible scenario for the exceptionally bright blue KNe identified in some sGRBs (Rossi et al. 2020; also see Tanvir et al. 2013;Jin et al. 2015;Gao et al. 2017;Troja et al. 2018;Gompertz et al. 2018;Jin et al. 2020 and the references therein).In this scenario the brightness of the red KN is determined by the ejecta mass, while the brightness of the blue KN is determined by energy deposition into the ejecta by the late jet (Hamidani et al. in preparation).
Previously, it has been suggested that AT2017gfo's early data can be explained with the prompt jet's cocoon emission (Piro & Kollmeier 2018).However, we have clarified that the prompt jet's cocoon, given its parameters, is not able to explain AT2017gfo's data (see grey lines Figures 4 and 5; also see arguments in Hamidani & Ioka 2023a,b).Instead, we showed that late jets (plateau in particular) with their cocoon present a viable alternative to explain the early AT2017gfo data and the blue KN (see blue lines Figures 4 and 5).
We also showed that the cocoon is bright in soft Xray 10 − 100 s after its breakout (see Figure 6).This confirmed previous claims in Hamidani & Ioka (2023b).We explained that such soft X-ray emission is within reach of the upcoming Einstein Probe mission.
We estimated the rate of detection of the cocoon emission by Einstein Probe, or as an EM counterpart to GWs in the O5 run of LIGO-VIRGO-KAGRA, and found that it can be detected on a yearly basis [see equations ( 14) and ( 15)].This detection rate is expected to increase by ∼ 2 − 3 orders of magnitude with future 3rd generation of GW detectors (Cosmic Explorer Reitze et al. 2019;and Einstein Telescope Punturo et al. 2010a,b;Maggiore et al. 2020) expected in the coming decade and with future wide field detectors such as LSST (Ivezić et al. 2019).
Traditionally BNS mergers have been tracked indirectly through the prompt emission of sGRBs.However, sGRBs jets are beamed, and considering that their typical opening angle is θ j ∼ 6 • (Rouco Escorial et al. 2023), only ∼ 1/200 sGRBs (or BNS mergers) is on-axis and can theoretically be detected.This is without considering that an unknown fraction sGRB-jets are choked and fail to produce a sGRB.Also, for instruments like Swift/BAT, this number is reduced by about one order of magnitude due to the limited field of view.This has resulted in a bias toward intrinsically bright, on-axis, and distant sGRBs, rendering them challenging to localize, follow up (e.g., in optical), and exploit.
GWs is another tool to detect BNS mergers.GWs offer great advantages, as they are less sensitive to the viewing angle, and give the possibility of detecting BNS even if the sGRB jet is off-axis.This feature favors the detection of off-axis nearby events.Closer distances for GW detected BNS events can make push these otherwise faint transients within the range of ground-based telescope facilities.GWs also can detect BNS source a few minutes before the actual merger.However, there are two disadvantages associated with GW detection of BNS: i) the spatial localization is poor (although Fermi GBM's localization can be poorer) requiring more GW detectors, and ii) the current limited detectors' horizon making BNS detections with GWs quite rare.
Here, we showed that the cocoon of sGRB-jets can be another channel to detect BNS mergers, in particular in soft X-ray (also in UV), with wide field and high cadence observations.The upcoming EP and ULTRASAT missions are ideal for detecting this emission (also the upcoming UVEX mission; Kulkarni et al. 2021).With EP, we estimate that ∼ 1/600 BNS mergers within ∼ 1 Gpc could be detected [see equation ( 14)].
The cocoon emission has its advantages for identifying BNS mergers, compared to sGRBs, as it is not biased to the jet; it may even accompany BNS mergers that are off-axis or with failed jets (also low-Γ jets; Lamb & Kobayashi 2016); and is good at detecting more local BNS mergers.The cocoon emission (e.g., with EP) has also its advantages compared to current GW detectors in terms of the horizon of the detection, exposure time, sky localization, and number of events [see equations ( 14) and ( 15)].Although off-axis afterglow emission (nonthermal) from structured sGRB-jet could contaminate this cooling emission (thermal), difference in timescales is expected.
Hence, the cocoon and its emission can be a useful tool.In the era of multi-messenger astronomy, understanding and tracking this cocoon emission, in coordination with sGRB and GW observations, could help make breakthroughs regarding the central engine of sGRBs, their jets, and the origin of heavy elements.

DATA AVAILABILITY
The data underlying this article will be shared on reasonable request to the corresponding author.

B. ANALYTIC MODEL FOR COCOON BREAKOUT AND ESCAPE
We follow the analytic model presented in Hamidani & Ioka (2023a) with the same approximations.The model is based on the parameter α, which is the ratio of the energy density of the cocoon (jet energy + kinetic energy of the expanding ejecta) over the energy density of the expanding ejecta (kinetic energy of the ejecta).α can be written as a function of the initial condition, breakout time t b , and cocoon volume (over the ejecta volume: V c /V e ), as: This parameter is crucial to determine the part of the cocoon that can escape the ejecta.Then, the cocoon part that is able to escape has the fractions of the mass, energy (kinetic+internal), and internal energy, which can be written as a function of α, as bellow [similarly to equations ( 44), ( 52), and (53) in Hamidani & Ioka (2023a)]: where f mix ∼ 2/3 (or ∼ 1 for failed jet case) is a mixing parameter for the shocked jet and ejecta.The cocoon mass, energy, and internal energy, can be found at the breakout time as Ve , and E c,i ≈ 3P c V c , respectively.It should be noted that the fraction of the escaped cocoon, in particularly in terms of mass, is very small (∼ 0.5−5%).After the breakout, and until t 1 ∼ 2 − 3t b , these quantities (all the cocoon's and escaped cocoon's) are expected to grow by a factor of f g ∼ 2 (due to the post-breakout cocoon growth).
Additionally as shown in Hamidani & Ioka (2023a), the escaped cocoon takes the shape of a cone with θ es c as its opening angle.A reasonable approximation of this angle is found (empirically) as determined by the solid angle of the cocoon ellipsoidal at its midpoint (at r b /2, at the breakout time; while accounting for the lateral expansion at t > t b ): The term f g accounts for the growth of the cocoon after the breakout (in terms of mass and energy; see Hamidani & Ioka 2023a).This is explained by the cocoon's volume growth (by a factor f g ), and more specifically by the lateral expansion of the cocoon (in r c , in the comoving frame of the ejecta) by a factor of f g (as

C. ANALYTIC MODEL FOR THE ESCAPED COCOON'S COOLING EMISSION
We follow the analytic model in Hamidani & Ioka (2023b), originally for the prompt jet case, with the exception that key approximations, t ≫ t b and r ≈ cβt (practical for the prompt jet case), are not ideal for late jets (as t b is much larger), and have to be bypassed.
For that, the formulations is generalized.In Figure 7 we show the laboratory time, radius, and observed time, with the merger as the reference (t, r, and t obs ; as used in Hamidani & Ioka 2023b).The alternative coordinates are also shown (in red), with the jet breakout time (or radius) as the new reference, as follows: For the cocoon breakout, the following picture is considered.As the escaped cocoon is originally located in the ellipsoidal cap of the pre-breakout cocoon, and as this cap is usually a small fraction of the cocoon (Hamidani & Ioka 2023a), we approximately consider the escaped cocoon's breakout take places in the polar axis of the ejecta (at the point r = r b and θ ∼ 0).Considering that the spread of the escaped cocoon can be approximated to a cone structure with the opening angle θ es c [see equation (B10)], the base of this cone has been set at the same point (r = r b and θ = 0 in the central engine frame).Note that, until the free-expansion phase is reached (t ≳ t 1 , introduced next), our analytic model cannot be applied (due to the singularity in around t ∼ t b ).
For the rest, the same procedure in Hamidani & Ioka (2023b) is followed.The escaped cocoon is assumed to accelerate outward after the breakout time, and reach the free expansion at t = t 1 ∼ 2t b (t = t 1 ∼ 1.5t b for the plateau jets where t b is long enough).The volume of the escaped cocoon up to a certain velocity β can be found as: where Ω = 4π(1 − cos θ es c ).The escaped cocoon is divided into two parts: relativistic (β out > β > β t ) where Γβ ∼ Γ, and nonrelativistic (β t > β > β m ) where Γβ ∼ β (with Γ out ≈ 10, β t ≈ 0.8, and β m ≈ 0.35), each with a different power-law distribution for the mass-density (indices are l ≈ 0 and m ≈ 8, respectively) and internal energy density (indices are −3 and 2, respectively; see Hamidani & Ioka 2023b).
Hence, the mass-densities can be found as: . These expression are equivalent to those in Hamidani & Ioka (2023b) with t substituted by ∆t.
Similarly, assuming grey opacity (κ = Const.)gives the observed time of a given shell Γ d β d as a function of its optical depth, as: τ r and τ nr are the optical depth.We use the same criteria as in Hamidani & Ioka (2023b) for photon diffusion from a thin shell: τ r = 1; and With the assumption that the free-expansion is reached at ∼ t 1 ∼ 2t b , adiabatic cooling starts after t 1 9 .One can find [using equation (C14)] that V (> β m , t > t 1 ) ∝ ∆t 3 and the coefficient for internal energy loss due to adiabatic cooling is V (β>βm,∆t1)   V (β>βm,∆t) Hence, the evolution of the internal energy can be found as: Then, the internal energy beyond a certain cocoon shell (moving with Γβ) and at a given time (∆t) can be found (for the relativistic and non-relativistic cocoon parts, respectively) as: where E j c,i,r (∆t) and E j c,i,nr (∆t) are the total internal energies (from jet heating) at a given time (∆t), in the relativistic, and non-relativistic parts of the cocoon, respectively.Here too, these expressions are equivalent to those in Hamidani & Ioka (2023b) (with t substituted by ∆t).
Next, using the time evolution of the diffusion velocity, the internal energy density [found using equations (C14) and (C18)] is used to estimate the radiated energy (in the form of blackbody photons); i.e., cooling luminosity.This is determined differentially (avoiding any one-zone approximations due to their inaccuracy), while taking into account the effect of the observer time and the relativistic beaming effect.The cooling bolometric luminosity can be found using: is the solid angle of the emission coming from the shell moving 9 This is consistent with numerical simulations in Hamidani & Ioka (2023a).Also, this is a similar effect to the geometrical effect discussed in Nakar & Sari (2010).One likely explanation is that, as the escaped cocoon fluid accelerates during its escape, and before it converges to a homologous structure at t 1 (free-expansion), internal shocks happen, between the unorganized escaping cocoon shells, enhancing the internal energy content.Therefore, the full effect of the adiabatic cooling is not seen until the free expansion is reached, as these interactions can no longer happen.As shown in Hamidani & Ioka (2023a), simulations show that internal energy is boosted by fg ∼ 2, and since t 1 ∼ 2t b , one can approximate the adiabatic expansion to start after t 1 .This approximation is not expected to affect the cocoon cooling emission, as photons are mostly released at t ≫ t 1 .with Γ d .The observed time can also be calculated as a function of the laboratory time as: ∆t obs = (1 − β)∆t. (C20) With equations (C16) and (C18), the full analytic expression of the jet-shock powered bolometric luminosity can be found as a function of ∆t obs (by varying the diffusion velocity Γ d β d ) as:  (C21) The photospheric radius r ph (and velocity Γ ph β ph ) and the observed temperature can be found using equations ( 10) and (11), after applying the criteria τ ph = 1 for the optical depth.All these expressions are similar to those found in Hamidani & Ioka (2023b), with t, r ph , and t obs , being substituted with ∆t, ∆r ph , and ∆t obs , respectively.The observed blackbody flux density can be found (for a given frequency ν) as

D. ANALYTIC MODEL FOR THE TRAPPED COCOON'S COOLING EMISSION
In addition to cooling emission from the escaped cocoon, here we present an analytic estimation of the cooling emission from the jet-shock-heated trapped cocoon (i.e., β < β m ).The same method is used to calculate the cooling emission from the escaped cocoon as in the following.Notable differences are the following: • The morphology of the trapped cocoon's shape, although better represented by that of an ellipsoidal at the breakout, is approximated to take the form Figure 7. Illustration of the relation between the laboratory time (since the merger time; with tm = 0) t lab ≡ t, the radius r (from the central engine's center), and the observed time t obs (for an observer at R obs , and since the observed merger time; t m,obs = 0) [black solid lines].Colored solid lines show the homologous expansion of the outer edge of the ejecta (magenta) and a given shell of the escaped cocoon with a velocity β (purple).Dashed lines show radiation (GWs in grey, breakout emission in red, and cooling emission from the shell moving with β in blue).The event of jet breakout is shown in red, in each of the three coordinates.The alternative coordinates used here with the event of jet breakout as their starting point (∆t, ∆r, and ∆t obs ) are shown with red solid lines.The choice of colors here is arbitrary with no connection to observed colors.
of a cone with a solid angle Ω c at much later times t ≫ t b , so that: The opening angle of the trapped cocoon is larger than the opening angle of the escaped cocoon (by about a factor ∼ 2).
• As the escaped cocoon mass is negligible compared to the overall mass of the cocoon, we approximate that mass-density of the trapped cocoon is comparable to the mass-density of the ejecta (see Figure 2 in Hamidani & Ioka 2023b).Hence, the trapped cocoon's mass-density profile is fitted with the same power-law index n, as the ejecta.
• When estimating the optical depth of the trapped cocoon, optical thickness of the outer escaped cocoon is neglected.
• Since the velocity of the trapped cocoon (and the ejecta) is significantly smaller, the non-relativistic limit β ≪ 1 is applied, allowing us to express the observed time and the laboratory time interchangeably; t obs ∼ t.
• Cooling emission from the trapped cocoon is approximated to start once the bolometric luminosity of the trapped cocoon dominates (over the escaped cocoon's).
The mass-density of the trapped cocoon can be found as Calculating the optical depth as τ ≈ βout β κρ(ct)dβ gives the following expression of time, at which the cooling emission from the β d shell is observed: The criteria for photon diffusion from the β d shell is (see Kisaka et al. 2015;Hamidani & Ioka 2023b): For the internal energy density of the trapped cocoon, taking into account that the pressure inside the trapped cocoon should roughly be constant, and after accounting for adiabatic internal energy loss and volume increase (in a homologously expanding system), one can find that P c (t) ∝ V with t 1 ∼ 2t b is the time from which the cocoon is considered to be expanding freely, and E tr c,i (t 1 , > β 0 ) = E c,i − E es c,i is the total internal energy in the trapped cocoon (accounting for the escaped internal energy).
In the non-relativistic limit, the cooling emission can be found as:

E. ANALYTIC MODEL FOR THE R-PROCESS POWERED KILONOVA
We present a simple analytic model for the r-process powered KN emission from the merger ejecta.The model shares several similarities to the one in Kisaka et al. (2015).The aim of this model is to estimate the early KN emission and give a relative comparison to the brightness of the cocoon emission.
Here too, the approximation of a diffusion shell (β d ), is used (see §C).The same mass-density and velocity profile is considered for the KN ejecta, although here, focusing on the early part of the KN (i.e., blue KN), peaking at ∼ 1 day (in terms of magnitudes), only the blue component M KN ∼ 0.02 M ⊙ is considered.
First, similarly to §D, the optical depth is calculated to find the diffusion velocity β d as a function of time (assuming the same grey opacity κ = 1 cm 2 g −1 ).The mass of the KN ejecta beyond the velocity β d can be found as: At later times (t ≫ t h ), the energy deposition (from r-process) in the optically thin KN ejecta is where k ∼ 1.3, and ε0 ∼ 1 × 10 18 erg g −1 s −1 is the heating rate at t h ∼ 0.1 s (or ∼ 2 × 10 10 erg g −1 s −1 ; see Wanajo et al. 2014 andIshizaki et al. 2021).
The bolometric luminosity of the KN is found as a function of time as the sum of two terms: where L KN (⩾ β d , t) is emission from the optically thin part of the KN ejecta, and L KN (< β d , t) is the cooling emission from the optically thick region due to the inward motion (in Lagrangian coordinate) of the diffusion shell.
Taking into account that only a fraction of this deposited energy is thermelized, the first term can be found as: where we consider the thermalisation efficiency around 1 day (around the peak of the blue KN) f tot ∼ 0.4 (Hotokezaka et al. 2016;Barnes et al. 2016).
For the second term, the released internal energy can be found using L KN (< β, t) = f . (E36) Note that here too, as β ≪ 1, the observed time and the laboratory time can be used interchangeably t obs ∼ t.

Figure 3 .
Figure 3. Cooling emission from the jet-shock-heated cocoon, for the model EX-W.Three cocoon components are shown: relativistic cocoon (dark purple), non-relativistic cocoon (light purple), and trapped cocoon (light blue) [using equations (C21) and (D29)].The vertical (black) line indicates the breakout time in the observer's frame [t b,obs = (1 − βm)t b ; see equation (C20)].Analytic estimation of the r-process powered KN emission is also shown (dotted grey line) [calculated using equation (E33)].The starting time of light curves is set as the start of the free-expansion phase t1 ∼ 2t b (see §C).The transition of emission from the nonrelativistic cocoon to the trapped cocoon is set to take place once the trapped cocoon's luminosity dominates over the non-relativistic cocoon's luminosity (here at ∼ 10 3 s; see §4.2).

Figure 4 .
Figure 4. Bolometric isotropic luminosity (top), photospheric four-velocity (middle), and observed temperature (bottom), for different jet models [prompt (grey), extended (green; left panels), and plateau (blue; right panels)] with different opening angle [narrow (thin), typical (medium), and wide (thick)], as a function of the observed time since the merger.Our analytic model of the r-process powered early KN (grey dotted line) is shown (see §E) for comparison, as well as the observations of the KN in AT2017gfo (grey circles;Villar et al. 2017;Waxman et al. 2018).The starting time of the cocoon emission is set as the start of the free-expansion phase t1 (see §C).For narrow and typical prompt jet models, T obs evolution at later times is nearly identical.

Figure 5 .
Figure 5. Apparent (left axis; for a distance D = 40 Mpc) and absolute (right axis) magnitudes of the cocoon emission from different jet phases [prompt (light grey), extended (green), and plateau (blue)] and different jet models [narrow (thin lines), typical (medium), and wide (thick)] as in Table 1].Top and bottom panels are for UVOT's UVM2 band, and g band; respectively.The expected magnitude of the r-process powered KN is shown (dotted grey line) using our analytic model in §E.Observation of the KN AT2017gfo are also shown (grey circles for detection, and triangles for upper limits; taken from Kasliwal et al. 2017; Drout et al. 2017; Villar et al. 2017 and the references therein).Limiting magnitude of the relevant detectors is indicated [horizontal black lines; UVOT: ∼ 22 for an exposure time of ∼ 1000 s (Roming et al. 2005); ULTRASAT: ∼ 22.4 for 900 s (Sagiv et al. 2014); ZTF: 20.8 for 30 s (Burdge et al. 2020); LSST: 24.9 for 30 s (Ivezić et al. 2019)].All models are shown in the frame of the source (no reddening applied) and AT2017gfo observations have been de-reddened [using Milky Way extinction model; Rv = 3.1 and E(B − V ) = 0.11] using(Brown et al. 2010; for UVM2 band) and(Schlegel et al. 1998; for g band).A jet with parameters similar to those of the typical plateau jet model (PL-T) appears to reproduce the blue kilonova (KN) of AT2017gfo effectively.There is a slight difference between the analytic models in the UV band and the g band, e.g., in comparison to AT2017gfo.This discrepancy, although the emission process is a blackbody in both cases, arises from our use of a simplified, non-wavelength-dependent mean grey opacity (κ = 1 cm 2 g −1 ; seeBanerjee et al. 2023).
achieved thanks to the following: Cray XC50 of the Center for Computational Astrophysics at the National Astronomical Observatory of Japan, and Cray XC40 at the Yukawa Institute Computer Facility.
C11) ∆r = r − r b ≈ cβ∆t, (C12) ∆t obs = t obs − t b,obs , (C13) with r b as the breakout radius (in the central engine frame), t b,obs = t b (1 − β m ) as the observed breakout time since the merger time (t m,obs = 0).