Update number of shells in SimulationState

tardis/model/base.py

@@ -305,14 +305,14 @@ def __init__(
                 )
             )
 
-        elif len(t_radiative) == self.no_of_shells + 1:
+        elif len(t_radiative) == self.no_of_raw_shells + 1:
             t_radiative = t_radiative[
                 self.geometry.v_inner_boundary_index
                 + 1 : self.geometry.v_outer_boundary_index
                 + 1
             ]
         else:
-            assert len(t_radiative) == self.no_of_shells
+            assert len(t_radiative) == self.no_of_raw_shells
 
         if dilution_factor is None:
             dilution_factor = 0.5 * (
@@ -446,7 +446,7 @@ def volume(self):
 
     @property
     def no_of_shells(self):
-        return self.geometry.no_of_shells
+        return self.geometry.no_of_shells_active
 
     @property
     def no_of_raw_shells(self):
@@ -479,7 +479,8 @@ def from_config(cls, config, atom_data=None):
             t_radiative = temperature
         elif config.plasma.initial_t_rad > 0 * u.K:
             t_radiative = (
-                np.ones(geometry.no_of_shells + 1) * config.plasma.initial_t_rad
+                np.ones(geometry.no_of_shells_active + 1)
+                * config.plasma.initial_t_rad
             )
         else:
             t_radiative = None
@@ -600,7 +601,7 @@ def from_csvy(cls, config, atom_data=None):
                 density_0, time_0, time_explosion
             )
 
-        no_of_shells = geometry.no_of_shells
+        no_of_shells = geometry.no_of_shells_active
 
         # TODO -- implement t_radiative
         # t_radiative = None
@@ -628,9 +629,7 @@ def from_csvy(cls, config, atom_data=None):
                 )
 
         elif config.plasma.initial_t_rad > 0 * u.K:
-            t_radiative = (
-                np.ones(geometry.no_of_shells) * config.plasma.initial_t_rad
-            )
+            t_radiative = np.ones(no_of_shells) * config.plasma.initial_t_rad
         else:
             t_radiative = None
 
