Dust evolution in a supernova interacting with the ISM

Supernovae (SN) explosions are thought to be an important source of dust in galaxies. At the same time strong shocks from SNe are known as an efficient mechanism of dust destruction via thermal and kinetic sputtering. A critically important question of how these two hypotheses of SNe activity control the dust budget in galaxies is still not quite clearly understood. In this paper we address this question within 3D multi-fluid hydrodynamical simulations, treating separately the SNe injected dust and the dust pre-existed in ambient interstellar gas. We focus primarily on how the injected and the pre-existing dust is destroyed by shock waves and hot gas in the SN bubble depending on the density of ambient gas. Within our model we estimate an upper limit of the SN-produced dust mass which can be supplied into interstellar medium. For a SN progenitor mass of 30 $M_\odot$ and the ejected dust mass $M_d=1~M_\odot$ we constrain the dust mass that can be delivered into the ISM as $\geq 0.13~M_\odot$, provided that the SN has injected large dust particles with $a\geq 0.1~\mu$m.


INTRODUCTION
Dust is an important constituent of the interstellar medium (ISM), playing key role in physical processes that determine its basic properties: chemical transformations, metal budget in gas phase, thermodynamic state, radiation transfer.Dust is known to convert stellar light and thermal energy of hot gas in infrared (IR) radiation.One of the intriguing questions concerns the dust mass budget in the ISM (see for recent discussion in Mattsson 2021;Kirchschlager et al. 2022;Péroux et al. 2023).The interrelation between the destruction and production dust rate in the ISM is still not quite well understood.Of particular concern is the processing the dust particles undergo behind the strong shock waves, > ∼ 150 km s −1 , penetrating the ISM.It is generally thought that dust particles experience efficient destruction from such shocks and from the hot gas behind them.The three processes dominate the dust destruction: the inertial and thermal sputtering (Barlow 1978;Draine & Salpeter 1979b,a;McKee 1989;Jones et al. 1994;Nath et al. 2008;Slavin et al. 2015;Priestley et al. 2022, and references therein), and the shattering in grain-grain collisions at higher densities (Borkowski & Dwek 1995;Jones et al. 1996;Slavin et al. 2004;Guillet et al. 2009;Bocchio et al. 2016;Kirchschlager et al. 2019).
Theoretical considerations show that the characteristic dust lifetime in the Milky Way ISM against sputtering is estimated to be < ∼ 3 × 10 8 yr (McKee 1989;Jones et al. 1994) to < ∼ 3 × 10 9 yr (Jones & Tielens 1994;Slavin et al. 2015), resulting in the decrease rate − < ∼ (0.1 − 0.01) ⊙ yr −1 ; more recent discussion can be found in (Bocchio et al. 2014;Ginolfi et al. 2018;Micelotta et al. 2018;Ferrara & Peroux 2021).On the other hand, the overall production rate from red giants winds and SNe explosions is + ∼ 10 −3 ⊙ yr −1 , where SFR ∼ 5 ⊙ yr −1 is assumed, indicating a severe disbalance between the dust destruction and its replenishment (Draine 2009;Bocchio et al. 2016).This discrepancy can be miti-gated when dust-to-gas decoupling under the action of shock waves and gravity for large particles is accounted for (Hopkins & Lee 2016;Mattsson et al. 2019a;Mattsson & Hedvall 2022).However, one to two orders of magnitude difference between the − and the + requires apparently a proportional amount of dust to be hidden of destructive SNe shocks, which seems unrealistic.An additional and apparently efficient dust mass supply can be connected with growth of grains in the ISM (Draine 1990;Chokshi et al. 1993;Dwek 1998), and more recently (Calura et al. 2008;Draine 2009;Mattsson 2011;Inoue 2011;Ginolfi et al. 2018;Heck et al. 2020, and references therin).Moreover, supersonic turbulence is shown to be an efficient mechanism that can stimulate formation of dust in the ISM of local galaxies and even in the early Universe, and as such can counteract efficient dust destruction by SNe (Hopkins & Lee 2016;Mattsson et al. 2019a,b;Mattsson 2020a,b;Li & Mattsson 2020;Commerçon et al. 2023).
The discrepancy between the production and destruction rate is seen in particular in galaxies at high redshifts as first pointed out by Todini & Ferrara (2001); Morgan & Edmunds (2003).The detection of dust in quasars at ∼ 5 (the universe's age < 1 Gyr) suggests the SNe II to be the dominant dust source in the early universe (Bertoldi et al. 2003;Maiolino et al. 2004;Beelen et al. 2006;Valiante et al. 2009Valiante et al. , 2011;;Dwek & Cherchneff 2011;Gall et al. 2011;Riechers et al. 2013).Moreover, further observations at intermediate redshifts ∼ 1 − 5 with the Herschel Space Observatory revealed a more generic problem -an apparent excess of dust in submillimeter and ultraluminous IR galaxies -the so-called 'dust budget crisis' (Michałowski et al. 2010;Dunne et al. 2011;Rowlands et al. 2014).
Direct IR observations in the local Universe also show that dust can be produced at early stages of the ejecta outflows as manifested in several nearby SN remnants (in particular, in Cas A and Crab, Dunne et al. 2003;Gomez et al. 2012;Arendt et al. 2014), in SN1987A (Indebetouw et al. 2014;Matsuura et al. 2015), and in the local group galaxy NGC 628 (Sugerman et al. 2006).More recent analysis of IR characteristics from SN1987A (Wesson & Bevan 2021) indicates that dust forms in SNe ejecta at latter stages > ∼ 1000 days.The latter may reflect the fact that the net dust product by SNe is environmentaly sensitive as demonstrated by Nozawa et al. (2006).
Recent measurements of abundance patterns in supernovae remnants (SNR) in the Milky Way (MW) and Large Magellanic Clouds (LMC) with strong non-radiative shocks, ≃ 350 − 2700 km s −1 , have also constrained the destroyed dust fraction by < ∼ 0.1 − 0.6, even in a rather dense environment (see Tables 1 and 3 in Zhu et al. 2019).It is therefore conceivable that older SNRs with weaker shock waves are less destructive than commonly thought.More recent observations of three SNRs in the MW do confirm such a conclusion (Priestley et al. 2021).
The total dust mass supplied by a SN in the ISM, and the dust size distribution sensitively depend not only on the progenitor mass, but on the density of ambient gas as well, because of the reverse shock from the interaction of the ejecta with the ambient gas.1D simulations by Nozawa et al. (2006Nozawa et al. ( , 2007)); Bianchi & Schneider (2007); Nath et al. (2008) and Bocchio et al. (2016, for more recent discussion) have shown that increase of the ambient gas density from 0.1 to 10 cm −3 can result in a drop of dust yield by one to two orders, particularly for higher progenitor masses.Equally important is the conclusion that destruction of small size dust by the reverse shock considerably flattens the size spectrum.More recently, theoretical analysis of dust destruction in the process of interaction of the ejecta with ambient ISM leads Slavin et al. (2020) to conclude that when relative motion of dust particles with respect to the gas is accounted, the heavier dust particles can penetrate the region affected by the reverse shock and can escape relatively intact into the ambient ISM gas.However, eventually a fraction of them is becoming destroyed in situ in the surrounding ambient gas.This result suggests that SNe inject mainly large dust particles into the ISM.This in turn can contribute to variations of the extinction law on smaller scales, unless dust particles are well mixed with dust in ambient gas.
In this paper we present results of 3D multi-fluid hydrodynamical simulations of how SN ejecta with dust particles having formed in the initial stages of its evolution propagates through the dusty interstellar gas.We concentrate on the question of how dust particles injected by the SN and those present in the surrounding ISM are destroyed by strong shock waves.The injected dust particles and the interstellar dust are treated as two different particle populations experiencing distinct evolutionary paths.Formation and growth of dust particles in dense cool regions of the shell surrounding the remnant are not included in our consideration.Section 2 sets up our model.Section 3 describes the obtained results, including dynamical aspects related to dust destruction and overall mass budget of the injected dust and the dust belonging to ambient ISM (refereed as the interstellar or pre-existing dust).Section 4 contains a general discussion and the summary.

MODEL DESCRIPTION
We consider the dynamics and destruction of dust particles in a supernova remnant during its 300 kyr-long interaction with a slightly inhomogeneous clumpy ISM.We use our gasdymanic code (Vasiliev et al. 2015(Vasiliev et al. , 2017) ) based on the unsplit total variation diminishing (TVD) approach that provides high-resolution capturing of shocks and prevents unphysical oscillations, and the Monotonic Upstream-Centered Scheme for Conservation Laws (MUSCL)-Hancock scheme with the Haarten-Lax-van Leer-Contact (HLLC) method (see e.g.Toro 2009) as approximate Riemann solver.This code has successfully passed the whole set of tests proposed in Klingenberg et al. (2007).In order to follow the dynamics of dust particles we have implemented the method similar to that proposed by Youdin & Johansen (2007), Mignone et al. (2019) and Moseley et al. (2023).A description and tests are given in Appendix.In this paper we only take destruction of dust particles by both thermal (in a hot gas) and kinetic (due to a relative motion between gas and grains) sputtering (Draine & Salpeter 1979a) into account.In order to separate different mechanisms that can contribute to dust processing, we do not consider here additional effects from possible growth of dust particles in denser regions of the remnant, and their fragmentation from shattering collisions.These processes operate on time scales of a few Myr (see, e.g., Hirashita & Kuo 2011;Mattsson 2020a), even in much denser environments (molecular clouds) > ∼ 10 2 cm −3 (Martínez- González et al. 2022).This is much longer than the ages of SNe remnants in a diffuse ISM with the shocks sufficiently strong ( > ∼ 150 km s −1 ) for dust destruction: < ∼ 100 kyr.The backward reaction of dust on to gas due to momentum transfer from dust particles is also accounted in order to ensure dynamical self-consistency.Generally dynamical effects from dust are minor, however in the post-shock domain, such as the postshock shell, the collisional coupling is weak, and dust particles move inertially with velocities ≥ 3 /4, is the shock velocity.As a result, their contribution to ram pressure 2 cannot be negligible until the dust velocity relaxes to the gas one |v − v | ≪ , being the sound speed (for further discussion see Sec. 3.4).In addition, the effects of backward reaction can be of importance in correct estimates of the kinetic sputtering for particles inertially entering the ambient gas at the interface between the bubble interior and the surrounding gas.
The effects of magnetic field are apparently of critical importance for dynamics of dust particles in SN remnants, particularly in terms of their destruction.It is known, that in presence of betatron acceleration of charged dust can enhance the efficacy of sputtering at radiative stages on an expanding SN remnant, however, only when dust particles are decoupled of gas (Shull 1977;Draine & Salpeter 1979a;Seab & Shull 1983;Slavin et al. 2004).As we show below (Sec.3.4), even massive particles, ≥ 0.1 m, lose collisional coupling only behind the shock front in outer regions of the hot bubble where radiation cooling is insignificant, whereas smaller ones remain collisionally linked to gas.In addition, as shown by Slavin et al. (2004), individual trajectories of dust grains are very sensitive to grains' sizes and charges, their chemical composition and shock velocity with gyration radii of 0.03 to 0.3 pc.Adequate numerical description on scales of SNe remnants requires an unprecedented resolution, and is worth separate consideration.
The gas distribution in the ambient (background) medium is set to be a slightly inhomogeneous with averaged number density and small-size low-amplitude density fluctuations ∼ 2 − 3 pc, / ∼ 0.1.The fluctuations are constructed by using the module pyFC 1 by Lewis & Austin (2002) which generates lognormal "fractal cubes" for gas density field.The averaged ambient density is equal to = 1 cm −3 as a fiducial value, we calculate also models with 0.1, 0.3, 3 and 10 cm −3 , the temperature in all models is 10 4 K. Initially the gas pressure is assumed uniform over the whole computational domain.The metallicity of the background gas is usually set to be solar [Z/H] = 0. We assume a dust-to-gas (DtG) mass ratio equal to = 10 −2− [Z/H] .We inject the mass and energy of a SN in a region of radius 0 = 3 pc, assuming commonly used values of mass and energy30 ⊙ and 10 51 erg.The energy is injected in thermal form.The injected masses of metals and dust are 10 ⊙ and 1 ⊙ , correspondingly.Dust particles are thought to start growing during the free expansion phase of a SN evolution on much earlier times, when the ejecta temperature drops to several thousand Kelvins (e.g.Dwek & Scalo 1980;Todini & Ferrara 2001), typically after > ∼ 300 days when the ejecta radius is of ∼ 0.03 pc.At this age the characteristic time of dust-gas collisional coupling (the stopping time, see Sec.A2) for the parameters within the ejecta is ∼ 10 6 −1 0.1 s, with 0.1 being the dust grain radius in 0.1 m.For 0 = 30 ⊙ , and 0 = 3 pc, the ejecta gas density ∼ 10 3 cm −3 , resulting in ∼ 10 3 0.1 s at gas temperatures typical for dust growth.The dust grains nucleate and grow in situ in the ejecta as first described by Dwek & Scalo (1980), and have developed further on in Todini & Ferrara (2001).This would suggest that the formed dust particles get mixed into the ejecta on a short time scale.At such conditions in the process of acceleration the relative velocity of dust grains and gas is small, and only a small fraction of dust particles experiences kinetic sputtering (Slavin et al. 2020, see also estimates in Sec.A4).The thermal sputtering remains apparently weak because of a low temperature in ejecta < 6000 K.However, a firm estimate of the surviving dust fraction requires more detailed consideration.In this sense the estimates of the contribution of dust supply into the ISM from SNe explosions presented in our model can be regarded as upper limits.When the ejecta meets the reverse shock and becomes thermalized at ∼ 300 yr, thermal sputtering increases and comes into play (Nozawa et al. 2007;Nath et al. 2008;Slavin et al. 2020).Overall, at ≥ 300 − 400 yr the destroyed dust mass fraction can vary from 0.3 to 0.8 depending on dynamics in the surrounding gas (Slavin et al. 2020).
The mass of dust is redistributed between many dust 'superparticles'.For instance, in case the number of dust 'superparticle' is 2 20 ∼ 10 6 , the mass of each dust 'superparticle' is ∼ 1 ⊙ /10 6 ∼ 10 −6 ⊙ for the monodisperse dust of total mass equal to = 1 ⊙ .In the equation of motion of dust we attribute to a single grain the mass and the size of physical dust particles, while in the momentum equation for the gas component the momentum from dust is treated as a sum of momentum from all dust grains that constitute the dust mass ∼ 10 −6 ⊙ allocated in a cell.
In this paper we assume dust grains to have equal sizes (monodisperse dust model) at = 0. We perform a set of independent runs for monodisperse dust populations with initial sizes equal to 0 = 0.03, 0.05, 0.1, 0.3, 0.5 m; 0 = 0.1 m is the fiducial value.During the evolution dust grains are destroyed and their sizes decrease depending on physical conditions in ambient gas; the minimum size is set to 0.01 m.We adopt the same initial distribution of sizes both for the dust pre-existed in ambient gas (referred hereinafter as the interstellar dust) and for the injected dust as well.We assume that the mass of metals produced due to sputtering is returned to the gaseous phase.Evolution of polydisperse dust will be described elsewhere.
In our models the spatial resolution is set to 0.5 pc, that is sufficient for adequate treatment of dynamics of a SN bubble in a medium with density ∼ 0.1 − 10 cm −3 .In particular, in a medium with = 1 cm −3 we follow the SN bubble evolution until late radiative phase: ∼ 10 ∼ 300 kyr, when its radius reaches around ≃ 40pc.Therefore, we set the size of the computational domain (96 pc) 3 , with the number of cells 192 3 .We set one dust 'superparticle' per a computational cell resulting in the total number of interstellar dust particles in the domain 192 3 ∼ 7 millions.The fiducial number of injected 'superparticles' per SN is 2 23 ∼ 8 millions.Our choice of both grid resolution and number of the injected particles is justified by our tests presented in Appendix B. Note that dust particles can escape their initial "mother" cells and move into neighbouring ones.As a result, depending on gas density the overall distribution of dust particles can be rather patchy and can lead to local deficiency of dust particles, also seen in averaged radial profiles, see in Sec.3.2.
Simulations are run with tabulated non-equilibrium cooling rates fitting the calculated ones for a gas that cools isochorically from 10 8 down to 10 K (Vasiliev 2011(Vasiliev , 2013)).The heating rate is assumed to be constant, with a value chosen such as to stabilize the radiative cooling of the ambient gas at = 10 4 K.

Dust radial distribution
Inside a SN bubble dust is efficiently destroyed by both the thermal and kinetic sputtering (Draine & Salpeter 1979a).The effect of the latter on interstellar dust particles in a unmagnetised case is around 15%, while for the injected ones it is even smaller because of a lower drift between the gas and dust 2 .Fig. 1 presents angle averaged radial profiles for gas density and temperature, and densities of both sorts of dust particles -the SN injected and the interstellar dust, for the fiducial grain size at 50 and 300 kyr after the SN has exploded.Because of initially introduced clumpy density distribution and partly because of growing Rayleigh-Taylor (RT) instabilities the averaged shock layers in ( ) and ( ) are wider than typical.At = 50 kyr the injected dust remains ∼ 7 pc deeper behind the shock and about 10% survive within the central ∼ 10 pc interior.The interstellar dust is swept by the shock and penetrates behind it by ∼ 5 pc inwards, but only in a ∼ 2 − 3 pc thick layer it survives within factor of 2, as seen from comparison of the magenta (with sputtering accounted) and the light-magenta (no sputtering is included) lines in upper paned of Fig. 1.At 300 kyr only a small fraction of injected dust survives within the bubble interior, around 10% if it penetrates ahead of the forward shock as shown by the line representing interstellar dust which has already crossed the shock, , ,SN , and partly experienced sputtering.They apparently experience acceleration episodes from the RT "tongues" that advance the shock front and penetrate outward 3 .Behind the front their density drops within a thin layer of < ∼ 5 pc, as can be observed in lower panel of Fig. 1.

Spatial dust distribution pattern around the shell
Dust particles are initially placed on the computational grid as described above in Sec. 2 and in Sec.A1.Dust particles can escape their original "mother" cells and move into neighbour ones, in those cases when the stopping time becomes longer than the crossing time for a grid cell sound , and the collisional coupling weakens, see Sec. 3.4.As a result, depending on gas density and number of particles the overall distribution of dust particles can become patchy, with a lack or deficiency of dust particles at certain regions of the bubble and in the averaged radial profiles.Fig. 2 shows the distribution of gas and dust density in and around the remnant at the radiative The angle averaged radial profiles of gas density ( , red line), gas temperature ( , green line), density of the injected dust ("i") particles ( , ) with sputtering (blue line) and without it (dark-blue line), density of the interstellar (pre-existed, "e") dust particles ( , ) with sputtering ( magenta dash line) and without it ( light-magenta dash line) and density of the interstellar dust particles behind the SN shock front ( , , ) with sputtering ( magenta solid line) and without it ( light-magenta solid line) at 50 and 300 kyrs (upper and lower panels, respectively).The initial grain size of the injected and interstellar particles is equal to the fiducial value: stage.The middle panel presents the positions of dust particles in a 2D slice (a single cell in thickness) at a certain time after explosion.Initially the interstellar dust particles are deposited to the computational grid homogeneously with a density of one particle per a cell in its centre, as seen in Fig. 2 in the unperturbed region.Evacuation of dust particles due to decoupling are clearly seen in several regions behind the front.Within the hydrodynamical description dust particles along with their physical characteristics mentioned in Sec.A5 are smeared out over the computational cell resulting in mitigation of this patchiness.This results in reducing the area with lacking dust seen in the right panel of Fig. 2, as compared with the pattern shown in the middle panel.The radial profiles shown in Fig. 1 include this effect, see also discussion in Sec.3.4.In this regard it is also worth noting that a relatively small number of "superparticles" -one per computational cell, may cause numerical artefacts in the immediate post-chock domain, where dust and gas decouple (see Fig. 8 in Sec.3.4).The smearing out procedure partly mitigates possible contaminations from such artefacts.
In the left panel we can see a tendency of gas density to vary quasiperiodically along the thin shell.Similar signs of quasi-periodicity can be found also in the distribution of dust "superparticles" in the middle panel.These signs indicate development of both the thermal instability and the Vishniac-Ostriker overstability (Ostriker & Cowie 1981;Vishniac 1983;White & Ostriker 1990) at their initial stages.It is difficult though to tell them apart at early times, whereas at later stages the thermal instability results in fragmentation of the shell and seems to hinder further development of the Vishniac-Ostriker overstability.However, the shell expansion at > ∼ 50 kyr follows the standard Oort law ∝ 1/4 , in which case the growth rate of the Vishniac-Ostriker overstability ≃ 0 (see Eq. 2.17d in Vishniac 1983).More recent on-puprose numerical study of the Vishniac-Ostriker overstability by Minière et al. (2018) shows that perturbations from the interior regions cause the instability to reappear at the radiative phase.In our case during the period from 50 to 300 kyr we do not observe signatures of the Vishniac-Ostriker overstability in the cooled shell while it becomes denser.The discrepancy can be due to differences in the initial perturbations introduced into the computational domain.However, it is worth noting that effects from possible growth of the Vishniac-Ostriker overstability cannot affect dust destruction.

Pre-existing dust
The interstellar grains in the swept up gas of the SN shell are destroyed efficiently during the early expansion < ∼ 50 kyr.This can be observed in the upper panel ( = 50 kyr) of Fig. 1, where the dust destruction is lower in the immediate post-shock region, whereas in the hotter innermost domain it remains considerable.The reason is that dust particles of radius > ∼ 0.1 m are destroyed on a timescale of ∼ 10 2 −1 0.1 kyr comparable to radiative cooling time cool ∼ 10 2 −1 kyr.The remnant enters the radiation dominated phase at ∼ 30 kyr, and as a consequence dust destruction rate decreases with decreasing temperature.Solid lines in Fig. 3 depict the mass of the interstellar dust included in the SN bubble.On longer timescale, sputtering rate falls below the necessary level, and asymptotically the swept up dust mass grows proportionally to ∝ 3 ∼ 3/4 .For comparison, the dot-dashed red line shows the dust mass density profiles for the case without dust destruction.

Injected dust
Dashed lines in Figure 3 show evolution of the total mass of dust populations of various initial sizes 0 at the time of injection = 0 by the SN.Small injected particles are destroyed quickly with characteristic sputtering time ∝ , such that grains of initial size 0.05 ≤ m are mostly sputtered away within 10 kyr after SN explosion.During the subsequent evolution the shock velocity and the gas temperature in the bubble and in the ejecta decrease, and consequently the efficiency of dust destruction also decreases, as can be judged from dust radial profiles on lower panel of Fig. 1.The apparent increase of the interstellar dust mass at later stages, > 100 kyr, is due to the swept up dust at the Oort expansion phase, ∝ 3 ∝ 3/4 .As a result, sputtering of injected dust particles of all sizes levels out at ∼ 100 kyr.In the end, the mass of surviving dust is ∼ 5% for particles of 0 = 0.1 m to ∼ 40% for 0 = 0.5 m.Before going into details about dust destruction on smaller scales, let us consider first more generic averaged behavioural features.This gives us an impression of a collective dynamics of the injected particles depending on their initial radii 0 .For this we calculate MNRAS 000, 1-?? ( 2015) Dust during the SN-ISM interaction 5 the averaged distance 0 ( ) of dust particles injected by the SN at = 0 with radius 0 .At > ∼ 100 kyr the SN remnant reaches the radiative stage and forms a relatively dense and thin shell.The injected dust particles, which are initially surrounded by the hot gas of the ejecta are destroyed mostly by thermal sputtering.The fraction of the remnant volume, that is occupied asymptotically at ≥ 100 kyr, depends on the particles' initial size.Fig. 4 illustrates the evolution of 0 ( ) for different 0 .At early times ( < ∼ 10 kyr), during the free expansion and and Sedov-Taylor phases, 0 ( ) shows a weak dependence on 0 .However, it changes further on: small particles, 0 < 0.1 m, are tightly coupled to the ejecta and the hot bubble interior, and remain within the inner 30-36 pc of the bubble.During the next ∼ 200 kyr such particles along with the hot  ejecta overtake the dense shell and remain inside the it.Heavier dust particles, 0 ≥ 0.1 m, reach either the post-shock layer, or even overcome it, as 0 = 0.5 m particles do.

Dust survival vs gas ambient density
The surviving fraction of the injected dust with a growing ambient density is shown by dashed lines in Fig. 5.As the ambient density increases the expansion velocity of the forward SN shock decreases, and as a consequence so does the expansion velocity of the ejecta, thus prolonging the exposure of the injected dust to the hostile en- vironment.Even though the SN shell enters the radiation dominated phase in a shorter time, it is sufficient to boost destruction of dust within the ejecta.
The interstellar dust is much less sensitive to the ambient density as shown in Fig. 5 by solid lines.This is connected with the fact that at gas temperature > 10 6 K for ∼ 0.1 m the sputtering time sp ∼ 10 5 −1 0.1 s, is nearly equal to the gas cooling time, cool ∼ 10 5 −1 s at ∼ 3 × 10 6 K, with ∼ ⊙ .For dust particles of ≥ 0.03 m the sputtering time is only ≤ 3 times shorter than cooling time.Below < 10 6 K the sputtering rate decreases as ∼ 3 , resulting in simultaneous reduction of dust destruction.As the gas cooling rate decreases later, at < ∼ 10 5 K, it is possible that the sputtering rate levels out before gas cools considerably.As the SN remnant becomes radiative at ≥ 50 −1 kyr, the independence of the survival yield of the pre-existing dust shown in Fig. 5 is consistent with this scenario.Moreover, an apparent increase of the surviving fraction of the interstellar dust can even be observed at ambient density > 0.2 cm −3 .This is related to the higher efficiency of radiative cooling at higher densities and faster transition of the shock velocity to the Oort law ∝ −3/4 .Similar results has been recently described by Kirchschlager et al. (2022).In this regard it is worth noting that the survival yield shown here is already established within the intermediate asymptotic at > ∼ 50 −1 kyr.

Time variation of dust sizes
Figure 6 shows the evolution of the grain size distributions by mass inside the SN bubble for the injected and interstellar particles of different initial sizes: 0 = 0.3, 0.1 and 0.05 m.The distribution function is defined as the mass fraction of grains in a given size bin Δ = / −1 as where is the bubble radius, i.e., the radius of the forward shock, logΔ = 0.165 with the number of bins being = 15 within e, 50 kyr e, 100 kyr e, 300 kyr i, 0 kyr i, 50 kyr i, 100 kyr i, 300 kyr Figure 6.The mass distribution functions of grain sizes for the injected (filled circles connected by dashed lines, "i") and the interstellar (open circles connected by solid lines, "e") monosized particles of initial (at = 0) sizes: 0 = 0.3, 0.1 and 0.05 m from top to bottom panels, respectively, in the SN bubble with age of 50, 100 and 300 kyr (red, green and blue lines).For 0 = 0.3 m (top panel) the distributions of the injected particles at = 100 and = 300 kyr coincide.For the injected dust the initial distribution is marked by the black solid triangle symbol.The average number density of the ambient gas is 1 cm −3 .0.01 − 0.3 m.Large grains are destroyed slower than the small ones proportionally to their radius sp ∝ .Consequently, within the first 50 kyr the injected particles with 0 = 0.3 m decrease in size by less than a factor of 2 and stay unchanged beyond = 100 kyr.This is consistent with the total mass of surviving particles shown in Fig. 3.
The interstellar particles that crossed the front are destroyed less MNRAS 000, 1-?? ( 2015) Dust during the SN-ISM interaction 7 Figure 7.The same as in Figure 6, but for particles of initial size 0 = 0.1 m and mumber densities of the ambient gas 3 cm −3 (upper panel) and 0.3 cm −3 (lower panel).
efficiently, because they spend in the postshock hostile region a shorter time, and only a negligible fraction pervade into the hot bubble and the ejecta.This is connected with the fact that at > 50 kyr the postshock gas becomes radiative, gas temperature falls below < ∼ 10 6 K and its density increases correspondingly.As mentioned above, the sputtering rate | |/ ∝ 3 decreases faster than the cooling rate | |/ ∝ −1 , which leads to a decrease in the dust destruction rate in the shell.The dust particles in their turn remains collisionally coupled to the gas, / ≃ 1, and after passing a few cells behind the front they remain in the denser and cooler shell gas.After ∼ 50 kyr only a minor fraction of the interstellar dust experinces sputtering.This explaines a slow growth of the surviving fraction of the interstellar dust.
The injected dust grains of smaller initial sizes 0.1 and 0.05 m are destroyed on a shorter time scale sp ∝ .The destruction is more efficient during the very initial phase < ∼ 50 kyr when the ejecta is denser and hotter.Within ∼ 50 − 100 kyr the maximum of ( ) shifts by factor 2 toward smaller radii.On longer timescales between 100 and 300 kyr the sputtering rate decreases slightly because of a lower temperature in the ejecta.
Increase of the ambient gas density enhances dust sputtering.The resultant net effect can be qualitatively described as a shift of the size distribution function ( ) ∼ ( / ) at a given time, as can be found from comparison of the distributions in Fig. 7 with the corresponding curves on Fig. 6. ∼ 30 − 40 kyr the SN remnant passes the initial period of the radiative phase, the borders practically coincide; the time is shown on the right -axis.The interstellar particles lose collisional coupling starting from entering the shell and several (up to 10 pc) deeper inward.Color coding for the profiles corresponds to different ages: 50 kyr (red), 100 kyr (green), 200 kyr (blue) and 300 kyr ( magenta).

Purity of the interstellar and injected dust
The overall dynamics of dust destruction in either cases -the interstellar and the injected particles, is predominantly determined by the interrelation between the collisional coupling (stopping) time , the sputtering time sp and the local sound-crossing time sound = Δ / .The stopping time is roughly ∼ / (Sec.A2), the thermal sputtering time at > ∼ 10 6 K is sp ∼ / , is the mass of a target particle (Draine 1995).As can be inferred from Fig. 1 the ratio sp / ∼ H / < ∼ 0.1 for > ∼ 10.At early times ≤ 50 kyr the interstellar pre-existing particles decouple from the gas shell and escape, since / > ∼ 1. this ratio decreases as dust passes through the shell / < ∼ 1 (Fig. 8).Thus, inside the shell, the dust connects to the postshock gas flow and are then destroyed.Only negligible fraction of them penetrate deeper into the hot bubble and ejecta.This can be seen from dust density radial profiles in Fig. 1, where the amount of surviving interstellar dust beneath the shell at > ∼ 50 kyr falls below 7%.During the evolution the injected and interstellar dust remain practically isolated: the former being locked in the ejecta and the hot bubble, the latter is swept up into the shell.Only a minor fraction of the interstellar dust can reach deep in the bubble interior under conditions of weak collisional coupling / ≫ 1. Mild signs of intermixing can be seen only in a thin layer deeply behind the shock front, as can be seen in the radial distributions of dust density in Fig. 1.
The tendency of the injected and the interstellar dust to stay in the remnant rather unmixed can also be seen in Fig. 9, where the radial profiles of their abundances are given at different times.Both the interstellar and the injected dust-to-gas ratios are normalized to the initial value 0 = 10 −2+[Z/H] of the interstellar dust.As the dust-to-gas ratio of injected dust is set = 1/30 (see in Sec. 2), its normalized dust-to-gas ratio in the lower panel of Fig. 9 for [Z/H] = −1 is higher than in the upper panel for [Z/H] = 0.At = 50 kyr, a considerable fraction of the injected dust in the thermalized bubble is already destroyed, ∼ 0.8.During the following expansion over the next 150 kyr, dust destruction becomes less efficient: at > 100 kyr around 8-9 % of the injected dust particles survive, at > ∼ 300 kyr the injected dust practically merges to the shell and the surviving dust fraction falls below ∼ 5%; no injected dust particles penetrate beyond the forward shock.The interstellar dust is efficiently destroyed behind the shock only at early stages < ∼ 100 kyr, while at later stages the post-shock gas cools radiatively and temperature falls below the sputtering threshold, as discussed above, Sec.3.3.2.In both high and low metallicity gas: [Z/H] = 0 and [Z/H] = −1, intermixing between the injected and the interstellar particles is seen in the rather narrow region inside the bubble.While the injected particles tend to occupy practically the entire bubble, the interstellar dust penetrate inwards until being mixed into expanding remnant gas, since | − | ≪ .As sp / < ∼ 0.3 only a negligible amount of interstellar dust can penetrate deep into the hot bubble interior, resulting in a sharp drop of their abundance seen in Fig. 9.

DISCUSSION AND CONCLUSION
The dust properties in the ISM can be spatially inhomogeneous due to contributions from the two dust populations: the population injected by SNe and the shock processed pre-existing one.Both depend on the local gas density and the local SN rate.For the ambient gas density ≥ 0.1 cm −3 only < ∼ 13% of injected dust particles with radii 0 ∼ 0.1 m do survive; at higher this fraction decreases as −1.2 .Larger particles are less sensitive to the ambient density: the surviving fraction of dust grains with 0 ∼ 0.3 m is more than 30% and weakly depends on ambient density ∝ −0.15 .
Destruction of the pre-existing dust particles in ambient gas is less severe.Estimates of the dust budget for a SN remnant at its radiative stage, > ∼ 30 −1/3 kyr, show that up to ∼ 80% of the swept up interstellar dust survives.Therefore, after ≥ 30 −1/3 kyr, the remnant shell contains the mass of the surviving swept up dust up to ≥ 0.8 , ∼ 32 2/5 ⊙ .On the other hand, the mass of the surviving injected dust in the remnant is only ≤ 0.13 , , which results in a ratio of injected-to-interstellar dust in the remnant at its later stages of , / , ≥ 0.004.This impurity is too low to be detected in the swept up shell.Possible observational signatures of differences between optical properties of the pre-existing and the shock processed dust can be recognized in the interface region separating the dense SN shell and its rarefied bubble.Spatial variations of the exinction law described by Fitzpatrick & Massa (1990) may be connected with contributions from such interfaces.Quite recently signatures of the two populations of dust with rather distinct spectral features are recognized in the Planck emission bands at high Galactic latitudes (Shiu et al. 2023, their Figs. 3, 4 & 7).
Dust supply from SN explosions in galaxies at < 5 plays apparently a minor role, particularly when the small dust particles with < 0.1 m are concerned.The survival percentage of > ∼ 0.1 m dust in high density environment > 1 cm −3 is less than 10%.This estimate is qualitatively consistent with early consideration by Nozawa et al. (2006Nozawa et al. ( , 2007)).The percentage of particles of small sizes is much lower than 1% in a > ∼ 1 cm −3 environment.Our results show that the survival fraction of dust particles of ∼ 0.1 m is , ∼ 10 −2 ¯ , − , is the ambient density in cm −3 , > 1, ¯ , is the dust mass injected by a single SN.This results in an estimate of dust supply rate due to SN where ¯ , is in 1 ⊙ , SFR is the SF rate in 1 ⊙ yr −1 , for the specific supernova rate ¯ ∼ 10 −2 ⊙ .For comparison, the typical SN dust destruction rate in the Milky Way is (see, e.g., in Draine 2009) This estimate seems to be applied particularly for dust production at the "cosmic noon" epoch, where the SNe injected dust is to be efficiently destroyed in the ejecta themselves at a higher density environment.SN explosions as the the primary dust-formation channel in galaxies can be important only for > 5 galaxies.Moreover, in most distant galaxies at > 10 only large (grey) dust can survive the extensive sputtering in SN shocks, because of their on average higher density environment ∝ (1 + ) 3 .Estimates of the global dust budget in the ISM of high-redshift galaxies have to take into account highly sensitive sputtering yields on the shock velocity and ambient gas density, particularly when small dust particles are considered.At such conditions, a general view that SN explosions in high-galaxies are exceptionally selective with respect to a preferential destruction of small size dust particles.At such circumstances possible growth of dust particles in situ in interstellar medium, as recently discussed (Hirashita & Kuo 2011;Mattsson 2020a;Hirashita & Chen 2023).

APPENDIX A: DYNAMICS OF DUST PARTICLES
Here we describe an implementation of the particle dynamics in our gasdynamic code (Vasiliev et al. 2015(Vasiliev et al. , 2017(Vasiliev et al. , 2019)).This code has successfully passed the whole set of tests proposed in Klingenberg et al. (2007).Several additional tests have been given in the appendix of Vasiliev et al. (2017).
Our description of dust dynamics basically follows to the method proposed by Youdin & Johansen (2007), several parts are similar to Mignone et al. (2019) and Moseley et al. (2023).We add several points allow to discriminate sorts of particles, to follow the evolution of dust (macroscopic) mass.We include destruction processes due to thermal and kinetic sputtering.Since we are going to study gas-dominated (by mass) flows in the ISM, we use explicit methods for solving the equations of dust dynamics.

A1 Properties
Dust particles trace the motion of grains in a gas.They are suffered by drag force from a gas, and due to it they transfer momentum and heat to a gas.Also they may be exposed to other external forces like gravity, radiation and so on.
Each dust particle is described by several features allow one to identify its evolution and physical properties, namely, colour, time of injection, sort, size and mass.Two former can be used to identify a source of each particle.Dust particles can consist of different material -a sort of particle, e.g, carbon or silicate.They can be of various sizes, which are distributed according to some spectrum or have a single size.In the former case the dust is initially polydisperse and one can choose several sizes distributed by some spectrum.In the latter it is monodisperse, but due to destruction processes it becomes polydisperse.These features are microscopic.
To follow the transport of dust mass (not of an individual grain) in a gaseous flow we introduce a macroscopic mass of a 'superparticle' or macroparticle.In this approach a particle is a conglomerate of microscopic grains.This value of mass (or supermass) is used to find the dust-to-gas mass ratio and so on.If it is necessary, any feature of dust particle can be added.can see a good coincidence between our numerical results and the above-mentioned analytic solution (depicted by black dashed lines).Only for the coarse grid ( = 100) and the shortest stopping time one can find small deviations from the analytic curve.

B2 SN evolution: a convergence
We follow the evolution of a single SN remnant with different spatial resolution for the grid and fixed number of dust particles injected by SN.The number of dust particles in the ambient gas (or preexisting dust particles) is proportional to the grid resolution, because we set one such particle per grid cell.Figure B2 presents the onedimensional (along line of sight crossed the center of SN bubble) profile of the dust density depositing to the grid for different spatial resolution of the grid: from 1 to 0.1875 pc.The slice is for a SN bubble with age of 20 kyr.The number of injected particles remains the same and equal to 8 millions.For comparison the gas density profile along the same line-of-sight and for the highest spatial resolution (0.1875 pc) is added.The dust density is multiplied by a factor of 100 for better presentation (the dust-to-gas density ratio in the ambient gas is initially set to 0.01).One can see that the dust distributions for the cell size lower than 0.5 pc are close.Therefore, we conclude that the grid resolution of 0.5 pc is sufficient to follow the dust dynamics during the SN bubble evolution.The gaps in profiles seen for high resolution are due to a limited number of particles in a cell.
Figure B3 shows the one-dimensional distribution of dust density for the runs with different number of dust particles: 8, 16 and 32 millions, but the fixed grid resolution equal to 0.5 pc.The slice is for a SN bubble with age of 20 kyr.It is seen that the dust density profiles for the runs with different number of particles are very close.The number of dust particles decreases within the bubble faster than the gas density, because of decoupling in the innermost regions of the bubble.Therefore, for our runs we adopt the spatial resolution equal to 0.5 pc, and the number of injected dust particles equal to Figure1.The angle averaged radial profiles of gas density ( , red line), gas temperature ( , green line), density of the injected dust ("i") particles ( , ) with sputtering (blue line) and without it (dark-blue line), density of the interstellar (pre-existed, "e") dust particles ( , ) with sputtering ( magenta dash line) and without it ( light-magenta dash line) and density of the interstellar dust particles behind the SN shock front ( , , ) with sputtering ( magenta solid line) and without it ( light-magenta solid line) at 50 and 300 kyrs (upper and lower panels, respectively).The initial grain size of the injected and interstellar particles is equal to the fiducial value:

Figure 2 .Figure 3 .
Figure 2. A part of the SN shell at the radiative stage: 2D slices of the gas density (left), positions of dust particles (middle), average dust density (right) at 50 kyrs.The initial grain size is equal to the fiducial value: 0 = 0.1 m.

Figure 4 .
Figure 4.The averaged distance 0 ( ) of ensembles of the dust particles from the point of their injection, i.e. the point of a SN explosion.The value of 0 ( ) is weighted to the current particle's supermass : 0 ( ) = /.The particles have been injected at = 0 with the initial radii from 0 = 0.03 m to 0.5 m as shown in the legend.Color lines mark 0 ( ) for an ensemble corresponding to a given 0 ; the thick grey line depicts the radius of the SN bubble.

Figure 5 .
Figure5.Dependence of survival fraction of the pre-existing (solid lines) and injected (dashed lines) dust at the asymptotic ≥ 100 kyr versus the ambient gas density; the horizontal thin solid line is the fraction without destruction.

Figure 8 .
Figure8.The radial profiles averaged over the solid angle for the size avegared ratio of the grain stopping time of the injected and interstellar dust to the local sound time sound , 0 = 0.1 m.Solid lines with open symbols are for the interstellar dust, dashed lines with filled symbols are for the injected particles.Shaded band represents evolution of the bubble shell: the right border corresponds to the forward front where the gas velocity jumps, and the left one shows the border where the gas temperature ≥ 10 5 K.At < ∼ 30 − 40 kyr the SN remnant passes the initial period of the radiative phase, the borders practically coincide; the time is shown on the right -axis.The interstellar particles lose collisional coupling starting from entering the shell and several (up to 10 pc) deeper inward.Color coding for the profiles corresponds to different ages: 50 kyr (red), 100 kyr (green), 200 kyr (blue) and 300 kyr ( magenta).

Figure 9 .
Figure9.The radial distribution of the normalized dust density (i.e. the dustto-gas ratio normalized to the initial background value 0 = 10 −2+[Z/H] ) at times 50, 100, 200 and 300 kyr (color lines).Solid lines show the interstellar dust: thin lines represent the case without dust destruction, thick lines illustrate the case when destruction is accounted.Dashed lines depict the injected dust: thin lines show dust without destruction, whereas the thick ones are for dust with destruction being accounted.Shaded area represents the evolution of the bubble shell as defined in Fig.8upper panel shows the model with ambient gas of metallicity [Z/H] = 0 and = 10 −2 , lower panel is for [Z/H] = −1 and correspondingly = 10 −3 , the injected dust mass is = 1 ⊙ in both cases.

Figure B1 .
Figure B1.The velocity of dust particles for stopping times = 10 5 , 2 × 10 5 , 3 × 10 5 and 4 × 10 5 yr (red solid lines from top to bottom, respectively) in a homogeneous gas flow for grid with number of cells along each spatial direction = 200.The green and blue lines show the velocity for = 10 5 yr and number of cells = 200 and 400.The black dash lines represent the analytic solution for stopping times = 10 5 , 2 × 10 5 , 3 × 10 5 and 4 × 10 5 yr (from top to bottom).

Figure B2 .Figure B3 .
Figure B2.The one-dimensional profile of the dust density depositing to the grid for different spatial resolution of the grid: the colour lines correspond to the cell size from 1 to 0.1875 pc.The slice is for a SN bubble with age of 20 kyr.The number of injected particles in all models is equal to 8 millions.The thick grey line depicts the gas density profile along the same line of sight and for the highest spatial resolution (0.1875 pc).The dust density is multiplied by a factor of 100.