In [56]:

    

%pylab 
from nugridpy import mesa as ms
import os


    

    




Using matplotlib backend: nbAgg
Populating the interactive namespace from numpy and matplotlib

In [85]:

    

# Before we start have a look at NuGrid Set1 model (Ritter+ 17)
ms.set_nugrid_path('/data/nugrid_apod2')
nsa=ms.star_log(mass=2,Z=0.0001)


    

    




nugrid_path = /data/nugrid_apod2
closest set is set1.5a (Z = 0.0001)
closest mass is 2.0
Using old star.logsa file ...
 reading ...100% 

Closing star.log  tool ...

In [86]:

    

# check we got the right model
nsa.header_attr
#nsa.get('surface_o16')[0]+nsa.get('surface_c12')[0]


    

    
Out[86]:





{'burn_min1': 50.0,
 'burn_min2': 1000.0,
 'c12_boundary_limit': 0.0001,
 'h1_boundary_limit': 0.0001,
 'he4_boundary_limit': 0.0001,
 'initial_mass': 2.0,
 'initial_z': 0.0001}

In [99]:

    

figure(101); shapes=['-',':','--','-.']; i=0
nsa.plot('model_number','log_LH',shape=shapes[mod(i,4)]); i+=1
nsa.plot('model_number','log_LHe',shape=shapes[mod(i,4)]); i+=1
title('Z='+str(nsa.header_attr['initial_z'])); ylim(2,8); xlim(0,40000); ylabel("log L shell")


    

    















    












    
Out[99]:





((0, 40000), <matplotlib.text.Text at 0x7fa514e0c7b8>)

In [101]:

    

nsa.kippenhahn(102,'model'); xlim(0,40000)


    

    















    












    
Out[101]:





(0, 40000)

In [58]:

    

from nugridpy import utils as ut


    

In [59]:

    

data_dir='/data/nugrid_apod2/scratch/MESA/Set_M2Zrange/'


    

In [60]:

    

!ls $data_dir


    

    




M2.00Z0.e-0  README.Set_template    plot_radius.py~    run_ssh.sh~
M2.00Z1.e-2  do_in_all_rundirs.sh   print_log_time.sh  vis2.py
M2.00Z1.e-3  do_in_all_rundirs.sh~  prof.py	       vis3.py
M2.00Z1.e-4  findCOgt1p15.py	    prof.py~	       vis3.py~
M2.00Z1.e-5  ipython_log.py	    rhoc-Tc.py	       vis_runs.py
M2.00Z1.e-6  kill_all.sh	    run.list	       work_template
M2.00Z1.e-7  mkidx.sh		    run.list.all
M2.00Z1.e-8  out.png		    run.list~
M2.00Z2.e-2  plot_radius.py	    run_ssh.sh

In [61]:

    

cases = [ case for case in os.listdir(data_dir) if  'M2.00' in case]
cases.sort()
cases[0], cases[-1] = cases[-1], cases[0] 
print(cases)


    

    




['M2.00Z2.e-2', 'M2.00Z1.e-2', 'M2.00Z1.e-3', 'M2.00Z1.e-4', 'M2.00Z1.e-5', 'M2.00Z1.e-6', 'M2.00Z1.e-7', 'M2.00Z1.e-8', 'M2.00Z0.e-0']

In [62]:

    

%%capture
sints=[]
for case in cases:
    sints.append(ms.star_log(os.path.join(data_dir,case,'LOGS')))


    

In [63]:

    

ZZ=[]
for s in sints:
    ZZ.append(s.header_attr['initial_z'])
print(ZZ  )


    

    




[0.02, 0.01, 0.001, 0.0001, 1e-05, 1e-06, 1e-07, 1e-08, 0.0]

In [64]:

    

figure(1)
for s in sints:
    s.hrd_new()


    

In [65]:

    

