Stellar Evolution Model Generation Script
A Python script designed for creating a grid of stellar models based on different parameters such as metallicity (Z), alpha-enhancement (A), and stellar mass (M). The script utilizes the PYSEP library for stellar physics calculations, handling input ranges for the desired parameters, and generating pre-main sequence models. It supports parallel processing to evolve the models and aims to build a grid of models to interpolate effective temperature (Teff) and luminosity (Lum) based on composition and mass. Suitable for astrophysical research and study in stellar evolution.
from pysep.api import stellarModel as sm from pysep.api.parallel import pStellarModel from pysep.misc.runDB.dbmanage import model_DB from pysep.dm import filetypes as pyFts from pysep.io.nml.control.defaults import solar from pysep.io.nml.control.defaults import controlLow from pysep.io.nml.control.defaults import controlHigh from pysep.io.nml.physics.defaults import phys1 from pysep.io.nml.physics.defaults import phys1Low from pysep.io.nml.physics.defaults import phys1High from pysep.opac.opal.defaults import GS98hz from pysep.opac.tops import open_and_parse from pysep.newpoly.popFileToPolyModel import generate_prems_models from pysep.prems import defaults as premsD import argparse import re import os import shutil from tqdm import tqdm import numpy as np def get_grid_points(Zarange, Aarange, Marange): # TODO Use proper X inference from BBN for a given Z Zvalues = np.arange(*Zarange) Avalues = np.arange(*Aarange) Mvalues = np.arange(*Marange) elements = np.empty(shape=(Zvalues.shape[0]*Avalues.shape[0]*Mvalues.shape[0], 3)) # TODO do idx better idx = 0 for iz, z in enumerate(Zvalues): for ia, a in enumerate(Avalues): for im, m in enumerate(Mvalues): elements[idx] = np.array([z, a, m]) idx += 1 return elements def build_prems(mass, abun, initpath): path = generate_prems_models(abun, mass, initpath, compName="TEMP", DDAGE=1000, fnameT="TEMP_%e%m.gs98") premsf = pyFts.premsf_file(os.path.join(initpath, path[mass])) return premsf def build_model(X,Z,A,M,outDir,abun,premsPath): prems = build_prems(M, abun, premsPath) if M < 0.7: useCTL = controlLow.copy() usePHYS = phys1Low.copy() elif 0.7 <= M < 2: useCTL = solar.copy() usePHYS = phys1.copy() else: useCTL = controlHigh.copy() usePHYS = phys1High.copy() useCTL[0]['descrip'] = f"Evolution of {M} star" useCTL[1]['descrip'] = "Goal: Build grid of models to interpolate Teff Lum" useCTL['lbnout'] = True useCTL['nbn'] = 1 usePHYS.add_key("ATIME", [0.001, 0.002, 0.5, 0.02, 0.3, 0.0015, 0.1, 0.02, 0.4, 0.02, 0.02]) useCTL[1]['nmodls'] = 3 useCTL[0]['cmixla'] = A useCTL[1]['cmixla'] = A tmpModel = sm(outDir,useCTL, usePHYS, GS98hz, prems) tmpModel.update_composition(X=X, Z=Z, ID=0) return tmpModel def main(args): # gridPoints should be iterable of three elemenet iterable [Z, alpha, M]... gridPoints = get_grid_points(args.Zarange, args.Aarange, args.Marange) models = list() abun = open_and_parse(args.chem) for gsZ, gsA, gsM in tqdm(gridPoints): model = build_model(args.X, gsZ, gsA, gsM, args.outDir, abun, args.premsPath) models.append(model) DB = model_DB("TeffLumModelGrid.db", args.outDir, mode='wr') nModels = DB.enroll_models_to_dir(models, pbar=True) pModels = pStellarModel(nModels) pModels.pEvolve() pModels.pStash(data=True) if __name__ == "__main__": parser = argparse.ArgumentParser(description="Run a grid of models to " "generate interpolation " "coefficients for getting teff" " and lum from composition " "and mass") parser.add_argument("Zarange", type=float, nargs=3, help="Z arange, high, low, step") parser.add_argument("Aarange", type=float, nargs=3, help="Alpha arange, high, low, step") parser.add_argument("Marange", type=float, nargs=3, help="M arange, high, low, step") parser.add_argument("chem", type=str, help="path to chemical abundance file") parser.add_argument("--X", type=float, help="X to use", default=0.735) parser.add_argument("--outDir", type=str, help="path to save models too", default="models/") parser.add_argument("--premsPath", type=str, help="path to save premain sequence models too", default="premsTMP/") args = parser.parse_args() main(args)