# Observable signatures of enhanced axion emission from protoneutron stars

Tobias Fischer<EMAIL_ADDRESS>Institute of Theoretical Physics,
University of Wrocław, 50-204 Wrocław, Poland.  Pierluca Carenza
<EMAIL_ADDRESS>Dipartimento Interateneo di Fisica “Michelangelo
Merlin”, Via Amendola 173, 70126 Bari, Italy.  Istituto Nazionale di Fisica
Nucleare–Sezione di Bari, Via Orabona 4, 70126 Bari, Italy.  Bryce Fore
<EMAIL_ADDRESS>Institute for Nuclear Theory, University of Washington, Seattle,
Washington 98195, USA.  Department of Physics, University of Washington,
Seattle, Washington 98195, USA.  Maurizio Giannotti<EMAIL_ADDRESS>Physical Sciences, Barry University, 11300 Northeast 2nd Avenue, Miami Shores,
FL 33161, USA.

Alessandro Mirizzi<EMAIL_ADDRESS>Dipartimento Interateneo di
Fisica “Michelangelo Merlin”, Via Amendola 173, 70126 Bari, Italy.  Istituto
Nazionale di Fisica Nucleare–Sezione di Bari, Via Orabona 4, 70126 Bari,
Italy.  Sanjay Reddy<EMAIL_ADDRESS>Institute for Nuclear Theory, University
of Washington, Seattle, Washington 98195, USA.  Department of Physics,
University of Washington, Seattle, Washington 98195, USA.

(12 August 2021; 11 October 2021; 9 November 2021)

###### Abstract

We perform general relativistic one-dimensional supernova (SN) simulations to
identify observable signatures of enhanced axion emission from the pion-
induced reaction $\pi^{-}+p\to n+a$ inside a newly born protoneutron star
(PNS). We focus on the early evolution after the onset of the supernova
explosion to predict the temporal and spectral features of the neutrino and
axion emission during the first 10 s. Pions are included as explicit new
degrees of freedom in hot and dense matter. Their thermal population and their
role in axion production are both determined consistently to include effects
due to their interactions with nucleons. For a wide range of ambient
conditions encountered inside a PNS, we find that the pion-induced axion
production dominates over nucleon-nucleon bremsstrahlung processes. By
consistently including the role of pions on the dense matter equation of state
and on the energy loss, our simulations predict robust discernible features of
neutrino and axion emission from a galactic supernova that can be observed in
terrestrial detectors. For axion couplings that are compatible with current
bounds, we find a significant suppression with time of the neutrino luminosity
during the first 10 s. This suggests that current bounds derived from the
neutrino signal from SN 1987A can be improved, and that future galactic
supernovae may provide significantly more stringent constraints.

DOI: 10.1103/PhysRevD.104.103012

##  I Introduction

Core-collapse supernovae (SN) are ideal high-energy nuclear and particle
physics laboratories, with baryon densities reaching a few times nuclear
saturation density at their interiors and temperatures of several tens of MeV,
as well as large isospin asymmetry Fischer:2017 . The latter is quantified by
the proton abundance ranging from extremely neutron-rich conditions with
$Y_{p}\simeq 0.05$ to isospin symmetric conditions with $Y_{p}=0.5$. The
central compact stellar object of a core-collapse SN is the protoneutron star
(PNS). It forms when the collapsing stellar core of a massive star bounces
back at densities slightly in excess of nuclear saturation density ($\rho_{\rm
sat}\simeq 2.5\times 10^{14}$ g cm-3). The later postbounce evolution of the
PNS is determined by continuous mass accretion, from the still-collapsing
outer layers of the progenitor star. When the SN explosion develops—because of
the revival of the stalled bounce shock—the nascent PNS will deleptonize via
the emission of neutrinos of all flavors on a timescale of several tens of
seconds Pons:1998mm  and later cool to form the final SN remnant neutron star.
While for the understanding of the SN explosion, in particular aided through
the neutrino-driven paradigm, multidimensional phenomena play a crucial role
Janka:2007 ; Mirizzi:2016 , the PNS deleptonization is well approximated
within spherically symmetric simulation studies Huedepohl:2010 ;
Fischer:2009af ; Roberts:2012b ; Roberts:2012c ; Fischer:2016a .

With the recent advances in the modeling of core-collapse SN, which now allow
the inclusion of muons as the degree of freedom and to treat six-species
Boltzmann neutrino transport schemes Bollig:2017 ; Fischer:2020d , it becomes
possible to study novel high-energy astroparticle physics phenomena. Besides
muons the next heavier degrees of freedom to consider are the pions. These
bosons exist as charge-neutral $\pi^{0}$ and as charged $\pi^{\pm}$. Pions are
one of the key messengers in the context of heavy-ion collisions (for a recent
comprehensive review, see Ref. Redlich:2018 ), which makes them one of the
best-studied particle degrees of freedom. The role of pions in neutron stars
has been studied for a long time, particularly in relation to the emergence of
a pion condensate. However, in the context of core-collapse SN, only the core-
collapse SN simulations of Ref. Wilson:1993rx  included pions in the equation
of state (EOS). Recently, the classical neutron star EOS of Akmal,
Pandharipande and Ravenhall, denoted as APR AP:1997 ; APR:1998 , was extended
to finite temperatures and arbitrary isospin asymmetry applying a Skyrme
functional model, prepared for astrophysical applications Schneider:2019 .
Pions are included through the phase transition to a neutral pion condensate
at supersaturation density. This results in a softening of the EOS. As was
pointed out in Ref. Fore:2020 , an even greater impact on the high-density EOS
emerges already below saturation density when considering the contributions
from strongly interacting pions, giving rise to the pion self-energy. The
resulting pion yields can be comparable to those of muons, depending on the
nuclear EOS. In this article we employ this formalism, based on the virial EOS
featuring pion-nucleon scattering phase shifts, and study the appearance of
pions in simulations of the PNS deleptonization.

The inclusion of pions in the EOS allows us to describe, for the first time,
the full impact of axions on the SN evolution in a self-consistent way. Axions
are a prediction of the most promising solution of the strong CP problem,
which is related to the absence of the expected CP violation in the strong
interactions Peccei:1977hh ; Peccei:1977ur ; Weinberg:1977ma ; Wilczek:1977pj
. The axion phenomenology is widely discussed in the literature (for recent
reviews, see Refs. DiVecchia:2019ejf ; DiLuzio:2020wdo ; Agrawal:2021dbo ). It
is essential to remark that axions are expected to couple to photons,
electrons, and nucleons, with model-dependent couplings diCortona:2015ldu ;
Borsanyi:2016ksw . These can be constrained through laboratory experiments
Irastorza:2018dyq ; Sikivie:2020zpn  and astrophysical considerations
Raffelt:2006cw ; Giannotti:2015kwo ; Giannotti:2017hny . In particular, core-
collapse SN have long been used to constrain the axion couplings to nucleons
Carenza:2019pxu , photons Payez:2014xsa ; Calore:2020tjw , and, recently, to
muons Bollig:2020 . In the present paper, we are particularly interested in
axions coupled to nucleons. In this case, the most widely studied mechanism
for axion production in a SN core is the nucleon-nucleon bremsstrahlung,
$N+N\to N+N+a$ Raffelt:1987yt ; Turner:1987by ; Burrows:1988ah ;
Burrows:1990pk ; Raffelt:1990yz ; Carenza:2019pxu . Because of the efficient
energy loss through axions, the above process contributes to the shortening of
the neutrino emission from a core-collapse SN. Quantitatively, this time
reduction is parameterized by axion-nucleon coupling Fischer:2016b . The
current bound on the axion-nucleon coupling has been revisited Carenza:2019pxu
, with improved nucleon-nucleon bremsstrahlung rates for the axion emissivity
taking into account leading-order medium modifications. Furthermore, in the
context of core-collapse SN, there is a second competing process for the
emission of axions, which has long been omitted, namely the emission of axions
from Compton pionic processes, $\pi^{-}+p\to n+a$ Turner:1988 ; Raffelt:1993ix
; Keil:1996ab . This process has been reviewed recently in the context of
strongly interacting pions under the conditions of hot neutron stars
Carenza:2020 . It leads to the conclusion that the previous assumption of the
dominating axion emission from nucleon-nucleon bremsstrahlung over those from
pions must be relaxed.

In order to quantify the impact of axion losses originating from pions, the
present study implements the associated axion emissivities Carenza:2020  in
simulations of the PNS deleptonization phase of a core-collapse SN explosion,
together with the updated axion emissivities from nucleon-nucleon
bremsstrahlung Carenza:2019pxu . Particular emphasis is devoted to the impact
of the shortening of the neutrino emission due to the associated enhanced
cooling, and the potentially observable signature at the future-planned Hyper-
Kamiokande and presently operating Super-Kamiokande water-Cherenkov detectors.

The paper is organized as follows. In Sec. II, we review our core-collapse SN
model with the necessary updates to include pions. In Sec. III, we introduce
the axion emission processes, specifically the improved nucleon-nucleon
bremsstrahlung rates and the pion rates. Their impact on the PNS
deleptonization will be discussed in Sec. IV. In Sec. V we present the
observable signatures of the axion emission and discuss the impact on the
neutrino signal. The paper closes with a summary in Sec. VI.

##  II Supernova model with pions

The spherically symmetric core-collapse SN model employed in this study,
AGILE-BOLTZTRAN, is based on general relativistic neutrino radiation
hydrodynamics Mezzacappa:1993gm ; Mezzacappa:1993gn ; Mezzacappa:1993gx ;
Liebendoerfer:2004 . It includes an adaptive baryon mass mesh refinement
Liebendoerfer:2002 ; Fischer:2009af , where for the present study we use 207
radial mass shell grid points. AGILE-BOLTZTRAN has been extended to handle
six-species Boltzmann neutrino transport Fischer:2020d  with the inclusion of
the muon abundance $Y_{\mu}$, as an additional degree of freedom, together
with a comprehensive set of muonic weak interactions in the collision integral
of the Boltzmann transport equation Guo:2020 ; Fischer:2020d . The list of
standard, nonmuonic weak processes used here can be found in Table I of Ref.
Fischer:2020a , with the charged-current weak rates, including the (inverse)
neutron decay channel, in the full-kinematics treatment as well as including
self-consistently contributions from weak magnetism. The latter two are
particularly important for the PNS deleptonization phase, as was discussed in
Ref. Fischer:2020a . The neutrino distributions and transport are discretized
with six momentum angles, $\cos\vartheta\in[-1,1]$, and 36 neutrino energy
bins, $E_{\nu}\in[0.5,300]$ MeV, following the setup of Bruenn Bruenn:1985en .