figure(2)
shapes=['-',':','--','-.']; i=0
for s in sints:
    s.plot('model_number','log_LH',shape=shapes[mod(i,4)],\
           legend='Z='+str(s.header_attr['initial_z'])); i+=1
    ylim(3,12); xlim(0,10000)


    

In [66]:

    

s.get('h1_boundary_mass')[-1]


    

    
Out[66]:





0.66470864254115014

In [67]:

    

mhfree = []
for s in sints:
    mhfree.append(s.get('h1_boundary_mass')[-1])


    

In [68]:

    

mhfree


    

    
Out[68]:





[0.61192592225791298,
 0.58814950288822265,
 0.57882809114290656,
 0.62697949181837875,
 0.63573186924937164,
 0.66042181784460163,
 0.67270867808645951,
 0.67629498240165953,
 0.66470864254115014]

In [70]:

    

figure(121)
plot(log10(array(ZZ)),mhfree,'o')
xlabel('$ \log Z$'); ylabel('$M_\mathrm{H-free}$')


    

    















    












    




/usr/local/lib/python3.5/dist-packages/ipykernel_launcher.py:2: RuntimeWarning: divide by zero encountered in log10
  






    
Out[70]:





<matplotlib.text.Text at 0x7fa5155392e8>

In [72]:

    

modstart=4600; s=sints[3]
s.kip_cont(ifig=15, modstart=modstart, modstop=5300, \
            t0_model=modstart, ixaxis='age', xres=2000, yres=2000, ylims=[0.58,0.62])
title('Z='+str(s.header_attr['initial_z']))


    

    




 creating color map burn ...100% 

 creating color map mix ...100% 

engenstyle was  full
mixstyle was  full

 finished preparing color map






    















    












    




plot versus age
plotting contours






    




/usr/local/lib/python3.5/dist-packages/matplotlib/contour.py:1518: UserWarning: Log scale: values of z <= 0 have been masked
  warnings.warn('Log scale: values of z <= 0 have been masked')






    




plotting patches
plotting abund boundaries






    
Out[72]:





<matplotlib.text.Text at 0x7fa5157ab198>

In [73]:

    

shapes=['-',':','--','-.']; i=5
for s in sints[3:]:
    close(i);figure(i)
    s.kippenhahn_CO(i,'model'); i+=1
    title('Z='+str(s.header_attr['initial_z']))
    close(i);figure(i)
    s.plot('model_number','log_LHe',shape=shapes[mod(i,3)]); i+=1
    s.plot('model_number','log_LH',shape=shapes[mod(i,3)]); i+=1
    ylim(3,12)
    title('Shell luminosity evolution, Z='+str(s.header_attr['initial_z']))


    

In [74]:

    

modstart=3090; modstop=4100; s=sints[4]; ifig=16; close(ifig)
s.kip_cont(ifig=ifig, modstart=modstart, modstop=modstop, \
            t0_model=modstart, ixaxis='age', xres=2000, yres=2000, ylims=[0.6,0.65])
title('Z='+str(s.header_attr['initial_z']))


    

    




 creating color map burn ...100% 

 creating color map mix ...100% 

engenstyle was  full
mixstyle was  full

 finished preparing color map






    















    












    




plot versus age
plotting contours






    




/usr/local/lib/python3.5/dist-packages/matplotlib/contour.py:1518: UserWarning: Log scale: values of z <= 0 have been masked
  warnings.warn('Log scale: values of z <= 0 have been masked')






    




plotting patches
plotting abund boundaries






    
Out[74]:





<matplotlib.text.Text at 0x7fa5149404a8>

In [102]:

    

figure(201); shapes=['-',':','--','-.']; i=0
s.plot('model_number','log_LH',shape=shapes[mod(i,4)]); i+=1
s.plot('model_number','log_LHe',shape=shapes[mod(i,4)]); i+=1
title('Z='+str(s.header_attr['initial_z'])); ylim(2,12); xlim(3300,3500); ylabel("log L shell")


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[102]:





