Search for Synchrotron Emission from Secondary Electrons of Proton–Proton Interactions in Galactic PeVatron Candidate HESS J1641–463

HESS J1641−463 is an unidentified gamma-ray source with a hard TeV gamma-ray spectrum, and thus it has been proposed to be a possible candidate for a cosmic-ray (CR) accelerator up to PeV energies (a PeVatron candidate). The source spatially coincides with the radio supernova remnant G338.5+0.1 but has not yet been fully explored in the X-ray band. We analyzed newly taken NuSTAR data, pointing at HESS J1641−463, with 82 ks effective exposure time. There is no apparent X-ray counterpart of HESS J1641−463, while nearby stellar cluster, Mercer 81, and stray-light X-rays are detected. Combined with the archival Chandra data, partially covering the source, we derived an upper limit of ∼6 × 10−13 erg cm−2 s−1 in 2–10 keV (∼3 × 10−13 erg cm−2 s−1 in 10–20 keV). If the gamma-ray emission is originated from the decay of π 0 mesons produced in interactions between CR protons and ambient materials, secondary electrons in the proton–proton interactions can potentially emit synchrotron photons in the X-ray band, which can be tested by our X-ray observations. Although the obtained X-ray upper limits cannot place a constraint on the primary proton spectrum, it will be possible with a future hard X-ray mission.


INTRODUCTION
The origin of high-energy cosmic rays (CRs) remains one of the most important topics in high-energy and astroparticle physics.CRs up to the knee (∼3 PeV) have been widely believed to be produced in supernova remnants (SNRs) in the Galaxy through the diffusive shock acceleration mechanism (e.g., Vink 2011).One of the issues in this SNR paradigm is that most of the SNRs detected in the gamma-ray range have spectra with exponential cutoff of at most a few tens of TeV, which failed to give a decisive answer to whether SNR shocks are capable of accelerating particles up to the knee -PeVatrons (e.g., Funk 2015).However, a remarkable discovery was made by the High Energy Stereoscopic System (H.E.S.S.) collaboration: Relatively hard spectra without any sign of cutoff were obtained from the Galactic center (HESS Collaboration et al. 2016) and the unidentified source HESS J1641−463 (Abramowski et al. 2014a), making them as the first candidates for PeVatron sources.Recently, (sub-) PeV gamma-ray astronomy has begun by extensive air shower arrays (EASA), such as HAWC, Tibet ASγ, and LHAASO.These observations are accelerating the study of PeVatrons.
In this paper, we would like to revive the interest in the primary PeVatron candidate, HESS J1641−4631 (Abramowski et al. 2014a).Since it is close to another gamma-ray bright source HESS J1640−465 (Abramowski et al. 2014b), the confusion had been preventing the detection of HESS J1641−463 until recently.Thanks to the improved H.E.S.S. point spread function, HESS J1641−463 was discovered in the tail of HESS J1640−465 with a slight spatial extension up to ∼3 ′ (Abramowski et al. 2014a; H. E. S. S. Collaboration et al. 2018).HESS J1641−463 spatially coincides with the SNR G338.5+0.1, which has roughly a circular morphology with a diameter of 5 ′ -9 ′ in the radio band (Whiteoak & Green 1996), a distance of ≈11 kpc (Kothes & Dougherty 2007), and an estimated age of 1.1-17 kyr (Abramowski et al. 2014a).HESS J1641−463 has one of the hardest gamma-ray spectra with Γ = 2.07 ± 0.11 stat ± 0.20 sys in 0.64-100 TeV without any indication of a cutoff, implying that the cutoff energy of accelerated particles is higher than 100 TeV at 99% C.L. (Abramowski et al. 2014a).Thus, HESS J1641−463 was proposed to be a primary candidate for a PeVatron.
The origin of gamma-ray emission in the HESS J1641−463 and HESS J1640−465 system remains unsettled.In the GeV gamma-ray band, there exist counterparts of HESS J1641−463 that has a point-like morphology and log-parabola spectrum and of HESS J1640−465 that is characterized by ∼3.′ 2 spatial extension and power-law spectrum with Γ = 1.7 (Mares et al. 2021).As discussed by Mares et al. (2021), there are a couple of scenarios for this complicated system, such as (1) hadronic gamma-ray radiation (SNR or illuminated molecular cloud), (2) pulsar, or (3) PWN for the emission of HESS J1641−463.None of them, however, are conclusive for HESS J1641−463, while HESS J1640−465 is likely originated from a PWN (e.g., Abdelmaguid et al. 2023).
The spectrum without a cutoff generally prefers a hadronic gamma-ray production over leptonic emission (i.e., inverse Compton (IC) scattering of the ambient photons by accelereated electrons) because the Klein-Nishina effect significantly reduces the gamma-ray flux and inevitably makes a spectral cutoff at ≳10 TeV.The hadronic gamma-ray emission is attributed to decay of neutral pions produced in interactions of CR protons with ambient protons.The same proton-proton interactions produce charged pions, which decay into muons and electrons/positrons in turn.These electrons and positrons are often referred to as secondary electrons.If the energy of the parent proton is sufficiently high, the secondary electrons are also energetic enough to emit X-ray photons via synchrotron radiation, which can be potentially measured with X-ray observations.Although there are some dedicated studies to predict the secondary synchrotron component (e.g., Aharonian 2004Aharonian , 2013;;Huang et al. 2020), the emission has not been detected to date.In this paper, we aim at searching for the secondary synchrotron emission in the X-ray data of HESS J1641−463 and investigating how the synchrotron radiation from the secondary electrons depends on physical parameters.If detected, we show that the X-ray observation would be a new diagnostic tool to determine the maximum energy of accelerated protons, while only lower limit (>100 TeV) is obtained by the current gamma-ray observation.
We present X-ray observations of HESS J1641−463 in Section 2. Section 3 describes modeling of the broadband spectral energy distribution (SED), including calculation of the synchrotron radiation from secondary electrons.The discussion and conclusion are presented in Sections 4 and 5, respectively.