AGILE-BOLTZTRAN has a flexible EOS module which can handle a variety of
different nuclear matter EOS Lattimer:1991nc ; Shen:1998gg ; Hempel:2009mc ;
Hempel:2011 ; Steiner:2013 . In the present work the DD2 relativistic mean-
field EOS with density-dependent nucleon-meson coupling is employed Typel:2005
; Typel:2009sy ; Typel:2013rza , together with the modified nuclear
statistical equilibrium EOS of Ref. Hempel:2009mc  for the description of
medium nuclear clusters, henceforth denoted as HS(DD2). It provides a good
description of nuclear saturation properties as well as neutron star
configurations, consistent with current maximum mass constraints from the
observations of massive pulsars of about 2 $M_{\odot}$ Antoniadis:2013 ;
Fonseca:2016 . Furthermore, electrons, positrons and photons are treated
following Ref. Timmes:1999 , while the muon EOS is based on a tabulation
Fischer:2020d .

The present work extends previous studies of axion emission in core-collapse
SN Fischer:2016b , which exclusively focused on nucleon-nucleon
bremsstrahlung. Here, we add a novel process of axion emission involving pions
Carenza:2020 . To study self-consistently the effects of this process in the
SN evolution, pions have to be included in the SN simulation. Pions have a
mass of $m_{\pi^{0,\pm}}=135(140)$ MeV and are, hence, only slightly more
massive than muons. Unlike muons, which are not in chemical equilibrium and,
hence. are exclusively produced from weak processes under SN conditions, pions
are in chemical equilibrium through strong interactions. Hence, it is not
required to solve an evolution equation for the pion abundances, as is the
case with muons Bollig:2017 ; Fischer:2020d . Instead, their abundances can be
calculated from the corresponding local thermal and chemical equilibrium
conditions for a gas of massive and relativistic bosons. This implies that the
pions have a chemical potential, $\mu_{\pi^{\pm}}=\mp\hat{\mu}$, where
$\hat{\mu}\equiv\mu_{n}-\mu_{p}$ is the chemical potential for negative charge
and $\mu_{\pi^{0}}=0$.

Figure 1: Momentum $p_{\pi}$-dependent pion self-energy for $\pi^{-}$,
$\Sigma_{\pi^{-}}$, calculated based on the nucleon-nucleon potential, Eqs.
(21)–(22) in Ref. Fore:2020 , evaluated at a temperature of $T=20$ MeV and at
fixed electron and muon abundances of $Y_{e}=0.15$ and $Y_{\mu}=0.04$, as well
as two different baryon densities, at saturation density (solid line) and at
one-half of saturation density (dashed line). Figure 2: EOS comparison,
showing the particle yields $Y_{i}$ of protons ($p$) and pions ($\pi^{-}$),
the total pressure $P$, and the charged chemical potential $\hat{\mu}$, as a
function of the baryon density $\rho$, in units of the saturation density
$\rho_{0}$, comparing the model with only thermal pions (blue dashed lines)
and with interacting pions (brown solid lines), for the nucleon-nucleon
potential (Eqs. (21) and (22) in Ref. Fore:2020, ), evaluated at a temperature
of $T=20$ MeV and fixed electron and muon abundances of $Y_{e}=0.15$ and
$Y_{\mu}=0.04$.

Furthermore, pions are subject to strong interactions. At low density and high
temperature, the virial expansion provides an efficient and model-independent
approach to include pion-nucleon interactions. Under these conditions, when
the pion and nucleon fugacities are small, the modification of the pion
density and thermodynamics due to their interactions with neutrons and protons
is adequately incorporated through the second virial coefficients,
$b_{2}^{n\pi}$ and $b_{2}^{p\pi}$, respectively. We employ the methods
described in Ref. Fore:2020  to calculate these coefficients using measured
pion-nucleon scattering phase shifts.

To describe reactions involving pions and nucleons, we will need in addition a
model for the pion dispersion relation in the hot and dense medium. We adopt
the simple pseudopotential model described in Ref. Fore:2020  and write the
pion energy as

$E_{\pi^{-}}=\sqrt{p^{2}+m_{\pi^{-}}^{2}}+\Sigma_{\pi^{-}}(p_{\pi})~{}\,,$ (1)