@@ -644,7 +643,7 @@ def from_csvy(cls, config, atom_data=None):
         if hasattr(csvy_model_config, "abundance"):
             abundances_section = csvy_model_config.abundance
             abundance, isotope_abundance = read_uniform_abundances(
-                abundances_section, geometry.no_of_shells
+                abundances_section, no_of_shells
             )
         else:
             index, abundance, isotope_abundance = parse_csv_abundances(

tardis/model/tests/test_csvy_model.py

@@ -1,6 +1,7 @@
 from pathlib import Path
 import numpy as np
 import pandas as pd
+import copy
 import numpy.testing as npt
 import tardis
 import os
@@ -28,15 +29,33 @@ def model_config_fnames(request, example_csvy_file_dir):
     return csvy_config_file, old_config_file
 
 
-def test_compare_models(model_config_fnames):
+@pytest.fixture(scope="module")
+def model_configs(model_config_fnames):
+    """Function to retrieve tardis Configuration objects for tests"""
+    csvy_config_file, old_config_file = model_config_fnames
+    csvy_config = Configuration.from_yaml(csvy_config_file)
+    old_config = Configuration.from_yaml(old_config_file)
+    return csvy_config, old_config
+
+
+@pytest.fixture(scope="module")
+def model_simulation_states(model_configs, kurucz_atomic_data):
+    """Function to retrieve tardis SimulationState objects for tests"""
+    csvy_config, old_config = model_configs
+    csvy_model = SimulationState.from_csvy(
+        csvy_config, atom_data=kurucz_atomic_data
+    )
+    config_model = SimulationState.from_config(
+        old_config, atom_data=kurucz_atomic_data
+    )
+    return csvy_model, config_model
+
+
+def test_compare_models(model_simulation_states):
     """Compare identical models produced by .from_config and
     .from_csvy to check that velocities, densities and abundances
     (pre and post decay) are the same"""
-    csvy_config_file, old_config_file = model_config_fnames
-    tardis_config = Configuration.from_yaml(csvy_config_file)
-    tardis_config_old = Configuration.from_yaml(old_config_file)
-    csvy_model = SimulationState.from_csvy(tardis_config)
-    config_model = SimulationState.from_config(tardis_config_old)
+    csvy_model, config_model = model_simulation_states
     csvy_model_props = csvy_model.get_properties().keys()
     config_model_props = config_model.get_properties().keys()
     npt.assert_array_equal(csvy_model_props, config_model_props)
@@ -76,12 +95,12 @@ def test_compare_models(model_config_fnames):
 
 
 @pytest.fixture(scope="module")
-def csvy_model_test_abundances(example_csvy_file_dir):
+def csvy_model_for_test(example_csvy_file_dir, kurucz_atomic_data):
     """Returns SimulationState to use to test abundances dataframes"""
     csvypath = example_csvy_file_dir / "csvy_model_to_test_abundances.yml"
     config = Configuration.from_yaml(csvypath)
-    csvy_model_test_abundances = SimulationState.from_csvy(config)
-    return csvy_model_test_abundances
+    csvy_model = SimulationState.from_csvy(config, atom_data=kurucz_atomic_data)
+    return config, csvy_model
 
 
 @pytest.fixture(scope="function")
@@ -110,16 +129,34 @@ def reference_input_dataframes():
     return reference_input_abundance, reference_input_isotopes
 
 
-def test_read_csvy_abundances(
-    csvy_model_test_abundances, reference_input_dataframes
+def test_dimensionality_after_update_v_inner_boundary(
+    csvy_model_for_test, kurucz_atomic_data
 ):
+    """Test that the dimensionality of SimulationState parameters after updating v_inner_boundary
+    in the context of csvy models specifically"""
+    config, csvy_model = csvy_model_for_test
+    new_config = copy.deepcopy(config)
+    new_config.model.v_inner_boundary = copy.deepcopy(csvy_model.velocity[1])
+    new_csvy_model = SimulationState.from_csvy(
+        new_config, atom_data=kurucz_atomic_data
+    )
+
+    assert new_csvy_model.no_of_raw_shells == csvy_model.no_of_raw_shells
+    assert new_csvy_model.no_of_shells == csvy_model.no_of_shells - 1
+    assert new_csvy_model.velocity.shape[0] == csvy_model.velocity.shape[0] - 1
+    assert new_csvy_model.density.shape[0] == csvy_model.density.shape[0] - 1
+    assert new_csvy_model.volume.shape[0] == csvy_model.volume.shape[0] - 1
+
+
+def test_read_csvy_abundances(csvy_model_for_test, reference_input_dataframes):
     """Test if model reads abundances and isotope abundances
     and constructs dataframes correctly before applying decay"""
     (
         reference_input_abundance,
         reference_input_isotopes,
     ) = reference_input_dataframes
 
+    _, csvy_model_test_abundances = csvy_model_for_test
     model_abundance_shape = csvy_model_test_abundances.raw_abundance.shape
     reference_input_shape = reference_input_abundance.shape
     assert model_abundance_shape == reference_input_shape
@@ -204,11 +241,10 @@ def N2(N1_0, lambda_1, lambda_2, t=4.0 * u.d):
     return reference_decayed_abundance
 
 
-def test_csvy_model_decay(
-    csvy_model_test_abundances, reference_decayed_abundance
-):
+def test_csvy_model_decay(csvy_model_for_test, reference_decayed_abundance):
     """Compare model abundance decay against decay calculations
     done by hand."""
+    _, csvy_model_test_abundances = csvy_model_for_test
     model_decayed_abundance_shape = csvy_model_test_abundances.abundance.shape
     reference_decayed_abundance_shape = reference_decayed_abundance.shape
     assert model_decayed_abundance_shape == reference_decayed_abundance_shape