<matplotlib.text.Text at 0x7fa514e918d0>

In [180]:

    

modstart=3090; modstop=4100; s=sints[4]; ifig=19; close(ifig)
s.kip_cont(ifig=ifig, modstart=modstart, modstop=modstop, \
             ixaxis='model_number', xres=2000, yres=2000, ylims=[0.6,0.65], xlims=[3300,3500])
title('Z='+str(s.header_attr['initial_z']))


    

    




 creating color map burn ...100% 

 creating color map mix ...100% 

engenstyle was  full
mixstyle was  full

 finished preparing color map






    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    




Closing profile tool ...
plot versus model number
plotting contours






    




/usr/local/lib/python3.5/dist-packages/matplotlib/contour.py:1518: UserWarning: Log scale: values of z <= 0 have been masked
  warnings.warn('Log scale: values of z <= 0 have been masked')






    




plotting patches
plotting abund boundaries






    
Out[180]:





<matplotlib.text.Text at 0x7fa50d662400>

In [185]:

    

ifig=20; close(ifig); figure(ifig)
plot(s.get('model_number'),log10(s.get('surface_n14')+s.get('surface_c12')+\
                           s.get('surface_o16')),'-',label='CNO')
plot(s.get('model_number'),log10(s.get('surface_c12')+s.get('surface_o16')),':',label='CO')
plot(s.get('model_number'),log10(s.get('surface_c12')),'-.',label='C')
title('Z='+str(s.header_attr['initial_z'])); legend(loc=0)


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[185]:





<matplotlib.legend.Legend at 0x7fa4d4ffce80>

In [77]:

    

# Z = 1e-6
modstart=2820; modstop=3850; s=sints[5]; ifig=17; close(ifig)
s.kip_cont(ifig=ifig, modstart=modstart, modstop=modstop, \
            t0_model=modstart, ixaxis='age', xres=2000, yres=2000, ylims=[0.62,0.66])
title('Z='+str(s.header_attr['initial_z']))


    

    




 creating color map burn ...100% 

 creating color map mix ...100% 

engenstyle was  full
mixstyle was  full

 finished preparing color map






    















    












    




plot versus age
plotting contours






    




/usr/local/lib/python3.5/dist-packages/matplotlib/contour.py:1518: UserWarning: Log scale: values of z <= 0 have been masked
  warnings.warn('Log scale: values of z <= 0 have been masked')






    




plotting patches
plotting abund boundaries






    
Out[77]:





<matplotlib.text.Text at 0x7fa514a0c668>

In [245]:

    

%%capture
s5s=[]; case='M2.00Z1.e-5'
for i in range(3):
    s5s.append(ms.star_log(os.path.join(data_dir,case,'LOGS'+str(i+1))))


    

In [269]:

    

figure(29)
for s in s5s:
    s.plot('model_number','log_LHe',shape=shapes[mod(i+3,3)], legend='$L_He$_'+str(i))
    s.plot('model_number','log_LH',shape=shapes[mod(i,3)],\
           legend='$L_H$_'+str(i)); i+=1
ylim(3,12) #; xlim(750332457.29442167, 750332458.28505421)


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[269]:





(3, 12)

In [246]:

    

figure(30)
for s in s5s:
    s.plot('star_age','log_LHe',shape=shapes[mod(i+3,3)], legend='$L_He$_'+str(i))
    s.plot('star_age','log_LH',shape=shapes[mod(i,3)],\
           legend='$L_H$_'+str(i)); i+=1
ylim(3,12); xlim(750332457.29442167, 750332458.28505421)


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[246]:





(750332457.2944217, 750332458.2850542)

In [247]:

    

s=s5s[-1]
t = s.get('star_age')
LH = s.get('log_LH'); LHe = s.get('log_LHe')


    

In [248]:

    

mod0=where(s.get('model_number')==15000)[0][0]
print(s.get('model_number')[mod0])


    

    




15000.0

In [273]:

    