where $\Sigma_{\pi^{-}}(p_{\pi}),$ is the momentum $p_{p}i$-dependent real
part of the pion self-energy. The strength of the pseudopotential is tuned to
ensure that the number densities obtained are consistent with the virial
expansionFore:2020 . The model includes the weakly repulsive s-wave
interaction and the strongly attractive p-wave interactions between pions and
nucleons. The latter has the effect of enhancing the $\pi^{-}$ and proton
abundances (see Fig. 2 in Ref. Fore:2020  and make a substantial contribution
to the EOS. Since the pion self-energy $\Sigma_{\pi^{-}}$ uses the knowledge
of the nucleon’s equilibrium Fermi-Dirac distribution functions, with the
associated nucleon chemical potentials, together with the pion-nucleon
potential—see expressions (21) and (22) in Ref. Fore:2020 —it becomes evident
that interaction contributions have to be consistent with the underlying
nuclear HS(DD2) EOS. The pion self-energies are shown in Fig. 1 for selected
conditions relevant for the PNS deleptonization. It becomes clear that their
contributions are indeed significant and cannot be neglected, neither for the
calculation of the pion contributions to the EOS nor for the emission of
axions stemming from pions, as will be discussed further below.

Note that the yields of $\pi^{+}$ are strongly suppressed due to the
contributions of the Boltzmann suppression term,
$Y_{\pi^{+}}\propto\exp\\{-\hat{\mu}/T\\}$, in comparison to
$Y_{\pi^{-}}\propto\exp\\{\hat{\mu}/T\\}$. A similar suppression is found for
$\pi^{0}$, in comparison to $\pi^{-}$. Furthermore, we have calculated the
self-energies for $\pi^{+}$ which are found to be similar to those of
$\pi^{-}$, and, hence, we find that at all conditions, the abundances of
$\pi^{+}$ and $\pi^{0}$ are indeed negligible.

With the inclusion of pions, the charge neutrality condition is modified. It
results in the following definition of the proton abundance:
$Y_{p}=Y_{e}+Y_{\mu}+Y_{\pi^{-}}-Y_{\pi^{+}}$, which is one of the three
independent nuclear EOS variables, in addition to temperature and baryon (or
rest mass) density. This results in a feedback to the baryon EOS, which is
illustrated in Fig. 2 for a selected temperature of $T=20$ MeV, an electron
fraction and a muon fraction of $Y_{e}=0.15$ and $Y_{\mu}=0.04$, lowering the
abundance of neutrons and increasing the abundance of protons. For the
calculation of the pion pressure and energy density, we used expressions (9)
and (10) of Ref. Fore:2020 .

The question about weak interactions involving pions, i.e. the $\pi^{-}$ decay
and its inverse, $\bar{\nu}_{\mu}+\mu^{-}\leftrightarrows\pi^{-}$ and
$\nu_{\mu}+\pi^{-}\leftrightarrows\mu^{-}$, has already been addressed in Ref.
Fore:2020 . Note that the latter is possible only in the presence of
interacting pions, i.e., when $\Sigma_{\pi^{-}}$ is large and negative to
ensure energy and momentum conservation. These reactions are competing not
only with the neutral current scattering processes, as was discussed in Ref.
Fore:2020 , but also with the muonic charged-current reactions as well as the
purely leptonic reaction channels Guo:2020 ; Fischer:2020d , giving additional
contributions to the $\nu_{\mu}$ and $\bar{\nu}_{\mu}$ opacity accordingly.

The conditions for these weak processes requires a finite abundance of
$\pi^{-}$. This is the case at densities in excess of about $\rho\simeq
10^{13}$ g cm-3. This, in turn, corresponds to temperatures fairly above
$T>10$ MeV during the SN evolution, where the muon (anti)neutrinos are trapped
entirely Fischer:2012a . As was found in Ref. Fore:2020 , these weak processes
would contribute mostly to low muon-(anti)neutrino energies. Hence, it may
leave an imprint on the long-term muon-neutrino signal during the later PNS
cooling phase, i.e. when the temperature drops at the PNS surface and the
neutrinospheres shift to higher densities continuously. Here, for the study of
the early PNS deleptonization phase for no more than about 10 s, the inclusion
of weak reactions with pions has a negligible impact on the neutrino emission.

Results from SN simulations which include muons and pions, the latter both
noninteracting and interacting, are presented in the Appendix A, investigating
the impact of muons and pions on the early PNS deleptonization phase. In
comparison with the reference case without muons and pions, no significant
deviation could be found for the models including muons and pions for the
first 20 s of the PNS deleptonization, independent of the treatment of pions.

The implementation of the axion losses in AGILE-BOLTZTRAN as a sink term in
the evolution equation for the internal energy, corresponding to the
conservation of energy, is discussed in Ref. Fischer:2016b . Here, the same
approach is employed and we add to the previously included axion emission from
nucleon-nucleon bremsstrahlung the updated rates as well as those from pions.
Freely streaming axions are assumed. The axion-emission rates will be
introduced and discussed in the next section.

Figure 3: Axion emissivities $Q_{\rm a}$, with respect to the rest mass
density for fixed electron and muon abundances of $Y_{e}=0.3$ and
$Y_{\mu}=10^{-4}$, respectively, as well as for three selected temperatures of
$T=10$ (left panel), 25 (middle panel), and 50 MeV (right panel), comparing
the bremsstrahlung rates $aNN$ (thin red dashed lines) and $aNN^{*}$ (thin
orange dash-dotted lines) as well as axions stemming from pions $a\pi$ (thin
green solid lines) and $a\pi^{*}$ (thick blue dashed lines).

##  III Axion-emission rates

At present, the process for axion emission in the context of core-collapse SN
has long been the nucleon-nucleon bremsstrahlung,

$N+N\longrightarrow N+N+a~{}.$ (2)

The associated rate expression for the axion emissivity is a complicated high-
dimensional momentum integral, involving the matrix elements with the yet
incompletely known axion-nucleon coupling constant $g_{aN}$ (c.f. expressions
(2.2) and (2.6) in Ref. Carenza:2019pxu, , and references therein). Previously
employed simplifications, such as the vacuum one-pion exchange potential to
describe the nuclear interactions and the assumption of massless pions, have
lead to a simplified and semianalytical expression for the axion emissivity
Brinkmann:1988vi ; Raffelt:1993ix , henceforth denoted as $aNN$. This
treatment was previously implemented in simulations of core-collapse SN
Fischer:2016b . It was found that $aNN$ leads to an overestimation of the
axion emission rate Carenza:2019pxu , which, hence, artificially enhances the
axion losses in simulations of core-collapse SN. This caveat has been overcome
with updated axion-emission rates derived in Ref. Carenza:2019pxu , henceforth
denoted as $aNN^{*}$. It includes a substantially improved description of the
nuclear interactions, i.e., nonvanishing mass for the exchange pions,
contributions from the two-pion exchange, effective in-medium nucleon masses,
and multiple nucleon scattering effects, as well as the proper phase-space
contributions including final-state Pauli blocking terms. The latter are
realized numerically via a 64-point Gauss-quadrature integration. Note that
the same nuclear model has been the foundation for the updated neutrino-pair
emission and absorption rates from $N$–$N$ bremsstrahlung Guo:2019 . Globally,
they provide a substantial reduction of the axion emissivity with increasing
density and temperature. The present paper implements the $aNN^{*}$ treatment
into the SN model AGILE-BOLTZTRAN. Together with the previously employed
simplistic axion emissivity, $aNN$, this will allow us to quantify the impact
of the medium modifications from the associated axion cooling contributions.

Another very efficient axion production mechanism is the pion Compton
scattering

$\pi^{-}+p\longrightarrow n+a~{}.$ (3)

Reactions involving $\pi^{+}$ and $\pi^{0}$ are also possible but they are
strongly suppressed compared to $\pi^{-}$ since the abundances of $\pi^{+}$
and $\pi^{0}$ are negligible, as was discussed in Sec. II. The caveat of
previous studies of axion emission from pions was the generally too low pion
abundance, also for $\pi^{-}$, to be compatible with axions from nucleon-
nucleon bremsstrahlung Keil:1996ab . However, the medium-modified suppression
of the bremsstrahlung rates leads to the revision of this argument
Carenza:2020 . The subsequent microscopic axion-emission rate allows us to
further support this argument, taking into account the treatment of
interacting pions Fore:2020 . It results not only in the substantial
enhancement of the pion abundances at high density but also in an axion
emissivity which is compatible to those of the bremsstrahlung processes. The
axion emissivities of Ref. Carenza:2020  are included into AGILE-BOLTZTRAN,
complementary to the bremsstrahlung rates.

Important here is the treatment of the interacting contributions to the pions
consistently with the underlying HS(DD2) EOS used in the simulations.
Therefore, the computation of the axion emissivity for the process (3)
implements the pion dispersion relation (1), including the momentum dependent
pion self-energies [see Eq. (1) and Fig. 1]. The axion emissivity through pion
processes is henceforth denoted as $a\pi$, in the case of vanishing pion
interactions, where we set $\Sigma_{\pi^{-}}=0$, and as $a\pi^{*}$ in the
general case of interacting pions.

Table 1: Summary of the different supernova simulations including the references to the various treatments for the calculation of the axion emissivity. Label |  |  |  |  | $N+N\rightarrow N+N+a$ |  |  |  |  | $\pi^{-}+p\rightarrow n+a$
---|---|---|---|---|---|---|---|---|---|---
Ref. run (Appendix A) |  |  |  |  | … |  |  |  |  | …
$aNN$ |  |  |  |  | Vacuum one-$\pi$ exchange, $m_{\pi}=0$ Brinkmann:1988vi ; Raffelt:1993ix ; Fischer:2016b  |  |  |  |  | …
$aNN^{*}$ |  |  |  |  | Improvements according to Ref. Carenza:2019pxu  |  |  |  |  | …
$aNN^{*}+a\pi$ |  |  |  |  | Improvements according to Ref. Carenza:2019pxu  |  |  |  |  | Rates according to Ref. Carenza:2020  with $\Sigma_{\pi}=0$
$aNN^{*}+a\pi^{*}$ |  |  |  |  | improvements according to Ref. Carenza:2019pxu  |  |  |  |  | Rates according to Ref. Carenza:2020 ,
|  |  |  |  |  |  |  |  |  | with $\Sigma_{\pi}$ according to Ref. Fore:2020

Figure 3 compares the axion emissivity, denoted as $Q_{\rm a}$, which is the
rate divided by the rest mass density and, hence, reflects the asymptotic
behavior at low density which is proportional to the number density of
targets, i.e., the abundance of $\pi^{-}$ for the axions stemming from process
(3) and $NN$ pairs for the bremsstrahlung processes (2). The latter is similar
as the low-density behavior of the neutrino pair production from nucleon-
nucleon bremsstrahlung processes (see Fig. 3 in Ref. Fischer:2016a ). Guided
by the previously revisited constraints for the axion-nucleon coupling
constant Carenza:2019pxu , we choose here $g_{ap}=1.2\times 10^{-9}$ and
$g_{an}=0$ for the calculations of the axion emissivities for the reactions
(2) and (3), since the same axion-nucleon coupling gives rise to the axion
emission for the bremsstrahlung processes and the pion conversion. We will
employ the same values for the coupling constants in the SN simulations, which
will be discussed below in Sec. IV. Here, we compare the different treatments
$aNN$ and $aNN^{*}$ for several temperatures, $T=10$ (left panel), 25 (middle
panel), and 50 MeV (right panel). From this analysis becomes evident not only
the strong temperature dependence of this axion-emission process but also the
substantial overestimation of the axion emissivity for the simplified
treatment $aNN$. It is important to note that the different axion emissivities
depend on the nuclear EOS, for which we use for the HS(DD2) EOS for all
calculations and SN simulations.

Figure 3 compares the bremsstrahlung rates with the axion emissivity from
Compton pionic processes, for $a\pi$ (thin solid green lines) and $a\pi^{*}$
(thick blue dashed lines). There are several observations: (i) At low
temperatures, there is no significant axion emission, neither from
bremsstrahlung nor stemming from pions; (ii) the overall low-density behavior
reflects the limiting dependence on the temperature of the pion abundance;
(iii) the match between $a\pi$ and $a\pi^{*}$ at low density as strongly
interacting pions divert into noninteracting pions; and (iv) the substantial
enhancement of the emissivity for $a\pi^{*}$ in comparison to $a\pi$ at a
density in excess of $\rho\simeq 10^{14}$ g cm-3, dominating the axion
emission from pions over those of bremsstrahlung, however, only at high
temperatures.

Figure 4: Radial profiles of selected quantities, fractions of electrons
($Y_{e}$), muons ($Y_{\mu}$), and pions ($Y_{\pi^{-}}$), temperature ($T$),
rest mass density ($\rho$) and the chemical potentials of electrons
($\mu_{e}$) and muons ($\mu_{\mu}$) as well as the nuclear charge chemical
potential ($\hat{\mu}$) at about 1 s postbounce during the early PNS
deleptonization phase. Figure 5: The same as Fig. 4 but at about 7 s
postbounce.

##  IV PNS deleptonization with axions

A comparison of the PNS deleptonization with muons alone and in addition with
pions can be found in the Appendix A, with an overall negligible impact on the
PNS structure and evolution of the neutrino luminosities and average energies.
Here, we report about simulations of the PNS deleptonization including axions.
We distinguish five different SN simulation setups: (i) the reference model
without axions, (ii) including axion emission from only nucleon-nucleon
bremsstrahlung at the vacuum one-pion exchange level employing the axion rates
of Ref. Fischer:2016b , denoted as $aNN$ and including all the improvements of
Ref. Carenza:2019pxu , denoted as $aNN^{*}$, then further (iii) $aNN^{*}$ plus
axions stemming from noninteracting and (iv) interacting pions, henceforth
denoted as $a\pi$ and $a\pi^{*}$, respectively. Table 1 summarizes these
different setups, including the labels used throughout the manuscript and the
references for the different treatments concerning the calculations of the
axion emissivity. Note that all simulations include muons and the associated
muonic weak interactions, according to Ref. Fischer:2020d .

Furthermore, having evolved our SN simulations consistently through all SN
phases, we were able to assess that the axion losses do not have a significant
impact prior to the SN explosion onset. Though axions begin to be produced
already slightly after core bounce from both channels (2) and (3), mainly
because the material is shock heated as the bounce shock propagates outward
reaching temperatures on the order of 10 MeV, the axion luminosities remain
negligible, on the order of at most $1$–$2\times 10^{50}$ erg s-1 for the
bremsstrahlung (depending only slightly on the treatment of the rates $aNN$ or
$aNN^{*}$) and $1$–$10\times 10^{49}$ erg s-1 for axions from pions, and are,
thus, significantly below those of the neutrinos. During the postbounce
evolution, the continuous increase of the central density and temperature, due
to the mass accretion onto the bounce shock from the still infalling material
of the stellar progenitor, leads to the continuous rise of the axion
luminosity. However, the axion luminosities never reach values compatible to
those of the neutrinos and remain below $10^{51}$ erg s-1 for both
bremsstrahlung and pion rates. Furthermore, the axion losses are dominated by
bremsstrahlung processes during the postbounce evolution. This is associated
with the still generally low abundance of pions obtained during the entire
postbounce evolution, similar to those of muons. The finding that axion losses
are negligible during the SN postbounce evolution prior to the SN explosion
onset is consistent with previous studies (cf. Ref. Fischer:2016b  and
references therein).

The situation changes after the SN explosion onset, i.e., when the SN shock is
accelerating to increasingly larger radii and when the central PNS enters the
deleptonization and later Helmholtz cooling phase.

Figure 6: Evolution of the $\bar{\nu}_{e}$ and axion
luminosities,$L_{\bar{\nu}_{e}}$ and $L_{\rm a}$ as well as the average
neutrino energies $\langle E_{\nu}\rangle$ for all neutrino species (here we
represent the heavy neutrino flavors collectively as $\nu_{x}=\nu_{\mu}$),
sampled in the comoving frame of reference at 500 km. For the axion
luminosities we distinguish axions from nucleon-nucleon bremsstrahlung (2)
($NN$) and axions stemming from pion ($\pi$) processes (3).

First, we consider the role of axion losses from bremsstrahlung processes (2).
With the axion treatment $aNN$, at about 1 s post bounce differences are small
compared to the reference simulation without axions (see Fig. 4). However,
already at about 5 s post bounce the axion emission causes the PNS
deleptonization evolution to accelerate in comparison to the reference case.
The result is an enhanced deleptonization and, hence, compression of the PNS,
featuring lower peak temperatures and a higher central density. However, the
comparison of these results with the updated axion emissivity $aNN^{*}$
clearly states a strong overestimation of the axion losses with respect to the
simplified bremsstrahlung rates $aNN$. The improved rates,$aNN^{*}$ leave a
substantially weaker impact on the PNS deleptonization, as illustrated in
Figs. 4 and 5 at two selected times during the PNS deleptonization. Note the
substantially slower deleptonization, given by a slower temperature drop, for
$aNN^{*}$ in comparison to $aNN$.

Consequently, the axion luminosity for $aNN$ is substantially higher than for
$aNN^{*}$ by nearly one order of magnitude, as illustrated in Fig. 6. The
associated accelerated PNS compression and cooling results in a shortening of
the neutrino emission timescale. In turn, the neutrino luminosities and
average energies of all flavors are reduced as soon as axion losses become
significant. For the sake of simplicity, in Fig. 6 we show only the neutrino
luminosity for $\bar{\nu}_{e}$, as those are the most relevant ones for
possible detection prospects through the inverse $\beta$ decay, e.g., at the
Super-Kamiokande detector Wu:2015 . The magnitude of the reduced neutrino
emission timescale depends on the magnitude of the axion emissivity; high
rates $aNN$ result in a more sever shortening which is substantially less
pronounced for the updated and reduced rates $aNN^{*}$. Note that the latter
neutrino fluxes and average energies are nearly indistinguishable from the
reference case without any axion losses, during the early PNS evolution.

Note further that since the nucleon-nucleon bremsstrahlung rates have a strong
temperature dependence, their impact starts to decrease as the PNS interior
cools (here, after about 5–10 s), the axion luminosity from bremsstrahlung
drops and its impact on the enhanced cooling becomes negligible. After that,
the PNS evolution is entirely governed by neutrino losses.

The inclusion of axion losses stemming from pions accelerates further the
cooling of the PNS interior. Since shortly after the SN explosion onset the
abundance of pions increases rapidly, the associated axion luminosity rises
too. This is illustrated in Figs. 4 and 6 for the non-interacting pions (solid
green lines), for which we set the pion self-energy to zero. Since the axion
luminosity from bremsstrahlung processes remains similar as for the simulation
without pions included, also the impact on the PNS structure remains
comparable (c.f. the central density at 7 s post bounce in Fig. 5). Only the
central temperature is somewhat lower, indicating the additional impact from
axion losses at the PNS interior from pions. Consequently, also the impact on
the neutrino luminosities and average energies is comparable to the setup
without axion losses stemming from non-interacting pions (see Fig. 6). The
situation changes substantially when considering interacting pions (blue
dashed lines in Figs. 4 and 5). The significantly higher abundance of
$\pi^{-}$ results in an enhanced axion luminosity, which, in turn, leads to an
accelerated PNS cooling. This is illustrated via the significantly lower
temperatures at the PNS interior already at about 7 s in Fig. 5.
Correspondingly, the faster PNS compression leads to higher central densities.
However, the faster temperature drop feeds back to a faster drop of the pion
abundance, which still exhibits a strong temperature dependence. Hence, also
the axion luminosities drop faster as illustrated in Fig. 6 (blue dashed
lines), and the impact from axion losses on the neutrino emission terminates
earlier than for the case with non-interacting pions (solid green lines).
Nevertheless, the lower neutrino luminosities obtained together with the
reduced average neutrino energies, in particular during the early PNS
deleptonization phase, at about 2–5 s, may leave interesting prospects for the
neutrino detection. In particular we find that the luminosity power law,
$L\propto t^{-1}\exp\\{-(t/\tau)^{\alpha}\\}$ Roberts:2021  with
characteristic cooling timescale $\tau$ and parameter $\alpha\simeq 2.5$,
changes substantially. This is illustrated when considering the quantity
$t\cdot L_{\nu}$ in Fig. 6, which indicates that when considering additional
losses, the neutrino losses cannot be considered constant. The latter is the
case only for the standard reference setup, during the early PNS
deleptonization phase up to about 10 s postbounce.

##  V Observable signatures

In this section we discuss the possible observable signatures of axion
emission in SN, especially focusing on the case of pionic processes. A first
study considering only bremsstrahlung processes was presented in Ref.
Payez:2014xsa , which we follow closely. In the following, we will show how
the modification of the SN neutrino signal, due to the emission of axions,
would affect the observable neutrino signal in large underground detectors.
For definiteness we focus on the future-generation Hyper-Kamiokande water-
Cherenkov detector, planned in Japan with a fiducial mass of 260 kton Hyper-
Kamiokande:2018ofw . When needed, we will also show, for comparison, results
from the current Super-Kamiokande neutrino detector, for which we assume a
mass of 32 kton.

The neutrino event rate, $N_{e}$, at Earth can be expressed symbolically as
follows Fogli:2004ff :

$N_{e}=F_{\nu}\otimes\sigma_{e}\otimes R_{e}\otimes\varepsilon~{},$ (4)

where the neutrino flux at Earth is convoluted with the interaction cross
section $\sigma_{e}$ in the detector for the production of an electron or a
positron, as well as the energy-resolution function $R_{e}$ of the detector
and the detection efficiency $\varepsilon$. The threshold of the experiment is
$E_{\rm th}=5$ MeV. For simplicity, in the following, we will neglect the
effect of energy-dependent features. Moreover, we assume $\varepsilon=1$ above
the threshold.

Concerning the detected neutrino fluxes, notably they are affected by peculiar
flavor conversion effects. Their characterization is far from being settled
due to peculiar effects in the SN interior, associated with neutrino-neutrino
interactions and with matter turbulence (for a review about SN neutrino
oscillation effects, see Ref. Mirizzi:2016 ). Furthermore, we are mostly
interested in the PNS deleptonization phase ($t\gtrsim 1$ s) after which the
spectra of all neutrino flavors become increasingly similar Fischer:2012a ,
which, in turn, reduced the impact of flavor conversions on the observable
signal. Therefore, for simplicity, we will neglect also these effects,
assuming unoscillated neutrino spectra at Earth.

Water-Cherenkov detectors are mostly sensitive to SN electron antineutrinos
through the inverse-beta decay process ${\bar{\nu}_{e}}+p\to n+e^{+}$, which
we characterize from Ref. Strumia:2003zx . In the following, we assume a
fiducial SN at a distance $d=10$ kpc.

Figure 7: Count rate in each log-time bin for $\bar{\nu}_{e}$ in the Hyper-
Kamiokande (left scale) and Super-Kamiokande detectors (right scale) due to
the inverse-beta interaction for the different scenarios discussed in Sec. IV.
The vertical error bars show the Poisson uncertainties on the counts in Hyper-
Kamiokande and the horizontal error bars show the bin widths.

In Fig. 7, we show the $\bar{\nu}_{e}$ event rate in the Hyper-Kamiokande
detector (left scale) due to the inverse-beta interaction for the different
scenarios discussed in sec. IV (see Fig. 6). For comparison, the right scale
refers to the event rate in the Super-Kamiokande detector. In particular,
following Ref. Roberts:2021 , we adopt a log scale for the time axis, choosing
ten equal-width bins per factor of 10 in time. Indeed, using a logarithmic
time axis, one can use the following algebraic relation: $dA/dt=dA/d\,{\rm
ln}t=(2.3)^{-1}dA/d\,{\rm log}_{10}t$. Therefore, Fig. 7 shows $\Delta
N=(2.3)^{-1}dN/d\,{\rm log}_{10}t$. We see that in the reference run without
axion emission (black line) the behavior of $\Delta N$ is rather constant with
time, for the interval considered here. Including axion emission from
bremsstrahlung process in the $aNN$ treatment would cause a fast drop already
at $t\gtrsim 2$ s, leading to a reduction by a factor of about 2 of the signal
at $t\sim 7$ s. The updated bremsstrahlung rate $aNN^{*}$ has a considerably
reduced impact and the complete run becomes nearly indistinguishable from the
reference one. Differences become visible only at $t\gtrsim 4$ s, reducing the
signal by a factor of about $1.2$ with respect to the reference case. If we
now include the pionic process, the drop of the neutrino event rate becomes
dramatic, being evident with respect to the reference case already at about
$t\gtrsim 1$ s. In the case of non-interacting pions $a\pi$ the final
reduction of the signal is comparable to what is obtained with the $aNN$
bremsstrahlung process. However, the drop of the signal is faster at early
times. Including interacting pions $a\pi^{*}$ features an even more dramatic
drop of the event rate at early times, resulting in a reduction by a factor of
3 at $t\sim 7$ s with respect to the reference case.

Figure 8: Positron spectra in the Hyper-Kamiokande (left scale) and Super-
Kamiokande detectors (right scale) in selected time ranges, between $1$–$2$ s
(upper panel) and $2$–$7$ s (bottom panel) for the models discussed in Sec.
IV.

Figure 8 shows the detectable positron spectra for the different models
considered, integrating the signal in two time windows, namely $t\in[1;2]~{}s$
(upper panel) and $t\in[2;7]~{}s$ (lower panel), respectively. We see that at
early times the effect of axion emission is marginal in the case of only
bremsstrahlung processes, regardless of whether the rate is corrected
($aNN^{*}$) or not ($aNN$). On the other hand, in the presence of pionic
processes, the suppression of the positron spectrum already shows a reduction
of a factor about $1.6$ at the peak with respect to the reference case, for
the interacting pion scenario $a\pi^{*}$. In the second, later time window the
effect of the suppression of the positron spectrum is more remarkable, except
in the case of $aNN^{*}$, which remains nearly indistinguishable in comparison
to the reference case. On the other hand, we find a reduction of the peak of
the spectrum of a factor of about 3 in the case of interacting pions
$a\pi^{*}$.

Figure 9: Evolution of the SN binding energy carried by neutrinos (upper
panel) and axions due to nucleon-nucleon bremsstrahlungs processes (middle
panel) and stemming form pions (lower panel), for the scenarios discussed in
Sec. IV.

A further signature of the axion emission would be the impact on the SN
binding energy, $E_{B}$, carried by neutrinos. Note that the axion emission
would lead to a reduction of the SN binding energy. This effect is illustrated
in Fig. 9, where we show the evolution of the energy carried away to infinity
by neutrinos of all flavors (top panel) and axions, distinguishing
bremsstrahlung processes (middle panels) and pionic processes (bottom panel),
for the different scenarios discussed above. Note that the early and sudden
rise of the energy carried by neutrinos is due to the contributions from the
stellar core collapse, bounce, and early postbounce phases prior to the SN
explosion onset, during which all simulations are identical and during which
we did not find any impact from the inclusion of axions.

We realize that in the case of bremsstrahlung processes the energy carried by
axions rises continuously, as indicated by the slowly increasing temperature
at the PNS center where most axions are being produced from bremsstrahlung
processes. When the temperature starts to decrease, the axion emissivity
decreases simultaneously, which is the case at around 5 s for $aNN$ and only
at around 10 s for $aNN^{*}$. The total energy carried away is similar for
$aNN$ and $aNN^{*}$; however, both of their magnitudes, $8\times 10^{52}$ and
$4\times 10^{52}$ erg, respectively, are much lower than that of neutrinos,
i.e., neutrinos dominate the losses during the PNS deleptonization phase. In
particular, the weak feedback for $aNN^{*}$ leaves nearly no difference
notable in the energy carried away by neutrinos, in comparison to the
reference case.

Contrarily, in the case of pionic processes a substantially larger fraction of
the binding energy, comparable to the energy carried by neutrinos, is carried
by axions already at early times. However, the associated rapid cooling
contributions lead to a fast reduction of the PNS interior temperature, which,
in turn, result in a faster drop of the axion emissivity. Hence, the process
of energy loss through axions stemming from pionic processes ceases much
earlier compared to the case of from bremsstrahlung processes. Consequently,
the impact on the total energy carried by neutrinos is small, a factor of
about 1.3 at 7 s. The total energy carried away by neutrinos and axions is
listed in Table 2, evaluated at about 20 s post bounce for the different
simulations. While for the axion losses the asymptotic values have been
reached within the first 7 s after bounce, for the neutrino losses we
extrapolate assuming monotonically decreasing neutrino luminosities of all
flavors.

Table 2: Total energy carried away by neutrinos, $L_{\nu}$, and axions,
$L_{\rm a}$, evaluated at a postbounce time of 20 s.

Model | $E_{\nu}$ | $E_{\rm a}^{NN}$ | $E_{\rm a}^{\pi}$
---|---|---|---
| $[10^{53}\,\,\,{\rm erg}]$ | $[10^{53}\,\,\,{\rm erg}]$ | $[10^{53}\,\,\,{\rm erg}]$
Ref. | 2.2 | … | …
$aNN$ | 1.8 | 0.88 | …
$aNN^{*}$ | 2.0 | 0.58 | …
$aNN^{*}+a\pi$ | 1.7 | 0.29 | 0.79
$aNN^{*}+a\pi^{*}$ | 1.5† | 0.06 | 0.81

†: Asymptotic values extrapolated assuming monotonically decreasing neutrino
luminosities for all flavors beyond 7 s postbounce evolution.

What is described above are all indirect signatures of the axion effects in a
SN. However, a direct detection of SN axions is perhaps also possible. To
explore this intriguing possibility, in Fig. 10 we show the axion spectrum,
integrated over the entire PNS profile, for the bremsstrahlung and pionic
processes, distinguishing non-interacting pions and strongly interacting
pions, evaluated at different postbounce times. Two features are clearly
visible in the pion-induced axion spectrum: (i) the energy threshold at the
pion mass and (ii) the energy peak at $\sim$150–200 MeV. Both features can be
understood from the kinematics of the process (3). Assuming equal mass of
protons and neutrons, it is evident that the minimal axion energy corresponds
to the minimal pion energy and, hence, with its mass. This is true also when
the pion interactions are turned on, and the pion dispersion relation is
corrected with $\Sigma_{\pi^{-}}$ [cf. Eq. (1) and Fig. 1]. The higher
threshold also moves the peak energy to higher values with respect to the
bremsstrahlung process.

Notice also the modified dispersion relation (1), with $\Sigma_{\pi^{-}}$
becoming increasingly more negative at large momenta, which, in turn, makes it
more expensive to increase the pion energy. Since the axion energy in the
reaction (3) is controlled by the pion energy (the nucleons are much heavier),
this reflects in a narrower axion spectrum when the pion interactions are
turned on ($\Sigma_{\pi^{-}}\neq 0$). This feature is clearly evident in the
spectra in Fig. 10 (cf. the green line, which shows the axion spectrum for
non-interacting pions, with the blue dashed line, for which the interactions
are turned on).

Figure 10: Normalized axion-emission spectra for nucleon-nucleon
bremsstrahlung processes (light red dash-dotted lines) and for axions stemming
from noninteracting pions (solid green lines) and interacting pions (thick
blue dashed lines), at different postbounce times 1 (top panel), at 5 (middle
panel) and at 7 s (bottom panel), corresponding to the SN simulations
discussed in Sec. IV.

It is also apparent that at early times, when the PNS is still hot and rather
dilute with large radii, the axion spectra can have extended tails up to 1
GeV. It was estimated in Ref. Carenza:2020  that for a close-by SN, i.e.,
$d\lesssim 1$ kpc, and a megaton-class water-Cherenkov underground detector,
this axion burst would be detectable through the processes $a+N\to N+\pi$,
giving an unexpected signal of ${\mathcal{O}}(10^{3})$ pions during a SN
emission. The detection of this pionic signal, together with a rapid drop of
the neutrino signal, would be strong evidence for the axion emission via
pionic processes.

##  VI Summary and conclusions

Results are reported from simulations of the PNS deleptonization phase of
core-collapse SN, launched from a progenitor of 18 $M_{\odot}$ Woosley:2002zz
, where in addition to the neutrinos, axions and the associated energy losses
are considered. The SN simulations are based on general relativistic neutrino
radiation hydrodynamics in spherical symmetry, employing six-species Boltzmann
neutrino transport. The restriction to spherical symmetry is justified for the
PNS deleptonization phase of a core-collapse SN, since the compact central
remnant PNS obeys spherical symmetry. The description of approximate
convection can be implemented in such spherically symmetric models
Roberts:2012c . All present simulations include muons as an additional degree
of freedom, and the associated muonic weak processes Guo:2020 ; Fischer:2020d
.

In a first step, updated nucleon-nucleon bremsstrahlung rates for the axion
emissivity are considered Carenza:2019pxu . These yield a substantial
reduction of the previously reported shortening of the neutrino emission
timescale Fischer:2016b , which adopted the simplified bremsstrahlung vacuum
rate expressions based on the one-pion exchange approximation. In a second
step, the inclusion of pions enables us to extend the previous studies
including axions stemming from pions. This novel channel has long been
omitted, due to the shortcoming of a generally negligible pion abundance.
Considering interacting pions, however, relaxes this argument. Then, the
inclusion of axion cooling from this channel significantly shortens the
timescale of neutrino emission as a direct feedback, in comparison to the case
when only bremsstrahlung processes are considered.

The resulting SN simulation neutrino signals are then further analyzed with
respect to the detection prospects at the future-planned Hyper-Kamiokande and
presently operating Super-Kamiokande water-Cherenkov detectors. We find that
the sudden drop of the count rate implies the presence of additional,
nonstandard energy losses, that could be related with axions originating from
pions due to Compton processes, if occurring during the early phase of
detection.

Furthermore, according to the recent study of Ref. Reed:2020gva , a
measurement of the binding energy of less than $1.5\times 10^{53}$ erg carried
by neutrinos, would be strong evidence for the presence of exotic energy
losses. It would be otherwise incompatible with SN models for different
neutron star mass and EOS. An experiment like Hyper-Kamiokande has the
potential to measure the total neutrino energy with an accuracy of few $\%$
GalloRosso:2017mdz , allowing us to diagnose with accuracy the presence of an
additional energy-loss.

Note that the uncertainty due to the elastic approximation, i.e., zero nucleon
recoil, implemented for the calculation of the emissivity of axions stemming
from pions, is unlikely to alter the findings of the present work. Even though
this gives rise to a sharp cut of the rate given by the pion rest mass, the
high-energy tail of the axion spectrum dominates the phase space integrals and
hence the axion emissivity.

It is important to note that all these results, as well as those reported in
Ref. Fischer:2020d , with a negligible impact of muons and associated muonic
weak processes during the SN bounce and early postbounce evolution, are based
on the HS(DD2) relativistic mean-field EOS. It represents a stiff nuclear
model at densities in excess of saturation density. These findings, e.g., the
negligible impact from the inclusion of muons and pions reported here, may
alter when adopting a softer nuclear EOS, i.e. when higher central densities
are encountered (see Refs. Bollig:2017 ; Bollig:2020 ; Steiner:2013  and the
left panel of Fig. 12 in the Appendix A). This yet incompletely understood
aspect remains to be explored in a more systematic fashion. Besides the
uncertainty of the EOS in terms of bulk properties, there are other aspects
that are likely to affect the signal proposed here, e.g., impact of the EOS on
the neutrino opacity, PNS convection, and fall back. All this remains to be
explored in a systematic fashion.

We note further that magnetic fields are not considered here, due to the
restriction to spherical symmetry. The magnetic field can be large at the
interior of PNS Braithwaite:2004 ; Takiwaki:2009 , which is usually associated
with magnetorotationally driven SN explosions LeBlanc:1970kg ; Winteler:2012 ,
giving rise to the formation of magnetars—an observationally known class of
highly magnetized neutron stars with surface magnetic field on the order of
$10^{15}$ G. This is still an active subject of research, in particular
regarding the development of jets in the context of core-collapse SNe
Moesta:2014 ; Moesta:2015 , especially the magnitude of the initial progenitor
magnetic field strength and the possible amplification during the SN
evolution. Furthermore, it is expected that strong magnetic fields affect the
emission of axions Kachelriess:1997 ; Maruyama:2018 . An interesting
possibility is the possible mixing of the axions with photons
RaffeltStodolsky:1988 , for which it has been demonstrated recently that the
axion mixing with the photon longitudinal mode is suppressed for light axions
due to the electron degeneracy CaputoCarenzaLucente:2021 . The only
possibility is the mixing with the photon transverse mode which requires the
photon mass to match the axion mass. Hence, in the present study we consider
light axions for which the mean free path is large compared to the PNS radius,
even in the presence of large magnetic fields. Revisiting the impact of
magnetic fields in detail would extend beyond the scope of the present study,
which identifies possibly observable signatures during the PNS
deleptonization, that originate from enhanced cooling contributions associated
with axion emission. We postpone the full study of magnetic field impact on SN
axions to a future work, which must include self-consistently the effect of
magnetic fields not only on the axion emission but also on the PNS structure
and the associated neutrino emission.

###### Acknowledgements.

T.F. acknowledges support from the Polish National Science Center (NCN) under
Grant No. 2020/37/B/ST9/00691. The work of P.C. and A.M. is partially
supported by the Italian Istituto Nazionale di Fisica Nucleare (INFN) through
the “Theoretical Astroparticle Physics” project and Research Grant No.
2017W4HA7S “NAT-NET: Neutrino and Astroparticle Theory Network” under the
program PRIN 2017 funded by the Italian Ministero dell’Università e della
Ricerca (MUR). B.F. acknowledges support from the SciDAC Grant No.
A18-0354-S002 (de-sc0018232). The work of S.R. was supported by the U.S.
Department of Energy under Grant No. DE-FG02-00ER41132. The work of M.G. was
supported by funding from a grant provided by the Fulbright U.S. Scholar
Program. M.G. thanks the Departamento de Física Teórica and the Centro de
Astropartículas y Física de Altas Energías (CAPA) of the Universidad de
Zaragoza for hospitality during the completion of this work. This work was
supported by the COST Actions No. CA16117 “ChETEC” and No. CA16214 “PHAROS”.
The SN simulations were performed at the Wroclaw Center for Scientific
Computing and Networking (WCSS) in Wroclaw (Poland).

Figure 11: Evolution of the neutrino luminosities and average energies sampled
in the comoving frame of reference at 500 km, comparing the simulations
without muons and without pions (black lines) and including muons but no pions
(gray lines), muons and noninteracting pions (blue lines) as well as muons and
interacting pions (brown lines). All simulations employ the six-species
neutrino transport scheme. Figure 12: Evolution of selected central quantities
during the PNS deleptonization of the reference simulation including muons and
associated muonic weak processes. Left panel: central temperature $T_{\rm
central}$ (thick solid line), in comparison to the maximum temperature ${\rm
max}(T)$ (thin solid line), and central density $\rho$ (dashed line). Middle
panel: electron and muon abundances, $Y_{e}$ (thin dash-dotted line) and
$Y_{\mu}$ (thin solid line), and the net neutrino abundances for $\nu_{e}$
(thick dash-dotted line) and $\nu_{\mu}$ (thick solid line). Right panel:
chemical potentials for muons, $\mu_{\mu}$ (thick solid line), electrons,
$\mu_{e}$ (thin dash-dotted line), and $\hat{\mu}$ (thin dashed line).

##  Appendix A PNS DELEPTONIZATION WITH MUONS AND PIONS

A comparison is presented here for the simulations of the PNS deleptonization
phase, comparing the following four different setups: (i) without muons and
without pions, henceforth denoted as (no $\mu$, no $\pi$), (ii) including
muons and associated muonic weak processes but without pions ($\mu$, no
$\pi$), (iii) with muons and noninteracting pions (no $\mu$ and $\pi$), and
(iv) with muons and interacting pions following Ref. Fore:2020  ($\mu$ and
$\pi^{*}$). The discussion about the implementation of the pions is given in
Sec. II, concerning both treatment for interacting and noninteracting pions.
All simulations are based on the six-species Boltzmann neutrino transport
scheme of AGILE-BOLTZTRAN Fischer:2020d . The standard set of weak reactions
used here is listed in Table I of Ref. Fischer:2020a , including the
corresponding references. Furthermore, we use the fully inelastic charged-
current electronic and muonic rates of Ref. Guo:2020 . The latter include
self-consistent contributions of weak magnetism.

The SN simulations are launched from the 18 $M_{\odot}$ progenitor from the
stellar evolution series of Ref. Woosley:2002zz . It has been evolved through
all SN phases self-consistently. Since the present focus is on the PNS
deleptonization, and since neutrino-driven SN explosions cannot be obtained in
spherical symmetry for this class of iron-core progenitors, we enhance the
charged-current electronic weak reaction rates in the gain layer, applying the
procedure introduced in Refs. Fischer:2009af ; Fischer:2012a . It results in
an explosion onset at about 350 ms postbounce, defined when the SN bounce
shock reaches a radius of about 1000 km. Note that at that stage the
enhancement of the charged-current weak rates is turned off and the standard
weak rates setup is used again. The remnant PNS contain about 1.54 $M_{\odot}$
of baryonic mass.

Figure 11 shows the resulting evolution of the neutrino luminosities, and Fig.
11 shows the average energies, sampled in the comoving frame of reference at a
radius of 500 km. The simulations of the early PNS deleptonization phase, up
to about 10 s postbounce, which include muons and associated muonic weak
reactions Guo:2020 ; Fischer:2020d , show a negligible deviation from the
reference setup without muons. This is attributed to two facts: (i) The early
PNS deleptonization phase is governed by neutrino losses at the PNS surface,
where the neutrinospheres for all flavors are located, still at moderately low
densities on the order of $<10^{13}$ g cm${-3}$, while muons are abundant at
higher densities, and (ii) the abundance of muons is generally too low to
leave a notable impact on the softening of the high-density EOS. We report
that only at later times, $>10$ s postbounce, do muons start to impact the PNS
deleptonization and cooling, by means of an accelerated PNS compression.

A second set of simulations focuses on the role of pions. The latter are in
thermal and chemical equilibrium. Note that only $\pi^{-}$ are considered;
$\pi^{+}$ and $\pi^{0}$ are strongly suppressed. Here we consider
noninteracting pions, which obey the standard thermal bosonic distribution,
and interacting pions following Ref. Fore:2020 . The latter gives rise to a
nontrivial pion dispersion relation featuring the pion self-energy. The
inclusion of noninteracting pions leaves a negligible impact on the early PNS
deleptonization phase, for the same reasons (i) and (ii) for muons
aforementioned and, moreover, because the abundance of noninteracting pions is
lower, by more than one order of magnitude, than the abundance of muons. The
situation changes only marginally when considering interacting pions, whose
abundance reaches, and even exceeds, those of the muons at the center of the
PNS. However, their overall impact on the PNS structure and evolution during
the early (about 1–10 s) PNS deleptonization remains negligible. An impact on
the neutrino fluxes and spectra could not be found.

It is interesting to note that toward later times during the PNS
deleptonization, on the order of 20–30 s, when the central temperature starts
to decrease the simulation enters the cooling phase (see the left panel in
Fig. 12). Relatedly, the central net neutrino abundance of the electron
flavors starts to decrease. Furthermore, note that only after that, do the
chemical potentials of muons approach the one of the electrons (see the right
panel in Fig. 12); i.e., weak equilibrium will be established only at post
bounce times greater than 30 s when $\mu_{\mu}\simeq\mu_{e}$.

## References

  * (1) T. Fischer, N.-U. Bastian, D. Blaschke, M. Cierniak, M. Hempel, T. Klähn, G. Martínez-Pinedo, W. G. Newton, G. Röpke, and S. Typel, The State of Matter in Simulations of Core-Collapse supernovae-Reflections and Recent Developments,  Publ. Astron. Soc. Austr. 34, e067 (2017) .
  * (2) J. Pons, S. Reddy, M. Prakash, J. Lattimer, and J. Miralles, Evolution of protoneutron stars,  Astrophys.J. 513, 780 (1999).
  * (3) H.-T. Janka, K. Langanke, A. Marek, G. Martínez-Pinedo, and B. Müller, Theory of core-collapse supernovae,  Phys. Rep. 442, 38 (2007).
  * (4) A. Mirizzi, I. Tamborra, H.-T. Janka, N. Saviano, K. Scholberg, R. Bollig, L. Hüdepohl, and S. Chakraborty, Supernova neutrinos: production, oscillations and detection,  Nuovo Cimento Rivista Serie 39, 1 (2016).
  * (5) L. Hüdepohl, B. Müller, H.-T. Janka, A. Marek, and G. G. Raffelt, Neutrino Signal of Electron-Capture Supernovae from Core Collapse to Cooling,  Phys. Rev. Lett. 104, 251101 (2010).
  * (6) T. Fischer, S. Whitehouse, A. Mezzacappa, F.-K. Thielemann, and M. Liebendörfer, Protoneutron star evolution and the neutrino driven wind in general relativistic neutrino radiation hydrodynamics simulations,  Astron. Astrophys. 517, A80 (2010).
  * (7) L. F. Roberts, A New Code for Proto-neutron Star Evolution,  Astrophys. J. 755, 126 (2012).
  * (8) L. F. Roberts, G. Shen, V. Cirigliano, J. A. Pons, S. Reddy, and S. E. Woosley, Protoneutron Star Cooling with Convection: The Effect of the Symmetry Energy,  Phys. Rev. Lett. 108, 061103 (2012).
  * (9) T. Fischer, The role of medium modifications for neutrino-pair processes from nucleon-nucleon bremsstrahlung. Impact on the protoneutron star deleptonization,  Astron. Astrophys. 593, A103 (2016).
  * (10) R. Bollig, H.-T. Janka, A. Lohs, G. Martínez-Pinedo, C. J. Horowitz, and T. Melson, Muon Creation in Supernova Matter Facilitates Neutrino-Driven Explosions,  Phys. Rev. Lett. 119, 242702 (2017).
  * (11) T. Fischer, G. Guo, G. Martínez-Pinedo, M. Liebendörfer, and A. Mezzacappa, Muonization of supernova matter,  Phys. Rev. D 102, 123001 (2020).
  * (12) A. Andronic, P. Braun-Munzinger, K. Redlich, and J. Stachel, Decoding the phase structure of QCD via particle production at high energy,  Nature (London) 561, 321 (2018).
  * (13) J. Wilson and R. Mayle, Report on the progress of supernova research by the Livermore group,  Phys.Rept. 227, 97 (1993).
  * (14) A. Akmal and V. R. Pandharipande, Spin-isospin structure and pion condensation in nucleon matter,  Phys. Rev. C 56, 2261 (1997).
  * (15) A. Akmal, V. R. Pandharipande, and D. G. Ravenhall, Equation of state of nucleon matter and neutron star structure,  Phys. Rev. C 58, 1804 (1998).
  * (16) A. S. Schneider, C. Constantinou, B. Muccioli, and M. Prakash, Akmal-Pandharipande-Ravenhall equation of state for simulations of supernovae, neutron stars, and binary mergers,  Phys. Rev. C 100, 025803 (2019).
  * (17) B. Fore and S. Reddy, Pions in hot dense matter and their astrophysical implications,  Phys. Rev. C 101, 035809 (2020).
  * (18) R. D. Peccei and H. R. Quinn, CP Conservation in the Presence of Instantons,  Phys. Rev. Lett. 38, 1440 (1977).
  * (19) R. D. Peccei and H. R. Quinn, Constraints Imposed by CP Conservation in the Presence of Instantons,  Phys. Rev. D 16, 1791 (1977).
  * (20) S. Weinberg, A New Light Boson?,  Phys. Rev. Lett. 40, 223 (1978).
  * (21) F. Wilczek, Problem of Strong $P$ and $T$ Invariance in the Presence of Instantons,  Phys. Rev. Lett. 40, 279 (1978).
  * (22) P. Di Vecchia, M. Giannotti, M. Lattanzi, and A. Lindner, Round Table on Axions and Axion-like Particles,  PoS Confinement2018, 034 (2019).
  * (23) L. Di Luzio, M. Giannotti, E. Nardi, and L. Visinelli, The landscape of QCD axion models,  Phys. Rept. 870, 1 (2020).
  * (24) P. Agrawal et al., Feebly-Interacting Particles:FIPs 2020 Workshop Report, arXiv:2102.12143 [hep-ph].
  * (25) G. Grilli di Cortona, E. Hardy, J. Pardo Vega, and G. Villadoro, The QCD axion, precisely,  JHEP 01, 034 (2016).
  * (26) S. Borsanyi et al., Calculation of the axion mass based on high-temperature lattice quantum chromodynamics,  Nature (London) 539, 69 (2016).
  * (27) I. G. Irastorza and J. Redondo, New experimental approaches in the search for axion-like particles,  Prog. Part. Nucl. Phys. 102, 89 (2018).
  * (28) P. Sikivie, Invisible Axion Search Methods,  Rev. Mod. Phys. 93, 015004 (2021).
  * (29) G. G. Raffelt, Astrophysical axion bounds,  Lect. Notes Phys. 741, 51 (2008).
  * (30) M. Giannotti, I. Irastorza, J. Redondo, and A. Ringwald, Cool WISPs for stellar cooling excesses, J. Cosmol. Astropart. Phys. 05, 057 (2016).
  * (31) M. Giannotti, I. G. Irastorza, J. Redondo, A. Ringwald, and K. Saikawa, Stellar Recipes for Axion Hunters, J. Cosmol. Astropart. Phys. 10, 010 (2017).
  * (32) P. Carenza, T. Fischer, M. Giannotti, G. Guo, G. Martínez-Pinedo, and A. Mirizzi, Improved axion emissivity from a supernova via nucleon-nucleon bremsstrahlung,  J. Cosmol. Astropart. Phys. 10, 016 (2019). [Erratum: J. Cosmol. Astropart. Phys.05, E01 (2020)].
  * (33) A. Payez, C. Evoli, T. Fischer, M. Giannotti, A. Mirizzi, and A. Ringwald, Revisiting the SN1987A gamma-ray limit on ultralight axion-like particles, J. Cosmol. Astropart. Phys. 02, 006 (2015).
  * (34) F. Calore, P. Carenza, M. Giannotti, J. Jaeckel, and A. Mirizzi, Bounds on axionlike particles from the diffuse supernova flux,  Phys. Rev. D 102, 123005 (2020).
  * (35) R. Bollig, W. DeRocco, P. W. Graham, and H.-T. Janka, Muons in Supernovae: Implications for the Axion-Muon Coupling,  Phys. Rev. Lett. 125, 051104 (2020).
  * (36) G. Raffelt and D. Seckel, Bounds on Exotic Particle Interactions from SN 1987a,  Phys. Rev. Lett. 60, 1793 (1988).
  * (37) M. S. Turner, Axions from SN 1987a,  Phys. Rev. Lett. 60, 1797 (1988).
  * (38) A. Burrows, M. S. Turner, and R. P. Brinkmann, Axions and SN 1987a,  Phys. Rev. D 39, 1020 (1989).
  * (39) A. Burrows, M. T. Ressell, and M. S. Turner, Axions and SN1987A: Axion trapping,  Phys. Rev. D 42, 3297 (1990).
  * (40) G. G. Raffelt, Astrophysical methods to constrain axions and other novel particle phenomena,  Phys. Rept. 198, 1 (1990).
  * (41) T. Fischer, S. Chakraborty, M. Giannotti, A. Mirizzi, A. Payez, and A. Ringwald, Probing axions with the neutrino signal from the next Galactic supernova,  Phys. Rev. D 94, 085012 (2016).
  * (42) M. S. Turner, Axions from SN1987A,  Phys. Rev. Lett. 60, 1797 (1988).
  * (43) G. Raffelt and D. Seckel, A selfconsistent approach to neutral current processes in supernova cores,  Phys. Rev. D 52, 1780 (1995).
  * (44) W. Keil, H. Janka, and E. Muller, Ledoux convection in protoneutron stars: A Clue to supernova nucleosynthesis?,  Astrophys. J. 473, L111 (1996).
  * (45) P. Carenza, B. Fore, M. Giannotti, A. Mirizzi, and S. Reddy, Enhanced Supernova Axion Emission and Its Implications,  Phys. Rev. Lett. 126, 071102 (2021).
  * (46) A. Mezzacappa and S. Bruenn, Type II supernovae and Boltzmann neutrino transport: The Infall phase,  Astrophys. J. 405, 637 (1993).
  * (47) A. Mezzacappa and S. Bruenn, A numerical method for solving the neutrino Boltzmann equation coupled to spherically symmetric stellar core collapse,  Astrophys. J. 405, 669 (1993).
  * (48) A. Mezzacappa and S. W. Bruenn, Stellar Core Collapse: A Boltzmann Treatment of Neutrino-Electron Scattering,  Astrophys. J. 410, 740 (1993).
  * (49) M. Liebendörfer, O. Messer, A. Mezzacappa, S. Bruenn, C. Cardall, et al., A Finite difference representation of neutrino radiation hydrodynamics for spherically symmetric general relativistic supernova simulations,  Astrophys. J. Suppl. 150, 263 (2004).
  * (50) M. Liebendörfer, S. Rosswog, and F.-K. Thielemann, An Adaptive grid, implicit code for spherically symmetric, general relativistic hydrodynamics in comoving coordinates,  Astrophys. J. Suppl. 141, 229 (2002).
  * (51) G. Guo, G. Martínez-Pinedo, A. Lohs, and T. Fischer, Charged-current muonic reactions in core-collapse supernovae,  Phys. Rev. D 102, 023037 (2020).
  * (52) T. Fischer, G. Guo, A. A. Dzhioev, G. Martínez-Pinedo, M.-R. Wu, A. Lohs, and Y.-Z. Qian, Neutrino signal from proto-neutron star evolution: Effects of opacities from charged-current-neutrino interactions and inverse neutron decay,  Phys. Rev. C 101, 025804 (2020).
  * (53) S. W. Bruenn, Stellar core collapse: Numerical model and infall epoch,  Astrophys. J. Suppl. 58, 771 (1985).
  * (54) J. M. Lattimer and F. Swesty, A Generalized equation of state for hot, dense matter,  Nucl.Phys. A535, 331 (1991).
  * (55) H. Shen, H. Toki, K. Oyamatsu, and K. Sumiyoshi, Relativistic equation of state of nuclear matter for supernova and neutron star,  Nucl.Phys. A637, 435 (1998).
  * (56) M. Hempel and J. Schaffner-Bielich, Statistical Model for a Complete Supernova Equation of State,  Nucl.Phys. A837, 210 (2010).
  * (57) M. Hempel, J. Schaffner-Bielich, S. Typel, and G. Röpke, Light clusters in nuclear matter: Excluded volume versus quantum many-body approaches,  Phys. Rev. C 84, 055804 (2011).
  * (58) A. W. Steiner, M. Hempel, and T. Fischer, Core-collapse Supernova Equations of State Based on Neutron Star Observations,  Astrophys. J. 774, 17 (2013).
  * (59) S. Typel, Relativistic model for nuclear matter and atomic nuclei with momentum-dependent self-energies,  Phys. Rev. C 71, 064301 (2005).
  * (60) S. Typel, G. Röpke, T. Klähn, D. Blaschke, and H. Wolter, Composition and thermodynamics of nuclear matter with light clusters,  Phys. Rev. C 81, 015803 (2010).
  * (61) S. Typel, M. Oertel, and T. Klähn, CompOSE CompStar online supernova equations of state harmonising the concert of nuclear physics and astrophysics compose.obspm.fr, Physics of Particles and Nuclei 46, 633 (2015).
  * (62) J. Antoniadis et. al., A Massive Pulsar in a Compact Relativistic Binary,  Science 340, 448 (2013).
  * (63) E. Fonseca et. al., The NANOGrav Nine-year Data Set: Mass and Geometric Measurements of Binary Millisecond Pulsars,  Astrophys. J. 832, 167 (2016).
  * (64) F. X. Timmes and D. Arnett, The Accuracy, Consistency, and Speed of Five Equations of State for Stellar Hydrodynamics,  Astrophys. J. Suppl. 125, 277 (1999).
  * (65) T. Fischer, G. Martínez-Pinedo, M. Hempel, and M. Liebendörfer, Neutrino spectra evolution during protoneutron star deleptonization,  Phys. Rev. D 85, 083003 (2012).
  * (66) R. P. Brinkmann and M. S. Turner, Numerical Rates for Nucleon-Nucleon Axion Bremsstrahlung,  Phys. Rev. D 38, 2338 (1988).
  * (67) G. Guo and G. Martínez-Pinedo, Chiral Effective Field Theory Description of Neutrino Nucleon-Nucleon Bremsstrahlung in Supernova Matter,  Astrophys. J. 887, 58 (2019).
  * (68) M.-R. Wu, Y.-Z. Qian, G. Martínez-Pinedo, T. Fischer, and L. Huther, Effects of neutrino oscillations on nucleosynthesis and neutrino signals for an 18 M⊙ supernova model,  Phys. Rev. D 91, 065016 (2015).
  * (69) S. W. Li, L. F. Roberts, and J. F. Beacom, Exciting prospects for detecting late-time neutrinos from core-collapse supernovae,  Phys. Rev. D 103, 023016 (2021).
  * (70) K. Abe et al. (Hyper-Kamiokande Collaboration), Hyper-Kamiokande design report, arXiv:1805.04163.
  * (71) G. L. Fogli, E. Lisi, A. Mirizzi, and D. Montanino, Probing supernova shock waves and neutrino flavor transitions in next-generation water-Cerenkov detectors, J. Cosmol. Astropart. Phys. 04, 002 (2005).
  * (72) A. Strumia and F. Vissani, Precise quasielastic neutrino/nucleon cross-section,  Phys. Lett. B 564, 42 (2003).
  * (73) S. Woosley, A. Heger, and T. Weaver, The evolution and explosion of massive stars,  Rev. Mod. Phys. 74, 1015 (2002).
  * (74) B. Reed and C. J. Horowitz, Total Energy in Supernova Neutrinos and the Tidal Deformability and Binding Energy of Neutron Stars,  Phys. Rev. D 102, 103011 (2020).
  * (75) A. Gallo Rosso, F. Vissani, and M. C. Volpe, What can we learn on supernova neutrino spectra with water Cherenkov detectors?, J. Cosmol. Astropart. Phys. 04, 040 (2018).
  * (76) J. Braithwaite and H. C. Spruit, A fossil origin for the magnetic field in A stars and white dwarfs,  Nature (London) 431, 819 (2004).
  * (77) T. Takiwaki, K. Kotake, and K. Sato, Special Relativistic Simulations of Magnetically Dominated Jets in Collapsing Massive Stars,  Astrophys. J. 691, 1360 (2009).
  * (78) J. M. LeBlanc and J. R. Wilson, A Numerical Example of the Collapse of a Rotating Magnetized Star,  Astrophys. J. 161, 541 (1970).
  * (79) C. Winteler, R. Käppeli, A. Perego, A. Arcones, N. Vasset, N. Nishimura, M. Liebendörfer, and F. K. Thielemann, Magnetorotationally Driven Supernovae as the Origin of Early Galaxy r-process Elements?,  Astrophys. J. 750, L22 (2012).
  * (80) P. Mösta, S. Richers, C. D. Ott, R. Haas, A. L. Piro, K. Boydstun, E. Abdikamalov, C. Reisswig, and E. Schnetter, Magnetorotational Core-collapse Supernovae in Three Dimensions,  Astrophys. J. 785, L29 (2014).
  * (81) P. Mösta, C. D. Ott, D. Radice, L. F. Roberts, E. Schnetter, and R. Haas, A large-scale dynamo and magnetoturbulence in rapidly rotating core-collapse supernovae,  Nature (London) 528, 376 (2015).
  * (82) M. Kachelriess, C. Wilke, and G. Wunner, Axion cyclotron emissivity of magnetized white dwarfs and neutron stars,  Phys. Rev. D 56, 1313 (1997).
  * (83) T. Maruyama, A. B. Balantekin, M.-K. Cheoun, T. Kajino, and G. J. Mathews, Axion production from Landau quantization in the strong magnetic field of magnetars,  Phys. Lett. B 779, 160 (2018).
  * (84) G. Raffelt and L. Stodolsky, Mixing of the photon with low-mass particles,  Phys. Rev. D 37, 1237 (1988).
  * (85) A. Caputo, P. Carenza, G. Lucente, E. Vitagliano, M. Giannotti, K. Kotake, T. Kuroda, and A. Mirizzi, Axion-like Particles from Hypernovae, arXiv e-prints (2021) arXiv:2104.05727, arXiv:2104.05727.


# Comment on the article ”Detailed study of the astrophysical direct capture
reaction 6Li$(p,\gamma)^{7}$Be in a potential model approach” by E. M.
Tursunov, S. A. Turakulov, and K. I. Tursunmakhatov

S. B. Dubovichenko  Fesenkov Astrophysical Institute, 050020 Almaty,
Kazakhstan  A.S. Tkachenko  Fesenkov Astrophysical Institute, 050020 Almaty,
Kazakhstan  R. Ya. Kezerashvili  New York City College of Technology, City
University of New York, Brooklyn, 11201 New York, USA  Graduate School and
University Center, City University of New York, 10016 New York, USA

###### Abstract

We explicitly present the comparison of the results for the astrophysical
$S-$factor and reaction rate for the 6Li($p,\gamma$)7Be capture process at
astrophysical energies, presented in Phys. Rev. Phys. Rev. C 105, 065806
(2022) and Phys. Rev. C 108. 065801 (2023) obtained within the famework of
potential models. We demonstrate that both potential model approaches describe
successfully the astropgisical $S-$factor and reaction rate similtaneusly and
reproduce the LUNA Collaboration [Phys. Rev. C 102, 052802(R) (2020)] results.

low and astrophysical energies, $p^{6}$Li system, thermonuclear reaction rate,
potential cluster model

††preprint: APS/123-QED

In the recent paper [1] there is a statement: ”Very recently, a detailed study
of the above 6Li($p,\gamma$)7Be direct capture process at astrophysical
energies was performed within the potential model [28]. Various versions of
the potential model have been suggested; however, none of them describe the
astrophysical $S-$factor and the reaction rates simultaneously. More
precisely, the temperature dependence of the reaction rates of the LUNA
Collaboration [25] was not reproduced within that model. Thus, the question of
whether a potential model can simultaneously describe the astrophysical
$S-$factor and the reaction rates remains open.”

Reference [28] is our article [2], and [25] is the LUNA Collaboration [3].
This statement is written without making any comparison between authors’
results and results [2] obtained within a potential model approach. The
statement motivated us to compare the astrophysical $S-$factor and reaction
rate obtained in Refs. [2] and [1]. To shed light on the contradiction between
this statement and our results [2], we present a comparison of calculations
for astrophysical factor and reaction rate reported in [2] and [1] along with
reported in literature experimental data.

Figure 1: (Color online) ($a$) Comparison of 6Li$(p,\gamma)^{7}$Be reaction
astrophysical $S-$factors. Experimental data are from $\medblacktriangleup$ –
[5], $\medblackcircle$ – [8], $\medblacksquare$ – [3], $\medblackdiamond$ –
[7], $\medsquare$ – [4], $\medtriangledown$ – [6], $\medblackstar$ – [9].
Results of calculations: red curve – our work [2]; pink curve – Ref. [1]; blue
curve – Ref. [3]; black and green curves – Ref. [11]. ($b$) Comparison of the
astrophysical reaction rates from [12, 3, 2, 1] in the range of 0.001 to 1
$T_{9}$, normalized to the NACRE rate [13]. The shaded area represents the
uncertainties from LUNA [3] (pale red band).

In Fig. 1$a$ we present the results of calculations for the astrophysical
$S-$factor in the framework of different models and all reported experimental
data, including the most recent Luna Collaboration [3] measurements. Our
astrophysical $S-$factor is given with a solid red curve, and results [1] are
presented with the pink curve. The $R-$matrix fit of the data from LUNA
Collaboration [3] is represented with the solid blue curve. A solid green
curve was obtained by Kiss et al. [11] using the weighted means of the ANCs
from the analysis of the 6Li(3He,$d$)7Be transfer reaction within the modified
two-body potential method (MTBPM). In addition, [11] contains the results for
the $S-$factor of the 6Li$(p,\gamma)^{7}$Be reaction calculated within the
MTBPM, using the values of ANCs obtained from the analysis of the experimental
astrophysical $S-$factors of the 6Li$(p,\gamma)^{7}$Be reaction [3]. This
result is given in Fig. 1$a$ with the solid black curve. One can see that both
[2] and [1] potential model approaches provide results within the error bar of
the LUNA data [3]. There is a slight difference between them, which is not
surprising, since in [1] slight modifications of the potentials from [2] were
used. The uncertainty of experimental data allows to select other options for
the potentials to describe the LUNA data [3]. Thus, the potential model [2]
reproduces LUNA Collaboration [3] data for the astrophysical $S-$factor.

Next, we compare the reaction rates from Refs. [2] and [1]. A comparison of
the reaction rates from [2] and [1] with the experimental data of the LUNA
Collaboration [3] and the results of the NACRE II compilation [12], normalized
to the NACRE rate [13] in the temperature range of $T_{9}=0.001-1$ is shown in
Fig. 1$b$. The shaded area represent the uncertainties from the LUNA
Collaboration [3]. The deviation between the adopted reaction rate obtained in
[3] and our calculation [2] in the range of $T_{9}=0.001-1$ does not exceed
5%, while the deviation of the results [1] in a temperature range of
$T_{9}=0.001-0.015$ is larger. The analytical approximations of the reaction
rates reported in [3], [2], and [1] are different, but the corresponding
reaction rates are within the uncertainties given by the LUNA Collaboration
[3].

Thus, both [2] and [1] single-channel potential model approaches reproduce not
only the absolute values of the reaction rates of Ref. [3] but also the
temperature dependence of the reaction rates. Consequently, both [2] and [1]
potential model approaches can describe simultaneously the $S-$factor and
reaction rate reported by the LUNA Collaboration [3] for the
6Li($p,\gamma$)7Be reaction. Therefore, it’s not clear how the results of [1]
differ from ours [2] and what is the novelty of [1]. In addition, it is
important to mention that the advantage of the potential model approach in [2]
relates to a realistic description of the $p^{6}$Li interaction in scattering
and bound-state channels and it does not need the use of any additional
fitting parameters for reproducing experimental data in contrast to [1].

## References

  * [1] E. M. Tursunov, S. A. Turakulov, K. I. Tursunmakhatov, Detailed study of the astrophysical direct capture reaction 6Li($p,\gamma$)7Be in a potential model approach. Phys. Rev. C 108. 065801 (2023). DOI: 10.1103/PhysRevC.108.065801
  * [2] S. B. Dubovichenko, A. S. Tkachenko, R. Ya. Kezerashvili, N. A. Burkova, and A. V. Dzhazairov-Kakhramanov 6Li($p,\gamma$)7Be reaction rate in the light of the new data of the Laboratory for Underground Nuclear Astrophysics. Phys. Rev. C 105, 065806 (2022). DOI: 10.1103/PhysRevC.105.065806
  * [3] D. Piatti, T. Chillery, R. Depalo, M. Aliotta, D. Bemmerer, A. Best, et al., Underground experimental study finds no evidence of low-energy resonance in the 6Li($p,\gamma$)7Be reaction. Phys. Rev. C 102, 052802(R) (2020).
  * [4] S. Bashkin and R. R. Carlson, Phys. Rev. 97, 1245 (1955).
  * [5] Z. E. Switkowski, J. C. Heggie, D. L. Kennedy, D. G. Sargood, F. C. Barker, and R. H. Spear, Nucl. Phys. A 331, 50 (1979).
  * [6] R. Ostojić, K. Subotić, and B. Stepančić, Nuovo Cimento A 76, 73 (1983).
  * [7] T. Paradellis, data taken from K. Arai, D. Baye, and P. Descouvemont, Nucl. Phys. A 699, 963 (2002).
  * [8] J. J. He, S. Z. Chen, C. E. Rolfs, S. W. Xu, J.Hu, X. W. Ma,M. Wiescher, et al., Phys. Lett. B 725, 287 (2013).
  * [9] A. Amar and N. Burtebayev, J. Nucl. Sci. 1, 15 (2014).
  * [10] [15] K. Arai, D. Baye, and P. Descouvemont, Nucl. Phys. A 699, 963 (2002).
  * [11] G. G. Kiss, M. La Cognata, R. Yarmukhamedov, et al., Phys. Rev. C 104, 015807 (2021).
  * [12] Y. Xu, K. Takahashi, S. Goriely, M. Arnould, M. Ohta, and H. Utsunomiya, Nucl. Phys. A 918, 61 (2013).
  * [13] C. Angulo, M. Arnould, M. Rayet, et al., Nucl. Phys. A 656, 3 (1999).