From 0a0c9a158ea3e0bd7c84e1461fdc15fd29141a10 Mon Sep 17 00:00:00 2001 From: karandesai-96 Date: Wed, 30 Mar 2016 08:28:00 +0530 Subject: [PATCH] Move structs to struct.py, init default values. - __init__ in RPacket and StorageModel resemble the same init methods in test_cmontecarlo.c - The way arguments are passed in paramterization is changed from individual params to a dict for better aesthetic sense. --- tardis/montecarlo/struct.py | 190 ++++++++++++++++++++ tardis/montecarlo/tests/test_cmontecarlo.py | 163 +++-------------- 2 files changed, 215 insertions(+), 138 deletions(-) create mode 100644 tardis/montecarlo/struct.py diff --git a/tardis/montecarlo/struct.py b/tardis/montecarlo/struct.py new file mode 100644 index 00000000000..87a1f976c78 --- /dev/null +++ b/tardis/montecarlo/struct.py @@ -0,0 +1,190 @@ +import os +from ctypes import * + +import numpy as np +from tardis import __path__ as path + +test_path = os.path.join(path[0], 'montecarlo', 'test_montecarlo.so') +tests = CDLL(test_path) + +cmontecarlo_filepath = os.path.join(path[0], 'montecarlo', 'montecarlo.so') +cmontecarlo_methods = CDLL(cmontecarlo_filepath) + + +class RPacket(Structure): + _fields_ = [ + ('nu', c_double), + ('mu', c_double), + ('energy', c_double), + ('r', c_double), + ('tau_event', c_double), + ('nu_line', c_double), + ('current_shell_id', c_int64), + ('next_line_id', c_int64), + ('last_line', c_int64), + ('close_line', c_int64), + ('recently_crossed_boundary', c_int64), + ('current_continuum_id', c_int64), + ('virtual_packet_flag', c_int64), + ('virtual_packet', c_int64), + ('d_line', c_double), + ('d_electron', c_double), + ('d_boundary', c_double), + ('d_cont', c_double), + ('next_shell_id', c_int64), + ('status', c_int), + ('id', c_int64), + ('chi_th', c_double), + ('chi_cont', c_double), + ('chi_ff', c_double), + ('chi_bf', c_double) + ] + + def __init__(self, **kwargs): + super(RPacket, self).__init__() + self.nu = 0.4 + self.mu = 0.3 + self.energy = 0.9 + self.r = 7.5e14 + self.tau_event = 2.9e13 + self.nu_line = 0.2 + self.current_shell_id = 0 + self.next_line_id = 1 + self.last_line = 0 + self.close_line = 0 + self.recently_crossed_boundary = 1 + self.current_continuum_id = 1 + self.virtual_packet_flag = 1 + self.virtual_packet = 0 + self.status = 0 + self.id = 0 + + for key, value in kwargs.iteritems(): + setattr(self, key, value) + + +class StorageModel(Structure): + _fields_ = [ + ('packet_nus', POINTER(c_double)), + ('packet_mus', POINTER(c_double)), + ('packet_energies', POINTER(c_double)), + ('output_nus', POINTER(c_double)), + ('output_energies', POINTER(c_double)), + ('last_interaction_in_nu', POINTER(c_double)), + ('last_line_interaction_in_id', POINTER(c_int64)), + ('last_line_interaction_out_id', POINTER(c_int64)), + ('last_line_interaction_shell_id', POINTER(c_int64)), + ('last_line_interaction_type', POINTER(c_int64)), + ('no_of_packets', c_int64), + ('no_of_shells', c_int64), + ('r_inner', POINTER(c_double)), + ('r_outer', POINTER(c_double)), + ('v_inner', POINTER(c_double)), + ('time_explosion', c_double), + ('inverse_time_explosion', c_double), + ('electron_densities', POINTER(c_double)), + ('inverse_electron_densities', POINTER(c_double)), + ('line_list_nu', POINTER(c_double)), + ('continuum_list_nu', POINTER(c_double)), + ('line_lists_tau_sobolevs', POINTER(c_double)), + ('line_lists_tau_sobolevs_nd', c_int64), + ('line_lists_j_blues', POINTER(c_double)), + ('line_lists_j_blues_nd', c_int64), + ('no_of_lines', c_int64), + ('no_of_edges', c_int64), + ('line_interaction_id', c_int64), + ('transition_probabilities', POINTER(c_double)), + ('transition_probabilities_nd', c_int64), + ('line2macro_level_upper', POINTER(c_int64)), + ('macro_block_references', POINTER(c_int64)), + ('transition_type', POINTER(c_int64)), + ('destination_level_id', POINTER(c_int64)), + ('transition_line_id', POINTER(c_int64)), + ('js', POINTER(c_double)), + ('nubars', POINTER(c_double)), + ('spectrum_start_nu', c_double), + ('spectrum_delta_nu', c_double), + ('spectrum_end_nu', c_double), + ('spectrum_virt_start_nu', c_double), + ('spectrum_virt_end_nu', c_double), + ('spectrum_virt_nu', POINTER(c_double)), + ('sigma_thomson', c_double), + ('inverse_sigma_thomson', c_double), + ('inner_boundary_albedo', c_double), + ('reflective_inner_boundary', c_int64), + ('current_packet_id', c_int64), + ('chi_bf_tmp_partial', POINTER(c_double)), + ('t_electrons', POINTER(c_double)), + ('l_pop', POINTER(c_double)), + ('l_pop_r', POINTER(c_double)), + ('cont_status', c_int), + ('virt_packet_nus', POINTER(c_double)), + ('virt_packet_energies', POINTER(c_double)), + ('virt_packet_last_interaction_in_nu', POINTER(c_double)), + ('virt_packet_last_interaction_type', POINTER(c_int64)), + ('virt_packet_last_line_interaction_in_id', POINTER(c_int64)), + ('virt_packet_last_line_interaction_out_id', POINTER(c_int64)), + ('virt_packet_count', c_int64), + ('virt_array_size', c_int64) + ] + + def __init__(self, **kwargs): + super(StorageModel, self).__init__() + + self.last_line_interaction_in_id = np.ctypeslib.as_ctypes(np.array([0] * 2)) + self.last_line_interaction_shell_id = np.ctypeslib.as_ctypes(np.array([0] * 2)) + self.last_line_interaction_type = np.ctypeslib.as_ctypes(np.array([2])) + + self.no_of_shells = 2 + + self.r_inner = np.ctypeslib.as_ctypes(np.array([6.912e14, 8.64e14])) + self.r_outer = np.ctypeslib.as_ctypes(np.array([8.64e14, 1.0368e15])) + + self.time_explosion = 5.2e7 + self.inverse_time_explosion = 1 / 5.2e7 + + # self.electron_densities = np.ctypeslib.as_ctypes(np.array([1.0e9] * 2)) + # self.inverse_electron_densities = np.ctypeslib.as_ctypes(np.array([1.0e-9] * 2)) + + self.line_list_nu = np.ctypeslib.as_ctypes(np.array([1.26318289e+16, 1.26318289e+16, + 1.23357675e+16, 1.23357675e+16, + 1.16961598e+16])) + + self.continuum_list_nu = np.ctypeslib.as_ctypes(np.array([1.e13] * 20000)) + + self.line_lists_tau_sobolevs = np.ctypeslib.as_ctypes(np.array([1.e-5] * 1000)) + # self.line_lists_j_blues = np.ctypeslib.as_ctypes(np.array([1.e-10] * 2)) + self.line_lists_j_blues_nd = 0 + + self.no_of_lines = 2 + self.no_of_edges = 100 + + self.line_interaction_id = 0 + # self.line2macro_level_upper = np.ctypeslib.as_ctypes(np.array([0] * 2)) + + # self.js = np.ctypeslib.as_ctypes(np.array([0.0] * 2)) + # self.nubars = np.ctypeslib.as_ctypes(np.array([0.0] * 2)) + + self.spectrum_start_nu = 1.e14 + self.spectrum_delta_nu = 293796608840.0 + self.spectrum_end_nu = 6.e15 + + self.spectrum_virt_start_nu = 1e14 + self.spectrum_virt_end_nu = 6e15 + self.spectrum_virt_nu = np.ctypeslib.as_ctypes(np.array([0.0] * 20000)) + + self.sigma_thomson = 6.652486e-25 + self.inverse_sigma_thomson = 1 / self.sigma_thomson + + self.inner_boundary_albedo = 0.0 + self.reflective_inner_boundary = 0 + + self.chi_bf_tmp_partial = np.ctypeslib.as_ctypes(np.array([160.0] * 20000)) + # self.t_electrons = np.ctypeslib.as_ctypes(np.array([0.0] * 2)) + + self.l_pop = np.ctypeslib.as_ctypes(np.array([2.0] * 20000)) + self.l_pop_r = np.ctypeslib.as_ctypes(np.array([3.0] * 20000)) + + for key, value in kwargs.iteritems(): + setattr(self, key, value) + diff --git a/tardis/montecarlo/tests/test_cmontecarlo.py b/tardis/montecarlo/tests/test_cmontecarlo.py index 38b736eb4a1..e9f199b9af5 100644 --- a/tardis/montecarlo/tests/test_cmontecarlo.py +++ b/tardis/montecarlo/tests/test_cmontecarlo.py @@ -1,163 +1,50 @@ -import os -import random -from ctypes import * - import pytest -import numpy as np -import numpy.ctypeslib as npclib -import numpy.testing as nptesting -from tardis import __path__ as path - -test_path = os.path.join(path[0], 'montecarlo', 'test_montecarlo.so') -tests = CDLL(test_path) - -cmontecarlo_filepath = os.path.join(path[0], 'montecarlo', 'montecarlo.so') -cmontecarlo_methods = CDLL(cmontecarlo_filepath) - - -class RPacket(Structure): - _fields_ = [ - ('nu', c_double), - ('mu', c_double), - ('energy', c_double), - ('r', c_double), - ('tau_event', c_double), - ('nu_line', c_double), - ('current_shell_id', c_int64), - ('next_line_id', c_int64), - ('last_line', c_int64), - ('close_line', c_int64), - ('recently_crossed_boundary', c_int64), - ('current_continuum_id', c_int64), - ('virtual_packet_flag', c_int64), - ('virtual_packet', c_int64), - ('d_line', c_double), - ('d_electron', c_double), - ('d_boundary', c_double), - ('d_cont', c_double), - ('next_shell_id', c_int64), - ('status', c_int), - ('id', c_int64), - ('chi_th', c_double), - ('chi_cont', c_double), - ('chi_ff', c_double), - ('chi_bf', c_double) - ] - - def __init__(self, **kwargs): - super(RPacket, self).__init__() - for key, value in kwargs.iteritems(): - setattr(self, key, value) - - -class StorageModel(Structure): - _fields_ = [ - ('packet_nus', POINTER(c_double)), - ('packet_mus', POINTER(c_double)), - ('packet_energies', POINTER(c_double)), - ('output_nus', POINTER(c_double)), - ('output_energies', POINTER(c_double)), - ('last_interaction_in_nu', POINTER(c_double)), - ('last_line_interaction_in_id', POINTER(c_int64)), - ('last_line_interaction_out_id', POINTER(c_int64)), - ('last_line_interaction_shell_id', POINTER(c_int64)), - ('last_line_interaction_type', POINTER(c_int64)), - ('no_of_packets', c_int64), - ('no_of_shells', c_int64), - ('r_inner', POINTER(c_double)), - ('r_outer', POINTER(c_double)), - ('v_inner', POINTER(c_double)), - ('time_explosion', c_double), - ('inverse_time_explosion', c_double), - ('electron_densities', POINTER(c_double)), - ('inverse_electron_densities', POINTER(c_double)), - ('line_list_nu', POINTER(c_double)), - ('continuum_list_nu', POINTER(c_double)), - ('line_lists_tau_sobolevs', POINTER(c_double)), - ('line_lists_tau_sobolevs_nd', c_int64), - ('line_lists_j_blues', POINTER(c_double)), - ('line_lists_j_blues_nd', c_int64), - ('no_of_lines', c_int64), - ('no_of_edges', c_int64), - ('line_interaction_id', c_int64), - ('transition_probabilities', POINTER(c_double)), - ('transition_probabilities_nd', c_int64), - ('line2macro_level_upper', POINTER(c_int64)), - ('macro_block_references', POINTER(c_int64)), - ('transition_type', POINTER(c_int64)), - ('destination_level_id', POINTER(c_int64)), - ('transition_line_id', POINTER(c_int64)), - ('js', POINTER(c_double)), - ('nubars', POINTER(c_double)), - ('spectrum_start_nu', c_double), - ('spectrum_delta_nu', c_double), - ('spectrum_end_nu', c_double), - ('spectrum_virt_start_nu', c_double), - ('spectrum_virt_end_nu', c_double), - ('spectrum_virt_nu', POINTER(c_double)), - ('sigma_thomson', c_double), - ('inverse_sigma_thomson', c_double), - ('inner_boundary_albedo', c_double), - ('reflective_inner_boundary', c_int64), - ('current_packet_id', c_int64), - ('chi_bf_tmp_partial', POINTER(c_double)), - ('t_electrons', POINTER(c_double)), - ('l_pop', POINTER(c_double)), - ('l_pop_r', POINTER(c_double)), - ('cont_status', c_int), - ('virt_packet_nus', POINTER(c_double)), - ('virt_packet_energies', POINTER(c_double)), - ('virt_packet_last_interaction_in_nu', POINTER(c_double)), - ('virt_packet_last_interaction_type', POINTER(c_int64)), - ('virt_packet_last_line_interaction_in_id', POINTER(c_int64)), - ('virt_packet_last_line_interaction_out_id', POINTER(c_int64)), - ('virt_packet_count', c_int64), - ('virt_array_size', c_int64) - ] - - def __init__(self, **kwargs): - super(StorageModel, self).__init__() - for key, value in kwargs.iteritems(): - setattr(self, key, value) +from tardis.montecarlo.struct import * @pytest.mark.parametrize( - ['mu', 'r', 'recently_crossed_boundary', 'r_inner', 'r_outer', 'expected'], - [(0.3, 7.5e14, 1, [6.912e14, 8.64e14], [8.64e14, 1.0368e15], 259376919351035.88), - (0.3, 7.5e14, 0, [6.912e14, 8.64e14], [8.64e14, 1.0368e15], 259376919351035.88), - (-0.3, 7.5e13, 0, [6.912e14, 8.64e14], [8.64e14, 1.0368e15], -838532664885601.1), - (-0.3, 7.5e14, 0, [6.912e14, 8.64e14], [8.64e14, 1.0368e15], -259376919351035.88)] + ['packet_params', 'expected'], + [({'mu': 0.3, + 'r': 7.5e14, + 'recently_crossed_boundary': 1}, 259376919351035.88), + ({'mu': 0.3, + 'r': 7.5e14, + 'recently_crossed_boundary': 0}, 259376919351035.88), + ({'mu': -0.3, + 'r': 7.5e13, + 'recently_crossed_boundary': 0}, -838532664885601.1), + ({'mu': -0.3, + 'r': 7.5e14, + 'recently_crossed_boundary': 0}, -259376919351035.88)] ) -def test_compute_distance2boundary(mu, r, recently_crossed_boundary, r_inner, r_outer, expected): - packet = RPacket(mu=mu, r=r, current_shell_id=0, next_shell_id=1, - recently_crossed_boundary=recently_crossed_boundary) - model = StorageModel(r_inner=npclib.as_ctypes(np.array(r_inner)), - r_outer=npclib.as_ctypes(np.array(r_outer))) +def test_compute_distance2boundary(packet_params, expected): + packet = RPacket(**packet_params) + model = StorageModel() cmontecarlo_methods.compute_distance2boundary.restype = c_double obtained = cmontecarlo_methods.compute_distance2boundary(byref(packet), byref(model)) - nptesting.assert_almost_equal(obtained, expected) + np.testing.assert_almost_equal(obtained, expected) def test_compute_distance2line(): distance_to_line = 7.792353908000001e+17 tests.test_compute_distance2line.restype = c_double - nptesting.assert_almost_equal(tests.test_compute_distance2line(), + np.testing.assert_almost_equal(tests.test_compute_distance2line(), distance_to_line) def test_compute_distance2continuum(): distance_to_electron = 4.359272608766106e+28 tests.test_compute_distance2continuum.restype = c_double - nptesting.assert_almost_equal(tests.test_compute_distance2continuum(), + np.testing.assert_almost_equal(tests.test_compute_distance2continuum(), distance_to_electron) def test_rpacket_doppler_factor(): doppler_factor = 0.9998556693818854 tests.test_rpacket_doppler_factor.restype = c_double - nptesting.assert_almost_equal(tests.test_rpacket_doppler_factor(), + np.testing.assert_almost_equal(tests.test_rpacket_doppler_factor(), doppler_factor) @@ -165,14 +52,14 @@ def test_rpacket_doppler_factor(): def test_move_packet(): doppler_factor = 0.9998556693818854 tests.test_move_packet.restype = c_double - nptesting.assert_almost_equal(tests.test_move_packet(), + np.testing.assert_almost_equal(tests.test_move_packet(), doppler_factor) def test_increment_j_blue_estimator(): j_blue = 1.1249855669381885 tests.test_increment_j_blue_estimator.restype = c_double - nptesting.assert_almost_equal(tests.test_increment_j_blue_estimator(), + np.testing.assert_almost_equal(tests.test_increment_j_blue_estimator(), j_blue) @@ -199,7 +86,7 @@ def test_montecarlo_thomson_scatter(): def test_calculate_chi_bf(): chi_bf = 1.0006697327643788 tests.test_calculate_chi_bf.restype = c_double - nptesting.assert_almost_equal(tests.test_calculate_chi_bf(), + np.testing.assert_almost_equal(tests.test_calculate_chi_bf(), chi_bf) @@ -212,7 +99,7 @@ def test_montecarlo_bound_free_scatter(): def test_bf_cross_section(): bf_cross_section = 0.0 tests.test_bf_cross_section.restype = c_double - nptesting.assert_almost_equal(tests.test_bf_cross_section(), + np.testing.assert_almost_equal(tests.test_bf_cross_section(), bf_cross_section)