21./3600


    

    
Out[273]:





0.005833333333333334

In [272]:

    

figure(28)
s=s5s[-1]
s.plot('model_number','log_LH',shape='o',legend='$L_H$_'); i+=1
s.plot('model_number','log_LHe',shape=shapes[1], legend='$L_He$_')
ylim(3,12)


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[272]:





(3, 12)

In [268]:

    

close(31); figure(31,figsize=(6, 4), dpi=125)
t0 = t[mod0]; tend = t[where(s.get('model_number')==32000)[0][0]]
plot(t-t[mod0],LH,shapes[0], label='$L_\mathrm{H}$')
plot(t-t[mod0],LHe,shapes[2], label='$L_\mathrm{He}$')
legend(loc=2)
ylim(3,11); xlim(0,tend-t0)
xlabel('t/yr'); ylabel('$\log \ L_\mathrm{shell}/L_\odot$')
savefig('TP_stellar_evolution_shell_L.pdf')


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    

In [139]:

    

print(k) for k in s.cols and


    

    
Out[139]:





{'burn_c': 230,
 'burn_n': 231,
 'burn_o': 232,
 'burn_qtop_1': 105,
 'burn_qtop_10': 123,
 'burn_qtop_11': 125,
 'burn_qtop_12': 127,
 'burn_qtop_13': 129,
 'burn_qtop_14': 131,
 'burn_qtop_15': 133,
 'burn_qtop_16': 135,
 'burn_qtop_17': 137,
 'burn_qtop_18': 139,
 'burn_qtop_19': 141,
 'burn_qtop_2': 107,
 'burn_qtop_20': 143,
 'burn_qtop_21': 145,
 'burn_qtop_22': 147,
 'burn_qtop_23': 149,
 'burn_qtop_24': 151,
 'burn_qtop_25': 153,
 'burn_qtop_26': 155,
 'burn_qtop_27': 157,
 'burn_qtop_28': 159,
 'burn_qtop_29': 161,
 'burn_qtop_3': 109,
 'burn_qtop_30': 163,
 'burn_qtop_31': 165,
 'burn_qtop_32': 167,
 'burn_qtop_33': 169,
 'burn_qtop_34': 171,
 'burn_qtop_35': 173,
 'burn_qtop_36': 175,
 'burn_qtop_37': 177,
 'burn_qtop_38': 179,
 'burn_qtop_39': 181,
 'burn_qtop_4': 111,
 'burn_qtop_40': 183,
 'burn_qtop_5': 113,
 'burn_qtop_6': 115,
 'burn_qtop_7': 117,
 'burn_qtop_8': 119,
 'burn_qtop_9': 121,
 'burn_type_1': 104,
 'burn_type_10': 122,
 'burn_type_11': 124,
 'burn_type_12': 126,
 'burn_type_13': 128,
 'burn_type_14': 130,
 'burn_type_15': 132,
 'burn_type_16': 134,
 'burn_type_17': 136,
 'burn_type_18': 138,
 'burn_type_19': 140,
 'burn_type_2': 106,
 'burn_type_20': 142,
 'burn_type_21': 144,
 'burn_type_22': 146,
 'burn_type_23': 148,
 'burn_type_24': 150,
 'burn_type_25': 152,
 'burn_type_26': 154,
 'burn_type_27': 156,
 'burn_type_28': 158,
 'burn_type_29': 160,
 'burn_type_3': 108,
 'burn_type_30': 162,
 'burn_type_31': 164,
 'burn_type_32': 166,
 'burn_type_33': 168,
 'burn_type_34': 170,
 'burn_type_35': 172,
 'burn_type_36': 174,
 'burn_type_37': 176,
 'burn_type_38': 178,
 'burn_type_39': 180,
 'burn_type_4': 110,
 'burn_type_40': 182,
 'burn_type_5': 112,
 'burn_type_6': 114,
 'burn_type_7': 116,
 'burn_type_8': 118,
 'burn_type_9': 120,
 'c12_boundary_lgRho': 195,
 'c12_boundary_lgT': 194,
 'c12_boundary_mass': 192,
 'c12_boundary_radius': 193,
 'center_c12': 213,
 'center_h1': 211,
 'center_he4': 212,
 'center_mg24': 217,
 'center_mu': 209,
 'center_n14': 214,
 'center_ne20': 216,
 'center_o16': 215,
 'center_si28': 218,
 'center_ye': 210,
 'cno': 228,
 'conv_mx1_bot': 9,
 'conv_mx1_top': 8,
 'conv_mx2_bot': 11,
 'conv_mx2_top': 10,
 'delta_mass': 4,
 'epsnuc_M_1': 96,
 'epsnuc_M_2': 97,
 'epsnuc_M_3': 98,
 'epsnuc_M_4': 99,
 'epsnuc_M_5': 100,
 'epsnuc_M_6': 101,
 'epsnuc_M_7': 102,
 'epsnuc_M_8': 103,
 'gravity': 205,
 'h1_boundary_lgRho': 187,
 'h1_boundary_lgT': 186,
 'h1_boundary_mass': 184,
 'h1_boundary_radius': 185,
 'he4_boundary_lgRho': 191,
 'he4_boundary_lgT': 190,
 'he4_boundary_mass': 188,
 'he4_boundary_radius': 189,
 'log_L': 202,
 'log_LH': 196,
 'log_LHe': 197,
 'log_LZ': 198,
 'log_Lnuc': 199,
 'log_R': 203,
 'log_Teff': 200,
 'log_abs_mdot': 5,
 'log_center_P': 208,
 'log_center_Rho': 207,
 'log_center_T': 206,
 'log_dt': 6,
 'log_g': 204,
 'luminosity': 201,
 'mix_qtop_1': 17,
 'mix_qtop_10': 35,
 'mix_qtop_11': 37,
 'mix_qtop_12': 39,
 'mix_qtop_13': 41,
 'mix_qtop_14': 43,
 'mix_qtop_15': 45,
 'mix_qtop_16': 47,
 'mix_qtop_17': 49,
 'mix_qtop_18': 51,
 'mix_qtop_19': 53,
 'mix_qtop_2': 19,
 'mix_qtop_20': 55,
 'mix_qtop_21': 57,
 'mix_qtop_22': 59,
 'mix_qtop_23': 61,
 'mix_qtop_24': 63,
 'mix_qtop_25': 65,
 'mix_qtop_26': 67,
 'mix_qtop_27': 69,
 'mix_qtop_28': 71,
 'mix_qtop_29': 73,
 'mix_qtop_3': 21,
 'mix_qtop_30': 75,
 'mix_qtop_31': 77,
 'mix_qtop_32': 79,
 'mix_qtop_33': 81,
 'mix_qtop_34': 83,
 'mix_qtop_35': 85,
 'mix_qtop_36': 87,
 'mix_qtop_37': 89,
 'mix_qtop_38': 91,
 'mix_qtop_39': 93,
 'mix_qtop_4': 23,
 'mix_qtop_40': 95,
 'mix_qtop_5': 25,
 'mix_qtop_6': 27,
 'mix_qtop_7': 29,
 'mix_qtop_8': 31,
 'mix_qtop_9': 33,
 'mix_type_1': 16,
 'mix_type_10': 34,
 'mix_type_11': 36,
 'mix_type_12': 38,
 'mix_type_13': 40,
 'mix_type_14': 42,
 'mix_type_15': 44,
 'mix_type_16': 46,
 'mix_type_17': 48,
 'mix_type_18': 50,
 'mix_type_19': 52,
 'mix_type_2': 18,
 'mix_type_20': 54,
 'mix_type_21': 56,
 'mix_type_22': 58,
 'mix_type_23': 60,
 'mix_type_24': 62,
 'mix_type_25': 64,
 'mix_type_26': 66,
 'mix_type_27': 68,
 'mix_type_28': 70,
 'mix_type_29': 72,
 'mix_type_3': 20,
 'mix_type_30': 74,
 'mix_type_31': 76,
 'mix_type_32': 78,
 'mix_type_33': 80,
 'mix_type_34': 82,
 'mix_type_35': 84,
 'mix_type_36': 86,
 'mix_type_37': 88,
 'mix_type_38': 90,
 'mix_type_39': 92,
 'mix_type_4': 22,
 'mix_type_40': 94,
 'mix_type_5': 24,
 'mix_type_6': 26,
 'mix_type_7': 28,
 'mix_type_8': 30,
 'mix_type_9': 32,
 'model_number': 1,
 'mx1_bot': 13,
 'mx1_top': 12,
 'mx2_bot': 15,
 'mx2_top': 14,
 'num_backups': 235,
 'num_retries': 234,
 'num_zones': 7,
 'pp': 227,
 'star_age': 2,
 'star_mass': 3,
 'surface_c12': 221,
 'surface_c13': 222,
 'surface_h1': 219,
 'surface_he4': 220,
 'surface_n14': 223,
 'surface_o16': 224,
 'total_mass_h1': 225,
 'total_mass_he4': 226,
 'tri_alfa': 229,
 'v_div_csound_surf': 233}

