From b5854bc6cea01d67d398115f8a551b3cbd10f204 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Tue, 26 Jul 2016 19:33:15 +0300 Subject: [PATCH 01/15] Change MontecarloRunner so it doesn't depend on a Configuration object. --- tardis/montecarlo/base.py | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index 9ce42a283e6..efa6df53da0 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -53,18 +53,18 @@ def _initialize_estimator_arrays(self, no_of_shells, tau_sobolev_shape): self.Edotlu_estimator = np.zeros(tau_sobolev_shape) - def _initialize_geometry_arrays(self, structure): + def _initialize_geometry_arrays(self, model): """ Generate the cgs like geometry arrays for the montecarlo part Parameters ---------- - structure: ~ConfigurationNameSpace + model : model.Radial1DModel """ - self.r_inner_cgs = structure.r_inner.to('cm').value - self.r_outer_cgs = structure.r_outer.to('cm').value - self.v_inner_cgs = structure.v_inner.to('cm/s').value + self.r_inner_cgs = model.r_inner.to('cm').value + self.r_outer_cgs = model.r_outer.to('cm').value + self.v_inner_cgs = model.v_inner.to('cm/s').value def _initialize_packets(self, T, no_of_packets, no_of_virtual_packets=None): nus, mus, energies = self.packet_source.create_packets(T, no_of_packets) @@ -123,10 +123,10 @@ def run(self, model, no_of_packets, no_of_virtual_packets=0, nthreads=1,last_run :return: """ self.time_of_simulation = model.time_of_simulation - self.volume = model.tardis_config.structure.volumes + self.volume = model.volume self._initialize_estimator_arrays(self.volume.shape[0], model.plasma.tau_sobolevs.shape) - self._initialize_geometry_arrays(model.tardis_config.structure) + self._initialize_geometry_arrays(model) self._initialize_packets(model.t_inner.value, no_of_packets) From 5f83ef539c5b1d94764cf1a75ae214e874e7bc1a Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Tue, 26 Jul 2016 19:41:31 +0300 Subject: [PATCH 02/15] Add a plasma parameter to MontecarloRunner.run The MontecarloRunner now needs a separate reference to the plasma since the restructured Radial1DModel doesn't have the plasma as an attribute anymore. --- tardis/montecarlo/base.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index efa6df53da0..add7107e297 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -110,7 +110,7 @@ def legacy_update_spectrum(self, no_of_virtual_packets): self.spectrum_virtual.update_luminosity( self.montecarlo_virtual_luminosity) - def run(self, model, no_of_packets, no_of_virtual_packets=0, nthreads=1,last_run=False): + def run(self, model, plasma, no_of_packets, no_of_virtual_packets=0, nthreads=1,last_run=False): """ Running the TARDIS simulation @@ -118,6 +118,7 @@ def run(self, model, no_of_packets, no_of_virtual_packets=0, nthreads=1,last_run ---------- :param model: + :param plasma: :param no_of_virtual_packets: :param nthreads: :return: @@ -125,7 +126,7 @@ def run(self, model, no_of_packets, no_of_virtual_packets=0, nthreads=1,last_run self.time_of_simulation = model.time_of_simulation self.volume = model.volume self._initialize_estimator_arrays(self.volume.shape[0], - model.plasma.tau_sobolevs.shape) + plasma.tau_sobolevs.shape) self._initialize_geometry_arrays(model) self._initialize_packets(model.t_inner.value, From 9f4385ec00740cc52be07d5c48b7bd7ba47f08a8 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Thu, 4 Aug 2016 02:54:40 +0300 Subject: [PATCH 03/15] Move and adapt calculate_j_blues from the old Radial1DModel to MontecarloRunner. --- tardis/montecarlo/base.py | 41 +++++++++++++++++++++++++++++++++++++++ 1 file changed, 41 insertions(+) diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index add7107e297..5a095720aa2 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -8,6 +8,7 @@ from spectrum import TARDISSpectrum from tardis.montecarlo import montecarlo, packet_source +from tardis.util import intensity_black_body from tardis.io.util import to_hdf import numpy as np @@ -284,6 +285,46 @@ def calculate_f_nu(self, frequency): def calculate_f_lambda(self, wavelength): pass + def calculate_j_blues(self, model, plasma, atom_data, + init_detailed_j_blues=False): + nus = atom_data.lines.nu.values + radiative_rates_type = plasma.radiative_rates_type + w_epsilon = plasma.w_epsilon + + if radiative_rates_type == 'blackbody': + logger.info('Calculating J_blues for radiative_rates_type=lte') + j_blues = intensity_black_body(nus[np.newaxis].T, model.t_rad.value) + return pd.DataFrame( + j_blues, index=atom_data.lines.index, + columns=np.arange(len(model.t_rad))) + + elif radiative_rates_type == 'dilute-blackbody' or init_detailed_j_blues: + logger.info('Calculating J_blues for radiative_rates_type=dilute-blackbody') + j_blues = model.w * intensity_black_body(nus[np.newaxis].T, model.t_rad.value) + return pd.DataFrame( + j_blues, index=atom_data.lines.index, + columns=np.arange(len(model.t_rad))) + + elif radiative_rates_type == 'detailed': + logger.info('Calculating J_blues for radiate_rates_type=detailed') + j_blues_norm_factor = (const.c.cgs * model.time_explosion / + (4 * np.pi * model.time_of_simulation * model.volume)) + j_blues = pd.DataFrame( + self.j_blue_estimator * + j_blues_norm_factor.value, + index=atom_data.lines.index, + columns=np.arange(len(model.t_rad))) + + for i in xrange(model.no_of_shells): + zero_j_blues = j_blues[i] == 0.0 + j_blues[i][zero_j_blues] = ( + w_epsilon * intensity_black_body( + atom_data.lines.nu[zero_j_blues].values, + model.t_rad.value[i])) + + else: + raise ValueError('radiative_rates_type type unknown - %s', radiative_rates_type) + def to_hdf(self, path_or_buf, path=''): """ Store the runner to an HDF structure. From fc8da950ee5b62a9be9d485027fa3bf44dcaaa88 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Thu, 4 Aug 2016 17:20:31 +0300 Subject: [PATCH 04/15] Add the time_of_simulation property to the MontecarloRunner. --- tardis/montecarlo/base.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index 5a095720aa2..6d7ad3ee31a 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -124,7 +124,7 @@ def run(self, model, plasma, no_of_packets, no_of_virtual_packets=0, nthreads=1, :param nthreads: :return: """ - self.time_of_simulation = model.time_of_simulation + self.time_of_simulation = self.calculate_time_of_simulation(model) self.volume = model.volume self._initialize_estimator_arrays(self.volume.shape[0], plasma.tau_sobolevs.shape) @@ -325,6 +325,13 @@ def calculate_j_blues(self, model, plasma, atom_data, else: raise ValueError('radiative_rates_type type unknown - %s', radiative_rates_type) + def calculate_time_of_simulation(self, model): + luminosity_inner = ( + 4 * np.pi * const.sigma_sb.cgs * model.r_inner[0] ** 2 * + model.t_inner ** 4).to('erg/s') + + return 1.0 * u.erg / luminosity_inner + def to_hdf(self, path_or_buf, path=''): """ Store the runner to an HDF structure. From f777f8fe33ea3c5e712bc6861d111d8c4016b2ea Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Fri, 5 Aug 2016 22:33:12 +0300 Subject: [PATCH 05/15] Add the plasma as an argument to montecarlo.montecarlo_radial1d. --- tardis/montecarlo/base.py | 2 +- tardis/montecarlo/montecarlo.pyx | 18 +++++++++--------- 2 files changed, 10 insertions(+), 10 deletions(-) diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index 6d7ad3ee31a..e661d9feb4b 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -134,7 +134,7 @@ def run(self, model, plasma, no_of_packets, no_of_virtual_packets=0, nthreads=1, no_of_packets) montecarlo.montecarlo_radial1d( - model, self, virtual_packet_flag=no_of_virtual_packets, + model, plasma, self, virtual_packet_flag=no_of_virtual_packets, nthreads=nthreads,last_run=last_run) # Workaround so that j_blue_estimator is in the right ordering # They are written as an array of dimension (no_of_shells, no_of_lines) diff --git a/tardis/montecarlo/montecarlo.pyx b/tardis/montecarlo/montecarlo.pyx index 8e412fde4b3..29244ee4748 100644 --- a/tardis/montecarlo/montecarlo.pyx +++ b/tardis/montecarlo/montecarlo.pyx @@ -93,7 +93,7 @@ cdef extern from "src/cmontecarlo.h": -cdef initialize_storage_model(model, runner, storage_model_t *storage): +cdef initialize_storage_model(model, plasma, runner, storage_model_t *storage): """ Initializing the storage struct. @@ -120,10 +120,10 @@ cdef initialize_storage_model(model, runner, storage_model_t *storage): storage.inverse_time_explosion = 1.0 / storage.time_explosion #electron density storage.electron_densities = PyArray_DATA( - model.plasma.electron_densities.values) + plasma.electron_densities.values) runner.inverse_electron_densities = ( - 1.0 / model.plasma.electron_densities.values) + 1.0 / plasma.electron_densities.values) storage.inverse_electron_densities = PyArray_DATA( runner.inverse_electron_densities) # Switch for continuum processes @@ -149,7 +149,7 @@ cdef initialize_storage_model(model, runner, storage_model_t *storage): storage.no_of_lines = model.atom_data.lines.nu.values.size storage.line_list_nu = PyArray_DATA(model.atom_data.lines.nu.values) runner.line_lists_tau_sobolevs = ( - model.plasma.tau_sobolevs.values.flatten(order='F') + plasma.tau_sobolevs.values.flatten(order='F') ) storage.line_lists_tau_sobolevs = PyArray_DATA( runner.line_lists_tau_sobolevs @@ -166,13 +166,13 @@ cdef initialize_storage_model(model, runner, storage_model_t *storage): # macro atom & downbranch if storage.line_interaction_id >= 1: runner.transition_probabilities = ( - model.plasma.transition_probabilities.values.flatten(order='F') + plasma.transition_probabilities.values.flatten(order='F') ) storage.transition_probabilities = PyArray_DATA( runner.transition_probabilities ) storage.transition_probabilities_nd = ( - model.plasma.transition_probabilities.values.shape[0]) + plasma.transition_probabilities.values.shape[0]) storage.line2macro_level_upper = PyArray_DATA( model.atom_data.lines_upper2macro_reference_idx) storage.macro_block_references = PyArray_DATA( @@ -218,10 +218,10 @@ cdef initialize_storage_model(model, runner, storage_model_t *storage): storage.reflective_inner_boundary = model.tardis_config.montecarlo.enable_reflective_inner_boundary storage.inner_boundary_albedo = model.tardis_config.montecarlo.inner_boundary_albedo # Data for continuum implementation - cdef np.ndarray[double, ndim=1] t_electrons = model.plasma.t_electrons + cdef np.ndarray[double, ndim=1] t_electrons = plasma.t_electrons storage.t_electrons = t_electrons.data -def montecarlo_radial1d(model, runner, int_type_t virtual_packet_flag=0, +def montecarlo_radial1d(model, plasma, runner, int_type_t virtual_packet_flag=0, int nthreads=4,last_run=False): """ Parameters @@ -253,7 +253,7 @@ def montecarlo_radial1d(model, runner, int_type_t virtual_packet_flag=0, cdef storage_model_t storage - initialize_storage_model(model, runner, &storage) + initialize_storage_model(model, plasma, runner, &storage) montecarlo_main_loop(&storage, virtual_packet_flag, nthreads, model.tardis_config.montecarlo.seed) From 097454c1677af9da941f00da61ff4fa809404e26 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Fri, 5 Aug 2016 22:39:47 +0300 Subject: [PATCH 06/15] Remove every Configuration object usage from the montecarlo.pyx --- tardis/montecarlo/base.py | 11 ++++++++++- tardis/montecarlo/montecarlo.pyx | 33 +++++++++++++++----------------- 2 files changed, 25 insertions(+), 19 deletions(-) diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index e661d9feb4b..ea85dd6bddd 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -29,9 +29,18 @@ class MontecarloRunner(object): t_rad_estimator_constant = ((np.pi**4 / (15 * 24 * zeta(5, 1))) * (const.h / const.k_B)).cgs.value - def __init__(self, seed, spectrum_frequency, distance=None): + def __init__(self, seed, spectrum_frequency, virtual_spectrum_range, + sigma_thomson, enable_reflective_inner_boundary, + inner_boundary_albedo, line_interaction_type, distance=None): + + self.seed = seed self.packet_source = packet_source.BlackBodySimpleSource(seed) self.spectrum_frequency = spectrum_frequency + self.virtual_spectrum_range = virtual_spectrum_range # TODO: linspace handling + self.sigma_thomson = sigma_thomson + self.enable_reflective_inner_boundary = enable_reflective_inner_boundary + self.inner_boundary_albedo = inner_boundary_albedo + self.line_interaction_type = line_interaction_type self.spectrum = TARDISSpectrum(spectrum_frequency, distance) self.spectrum_virtual = TARDISSpectrum(spectrum_frequency, distance) self.spectrum_reabsorbed = TARDISSpectrum(spectrum_frequency, distance) diff --git a/tardis/montecarlo/montecarlo.pyx b/tardis/montecarlo/montecarlo.pyx index 29244ee4748..d58b0548554 100644 --- a/tardis/montecarlo/montecarlo.pyx +++ b/tardis/montecarlo/montecarlo.pyx @@ -105,8 +105,7 @@ cdef initialize_storage_model(model, plasma, runner, storage_model_t *storage): storage.packet_energies = PyArray_DATA(runner.input_energy) # Setup of structure - structure = model.tardis_config.structure - storage.no_of_shells = structure.no_of_shells + storage.no_of_shells = model.no_of_shells storage.r_inner = PyArray_DATA(runner.r_inner_cgs) @@ -115,8 +114,7 @@ cdef initialize_storage_model(model, plasma, runner, storage_model_t *storage): # Setup the rest # times - storage.time_explosion = model.tardis_config.supernova.time_explosion.to( - 's').value + storage.time_explosion = model.time_explosion.to('s').value storage.inverse_time_explosion = 1.0 / storage.time_explosion #electron density storage.electron_densities = PyArray_DATA( @@ -161,7 +159,7 @@ cdef initialize_storage_model(model, plasma, runner, storage_model_t *storage): runner.Edotlu_estimator) storage.line_interaction_id = runner.get_line_interaction_id( - model.tardis_config.plasma.line_interaction_type) + runner.line_interaction_type) # macro atom & downbranch if storage.line_interaction_id >= 1: @@ -203,20 +201,20 @@ cdef initialize_storage_model(model, plasma, runner, storage_model_t *storage): storage.js = PyArray_DATA(runner.j_estimator) storage.nubars = PyArray_DATA(runner.nu_bar_estimator) - storage.spectrum_start_nu = model.tardis_config.spectrum.frequency.value.min() - storage.spectrum_end_nu = model.tardis_config.spectrum.frequency.value.max() - storage.spectrum_virt_start_nu = model.tardis_config.montecarlo.virtual_spectrum_range.end.to('Hz', units.spectral()).value - storage.spectrum_virt_end_nu = model.tardis_config.montecarlo.virtual_spectrum_range.start.to('Hz', units.spectral()).value - storage.spectrum_delta_nu = model.tardis_config.spectrum.frequency.value[1] - model.tardis_config.spectrum.frequency.value[0] + storage.spectrum_start_nu = runner.spectrum_frequency.value.min() + storage.spectrum_end_nu = runner.spectrum_frequency.value.max() + # TODO: Linspace handling for virtual_spectrum_range + storage.spectrum_virt_start_nu = runner.virtual_spectrum_range.end.to('Hz', units.spectral()).value + storage.spectrum_virt_end_nu = runner.virtual_spectrum_range.start.to('Hz', units.spectral()).value + storage.spectrum_delta_nu = runner.spectrum_frequency.value[1] - runner.spectrum_frequency.value[0] storage.spectrum_virt_nu = PyArray_DATA( runner.legacy_montecarlo_virtual_luminosity) - storage.sigma_thomson = ( - model.tardis_config.montecarlo.sigma_thomson.cgs.value) + storage.sigma_thomson = runner.sigma_thomson.cgs.value storage.inverse_sigma_thomson = 1.0 / storage.sigma_thomson - storage.reflective_inner_boundary = model.tardis_config.montecarlo.enable_reflective_inner_boundary - storage.inner_boundary_albedo = model.tardis_config.montecarlo.inner_boundary_albedo + storage.reflective_inner_boundary = runner.enable_reflective_inner_boundary + storage.inner_boundary_albedo = runner.inner_boundary_albedo # Data for continuum implementation cdef np.ndarray[double, ndim=1] t_electrons = plasma.t_electrons storage.t_electrons = t_electrons.data @@ -255,8 +253,7 @@ def montecarlo_radial1d(model, plasma, runner, int_type_t virtual_packet_flag=0, initialize_storage_model(model, plasma, runner, &storage) - montecarlo_main_loop(&storage, virtual_packet_flag, nthreads, - model.tardis_config.montecarlo.seed) + montecarlo_main_loop(&storage, virtual_packet_flag, nthreads, runner.seed) cdef np.ndarray[double, ndim=1] virt_packet_nus = np.zeros(storage.virt_packet_count, dtype=np.float64) cdef np.ndarray[double, ndim=1] virt_packet_energies = np.zeros(storage.virt_packet_count, dtype=np.float64) cdef np.ndarray[double, ndim=1] virt_packet_last_interaction_in_nu = np.zeros(storage.virt_packet_count, dtype=np.float64) @@ -298,8 +295,8 @@ def montecarlo_radial1d(model, plasma, runner, int_type_t virtual_packet_flag=0, def postprocess(model, runner): Edotlu_norm_factor = (1 / - (model.time_of_simulation * model.tardis_config.structure.volumes)) + (runner.time_of_simulation * model.volume)) exptau = 1 - np.exp(- runner.line_lists_tau_sobolevs.reshape(-1, runner.j_estimator.shape[0]) ) - model.Edotlu = Edotlu_norm_factor*exptau*runner.Edotlu_estimator + runner.Edotlu = Edotlu_norm_factor*exptau*runner.Edotlu_estimator From 894510cd9fe8288983e3b54a549a43fd92b7d400 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Fri, 5 Aug 2016 22:47:44 +0300 Subject: [PATCH 07/15] Add a from_config classmethod to the MontecarloRunner. --- tardis/montecarlo/base.py | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index ea85dd6bddd..13960f61b9d 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -366,3 +366,14 @@ def to_hdf(self, path_or_buf, path=''): 'spectrum_virtual') self.spectrum_reabsorbed.to_hdf(path_or_buf, runner_path, 'spectrum_reabsorbed') + + @classmethod + def from_config(cls, config): + return cls(seed=config.montecarlo.seed, + spectrum_frequency=config.spectrum.frequency, + virtual_spectrum_range=config.montecarlo.virtual_spectrum_range, + sigma_thomson=config.montecarlo.sigma_thomson, + enable_reflective_inner_boundary=config.montecarlo.enable_reflective_inner_boundary, + inner_boundary_albedo=config.montecarlo.inner_boundary_albedo, + line_interaction_type=config.plasma.line_interaction_type, + distance=config.supernova.get('distance', None)) \ No newline at end of file From f69c995df8e1827fcd5b46a4c8f91a4bf36c3f56 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Mon, 8 Aug 2016 13:30:08 +0300 Subject: [PATCH 08/15] montecarlo.pyx: Access atomic_data from the plasma, instead of model. --- tardis/montecarlo/montecarlo.pyx | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/tardis/montecarlo/montecarlo.pyx b/tardis/montecarlo/montecarlo.pyx index d58b0548554..682e3c8859f 100644 --- a/tardis/montecarlo/montecarlo.pyx +++ b/tardis/montecarlo/montecarlo.pyx @@ -144,8 +144,8 @@ cdef initialize_storage_model(model, plasma, runner, storage_model_t *storage): storage.l_pop_r = l_pop_r.data # Line lists - storage.no_of_lines = model.atom_data.lines.nu.values.size - storage.line_list_nu = PyArray_DATA(model.atom_data.lines.nu.values) + storage.no_of_lines = plasma.atomic_data.lines.nu.values.size + storage.line_list_nu = PyArray_DATA(plasma.atomic_data.lines.nu.values) runner.line_lists_tau_sobolevs = ( plasma.tau_sobolevs.values.flatten(order='F') ) @@ -172,17 +172,17 @@ cdef initialize_storage_model(model, plasma, runner, storage_model_t *storage): storage.transition_probabilities_nd = ( plasma.transition_probabilities.values.shape[0]) storage.line2macro_level_upper = PyArray_DATA( - model.atom_data.lines_upper2macro_reference_idx) + plasma.atomic_data.lines_upper2macro_reference_idx) storage.macro_block_references = PyArray_DATA( - model.atom_data.macro_atom_references['block_references'].values) + plasma.atomic_data.macro_atom_references['block_references'].values) storage.transition_type = PyArray_DATA( - model.atom_data.macro_atom_data['transition_type'].values) + plasma.atomic_data.macro_atom_data['transition_type'].values) # Destination level is not needed and/or generated for downbranch storage.destination_level_id = PyArray_DATA( - model.atom_data.macro_atom_data['destination_level_idx'].values) + plasma.atomic_data.macro_atom_data['destination_level_idx'].values) storage.transition_line_id = PyArray_DATA( - model.atom_data.macro_atom_data['lines_idx'].values) + plasma.atomic_data.macro_atom_data['lines_idx'].values) storage.output_nus = PyArray_DATA(runner._output_nu) storage.output_energies = PyArray_DATA(runner._output_energy) From 0e7cf7dbe97432d1ce28cf76aa7a742961050ff3 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Tue, 9 Aug 2016 19:24:53 +0300 Subject: [PATCH 09/15] Revert "Move and adapt calculate_j_blues from the old Radial1DModel to MontecarloRunner." This reverts commit 9f4385ec00740cc52be07d5c48b7bd7ba47f08a8, because it was decided that the j_blues will become plasma properties. --- tardis/montecarlo/base.py | 48 --------------------------------------- 1 file changed, 48 deletions(-) diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index 13960f61b9d..6be054e55ce 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -8,7 +8,6 @@ from spectrum import TARDISSpectrum from tardis.montecarlo import montecarlo, packet_source -from tardis.util import intensity_black_body from tardis.io.util import to_hdf import numpy as np @@ -294,53 +293,6 @@ def calculate_f_nu(self, frequency): def calculate_f_lambda(self, wavelength): pass - def calculate_j_blues(self, model, plasma, atom_data, - init_detailed_j_blues=False): - nus = atom_data.lines.nu.values - radiative_rates_type = plasma.radiative_rates_type - w_epsilon = plasma.w_epsilon - - if radiative_rates_type == 'blackbody': - logger.info('Calculating J_blues for radiative_rates_type=lte') - j_blues = intensity_black_body(nus[np.newaxis].T, model.t_rad.value) - return pd.DataFrame( - j_blues, index=atom_data.lines.index, - columns=np.arange(len(model.t_rad))) - - elif radiative_rates_type == 'dilute-blackbody' or init_detailed_j_blues: - logger.info('Calculating J_blues for radiative_rates_type=dilute-blackbody') - j_blues = model.w * intensity_black_body(nus[np.newaxis].T, model.t_rad.value) - return pd.DataFrame( - j_blues, index=atom_data.lines.index, - columns=np.arange(len(model.t_rad))) - - elif radiative_rates_type == 'detailed': - logger.info('Calculating J_blues for radiate_rates_type=detailed') - j_blues_norm_factor = (const.c.cgs * model.time_explosion / - (4 * np.pi * model.time_of_simulation * model.volume)) - j_blues = pd.DataFrame( - self.j_blue_estimator * - j_blues_norm_factor.value, - index=atom_data.lines.index, - columns=np.arange(len(model.t_rad))) - - for i in xrange(model.no_of_shells): - zero_j_blues = j_blues[i] == 0.0 - j_blues[i][zero_j_blues] = ( - w_epsilon * intensity_black_body( - atom_data.lines.nu[zero_j_blues].values, - model.t_rad.value[i])) - - else: - raise ValueError('radiative_rates_type type unknown - %s', radiative_rates_type) - - def calculate_time_of_simulation(self, model): - luminosity_inner = ( - 4 * np.pi * const.sigma_sb.cgs * model.r_inner[0] ** 2 * - model.t_inner ** 4).to('erg/s') - - return 1.0 * u.erg / luminosity_inner - def to_hdf(self, path_or_buf, path=''): """ Store the runner to an HDF structure. From f410a6073783bab145508435dc1364c8eb3c0329 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Tue, 9 Aug 2016 22:48:55 +0300 Subject: [PATCH 10/15] Add a calculate_time_of_simulation and a calculate_luminosity_inner to the MontecarloRunner. --- tardis/montecarlo/base.py | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index 6be054e55ce..2e312a2a9c9 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -287,6 +287,13 @@ def calculate_radiationfield_properties(self): return t_rad * u.K, w + def calculate_luminosity_inner(self, model): + return (4 * np.pi * const.sigma_sb.cgs * + model.r_inner[0] ** 2 * model.t_inner ** 4).to('erg/s') + + def calculate_time_of_simulation(self, model): + return (1.0 * u.erg / self.calculate_luminosity_inner(model)) + def calculate_f_nu(self, frequency): pass From 3895e14e20fa1d1ba28e4bc53f74773f2ca771d1 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Sun, 14 Aug 2016 00:13:26 +0300 Subject: [PATCH 11/15] Move the sigma_thomson calculation from config_reader to montecarlo. --- tardis/io/config_reader.py | 16 ---------------- tardis/montecarlo/base.py | 10 +++++++++- 2 files changed, 9 insertions(+), 17 deletions(-) diff --git a/tardis/io/config_reader.py b/tardis/io/config_reader.py index 6fc603d5e9f..960ff5ce5d3 100644 --- a/tardis/io/config_reader.py +++ b/tardis/io/config_reader.py @@ -943,22 +943,6 @@ def from_config_dict(cls, config_dict, atom_data=None, test_parser=False, else: plasma_section['t_rads'] = None - if plasma_section['disable_electron_scattering'] is False: - logger.debug("Electron scattering switched on") - validated_config_dict['montecarlo']['sigma_thomson'] = 6.652486e-25 / (u.cm ** 2) - else: - logger.warn('Disabling electron scattering - this is not physical') - validated_config_dict['montecarlo']['sigma_thomson'] = 1e-200 / (u.cm ** 2) - - if plasma_section['helium_treatment'] == 'recomb-nlte': - validated_config_dict['plasma']['helium_treatment'] == 'recomb-nlte' - else: - validated_config_dict['plasma']['helium_treatment'] == 'dilute-lte' - - - - - ##### NLTE subsection of Plasma start nlte_validated_config_dict = {} nlte_species = [] diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index 2e312a2a9c9..279eb36f721 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -328,10 +328,18 @@ def to_hdf(self, path_or_buf, path=''): @classmethod def from_config(cls, config): + if config.plasma.disable_electron_scattering: + logger.warn('Disabling electron scattering - this is not physical') + sigma_thomson = 1e-200 / (u.cm ** 2) + else: + logger.debug("Electron scattering switched on") + sigma_thomson = 6.652486e-25 / (u.cm ** 2) + + return cls(seed=config.montecarlo.seed, spectrum_frequency=config.spectrum.frequency, virtual_spectrum_range=config.montecarlo.virtual_spectrum_range, - sigma_thomson=config.montecarlo.sigma_thomson, + sigma_thomson=sigma_thomson, enable_reflective_inner_boundary=config.montecarlo.enable_reflective_inner_boundary, inner_boundary_albedo=config.montecarlo.inner_boundary_albedo, line_interaction_type=config.plasma.line_interaction_type, From 52071c8b2e163c7f8f14e968b7ec16ea072c7439 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Sun, 14 Aug 2016 00:22:14 +0300 Subject: [PATCH 12/15] Move the spectrum frequency calculation from config_reader to montecarlo. --- tardis/io/config_reader.py | 33 --------------------------- tardis/io/tests/test_config_reader.py | 8 ------- tardis/montecarlo/base.py | 8 +++++-- 3 files changed, 6 insertions(+), 43 deletions(-) diff --git a/tardis/io/config_reader.py b/tardis/io/config_reader.py index 960ff5ce5d3..c7f2aaf9fcf 100644 --- a/tardis/io/config_reader.py +++ b/tardis/io/config_reader.py @@ -480,35 +480,6 @@ def parse_supernova_section(supernova_dict): return config_dict -def parse_spectrum_list2dict(spectrum_list): - """ - Parse the spectrum list [start, stop, num] to a list - """ - if 'start' in spectrum_list and 'stop' in spectrum_list \ - and 'num' in spectrum_list: - spectrum_list = [spectrum_list['start'], spectrum_list['stop'], - spectrum_list['num']] - if spectrum_list[0].unit.physical_type != 'length' and \ - spectrum_list[1].unit.physical_type != 'length': - raise ValueError('start and end of spectrum need to be a length') - - - spectrum_config_dict = {} - spectrum_config_dict['start'] = spectrum_list[0] - spectrum_config_dict['end'] = spectrum_list[1] - spectrum_config_dict['bins'] = spectrum_list[2] - - spectrum_frequency = quantity_linspace( - spectrum_config_dict['end'].to('Hz', u.spectral()), - spectrum_config_dict['start'].to('Hz', u.spectral()), - num=spectrum_config_dict['bins'] + 1) - - spectrum_config_dict['frequency'] = spectrum_frequency - - return spectrum_config_dict - - - def parse_convergence_section(convergence_section_dict): """ Parse the convergence section dictionary @@ -1010,10 +981,6 @@ def from_config_dict(cls, config_dict, atom_data=None, test_parser=False, ###### END of convergence section reading - - validated_config_dict['spectrum'] = parse_spectrum_list2dict( - validated_config_dict['spectrum']) - return cls(validated_config_dict, atom_data) diff --git a/tardis/io/tests/test_config_reader.py b/tardis/io/tests/test_config_reader.py index 073f8b39a68..56aafe9ddbf 100644 --- a/tardis/io/tests/test_config_reader.py +++ b/tardis/io/tests/test_config_reader.py @@ -31,14 +31,6 @@ def test_quantity_linspace(): assert len(quantity_linspace) == 1000 -def test_spectrum_list2_dict(): - spectrum_dict = config_reader.parse_spectrum_list2dict( - [200*u.angstrom, 10000 * u.angstrom, 100]) - assert_almost_equal(spectrum_dict['start'].to(u.angstrom).value, 200) - assert_almost_equal(spectrum_dict['end'].to(u.angstrom).value, 10000) - assert_almost_equal(spectrum_dict['bins'], 100) - - def test_convergence_section_parser(): test_convergence_section = {'type': 'damped', 'lock_t_inner_cyles': 1, diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index 279eb36f721..feda90244f3 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -7,6 +7,7 @@ from scipy.special import zeta from spectrum import TARDISSpectrum +from tardis.util import quantity_linspace from tardis.montecarlo import montecarlo, packet_source from tardis.io.util import to_hdf @@ -335,9 +336,12 @@ def from_config(cls, config): logger.debug("Electron scattering switched on") sigma_thomson = 6.652486e-25 / (u.cm ** 2) - + spectrum_frequency = quantity_linspace( + config.spectrum.stop.to('Hz', u.spectral()), + config.spectrum.start.to('Hz', u.spectral()), + num=config.spectrum.num + 1) return cls(seed=config.montecarlo.seed, - spectrum_frequency=config.spectrum.frequency, + spectrum_frequency=spectrum_frequency, virtual_spectrum_range=config.montecarlo.virtual_spectrum_range, sigma_thomson=sigma_thomson, enable_reflective_inner_boundary=config.montecarlo.enable_reflective_inner_boundary, From 790737ca92ce178b3f54bf91233aee513ad18a7c Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Sun, 14 Aug 2016 00:50:26 +0300 Subject: [PATCH 13/15] Remove supernova distance checking from config_dict --- tardis/io/config_reader.py | 5 ----- 1 file changed, 5 deletions(-) diff --git a/tardis/io/config_reader.py b/tardis/io/config_reader.py index c7f2aaf9fcf..ca958c77e54 100644 --- a/tardis/io/config_reader.py +++ b/tardis/io/config_reader.py @@ -460,11 +460,6 @@ def parse_supernova_section(supernova_dict): config_dict['time_explosion'] = parse_quantity(supernova_dict['time_explosion']).to('s') - if 'distance' in supernova_dict: - config_dict['distance'] = parse_quantity(supernova_dict['distance']) - else: - config_dict['distance'] = None - if 'luminosity_wavelength_start' in supernova_dict: config_dict['luminosity_nu_end'] = parse_quantity(supernova_dict['luminosity_wavelength_start']). \ to('Hz', u.spectral()) From 931751efa5bc0462822767a61fe7378c6d0936bc Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Sun, 14 Aug 2016 01:06:23 +0300 Subject: [PATCH 14/15] Don't rename virtual_spectrum_range's stop and num to end and samples. --- tardis/io/config_reader.py | 8 -------- tardis/montecarlo/base.py | 2 ++ tardis/montecarlo/montecarlo.pyx | 2 +- 3 files changed, 3 insertions(+), 9 deletions(-) diff --git a/tardis/io/config_reader.py b/tardis/io/config_reader.py index ca958c77e54..77a72d6c32d 100644 --- a/tardis/io/config_reader.py +++ b/tardis/io/config_reader.py @@ -965,14 +965,6 @@ def from_config_dict(cls, config_dict, atom_data=None, test_parser=False, black_body_section['stop'] montecarlo_section['black_body_sampling']['samples'] = \ black_body_section['num'] - virtual_spectrum_section = montecarlo_section['virtual_spectrum_range'] - montecarlo_section['virtual_spectrum_range'] = {} - montecarlo_section['virtual_spectrum_range']['start'] = \ - virtual_spectrum_section['start'] - montecarlo_section['virtual_spectrum_range']['end'] = \ - virtual_spectrum_section['stop'] - montecarlo_section['virtual_spectrum_range']['samples'] = \ - virtual_spectrum_section['num'] ###### END of convergence section reading diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index feda90244f3..00fbcd6867e 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -340,8 +340,10 @@ def from_config(cls, config): config.spectrum.stop.to('Hz', u.spectral()), config.spectrum.start.to('Hz', u.spectral()), num=config.spectrum.num + 1) + return cls(seed=config.montecarlo.seed, spectrum_frequency=spectrum_frequency, + # TODO: Linspace handling for virtual_spectrum_range virtual_spectrum_range=config.montecarlo.virtual_spectrum_range, sigma_thomson=sigma_thomson, enable_reflective_inner_boundary=config.montecarlo.enable_reflective_inner_boundary, diff --git a/tardis/montecarlo/montecarlo.pyx b/tardis/montecarlo/montecarlo.pyx index 682e3c8859f..82f060b6af5 100644 --- a/tardis/montecarlo/montecarlo.pyx +++ b/tardis/montecarlo/montecarlo.pyx @@ -204,7 +204,7 @@ cdef initialize_storage_model(model, plasma, runner, storage_model_t *storage): storage.spectrum_start_nu = runner.spectrum_frequency.value.min() storage.spectrum_end_nu = runner.spectrum_frequency.value.max() # TODO: Linspace handling for virtual_spectrum_range - storage.spectrum_virt_start_nu = runner.virtual_spectrum_range.end.to('Hz', units.spectral()).value + storage.spectrum_virt_start_nu = runner.virtual_spectrum_range.stop.to('Hz', units.spectral()).value storage.spectrum_virt_end_nu = runner.virtual_spectrum_range.start.to('Hz', units.spectral()).value storage.spectrum_delta_nu = runner.spectrum_frequency.value[1] - runner.spectrum_frequency.value[0] From dde78a86d4b97e65a1faa4efad608ce65133c010 Mon Sep 17 00:00:00 2001 From: Fotis Tsamis Date: Tue, 16 Aug 2016 22:33:42 +0300 Subject: [PATCH 15/15] Fix a bug when initializing the logger for montecarlo/base.py --- tardis/montecarlo/base.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tardis/montecarlo/base.py b/tardis/montecarlo/base.py index 00fbcd6867e..c870834672c 100644 --- a/tardis/montecarlo/base.py +++ b/tardis/montecarlo/base.py @@ -14,7 +14,7 @@ import numpy as np import pandas as pd -logger = logging.getLevelName(__name__) +logger = logging.getLogger(__name__) class MontecarloRunner(object): """