Issue 
A&A
Volume 642, October 2020



Article Number  A41  
Number of page(s)  9  
Section  Numerical methods and codes  
DOI  https://doi.org/10.1051/00046361/202038151  
Published online  02 October 2020 
Correlated energy uncertainties in reaction rate calculations
^{1}
North Carolina State University, Raleigh, NC 27695, USA
email: richard_longland@ncsu.edu
^{2}
Triangle Universities Nuclear Laboratory, Durham, NC 27708, USA
^{3}
Université ParisSaclay, CNRS/IN2P3, IJCLab, 91405 Orsay, France
Received:
11
April
2020
Accepted:
7
July
2020
Context. Monte Carlo methods can be used to evaluate the uncertainty of a reaction rate that arises from many uncertain nuclear inputs. However, until now no attempt has been made to find the effect of correlated energy uncertainties in input resonance parameters.
Aims. Our goal is to investigate the impact of correlated resonance energy uncertainties on reaction rates.
Methods. Using a combination of numerical and Monte Carlo variation of resonance energies, the effect of correlations are investigated. Five reactions are considered: two fictional, illustrative cases and three reactions whose rates are of current interest.
Results. The effect of correlations in resonance energies depends on the specific reaction cross section and temperatures considered. When several resonances contribute equally to a reaction rate, and when they are located on either side of the Gamow peak, correlations between their energies dilute their effect on reaction rate uncertainties. If they are both located above or below the maximum of the Gamow peak, however, correlations between their resonance energies can increase the reaction rate uncertainties. This effect can be hard to predict for complex reactions with wide and narrow resonances contributing to the reaction rate.
Key words: methods: numerical / methods: statistical / nuclear reactions / nucleosynthesis / abundances
© ESO 2020
1. Introduction
Thermonuclear reaction rates dictate energy generation and elemental synthesis in stars and stellar explosions. They are a key physical input to computational stellar models that attempt to explain astrophysical phenomena in conjunction with observational data. It is for this reason that reaction rates must be well known and, moreover, their uncertainties must be well understood if comparisons between stellar models and observations are to be made reliably. Reaction rate uncertainties were first addressed by the Nuclear Astrophysics Compilation of REaction rates (NACRE) collaboration in an evaluation of 86 reactions on A = 1 − 30 nuclei (Angulo et al. 1999). In that time period, more sophisticated techniques were developed to quantify reaction rate uncertainties (Thompson & Iliadis 1999; Iliadis et al. 1999). Those techniques were used in the reaction rate evaluation of Iliadis et al. (2001). While that method was developed to accurately propagate cross section uncertainties to nuclear reaction rates, some approximations were necessary. For example, the propagation of energy uncertainties required that they be “very small” (Thompson & Iliadis 1999). The method was also not able to propagate the uncertainty in rates dominated by broad resonances.
Over the last decade, advances have been made in characterising thermonuclear reaction rate uncertainties using statistically rigorous practices. Longland et al. (2010) introduced a Monte Carlo method for propagating the uncertainties in resonant reaction cross sections to reaction rates by carefully assigning probability density distributions to individual experimental inputs. Similar methods have been developed in, for example, Rauscher (2012). The method of Longland et al. (2010) was further improved by Sallaska et al. (2013) and extended to account for uncertain spinparities (Mohr et al. 2014). The methods have been applied to nucleosynthesis in a number of environments, for example, Longland (2012), Rauscher et al. (2016), and Cescutti et al. (2018), including Monte Carlo variations in full hydrodynamic models (Fields et al. 2016, 2018). Correlations between input parameters has recently been addressed by Longland (2017), who introduced a formalism for accounting for resonance parameters that are correlated through some standard normalisation point. That study does not, however, address correlated resonance energies. Resonance energies can be strongly correlated, particularly for unobserved resonances located according to excitation energies in the compound nucleus. If the Qvalue of the reaction is poorly known, those resonance energies are also equally unknown and strongly correlated. This is particularly important for radioactive nuclei where their resonances are rarely measured directly. Another common case arises from resonance energies determined from an excitation function measurement using stable beams with moderately wellknown beam energy. If the beam energy is selected using an analysing or switching magnet, the resonance energies are highly correlated with uncertainties dictated by the calibration of said magnet.
In this paper, we extend the formalism laid out in Longland (2017) to account for correlations in resonance energies. A brief overview of the reaction rate formalism is outlined in Sect. 2, particularly as it pertains to resonance energies. The Monte Carlo reaction rate uncertainty propagation method is reviewed in Sect. 3 and the method for computing correlated resonance energies is developed. The method is tested with illustrative examples in Sect. 4 and real reactions in Sect. 5. A summary of our findings and recommendations is given in Sect. 6.
2. Reaction rate formalism
What follows is a brief overview of astrophysical reaction rates. For more details the reader is encouraged to refer to Rolfs & Rodney (1988), Angulo et al. (1999), Longland et al. (2010), Sallaska et al. (2013), Iliadis (2015). The reaction rate per particle pair, ⟨σv⟩, is defined as
$$\begin{array}{c}\hfill \u27e8\sigma v\u27e9={\left(\frac{8}{\pi \mu}\right)}^{1/2}\frac{1}{{(kT)}^{3/2}}{\displaystyle {\int}_{0}^{\infty}E\sigma (E){e}^{E/kT}dE,}\end{array}$$(1)
where μ is the reduced mass of the system, k is the Boltzmann constant, T is the temperature at which the reaction rate is being calculated, and σ(E) is the energydependent cross section of the reaction.
The cross section can be parameterised by removing the swave Coulomb barrier tunnelling cross section as
$$\begin{array}{c}\hfill \sigma (E)=\frac{1}{E}S(E){e}^{2\pi \eta},\end{array}$$(2)
where η is the Sommerfeldt parameter defined by
$$\begin{array}{c}\hfill 2\pi \eta =0.98951{Z}_{0}{Z}_{1}\sqrt{\frac{\mu}{E}}.\end{array}$$(3)
Here Z_{0} and Z_{1} are the atomic numbers of the interacting particles. The quantity S(E) is the socalled astrophysical Sfactor, which contains any details of the cross section that are not accounted for by simple swave Coulomb barrier scattering. We see from combining Eqs. (2) and (1) that the reaction rate is defined as
$$\begin{array}{c}\hfill \u27e8\sigma v\u27e9={\left(\frac{8}{\pi \mu}\right)}^{1/2}\frac{1}{{(kT)}^{3/2}}{\displaystyle {\int}_{0}^{\infty}S(E){e}^{2\pi \eta}{e}^{E/kT}\mathrm{d}E.}\end{array}$$(4)
The product of the two exponentials: the Gamow factor, e^{−2πη} and the Boltzmann factor, e^{−E/kT} approximates the energy range over which the astrophysical Sfactor should be known. This product is the “Gamow peak”, and the location of the maximum is at E_{0}. This maximum and the width of the Gamow peak are defined by:
$$\begin{array}{cc}& {E}_{0}={\left[{\left(\frac{\pi}{\u0127}\right)}^{2}{({Z}_{0}{Z}_{1}{e}^{2})}^{2}\left(\frac{\mu}{2}\right){(kT)}^{2}\right]}^{1/3}=0.1220{\left({Z}_{0}^{2}{Z}_{1}^{2}\mu {T}_{9}^{2}\right)}^{1/3}\hfill \end{array}$$(5)
$$\begin{array}{cc}& \mathrm{\Delta}=\frac{4}{\sqrt{3}}\sqrt{{E}_{0}kT}=0.2368{\left({Z}_{0}^{2}{Z}_{1}^{2}\mu {T}_{9}^{5}\right)}^{1/6},\hfill \end{array}$$(6)
where Z_{0} and Z_{1} are the atomic numbers of the interacting nuclei, e is the elementary electric charge, and T_{9} is the temperature in 10^{9} K.
Many reactions of astrophysical importance proceed through nuclear resonances, populating compound nuclear states that subsequently decay. We only consider twobody reactions in the following. For a single, isolated resonance, the cross section in Eq. (1) can be replaced by
$$\begin{array}{c}\hfill \sigma (E)=\frac{2J+1}{(2{J}_{1}+1)(2{J}_{2}+1)}\frac{\pi}{{k}^{2}}\frac{{\mathrm{\Gamma}}_{a}(E){\mathrm{\Gamma}}_{b}(E)}{{(E{E}_{r})}^{2}+{\mathrm{\Gamma}}^{2}/4}.\end{array}$$(7)
Here, J, J_{1}, and J_{2} are the spins of the resonance, target, and projectile particles, respectively. This angular momentum term is often denoted by the symbol, ω. Γ_{a}(E) and Γ_{b}(E) are energydependent quantities describing the entrance and exit partial widths of the state in question. For example, for a (p, γ) reaction, Γ_{a}(E) corresponds to the proton partial width and Γ_{b}(E) is the γray partial width. Γ corresponds to the total width of the state, and E_{r} is the resonance energy. For a charged particle, Γ_{a}(E) is defined as
$$\begin{array}{c}\hfill {\mathrm{\Gamma}}_{a}(E)=2\frac{{\u0127}^{2}}{\mu {R}^{2}}{P}_{\ell}(E)\phantom{\rule{0.166667em}{0ex}}{C}^{2}S\phantom{\rule{0.166667em}{0ex}}{\theta}_{a}^{2},\end{array}$$(8)
where P_{ℓ}(E) is the energydependent penetration factor describing the probability of the particles tunnelling though the Coulomb barrier, C^{2}S is the product of the Isospin ClebshGordan coefficient for the interacting particles and spectroscopic factor. This latter quantity describes how well an excited state in the compound nucleus can be described by a singleparticle state. ${\theta}_{a}^{2}$ is the dimensionless singleparticle reduced width, which can be calculated theoretically from the particle’s wavefunction (see Iliadis 1997). We consider this latter quantity to be unity, here, since we are only considering relative effects under energy variations. The channel radius, R is calculated as $R={R}_{0}{A}_{p}^{1/3}+{A}_{t}^{1/3}$. This choice does not have a large effect on the calculations provided it is used consistently throughout.
If the resonance is sufficiently narrow such that its partial widths do not vary significantly over its width, we can assume that Γ_{a} and Γ_{b} are constant. In that case, the integral in Eq. (1) can be evaluated algebraically. Now the resonance cross section is replaced by a single quantity: the resonance strength, ωγ, which is defined by
$$\begin{array}{c}\hfill \omega \gamma =\omega \frac{{\mathrm{\Gamma}}_{p}{\mathrm{\Gamma}}_{\gamma}}{{\mathrm{\Gamma}}_{p}+{\mathrm{\Gamma}}_{\gamma}}.\end{array}$$(9)
It is important to recognise that the partial widths, Γ_{p} and Γ_{γ}, in Eq. (9) depend on resonance energies (through Eq. (8) for charged particles). Any energy shifts must be carefully propagated through these equations to determine shifted resonance strengths. Once this procedure is followed, Eq. (1) can be replaced with
$$\begin{array}{c}\hfill \u27e8\sigma v\u27e9={\left(\frac{2\pi}{\mu kT}\right)}^{3/2}{\u0127}^{2}{\displaystyle \sum _{i}\omega {\gamma}_{i}{e}^{{E}_{r}/kT}.}\end{array}$$(10)
This paper focuses on reaction rates dominated by resonances in the absence of interference. Addressing interfering resonance reaction rates that are best described by Rmatrix or other complex models is left for future work.
3. Monte Carlo reaction rates
3.1. Formalism
In order to investigate the influence that correlated resonance energy uncertainties have on reaction rates, the Monte Carlo reaction rate method first described in Longland et al. (2010) was used as a starting point. The general strategy laid out in that study was to first assign a probability density distribution to every uncertain input parameter to the reaction rate calculation. For the case of resonant cross sections, these uncertain parameters are the resonance energies, E_{r}, partial widths, Γ_{a}, and resonance strengths, ωγ, in Eqs. (7), (9), and (10). The strategy was extended to allow for uncertain spinparities in Mohr et al. (2014).
Once probability density distributions have been obtained for uncertain input parameters, sample parameters are randomly chosen from the distributions assuming that all parameters are independent (the validity of this assumption is discussed below). The reaction rate calculated from the sampled parameters represents a single sample rate. This procedure is repeated many times (10 000 is preferred but at least 3000 samples was found to produce stable results) to obtain a distribution of reaction rates that can be summarised with a reaction rate probability density distribution. Longland et al. (2010) found that the reaction rate probability density can often be summarised with a lognormal distribution with shape parameters μ and σ. The recommended rate is then given by:
$$\begin{array}{c}\hfill {\u27e8\sigma v\u27e9}_{\phantom{\rule{0.333333em}{0ex}}\text{rec.}}={e}^{\mu (T)}.\end{array}$$(11)
The “low” and “high” rates given by the 1 − σ uncertainties are found using
$$\begin{array}{c}\hfill {\u27e8\sigma v\u27e9}_{l}={e}^{\mu (T)}{e}^{\sigma (T)}\phantom{\rule{2em}{0ex}}{\u27e8\sigma v\u27e9}_{h}={e}^{\mu (T)}{e}^{+\sigma (T)}.\end{array}$$(12)
The procedure described above was found to be a flexible method for estimating the reaction rate uncertainties. Provided energy uncertainties are correctly propagated into partial width uncertainties, it works equally well for narrow and wide resonant reaction rates with large uncertainties. However, the method did not account for the case in which there are correlations between parameters.
3.2. Correlated energies
Recently, efforts were made to include the effects of correlated partial width and resonance strength uncertainties in the Monte Carlo method described above (Longland 2017). This case of correlated widths and resonance strengths needed to be investigated because resonances are usually normalised to some standard, well known resonance. Given that experimental reality, the assumption of independent variables made in Longland et al. (2010) is not strictly accurate.
Consider the following example: a partial width for resonance j, Γ_{j}, is correlated with a reference resonance, r. These partial widths carry factor uncertainties, f.u._{j} ≡ σ_{Γ, j}/Γ_{j} and f.u._{r} ≡ σ_{Γ, r}/Γ_{r}, respectively where σ_{Γ, j} is the uncertainty in Γ_{j}, etc. A single correlation parameter, ρ_{j}, can be used to describe the magnitude of their correlation:
$$\begin{array}{c}\hfill {\rho}_{j}=\frac{{\sigma}_{\mathrm{\Gamma},r}}{{\mathrm{\Gamma}}_{r}}\frac{{\mathrm{\Gamma}}_{j}}{{\sigma}_{\mathrm{\Gamma},j}}\equiv \frac{f.u{.}_{r}}{f.u{.}_{j}}.\end{array}$$(13)
During the Monte Carlo procedure, the following steps are taken: first, random, uncorrelated samples are produced for the reference resonance and resonance j. These samples, denoted x_{i} and y_{j, i}, are drawn from standard normal distributions (i.e., with a mean of zero and standard deviation of 1). Longland (2017) showed that a simple, 2step procedure can then be used to calculate correlated partial widths. Second, the samples are correlated using
$$\begin{array}{c}\hfill {y}_{j,i}^{\prime}={\rho}_{j}{x}_{i}+\sqrt{1{\rho}_{j}^{2}}{y}_{j,i},\end{array}$$(14)
where ${y}_{j,i}^{\prime}$ are the correlated samples for resonance j. Finally, the partial width samples can be computed with the knowledge that partial widths should be lognormally distributed (see Longland et al. 2010):
$$\begin{array}{c}\hfill {\mathrm{\Gamma}}_{j,i}={\mathrm{\Gamma}}_{j,\phantom{\rule{0.333333em}{0ex}}\text{rec.}}\phantom{\rule{0.166667em}{0ex}}{(f.u.)}^{{y}_{j,i}^{\prime}}.\end{array}$$(15)
Here, the recommended partial width for resonance j is Γ_{j, rec.}.
Note that the probability density distributions defined by Eq. (15) are a lognormal distribution. For small uncertainties they resemble Gaussian distributions, but are not defined for negative values. They therefore describe physical parameters such as partial widths and resonance strengths appropriately. They are also well motivated by the central limit theorem, as discussed in detail in Longland et al. (2010).
In the present study we follow the same strategy as above, but must consider two modifications: (i) energy uncertainties must be propagated through partial widths. For the charged particle reactions considered here, this is accomplished using Eq. (8); and (ii) energy uncertainties are not lognormal. Often large resonance energy uncertainties arise because of large uncertainties in reaction Qvalues. For example, imagine a lowenergy resonance with a large energy uncertainty that has an appreciable probability of being a subthreshold resonance. Modification (i) is already accounted for in the RatesMC code implementation of the Monte Carlo methods outlined in Longland et al. (2010)^{1} The second case requires modifications to Eqs. (13) and (15):
$$\begin{array}{cc}& {\rho}_{j}=\frac{min({\sigma}_{E})}{{\sigma}_{E,j}}.\hfill \end{array}$$(16)
$$\begin{array}{cc}& {E}_{j,i}={E}_{j,\phantom{\rule{0.333333em}{0ex}}\text{rec.}}+{y}_{j,i}^{\prime}{\sigma}_{E,j}.\hfill \end{array}$$(17)
The correlated standard normal samples, ${y}_{j,i}^{\prime}$ are computed using Eq. (14). The purpose of Eq. (16) is similar to that of Eq. (13). It ensures that correlations between resonance energies do not exceed their experimental limits. For example, consider a reaction in which the majority of the resonance energy uncertainty arises from an uncertain reaction Qvalue, σ_{Q}. Those resonance energies will be correlated with ρ ≈ 1. Now assume another resonance in that reaction populates an excited state that also has an uncertain excitation energy, σ_{Ex}. The energy of that resonance will have a larger overall uncertainty (given by ${\sigma}^{2}={\sigma}_{Q}^{2}+{\sigma}_{Ex}^{2}$) and will not be fully correlated with the other resonance energies. For this example resonance, ρ < 1.
Finally, cases in which resonance energies are determined through different means must be considered. In most cases, high energy resonance energies are determined through direct measurement of the reaction cross section. Resonance energies can be determined by measuring a yield curve, for example. Lowenergy resonances in the same reaction may be determined from excitation energies and an uncertain Qvalue. Therefore, the ability to enable or disable energy correlations for a resonance must be available. The RatesMC code has been modified to allow this.
4. Test cases
To investigate the effect of correlated energy uncertainties on reaction rates, two fictional reactions designed to probe low and highdensity resonance regimes are first investigated. The effect on actual reaction rates will be detailed in Sect. 5.
4.1. Low resonance density
The first test case considered is shown in Fig. 1. The physical system consists of a fictional reaction whose cross section consists entirely of three isolated, narrow resonances. Thus, Eq. (10) can be used to calculate the reaction rate. The partial widths are calculated using Eq. (8) by assuming R_{0} = 1.25 fm and a (p, γ) reaction occurring on an isotope with an atomic number of 12 and mass of 23. The resonance parameters chosen are shown in Table 1.
Fig. 1. Test cases for our investigation of the impact of energy uncertainty correlations on nuclear reaction rate uncertainties. Shown as a solid line is the Gamow peak in arbitrary units calculated from Eq. (10) at three temperatures: 0.08 GK, 0.13 GK, and 0.2 GK in panels a, b, and c respectively. The reaction rate for each of the three narrow resonances described in Table 1 are shown as coloured points. The ydirection is scaled arbitrarily for clarity to highlight which resonances contribute most to the reaction rate and where they are located in comparison to the Gamow peak. This information is also displayed by the bar on the right of each panel. For example, in panel a, only one resonance – the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=100$ keV resonance shown in green – contributes to the reaction rate. In panel b, two resonances – ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=100$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV in green and orange – contribute approximately 50% each to the reaction rate at 0.13 GK. 
Parameters of the three resonances considered in our fictional test case.
These resonances contribute different amounts to the reaction rate depending on the temperature. The lowenergy chargedparticle reaction resonances most important at low temperatures can be well predicted by the Gamow peak defined in Eq. (4). Indeed, this is the case as shown in Fig. 1, where the three resonance contributions to the total reaction rate follow closely the progression of the Gamow peak as temperatures increase. Higher energy resonances, though, will not follow this pattern. Their cross sections become constrained by γray partial widths, which does not exhibit the Coulomb barrier energy dependence (Newton et al. 2007).
First, consider the reaction rate at T = 0.13 GK. At this temperature, the two resonances at ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=100$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV contribute approximately equally to the total reaction rate as shown by the bar on the right of the centre panel. They are situated either side of the maximum of the Gamow peak. To investigate the effect that their energy uncertainties have on the reaction rate, their resonance energies are varied over a given range. For each trial resonance energy, the proton partial width is recalculated using Eq. (8) and the parameters in Table 1. Using this, the reaction rate is determined using Eqs. (9) and (10). These resonance energies are varied using two schemes: (i) correlated energies, so any increase in the energy of one resonance corresponds to an equal increase in the other, and (ii) anticorrelated energies, in which any energy increase in one resonance energy corresponds to an equal magnitude decrease in the other. The resonance energies are varied by ±60 keV this way. The variations affect the individual contributions to the rate as well as the total rate. These are shown in Fig. 2.
Fig. 2. Effects of correlated and anticorrelated energy variations to a pair of resonances on the calculated reaction rate at 0.13 GK. The green and orange dashed lines show the individual resonance contributions from the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=100$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV resonances to the total reaction rate shown by the black line. The energy shift is defined as the shift of the lower energy resonance. All are normalised to their recommended values at ΔE = 0. 
Figure 2 shows that larger reaction rate variations are expected if resonance energies are anticorrelated. To understand this effect, consider first the correlated energy case in the lefthand panel as well as the middle panel of Fig. 1. As the resonances both increase in energy, the one at ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=100$ keV shifts closer to the maximum of the Gamow peak, thus increasing its contribution to the total rate. Conversely, the resonance at ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV moves away from the maximum of the Gamow peak and decreases its contribution. The net effect is that the total rate does increase, but the magnitude is weakened by the opposite contributions of the two resonances. As the resonance energies decrease, a similar effect is apparent: the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=100$ keV resonance contributes less while the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV resonance contributes more, resulting again in an increase in reaction rate that is weakened by the opposite contributions.
In the case of anticorrelated resonance energies in the righthand panel of Fig. 2, the resonances contributions work in tandem. When the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=100$ keV resonance energy is decreased, it moves away from the maximum of the Gamow peak to lower energies, while the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV resonance also moves away, but to higher energies. Thus the contributions of both resonances decrease, resulting in a reduced total reaction rate. Similarly, as one moves towards the maximum of the Gamow peak, so will the other, resulting in a strengthened increase in the total rate. Anticorrelated energy uncertainties in this case result in an increased reaction rate uncertainty.
The anticorrelated resonance energy example discussed above will rarely occur in experimental resonance energy measurements. However, if the resonance energies are treated as completely uncorrelated during the Monte Carlo procedure, their relative variations will be somewhere between the fully correlated and fully anticorrelated cases. Thus, taking into account the effect illustrated in Fig. 2 we expect larger uncertainties for uncorrelated resonance energies than correlated energies. This is indeed the case, as shown in Fig. 3, which shows the reaction rate probability density distributions for these two cases. In grey, a broad, approximately Gaussian peak represents the probability density distribution for the uncorrelated energy case. In the correlated case, the effect shown in Fig. 2 is clear in that the probability distribution is not only narrower, but is also highly skewed. Since we chose a temperature at which the resonances are either side of the Gamow peak and contribute approximately equally to the total rate, the reaction rate can only increase as their energies are varied.
Fig. 3. Reaction rate probability distributions from the Monte Carlo variation of resonance energies at 0.13 GK. Shown in grey is the uncorrelated case, in which the resonance energies are allowed to vary independently. In blue is the case in which the resonance energies are fully correlated. The distribution becomes narrower and highly skewed in this case owing to the effect illustrated in Fig. 2. 
This effect is not universal. In the example above, a specific temperature was chosen to correspond to two resonances either side of the Gamow peak. To investigate more possibilities, the third panel in Fig. 1 is illustrative. In this case, the resonances at ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=300$ keV contribute about 75% and 25% to the reaction rate, respectively. They are also both located above the maximum of the Gamow peak, so as we shift their energies in a correlated manner, they will both shift towards or away from the peak in unison. The effect of their variations on the total reaction rate is shown in Fig. 4. In this case, the opposite effect to that observed in Fig. 2 is apparent. If the resonances are correlated, they move together to reduce or increase the reaction rate. If they are anticorrelated, their contributions essentially cancel out to produce very little variation in the total reaction rate.
Fig. 4. Effects of correlated and anticorrelated resonance energies on the reaction rate at T = 0.2 GK. The red and blue dashed lines correspond to the contributions from the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=300$ keV resonances, respectively. See Fig. 2 for details. 
The Monte Carlo reaction rate comparison at T = 0.2 GK is shown in Fig. 5. In this case, we see the opposite effect to the example at T = 0.13 GK. If energy correlations are taken into account, the reaction rate uncertainties increase. Clearly, these effects are hard to predict, particularly when large resonance energies are concerned. However, using Monte Carlo uncertainty propagation, we are able to account for the effect of energy correlations on the reaction rate uncertainties.
Fig. 5. Reaction rate probability distributions from the Monte Carlo variation of resonance energies at 0.2 GK. See Fig. 3 for details. The probability distribution becomes wider when correlated energies are considered in this case owing to the effect illustrated in Fig. 4. 
4.2. High resonance density
Following the same procedure for a high density of resonances produces results that are easier to predict. In this case, the reaction rate uncertainty when correlated energy uncertainties are taken into account reliably decreases. Now, resonance are placed at energies between ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=50$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=400$ keV with a spacing of 20 keV.
The reaction rate uncertainty is shown in Fig. 6 for T = 0.13 GK. The uncertainty for correlated energies (blue) clearly decreases in comparison to the uncorrelated case (grey). As temperature increases we find that the effect of correlations decreases because more resonances contribute to the rate.
Fig. 6. Reaction rate probability distributions from the Monte Carlo variation of resonance energies at 0.13 GK for the high resonance density case. See Fig. 3 for details. 
4.3. Discussion
The results discussed above are difficult to predict a priori. However, a conservative estimate of reaction rate uncertainties would only be concerned with the case in which resonance energy correlations increase the reaction rate uncertainties. This case is shown in Fig. 5. For these cases to occur, the resonances should be located on the same side of the Gamow peak (i.e., as their energy increases or decreases, correlations would cause them all to increase or decrease their contributions to the total rate in unison). This requires more than one resonance contributing to the rate, but the resonance density cannot be too high, else resonances would be distributed throughout the Gamow peak. For example, it would not be realistic to expect a case where a high resonance density is found on one side of the Gamow peak, but no resonances on the other. Furthermore, the resonance energy uncertainty in contributing resonances should be large in comparison to the Gamow peak’s width, defined by Eq. (6).
5. Physical cases
Now that the general behaviour of how reaction rate uncertainties change when resonance energy correlations are taken into account, some physically realistic cases will be considered. These are the ^{35}Ar(p, γ)^{36}K reaction; the ^{39}Ca(p, γ)^{40}Sc reaction; and the ^{13}N(α, p)^{16}O reaction. They span a range of resonance densities and represent cases where the resonance energies are known to be uncertain and correlated.
5.1. The ^{35}Ar(p, γ)^{36}K reaction
The ^{35}Ar(p, γ)^{36}K reaction is a key reaction in explosive hydrogen burning. In Xray bursts this reaction is expected to occur faster than its competing βdecay, but in novae (i.e., lower temperatures) the rate is less well known. Its effect on the nucleosynthesis of heavier elements is not well understood (Glasner & Truran 2009). The reaction rate was evaluated in Iliadis et al. (1999). At that time, the Qvalue of this reaction was poorly known (Audi & Wapstra 1995), leading to large, correlated resonance energy uncertainties. Since that time the excitation energy uncertainties in ^{36}K have been dramatically reduced by an order of magnitude by Wrede et al. (2010). However, the resonance energies are still expected to be correlated and the resonance density of this reaction is very low with just four known resonances below 1 MeV. For these reasons it is an ideal case with which to investigate correlated resonance energies in the Monte Carlo framework. The resonance parameters from Wrede et al. (2010) are listed in Table 2. Note that we have assigned very small uncertainties (1%) to the partial widths so that the resonance energy effects can be clearly identified. Separate calculations confirm that the uncertainties due to other sources sum quadratically, as expected.
Resonance parameters for the ^{35}Ar(p, γ)^{36}K reaction taken from Wrede et al. (2010).
The reaction rate uncertainties for the ^{35}Ar(p, γ)^{36}K reaction assuming the resonance parameters shown in Table 2 are shown in Fig. 7. The coloured contour represents the reaction rate uncertainties arising from correlated energy uncertainties, with thick and thin lines representing the 1σ and 2σ uncertainty bands, respectively. These rates have been normalised to the recommended (median) rate, which is shown by a horizontal line at unity. The blue lines show the reaction rate uncertainties when resonance energy correlations are not taken into account. The thick blue line represents the median rate, and the thin dashed blue lines represent the 1σ uncertainties. This figure shows that over most of the temperature range, energy correlations do not strongly affect the reaction rate uncertainties. They are only slightly smaller when taking resonance energy correlations into account. This is mostly due to the fact that below 200 MK and above 1 GK, only one resonance is contributing and the effect of energy correlation is almost inexistent. In between these temperatures two resonances contribute to the reaction rate. At 400 MK, for example, the resonances are located either side of the maximum of the Gamow peak. In this case the effect of correlations is small, but in line with the case described in Sect. 4.1: as the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=260$ keV resonance moves to a lower energy, for example, the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=623$ keV resonance also moves to a lower energy, thus reducing the impact of resonance energies on the reaction rate.
Fig. 7. Reaction rate uncertainties for the ^{35}Ar(p, γ)^{36}K reaction assuming the resonance parameters shown in Table 2. Recall that only resonance energies are taken into account. The rate has been normalised to the median rate, which is shown by a dashed line at unity. The thick and thin black lines represent the 1σ and 2σ uncertainties in the correlated resonance energy calculation. The blue lines show the uncorrelated case, with the thick line representing the median rate (again, normalised to the recommended correlated energy rate), and dashed lines showing the 1σ uncertainties. At T = 400 MK, the effect of correlations slightly decreases the reaction rate uncertainty. 
These calculations were also performed assuming the (obsolete) resonance parameters reported in Iliadis et al. (1999). Although the total rate uncertainty is much larger owing to the larger energy uncertainties, the correlations between those energies have the same, minor, effect on the uncertainty as outlined above.
5.2. The ^{39}Ca(p, γ)^{40}Sc reaction
The ^{39}Ca(p, γ)^{40}Sc reaction is also important in explosive nucleosynthesis. It influences the endpoint of the rpprocess in Xray bursts, and has also been evaluated in Iliadis et al. (1999). The resonance density is similar to the ^{35}Ar(p, γ)^{36}K reaction, but in this case the Qvalue is better known. The resonance energy uncertainties are just 56 keV, as shown in Table 3.
Resonance parameters for the ^{39}Ca(p, γ)^{40}Sc reaction taken from Iliadis et al. (1999).
In this case, the predictions in Sect. 4.3 indicate that there should, indeed, be an effect of resonance energy correlations on the reaction rate. The average resonance energy separation is 300 keV compared with Δ = 200 keV at 400 MK. In contrast to the ^{35}Ar(p, γ)^{36}K reaction, though, there is a temperature range at which both resonances at ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=223$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=353$ keV lay on the lowenergy side of the Gamow peak. Thus we expect an effect of correlated energies on the reaction rate.
Figure 8 shows that correlations between resonance energy uncertainties do, indeed, affect the reaction rate uncertainty strongly at 300–400 MK where both resonances at ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=233$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=353$ keV contribute to the reaction rate and are both are on the lowenergy side of the Gamow peak. In this particular scenario, as the resonances comove to lower energies, both contribute to a lower reaction rate. Conversely if they both comove to higher resonance energy they both contribute to a higher reaction rate. In the uncorrelated energy case, this scenario is more rare, thus the rate uncertainty is smaller. Note that over most of the temperature range, 0.2 < T_{9} < 1, the rate uncertainties are considerable. At those temperatures, only one or two resonances ever contributes towards the rate, and the uncertainty is arising from the strong Coulomb barrier energy dependence in Eq. (8).
Fig. 8. Reaction rate uncertainties for the ^{39}Ca(p, γ)^{40}Sc reaction assuming the resonance parameters shown in Table 3. Only uncertainties in the resonance energies are accounted for. See Fig. 7 for description. 
5.3. The ^{13}N(α, p)^{16}O reaction
The ^{13}N(α, p)^{16}O reaction affects nitrogen production in supernova explosions as the shockwave passes through the outer regions of the exploding star (Pignatari et al. 2013). That material can eventually go on to form presolar grains, whose isotopic nitrogen ratios provide a precise test of astrophysical models (Zinner 2014). The ^{13}N(α, p)^{16}O reaction rate should be known, therefore, to a high precision. The rate was recently evaluated in Meyer et al. (2020). In this case the ^{13}N+α threshold (S_{α+13N} = 5817.4 (4) keV) is accurately known from Tilley et al. (1993). However the energy of the resonances in the compound nucleus ^{17}F suffers from systematic uncertainty of several tens of keV which introduces a strong correlation between resonance energies. This uncertainty originates from a possible error in the calibration of one of the magnets used during the measurement of the excitation functions of the ^{16}O(p, p)^{16}O reaction by Salisbury et al. (1962), Salisbury & Richards (1962), and the ^{16}O(p, p′)^{16}O and ^{16}O(p, α)^{13}13 reactions by Dangle et al. (1964). The properties of the most influential resonances with αparticle partial width determined either from direct measurements or mirror symmetry considerations are taken from Meyer et al. (2020) and summarised in Table 4.
Resonance parameters for the ^{13}N(α, p)^{16}O reaction taken from Meyer et al. (2020).
In order to emphasise the effect of the (un)correlated uncertainties on the resonance energy a very small uncertainty (1%) has been assigned to the partial widths and the tentative spin and parity have been considered as firmed assignment. In this particular case, the resonances have very large (factor of 2.5) uncertainties, which completely dominates the energy uncertainties under investigation here.
Even though some resonances have large total widths their number is relatively small and they can be considered as isolated. The case of low resonance density discussed in Sect. 4.1 should then apply and an effect of correlated energies on the reaction rate is then expected. This is indeed the case as shown in Fig. 9 where the ^{13}N(α, p)^{16}O reaction rate uncertainties are presented.
Fig. 9. Reaction rate uncertainties for the ^{13}N(α, p)^{16}O reaction assuming the resonance parameters shown in Table 4. Only uncertainties in the resonance energies are accounted for. See Fig. 7 for description. 
Let us first consider a typical temperature of 1 GK which corresponds to a Gamow peak with a maximum of about 1 MeV and width of Δ ≈ 700 keV from Eq. (6). In this case the reaction rate is dominated by the two resonances at ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=741$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=1213$ keV which are siting on each side of the Gamow peak’s maximum. As these two resonances comove when their energy uncertainties are correlated, one resonance will have an increased contribution to the reaction rate while the other one’s contribution will decrease. A smaller reaction rate uncertainty is therefore observed for the correlated case (black solid line in Fig. 9) with respect to the uncorrelated case (dashed blue line in Fig. 9) between 0.8 GK and 2 GK, which is in line with the findings presented in Fig 3.
The opposite behaviour is observed for temperatures lower than 400 MK and greater than 6 GK where the uncorrelated case gives smaller reaction rate uncertainties than in the correlated case. For these temperatures all resonances reported in Table 4 are on one side of the Gamow peak, for example, at higher energies for 400 MK and lower energies for 6 GK. The reaction rate therefore spans a larger range in the correlated case inducing a greater rate uncertainty than in the uncorrelated case as presented in Fig. 5. At temperatures below about 600 MK, the rate uncertainties are large. This is because at these temperatures, the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=741$ keV resonance dominates the rate and is located on the lowenergy tail of the Gamow peak. Any variation in the resonance energy has a large effect on the rate through Eq. (8).
6. Summary and conclusions
Monte Carlo methods can be a powerful tool for computing statistically rigorous uncertainties of thermonuclear reaction rates. While the methods have been in use for some time, no effort had previously been made to account for correlations between resonance energies. These effects are particularly important for radioactive nuclei where resonances are not often directly measured.
In this paper, we expanded on the correlation scheme developed in Longland (2017) to allow for correlations between resonance energies. We found that the effects are not necessarily easy to predict. Reactions rates dominated by many resonances are not strongly affected by correlations, whereas those dominated by only a single resonance at astrophysically important temperatures are also not significantly affected. The cases that matter most are those where multiple resonances contributed to the reaction rate. This effect is enhanced if they are on the same side of the Gamow peak’s maximum value.
Correlations between resonance parameters can be an important effect in thermonuclear reaction rate calculations. The correlation of resonance energies was previously unexplored, which has now been accounted for in this work. Since the effects of these correlations are rather unpredictable, we recommend that any reaction rate uncertainty calculation be carefully checked to ensure corrections due to resonance energy correlations do not significantly affect the results. This will be of particular importance for reactions on isotopes far from stability, where the energies of excited states can carry large, correlated uncertainties because they are determined from uncertain reaction Qvalues.
Acknowledgments
This material is based partly upon work supported by the US Department of Energy, Office of Science, Office of Nuclear Physics, under Award Number DESC0017799 and under Contract No. DEFG0297ER41041.
References
 Angulo, C., Arnould, M., Rayet, M., et al. 1999, Nucl. Phys. A, 656, 3 [NASA ADS] [CrossRef] [Google Scholar]
 Audi, G., & Wapstra, A. H. 1995, Nucl. Phys. A, 595, 409 [NASA ADS] [CrossRef] [Google Scholar]
 Cescutti, G., Hirschi, R., Nishimura, N., et al. 2018, MNRAS, 478, 4101 [NASA ADS] [CrossRef] [Google Scholar]
 Dangle, R. L., Oppliger, L. D., & Hardie, G. 1964, Phys. Rev., 133, 647 [Google Scholar]
 Fields, C. E., Farmer, R., Petermann, I., Iliadis, C., & Timmes, F. X. 2016, ApJ, 823, 46 [NASA ADS] [CrossRef] [Google Scholar]
 Fields, C. E., Timmes, F. X., Farmer, R., et al. 2018, ApJS, 234, 19 [NASA ADS] [CrossRef] [Google Scholar]
 Glasner, S. A., & Truran, J. W. 2009, ApJ, 692, L58 [NASA ADS] [CrossRef] [Google Scholar]
 Iliadis, C. 1997, Nucl. Phys. A, 618, 166 [NASA ADS] [CrossRef] [Google Scholar]
 Iliadis, C. 2015, Nuclear Physics of Stars (John Wiley & Sons) [Google Scholar]
 Iliadis, C., Endt, P. M., Prantzos, N., & Thompson, W. J. 1999, ApJ, 524, 434 [CrossRef] [Google Scholar]
 Iliadis, C., D’Auria, J. M., Starrfield, S., Thompson, W. J., & Wiescher, M. 2001, ApJS, 134, 151 [NASA ADS] [CrossRef] [Google Scholar]
 Longland, R. 2012, A&A, 548, A30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Longland, R. 2017, A&A, 604, A34 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Longland, R., Iliadis, C., Champagne, A. E., et al. 2010, Nucl. Phys. A, 841, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Meyer, A., de Séréville, N., Laird, A. M., et al. 2020, Phys. Rev. C, 102, 035803 [CrossRef] [Google Scholar]
 Mohr, P., Longland, R., & Iliadis, C. 2014, Phys. Rev. C, 90, 065806 [NASA ADS] [CrossRef] [Google Scholar]
 Newton, J. R., Iliadis, C., Champagne, A. E., et al. 2007, Phys. Rev. C, 75, 045801 [CrossRef] [Google Scholar]
 Pignatari, M., Wiescher, M., Timmes, F. X., et al. 2013, ApJ, 767, L22 [NASA ADS] [CrossRef] [Google Scholar]
 Rauscher, T. 2012, ApJS, 201, 26 [NASA ADS] [CrossRef] [Google Scholar]
 Rauscher, T., Nishimura, N., Hirschi, R., et al. 2016, MNRAS, 463, 4153 [NASA ADS] [CrossRef] [Google Scholar]
 Rolfs, C. E., & Rodney, W. S. 1988, Cauldrons in the Cosmos: Nuclear Astrophysics (University of Chicago Press) [Google Scholar]
 Salisbury, S. R., & Richards, H. T. 1962, Phys. Rev., 126, 2147 [CrossRef] [Google Scholar]
 Salisbury, S. R., Hardie, G., Oppliger, L., & Dangle, R. 1962, Phys. Rev., 126, 2143 [CrossRef] [Google Scholar]
 Sallaska, A. L., Iliadis, C., Champange, A. E., et al. 2013, ApJS, 207, 18 [NASA ADS] [CrossRef] [Google Scholar]
 Schulz, N., Alford, W. P., & Jamshidi, A. 1971, Nucl. Phys. A, 162, 349 [NASA ADS] [CrossRef] [Google Scholar]
 Thompson, W. J., & Iliadis, C. 1999, Nucl. Phys. A, 647, 259 [NASA ADS] [CrossRef] [Google Scholar]
 Tilley, D. R., Weller, H. R., & Cheves, C. M. 1993, Nucl. Phys. A, 564, 1 [CrossRef] [Google Scholar]
 Wrede, C., Clark, J. A., Deibel, C. M., et al. 2010, Phys. Rev. C, 82, 035805 [CrossRef] [Google Scholar]
 Zinner, E. 2014, Presolar Grains (Elsevier), 1, 181 [Google Scholar]
