Source code for tardis.model.base

import logging
import os
from pathlib import Path

import numpy as np
from astropy import units as u

from tardis.io.configuration.config_reader import Configuration
from tardis.io.configuration.config_validator import validate_dict
from tardis.io.model.parse_composition_configuration import (
    parse_composition_from_config,
    parse_composition_from_csvy,
)
from tardis.io.model.parse_geometry_configuration import (
    parse_geometry_from_config,
    parse_geometry_from_csvy,
)
from tardis.io.model.parse_packet_source_configuration import (
    parse_packet_source_from_config,
)
from tardis.io.model.parse_radiation_field_configuration import (
    parse_radiation_field_state_from_config,
    parse_radiation_field_state_from_csvy,
)
from tardis.io.model.readers.csvy import (
    load_csvy,
)
from tardis.io.util import HDFWriterMixin
from tardis.util.base import is_valid_nuclide_or_elem

logger = logging.getLogger(__name__)


SCHEMA_DIR = (
    Path(__file__).parent / ".." / "io" / "configuration" / "schemas"
).resolve()


[docs] class SimulationState(HDFWriterMixin): """ An object that hold information about the individual shells. Parameters ---------- velocity : np.ndarray An array with n+1 (for n shells) velocities "cut" to the provided boundaries .. note:: To access the entire, "uncut", velocity array, use `raw_velocity` abundance : pd.DataFrame time_explosion : astropy.units.Quantity Time since explosion t_inner : astropy.units.Quantity elemental_mass: pd.Series luminosity_requested : astropy.units.quantity.Quantity t_radiative : astropy.units.Quantity Radiative temperature for the shells dilution_factor : np.ndarray If None, the dilution_factor will be initialized with the geometric dilution factor. v_boundary_inner : astropy.units.Quantity v_boundary_outer : astropy.units.Quantity raw_velocity : np.ndarray The complete array of the velocities, without being cut by `v_boundary_inner` and `v_boundary_outer` electron_densities : astropy.units.quantity.Quantity Attributes ---------- w : numpy.ndarray Shortcut for `dilution_factor` t_rad : astropy.units.quantity.Quantity Shortcut for `t_radiative` radius : astropy.units.quantity.Quantity r_inner : astropy.units.quantity.Quantity r_outer : astropy.units.quantity.Quantity r_middle : astropy.units.quantity.Quantity v_inner : astropy.units.quantity.Quantity v_outer : astropy.units.quantity.Quantity v_middle : astropy.units.quantity.Quantity density : astropy.units.quantity.Quantity volume : astropy.units.quantity.Quantity no_of_shells : int The number of shells as formed by `v_boundary_inner` and `v_boundary_outer` no_of_raw_shells : int """ hdf_properties = [ "t_inner", "dilution_factor", "t_radiative", "v_inner", "v_outer", "density", "r_inner", "time_explosion", "abundance", ] hdf_name = "simulation_state" def __init__( self, geometry, composition, radiation_field_state, time_explosion, packet_source, electron_densities=None, ): self.geometry = geometry self.composition = composition self.radiation_field_state = radiation_field_state self.time_explosion = time_explosion self._electron_densities = electron_densities self.packet_source = packet_source @property def blackbody_packet_source(self): return self.packet_source @property def t_inner(self): return self.packet_source.temperature @t_inner.setter def t_inner(self, value): self.packet_source.temperature = value @property def dilution_factor(self): return self.radiation_field_state.dilution_factor[ self.geometry.v_inner_boundary_index : self.geometry.v_outer_boundary_index ] @dilution_factor.setter def dilution_factor(self, new_dilution_factor): if len(new_dilution_factor) == self.no_of_shells: self.radiation_field_state.dilution_factor[ self.geometry.v_inner_boundary_index : self.geometry.v_outer_boundary_index ] = new_dilution_factor else: raise ValueError( "Trying to set dilution_factor for unmatching number" "of shells." ) @property def t_radiative(self): return self.radiation_field_state.temperature[ self.geometry.v_inner_boundary_index : self.geometry.v_outer_boundary_index ] @t_radiative.setter def t_radiative(self, new_t_radiative): if len(new_t_radiative) == self.no_of_shells: self.radiation_field_state.temperature[ self.geometry.v_inner_boundary_index : self.geometry.v_outer_boundary_index ] = new_t_radiative else: raise ValueError( "Trying to set t_radiative for different number of shells." ) @property def elemental_number_density(self): elemental_number_density = ( ( self.composition.elemental_mass_fraction * self.composition.density ) .divide(self.composition.element_masses, axis=0) .dropna() ) elemental_number_density = elemental_number_density.iloc[ :, self.geometry.v_inner_boundary_index : self.geometry.v_outer_boundary_index, ] elemental_number_density.columns = range( len(elemental_number_density.columns) ) return elemental_number_density @property def isotopic_number_density(self): isotopic_number_density = ( self.composition.isotopic_mass_fraction * self.composition.density ).divide( self.composition.isotope_masses.loc[ self.composition.isotopic_mass_fraction.index ] * u.u.to(u.g), axis=0, ) isotopic_number_density = isotopic_number_density.iloc[ :, self.geometry.v_inner_boundary_index : self.geometry.v_outer_boundary_index, ] isotopic_number_density.columns = range( len(isotopic_number_density.columns) ) return isotopic_number_density @property def radius(self): return self.time_explosion * self.velocity @property def v_boundary_inner(self): return self.geometry.v_inner_boundary @property def v_inner_boundary(self): return self.v_boundary_inner @property def v_boundary_outer(self): return self.geometry.v_outer_boundary @property def r_inner(self): return self.geometry.r_inner_active @property def r_outer(self): return self.geometry.r_outer_active @property def r_middle(self): return 0.5 * self.r_inner + 0.5 * self.r_outer @property def velocity(self): velocity = self.geometry.v_outer_active.copy() return velocity.insert(0, self.geometry.v_inner_active[0]) @property def v_inner(self): return self.geometry.v_inner_active @property def v_outer(self): return self.geometry.v_outer_active @property def v_middle(self): return 0.5 * self.v_inner + 0.5 * self.v_outer @property def density(self): return self.composition.density[ self.geometry.v_inner_boundary_index : self.geometry.v_outer_boundary_index ] @property def abundance(self): elemental_mass_fraction = ( self.composition.elemental_mass_fraction.copy() ) elemental_mass_fraction = elemental_mass_fraction.iloc[ :, self.geometry.v_inner_boundary_index : self.geometry.v_outer_boundary_index, ] elemental_mass_fraction.columns = range( len(elemental_mass_fraction.columns) ) return elemental_mass_fraction @property def volume(self): return ((4.0 / 3) * np.pi * (self.r_outer**3 - self.r_inner**3)).cgs @property def no_of_shells(self): return self.geometry.no_of_shells_active @property def no_of_raw_shells(self): return self.geometry.no_of_shells
[docs] @classmethod def from_config(cls, config, atom_data, legacy_mode_enabled=False): """ Create a new SimulationState instance from a Configuration object. Parameters ---------- config : tardis.io.config_reader.Configuration atom_data : tardis.io.AtomData Returns ------- SimulationState """ time_explosion = config.supernova.time_explosion.cgs geometry = parse_geometry_from_config(config, time_explosion) composition, electron_densities = parse_composition_from_config( atom_data, config, time_explosion, geometry ) packet_source = parse_packet_source_from_config( config, geometry, legacy_mode_enabled ) radiation_field_state = parse_radiation_field_state_from_config( config, geometry, dilution_factor=None, packet_source=packet_source, ) return cls( geometry=geometry, composition=composition, radiation_field_state=radiation_field_state, time_explosion=time_explosion, packet_source=packet_source, electron_densities=electron_densities, )
[docs] @classmethod def from_csvy(cls, config, atom_data=None, legacy_mode_enabled=False): """ Create a new SimulationState instance from a Configuration object. Parameters ---------- config : tardis.io.config_reader.Configuration atom_data : tardis.io.AtomData Returns ------- SimulationState """ CSVY_SUPPORTED_COLUMNS = { "velocity", "density", "t_rad", "dilution_factor", } if os.path.isabs(config.csvy_model): csvy_model_fname = config.csvy_model else: csvy_model_fname = os.path.join( config.config_dirname, config.csvy_model ) csvy_model_config, csvy_model_data = load_csvy(csvy_model_fname) csvy_schema_fname = SCHEMA_DIR / "csvy_model.yml" csvy_model_config = Configuration( validate_dict(csvy_model_config, schemapath=csvy_schema_fname) ) if hasattr(csvy_model_data, "columns"): abund_names = set( [ name for name in csvy_model_data.columns if is_valid_nuclide_or_elem(name) ] ) unsupported_columns = ( set(csvy_model_data.columns) - abund_names - CSVY_SUPPORTED_COLUMNS ) field_names = set( [field["name"] for field in csvy_model_config.datatype.fields] ) assert ( set(csvy_model_data.columns) - field_names == set() ), "CSVY columns exist without field descriptions" assert ( field_names - set(csvy_model_data.columns) == set() ), "CSVY field descriptions exist without corresponding csv data" if unsupported_columns != set(): logger.warning( "The following columns are " "specified in the csvy model file," f" but are IGNORED by TARDIS: {str(unsupported_columns)}" ) time_explosion = config.supernova.time_explosion.cgs electron_densities = None geometry = parse_geometry_from_csvy( config, csvy_model_config, csvy_model_data, time_explosion ) composition = parse_composition_from_csvy( atom_data, csvy_model_config, csvy_model_data, time_explosion, geometry, ) packet_source = parse_packet_source_from_config( config, geometry, legacy_mode_enabled ) radiation_field_state = parse_radiation_field_state_from_csvy( config, csvy_model_config, csvy_model_data, geometry, packet_source ) return cls( geometry=geometry, composition=composition, time_explosion=time_explosion, radiation_field_state=radiation_field_state, packet_source=packet_source, electron_densities=electron_densities, )