X-RAY ANALYSIS
This section describes X-ray observations and analysis of the gamma-ray source HESS J1641−463.In the soft X-ray band (≲ 10 keV), on-axis observations pointing at this source have not been performed.There are archival observations with Chandra (Obs IDs 11008 and 12508) and XMM-Newton (Obs ID 0302560201), which partially cover the gamma-ray extent of HESS J1641−463.We made use of these Chandra data in this paper.2In addition to the soft X-ray data, we also used newly retrieved hard X-ray data by NuSTAR (Obs ID 40401004002).Table 1 summarizes information on these X-ray data analyzed in this paper.
The Chandra data were processed using CALDB version 4.9.6 in CIAO version 4.14.The NuSTAR data were calibrated and screened by using nupipeline of NuSTAR Data Analysis Software (NuSTARDAS version 2.1.1 with CALDB version 20220118) included in HEAsoft version 6.29.We used the strictest mode (i.e., SAAMODE=STRICT and TENTACLE=YES cut) for screening the NuSTAR data.We present the analysis and results of X-ray images and spectra in Sections 2.1 and 2.2, respectively.

Image
We used fluximage to produce Chandra flux images.For NuSTAR images, we generated count maps and exposure maps by Xselect and nuexpomap, respectively, and then produced flux images by dividing the count maps by the exposure maps.Background was not taken into account for the X-ray flux images.
Figure 1 shows flux images with Chandra in the 0.5-7 keV band, and Figure 2 shows flux images with NuSTAR in the 3-20 and 20-79 keV bands.There is no apparent counterpart of HESS J1641−463 (see Figures 1 and 2).Another source, stellar cluster candidate Mercer 81, clearly appeared in both the Chandra and NuSTAR images.The Chandra images contain the other point-like sources, reported by Abramowski et al. (2014a).
NuSTAR data are sometimes contaminated by nearby bright sources located outside the field of view (FoV), known as stray-light contamination.As shown in Figure 2, unfortunately our observation was heavily contaminated by straylight sources, since the target region is crowded in the Galactic plane with a low latitude of b = 0. • 1.The FPMA image included stray light from AX J1631.9−4752 and Big Dipper, while the FPMB image was contaminated by stray light from 2MASS J16480656−4512068 (see the panels (a) and (d) in Figure 2).Because the TeV gamma-ray emitting region was entirely overwrapped with the stray-light component in the FPMB image, it would be better to use local background in the same FoV for estimation of the background rather than using Nuskybgd (Wik et al. 2014).It should be noted that we checked consistency between the two background methods.For the spectral analysis in Section 2.2, the stray-light region of AX J1631.9−4752 in the FPMA data was excluded, and we extracted the background in the same region contaminated by stray-light 2MASS J16480656−4512068 in the FPMB data.

