Integrating experimental model state

tardis/io/tests/data/non_uniform_isotope_abundance.dat

@@ -0,0 +1,4 @@
+Index H He Ni56 Ni58
+
+0 0.0 0.91 0.06 0.03
+1 0.4 0.51 0.05 0.04

tardis/io/tests/data/tardis_configv1_isotope_iabund.yml

@@ -0,0 +1,45 @@
+tardis_config_version: v1.0
+
+supernova:
+  luminosity_requested: 2.8e9 solLum
+  time_explosion: 13 day
+
+atom_data: kurucz_atom_pure_simple.h5
+
+model:
+  structure:
+    type: specific
+    velocity:
+      start: 1.1e4 km/s
+      stop: 2.0e4 km/s
+      num: 2
+    density:
+      type: branch85_w7
+  abundances:
+    type: file
+    filename: non_uniform_isotope_abundance.dat
+    filetype: custom_composition
+
+plasma:
+  ionization: lte
+  excitation: lte
+  radiative_rates_type: dilute-blackbody
+  line_interaction_type: macroatom
+
+montecarlo:
+  seed: 23111963
+  no_of_packets: 2.0e+5
+  iterations: 5
+  last_no_of_packets: 5.0e+5
+  no_of_virtual_packets: 5
+  convergence_strategy:
+    type: damped
+    damping_constant: 0.5
+    threshold: 0.05
+    lock_t_inner_cycles: 1
+    t_inner_update_exponent: -0.5
+
+spectrum:
+  start: 500 angstrom
+  stop: 20000 angstrom
+  num: 10000

tardis/model/base.py

@@ -8,7 +8,6 @@
 from radioactivedecay.utils import Z_DICT
 
 from tardis.util.base import quantity_linspace, is_valid_nuclide_or_elem
-from tardis.io.atom_data import AtomData
 from tardis.io.parsers.csvy import load_csvy
 from tardis.io.model_reader import (
     read_density_file,
@@ -119,13 +118,15 @@ class Composition:
     density : astropy.units.quantity.Quantity
         An array of densities for each shell.
     isotopic_mass_fraction : pd.DataFrame
-    atomic_mass : pandas.core.series.Series
+    atomic_mass : pd.DataFrame
     atomic_mass_unit: astropy.units.Unit
 
     Attributes
     ----------
-    number_density : pd.DataFrame
-        Number density of each isotope in each shell.
+    atomic_mass : pd.DataFrame
+        Atomic mass of elements calculated for each shell.
+    elemental_number_density : pd.DataFrame
+        Number density of each element in each shell.
     """
 
     def __init__(
@@ -184,6 +185,8 @@ def __init__(self, composition, geometry, time_explosion):
 
     @property
     def mass(self):
+        """Mass calculated using the formula:
+        mass_fraction * density * volume"""
         return (
             self.composition.elemental_mass_fraction
             * self.composition.density
@@ -192,6 +195,8 @@ def mass(self):
 
     @property
     def number(self):
+        """Number calculated using the formula:
+        mass / atomic_mass"""
         if self.composition.atomic_mass is None:
             raise AttributeError(
                 "ModelState was not provided elemental masses."
@@ -215,6 +220,7 @@ class Radial1DModel(HDFWriterMixin):
     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
@@ -669,6 +675,7 @@ def from_config(cls, config, atom_data=None):
         Parameters
         ----------
         config : tardis.io.config_reader.Configuration
+        atom_data : tardis.io.AtomData
 
         Returns
         -------
@@ -787,6 +794,7 @@ def from_csvy(cls, config, atom_data=None):
         Parameters
         ----------
         config : tardis.io.config_reader.Configuration
+        atom_data : tardis.io.AtomData
 
         Returns
         -------

tardis/model/tests/test_base.py

@@ -388,6 +388,61 @@ def test_model_state_number(simulation_verysimple, index, expected):
     assert_almost_equal((model_state.number).loc[index], expected, decimal=-47)
 
 
+@pytest.fixture
+def non_uniform_model_state(atomic_dataset):
+    filename = "tardis_configv1_isotope_iabund.yml"
+    config = Configuration.from_yaml(data_path(filename))
+    atom_data = atomic_dataset
+    model = Radial1DModel.from_config(config, atom_data=atom_data)
+    return model.model_state
+
+
+@pytest.mark.parametrize(
+    ("index", "expected"),
+    [
+        ((1, 0), 1.67378172e-24),
+        ((28, 0), 9.51707707e-23),
+        ((28, 1), 9.54725917e-23),
+    ],
+)
+def test_radial_1d_model_atomic_mass(non_uniform_model_state, index, expected):
+    atomic_mass = non_uniform_model_state.composition.atomic_mass
+
+    assert_almost_equal(
+        atomic_mass.loc[index],
+        expected,
+        decimal=30,
+    )
+
+
+class TestModelStateFromNonUniformAbundances:
+    @pytest.fixture
+    def model_state(self, non_uniform_model_state):
+        return non_uniform_model_state
+
+    def test_atomic_mass(self, model_state):
+        atomic_mass = model_state.composition.atomic_mass
+        assert_almost_equal(atomic_mass.loc[(1, 0)], 1.67378172e-24, decimal=30)
+        assert_almost_equal(
+            atomic_mass.loc[(28, 0)], 9.51707707e-23, decimal=30
+        )
+        assert_almost_equal(
+            atomic_mass.loc[(28, 1)], 9.54725917e-23, decimal=30
+        )
+
+    def test_elemental_number_density(self, model_state):
+        number = model_state.composition.elemental_number_density
+        assert_almost_equal(number.loc[(1, 0)], 0)
+        assert_almost_equal(number.loc[(28, 0)], 10825427.035, decimal=2)
+        assert_almost_equal(number.loc[(28, 1)], 1640838.763, decimal=2)
+
+    def test_number(self, model_state):
+        number = model_state.number
+        assert_almost_equal(number.loc[(1, 0)], 0)
+        assert_almost_equal(number.loc[(28, 0)], 1.53753476e53, decimal=-47)
+        assert_almost_equal(number.loc[(28, 1)], 4.16462779e52, decimal=-47)
+
+
 ###
 # Save and Load
 ###