In [254]:

    

[ print(k) for k in s.cols.keys() if 'bound' in k]


    

    




c12_boundary_lgRho
h1_boundary_radius
h1_boundary_lgT
he4_boundary_radius
h1_boundary_lgRho
c12_boundary_mass
he4_boundary_lgT
c12_boundary_radius
he4_boundary_mass
c12_boundary_lgT
h1_boundary_mass
he4_boundary_lgRho






    
Out[254]:





[None, None, None, None, None, None, None, None, None, None, None, None]

In [252]:

    

figure(301);s5s[-1].plot('model_number','h1_boundary_mass')


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    

In [258]:

    

figure(309)  #s5s[-1].plot('model_number','h1_boundary_radius')
plot(s5s[-1].get('model_number'),\
     s5s[-1].get('h1_boundary_radius')*ast.rsun_cm/1.e5,'o')
title('radius of H-free core immediately before and during HIF')
ylabel('radius $M_\mathrm{H-free}$ / [Mm]') 
xlabel('model number')


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[258]:





<matplotlib.text.Text at 0x7fa50d101c18>

In [263]:

    

#amount of mass by which M_Hfree decreases at the first time step of H-ingestion
# 0.0060725578854978757 Msun

# delta_R by which the H-free core is reduced during the one time step when the HIF starts:
delta_R = 1000*(30.5 - 21.5)
print('delta_R H-free: ','%7.1fkm'%(delta_R))
print('delta_R H-free: ','%7.1fHp'%(delta_R/Hp))
# model number after which HIF happens: 
mod_phif  = 29420


    

    




