# Transient Gamma Rays from the 2021 Outburst of the Recurrent Nova RS Ophiuchi: The Effect of Gamma-Ray Absorption

Vo Hong Minh Phan, Pierre Cristofari, Enrico Peretti, Vincent Tatischeff, Andrea Ciardi

► To cite this version:

Vo Hong Minh Phan, Pierre Cristofari, Enrico Peretti, Vincent Tatischeff, Andrea Ciardi. Transient Gamma Rays from the 2021 Outburst of the Recurrent Nova RS Ophiuchi: The Effect of Gamma-Ray Absorption. The Astrophysical Journal Letters, 2025, 990 (1), pp.L17. 10.3847/2041-8213/ade4b6 hal-05248130

HAL Id: hal-05248130

https://hal.science/hal-05248130v1

Submitted on 11 Sep 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Transient Gamma Rays from the 2021 Outburst of the Recurrent Nova RS Ophiuchi: The Effect of Gamma-Ray Absorption

Vo Hong Minh Phan¹, Pierre Cristofari², Enrico Peretti³,⁴, Vincent Tatischeff⁵, and Andrea Ciardi¹

¹ Sorbonne Université, Observatoire de Paris, PSL Research University, LUX, CNRS, 75005 Paris, France; vhmphan@obspm.fr

² LUX³ Observatoire de Paris, Université PSL, Sorbonne Université, CNRS, 92190 Meudon, France

³ INAF—Astrophysical Observatory of Arcetri, Largo E. Fermi 5, 50125 Florence, Italy

⁴ APC, Université Paris Diderot, CNRS/IN2P3, CEA/Irfu, Observatoire de Paris, Sorbonne Paris Cité, France

⁵ Université Paris-Saclay, CNRS/IN2P3, IJCLab, 91405 Orsay, France

Received 2025 March 18; revised 2025 June 10; accepted 2025 June 13; published 2025 August 28

## Abstract

In 2021, RS Ophiuchi was the first nova to be detected in the very-high-energy (TeV) gamma-ray domain, directly testifying to efficient acceleration of charged particles up to at least the TeV range at the nova shock. Surprisingly, the TeV gamma-ray signal peaks ~2 days after the GeV signal, and the origin of this delay has still not been clearly understood. We investigate the possibility that this delay is due to the effect of gamma-ray absorption resulting from interactions between gamma rays and optical photons copiously emitted during the outburst. We model particle acceleration at a nova shock to obtain the gamma-ray emission produced in interactions between the accelerated particles and the shocked gas. The effect of gamma-ray absorption is then included in detail using the radiative transfer equation. We find that this can naturally account for the delay between the peaks of GeV and TeV gamma-ray lightcurves. This result emphasizes the importance of gamma-ray absorption for interpreting gamma-ray observations of novae in the TeV range, which, in turn, demonstrates the necessity of a multiwavelength view for unraveling the underlying physics of particle acceleration in these systems.

*Unified Astronomy Thesaurus concepts: Recurrent novae (1366); Gamma-ray transient sources (1853); Shocks (2086)*

## 1. Introduction

Novae are transients observed in many different wavelengths from radio to X-rays and gamma rays. In classical and recurrent novae, the increased brightness results from thermonuclear explosion on the surface of white dwarfs (WDs). These explosions lead to formation of shocks, which can accelerate particles up to TeV energies or beyond in certain cases (B. D. Metzger et al. 2015, 2016). The presence of these high-energy particles indeed means that nonthermal emissions, especially gamma rays induced by leptonic processes (e.g., bremsstrahlung or inverse Compton scattering) and hadronic processes (e.g., decay of neutral pions produced in proton-proton interactions) are expected from these sources.

Many novae were, in fact, observed in GeV gamma rays by the Fermi Large Area Telescope (LAT) gamma-ray telescope (see, e.g., L. Chomiuk et al. 2021, for a review). Recently, both HESS and MAGIC gamma-ray telescopes reported the detection of TeV gamma rays from the 2021 outburst of the recurrent nova RS Ophiuchi (RS Oph), making this system the first nova ever to be observed in the TeV gamma-ray domain (V. A. Acciari et al. 2022; F. Aharonian et al. 2022). RS Oph consists of a WD embedded in the wind of a red giant (RG). The WD continuously accretes materials from the RG wind, and a nova explosion occurs when there is a thermonuclear runaway due to hydrogen thermonuclear burning in degenerate conditions of the WD's outer layers. In the case of RS Oph, the nova explosions occur approximately every 10 or 20 yr (B. E. Schaefer 2010). The most recent explosion was in

2021, with the shock speed inferred by X-ray and radio observations to be around a few thousand km s⁻¹ (C. C. Cheung et al. 2022; R. Pandey et al. 2022), which is of the same order of magnitude as in the previous explosion in 2006 (M. F. Bode et al. 2006; R. Das et al. 2006; J. L. Sokoloski et al. 2006). The relatively high shock speed, the strong magnetic field strength (believed to be of order of a few Gauss close to the RG), and their detection in gamma rays from GeV to TeV energy make novae ideal sources for the study of particle acceleration (see e.g., P. Martin & G. Dubus 2013; P. Martin et al. 2018).