Spectrum
This section presents the spectral analysis of the Chandra and NuSTAR data.We used specextract and nuproducts with 'extended = yes', respectively, to produce Chandra and NuSTAR spectra.Xspec (Arnaud 1996, version 12.12.0)was used in spectral fitting.Errors are shown in 1σ confidence level.Source and background regions are shown in Figure 1 and Figure 2.
Besides HESS J1641−463, we first conducted spectral analysis of Mercer 81, which was clearly detected by both Chandra and NuSTAR, and thus it is helpful to check consistency between the Chandra and NuSTAR analyses.For the Chandra data, the source spectrum was extracted from a circle with a center at (RA, Dec) = (16 h 40 m 29.s 4, −46 • 23 ′ 32.′′ 1) and a radius of 0. ′ 9, while the background region was taken from an annulus with 1 ′ -2 ′ .We used only the data with ObsID 11008, which fully covered the entire emission of Mercer 81.For the NuSTAR data, the source region was the same as the Chandra data, while the background was extracted from an annulus with 1. ′ 5-2.′ 5.The obtained spectra of Mercer 81 are shown in Figure 3.We fit the spectra with two models; nonthermal (i.e., power law plus Gaussian that accounts for ∼6 keV line emission) and thermal (apec) models (Table 2).For both models, absorption was taken into account by using the TBabs model (Wilms et al. 2000).When fitting with the nonthermal model, we obtained the best-fit parameters of N H = (17 ± 1) × 10 22 cm −2 , Γ = 3.9 ± 0.2, and Gaussian line with center energy of 6.66 ± 0.03 keV, σ = (9.3± 3.7) × 10 −2 keV, and normalization of (3.8 ± 0.6) × 10 −6 photons cm −2 s −1 .The fitting result with 3.9 ±0.2 2.1±0.2F2−10 (10 −13 erg cm −2 s −1 ) 2.9 +0.2 −0.5 2.9 +0.2 −0.4 Gaussian line energy (keV) 6.66 ± 0.03 -Gaussian line width (keV) the thermal model is as follows: N H = (14 ± 1) × 10 22 cm −2 , kT = 2.1 ± 0.2 keV, and the abundance of 1.0 ± 0.2.The thermal model gave a smaller value of χ 2 : χ 2 was 267 (d.o.f.229) and 258 (d.o.f.231) for the nonthermal and thermal models, respectively.We found that the normalization of the Chandra and NuSTAR spectra was in good agreement, as shown in Figure 3, thus our estimation of background (local background method) for the NuSTAR data is reliable to some extent.
The source region of HESS J1641−463 was set to a circle with a center at (RA, Dec) = (16 h 41 m 2. s 1, −46 • 18 ′ 13. ′′ 0) and a radius of 3 ′ , which is comparable to the TeV gamma-ray extent (Abramowski et al. 2014a).For the Chandra data, the background was extracted from an annulus region of 3. ′ 3-4.′ 5, centered at the aforementioned position of HESS J1641−463.For the NuSTAR data, the background was taken from an annulus region of 3. ′ 5-4.′ 1 for both FPMA and FPMB.We show the spectra of HESS J1641−463 in Figure 4.In the case of a partial coverage, we scaled the obtained normalization to the TeV gamma-ray size of r = 3 ′ on the assumption of a uniform distribution.This normalization scaling was applied to the following data: the Chandra data with ObsID 11008, where the region outside  the FoV (46% of the full coverage) was masked, and the NuSTAR FPMA data, where the region contaminated by stray-light (36%) was excluded.Since there was no apparent X-ray signal from the HESS J1641−463 region, we assumed a fitting model of an absorbed power law and fixed a column density to the Galactic H I value of N H = 2 × 10 22 cm −2 (Dickey & Lockman 1990;Kalberla et al. 2005;Mares et al. 2021) and a photon index to Γ = 2. On this assumption, 2σ upper limits were derived (Table 3).We checked how the choice of the N H value depends the flux: the normalization of the NuSTAR spectrum, which is not largely affected by the absorption, remains within 50% in a range of N H from 0.7 × 10 22 cm −2 to 20 ×10 22 cm −2 , while that of Chandra is consistent within a factor of 3 in the same range.The flux with Chandra ObsID 12508 in the 2-10 keV energy band was 6.3 × 10 −13 erg cm −2 s −1 , which is reconciled with that with NuSTAR of (6.2-7.1)×10−13 erg cm −2 s −1 .The scaled flux obtained with Chandra ObsID 11008 is 2.8 × 10 −13 erg cm −2 s −1 , which is in agreement with 2.31 × 10 −13 erg cm −2 s −1 in the literature (Mares et al. 2021).It should be noted that these values are all consistent within the systematic uncertainty caused by the choice of N H .
We evaluated X-ray detection significance of HESS J1641−463 as follows.Background-subtracted count rates are (1.26 ± 0.48) × 10 −2 cts s −1 , (1.87 ± 1.37) × 10 −3 cts s −1 , and (6.19 ± 2.16) × 10 −3 cts s −1 for Chandra, NuSTAR FPMA, and FPMB, respectively.This implies 1.4-2.9σsignificance above statistical fluctuation.We also found that χ 2 reduced from 501 (d.o.f.435) to 457 (d.o.f.434) when fitting the Chandra + NuSTAR spectra with the normalization being fixed to 0 and free, respectively.This gave us F-test statistic value of 15.7 and probability of 1.38 × 10 −4 , corresponding to a 3.5σ detection.Because the detection significance was as low as ∼ 3σ and there was no apparent counterpart seen in the images, we conservatively adopt the upper limits in this paper.