delta_R H-free:   9000.0km
delta_R H-free:      2.6Hp

In [265]:

    

5200/Hp


    

    
Out[265]:





1.4917849083879469

In [253]:

    

figure(302);s.plot('model_number','log_dt');xlim((29349, 29448))


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[253]:





(29349, 29448)

In [266]:

    

# time step just before the HIF:
dt_evol = 10**-3.945 * 3.14e7/3600 # in hours
print("Time step just before HIF: ","%4.2fhr"%dt_evol)
dt_evol = 10**-5.27 * 3.14e7/60 # in min
print("Third time step into HIF: ","%4.2fmin"%dt_evol)
dt_evol = 10**-6.175 * 3.14e7 # in sec
print("Time step in developed HIF: ","%4.2fs"%dt_evol)


    

    




Time step just before HIF:  0.99hr
Third time step into HIF:  2.81min
Time step in developed HIF:  20.99s

In [175]:

    

case='M2.00Z1.e-5'   ; print(mod_phif) 
pphif=ms.mesa_profile(os.path.join(data_dir,case,'LOGS3'),mod_phif)


    

    




29420
983 in profiles.index file ...
Found and load nearest profile for cycle 29400
reading /data/nugrid_apod2/scratch/MESA/Set_M2Zrange/M2.00Z1.e-5/LOGS3/log892.data ...
 reading ...100% 

