Stellar model grids

Background and EEPs

Stellar model grids are typically constructed as a set of evolutionary tracks, where models of stellar evolution are run on grids of initial mass and metallicity, often with some other physical parameter varied as well (e.g., rotation, helium fraction, \(\alpha\)-abundance, etc.). Each of these evolutionary tracks predicts various physical properties (temperature, luminosity, etc.) of a star with given initial mass and metallicity, as a function of age.

It is also often of interest to re-organize these evolution track grids into “isochrones”—sets of stars at a range of masses, all with the same age. As described in this reference, in order to construct these isochrones, the time axis of each evolution track gets mapped into a new coordinate, called “equivalent evolutionary phase,” or EEP. The principle of the EEPs is to first identify physically significant stages in stellar evolution, and then subdivide each of these stages into a number of equal steps. This adaptive sampling enables accurate interpolation between evolution tracks even at ages when stars are evolving quickly, in the post-main sequence phases.

Previous versions of isochrones relied directly on these precomputed isochrone grids and interpolated between grid points in (mass, age, feh) space. This returned inaccurate results for post-MS stages of stellar evolution, and thus was not reliable for modeling evolved stars. However, beginning with v2.0, isochrones now implements all interpolation using EEPs. In addition, it provides direct access to the evolution track grids, in addition to precomputed isochrone grids. Note that version 2.0 includes only the MIST models; future updates will include more (e.g. PARSEC, YAPSI).

Model Grid Objects and Interpolation

Isochrones provides a simple and direct interface to full grids of stellar models. Upon first access, the grids are downloaded in original form, reorganized, and written to disk in binary format in order to load quickly with subsequent access. The grids are loaded as pandas dataframes with multi-level indexing that reflects the structure of the grids: evolution track grids are indexed by metallicity, initial mass, and EEP; and isochone grids by metallicity, age, and EEP.

[1]:
from isochrones.mist import MISTEvolutionTrackGrid, MISTIsochroneGrid

track_grid = MISTEvolutionTrackGrid()
track_grid.df.head()  # just show first few rows
[1]:
nu_max logg eep initial_mass radius logTeff mass density Mbol phase feh Teff logL delta_nu interpolated star_age age dt_deep
initial_feh initial_mass EEP
-4.0 0.1 1 143.524548 3.033277 1.0 0.1 1.593804 3.620834 0.1 0.034823 5.132871 -1.0 -3.978406 4176.707371 -0.157148 21.776686 False 13343.289397 4.125263 0.026168
2 145.419039 3.038935 2.0 0.1 1.583455 3.620769 0.1 0.035510 5.147664 -1.0 -3.978406 4176.085183 -0.163066 21.993078 False 14171.978264 4.151430 0.026121
3 147.409881 3.044805 3.0 0.1 1.572790 3.620702 0.1 0.036237 5.163015 -1.0 -3.978406 4175.435381 -0.169206 22.219791 False 15048.910447 4.177505 0.026016
4 149.499346 3.050886 4.0 0.1 1.561817 3.620631 0.1 0.037006 5.178922 -1.0 -3.978406 4174.757681 -0.175569 22.457004 False 15975.827275 4.203463 0.025996
5 151.703570 3.057203 5.0 0.1 1.550499 3.620558 0.1 0.037823 5.195452 -1.0 -3.978406 4174.049081 -0.182181 22.706349 False 16962.744747 4.229496 0.025996
[2]:
iso_grid = MISTIsochroneGrid()
iso_grid.df.head()  # just show first few rows
[2]:
eep age feh mass initial_mass radius density logTeff Teff logg logL Mbol delta_nu nu_max phase dm_deep
log10_isochrone_age_yr feh EEP
5.0 -4.0 35 35 5.0 -3.978406 0.100000 0.100000 1.106082 0.104184 3.617011 4140.105252 3.350571 -0.489734 5.964335 37.987066 299.346079 -1.0 0.002885
36 36 5.0 -3.978406 0.102885 0.102885 1.122675 0.102507 3.618039 4149.909661 3.347798 -0.472691 5.921728 37.739176 298.570836 -1.0 0.003573
37 37 5.0 -3.978406 0.107147 0.107147 1.147702 0.099921 3.619556 4164.436984 3.343658 -0.447471 5.858678 37.345115 297.180748 -1.0 0.004247
38 38 5.0 -3.978406 0.111379 0.111379 1.173015 0.097287 3.621062 4178.903372 3.339612 -0.422498 5.796244 36.923615 295.526946 -1.0 0.004217
39 39 5.0 -3.978406 0.115581 0.115581 1.198615 0.094627 3.622555 4193.289262 3.335660 -0.397776 5.734440 36.473151 293.589960 -1.0 0.004189