Note-
Size is in units of arcmin 2 .Flux is 2σ upper limit and in units of 10 −13 erg cm −2 s −1 .

BROADBAND SED MODELING
Figure 5 shows the broadband SED of HESS J1641−463, including the X-ray upper limits in this paper as well as the radio and gamma-ray flux in the literature.SNR G338.5+0.1, which is a counterpart candidate of HESS J1641−463, was identified as a radio SNR by Shaver & Goss (1970) in an investigation of extended radio sources with the Parkes telescope at 5000 MHz and the Molonglo radio telescope at 408 MHz.However, because of the complex region in the vicinity, the association is uncertain, and the SNR itself has an undefined shape, preventing precise determination of its radio flux and size.It was reported that G338.5+0.1 has a diameter of 5 ′ -9 ′ and ∼12 Jy at 1 GHz (Green 2019;Whiteoak & Green 1996).Although the radio flux of the entire SNR is roughly 12 Jy at 1 GHz, Mares et al. (2021) calculated the radio emission directly associated with the gamma-ray emitting region of HESS J1641−463 (i.e., the size in FWHM of 3 ′ ) based on the Southern Galactic Plane Survey (SGPS) data (Haverkorn et al. 2006), resulting in an upper limit of ∼2 Jy at 1.4 GHz.The former value is used in this paper, otherwise mentioned.Note that the usage of the radio flux affects only a leptonic modeling, not a hadronic component that is our main point in this paper.
GeV gamma rays are detected with Fermi-LAT in the region of HESS J1641−463 (Mares et al. 2021).The GeV gamma-ray source is point-like and has a spectrum with a curvature, implying that it would likely be a pulsar.Since there is no counterpart of the source at the other wavelength, we do not include the GeV gamma-ray emission in our models.
In the H.E.S.S. Galactic plane survey (HGPS), the obtained spectrum of HESS J1641−463 was steeper than the one from the analysis dedicated for HESS J1641−463 by Abramowski et al. (2014a): the photon index was 2.47±0.11(H.E. S. S. Collaboration et al. 2018).See Figure 5 for comparison of these two spectra.In this paper, we mainly use the spectrum of Abramowski et al. (2014a).
The upper limits of the X-ray flux obtained in Section 2, combined with the gamma-ray and radio flux in the literature, enable us to model the broadband SED (Figure 5) in order to interpret the origin of the very-high-energy gamma-ray emission.In the case of a leptonic origin, IC scattering of high-energy electrons is accounted for the gammaray radiation, and synchrotron radiation from the same electrons produce radio and X-ray photons.In a hadronic scenario, decays of π 0 produced in proton-proton interactions are applied to the gamma-ray data.Other channels in the proton-proton interactions yield secondary electrons and positrons, which are potentially energetic enough to produce synchrotron radiation in the X-ray energy band.We summarize derivation of the secondary electrons in Section 3.1 and apply it to the observation in Section 3.2.

