Skip to content

Commit

Permalink
Initial application of the configuration object (incomplete)
Browse files Browse the repository at this point in the history
  • Loading branch information
andrewfullard committed Feb 28, 2024
1 parent 434cef5 commit 8acdcf0
Show file tree
Hide file tree
Showing 23 changed files with 174 additions and 120 deletions.
1 change: 1 addition & 0 deletions tardis/io/model_reader.py
Original file line number Diff line number Diff line change
Expand Up @@ -741,6 +741,7 @@ def transport_from_hdf(fname):
nthreads=d["nthreads"],
enable_virtual_packet_logging=d["virt_logging"],
use_gpu=d["use_gpu"],
montecarlo_configuration=d["montecarlo_configuration"],
)

new_transport.Edotlu_estimator = d["Edotlu_estimator"]
Expand Down
12 changes: 7 additions & 5 deletions tardis/montecarlo/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,6 @@
initialize_estimator_statistics,
)
from tardis.montecarlo.montecarlo_configuration import (
configuration_initialize,
MonteCarloConfiguration,
)
from tardis.montecarlo.montecarlo_numba import (
Expand Down Expand Up @@ -126,7 +125,10 @@ def initialize_transport_state(

geometry_state = simulation_state.geometry.to_numba()
opacity_state = opacity_state_initialize(
plasma, self.line_interaction_type
plasma,
self.line_interaction_type,
self.montecarlo_configuration.DISABLE_LINE_SCATTERING,
self.montecarlo_configuration.CONTINUUM_PROCESSES_ENABLED,
)
transport_state = MonteCarloTransportState(
packet_collection,
Expand All @@ -141,8 +143,8 @@ def initialize_transport_state(
transport_state._integrator = FormalIntegrator(
simulation_state, plasma, self
)
configuration_initialize(
self.montecarlo_configuration, self, no_of_virtual_packets
self.montecarlo_configuration.configuration_initialize(
transport=self, number_of_vpackets=no_of_virtual_packets
)

return transport_state
Expand Down Expand Up @@ -305,7 +307,7 @@ def from_config(
valid values are 'GPU', 'CPU', and 'Automatic'."""
)

montecarlo_configuration = MonteCarloConfiguration
montecarlo_configuration = MonteCarloConfiguration()

montecarlo_configuration.DISABLE_LINE_SCATTERING = (
config.plasma.disable_line_scattering
Expand Down
6 changes: 5 additions & 1 deletion tardis/montecarlo/montecarlo_numba/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -88,7 +88,11 @@ def montecarlo_main_loop(
montecarlo_configuration.TEMPORARY_V_PACKET_BINS,
)
)
rpacket_trackers.append(RPacketTracker())
rpacket_trackers.append(
RPacketTracker(
montecarlo_configuration.INITIAL_TRACKING_ARRAY_LENGTH
)
)

# Get the ID of the main thread and the number of threads
main_thread_id = get_thread_id()
Expand Down
11 changes: 6 additions & 5 deletions tardis/montecarlo/montecarlo_numba/formal_integral.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,9 +11,6 @@


from tardis.montecarlo.montecarlo_numba.numba_config import SIGMA_THOMSON
from tardis.montecarlo import (
montecarlo_configuration as montecarlo_configuration,
)
from tardis.montecarlo.montecarlo_numba import njit_dict, njit_dict_no_parallel
from tardis.montecarlo.montecarlo_numba.numba_interface import (
opacity_state_initialize,
Expand Down Expand Up @@ -284,9 +281,13 @@ def __init__(self, simulation_state, plasma, transport, points=1000):
self.simulation_state = simulation_state
self.transport = transport
self.points = points
self.montecarlo_configuration = self.transport.montecarlo_configuration
if plasma:
self.plasma = opacity_state_initialize(
plasma, transport.line_interaction_type
plasma,
transport.line_interaction_type,
self.montecarlo_configuration.DISABLE_LINE_SCATTERING,
self.montecarlo_configuration.CONTINUUM_PROCESSES_ENABLED,
)
self.atomic_data = plasma.atomic_data
self.original_plasma = plasma
Expand Down Expand Up @@ -361,7 +362,7 @@ def raise_or_return(message):
'and line_interaction_type == "macroatom"'
)

if montecarlo_configuration.CONTINUUM_PROCESSES_ENABLED:
if self.montecarlo_configuration.CONTINUUM_PROCESSES_ENABLED:
return raise_or_return(
"The FormalIntegrator currently does not work for "
"continuum interactions."
Expand Down
99 changes: 71 additions & 28 deletions tardis/montecarlo/montecarlo_numba/interaction.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,9 +2,6 @@
from numba import njit

from tardis import constants as const
from tardis.montecarlo import (
montecarlo_configuration as montecarlo_configuration,
)
from tardis.montecarlo.montecarlo_numba import njit_dict_no_parallel
from tardis.montecarlo.montecarlo_numba.macro_atom import (
MacroAtomTransitionType,
Expand Down Expand Up @@ -150,6 +147,8 @@ def continuum_event(
chi_ff,
chi_bf_contributions,
current_continua,
continuum_processes_enabled,
enable_full_relativity,
):
"""
continuum event handler - activate the macroatom and run the handler
Expand All @@ -162,12 +161,12 @@ def continuum_event(
continuum : tardis.montecarlo.montecarlo_numba.numba_interface.Continuum
"""
old_doppler_factor = get_doppler_factor(
r_packet.r, r_packet.mu, time_explosion
r_packet.r, r_packet.mu, time_explosion, enable_full_relativity
)

r_packet.mu = get_random_mu()
inverse_doppler_factor = get_inverse_doppler_factor(
r_packet.r, r_packet.mu, time_explosion
r_packet.r, r_packet.mu, time_explosion, enable_full_relativity
)
comov_energy = r_packet.energy * old_doppler_factor
comov_nu = (
Expand All @@ -186,13 +185,23 @@ def continuum_event(
)

macro_atom_event(
destination_level_idx, r_packet, time_explosion, opacity_state
destination_level_idx,
r_packet,
time_explosion,
opacity_state,
continuum_processes_enabled,
enable_full_relativity,
)


@njit(**njit_dict_no_parallel)
def macro_atom_event(
destination_level_idx, r_packet, time_explosion, opacity_state
destination_level_idx,
r_packet,
time_explosion,
opacity_state,
continuum_processes_enabled,
enable_full_relativity,
):
"""
Macroatom event handler - run the macroatom and handle the result
Expand All @@ -209,32 +218,45 @@ def macro_atom_event(
)

if (
montecarlo_configuration.CONTINUUM_PROCESSES_ENABLED
continuum_processes_enabled
and transition_type == MacroAtomTransitionType.FF_EMISSION
):
free_free_emission(r_packet, time_explosion, opacity_state)
free_free_emission(
r_packet, time_explosion, opacity_state, enable_full_relativity
)

elif (
montecarlo_configuration.CONTINUUM_PROCESSES_ENABLED
continuum_processes_enabled
and transition_type == MacroAtomTransitionType.BF_EMISSION
):
bound_free_emission(
r_packet, time_explosion, opacity_state, transition_id
r_packet,
time_explosion,
opacity_state,
transition_id,
enable_full_relativity,
)
elif (
montecarlo_configuration.CONTINUUM_PROCESSES_ENABLED
continuum_processes_enabled
and transition_type == MacroAtomTransitionType.BF_COOLING
):
bf_cooling(r_packet, time_explosion, opacity_state)

elif (
montecarlo_configuration.CONTINUUM_PROCESSES_ENABLED
continuum_processes_enabled
and transition_type == MacroAtomTransitionType.ADIABATIC_COOLING
):
adiabatic_cooling(r_packet)

elif transition_type == MacroAtomTransitionType.BB_EMISSION:
line_emission(r_packet, transition_id, time_explosion, opacity_state)
line_emission(
r_packet,
transition_id,
time_explosion,
opacity_state,
continuum_processes_enabled,
enable_full_relativity,
)
else:
raise Exception("No Interaction Found!")

Expand Down Expand Up @@ -295,7 +317,9 @@ def get_current_line_id(nu, line_list):


@njit(**njit_dict_no_parallel)
def free_free_emission(r_packet, time_explosion, opacity_state):
def free_free_emission(
r_packet, time_explosion, opacity_state, full_relativity
):
"""
Free-Free emission - set the frequency from electron-ion interaction
Expand All @@ -313,14 +337,16 @@ def free_free_emission(r_packet, time_explosion, opacity_state):
current_line_id = get_current_line_id(comov_nu, opacity_state.line_list_nu)
r_packet.next_line_id = current_line_id

if montecarlo_configuration.ENABLE_FULL_RELATIVITY:
if full_relativity:
r_packet.mu = angle_aberration_CMF_to_LF(
r_packet, time_explosion, r_packet.mu
)


@njit(**njit_dict_no_parallel)
def bound_free_emission(r_packet, time_explosion, opacity_state, continuum_id):
def bound_free_emission(
r_packet, time_explosion, opacity_state, continuum_id, full_relativity
):
"""
Bound-Free emission - set the frequency from photo-ionization
Expand All @@ -342,14 +368,14 @@ def bound_free_emission(r_packet, time_explosion, opacity_state, continuum_id):
current_line_id = get_current_line_id(comov_nu, opacity_state.line_list_nu)
r_packet.next_line_id = current_line_id

if montecarlo_configuration.ENABLE_FULL_RELATIVITY:
if full_relativity:
r_packet.mu = angle_aberration_CMF_to_LF(
r_packet, time_explosion, r_packet.mu
)


@njit(**njit_dict_no_parallel)
def thomson_scatter(r_packet, time_explosion):
def thomson_scatter(r_packet, time_explosion, enable_full_relativity):
"""
Thomson scattering — no longer line scattering
\n1) get the doppler factor at that position with the old angle
Expand All @@ -364,7 +390,7 @@ def thomson_scatter(r_packet, time_explosion):
time since explosion in seconds
"""
old_doppler_factor = get_doppler_factor(
r_packet.r, r_packet.mu, time_explosion
r_packet.r, r_packet.mu, time_explosion, enable_full_relativity
)
comov_nu = r_packet.nu * old_doppler_factor
comov_energy = r_packet.energy * old_doppler_factor
Expand All @@ -375,18 +401,23 @@ def thomson_scatter(r_packet, time_explosion):

r_packet.nu = comov_nu * inverse_new_doppler_factor
r_packet.energy = comov_energy * inverse_new_doppler_factor
if montecarlo_configuration.ENABLE_FULL_RELATIVITY:
if enable_full_relativity:
r_packet.mu = angle_aberration_CMF_to_LF(
r_packet, time_explosion, r_packet.mu
)
temp_doppler_factor = get_doppler_factor(
r_packet.r, r_packet.mu, time_explosion
r_packet.r, r_packet.mu, time_explosion, enable_full_relativity
)


@njit(**njit_dict_no_parallel)
def line_scatter(
r_packet, time_explosion, line_interaction_type, opacity_state
r_packet,
time_explosion,
line_interaction_type,
opacity_state,
continuum_processes_enabled,
enable_full_relativity,
):
"""
Line scatter function that handles the scattering itself, including new angle drawn, and calculating nu out using macro atom
Expand All @@ -399,7 +430,7 @@ def line_scatter(
opacity_state : tardis.montecarlo.montecarlo_numba.numba_interface.OpacityState
"""
old_doppler_factor = get_doppler_factor(
r_packet.r, r_packet.mu, time_explosion
r_packet.r, r_packet.mu, time_explosion, enable_full_relativity
)
r_packet.mu = get_random_mu()

Expand All @@ -412,7 +443,12 @@ def line_scatter(

if line_interaction_type == LineInteractionType.SCATTER:
line_emission(
r_packet, r_packet.next_line_id, time_explosion, opacity_state
r_packet,
r_packet.next_line_id,
time_explosion,
opacity_state,
continuum_processes_enabled,
enable_full_relativity,
)
else: # includes both macro atom and downbranch - encoded in the transition probabilities
comov_nu = r_packet.nu * old_doppler_factor # Is this necessary?
Expand All @@ -421,12 +457,19 @@ def line_scatter(
r_packet.next_line_id
]
macro_atom_event(
activation_level_id, r_packet, time_explosion, opacity_state
activation_level_id,
r_packet,
time_explosion,
opacity_state,
continuum_processes_enabled,
enable_full_relativity,
)


@njit(**njit_dict_no_parallel)
def line_emission(r_packet, emission_line_id, time_explosion, opacity_state):
def line_emission(
r_packet, emission_line_id, time_explosion, opacity_state, full_relativity
):
"""
Sets the frequency of the RPacket properly given the emission channel
Expand All @@ -451,7 +494,7 @@ def line_emission(r_packet, emission_line_id, time_explosion, opacity_state):
r_packet.next_line_id = emission_line_id + 1
nu_line = opacity_state.line_list_nu[emission_line_id]

if montecarlo_configuration.ENABLE_FULL_RELATIVITY:
if full_relativity:
r_packet.mu = angle_aberration_CMF_to_LF(
r_packet, time_explosion, r_packet.mu
)
Loading

0 comments on commit 8acdcf0

Please sign in to comment.