This generally contains only a subset of the original columns provided by the underlying grid, with standardized names. There are also additional computed columns, such as stellar radius and density. The full, original grids, can be found with the .df_orig attribute if desired:

[3]:
iso_grid.df_orig.head()  # just show first few rows
[3]:
EEP log10_isochrone_age_yr initial_mass star_mass star_mdot he_core_mass c_core_mass o_core_mass log_L log_L_div_Ledd ... nu_max acoustic_cutoff max_conv_vel_div_csound max_gradT_div_grada gradT_excess_alpha min_Pgas_div_P max_L_rad_div_Ledd e_thermal phase feh
log10_isochrone_age_yr feh EEP
5.0 -4.0 35 35 5.0 0.100000 0.100000 -1.455094e-13 0.0 0.0 0.0 -0.489734 -4.167035 ... 299.346079 2233.536029 0.127243 1.095544 0.0 0.999989 0.000016 3.002314e+46 -1.0 -4.0
36 36 5.0 0.102885 0.102885 -1.562027e-13 0.0 0.0 0.0 -0.472691 -4.129623 ... 298.570836 2228.014832 0.128938 1.101114 0.0 0.999989 0.000017 3.106838e+46 -1.0 -4.0
37 37 5.0 0.107147 0.107147 -1.707298e-13 0.0 0.0 0.0 -0.447471 -4.073591 ... 297.180748 2218.440338 0.130528 1.109114 0.0 0.999988 0.000019 3.264230e+46 -1.0 -4.0
38 38 5.0 0.111379 0.111379 -1.836256e-13 0.0 0.0 0.0 -0.422498 -4.017245 ... 295.526946 2207.403678 0.132657 1.116760 0.0 0.999987 0.000021 3.424460e+46 -1.0 -4.0
39 39 5.0 0.115581 0.115581 -1.949639e-13 0.0 0.0 0.0 -0.397776 -3.960633 ... 293.589960 2194.776391 0.134294 1.124050 0.0 0.999986 0.000022 3.587613e+46 -1.0 -4.0

5 rows × 80 columns

[4]:
iso_grid.df_orig.columns
[4]:
Index(['EEP', 'log10_isochrone_age_yr', 'initial_mass', 'star_mass',
       'star_mdot', 'he_core_mass', 'c_core_mass', 'o_core_mass', 'log_L',
       'log_L_div_Ledd', 'log_LH', 'log_LHe', 'log_LZ', 'log_Teff',
       'log_abs_Lgrav', 'log_R', 'log_g', 'log_surf_z', 'surf_avg_omega',
       'surf_avg_v_rot', 'surf_num_c12_div_num_o16', 'v_wind_Km_per_s',
       'surf_avg_omega_crit', 'surf_avg_omega_div_omega_crit',
       'surf_avg_v_crit', 'surf_avg_v_div_v_crit', 'surf_avg_Lrad_div_Ledd',
       'v_div_csound_surf', 'surface_h1', 'surface_he3', 'surface_he4',
       'surface_li7', 'surface_be9', 'surface_b11', 'surface_c12',
       'surface_c13', 'surface_n14', 'surface_o16', 'surface_f19',
       'surface_ne20', 'surface_na23', 'surface_mg24', 'surface_si28',
       'surface_s32', 'surface_ca40', 'surface_ti48', 'surface_fe56',
       'log_center_T', 'log_center_Rho', 'center_degeneracy', 'center_omega',
       'center_gamma', 'mass_conv_core', 'center_h1', 'center_he4',
       'center_c12', 'center_n14', 'center_o16', 'center_ne20', 'center_mg24',
       'center_si28', 'pp', 'cno', 'tri_alfa', 'burn_c', 'burn_n', 'burn_o',
       'c12_c12', 'delta_nu', 'delta_Pg', 'nu_max', 'acoustic_cutoff',
       'max_conv_vel_div_csound', 'max_gradT_div_grada', 'gradT_excess_alpha',
       'min_Pgas_div_P', 'max_L_rad_div_Ledd', 'e_thermal', 'phase', 'feh'],
      dtype='object')