Secondary electrons in proton-proton interaction
When CR protons collide with protons in ambient medium (proton-proton interaction), pions, π 0 and π ± , are produced.They decay into the following products: (1) π 0 decay yields gamma-ray photons (Equation 1).Electrons and positrons (e ± ) are also produced from π ± through Equations 2 and 3, and these electrons and positrons are referred to as secondary electrons hereafter.The secondary electrons have ∼10% energy of the primary protons.The cross section of the π 0 -decay photons, , is well formulated in the literature (e.g., Kelner et al. 2006;Kafexhiu et al. 2014), and Geant4-based cross section in Kafexhiu et al. ( 2014) is adopted for in this paper.We make use of a library of LibppGam3 for calculating the gamma-ray cross section.Kelner et al. (2006) derived the spectral ratio of the products of proton-proton interactions, F i (E i /E p , E p ), where the subscript i indicates gamma rays, electrons, or neutrinos and E p is the energy of the primary proton.The spectrum for the secondary electrons, F e (E e /E p , E p ), is given by Equation ( 62) in Kelner et al. (2006).Thus, the cross section of the secondary electrons is described by Figure 6 compares the cross section of the secondary electrons by Equation 4with that of the gamma rays in Kafexhiu et al. (2014).Following Equation (25) in Kafexhiu et al. (2014), the energy distribution of the secondary electrons can be expressed as Here, n is the ambient density, c the speed of light, and dNp dEp the spectrum of protons.Technically, we made use of a library of GAMERA4 to solve Equation 6.We multiplied the cross section by the so-called nuclear enhancement factor ϵ(E p ) to account for proton-nucleus interactions (Kafexhiu et al. 2014).In this paper, the nuclear enhancement factor was taken from Sihver et al. (1993).
The secondary electrons, injected with Equation 5, suffer radiation cooling.Therefore, we took the cooling effect into consideration in the following equation of time evolution: where N e is the spectrum of electrons, and b(E) indicates the energy loss rate.Although b(E) consists of synchrotron radiation, IC scattering, Bremsstrahlung, and Coulomb collision, energy loss due to synchrotron radiation is dominant for >1 GeV electrons for a parameter range we consider in the paper.In the case of synchrotron radiation, the break energy is derived by equating the radiation timescale to the age of the source (T ), where B indicates the strength of the magnetic field.Figure 7 shows an example of the model using parameters summarized in Table 4.We assumed a power-law spectrum with an exponential cutoff for primary protons; where K is the normalization, E ref (=1 TeV) is the reference energy, s p is the spectral index, and E c is the cutoff energy.
In Figure 7 (right), we show an example of our radiation models, overlaid with the observation data of the TeV gamma ray (Abramowski et al. 2014a), our X-ray upper limits, and the radio flux (Green 2019;Mares et al. 2021).We calculated radiation spectra of π 0 -decay gamma rays from the primary protons and synchrotron X-rays from the secondary electrons by using Naima5 (Zabalza 2015).It should be noted that the IC component from the secondary electrons is smaller than the pion-decay emission by two orders of magnitude.

