On Trapped Modes In Variable White Dwarfs As Probes Of The $^{12}$C($\alpha, \gamma$)$^{16}$O Reaction Rate

We seek signatures of the current experimental $^{12}$C$(\alpha,\gamma)^{16}$O reaction rate probability distribution function in the pulsation periods of carbon-oxygen white dwarf models. We find that adiabatic g-modes trapped by the interior carbon-rich layer offer potentially useful signatures of this reaction rate probability distribution function. Probing the carbon-rich region is relevant because it forms during the evolution of low-mass stars under radiative helium burning conditions, mitigating the impact of convective mixing processes. We make direct quantitative connections between the pulsation periods of the identified trapped g-modes in variable WD models and the current experimental $^{12}$C$(\alpha,\gamma)^{16}$O reaction rate probability distribution function. We find an average spread in relative period shifts of $\Delta P/P \simeq \pm$ 2\% for the identified trapped g-modes over the $\pm$ 3$\sigma$ uncertainty in the $^{12}$C$(\alpha,\gamma)^{16}$O reaction rate probability distribution function -- across the effective temperature range of observed DAV and DBV white dwarfs and for different white dwarf masses, helium shell masses, and hydrogen shell masses. The g-mode pulsation periods of observed white dwarfs are typically given to 6-7 significant figures of precision. This suggests that an astrophysical constraint on the $^{12}$C$(\alpha,\gamma)^{16}$O reaction rate could, in principle, be extractable from the period spectrum of observed variable white dwarfs.