All Tables
Resonance parameters for the ^{35}Ar(p, γ)^{36}K reaction taken from Wrede et al. (2010).
Resonance parameters for the ^{39}Ca(p, γ)^{40}Sc reaction taken from Iliadis et al. (1999).
Resonance parameters for the ^{13}N(α, p)^{16}O reaction taken from Meyer et al. (2020).
All Figures
Fig. 1. Test cases for our investigation of the impact of energy uncertainty correlations on nuclear reaction rate uncertainties. Shown as a solid line is the Gamow peak in arbitrary units calculated from Eq. (10) at three temperatures: 0.08 GK, 0.13 GK, and 0.2 GK in panels a, b, and c respectively. The reaction rate for each of the three narrow resonances described in Table 1 are shown as coloured points. The ydirection is scaled arbitrarily for clarity to highlight which resonances contribute most to the reaction rate and where they are located in comparison to the Gamow peak. This information is also displayed by the bar on the right of each panel. For example, in panel a, only one resonance – the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=100$ keV resonance shown in green – contributes to the reaction rate. In panel b, two resonances – ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=100$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV in green and orange – contribute approximately 50% each to the reaction rate at 0.13 GK. 

In the text 
Fig. 2. Effects of correlated and anticorrelated energy variations to a pair of resonances on the calculated reaction rate at 0.13 GK. The green and orange dashed lines show the individual resonance contributions from the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=100$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV resonances to the total reaction rate shown by the black line. The energy shift is defined as the shift of the lower energy resonance. All are normalised to their recommended values at ΔE = 0. 