In [221]:

    

# look for variables containing velocity
[ print(k) for k in pphif.cols.keys() if 'pres' in k]


    

    




pressure
pressure_scale_height






    
Out[221]:





[None, None]

In [188]:

    

figure(303)
pphif.plot('mass','conv_vel_div_csound',logy=True)
xlim(0.6,0.65)


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[188]:





(0.6, 0.65)

In [215]:

    

figure(304)
pphif.plot('mass','log_conv_vel',shape='--')
xlim(0.6,0.65);ylim(3,6)


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[215]:





(3, 6)

In [238]:

    

# range for Ma numbers in convection zone just before HIF
sef = '%10.3e'
print(sef%(10**-2.8),sef%(10**-2.895))


    

    




 1.585e-03  1.274e-03

In [222]:

    

figure(306)
pphif.plot('mass','pressure_scale_height',logy=True)
xlim(0.6,0.65)


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[222]:





(0.6, 0.65)

In [234]:

    

# pressure scale height in convection zone
Hp = 1.e-5*ast.rsun_cm*10**-2.3     # km
print("pressure scale height: ","%6.1fkm"%Hp)


    

    




pressure scale height:  3485.8km

In [216]:

    

figure(305)
#pphif.plot('radius','log_conv_vel')
plot(pphif.get('radius')*ast.rsun_cm/1.e5,1.e-5*10**pphif.get('log_conv_vel'))
xlim(5900, 38800); ylim(0, 2)
xlabel('radius / [Mm]'); ylabel('$v_\mathrm{conv}$ / [km/s]')
title('convective velocity in PDCZ')


    

    




/usr/local/lib/python3.5/dist-packages/matplotlib/pyplot.py:524: RuntimeWarning: More than 20 figures have been opened. Figures created through the pyplot interface (`matplotlib.pyplot.figure`) are retained until explicitly closed and may consume too much memory. (To control this warning, see the rcParam `figure.max_open_warning`).
  max_open_warning, RuntimeWarning)






    















    












    
Out[216]:





<matplotlib.text.Text at 0x7fa4d48b7f98>

In [230]:

    

delta_R_conv = (30.38-11.25)*1e3 # km
v_conv = 1.5                     # km/s
t_conv = delta_R_conv/v_conv
print('convection zone crossing time: ',t_conv/3600.,"hr")
t_conv = 1.7*Hp/v_conv
print('convective turnover: ',t_conv/3600.,"hr")


    

    




convection zone crossing time:  3.5425925925925927 hr
convective turnover:  1.0973680105170474 hr

In [283]:

    

R_c = 15.e8      # cm
M_PDCZ = 2e-2   # Msun
M_c = 0.61      # Msun
Lh = 10**10.5   # Lsun
t_conv = 3600   # s

Lh*ast.lsun_erg_s*t_conv*R_c/(2.*ast.grav_const*M_PDCZ*M_c*ast.msun_g**2)


    

    
Out[283]:





0.1017687336739877

In [278]:

    

E_velconv = 0.5*M_PDCZ*ast.msun_g*(v_conv * 1.e5)**2


    

In [279]:

    

Lh*ast.lsun_erg_s*t_conv/E_velconv


    

    
Out[279]:





976520.7530952818

In [280]:

    

Lh*ast.lsun_erg_s*t_conv


    

    
Out[280]:





4.3703942174591065e+47

In [282]:

    

10**5.9*ast.lsun_erg_s*t_conv/E_velconv


    

    
Out[282]:





24.52909229787557

In [ ]: