From 45fa45ba0ae4c31e66d1c1da37a69a5f0a886f8b Mon Sep 17 00:00:00 2001 From: Joshua Shields <54691495+jvshields@users.noreply.github.com> Date: Fri, 18 Oct 2024 16:04:11 -0400 Subject: [PATCH] add microturbulence to simulation (#222) * add microturbulence to simulation * fix broadening test * explanation of 0.0 in test --- stardis/config_schema.yml | 6 +++++- stardis/io/base.py | 2 ++ stardis/io/model/marcs.py | 8 ++++++-- stardis/model/base.py | 5 ++++- .../opacities/opacities_solvers/broadening.py | 15 +++++++++++---- .../opacities_solvers/tests/test_broadening.py | 4 ++-- 6 files changed, 30 insertions(+), 10 deletions(-) diff --git a/stardis/config_schema.yml b/stardis/config_schema.yml index f20d4b6b..4385ae3f 100644 --- a/stardis/config_schema.yml +++ b/stardis/config_schema.yml @@ -109,7 +109,11 @@ properties: - van_der_waals - radiation default: [] - description: Types of broadening to include. + description: Types of broadening to include. + disable_microturbulence: + type: boolean + default: false + description: Whether to disable microturbulence. broadening_range: type: [quantity, "null"] default: null diff --git a/stardis/io/base.py b/stardis/io/base.py index eed7260e..8716314b 100644 --- a/stardis/io/base.py +++ b/stardis/io/base.py @@ -75,6 +75,8 @@ def parse_config_to_model(config_fname, add_config_dict): stellar_model = raw_marcs_model.to_stellar_model( adata, final_atomic_number=config.model.final_atomic_number ) + if config.opacity.line.disable_microturbulence: + stellar_model.microturbulence = stellar_model.microturbulence * 0.0 elif config.model.type == "mesa": raw_mesa_model = read_mesa_model(Path(config.model.fname)) diff --git a/stardis/io/model/marcs.py b/stardis/io/model/marcs.py index 3ed94e0d..609a8aee 100644 --- a/stardis/io/model/marcs.py +++ b/stardis/io/model/marcs.py @@ -143,8 +143,12 @@ def to_stellar_model(self, atom_data, final_atomic_number=118): temperatures = ( self.data.t.values[::-1] * u.K ) # Flip data to move from innermost stellar point to surface - # First two none values are old fv_geometry and abundances which are replaced by the new structures. - return StellarModel(temperatures, marcs_geometry, marcs_composition) + return StellarModel( + temperatures, + marcs_geometry, + marcs_composition, + microturbulence=self.metadata["microturbulence"], + ) def read_marcs_metadata(fpath, gzipped=True): diff --git a/stardis/model/base.py b/stardis/model/base.py index 92a1b088..9d6f7f35 100644 --- a/stardis/model/base.py +++ b/stardis/model/base.py @@ -21,14 +21,17 @@ class StellarModel(HDFWriterMixin): Composition of the model. Includes density and atomic mass fractions. no_of_depth_points : int Class attribute to be easily accessible for initializing arrays that need to match the shape of the model. + microturbulence : float + Microturbulence in km/s. """ hdf_properties = ["temperatures", "geometry", "composition"] - def __init__(self, temperatures, geometry, composition): + def __init__(self, temperatures, geometry, composition, microturbulence=0.0): self.temperatures = temperatures self.geometry = geometry self.composition = composition + self.microturbulence = microturbulence @property def no_of_depth_points(self): diff --git a/stardis/radiation_field/opacities/opacities_solvers/broadening.py b/stardis/radiation_field/opacities/opacities_solvers/broadening.py index 059badf9..0e48dfb0 100644 --- a/stardis/radiation_field/opacities/opacities_solvers/broadening.py +++ b/stardis/radiation_field/opacities/opacities_solvers/broadening.py @@ -24,7 +24,7 @@ @numba.njit -def _calc_doppler_width(nu_line, temperature, atomic_mass): +def _calc_doppler_width(nu_line, temperature, atomic_mass, microturbulence): """ Calculates doppler width. https://ui.adsabs.harvard.edu/abs/2003rtsa.book.....R/ @@ -51,13 +51,18 @@ def _calc_doppler_width(nu_line, temperature, atomic_mass): return ( nu_line / SPEED_OF_LIGHT - * math.sqrt(2.0 * BOLTZMANN_CONSTANT * temperature / atomic_mass) + * ( + math.sqrt( + 2.0 * BOLTZMANN_CONSTANT * temperature / atomic_mass + + microturbulence**2 + ) + ) ) @numba.vectorize(nopython=True) -def calc_doppler_width(nu_line, temperature, atomic_mass): - return _calc_doppler_width(nu_line, temperature, atomic_mass) +def calc_doppler_width(nu_line, temperature, atomic_mass, microturbulence=0.0): + return _calc_doppler_width(nu_line, temperature, atomic_mass, microturbulence) @cuda.jit @@ -699,6 +704,7 @@ def calculate_broadening( stellar_model.composition.nuclide_masses.loc[lines.atomic_number].values[ :, np.newaxis ], + stellar_model.microturbulence.cgs.value, ) return gammas, doppler_widths @@ -727,6 +733,7 @@ def calculate_molecule_broadening( lines.nu.values[:, np.newaxis], stellar_model.temperatures.value, molecule_masses, + stellar_model.microturbulence.cgs.value, ) return gammas, doppler_widths diff --git a/stardis/radiation_field/opacities/opacities_solvers/tests/test_broadening.py b/stardis/radiation_field/opacities/opacities_solvers/tests/test_broadening.py index c97805ce..4781027a 100644 --- a/stardis/radiation_field/opacities/opacities_solvers/tests/test_broadening.py +++ b/stardis/radiation_field/opacities/opacities_solvers/tests/test_broadening.py @@ -65,9 +65,10 @@ def test_calc_doppler_width_sample_values( calc_doppler_width_sample_values_input_nu_line, calc_doppler_width_sample_values_input_temperature, calc_doppler_width_sample_values_input_atomic_mass, + 0.0, ), calc_doppler_width_sample_values_expected_result, - ) + ) # No microturbulence for legacy reasons - 0.0 @pytest.mark.skipif( @@ -519,7 +520,6 @@ def test_calc_gamma_van_der_waals_sample_values( calc_gamma_van_der_waals_sample_values_input_h_density, calc_gamma_van_der_waals_sample_values_expected_result, ): - assert np.allclose( calc_gamma_van_der_waals( calc_gamma_van_der_waals_sample_values_input_ion_number,