Interestingly, gamma-ray observations revealed a time delay between the peaks of the GeV and TeV lightcurves. It has been suggested in some studies that the standard acceleration model, where the injection spectrum of accelerated particles around the shock follows a simple power law in momentum with an exponential cutoff around the maximum particle energy, cannot explain this feature, and thus, such a time delay can potentially provide new insights into the process of shock acceleration. In fact, many of the previous works have attempted to explain this feature using a more gradual exponential cutoff for the shock injection spectrum (F. Aharonian et al. 2022; J.-H. Zheng et al. 2022). In other approaches, the typical shock injection spectrum is kept, but multiple populations of accelerated protons have to be introduced. For example, R. Diesing et al. (2023) suggest the presence of two shocks where the slower one dominates the GeV gamma-ray emissions in the first 1 or 2 days of the outburst, and the faster shock accounts for the TeV emissions at a later time.

Since these models have to introduce more parameters than a typical shock acceleration model in order to obtain good fits for the delay between the peaks of the GeV and TeV lightcurves, we would like to offer an alternative explanation, which relies on the effect of gamma-ray absorption. Given that

<figure><img src="image_3.png" /></figure>

Original content from this work may be used under the terms of the [Creative Commons Attribution 4.0 licence](https://creativecommons.org/licenses/by/4.0/). Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
the nova is also very bright in optical light, TeV gamma rays can be absorbed as they interact with optical photons. In fact, the role of gamma-ray absorption has previously been studied, most notably by V. A. Acciari et al. (2022) and R. Diesing et al. (2023), who concluded that this effect should be negligible based on their approximate estimates of the gamma-ray opacity. However, none of the abovementioned works have performed a thorough analysis with the radiative transfer of gamma rays taken into account to verify if modifications of the standard particle acceleration model are actually required. The aim of this work is to fill in this gap.

We study gamma-ray emission from the 2021 outburst of RS Oph using a model of particle acceleration in a single nova shock. We choose $t=0$ day at 1 day before the optical peak of the 2021 outburst, as in F. Aharonian et al. (2022), which corresponds to Modified Julian Date 59434.25 (or 2021 August 8.25 in Coordinated Universal Time). The Letter is organized as follows. We first discuss the dynamics of nova shock and the corresponding acceleration process of non-thermal particles. The gamma-ray emission from accelerated hadrons interacting with the compressed material of the RG wind is then modeled. Finally, we include the gamma-ray absorption due to interactions between gamma rays and optical photons using the radiative transfer equation. We find that the gamma-ray absorption can naturally account for the time delay of ~2 days between the GeV and TeV lightcurves.

## 2. Particle Acceleration in Nova Shocks

### 2.1. Dynamics of Nova Shocks

The gas mass density profile at a distance $r$ from the WD can be modeled as follows (F. Aharonian et al. 2022):

$$ \rho(r) \approx \frac{\dot{M}}{4\pi v_{\text{wind}}(r_{\text{orb}}^2 + r^2)}, \qquad (1) $$

where $\dot{M}$ is the mass-loss rate of the RG, $v_{\text{wind}}$ is the wind speed, and $r_{\text{orb}}$ is the orbital radius of the RG. Recent analyses of Hα emission line and Na I D absorption lines indicate a mass-loss rate of about $5 \times 10^{-7} M_{\odot} \text{yr}^{-1}$ (R. A. Booth et al. 2016), which, for the typical wind speed for the RG between 10 and 30 km s$^{-1}$ (H. J. G. L. M. Lamers & J. P. Cassinelli 1999), gives the ratio $\dot{M}/v_{\text{wind}}$ of about a few times $10^{13} \text{ g cm}^{-1}$. This ratio determines the normalization of the gamma flux, as we shall see later, but it is degenerate with several other parameters. Thus, for simplicity, we will fix this ratio to $\dot{M}/v_{\text{wind}} = 2 \times 10^{13} \text{ g cm}^{-1}$ for the fit of gamma-ray data similar to the value used in modeling multiwavelength data in the previous outburst by V. Tatischeff & M. Hernanz (2007). As for the orbital radius, it is relatively well constrained to the value of about 1.48 au with spectroscopic optical observations (E. Brandi et al. 2009).

The nova shock speed $v_{\text{sh}}$ is expected to evolve as follows:

$$ v_{\text{sh}}(t) = \begin{cases} v_{\text{sh},0} & t \le t_r, \\ v_{\text{sh},0} \left(\frac{t}{t_r}\right)^{-\alpha} & t > t_r, \end{cases} \qquad (2) $$

where initially the shock expands with a constant speed in the free expansion phase, and then, at $t \simeq t_r$, it transitions into the radiative phase where $v_{\text{sh}} \sim t^{-\alpha}$ with $\alpha \sim 0.5$. It has been suggested also that there might exist a Sedov–Taylor phase (or

adiabatic phase) where $v_{\text{sh}} \sim t^{-1/3}$ (M. F. Bode & F. D. Kahn 1985; R. Das et al. 2006). Optical and X-ray observations, however, indicate that this phase may be too brief to be detected. We note also that, in the above equations for shock evolution, the eruption time has been implicitly assumed to happen at $t_{\text{er}}=0$ day. We checked, however, that the fit to gamma-ray data is not very sensitive to the exact value of $t_{\text{er}}$ as long as $t_{\text{er}} \lesssim 0.2$ day.

Shock speeds inferred from postshock gas temperature measurements using X-ray observations seem to indicate $v_{\text{sh},0} \simeq 2470 \text{ km s}^{-1}$, $\alpha \simeq 0.43$, and $t_r \simeq 6$ days (C. C. Cheung et al. 2022; see also M. Orio et al. 2023). Spectroscopic data from optical observations, on the other hand, are well fitted for $v_{\text{sh},0} \gtrsim 4000 \text{ km s}^{-1}$, $\alpha \simeq 0.6$, and $t_r \simeq 4$ days (R. Pandey et al. 2022). In fact, V. Tatischeff & M. Hernanz (2007) have pointed out that back-reactions of particles accelerated in nova shocks can modify relations between shock speeds and postshock gas temperatures and explain the differences between shock properties inferred by X-ray and optical observations. In the following, we will treat $v_0$, $t_r$, and $\alpha$ as free parameters that can be obtained by fitting simultaneously the shock speed profile inferred from optical and X-ray observations and the gamma-ray emissions.

### 2.2. Maximum Proton Energy

The acceleration rate of protons in shocks can be estimated as (see e.g., V. Tatischeff & M. Hernanz 2007)

$$ \dot{E}_{\text{acc}}(t) \approx \frac{qB_2(t)v_{\text{sh}}^2(t)}{\eta c}, \qquad (3) $$

where $q$ is the charge of proton, $B_2$ is the magnetic field strength downstream of the shock, $c$ is the speed of light, and $1/\eta$ is the acceleration efficiency. It is generally suggested that $\eta \gtrsim 2\pi$, but the exact value is not well known. Since $\eta$ is degenerate with $B_2$, we can fix $\eta = 10\pi$ for simplicity (roughly similar to the acceleration efficiency assumed by F. Aharonian et al. 2022). We will assume $B_2$ to be the compressed background magnetic field around the RG as in F. Aharonian et al. (2022) with the following profile:

$$ B_2 = B_{2,0} \left( \frac{\sqrt{r_{\text{orb}}^2 + R_{\text{sh}}^2}}{0.35 \text{ au}} \right)^{-2}, \qquad (4) $$

where $B_{2,0}$ ranges from 1 to 10 G and $R_{\text{sh}}(t) = \int_0^t v_{\text{sh}}(t')dt'$ is the shock radius.

Particles are also losing energy over time due to adiabatic expansion of the shock. The energy loss rate can be modeled as follows (D. Caprioli et al. 2010; P. Cristofari et al. 2020):

$$ \dot{E}_{\text{ad}}(E, t) \approx \frac{p v}{(\rho v_{\text{sh}}^2(t))^{\frac{1}{3}\gamma_{\text{gas}}}} \frac{d}{dt} \left[ (\rho v_{\text{sh}}^2(t))^{\frac{1}{3}\gamma_{\text{gas}}} \right], \qquad (5) $$

where $\gamma_{\text{gas}} = 5/3$ is the adiabatic index of the gas.

The maximum kinetic energy of protons at any time $t$ can then be evaluated by solving

$$ \frac{dE_{\text{max}}}{dt} = \dot{E}_{\text{acc}}(t) + \dot{E}_{\text{ad}}(E_{\text{max}}, t). \qquad (6) $$
Once the nova shock dynamics and the maximum energy of accelerated particles are defined, we can proceed to discuss the spectrum of these particles accumulated around the shock.

## 2.3. Cumulative Proton Spectrum

The injection spectrum of protons around the shocks is typically assumed to be a power law in momentum $p$ with an exponential cutoff around $p_{max}$, which can be estimated from $E_{max}$ with Equation (6). The spectrum is also normalized such that, at any given time, the pressure of accelerated protons is a fraction $\xi_p$ of the shock ram pressure (P. Cristofari et al. 2021), which gives

$$f_p(E, t) = \frac{3\xi_p \rho v_{\text{sh}}^2}{m_p^2 c^4 \beta I(\delta)} \left(\frac{p}{m_p c}\right)^{2-\delta} \exp\left(-\frac{p}{p_{\text{max}}}\right). \quad (7)$$

Here, $m_p$ is the proton mass, $\rho$ and $v_{sh}$ are estimated at the shock radius, and $I(\delta) = \int_{x_{min}}^{x_{max}} dx x^{4-\delta} e^{-x/x_{max}} / \sqrt{1 + x^2}$ with $x = p/(m_p c)$. The cumulative spectrum of protons accelerated at the nova shock can then be obtained by solving the following transport equation (see e.g., F. Aharonian et al. 2022):

$$\frac{\partial N_p(E, t)}{\partial t} + \frac{\partial}{\partial E}[\dot{E}_{\text{ad}}(E, t)N_p(E, t)] = \frac{4\pi R_{\text{sh}}^2(t)v_{\text{sh}}(t)}{r_c} f_p(E, t), \quad (8)$$

where $r_c$ is the shock compression ratio. We will assume $r_c \approx 4$, as commonly adopted for strong shocks (P. Cristofari et al. 2021).

# 3. Gamma-Ray Emissions from Nova Shocks

## 3.1. Gamma-Ray Production

Interactions between accelerated protons and ambient protons (materials swept up by the nova shock) can lead to the production of neutral pions, which decay quickly into gamma rays. These gamma rays can be partially absorbed by interactions with optical photons produced from the outburst (see the next subsection). However, previous works (e.g., R. Diesing et al. 2023), based on their approximate estimates of the gamma-ray opacities, suggest that gamma-ray absorption is negligible, such that the expected gamma-ray flux from the nova can be estimated as follows:

$$\phi_0(E_\gamma, t) \approx \frac{\rho_{\text{sh}}(t)}{4\pi d_N^2 m_p} \int_{E_{\text{min}}}^{E_{\text{max}}} dE N_p(E, t) \times v_p \epsilon_n(E) \frac{d\sigma_{\text{pp}}(E, E_\gamma)}{dE_\gamma}, \quad (9)$$

where $\rho_{sh}(t) \approx 4\rho(r=R_{sh}(t))$ is the mass density of shocked gas, $d_N$ is the distance of the nova to Earth, $\epsilon_n(E)$ is the nuclear enhancement factor to take into account gamma rays from accelerated nuclei heavier than proton, and $d\sigma_{pp}(E, E_\gamma)/dE_\gamma$ is the differential cross section for gamma-ray production in proton-proton interaction (E. Kafexhiu et al. 2014).

We note that electrons are also expected to be accelerated from nova shocks (with not well-constrained acceleration efficiency) and should be taken into account, e.g., to interpret

radio emissions or better constrain the magnetic field strength from this system (see, e.g., A. De Sarkar et al. 2023). These particles can also produce leptonic gamma rays via bremsstrahlung or inverse Compton scattering. The focus of this work is, however, to better understand the gamma-ray emissions, and we found, as explained below, that hardonic gamma rays alone are sufficient for a good fit to the data when the effect of gamma-ray absorption is taken into account. More thorough investigations with leptonic emissions will be examined in future works.

## 3.2. Gamma-Ray Absorption

### 3.2.1. Distribution of Optical Photons

As mentioned above, since nova outbursts can be very bright in optical light, the medium around these novae may be opaque to gamma rays of energy above about 200 GeV in the first days of the outbursts. Optical photons from the outburst can be modeled with a blackbody distribution, which can be described with the temperature $T_{opt}$ and the energy density $u(r, t)$. If we assume that optical photons are homogeneously and isotropically emitted within a photosphere of radius $R_{ph}$, it can be shown that the energy density $u(r, t)$ of optical photons follows

$$u(r, t) = \frac{L_{\text{opt}}(t)}{4\pi c} \frac{3r}{2R_{\text{ph}}^3} g_{\text{opt}}\left(\frac{r}{R_{\text{ph}}}\right), \quad (10)$$

with

$$g_{\text{opt}}(x) = \frac{1}{x} + \left(\frac{1}{x^2} - 1\right) \ln\left(\sqrt{|\frac{x+1}{x-1}|}\right). \quad (11)$$

Note that we can verify from the above equation that $u(r, t) \sim L_{opt}(t)/(4\pi r^2 c)$ in the limit where $r \gg R_{ph}$, which is commonly adopted for rough estimates of the gamma-ray absorption in previous works (see, e.g., R. Diesing et al. 2023).

For RS Oph 2021, we can adopt $T_{opt} \approx 1.1 \times 10^4$ K (C. C. Cheung et al. 2022), $R_{ph} \approx 200 R_\odot$ (V. A. Acciari et al. 2022), and the optical lightcurve $L_{opt}(t)$ can be inferred from the observed optical spectra. In fact, C. C. Cheung et al. (2022) have already performed the analyses for the optical spectra taking into account the interstellar absorption to derive the optical luminosity over time for the 2021 outburst of RS Oph. These analyses assume, however, that the distance from RS Oph to Earth is about 1.6 kpc. This distance, which was first estimated using the H I absorption line (R. M. Hjellming et al. 1986), is quite commonly adopted in the literature. However, such a value means that the expected accretion rate might be much smaller than required for the typical recurrence period of RS Oph (B. E. Schaefer 2009). More recent estimates of the distance, e.g., using high-resolution radio observations together with the measurements of the shock speed profile (M. P. Rupen et al. 2008) or the parallax distance provided by Gaia (Gaia Collaboration et al. 2021), lead to values around 2.45 or 2.68 kpc (see V. A. Acciari et al. 2022 for a more complete discussion).

In the following, we will follow V. A. Acciari et al. (2022) and choose $d_N = 2.45$ kpc for the fit of the gamma-ray emission. In order to precisely obtain the optical luminosity over time at $d = 2.45$ kpc, we need to reanalyze the optical spectra and correct for interstellar absorption of the source emission at this distance. Note however that RS Oph is at the
latitude of 10.37, meaning that the correction factor for interstellar absorption should be rather similar for $d_N = 1.6$ and 2.45 kpc. For this reason, the optical lightcurve can simply be rescaled with distance for $d_N > 1.6$ kpc. We then choose the optical lightcurve of the following form for $t \gtrsim 1$ day:

$$L_{\text{opt}}(t) \approx 7.8 \times 10^{38} \left( \frac{d_N}{1.6 \text{ kpc}} \right)^2 \left( \frac{t}{1 \text{ day}} \right)^{-1} \text{erg s}^{-1}, \quad (12)$$

where the normalization and the scaling with time have been fitted to optical luminosity data derived by C. C. Cheung et al. 2022 for $d_N = 1.6$ kpc (see Figure 5 in Appendix A).

### 3.2.2. Radiative Transfer of Gamma Rays

The effect of absorption can be taken into account by solving the radiative transfer equation for gamma rays. We have assumed that gamma rays are mostly produced in a thin layer around the nova shock, which allows us to express the gamma-ray flux in the following form:

$$\phi(E_{\gamma}, t) = \frac{\phi_0(E_{\gamma}, t)}{2} [e^{-\tau_1(E_{\gamma}, t)} + e^{-\tau_2(E_{\gamma}, t)}], \quad (13)$$

where the exact form of the opacities $\tau_1(E_{\gamma}, t)$ and $\tau_2(E_{\gamma}, t)$ is derived in more detail in Appendix A.

In this representation, we can interpret roughly that half of the gamma rays propagate through the region with $r \le R_{\text{sh}}$ and experience absorption as represented by the average opacity $\tau_1(E_{\gamma}, t)$. The other half escape directly into the region $r \ge R_{\text{sh}}$ and are attenuated with the average opacity $\tau_2(E_{\gamma}, t)$. In fact, $\tau_2(E_{\gamma}, t)$ is similar to the rough estimate of the gamma-ray opacity provided in many of the previous works (e.g., R. Diesing et al. 2023). However, the exact value of $\tau_2(E_{\gamma}, t)$ can vary by a factor of a few depending on the chosen value of $d_N$ (because $L_{\text{opt}} \sim d_N^2$; see Equation (12)), which can be significant for the gamma-ray lightcurve as the effect of absorption is exponential. More importantly, when the full radiative transfer is taken into account, we can see that half of the gamma rays have to pass through the region $r \le R_{\text{sh}}$ and are more strongly absorbed ($\tau_1(E_{\gamma}, t) > \tau_2(E_{\gamma}, t)$) since optical photon density is higher close to the WD.

### 3.3. Results for the 2021 Outburst of RS Oph

In this subsection, we present the results of our modeling of the transient gamma-ray emission from a nova with gamma-ray absorption applied for the 2021 outburst of RS Oph. All the parameters of the model are summarized in Table 1 in Appendix B. Some parameters are fixed to values obtained in previous works since they are already well constrained with other observations.

In Figure 1, we present the gamma-ray spectra assuming parameters that best fit the Fermi-LAT, HESS, and MAGIC data (V. A. Acciari et al. 2022; F. Aharonian et al. 2022). We found that the effect of absorption is most significant between $t = 1$ day and $t = 2$ days for gamma rays of energy $\gtrsim 200$ GeV. The gamma-ray spectrum $\gtrsim 200$ GeV at $t \approx 1.6$ days (corresponding to night 1 in HESS convention) is significantly suppressed such that it matches well with data without the need for modifications of the injection spectrum to suppress the amount of accelerated TeV protons, as commonly adopted in previous works. Since proton-proton interactions also produce neutrinos, a potentially observational consequence

<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /><figcaption>Figure 1. Gamma-ray spectra compared to data from F. Aharonian et al. (2022) and V. A. Acciari et al. (2022). Solid and dashed lines are for the case with and without gamma-ray absorption, respectively.</figcaption></figure>

of our model would be that one might be able to detect TeV neutrinos from novae like RS Oph in the first days of their outbursts even though TeV gamma rays might be undetectable in the same period. However, the predicted flux of TeV neutrinos from RS Oph would still be too low to be detected by current neutrino observatories, such that this phenomenon might be tested only with future instruments.
<figure><img src="image_7.png" /><figcaption>**Figure 2.** Shock speed evolution in time (dashed line) compared to data inferred from $H\alpha$ and $H\beta$ emission lines (red and green filled circles) by R. Pandey et al. (2022) and from X-ray data (yellow filled circles) by M. Orio et al. (2023).</figcaption></figure>

It is worth mentioning again that the shock speed can also be measured using infrared and X-ray observations (with rather large uncertainty). The shock speed measurements derived from the $H\alpha$ and $H\beta$ emission line profile modeling (R. Pandey et al. 2022) and from the postshock temperature<fnref n="6" /> inferred from observations with the NICER X-ray telescope (M. Orio et al. 2023) are collected in Figure 2. The shock speed evolution that best fits both these measurements and the gamma-ray data from Fermi-LAT and HESS is also presented.

We show also the integrated gamma-ray flux over time in Figure 3. The gamma-ray spectra were integrated in the ranges 0.1–100 and 250–2500 GeV to obtain lightcurves that can be compared with the GeV and TeV lightcurves measured with Fermi-LAT and HESS. The GeV and TeV lightcurves are presented without (dashed line) and with (solid line) gamma-ray absorption. In the case without absorption, the predicted GeV and TeV lightcurves peak at the same time. This is because the maximum energy of accelerated particles is reached roughly after $t \simeq 1$ day, which makes the spectral shape of accelerated particles rather similar for $t \gtrsim 1$ day. The lightcurves in both GeV and TeV energy ranges are, therefore, determined mostly by the evolution of the shock and have rather similar time dependence at a later time (after reaching the peaks) for the case without absorption. It is clear however that we overestimate the emission for TeV gamma rays at early time in this case.

A better fit of the TeV lightcurve is obtained when gamma-ray absorption is correctly taken into account. Figure 3 shows that the gamma-ray absorption is strongest at $t \simeq 1$ day, which is around the maximum of the optical lightcurve. Note also that the TeV lightcurve in this case is plotted differently (dotted red line) for $t \lesssim 1$ day to indicate that the shape of the TeV lightcurve depends strongly on our interpolation of the optical luminosity for −2 days $\lesssim t \lesssim 1$ day (data on the optical luminosity of RS Oph derived by C. C. Cheung et al. 2022 are not available in this period). We can also see that the effect of gamma-ray absorption naturally leads to a delay of about 1 or 2 days between the peaks of the GeV and TeV lightcurves.

<figure><img src="image_8.png" /><figcaption>**Figure 3.** Integrated gamma-ray flux over time compared to data from F. Aharonian et al. (2022). The lightcurves for Fermi-LAT and HESS (green and red line, respectively) are estimated by integrating the gamma-ray flux over the energy ranges 0.1–100 and 250–2500 GeV, respectively. Solid and dashed lines are correspondingly for the case with and without gamma-ray absorption. The TeV lightcurve with absorption for $t < 1$ day is plotted in dotted line to highlight the uncertainty of the model in this period of time where no data are available (see the discussion in Section 3.3).</figcaption></figure>

## 4. Summary and Prospects

We study the transient gamma-ray emissions from the most recent outburst of RS Oph in 2021 with a focus on interpreting the delay between the peaks of GeV and TeV lightcurves. Using a radiative transfer model for gamma rays produced in interactions between particles accelerated around the nova shock and the RG wind materials, we have shown gamma rays passing through the shock downstream are strongly absorbed such that half of the gamma rays with $E_γ \gtrsim 200$ GeV are almost completely absorbed in the first 1 or 2 days of the outburst. The remaining half escaping directly into the upstream region of the shock are absorbed with an opacity denoted as $\tau_2$ (see Equation (13)) comparable to the approximate opacity estimated by, e.g., R. Diesing et al. (2023). This means that the gamma-ray flux, for $t \lesssim 1$ day and $E_γ \gtrsim 200$ GeV, with and without absorption should be different by a factor $\sim e^{-\tau_2/2}$, which can be between 6 and more than 10 depending on the exact value of the nova distance $d_N = 1.4$ or 2.45 kpc.

More importantly, the effect of gamma-ray absorption naturally accounts for the time delay between the peaks of GeV and TeV gamma-ray lightcurves such that the modifications of the shock injection spectrum or the presence of multiple accelerated proton populations are no longer necessary for explaining the data (F. Aharonian et al. 2022; R. Diesing et al. 2023). This result highlights the important role of gamma-ray absorption for better understanding the TeV gamma-ray emission from novae, which is also crucial for deciphering the underlying process of acceleration, especially around the maximum energy of shock-accelerated particles.

## Acknowledgments

V.H.M.P. acknowledges support from the Initiative Physique des Infinis (IPI), a research training program of the Idex SUPER at Sorbonne Université. P.C. acknowledges support from the PSL Starting Grant GALAPAGOS. E.P. acknowledges support from INAF through “Assegni di ricerca per progetti di ricerca relativi a CTA e precursori.”

<footnote n="6">⁶ We convert the postshock temperature from Figure 6 of M. Orio et al. (2023) into shock speed using Equation (3) of M. F. Bode et al. (2006).</footnote>
## Appendix A
### Gamma-Ray Opacities

In order to take into account the effect of gamma-ray absorption, we will first find the gamma-ray intensity along each line of sight $I(E_γ, θ, s \approx d_N, t)$ by solving the radiative transfer equation

$$ \frac{dI(E_γ, θ, s, t)}{ds} = -I(E_γ, θ, s, t)η_{abs}(E_γ, r(θ, s), t) + I_0(E_γ, t)δ(r(θ, s) - R_{sh}), \quad (A1) $$

where $r(θ, s)$, $η_{abs}(E_γ, r, t)$, and $I_0(E_γ, t)$ are, respectively, the distance from the WD, the absorption coefficient, and the gamma-ray intensity around the shock, which could be estimated as follows:

$$ r(θ, s) = \sqrt{s^2 - 2s\sqrt{R_{sh}^2 - d^2 \sin^2 θ} + R_{sh}^2}, \quad (A2) $$

$$ η_{abs}(E_γ, r, t) = \int_0^∞ dE_{ph} f_{opt}(E_{ph}, r, t) σ_{γγ}(E_γ, E_{ph}), \quad (A3) $$

$$ I_0(E_γ, t) = \frac{\phi_0(E_γ, t)}{4\pi} \left(\frac{d_N}{R_{sh}}\right)^2. \quad (A4) $$

Note that we have also introduced the differential number density of optical photons $f_{opt}(E_{ph}, r, t)$, which, as mentioned above, can be modeled with a blackbody distribution as follows:

$$ f_{opt}(E_{ph}, r, t) = u_{opt}(r, t) \times \frac{15 \left( \frac{E_{ph}}{k_B T_{opt}} \right)^2}{\pi^4 (k_B T_{opt})^2 \left[ \exp \left( \frac{E_{ph}}{k_B T_{opt}} \right) - 1 \right]}, \quad (A5) $$

where $u(r, t)$ is the optical photon energy density (see Equation (10)), which can be derived from the optical lightcurve $L_{opt}(t)$ (the optical luminosity as a function of time; $k_B$ is the Boltzmann constant, and $T_{opt}$ the blackbody temperature. Please refer to Figure (4) for a schematic view of the system and Figure (5) for the optical lightcurve of the 2021 outburst of RS Oph.

<figure><img src="image_9.png" /><figcaption>Figure 4. Schematic view of the system.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Figure 5. Left: luminosity function of optical photons (see Equation (12)) compared to data collected in C. C. Cheung et al. (2022) assuming $d_N = 1.6$ kpc. Right: opacities over time estimated from Equations (A9) and (A10) for $E_γ = 1$ TeV.</figcaption></figure>
For each line of sight, we can solve the radiative transfer equation (Equation (A1)) to show that $I(E_γ, θ, s \approx d_N, t)$ always has two contributions corresponding to gamma rays coming from $s = s_1$ and $s = s_2$, which are the two solutions of $r(θ, s) = R_{sh}$ (see Figure 4). We have implicitly chosen the coordinate system such that for all lines of sight $s_1 = 0$ and $s_2 = 2\sqrt{R_{sh} - d_N^2} \sin^2 θ$. Given this choice of the coordinate system, we can write down the solution of the radiative transfer equation as follows:

$$I(E_γ, θ, s \approx d_N, t) = \frac{R_{sh} I_0(E_γ, t)}{\sqrt{R_{sh}^2 - d^2 \sin^2 θ}} \times \left[ \exp\left(-\int_{s_1}^{d_N} ds \, \eta_{abs}(E_γ, r(θ, s), t)\right) + \exp\left(-\int_{s_2}^{d_N} ds \, \eta_{abs}(E_γ, r(θ, s), t)\right) \right]. \quad (A6)$$

We can now derive the gamma-ray flux with absorption by integrating the intensity along each line of sight over the solid angle covering the emission region:

$$\phi(E_γ, t) = 2\pi \int_{0}^{R_{sh}/d_N} d\theta \sin\theta I(E_γ, θ, s = d_N, t). \quad (A7)$$

The factor $2\pi$ comes from integrating over the azimuthal angle. Note also that we have made the approximation $\arcsin(R_{sh}/d_N) \approx R_{sh}/d_N$ since $d_N \gg R_{sh}$. Also, we can show in this way that the gamma-ray flux becomes $\phi_0(E_γ, t)$ when the absorption coefficients are set to zero. For ease of

discussion, we will express the gamma-ray flux as follows:

$$\phi(E_γ, t) = \frac{\phi_0(E_γ, t)}{2} [e^{-\tau_1(E_γ, t)} + e^{-\tau_2(E_γ, t)}], \quad (A8)$$

where

$$e^{-\tau_1(E_γ, t)} = \left(\frac{d_N}{R_{sh}}\right)^2 \int_0^{R_{sh}/d_N} d\theta \frac{R_{sh} \sin\theta}{\sqrt{R_{sh}^2 - d^2 \sin^2\theta}} \exp \left( -\int_0^{d_N} ds \, \eta_{abs}(E_γ, r(θ, s), t) \right), \quad (A9)$$

and

$$e^{-\tau_2(E_γ, t)} = \left(\frac{d_N}{R_{sh}}\right)^2 \int_0^{R_{sh}/d_N} d\theta \frac{R_{sh} \sin\theta}{\sqrt{R_{sh}^2 - d^2 \sin^2\theta}} \exp \left( -\int_0^{d_N} ds \, \eta_{abs}(E_γ, r(θ, s), t) \right). \quad (A10)$$

We also present estimates of the opacities for $E_γ = 1$ TeV in Figure 5. It is clear from this plot and also from the above equations that gamma rays that have to pass through the shock downstream are more strongly absorbed, and thus, $\tau_1(E_γ, t) > \tau_2(E_γ, t)$ (see also the discussions in Section 3.2.2).

## Appendix B
### Fit Parameters

We provide in Table 1 of this appendix a summary for all the fit parameters for gamma-ray emissions from the 2021 outburst of RS Oph.

<table><caption>**Table 1**
Parameters for Fitting the Gamma-Ray Data from the 2021 Outburst of RS Oph</caption><thead><tr><th>Parameters</th><th>Descriptions</th><th>Comments</th><th>Values</th></tr></thead><tbody><tr><td>t<sub>er</sub></td><td>Eruption time</td><td>Fixed parameters</td><td>0 day</td></tr><tr><td>r<sub>orb</sub></td><td>Orbital radius of RG (E. Brandi et al. 2009)</td><td>...</td><td>1.48 au</td></tr><tr><td>d<sub>N</sub></td><td>Distance from RS Oph to Earth (M. P. Rupen et al. 2008)</td><td>...</td><td>2.45 kpc</td></tr><tr><td>T<sub>opt</sub></td><td>Temperature of optical photons (C. C. Cheung et al. 2022)</td><td>...</td><td>1.1 × 10<sup>4</sup> K</td></tr><tr><td>R<sub>ph</sub></td><td>Radius of the photosphere (V. A. Acciari et al. 2022)</td><td>...</td><td>200 R<sub>☉</sub></td></tr><tr><td>M/v<sub>wind</sub></td><td>Mass-loss rate of RG over wind speed of RG (V. Tatischeff & M. Hernanz 2007)</td><td>...</td><td>2 × 10<sup>13</sup> g cm<sup>-1</sup></td></tr><tr><td>ξ<sub>p</sub></td><td>Fraction of shock ram pressure converted into pressure of accelerated protons</td><td>Fitted parameters</td><td>0.14</td></tr><tr><td>v<sub>sh,0</sub></td><td>Initial speed of the shock</td><td>...</td><td>3500 km s<sup>-1</sup></td></tr><tr><td>t<sub>r</sub></td><td>Transition time of the shock speed (Equation (2))</td><td>...</td><td>2.2 days</td></tr><tr><td>α</td><td>Power-law index of the shock speed (Equation (2))</td><td>...</td><td>0.43</td></tr><tr><td>δ</td><td>Spectral index of the injection spectrum</td><td>...</td><td>4.2</td></tr><tr><td>B<sub>2,0</sub></td><td>Magnetic field strength close to RG (Equation (4))</td><td>...</td><td>6.5 G</td></tr></tbody></table>
ORCID iDs

Vo Hong Minh Phan https://orcid.org/0000-0002-5611-095X

Pierre Cristofari https://orcid.org/0000-0002-5856-7662

Enrico Peretti https://orcid.org/0000-0003-0543-0467

Vincent Tatischeff https://orcid.org/0009-0004-4190-6068

Andrea Ciardi https://orcid.org/0000-0001-5865-6772

References

Acciari, V. A., Ansoldi, S., Antonelli, L. A., et al. 2022, NatAs, 6, 689

Aharonian, F., Ait Benkhali, F., Angüner, E. O., et al. 2022, Sci, 376, 77

Bode, M. F., & Kahn, F. D. 1985, MNRAS, 217, 205

Bode, M. F., O'Brien, T. J., Osborne, J. P., et al. 2006, ApJ, 652, 629

Booth, R. A., Mohamed, S., & Podsiadlowski, P. 2016, MNRAS, 457, 822

Brandi, E., Quiroga, C., Mikołajewska, J., Ferrer, O. E., & Garcia, L. G. 2009, A&A, 497, 815

Caprioli, D., Amato, E., & Blasi, P. 2010, APh, 33, 160

Cheung, C. C., Johnson, T. J., Jean, P., et al. 2022, ApJ, 935, 44

Chomiuk, L., Metzger, B. D., & Shen, K. J. 2021, ARA&A, 59, 391

Cristofari, P., Blasi, P., & Amato, E. 2020, APh, 123, 102492

Cristofari, P., Blasi, P., & Caprioli, D. 2021, A&A, 650, A62

Das, R., Banerjee, D. P. K., & Ashok, N. M. 2006, ApJL, 653, L141

De Sarkar, A., Nayana, A. J., Roy, N., Razzaque, S., & Anupama, G. C. 2023, ApJ, 951, 62

Diesing, R., Metzger, B. D., Aydi, E., et al. 2023, ApJ, 947, 70

Gaia Collaboration, Brown, A. G. A., Vallenari, A., et al. 2021, A&A, 649, A1

Hjellming, R. M., van Gorkom, J. H., Taylor, A. R., et al. 1986, ApJL, 305, L71

Kafexhiu, E., Aharonian, F., Taylor, A. M., & Vila, G. S. 2014, PhRvD, 90, 123014

Lamers, H. J. G. L. M., & Cassinelli, J. P. 1999, Introduction to Stellar Winds (Cambridge: Cambridge Univ. Press), 452

Martin, P., & Dubus, G. 2013, A&A, 551, A37

Martin, P., Dubus, G., Jean, P., Tatischeff, V., & Dosne, C. 2018, A&A, 612, A38

Metzger, B. D., Caprioli, D., Vurm, I., et al. 2016, MNRAS, 457, 1786

Metzger, B. D., Finzell, T., Vurm, I., et al. 2015, MNRAS, 450, 2739

Orio, M., Gendreau, K., Giese, M., et al. 2023, ApJ, 955, 37

Pandey, R., Habtie, G. R., Bandyopadhyay, R., et al. 2022, MNRAS, 515, 4655

Rupen, M. P., Mioduszewski, A. J., & Sokoloski, J. L. 2008, ApJ, 688, 559

Schaefer, B. E. 2009, ApJ, 697, 721

Schaefer, B. E. 2010, ApJS, 187, 275

Sokoloski, J. L., Luna, G. J. M., Mukai, K., & Kenyon, S. J. 2006, Natur, 442, 276

Tatischeff, V., & Hernanz, M. 2007, ApJL, 663, L101

Zheng, J.-H., Huang, Y.-Y., Zhang, Z.-L., et al. 2022, PhRvD, 106, 103011