In the text 
Fig. 3. Reaction rate probability distributions from the Monte Carlo variation of resonance energies at 0.13 GK. Shown in grey is the uncorrelated case, in which the resonance energies are allowed to vary independently. In blue is the case in which the resonance energies are fully correlated. The distribution becomes narrower and highly skewed in this case owing to the effect illustrated in Fig. 2. 

In the text 
Fig. 4. Effects of correlated and anticorrelated resonance energies on the reaction rate at T = 0.2 GK. The red and blue dashed lines correspond to the contributions from the ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=250$ keV and ${E}_{r}^{\phantom{\rule{0.333333em}{0ex}}\text{c.m.}}=300$ keV resonances, respectively. See Fig. 2 for details. 

In the text 
Fig. 5. Reaction rate probability distributions from the Monte Carlo variation of resonance energies at 0.2 GK. See Fig. 3 for details. The probability distribution becomes wider when correlated energies are considered in this case owing to the effect illustrated in Fig. 4. 

In the text 
Fig. 6. Reaction rate probability distributions from the Monte Carlo variation of resonance energies at 0.13 GK for the high resonance density case. See Fig. 3 for details. 

In the text 
Fig. 7. Reaction rate uncertainties for the ^{35}Ar(p, γ)^{36}K reaction assuming the resonance parameters shown in Table 2. Recall that only resonance energies are taken into account. The rate has been normalised to the median rate, which is shown by a dashed line at unity. The thick and thin black lines represent the 1σ and 2σ uncertainties in the correlated resonance energy calculation. The blue lines show the uncorrelated case, with the thick line representing the median rate (again, normalised to the recommended correlated energy rate), and dashed lines showing the 1σ uncertainties. At T = 400 MK, the effect of correlations slightly decreases the reaction rate uncertainty. 

In the text 
Fig. 8. Reaction rate uncertainties for the ^{39}Ca(p, γ)^{40}Sc reaction assuming the resonance parameters shown in Table 3. Only uncertainties in the resonance energies are accounted for. See Fig. 7 for description. 

In the text 
Fig. 9. Reaction rate uncertainties for the ^{13}N(α, p)^{16}O reaction assuming the resonance parameters shown in Table 4. Only uncertainties in the resonance energies are accounted for. See Fig. 7 for description. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.