Any property (or properties) of these grids can be interpolated to any value of the index parameters via the .interp method:

[5]:
track_grid.interp([-0.12, 1.01, 353.1], ['mass', 'radius', 'logg', 'Teff'])
[5]:
array([1.00983180e+00, 1.04691913e+00, 4.40266419e+00, 6.03383320e+03])

Similarly, the .interp_orig method interpolates any of the original columns by name:

[6]:
track_grid.interp_orig([-0.12, 1.01, 353.1], ['v_wind_Km_per_s'])
[6]:
array([2.87408918e-05])

Note that these interpolations are fast—30-40x faster than the equivalent interpolation in scipy, for evaluating at a single point:

[7]:
from scipy.interpolate import RegularGridInterpolator

grid = track_grid.interp.grid[:, :, :, 4]  # subgrid corresponding to radius
interp = RegularGridInterpolator(track_grid.interp.index_columns, grid)
assert track_grid.interp([-0.12, 1.01, 353.1], ['radius']) == interp([-0.12, 1.01, 353.1])
[8]:
%timeit interp([-0.12, 1.01, 353.1])
%timeit track_grid.interp([-0.12, 1.01, 353.1], ['radius'])
The slowest run took 14.68 times longer than the fastest. This could mean that an intermediate result is being cached.
100 loops, best of 3: 1.13 ms per loop
The slowest run took 5.04 times longer than the fastest. This could mean that an intermediate result is being cached.
10000 loops, best of 3: 12.5 µs per loop

In order to select a subset of these grids, you can use pandas multi-index magic:

[9]:
iso_grid.df.xs((9.0, 0.0), level=(0, 1)).head()  # just show first few rows
[9]:
eep age feh mass initial_mass radius density logTeff Teff logg logL Mbol delta_nu nu_max phase dm_deep
EEP
193 193 9.0 0.042799 0.100000 0.100000 0.126216 70.115891 3.460248 2885.680821 5.235913 -3.002129 12.245322 1045.120425 27533.135930 -1.0 0.003449
194 194 9.0 0.042805 0.103449 0.103449 0.129201 67.622476 3.462649 2901.679870 5.227759 -2.972222 12.170556 1025.261668 27025.521973 -1.0 0.004051
195 195 9.0 0.042814 0.108103 0.108103 0.133343 64.282466 3.465890 2923.411541 5.216743 -2.931855 12.069637 998.429682 26339.315728 -1.0 0.004742
196 196 9.0 0.042824 0.112932 0.112932 0.137791 60.858291 3.469256 2946.160133 5.205249 -2.889887 11.964717 970.463953 25622.963627 -1.0 0.004720
197 197 9.0 0.042834 0.117543 0.117544 0.142177 57.660112 3.472471 2968.048014 5.194272 -2.849811 11.864529 943.753111 24938.691941 -1.0 0.004735

Example visualization

Just for fun, let’s plot a few isochrones:

[10]:
import hvplot.pandas

# Select two isochrones from the grid
iso_df1 = iso_grid.df.xs((9.0, 0.0), level=(0, 1))
iso_df2 = iso_grid.df.xs((9.5, 0.0), level=(0, 1))

options = dict(invert_xaxis=True, legend_position='bottom_left')

# Isn't hvplot/holoviews great?
plot1 = iso_df1.hvplot.line('logTeff', 'logL', label='Log(age) = 9.0')
plot2 = iso_df2.hvplot.line('logTeff', 'logL', label='Log(age) = 9.5')
(plot1 * plot2).options(**options)