22-NuclearEnergy (template)

import numpy as np, copy
import matplotlib.pyplot as plt
%matplotlib inline
from astropy import constants as const
import astropy.units as u
from astropy.units import cds
cds.enable() 

<astropy.units.core._UnitContext at 0x105f87370>

1. At home: pp vs CNO rates

TODO

1. Make a plot of \(\log(\epsilon)\) as a function of \(\log(T)\) for the PP-I chain and the CNO cycle, based on the expression given in class (from Kip 18.5.1). Use the temperature range already given in the code below.

2. In both cases (PP-I and CNO), by how much does the \(\epsilon\) changes for a change of temperature of one order of magnitude?

3. If you were to reproduce \(\epsilon\) in the range of temperature below by a power law of the form \(\epsilon \propto T^\nu\), what would the index \(\nu\) of the power law be? (Graph a curve to find your answer.)

fig, ax = plt.subplots(1,1)
ax.set_ylim(-5,6)

logT = np.linspace(6,8,50)


#########################
# At home

TODO: Write a caption for your graph, and answer the questions asked.

2. Let’s look at the energy generation in our MESA models

Reading in the models

First, use the cell below to read in the two MESA models from the convection notebook for the Sun (in variable Msun) and for a 10 \(M_\odot\) star (in variable M10).

def read_model(file):

    return np.genfromtxt(file, skip_header=5, names=True)

Msun = read_model('https://raw.githubusercontent.com/veropetit/PHYS633-F2024/main/Book/L19-Convection/19-Sun-profile8.data')
M10  = read_model('https://raw.githubusercontent.com/veropetit/PHYS633-F2024/main/Book/L19-Convection/19-M10-profile8.data')

1. In class: we will look at the energy generation in the Solar model.

In MESA, the column named eps_nuc gives the total power per gram coming from all of the nuclear reactions. The columns pp and cno give the power/g coming from the PP (I, II, and III) chains and the CNO cycle, respectively.

We will make a gapph of \(\epsilon_\mathrm{nuc}\), \(\epsilon_\mathrm{pp}\) and \(\epsilon_\mathrm{cno}\) as a function of radius, and then compare the values in MESA with the values we get from the formula given in the textbook.

At home:

TODO

1. On the left graph, add curves for the \(\epsilon_\mathrm{nuc}\), \(\epsilon_\mathrm{pp}\) and \(\epsilon_\mathrm{cno}\) given by the MESA model.

2. Make a similar plot for the 10 \(M_\odot\) model on the right graph.

3. In the interpretation box below, comment on the similarities/differences between the two models.

4. Use the analytical expression for \(\epsilon_\mathrm{pp-1}\) and \(\epsilon_\mathrm{cno}\) from Part 1 to calculate your expected values for \(\epsilon_\mathrm{pp-1}\) and \(\epsilon_\mathrm{cno}\) in both models, and add these curves to your plots. In the MESA model, you will some useful columns such as:

   • logRho for the logarithm10 of the density in g/cm2

   • logT for the logarithm10 of the temperature in K

   • h1 for the mass fraction of hydrogen

   • c12, n14, o16,o18 for the mass fraction of important isotopes of carbon, nitrogen, and oxygen.

5. In the interpretation box below, compare your calculation in #4 to the real \(\epsilon\) in the models. If there are any differences, comment on what could explain these discrepencies.

fig, ax = plt.subplots(1,2, figsize=(12,5.5))

ax[0].set_title('Sun')
ax[1].set_title(r'10 M$_\odot$')

ax[0].set_xlabel('$r$ ($R_\odot$)')
ax[1].set_xlabel('$r$ ($R_\odot$)')

ax[0].set_ylabel('erg/g')
ax[1].set_ylabel('erg/g')

plt.tight_layout()


#########################
# In class




#########################
# At home

TODO Write a caption for your graph, and answer the questions asked.