Application to HESS J1641-463
We adopt the hadronic model to the SED of HESS J1641−463.The hadronic interpretation is supported by the presence of a dense molecular cloud with n ∼ 100 cm −3 observed by NANTEN (Abramowski et al. 2014a).Computation of our model requires parameters of n, d (distance), T, B, s p , E c,p , and W p .Although these parameters are poorly constrained, we defined a baseline model as (n, d, T, B, s p , E c,p ) = (100 cm −3 , 11 kpc, 5 kyr, 300 µG, 2.0, 1 PeV), listed in Table 4. From the literature, we adopted the ambient density of 100 cm −3 , the distance of 11 kpc, and the age of 5 kyr (Dickey & Lockman 1990;Kalberla et al. 2005;Abramowski et al. 2014a).We set B = 300 µG, assuming that the gamma rays are produced in a molecular cloud where the magnetic field could be from tens to hundreds of µG (see Figure 8 for dependency on these values).s p = 2 and E c,p = 1 PeV are obtained to reproduce the TeV gamma-ray spectrum (Abramowski et al. 2014a).First, we fit the TeV gamma-ray data with the π 0 -decay model to determine the normalization of protons, K in Equation 8, which can be transferred to the total energy of protons (W p ). W p is 1.5×10 48 erg for the baseline model.Then, we calculated the secondary synchrotron component with the parameters in Table 4.The result of the baseline model is shown in Figure 7.  4. Left: the injection spectrum of the protons is shown with the blue line, the injection spectrum of the secondary electrons (Equation 5 multiplied by the age T ) is indicated by the pink line, and the spectrum of the cooled secondary electrons is illustrated by the red line.Right: See Figure 5 for the observational data.In order to explore how the secondary synchrotron emission depends on the parameters, we performed the aforementioned computations with the parameters different from the baseline model.We changed T in a range of 1-11 kyr, B in 100-2000 µG, s p in 1.5-2.5, and E c,p in 100-3000 TeV.The results are shown in Figure 8.The age, T , has an effect only on the break energy of the electrons (i.e., the energy where the age and the time scale of sycnhrotron cooling are balanced in Equation 7) and does not largely affect the model in the X-ray band (Figure 8a).The magnetic field, B, slightly modifies the synchrotron model (Figure 8b).However, the choice of B does not largely change the normalization of the synchrotron radiation, since large B makes the energy loss large and the electron spectrum steep, but, at the same time, makes the synchrotron power high, canceling out the steep spectrum of electrons.The spectral shape of the primary protons, s p and E c,p , largely change the X-ray spectrum (Figure 8c and d).It is worth noting that n and W p (or K) are degenerated, and they do not affect the X-ray spectrum.The hadronic gamma rays trace a molecular cloud (or dense gas) illuminated by cosmic rays.There are two cases proposed: (1) SNR-cloud interaction -particles are accelerated at the SNR shock and interacts with the molecular clouds in the immediate vicinity of the SNR and (2) CR-cloud interaction -accelerated particles have escaped the acceleration site, diffusively propagated, and hit nearby molecular clouds.In the former case, the higher-energy CRs can penetrate into clouds deeply, resulting in the relatively hard CR spectrum (e.g., Inoue et al. 2012;Gabici & Aharonian 2014;Celli et al. 2019;Higurashi et al. 2020).The SNR shock at such dense environment would produce thermal X-rays.In the CR-cloud scenario, the spectrum of the escaping particles generally depends on assumptions, such as a time from the escape and a distance from the source and the cloud (e.g., Aharonian & Atoyan 1996;Gabici & Aharonian 2007;Ohira et al. 2010).Particularly, a self-consistent calculation for the HESS J1641−463/HESS J1604−465 system was conducted in detail by Tang et al. (2015), showing that the distance affects the low-energy part of the spectrum of the escaping particles.
Since all of the models are lower than the X-ray upper limit obtained in this paper, we cannot place tight constraints by the X-ray observations at this time.As illustrated in Figure 8, we found that the synchrotron radiation from the secondary electrons is variable in the X-ray energy band, depending strongly on s p and E c,p and slightly on B. On the other words, a precise measurement of the secondary component can constrain the proton spectrum, which will be possible with the future missions (Section 4.2).
The hard TeV gamma-ray spectrum indicates particle acceleration up to the multi-TeV range, E c,p >100 TeV (Abramowski et al. 2014a, and Figure 8).Given the SNR age of ∼5 kyr, the shock velocity should slow down to tens of km s −1 , making it unlikely to produce such energetic particles at present.Thus, it implies the acceleration had been achieved in the past in the SNR-cloud scenario.The CR-cloud scenario, on the other hand, does not require such on-site and present acceleration.In both cases, these high energy CRs can be hardly confined in the SNR shell and escape the shock first.
Our modeling revealed that the total proton energy, W p (E > 1 TeV), is roughly ∼ 10 48 erg on the assumption of n = 100 cm −3 .Typically, an SN explosion releases an energy of E SN ∼ 10 51 erg, and a fraction of the energy is converted into acceleration of CRs (e.g., E CR ∼ 10 50 erg, assuming ∼10% for the conversion factor).W p ∼ 10 48 erg in HESS J1641−463 is smaller than the typical value for SNRs.In the SNR-cloud scenario, this may imply the conversion fraction from E SN to E CR is lower for this source.If this SNR is indeed interacting with a molecular cloud, parts of it may be radiative, in which case a lot of energy could have been radiated away already.Alternatively, the small proton energy might support the escape scenario.
It should be noted that the lack of bright thermal X-ray emission, which is expected for such SNR with dense gas (n ∼ 100 cm −3 ), may disfavor the SNR scenario.However, the SNR with the age of ∼5 kyr might be too old to heat up the gas to keV because of the decelerated shock speed.On the other hand, in the escape (CR-cloud) scenario, the accelerated particles have already left the source, which naturally explains the absence of the thermal emission.The existence of the thermal emission can be verified by deep observation in the soft X-ray band in the future.
Here, we note the result when we use the steep spectrum obtained in HGPS (H.E. S. S. Collaboration et al. 2018).With n being fixed to 100 cm −3 , we derived s p = 2.44 ± 0.13, E p >200 TeV, and W p = (1.4 ± 0.17) × 10 48 erg.Using these best-fit parameters of the protons and the others of the baseline model, the synchrotron radiation from the secondary electrons became about one order of magnitude lower than the model optimized for the harder one in Abramowski et al. (2014a).

