%pylab ipympl
from nugridpy import mesa as ms
data_dir = '/data/ASDR/CSA/Stellar_models/CBM_M25/template/LOGS'
s24_templ = ms.star_log(data_dir)
close(1)
s24_templ.kip_cont(ifig=1,ixaxis='log_time_left',fsize=12,xlims=[6.9,-3],ylims=[0,10],xres=3000,yres=3000)
If you have never used an ipython notebook, then here are the few basic rules you need to know:
Return
key to create a newlineShift
-Return
to execute the block of code%pylab ipympl
from nugridpy import mesa as ms
from nugridpy import nugridse as mp
# default data location (try mirror if default
# is not available)
data_dir="/data/ASDR/NuGrid"
# mirror NuGrid data location - uncomment both of the following lines
# ![ ! -h /user/data ] && ln -s /data/nugrid_data /user/data
# data_dir = '/user'
ms.set_nugrid_path(data_dir)
mp.set_nugrid_path(data_dir)
By default MESA is putting out two types of data. History data provides the time evolution of scalar quantities, one per time step. This data can be accessed with the mesa.star_log
(or mesa.history_data
which is the same) class.
MESA also outputs profile data at select time steps. Profiles are available via the mesa_profile
class.
Initialise the 2 solar-mass Z=0.02 MESA stellar evolution model from set1.2 using the seeker method:
s=ms.star_log(mass=2,Z=0.02) # or ms.history_data which is just an alias
Let's now see what the model looks like by making some Kippenhahn diagrams
ifig=111;close(ifig);figure(ifig)
s.kippenhahn_CO(111,'model')
What's happening at the centre of the star? ...
ifig=101;close(ifig);figure(ifig)
s.tcrhoc()
axis([0,7,7,8.5])
... and at the surface?
ifig=102;close(ifig);figure(ifig)
s.hrd_new()
legend(loc='lower right').draw_frame(False)
ifig=107;close(ifig);figure(ifig)
s.kip_cont(ifig=ifig, boundaries=True,engenlevels=[1e2,1e6,1e12])
#s.kip_cont?
The mesa.mesa_profile
class provides access to the available profile data. We can again use the seeker method to just find the righ model by initial mass and metallicity. An additional argument num=
specifies the model number. This behaviour can be changed, see the doc string for details. Profiles are not available for all models, and the method will report which nearby model it has found.
prof=ms.mesa_profile(mass=2,Z=0.02,num=4500)
ifig=108;close(ifig);figure(ifig)
prof.plot('mass','c12',shape='-')
# available columns:
prof.cols
# get individual columns
prof.get('logT')
# get all header attributes
prof.header_attr
# This nice abundance panel plot needs some dependency fixing inside NuGridPy :-(
# ms.abu_profiles(prof)
Now that we have an overview of the stellar evolution of this $2M_\odot$ star let's have a closer look at the nucleosynthesis.
mppnp
post-processing data¶The NuGrid mppnp
code would post-process detailed MESA hdf5 output. This post-processing data is available via the NuGridPy nugridse
class.
Initialise the 2 solar-mass Z=0.02 NuGrid nucleosynthesis data from set1.2 using the seeker method:
pt=mp.se(mass=2,Z=0.02)
Each of the se file sets (in fact each of the dozens of hdf5 files that make up the data set for one mass/metallicty combination, or stellar evolution track) has three types of data contained in them:
data type access | content |
---|---|
pt.se.hattrs |
a header section that holds the header attributes, including units in the form of factors so that if applied with the quantities the result is in cgs units |
pt.se.cattrs |
for each cycle (or time step) the cycle attributes are a number of scalar global quantities, such as total mass or star age |
pt.se.dcols |
again, for each time step these are the vector quantities available, i.e. the data table columns; one of the data columns, iso_massf is in fact a matrix where the matrix columns are different species, i.e. a radial vector of species vectors |
pt.se.hattrs
pt.se.cattrs
pt.se.dcols
pt.get('temperature_unit')/1.e9
ifig=122;close(ifig);figure(ifig)
pt.plot('mass','temperature',fname=1000,shape=':')
ylabel('$T / 10^9 \mathrm{K}$')
xlabel('$m_r/M_\odot$')
# pt.abu_profile?
and the nuclides that are being created or destroyed in this part of the star:
species=['H-1','He-4','C-12','C-13','N-14','O-16']
ifig=108;close(ifig);figure(ifig)
pt.abu_profile(isos=species, ifig=ifig, fname=18000, logy=True)
ylim(-7,0)
xlim(0.552,0.556)
species=['H-1','C-12','C-13','N-14','Fe-56','Sr-86','Ba-138','Pb-206']
ifig=121;close(ifig);figure(ifig)
pt.abu_profile(isos=species, ifig=ifig, fname=45500, logy=True, colourblind=True)
ylim(-9,0)
xlim(0.603,0.6033)
title("Formation of the $^\mathsf{13}\mathsf{C}$ pocket: the partial H-$^\mathsf{12}\mathsf{C}$ zone")
ifig=124;close(ifig);figure(ifig)
pt.abu_profile(isos=species, ifig=ifig, fname=47222, logy=True, colourblind=False)
ylim(-9,0)
xlim(0.60312,0.6032)
title("Final s-process formed in $^\mathsf{13}\mathsf{C}$ pocket: end of intershell period")
# help
#pt.iso_abund?
ifig=123;close(ifig);figure(ifig)
pt.iso_abund(47221, stable=True,amass_range=[50,210], mass_range=[0.60312,0.6032], ylim=[-9, -2])
We can also see the impact on the isotopic abundance chart for that part of the star:
# pt.abu_chart?
#ifig=1233;close(ifig);figure(ifig)
pt.abu_chart(47220,mass_range=[0.60312,0.6032], plotaxis=[0, 80, 0, 60],\
ilabel=False,imlabel=False,boxstable=False)
You can extract the trajectory and initial abundances from this model at the key mass coordinate for use in a 1-zone PPN simulation in order to study, for example, the impact of a new rate
pt.abund_at_masscoordinate(26100,0.57685,online=True)