INTRODUCTION
Tens of thousands of nuclear reactions can participate during the evolution of a star, depending on the environmental conditions. Only a few of these reactions have a strong impact on the overall chemical evolution of the elements, with significant consequences for the chemistry and the subsequent molecular evolution of baryonic matter. In particular, the helium burning 12 C(α,γ) 16 O reaction plays a major role in the energy production and nucleosynthesis of stars (e.g., Iben 1967;Fowler 1984a;Arnett 1996;Iliadis 2015;deBoer et al. 2017) and thus influences the 12 C/ 16 O ratio in the universe.
The difficulty in measuring the 12 C(α,γ) 16 O rate in nuclear experiments is due to the small cross section of the 12 C(α,γ) 16 O reaction at temperatures relevant for helium burning in stars. Nuclear experiments in terrestrial laboratories provide data for energies as low as about 2 MeV, with extrapolation to stellar conditions at kT 20 keV. At stellar conditions, two partial waves contribute, denoted by their spectroscopic E1 and E2 amplitudes in reaction and scattering theory (Fowler 1984b). The challenge is measurement of the low energy angular distributions of the 12 C(α,γ) 16 O reaction, from which the E1 and E2 cross sections are extracted.
Decreasing the uncertainty in the 12 C(α,γ) 16 O reaction rate from low-energy nuclear experiments has markedly improved in recent years. For example, to obtain a comprehensive evaluation deBoer et al. (2017) considered the entirety of existing experimental data, aggregating about 60 years of experimental data consisting of more than 50 independent experimental studies. More than 10,000 data points were then incorporated into a complete multi-channel phenomenological R-matrix analysis. A main result was the characterization of the uncertainty in the reaction rate, which was accomplished through a Monte-Carlo uncertainty analysis of the data and the extrapolation to low-energy using the R-matrix model. After finding an approximately Gaussian underlying probability distribution for the rate, there was statistical significance with the 1σ uncertainty of the reaction rate. A goal of forthcoming experiments is to further reduce the uncertainty in the 12 C(α,γ) 16 O reaction rate (e.g., Smith et al. 2021).
Partnering with this laboratory astrophysics quest are astrophysical constraints on the 12 C(α,γ) 16 O reaction rate. For example, Weaver & Woosley (1993) find a permissible range of the reaction rate by requiring the integrated nucleosynthesis yields, from a set of massive star explosion models over plausible initial mass functions, agree with the observed solar abundances for the intermediate mass isotopes. They found this range was insensitive to the assumed slope of the initial mass function within observational limits, and relatively insensitive to some details of convective boundary mixing.
As another example of an astrophysical constraint, models for the evolution of single stars predicts the existence of a gap in the black hole mass distribution for high mass stars due to the high temperature effects of electron-positron pair production (Heger & Woosley 2002). The location of the black hole mass gap is generally robust with respect to model uncertainties (Takahashi 2018;Farmer et al. 2019;Marchant & Moriya 2020), but depends sensitively on the uncertain 12 C(α,γ) 16 O reaction rate (Farmer et al. 2020). The location of the black hole mass gap, probed though LIGO/Virgo/Karga (LIGO Scientific Collaboration et al. 2015;Acernese et al. 2015;Akutsu et al. 2021) gravitational wave determinations of the masses and spins of merging binary black holes, thus allows a constraint on the 12 C(α,γ) 16 O reaction rate (Farmer et al. 2020;Woosley & Heger 2021;Mehta et al. 2022).
Seismology of hydrogen dominated atmosphere (DA class) and helium dominated atmosphere (DB class) carbon-oxygen white dwarfs (CO WDs) has also been used to place constraints on the 12 C(α,γ) 16 O reaction rate (Metcalfe et al. 2001(Metcalfe et al. , 2002. In this approach the central 16 O abundance and the location of the oxygento-carbon transition within the CO core are used as fitting parameters when minimizing the difference between the observed pulsations periods of specific WDs and the pulsation periods of WD models. The derived parameters can then imply a constraint on the 12 C(α,γ) 16 O reaction rate, but are sensitive to model choices (Metcalfe 2003), diffusion (Fontaine & Brassard 2002), and convective boundary mixing (Straniero et al. 2003).
De Gerónimo et al. (2017) analyzed two CO WD models with masses 0.548 and 0.837 M derived from evolutionary calculations from the zero-age main sequence (ZAMS). They considered models that varied the number of thermal pulses, the amount of overshooting, and the 12 C(α,γ) 16 O reaction rate within the uncertainties known at the time. They found that independent variations of these quantities produced significant changes in the resulting DA WD chemical profiles and the pulsation period spectrum. Pepper et al. (2022) calculated evolutionary models with initial masses in the range of 0.90≤ M/M ≤3.05. They considered different 12 C(α,γ) 16 O reaction rates within the uncertainties known at the time. As expected, they found no changes in the evolution prior to the core He-burning stage. However, they found that the subsequent stages of evolution produced differences in the convective He core mass, the number of thermal pulses during the asymptotic giant branch phase of evolution, and broad trends in the chemical profiles.
The main novelty of this article is a new search for potential signatures of the current experimental 12 C(α, γ) 16 O reaction rate probability distribution function in the pulsation periods of CO WD models. Section 2 describes our models, Section 3 describes our results, Section 4 discusses sensitivities, and Section 5 summarizes the results of our new search.  2017) tabulated reaction rate for 12 C(α, γ) 16 O to a finer temperature grid, from 52 to 2015 temperature points, to ensure that no temperature step results in variations in the rate of more than a factor of two. The recalculations also provided the formal ± 3 σ uncertainties on the experimental reaction rate probability distribution function in steps of 0.5σ. σ 0 is the median rate consistent with an astrophysical S-factor of S(300 keV) = 140 keV b with a ±1 σ = 21 keV b uncertainty. By exploring ±3 σ we effectively explore the range S(300 keV) = (77,203) keV b, where positive and negative σ indicate a stronger and weaker rate than the median value, respectively. This probability distribution function is shown in Figure 1 over the ± 3 σ region. Green bands represent positive σ i while the gray bands represent negative σ i . The 13 individual σ curves are depicted by respective dotted lines, with the ±1, 2, 3σ lines labeled. The blue region shows the temperature range spanned by core and shell He burning in our evolutionary 2.1 M ZAMS mass models. These reaction rate files are available at doi: 10.5281/zenodo.6779983.

MESA, wd builder, and GYRE
We use MESA version r15140 (Paxton et al. 2011(Paxton et al. , 2013(Paxton et al. , 2015(Paxton et al. , 2018(Paxton et al. , 2019 to evolve 2.1 M , Z = 0.0151 metallicity, non-rotating models from the ZAMS to the top of the WD cooling track in the Hertzsprung-Russell diagram. This results in 0.56 M DA WDs with 0.01 M He-rich mantles and 10 −3.5 M H-rich envelopes. One such evolutionary model is run for each 0.5σ step in the 12 C(α, γ) 16 O reaction rate probability distribution function. We refer to this set of models as the evolutionary DAVs. Each model used 20,000 cells with a 30 isotope nuclear reaction network and time resolution settings that consumed 3 months of wall-clock time on 16 cores to complete 550,000 time steps. Our models are similar to the lower resolution models used in Timmes et al. (2018) and Chidester et al. (2021).
We also use wd builder with MESA version r15140 to build ab initio WD models. By ab initio we mean calculations that begin with a hot WD model and an assumed chemical stratification, as opposed to a hot WD model that is the result of a stellar evolution calculation. The imposed initial 1 H, 4 He, 12 C, 14 N, 16 O, 22 Ne, and 56 Fe mass fraction profiles are taken from the evolutionary DAV models after the first thermal pulse on the asymptotic giant branch, defined by the first time the luminosity L > 10 4 L . The H envelope is then thinned to mimic the H envelope thickness of the evolutionary DAVs, the mass fraction profiles are smoothed at chemical transitions, and mass location where 16 O and 12 C exchange dominance is taken to the average mass location for all σ i . This initial conditions procedure is done for DAV and DBV wd builder models. The DBV wd builder models are then stripped of their H envelope.
We use ab initio WD models because they allow a more rapid exploration of the different WD classes, CO WD masses, and envelope masses needed to preliminary assess the robustness of our results. A potential disadvantage, or advantage, of ab initio WD models is that the imposed initial chemical stratification may not be attainable by a stellar model evolved from the ZAMS.
We use GYRE release 6.0 (Townsend & Teitler 2013;Townsend et al. 2018) embedded in MESA version r15140 to calculate the adiabatic pulsation properties as a WD model evolves. All the MESA + GYRE models begin from the top of the WD cooling track, terminate as cool WDs, and include the effects of element diffusion. Details of the MESA and wd builder models are in the files to reproduce our results at doi: 10.5281/zenodo.6779983.
We end this section by pointing out that we are not advocating for any specific model, or any specific setting used by an evolutionary model (e.g., convective mixing parameters). Rather our goal is to find, if there exists, potential signatures of the current experimental 12 C(α, γ) 16 O reaction rate probability distribution function in the pulsation periods of variable WD models. If such signatures exist and appear across the model space (including other researcher's variable WD models), then future uncertainty quantification studies could explore the impact of specific settings. Figure 1 shows the mass fraction profiles of the evolutionary DAV models for the 13 σ i reaction rates. We describe the labeled regions and chemical transitions. Region R1 extends from the center to the transition between 16 O and 12 C in the core, henceforth the O→C transition, encompassing the innermost 0.3 M . The reaction rate uncertainties have a large impact in R1, with the central 16 O mass fraction ranging from 0.45 for σ = -3.0 to 0.77 for σ = 3.0 in a regular pattern. The 16 O and 12 C mass fraction profiles are flat because this region forms during convective core He burning (see discussion of Figure 2). Deconvolving how much 16 O and 12 C is due to the 12 C(α, γ) 16 O reaction rate and how much is due to convective mixing processes, which has its own uncertainties, appears challenging (Straniero et al. 2003;De Gerónimo et al. 2017).

Evolutionary DAVs
Region R2 extends from the O→C transition to the transition between 12 C and 4 He, henceforth the C→He transition, encompassing 0.25 M . The reaction rate uncertainties also have a large impact in R2, with the 12 C and 16 O mass fraction profiles showing a regular pattern with σ i . This carbon-rich region forms during thick-shell He burning under radiative conditions (see Figure 2), mitigating the impact of convective mixing processes. The significant spread in the 16 O and 12 C mass fraction profiles is set by the 12 C(α, γ) 16 O reaction rate probability distribution function. The broad trend of the 12 C mass fraction increasing and the 16 O mass fraction decreasing with enclosed mass reflects the decreasing temperature and density in the radial direction during the formation of R2.
Region R3 extends from the C→He transition to the He→H transition. In the He-rich mantle of R3, the more  Figure 2 shows the evolution of the σ = 0 model near core He depletion and the onset of thick-shell He burning. Region R1 in the final WD model is formed during convective core He burning. The extent of the convective core does not intrude into region R2. Core convection ceases at core He depletion, marked by the red line. The O-rich R1 region then becomes radiative for the remainder of the evolution. Region R2 in the final WD model is formed after core He depletion during radiative thickshell He burning. The resulting C-rich R2 region then remains radiative for the remainder of the evolution.
It is useful to reference features with respect to mass or radius. Figure 3 thus shows the mass-radius relation for the σ = 0 model with regions and transitions labeled.
Some pulsation modes are more informative of these four regions than others. This can be due to a resonance, or near-resonance, between the mode's radial wavelength and thickness of one or more of the composition layers (Brassard et al. 1991). As these modes traverse a composition gradient within a local resonance region, they are partially reflected and become "trapped" within the local layer (Winget et al. 1981). Such modes are identified by showing a local minima in a kinetic energy diagram. We choose to refer to any mode dis-playing a local minima in the kinetic energy diagram as a trapped mode, regardless of the resonant region's location (e.g., upper or deeper layers). Modes trapped by the upper layers can reveal insights about WD envelopes (Kawaler & Weiss 1990;Brassard et al. 1992b;Kawaler 1995;Costa et al. 2008). Modes trapped deeper in the WD can be sensitive to the different regions, and thus reveal insights on the interior chemical profiles (Brassard et al. 1992a;Córsico et al. 2002;Giammichele et al. 2017). In particular, g-modes trapped by R2 can probe the current experimental 12 C(α, γ) reaction rate probability distribution function.
The Brunt-Väisälä frequency N is a characteristic frequency for pulsations, specifically the frequency of oscillation about an equilibrium position under gravity: where g is the gravitational acceleration, ρ is the mass density, P is the pressure, T is the temperature, χ T is the temperature exponent ∂(lnP )/∂(lnρ)| T,µ I , χ ρ is the density exponent ∂(lnP )/∂(lnT )| ρ,µ I , ∇ ad is the adiabatic temperature gradient, ∇ T is the actual temperature gradient, and B is the Ledoux term that accounts for composition gradients (e.g., Hansen & Kawaler 1994;Fontaine & Brassard 2008). The implementation of Equation 1 in MESA is described in Paxton et al. (2013). Figure 4 shows the Brunt-Väisälä frequency and the Ledoux B term profiles as a function of fractional radius after the 13 σ i evolutionary DAVs have cooled to T eff = 11,500 K. The composition gradients in the O→C and C→He transition regions, highlighted by the blue bands, induce bumps in the Ledoux B term profile and thus bumps in the N 2 profile. The first peak is located at the O→C transition and magnified in the inset plot. As σ i increases from -3.0 to 3.0, the location and magnitude of B increase with radius in a near-regular pattern.
The kinetic energy E kin of each g-mode can be expressed as (Unno et al. 1989;Córsico et al. 2002) where n is the radial order, the spherical harmonic degree, G the gravitational constant, M the stellar mass, R the stellar radius, ω 2 n = f 2 n (GM/R 3 ) −1 the dimensionless eigenfrequency, f n = 2π/P n the frequency, P n the period, C 1 = x 3 (M/M r ) the scaled density, r the radial distance from the center, M r the stellar mass enclosed at radius r, x = r/R the scaled radius, and y 1 and y 2 are the dimensionless Dziembowski (1971)   We numerically integrate Equation 2 for g-modes of radial orders n = 1-19 and harmonic degree = 1. Unless otherwise specified, all g-modes considered have = 1. Figure 5 shows the resulting E kin and period spacing ∆P n = P n+1 − P n diagrams for the evolutionary DAVs when they have cooled to T eff = 11,500 K. Trapped modes are identified by minima in E kin and by a corresponding minima in the period spacing (e.g., Winget et al. 1981;Brassard et al. 1991). Figure 5 highlights two strong local minima, the g 5 and g 10 trapped modes.
The frequency of an adiabatic g-mode is the area under a curve known as a weight function dζ/dr where, following Kawaler et al. (1985), dζ dr = [C(y, r) + N (y, r) + G(y, r)]ρr 2 r=R r=0 T (y, r)ρr 2 dr (4) and C(y, r) varies with the Lamb frequency, N (y, r) contains the Brunt-Väisälä frequency, G(y, r) involves the gravitational eigenfunctions, T (y, r) is proportional to the kinetic energy density, and y = (y 1 , y 2 , y 3 , y 4 ) are the Dziembowski (1971) variables. Figure 6 shows the weight functions of the g 4 , g 5 g 6 , and g 10 modes. The g 5 and g 10 trapped modes have larger weight functions in regions R1 and R2 compared to the g 4 and g 6 non-trapped modes. The frequency of the trapped modes is thus more strongly weighted by regions R1 and R2. The g 10 weight function is more equally distributed than the g 5 weight function, and the peak at the O→C transition is about the same height as the other peaks. Both factors decrease the ability of g 10 to probe the O→C transition and region R2. In contrast, the peak in the g 5 weight function at the O→C transition is the largest peak. About 35% of the g 5 frequency comes from the O→C transition, and 67% from the O→C transition and R2. Taken together, these weight functions indicate that R2 contributes 30-40% and R1 contributes 20% to the periods of g 5 and g 10 . This suggests that applying our results, conducted over the current experimental ± 3σ uncertainty in the 12 C(α, γ) 16 O reaction rate probability distribution function, to the derived pulsation spectrum of observed variable WD may face challenges obtaining precision constraints on the reaction rate.   Resonance regions can also be identified by checking when the width of a chemical stratification region is equal to an integer number of radial wavelengths, λ r = 2π/k r where the wavenumber is and S is the Lamb frequency. Figure 7 shows integer multiples q of λ r for the g 5 and g 10 modes versus radius for σ 0 . When q · λ r is near a black segment, the g-mode resonates, or nearly resonates, with the region's width. The 2·λ r curve for g 5 lies close to the width of R2, identifying R2 as the resonant cavity for the g 5 mode. For the g 10 mode the 4·λ r curve is close to the widths of R1 and R2. We also verified 6·λ r is close to the width of R3. Larger q values for g 10 may resonate with R4. As g 10 resonates with multiple regions, the mode is not uniquely trapped. This is commensurate with the more uniform distribution of peaks in the g 10 weight function of Figure 6. Figure 8 shows the history of the relative period differences for g 5 , g 6 , and g 10 as the evolutionary DAVs cool.
The g 5 trapped mode shows the most distinctive trend in the period with σ i out of every g-mode in the range 1 ≤ n ≤ 20. The relative period differences span ±2%, with positive σ i yielding shorter periods and negative σ i yielding longer periods. The g 10 trapped mode is not as distinctive, showing a smaller spread in periods with σ i . The g 6 non-trapped mode shows little spread with σ i and no distinctive trend with period. While Figure 8 highlights g 6 and g 10 , we emphasize that no other gmode within 1 ≤ n ≤ 20 shows any distinctive period pattern with σ i . Figure 8 also shows the relative period differences for g 5 with respect to σ i at T eff = 11,500 K. The relative period differences increase from −2.0% for σ=−3.0 to 3.0% for σ = 3.0. The relationship is non-monotonic due to variations in the location of the C→He transition, which impacts the width of R2 where the g 5 trapped mode resonates. Minimizing variations in this location may increase the monotonicity of this relationship.

Variations
We have presented evidence that adiabatic g-modes trapped by the radiatively formed carbon-rich layer in CO WD models offer potentially useful probes of the current experimental 12 C(α, γ) 16 O reaction rate probability distribution function. In this section we give a preliminary assessment of the robustness of this result by sampling different WD classes, masses, and shell masses. Each sequence contains 13 σ i models. For each sequence, we verify the existence of a trapped or partially trapped mode that probes region R2 as indicated by the sequences's kinetic energy diagram and weight functions. We also confirm that the sequences's R2 trapped g-mode gives the most distinct ∆P/P 0 versus σ i relationship. Figure 9 shows the results of this survey. The x-axis is the sequence and the y-axis is relative period difference ∆P/P 0 for the labeled R2 trapped g-mode. Proceeding from left to right, the first sequence is the evolutionary DAVs which are analyzed in detail above.
The next three sequences explore the impact of the WD mass using wd builder DAV models of 0.52, 0.56, and 0.73 M , respectively. All three sequences show g 6 is the most distinctive adiabatic trapped g-mode that probes R2. The differences between the 0.56 M evolutionary DAVs and the 0.56 M wd builder DAVs is due to their different composition profiles (see Section 2). This suggests that the most distinctive trapped g-mode depends on model details. =-3.0 =-2.5 =-2.0 =-1.5 =-1.0 =-0.5 =0.5 =1.0 =1.5 =2.0 =2.5 =3.0 Figure 9. Relative period differences for model sequences of varying WD masses, shell masses, and classes. The wd builder DAVs were measured at T eff = 11,500 K for the wd builder DBVs were measured at T eff = 25,000 K. Labeled are the ± 3 σ uncertainties on the experimental 12 C(α, γ) 16 O reaction rate in steps of 0.5σ, the trapped g-mode that most distinctly probes the radiatively-formed, carbon-rich region R2, and the mean relative period difference for each sequence.
∆P/P 0 +1 −2 % range than the evolutionary DAVs, and show non-monotonic spacings with σ i . The 0.73 M sequence spans ∆P/P 0 +2 −4 % and shows a monotonic spacing with σ i . This suggests that more massive WDs may give stronger signals with this method.
The fifth sequence shows the impact of a H envelope that is 10 times thinner than the 0.56 M wd builder DAVs. One dipole = 1 and two quadruple = 2 trapped g-modes are shown. The g 2 trapped mode with = 1 shows a regular pattern with σ i over the range ∆P/P 0 +1 −2.5 %. The g 2 trapped mode with = 2 is similar to its = 1 counterpart, and the g 7 trapped mode with = 2 shows a larger ∆P/P 0 +1.5 −2.5 % range. This sequence suggests that harmonic degrees = 1 and/or = 2 can have R2 trapped modes that distinguish σ i . The next sequence shows the wd builder DBVs. A regular pattern with σ i for the g 3 trapped mode emerges, and the sequence spans ∆P/P 0 +0.5 −1 %. The σ = 2.0 point in this sequence is an outlier that we cannot explain. The last sequence shows the impact of a He envelope that is 10 times thinner than the 0.56 M wd builder DBVs. A regular pattern with σ i persists at the ∆P/P 0 +0.5 −0.5 % level for the g 2 trapped mode. For all the sequences, the positive σ i have shorter periods than the σ = 0 model. Following Chidester et al. (2021), the main contributors to period differences are changes in the local pressure scale height H, mean molecular weight µ i , temperature, density exponent χ ρ , and first adiabatic index Γ 1 . Changes in these parameters compete with one another to drive the overall period difference. Models with larger mass fractions of 22 Ne give shorter periods due to a smaller H. As increasing the 12 C(α, γ) 16 O increases the 16 O content in R2, smaller H drive shorter periods. Similar logic applies to negative σ i models showing longer periods.

SENSITIVITIES
There are many potential sensitivities that we have not investigated. We highlight four and the edges of R2 that they affect.

Width of the O→C transition
Mode trapping by R2 may depend on the width of the O→C transition. Increasing the O→C transition width decreases the µ i gradient, which in turn may decrease the trapping ability of R2. The width of the O→C transition in our models is relatively narrow, in agreement with the MESA models of Pepper et al. (2022), while the O→C transition in Córsico et al. (2002) is considerably wider. Córsico et al. (2002) also find that their adia-batic g 5 mode depends on the O→C transition and R2. However, their g 5 mode is a local maximum in their E kin diagram, meaning it is not a trapped mode. With mode trapping suppressed for a wide enough O→C transition, the kinetic energy from trapped modes in R2 is released. As Equation 2 is weighted with the density, a non-trapped g-mode sensitive to R2 may appear as a local maximum in E kin , as found by Córsico et al. (2002). Relative extrema in E kin , not just minima, might find gmodes that probe R2 and thus also reveal inferences on σ i . Region R2 still exists, because it forms under thickshell radiative He-burning conditions. We caution that its trapping abilities, and thus our relative period shifts may change with the width of the O→C transition.

3α reaction rate probability distribution function
The 3α process, the fusion of three 4 He nuclei into one 12 C nucleus, impacts the innermost and outermost edges of R2. The 3α process is followed by the subsequent α-capture reaction 12 C(α, γ) 16 O. The final mass fractions of 12 C and 16 O, under radiative burning conditions, is determined by the competition between the 3α and 12 C(α, γ) 16 O reaction rates at a given temperature. The feeding of 12 C, driven by the 3α process, occurs early in the evolution when the mass fraction of 12 C is low and 4 He is high (the carbon bump at the outermost edge of R2). Oxygen production occurs later by α-capture on the freshly produced 12 C (the O→C transition at the innermost edge of R2). Current estimates of the uncertainty in the 3α reaction rate are 30% over the regions of typical astrophysical interest (Kibédi et al. 2020). However, at lower temperatures ( 0.1 GK), the uncertainty is likely much larger, because other reaction mechanisms become significant (Suno et al. 2016). Region R2 again still exists, and we again caution that the period shifts we find from the 12 C(α, γ) 16 O probability distribution may change when a 3α reaction rate probability distribution function is considered.

Mixing during core He-depletion
Mixing in low-and intermediate-mass stars during core He-burning is particularly challenging to model (Salaris & Cassisi 2017;Jermyn et al. 2022) and impacts the innermost edge of R2. The radiative gradient profile within a He core convection region develops a local minimum at some point during its evolution (e.g., see Figure 6 of Paxton et al. 2018). With further outward propagation of the convective boundary, or the action of overshooting, the mixing of fresh He into the core can lower the radiative gradient throughout the core to such an extent that it equals the adiabatic gradient at the local minimum of the radiative gradient. When this happens, the convective region interior to the minimum becomes decoupled from the region exterior to the minimum: the convective core splits (Eggleton 1972). In addition, even small amounts of He added to the convective core enhances the rate of energy production by the 12 C(α, γ) 16 O reaction. The resulting increase in the radiative gradient can lead to rapid growth in the convective He core boundary (a "breathing pulse"). The enhanced nuclear burning also increases the central 16 O mass fraction. A concensus on breathing pulses being physical or numerical has not yet been reached (Caputo et al. 1989;Cassisi et al. 2003;Farmer et al. 2016;Constantino et al. 2017;Paxton et al. 2019). Region R2 persists, and we caution that the absolute period shifts we find may change with different treatments of core He-depletion in evolutionary models.
De Gerónimo et al. (2017) found that overshooting during the core He burning leaves imprints on the Brunt-Väisälä frequency that result in absolute period differences of 2-5 s on average, relative to models with zero overshooting. Models that included overshooting had a larger central 16 O mass fraction and an extended R1. They found that these results are, on average, independent of the ZAMS mass. De Gerónimo et al.  Kunz et al. (2002) reaction rate is not the same as adopting a modern reaction rate sourced from a probability distribution function (Mehta et al. 2022). We suspect the De Gerónimo et al. (2017) result is partially due to using a 12 C(α, γ) 16 O reaction rate that is 10% larger than their base reaction rate, as the highest rate considered. In contrast our 12 C(α, γ) 16 O rate probability distribution function spans 0.5-1.5 times our nominal σ = 0 reaction rate (see Figure 1). Therefore, we find larger overall average period differences from the 12 C(α, γ) 16 O reaction rate. Future uncertainty quantification studies could explore a potential coupling between simultaneous variations in overshooting and the adopted 12 C(α, γ) 16 O reaction rate. It is possible that a strong coupling could alter the R2 trapped mode properties of our models.

Number of Thermal Pulses
The thermal pulse phase of evolution impacts the outermost edge of R2. For the case of a fixed number of thermal pulses, De Gerónimo et al. (2017) found period differences of 5-10 s for their 0.548 M model and 2-3 s for their 0.837 M model. This effect is mainly due to the C→He transition being less pronounced in their more massive WD models. The impact of the thermal pulses in our evolutionary models is shown in Figure 1. Each of our 13 σ i models experienced 14 thermal pulses, with the onset of each thermal pulse defined by the photon luminosity exceeding 10 4 L . In contrast, our wd builder models were inferred from the chemical profiles at the first thermal pulse, fixing the number of thermal pulses for those sequences. The sensitivity of our results to the number of thermal pulses can thus be estimated by comparing the g 5 trapped mode periods of our 0.56 M evolutionary and wd builder models.
We find a standard deviation of 2.3 for the thermal pulses and a standard deviation of 4.5 for the 12 C(α, γ) 16 O. This suggests that variations from the number of thermal pulses in our models is smaller than the variations from the 12 C(α, γ) 16 O reaction rate probability distribution function. In contrast, De Gerónimo et al. (2017) found larger variations from thermal pulses than the 12 C(α, γ) 16 O reaction rate. We speculate this difference is again due to our larger span of 12 C(α, γ) 16 O reaction rates. Our 14:1 thermal pulse ratio is also larger than the 10:3 thermal pulse ratio of De Gerónimo et al. (2017), which may strengthen our result that variations from the 12 C(α, γ) 16 O reaction rate probability distribution function have a larger impact than the number of thermal pulses. In addition, Pepper et al. (2022) found the number of thermal pulses is dependent on the 12 C(α, γ) 16 O reaction rate, which Figure 1 confirms. Also in agreement with Pepper et al. (2022), we find that smaller σ i increases the number of thermal pulses as smaller reaction rates have larger He-shell masses at the onset of each pulse.

SUMMARY
We conducted a new search for signatures of the current experimental 12 C(α, γ) 16 O reaction rate probability distribution function in the pulsation periods of CO WD models. We found that adiabatic g-modes trapped or partially trapped by the interior C-rich layer (region R2 in Figure 1) offer the most direct probe of this reaction rate because this region forms under radiative He burning conditions, mitigating the impact of uncertainties from convective mixing processes. We found an average spread in relative period shifts of ∆P/P ± 2% for the identified trapped g-modes over the experimental ± 3σ probability distribution function of the 12 C(α, γ) 16 O reaction rate. We found the effect persists across the observed T eff window of DAV and DBV WDs, and for different WD masses and smaller H/He shell masses. Figures 1, 8, and 9 make the first direct quantitative con-nection between the pulsation periods of variable WD models and the current, experimental 12 C(α, γ) 16 O reaction rate probability distribution function.
The C-rich layer is a "sweet spot" in DAV and DBV WD models for probing the 12 C(α, γ) 16 O reaction rate probability distribution function. Figure 9 suggests a corresponding "sweet spot" of g-modes with radial orders 2 n 7 can investigate R2. This suggestion is complemented by an analysis from Córsico et al. (2002), who found that all g-mode periods 500-600 s were trapped (or nearly trapped) in the H-rich envelope. They found the weight functions of those modes were low in amplitude, and similar to one another, and concluded that mode trapping vanishes for long periods (higher radial orders). Thus, inferences from trapped modes such as R2 are limited to g-modes with periods 500-600 s. This confines our models to have g-modes n 10.
In every model sequence explored, the g-mode that best distinguishes the 12 C(α, γ) 16 O reaction rate follows five specific patterns: 1. The g-mode is trapped, as confirmed by its local minimum in E kin .
2. The g-mode best resonates with the R2 region; it has a q·λ r matching best with the R2 region width, and gives a weight function substantially weighted in the R2 region.
3. The g-mode is within a radial order "sweet spot" of 2 n 7.
4. The g-mode period is shorter for positive σ i , and longer for negative σ i .
5. The g-mode shifts are within the detectable range.
The signatures persist because R2 forms under radiative helium burning conditions, but could be sensitive to the couplings with other uncertainties (see Section 4). Out of the entire g-mode spectrum, only one, the R2 trapped g-mode, consistently showed an identifiable trace to the 12 C(α, γ) 16 O reaction rate probability distribution function. Moreover, the R2 trapped mode consistently followed the 5 patterns listed above, irrespective of model type, WD class, WD mass, and envelope thickness. Thus, the trapped R2 g-mode signature found is the most direct information at tracing the 12 C(α, γ) 16 O reaction rate probablity distribution function using WD seismology.
The g-mode periods of observed variable WD are derived from a Fourier analysis of the photometric light curves and are typically given to 6−7 significant figures of precision (e.g., Duan et al. 2021). Usually WD composition profile templates are fit to the extracted g-mode period spectrum and other observed constraints (e.g., T eff , log g) of a specific WD. The root-mean-square residuals to the 150−400 s low-order g-mode periods are typically in the range σ rms 0.3 s (e.g., Bischoff-Kim et al. 2014), for a fit precision of σ rms /P 0.3%. Lower root-mean-square residuals using ab initio WD models are possible (Charpinet et al. 2019;Giammichele et al. 2022), although see De Gerónimo et al. (2019).
Our finding of relative period shifts of ∆P/P ± 2.0% suggests that an astrophysical constraint on the 12 C(α, γ) 16 O reaction rate probability distribution function could, in principle, be extractable from the derived period spectrum of observed variable WD. Our results can inform future inferences, including those from machine-learning (e.g., Bellinger et al. 2016), on the interior mass fraction profiles and the reaction rates that produce those chemical profiles.