{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Radioactive Decay Energy\n", "\n", "Within the ejecta of a supernova, the $\\gamma$-rays largely come from the decay of $^{56}Ni$ into $^{56}Co$, which releases a significant amount of energy. \n", "\n", "When $^{56}Ni$ decays into $^{56}Co$ it can release a $\\gamma$-ray at several different transition levels. Each transition level has an energy and an associated probability out of 100 decays. For example, the transition from Energy level 9 to Energy level 7 has an energy of 0.270 Mev and a probability of 36.5 out of 100 decays. To find the total energy per decay you multipliy each energy with its associated probability and add them all up." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "/home/ceceliapowers/software/tardis/tardis/__init__.py:23: UserWarning: Astropy is already imported externally. Astropy should be imported after TARDIS.\n", " warnings.warn(\n" ] }, { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "0fdee0253d604567880aec1ca9bdd8d4", "version_major": 2, "version_minor": 0 }, "text/plain": [ "Iterations: 0/? [00:00" ] }, "execution_count": 2, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# energies of each transition\n", "t_energies = np.array([0.270, 0.750, 0.480, 1.56, 0.812, 0.158]) * u.MeV\n", "# probabilities of each transition\n", "t_prob = np.array([.365, .495, .366, .140, .860, 1.00])\n", "\n", "energy_per_decay = sum(t_energies * t_prob)\n", "energy_per_decay" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "From the above cell, we get the energy per transition of 1.72 MeV. Note that this comes from a simplified scheme of energies and a more accurate total energy per $^{56}Ni$ decay is 1.75 MeV. [] " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The $^{56}Co$ produced from the decay of $^{56}Ni$ is also radioactive and will decay into $^{56}Fe$ and release more $\\gamma$-rays, however this decay is more complicated than the decay of $^{56}Ni$. Whereas $^{56}Ni$ only decays through electron capture, $^{56}Co$ can decay either by electron capture, which occurs for 81 out of 100 cases, or through positron decay, which occurs for 19 out of 100 cases.\n", "\n", "Positron decay produces positrons with a given kinetic energy, that will eventually annihilate with electrons to produce two 0.511 MeV $\\gamma$-rays. The scheme of decays for $^{56}Co$ is slightly more complicated than the $^{56}Ni$ scheme, but to find the total energy per decay, you follow the same process. The total energy per decay from $\\gamma$-rays is 3.61 MeV and the total kinetic energy of positrons is 0.12 MeV" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ " The total rate of energy production for a mass of $^{56}Ni$ at a given time is given by the following equation:\n", "\n", "$$\\epsilon = \\frac{M_\\odot}{56m_{u}}\\frac{1}{\\tau_{\\text{Co}}-\\tau_{\\text{Ni}}}[[Q_{\\text{Ni}}(\\frac{\\tau_{\\text{Co}}}{\\tau_{\\text{Ni}}}-1)-Q_{\\text{Co}}]\\exp(-t/\\tau_{\\text{Ni}})+Q_{\\text{Co}}\\exp(-t/\\tau_{\\text{Co}})]\\frac{M_{Ni0}}{M_\\odot}$$\n", "\n", "$M_\\odot$ is a solar mass. $56_{u}$ is 56 atomic mass units. \n", "\n", "$\\tau_{Ni}$ is the lifetime of $^{56}Ni$ which is 8.80 days and $\\tau_{Co}$ is the lifetime of $^{56}Co$ which is 111.3 days. \n", "\n", "$Q_{\\text{Ni}}$ is the energy per decay of $^{56}Ni$ which is 1.75 MeV and $Q_{\\text{Co}}$ is the sum of the energy per decay of $^{56}Co$ from $\\gamma$-rays and the kinetic energy from positrons which is 3.73 MeV\n", "\n", "If we plug these values into the equation we get the equation:\n", "\n", "$$\\epsilon = (6.45e43\\exp(-t/8.8)+1.45e43\\exp(-t/111.3))\\frac{M_{Ni0}}{M_\\odot} erg/s$$" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The total energy production rate for 1 solar mass of 56Ni after 10 days is: 3.40e+43 erg / s\n" ] } ], "source": [ "#time in days\n", "time = 10 * u.day\n", "#mass of Ni56 in solar masses\n", "m_Ni56 = 1 * const.M_sun\n", "\n", "energy_production_rate = (6.45e43 * np.exp(-time/(8.8*u.day)) + 1.45e43 * np.exp(-time/(111.3*u.day))) * (m_Ni56/ const.M_sun) * u.erg /u.s\n", "\n", "print(f\"The total energy production rate for 1 solar mass of 56Ni after 10 days is: {energy_production_rate:.2e}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The total total energy production integrated over time is given by:\n", "$$\n", "\\mathcal{E} = \\mathcal{E}_{Ni} + \\mathcal{E}_{Co} = 1.885e50 \\frac{M_{Ni0}}{M_\\odot} ergs\n", "$$\n", "where\n", "$$\n", "\\mathcal{E}_{Ni} = 6.22e49 \\frac{M_{Ni0}}{M_\\odot} ergs\n", "$$\n", "$$\n", "\\mathcal{E}_{Co} = 1.26e50 \\frac{M_{Ni0}}{M_\\odot} erg\n", "$$" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The total energy production for 1 solar mass of 56Ni is: 1.885e+50 erg\n" ] } ], "source": [ "total_energy_production = 1.885e50 * (m_Ni56/const.M_sun) * u.erg\n", "print(f\"The total energy production for 1 solar mass of 56Ni is: {total_energy_production}\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next we show how we can find the same result from functions in TARDIS-HE\n", "\n", "First we create a DataFrame with isotope abundances. This DataFrame typically comes from a model, but for our purposes we create an ejecta with 5 shells, made entirely of $^{56}Ni$." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
01234
atomic_numbermass_number
285611111
\n", "
" ], "text/plain": [ " 0 1 2 3 4\n", "atomic_number mass_number \n", "28 56 1 1 1 1 1" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "mass_fractions = {0:1,\n", " 1:1,\n", " 2:1,\n", " 3:1,\n", " 4:1}\n", "\n", "raw_isotope_abundance = pd.DataFrame(mass_fractions, index=pd.MultiIndex.from_tuples([(28, 56)], names=[\"atomic_number\", \"mass_number\"]))\n", "raw_isotope_abundance" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Then we create an array of shell masses. We want to find the total energy production of 1 solar mass of $^{56}Ni$, so we give each of the five shells a mass of $\\frac{1}{5}M_\\odot$" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$[3.9768197 \\times 10^{32},~3.9768197 \\times 10^{32},~3.9768197 \\times 10^{32},~3.9768197 \\times 10^{32},~3.9768197 \\times 10^{32}] \\; \\mathrm{g}$" ], "text/plain": [ "" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "shell_masses = np.ones(5)*0.2*const.M_sun.to(u.g)\n", "shell_masses" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we can use TARDIS-HE functions to get a DataFrame with information on the decays of this model. The first step is to create a dictionary which has the mass of each isotope in each shell." ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{0: {'Ni56': 3.976819741396102e+32},\n", " 1: {'Ni56': 3.976819741396102e+32},\n", " 2: {'Ni56': 3.976819741396102e+32},\n", " 3: {'Ni56': 3.976819741396102e+32},\n", " 4: {'Ni56': 3.976819741396102e+32}}" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "isotope_dict = create_isotope_dicts(raw_isotope_abundance, shell_masses)\n", "isotope_dict" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next, we create a dictionary that holds the inventory of each isotope given in the [`radioactivedecay`](https://radioactivedecay.github.io/) package." ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{0: Inventory activities (Bq): {'Ni-56': 5.653446204017355e+48}, decay dataset: icrp107_ame2020_nubase2020,\n", " 1: Inventory activities (Bq): {'Ni-56': 5.653446204017355e+48}, decay dataset: icrp107_ame2020_nubase2020,\n", " 2: Inventory activities (Bq): {'Ni-56': 5.653446204017355e+48}, decay dataset: icrp107_ame2020_nubase2020,\n", " 3: Inventory activities (Bq): {'Ni-56': 5.653446204017355e+48}, decay dataset: icrp107_ame2020_nubase2020,\n", " 4: Inventory activities (Bq): {'Ni-56': 5.653446204017355e+48}, decay dataset: icrp107_ame2020_nubase2020}" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "inventories = create_inventories_dict(isotope_dict)\n", "inventories" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we calculate the number of decays for each isotope in each shell and create a DataFrame to hold this information." ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " number_of_decays\n", "shell_number isotope \n", "0 Ni56 4.281026e+54\n", " Co56 4.281026e+54\n", "1 Ni56 4.281026e+54\n", " Co56 4.281026e+54\n", "2 Ni56 4.281026e+54\n", " Co56 4.281026e+54\n", "3 Ni56 4.281026e+54\n", " Co56 4.281026e+54\n", "4 Ni56 4.281026e+54\n", " Co56 4.281026e+54\n" ] } ], "source": [ "cumulative_decay_df = calculate_total_decays(inventories, time_delta= np.inf)\n", "print(cumulative_decay_df)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next we use the [kurucz_cd23_chianti_H_He.h5](https://github.com/tardis-sn/tardis-regression-data/raw/main/atom_data/kurucz_cd23_chianti_H_He.h5) dataset to get the decay radiation data for each isotope." ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "WARNING:tardis.io.atom_data.atom_web_download:Atomic Data kurucz_cd23_chianti_H_He already exists in /home/ceceliapowers/Downloads/tardis-data/kurucz_cd23_chianti_H_He.h5. Will not download - override with force_download=True.\n", "INFO:tardis.io.atom_data.util:\n", "\tAtom Data kurucz_cd23_chianti_H_He.h5 not found in local path.\n", "\tExists in TARDIS Data repo /home/ceceliapowers/Downloads/tardis-data/kurucz_cd23_chianti_H_He.h5\n", "INFO:tardis.io.atom_data.base:Reading Atom Data with: UUID = 6f7b09e887a311e7a06b246e96350010 MD5 = 864f1753714343c41f99cb065710cace \n", "INFO:tardis.io.atom_data.base:Non provided Atomic Data: synpp_refs, photoionization_data, yg_data, two_photon_data, linelist_atoms, linelist_molecules\n" ] }, { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
AElementZNParent E(level)UncertaintyJPiMetastableDecay ModeDecay Mode Value...RadiationRad subtypeRad EnergyUncertainty.1EP EnergyUncertainty.2Rad IntensityUncertainty.3DoseUncertainty.4
Isotope
Nn11Nn010.0NaN1/2+FalseB-100.0...bmNaN301.3700NaN782.34710.0100.0NaN0.301000NaN
Nn11Nn010.0NaN1/2+FalseB-100.0...bm avNaN301.3700NaNNaNNaN100.0NaN0.301000NaN
H33H120.0NaN1/2+FalseB-100.0...bmNaN5.681712.018.5913.0100.0NaN0.005682NaN
H33H120.0NaN1/2+FalseB-100.0...bm avNaN5.6820NaNNaNNaN100.0NaN0.005680NaN
He66He240NaN0+FalseB-100.0...bmNaN1567.620054.03507.80011.0100.0NaN1.567600NaN
..................................................................
Lv292292Lv1161760NaN0+FalseA100.0...NaNNaNNaNNaNNaNNaNNaNNaNNaNNaN
Lv293293Lv1161770NaNNaNFalseA100.0...NaNNaNNaNNaNNaNNaNNaNNaNNaNNaN
Ts293293Ts1171760NaNNaNFalseA100.0...NaNNaNNaNNaNNaNNaNNaNNaNNaNNaN
Og294294Og1181760NaN0+FalseA100.0...NaNNaNNaNNaNNaNNaNNaNNaNNaNNaN
Ts294294Ts1171770NaNNaNFalseA100.0...NaNNaNNaNNaNNaNNaNNaNNaNNaNNaN
\n", "

242295 rows × 41 columns

\n", "
" ], "text/plain": [ " A Element Z N Parent E(level) Uncertainty JPi Metastable \\\n", "Isotope \n", "Nn1 1 Nn 0 1 0.0 NaN 1/2+ False \n", "Nn1 1 Nn 0 1 0.0 NaN 1/2+ False \n", "H3 3 H 1 2 0.0 NaN 1/2+ False \n", "H3 3 H 1 2 0.0 NaN 1/2+ False \n", "He6 6 He 2 4 0 NaN 0+ False \n", "... ... ... ... ... ... ... ... ... \n", "Lv292 292 Lv 116 176 0 NaN 0+ False \n", "Lv293 293 Lv 116 177 0 NaN NaN False \n", "Ts293 293 Ts 117 176 0 NaN NaN False \n", "Og294 294 Og 118 176 0 NaN 0+ False \n", "Ts294 294 Ts 117 177 0 NaN NaN False \n", "\n", " Decay Mode Decay Mode Value ... Radiation Rad subtype Rad Energy \\\n", "Isotope ... \n", "Nn1 B- 100.0 ... bm NaN 301.3700 \n", "Nn1 B- 100.0 ... bm av NaN 301.3700 \n", "H3 B- 100.0 ... bm NaN 5.6817 \n", "H3 B- 100.0 ... bm av NaN 5.6820 \n", "He6 B- 100.0 ... bm NaN 1567.6200 \n", "... ... ... ... ... ... ... \n", "Lv292 A 100.0 ... NaN NaN NaN \n", "Lv293 A 100.0 ... NaN NaN NaN \n", "Ts293 A 100.0 ... NaN NaN NaN \n", "Og294 A 100.0 ... NaN NaN NaN \n", "Ts294 A 100.0 ... NaN NaN NaN \n", "\n", " Uncertainty.1 EP Energy Uncertainty.2 Rad Intensity Uncertainty.3 \\\n", "Isotope \n", "Nn1 NaN 782.347 10.0 100.0 NaN \n", "Nn1 NaN NaN NaN 100.0 NaN \n", "H3 12.0 18.591 3.0 100.0 NaN \n", "H3 NaN NaN NaN 100.0 NaN \n", "He6 54.0 3507.800 11.0 100.0 NaN \n", "... ... ... ... ... ... \n", "Lv292 NaN NaN NaN NaN NaN \n", "Lv293 NaN NaN NaN NaN NaN \n", "Ts293 NaN NaN NaN NaN NaN \n", "Og294 NaN NaN NaN NaN NaN \n", "Ts294 NaN NaN NaN NaN NaN \n", "\n", " Dose Uncertainty.4 \n", "Isotope \n", "Nn1 0.301000 NaN \n", "Nn1 0.301000 NaN \n", "H3 0.005682 NaN \n", "H3 0.005680 NaN \n", "He6 1.567600 NaN \n", "... ... ... \n", "Lv292 NaN NaN \n", "Lv293 NaN NaN \n", "Ts293 NaN NaN \n", "Og294 NaN NaN \n", "Ts294 NaN NaN \n", "\n", "[242295 rows x 41 columns]" ] }, "execution_count": 10, "metadata": {}, "output_type": "execute_result" } ], "source": [ "download_atom_data('kurucz_cd23_chianti_H_He')\n", "atom_data_file = 'kurucz_cd23_chianti_H_He.h5'\n", "atom_data = AtomData.from_hdf(atom_data_file)\n", "gamma_ray_lines = atom_data.decay_radiation_data\n", "gamma_ray_lines" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Finally, we can create a DataFrame with that holds all the information about the decays of each isotope in each shell." ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "data": { "text/html": [ "
\n", "\n", "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
number_of_decaysdecay_moderadiationradiation_energy_keVradiation_intensityenergy_per_channel_keVdecay_energy_keVdecay_energy_erg
shell_numberisotope
0Ni564.281026e+54ECbp65.0000.000600.0003901.669600e+512.674994e+42
Ni564.281026e+54ECbp408.0000.000030.0001225.239975e+508.395366e+41
Ni564.281026e+54ECbp478.0000.000030.0001436.138991e+509.835748e+41
Ni564.281026e+54ECbp av100.0000.000700.0007002.996718e+514.801271e+42
Ni564.281026e+54ECe6.07054.500003.3081501.416227e+552.269047e+46
..............................
4Co564.281026e+54ECg3369.8600.010100.3403561.457072e+542.334487e+45
Co564.281026e+54ECg3451.2320.9490032.7521921.402130e+562.246459e+47
Co564.281026e+54ECg3548.0500.195506.9364382.969507e+554.757674e+46
Co564.281026e+54ECg3600.8000.016700.6013342.574325e+544.124523e+45
Co564.281026e+54ECg3611.5300.008600.3105921.329651e+542.130335e+45
\n", "

485 rows × 8 columns

\n", "
" ], "text/plain": [ " number_of_decays decay_mode radiation \\\n", "shell_number isotope \n", "0 Ni56 4.281026e+54 EC bp \n", " Ni56 4.281026e+54 EC bp \n", " Ni56 4.281026e+54 EC bp \n", " Ni56 4.281026e+54 EC bp av \n", " Ni56 4.281026e+54 EC e \n", "... ... ... ... \n", "4 Co56 4.281026e+54 EC g \n", " Co56 4.281026e+54 EC g \n", " Co56 4.281026e+54 EC g \n", " Co56 4.281026e+54 EC g \n", " Co56 4.281026e+54 EC g \n", "\n", " radiation_energy_keV radiation_intensity \\\n", "shell_number isotope \n", "0 Ni56 65.000 0.00060 \n", " Ni56 408.000 0.00003 \n", " Ni56 478.000 0.00003 \n", " Ni56 100.000 0.00070 \n", " Ni56 6.070 54.50000 \n", "... ... ... \n", "4 Co56 3369.860 0.01010 \n", " Co56 3451.232 0.94900 \n", " Co56 3548.050 0.19550 \n", " Co56 3600.800 0.01670 \n", " Co56 3611.530 0.00860 \n", "\n", " energy_per_channel_keV decay_energy_keV \\\n", "shell_number isotope \n", "0 Ni56 0.000390 1.669600e+51 \n", " Ni56 0.000122 5.239975e+50 \n", " Ni56 0.000143 6.138991e+50 \n", " Ni56 0.000700 2.996718e+51 \n", " Ni56 3.308150 1.416227e+55 \n", "... ... ... \n", "4 Co56 0.340356 1.457072e+54 \n", " Co56 32.752192 1.402130e+56 \n", " Co56 6.936438 2.969507e+55 \n", " Co56 0.601334 2.574325e+54 \n", " Co56 0.310592 1.329651e+54 \n", "\n", " decay_energy_erg \n", "shell_number isotope \n", "0 Ni56 2.674994e+42 \n", " Ni56 8.395366e+41 \n", " Ni56 9.835748e+41 \n", " Ni56 4.801271e+42 \n", " Ni56 2.269047e+46 \n", "... ... \n", "4 Co56 2.334487e+45 \n", " Co56 2.246459e+47 \n", " Co56 4.757674e+46 \n", " Co56 4.124523e+45 \n", " Co56 2.130335e+45 \n", "\n", "[485 rows x 8 columns]" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "isotope_decay_df = create_isotope_decay_df(cumulative_decay_df, gamma_ray_lines)\n", "isotope_decay_df" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we can use the data from the isotope decay DataFrame to find the energy production rate and total decay energy for the model and compare our results to the simplified analytical solution above." ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "def percent_error(observed, expected):\n", " return np.abs(expected - observed)/expected * 100" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "First we need to find the lifetimes of each isotope using the `get_taus` function." ] }, { "cell_type": "code", "execution_count": 22, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The lifetime of 56Ni is: 8.764372373400452 d\n", "The lifetime of 56Co is: 111.41933800785463 d\n", "The percent error for 56Ni lifetime is: 0.40485939317669306 %\n", "The percent error for 56Co lifetime is: 0.1072219297885327 %\n" ] } ], "source": [ "tau_dict = get_taus(raw_isotope_abundance)\n", "\n", "tau_co56 = tau_dict[0]['Co-56'] * u.s\n", "tau_ni56 = tau_dict[0]['Ni-56'] * u.s\n", "tau_ni56_days = tau_ni56.to(u.day)\n", "tau_co56_days = tau_co56.to(u.day)\n", "print(f\"The lifetime of 56Ni is: {tau_ni56_days}\")\n", "print(f\"The lifetime of 56Co is: {tau_co56_days}\")\n", "print(f\"The percent error for 56Ni lifetime is: {percent_error(tau_ni56_days.value, 8.8)} %\")\n", "print(f\"The percent error for 56Co lifetime is: {percent_error(tau_co56_days.value, 111.3)} %\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next we find the energy per decay or Q value for both $^{56}Ni$ and $^{56}Co$ using the data in our isotope decay DataFrame, following a similar method as above. Here the \"probabilities\" are found in the `radiation_intensity` column and are in the form of number of decays per 100 decays so we have to divide by 100. The energy for each transition is found in the `radiation_energy_keV` column." ] }, { "cell_type": "code", "execution_count": 20, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The energy per decay for Ni56 is: 1.7249433621 MeV\n", "The energy per decay for Co56 is: 3.8521702007573437 MeV\n", "The percent error for Ni56 energy per decay is: 1.4318078799999943 %\n", "The percent error for Co56 energy per decay is: 3.2753405028778486 %\n" ] } ], "source": [ "# For Ni56\n", "ni56 = isotope_decay_df.xs('Ni56', level='isotope')\n", "ni56_energy_per_decay_MeV = (ni56['radiation_intensity']/100 * (ni56['radiation_energy_keV']*u.keV.to(\"MeV\"))).sum()/5 * u.MeV\n", "\n", "# For Co56\n", "co56 = isotope_decay_df.xs('Co56', level='isotope')\n", "co56_energy_per_decay_MeV = (co56['radiation_intensity']/100 * (co56['radiation_energy_keV']* u.keV.to(\"MeV\"))).sum()/5 * u.MeV\n", "\n", "print(f\"The energy per decay for Ni56 is: {ni56_energy_per_decay_MeV}\")\n", "print(f\"The energy per decay for Co56 is: {co56_energy_per_decay_MeV}\")\n", "print(f\"The percent error for Ni56 energy per decay is: {percent_error(ni56_energy_per_decay_MeV.value, 1.75)} %\")\n", "print(f\"The percent error for Co56 energy per decay is: {percent_error(co56_energy_per_decay_MeV.value, 3.73)} %\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Next we can plug the $\\tau$ values and Q values into the energy production rate equation to get the coefficients in the second simplified energy production rate equation." ] }, { "cell_type": "code", "execution_count": 24, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The coefficient for Ni56 is: 6.280298471107531e+43 erg / s\n", "The coefficient for Co56 is: 1.4795292955420012e+43 erg / s\n", "The percent error for Ni56 coefficient is: 2.631031455697185 %\n", "The percent error for Co56 coefficient is: 2.0365031408276604 %\n" ] } ], "source": [ "outside_factor = (1/ (tau_co56 - tau_ni56)) * ((const.M_sun.to('g'))/(56 * 1.67e-24 * u.g))\n", "ni56_coeff = outside_factor * (ni56_energy_per_decay_MeV.to(\"erg\") * ((tau_co56/tau_ni56) - 1) - co56_energy_per_decay_MeV.to(\"erg\"))\n", "co56_coeff = outside_factor * co56_energy_per_decay_MeV.to(\"erg\")\n", "\n", "print(f\"The coefficient for Ni56 is: {ni56_coeff}\")\n", "print(f\"The coefficient for Co56 is: {co56_coeff}\")\n", "print(f\"The percent error for Ni56 coefficient is: {percent_error(ni56_coeff.value, 6.45e43)} %\")\n", "print(f\"The percent error for Co56 coefficient is: {percent_error(co56_coeff.value, 1.45e43)} %\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Then we can use these coefficients to find the energy production rate after 10 days." ] }, { "cell_type": "code", "execution_count": 26, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The energy production rate from TARDIS-HE is: 3.3591101069090006e+43 erg / s\n", "The percent error for energy production rate from TARDIS-HE is: 1.0789985177098882 %\n" ] } ], "source": [ "energy_production_rate_from_HE = ni56_coeff * np.exp(-time/tau_ni56_days) + co56_coeff * np.exp(-time/tau_co56_days) * (m_Ni56/ const.M_sun)\n", "\n", "print(f\"The energy production rate from TARDIS-HE is: {energy_production_rate_from_HE}\")\n", "print(f\"The percent error for energy production rate from TARDIS-HE is: {percent_error(energy_production_rate_from_HE.value, energy_production_rate.value)} %\")" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can also sum the decay energy column of the DataFrame to find the total energy produced." ] }, { "cell_type": "code", "execution_count": 28, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The total energy production from TARDIS-HE is: 1.9126597273910088e+50 erg\n", "The percent error for total energy production from TARDIS-HE is: 1.4673595432896005 %\n" ] } ], "source": [ "total_energy_production_HE = sum(isotope_decay_df[\"decay_energy_erg\"]) * u.erg\n", "\n", "print(f\"The total energy production from TARDIS-HE is: {total_energy_production_HE}\")\n", "print(f\"The percent error for total energy production from TARDIS-HE is: {percent_error(total_energy_production_HE.value, total_energy_production.value)} %\")" ] } ], "metadata": { "kernelspec": { "display_name": "tardis", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.12.4" } }, "nbformat": 4, "nbformat_minor": 2 }