Other scenarios
Leptonic scenario -The leptonic-dominated model consists of IC scattering and synchrotron radiation from primary CR electrons.Assuming a cutoff power-law model of electrons (Equation 8) and typical values for the photon fields of CMB, NIR, and FIR (i.e., temperatures of 2.72 K, 30 K, and 3000 K, and energy densities of 0.261, 0.5, and 1 eV cm −3 , respectively (Zabalza 2015)), fitting the multiwavelength SED results in the spectral index (s e ) of 2.9 ± 0.1, the magnetic field of 2.5 +0.8 −0.5 µG, and the poorly determined cutoff energy of E c,e > 100 TeV.Adopting the upper limit of the radio flux derived from the SGPS data (Haverkorn et al. 2006;Mares et al. 2021), s e and E c,e cannot be well constrained while B should be less than 2 µG to surpress the synchrotron radiation.
Because the TeV gamma-ray flux is compatible with the X-ray upper limits, the magnetic field should be less than ∼3 µG in the leptonic case.This B-field value is comparable to the typical value of the interstellar magnetic field, implying that no shock compression and/or amplification are required and that the leptonic case seems very unlikely.Such a small value of the magnetic field may imply the TeV gamma-ray origin is a pulsar halo.The possible presence of the pulsar is also supported from the curved spectrum from the GeV gamma-ray point-like source, although there are no counterparts of the pulsar at the other wavelengths.Furthermore, the electron total energy (W e for E >1 TeV) is (3.6 ± 0.6) × 10 47 erg.If the electron to proton ratio, K ep , is 0.01, the proton energy would be ∼ 10 49 erg, which is larger then the hadronic case.For the leptonic component to be dominant, the ambient density should be <10 cm −3 , which conflicts with the observed value of 100 cm −3 .In conclusion, the leptonic-dominated case is disfavored, although it is not excluded.
Mixed scenario (one-zone) -Here, we consider adding primary electrons that have the same spectral model (Equation 8) with the protons, as one-zone model.Radiation cooling of electrons is taken into account.In the ranges of E p = 100-1000 TeV, B = 10-1000 µG, and s =2-2.5, no parameter set was found to be in agreements with all the observations: (1) K ep should be < 0.01 for not exceeding X-ray ULs, but then the synchrotron flux becomes much lower than the radio data, and (2) if we determine the normalization of the electrons to reproduce the radio flux, the model exceeds the X-ray upper limits or the Fermi-LAT data.E c,p (proton maximum energy) should be ≳100 TeV in order to reproduce the hard spectrum in the TeV gamma-ray energy band, resulting in higher cutoff energy in the spectrum of the primary synchrotron radiation, even though cooling is taken into consideration.In conclusion, one-zone model is unlikely acceptable.Adopting the upper limit of the radio flux derived from the SGPS data (Haverkorn et al. 2006;Mares et al. 2021) does not largely change the aforementioned result.
Mixed scenario (two-zone) -In addition to the protons and the secondary electrons fixed to the baseline model (Table 4), we add primary electrons with different spectral shape (i.e., the primary electron index of s e ̸ = s p , and the primary electron cutoff energy of E c,e ̸ = E c,p ).The magnetic field for the primary electrons to emit synchrotron photons (labeled as B e hereafter) is assumed to be independent of that of the secondary electrons (i.e., B e ̸ = B).After the synchrotron radiation from the primary electron is normalized to the radio data, we searched for parameters (E c,e and B e ) with which the total flux of all the radiation components does not exceed the observations of the X-ray upper limits and the GeV emission from the point-like source (Mares et al. 2021).Assuming s e = 2.5, B e should be larger than 4 µG to suppress the IC component and E c,e should be less than ∼10 TeV not to exceed the X-ray upper limit.If we adopt for K ep =0.01, E c,e and B e are constrained to be ≲4 TeV and 20-40 µG, respectively.

Future prospects
In the hadronic scenario, our models (Figure 8) suggest that the flux of the secondary synchrotron component is > 1.8 × 10 −14 erg cm −2 s −1 in 2-10 keV and > 4.9 × 10 −15 erg cm −2 s −1 in 10-50 keV for E c,p > 100 TeV.This flux level is sufficiently detectable with X-ray telescopes.In order to avoid a contamination of the primary synchrotron radiation which generally has a cutoff at the soft X-ray band, hard X-ray detectors would be needed.Proposed missions, such as HEX-P6 (Madsen et al. 2019(Madsen et al. , 2023;;Mori et al. 2023;Reynolds et al. 2023) or FORCE (Mori et al. 2016) will be able to measure the emission and determine the proton maximum energy by the X-ray observation for the first time.Long-exposure observations with the existing X-ray telescopes would also be helpful to search for the secondary synchrotoron radiation, and we are planning an XMM-Newton observation.
Besides the synchrotron emission from the secondary electrons, neutrino is also a smoking gun to evidence hadronic acceleration.Although the neutrino amount from this single source would be much smaller than the detection threshold of the existing neutrino observatories such as IceCube, neutrinos accumulated by the hadronic PeVatron candidates may have a contribution to the neutrino diffuse flux detected from the Galactic plane (ICECUBE COLLABORATION 2023) as unresolved sources.

CONCLUSION
We analyzed newly taken NuSTAR data and the archival Chandra data of the unidentified TeV gamma-ray source, HESS J1641−463, to search for the synchrotron radiation from the secondary electrons in the hadronic interpretation.No X-ray counterpart was detected in the gamma-ray extension, and we obtained 2σ upper limits of ∼ 6 × 10 −13 erg cm −2 s −1 in the 2-10 keV band and ∼ 3 × 10 −13 erg cm −2 s −1 in the 10-20 keV band.Combined with the published radio and gamma-ray spectra, we performed modeling of the broadband SED.Although the X-ray upper limit could not constrain the hadronic model, we found the flux of the secondary synchrotron component is larger than ∼ 10 −14 erg cm −2 s −1 , which would be detectable with a deep X-ray observation with existing telescopes or with future mission.

Figure 2 .
Figure 2. Left (a and d): 3-20 keV count maps of NuSTAR, shown with stray-light sources, the contour of H.E.S.S. significance map (Abramowski et al. 2014a), and a yellow cross of the best-fit position of the GeV gamma-ray source (Mares et al. 2021).Middle (b and e): Flux images with NuSTAR in 3-20 keV.Right (c and f): Flux images with NuSTAR in 20-79 keV.The FPMA and FPMB images are respectively shown in the top and bottom panels.Source and background regions for spectral analysis are indicated by green solid and green dashed lines, respectively.

Figure 3 .Figure 4 .
Figure 3. Chandra (shown in green) and NuSTAR (FPMA in black and FPMB in red) spectra of Mercer 81, fitted with the thermal model (left) and the nonthermal model (right).The fitting results are shown in Table2.

Figure 8 .
Figure 8. Models of Table 4, but changing the age T (a), the magnetic field B (b), the spectral index sp (c), and the maximum energy Ec,p (d).

Table 1 .
Dataset of X-ray observations

Table 2 .
Fitting results of Mercer 81

Table 4 .
Baseline model parameter