Skip to content

Commit

Permalink
Browse files Browse the repository at this point in the history
  • Loading branch information
1e9abhi1e10 committed Jan 4, 2023
1 parent 3dfff4b commit 97047ce
Show file tree
Hide file tree
Showing 11 changed files with 41 additions and 43 deletions.
1 change: 1 addition & 0 deletions .mailmap
Original file line number Diff line number Diff line change
@@ -1,3 +1,4 @@
Abhishek Patidar <1e9abhi1e10@gmail.com>
Adam Suban-Loewen <adam.subanloewen@gmail.com>

Alexander Holas <alexander.holas@h-its.org>
Expand Down
6 changes: 3 additions & 3 deletions tardis/analysis/opacities.py
Original file line number Diff line number Diff line change
Expand Up @@ -372,13 +372,13 @@ def _calc_planck_mean_opacity(self):

for i in range(self.nshells):
delta_nu = self.nu_bins[1:] - self.nu_bins[:-1]
T = self.mdl.plasma.t_rad[i]
bb_nu = Blackbody(T)
temperature = self.mdl.plasma.t_rad[i]
bb_nu = Blackbody(temperature)

tmp = (
bb_nu(self.nu_bins[:-1]) * delta_nu * self.kappa_tot[:, 0]
).sum()
tmp /= (bb_nu(self.nu_bins[:-1], T) * delta_nu).sum()
tmp /= (bb_nu(self.nu_bins[:-1], temperature) * delta_nu).sum()

kappa_planck_mean[i] = tmp

Expand Down
6 changes: 3 additions & 3 deletions tardis/montecarlo/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -208,7 +208,7 @@ def _initialize_geometry_arrays(self, model):
self.v_outer_cgs = model.v_outer.to("cm/s").value

def _initialize_packets(
self, T, no_of_packets, iteration, radius, time_explosion
self, temperature, no_of_packets, iteration, radius, time_explosion
):
# the iteration is added each time to preserve randomness
# across different simulations with the same temperature,
Expand All @@ -220,11 +220,11 @@ def _initialize_packets(
seeds = rng.choice(MAX_SEED_VAL, no_of_packets, replace=True)
if not self.enable_full_relativity:
radii, nus, mus, energies = self.packet_source.create_packets(
T, no_of_packets, rng, radius
temperature, no_of_packets, rng, radius
)
else:
radii, nus, mus, energies = self.packet_source.create_packets(
T, no_of_packets, rng, radius, time_explosion
temperature, no_of_packets, rng, radius, time_explosion
)

mc_config_module.packet_seeds = seeds
Expand Down
6 changes: 3 additions & 3 deletions tardis/montecarlo/montecarlo_numba/formal_integral.py
Original file line number Diff line number Diff line change
Expand Up @@ -743,12 +743,12 @@ def trapezoid_integration(array, h):


@njit(**njit_dict_no_parallel)
def intensity_black_body(nu, T):
def intensity_black_body(nu, temperature):
"""Get the black body intensity at frequency nu
and temperature T"""
and temperature """
if nu == 0:
return np.nan # to avoid ZeroDivisionError
beta_rad = 1 / (KB_CGS * T)
beta_rad = 1 / (KB_CGS * temperature)
coefficient = 2 * H_CGS * C_INV * C_INV
return coefficient * nu * nu * nu / (np.exp(H_CGS * nu * beta_rad) - 1)

Expand Down
8 changes: 4 additions & 4 deletions tardis/montecarlo/montecarlo_numba/formal_integral_cuda.py
Original file line number Diff line number Diff line change
Expand Up @@ -475,23 +475,23 @@ def trapezoid_integration_cuda(arr, dx):


@cuda.jit(device=True)
def intensity_black_body_cuda(nu, T):
def intensity_black_body_cuda(nu, temperature):
"""
Get the black body intensity at frequency nu
and temperature T
and temperature
Parameters
----------
nu : float64
T : float64
temperature : float64
Returns
-------
float64
"""
if nu == 0:
return np.nan # to avoid ZeroDivisionError
beta_rad = 1 / (KB_CGS * T)
beta_rad = 1 / (KB_CGS * temperature)
coefficient = 2 * H_CGS * C_INV * C_INV
return coefficient * nu * nu * nu / (math.exp(H_CGS * nu * beta_rad) - 1)

Expand Down
4 changes: 2 additions & 2 deletions tardis/montecarlo/montecarlo_numba/interaction.py
Original file line number Diff line number Diff line change
Expand Up @@ -116,9 +116,9 @@ def sample_nu_free_free(numba_plasma, shell):
Frequency of the free-free emission process
"""
T = numba_plasma.t_electrons[shell]
temperature = numba_plasma.t_electrons[shell]
zrand = np.random.random()
return -K_B * T / H * np.log(zrand)
return -K_B * temperature / H * np.log(zrand)


@njit(**njit_dict_no_parallel)
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -32,36 +32,36 @@
not GPUs_available, reason="No GPU is available to test CUDA function"
)
@pytest.mark.parametrize(
["nu", "T"],
["nu", "temperature"],
[
(1e14, 1e4),
(0, 1),
(1, 1),
],
)
def test_intensity_black_body_cuda(nu, T):
def test_intensity_black_body_cuda(nu, temperature):
"""
Initializes the test of the cuda version
against the numba implementation of the
intensity_black_body to 15 decimals. This
is done as both results have 15 digits of precision.
"""
actual = np.zeros(3)
black_body_caller[1, 3](nu, T, actual)
black_body_caller[1, 3](nu, temperature, actual)

expected = formal_integral_numba.intensity_black_body(nu, T)
expected = formal_integral_numba.intensity_black_body(nu, temperature)

ntest.assert_allclose(actual, expected, rtol=1e-14)


@cuda.jit
def black_body_caller(nu, T, actual):
def black_body_caller(nu, temperature, actual):
"""
This calls the CUDA function and fills out
the array
"""
x = cuda.grid(1)
actual[x] = formal_integral_cuda.intensity_black_body_cuda(nu, T)
actual[x] = formal_integral_cuda.intensity_black_body_cuda(nu, temperature)


@pytest.mark.skipif(
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -12,18 +12,18 @@


@pytest.mark.parametrize(
["nu", "T"],
["nu", "temperature"],
[
(1e14, 1e4),
(0, 1),
(1, 1),
],
)
def test_intensity_black_body(nu, T):
def test_intensity_black_body(nu, temperature):
func = formal_integral.intensity_black_body
actual = func(nu, T)
actual = func(nu, temperature)
print(actual, type(actual))
expected = intensity_black_body(nu, T)
expected = intensity_black_body(nu, temperature)
ntest.assert_almost_equal(actual, expected)


Expand Down
21 changes: 9 additions & 12 deletions tardis/montecarlo/packet_source.py
Original file line number Diff line number Diff line change
Expand Up @@ -53,7 +53,7 @@ def create_uniform_packet_energies(no_of_packets, rng):
return np.ones(no_of_packets) / no_of_packets

@staticmethod
def create_blackbody_packet_nus(T, no_of_packets, rng, l_samples=1000):
def create_blackbody_packet_nus(temperature, no_of_packets, rng, l_samples=1000):
"""
Create packet :math:`\\nu` distributed using the algorithm described in
Bjorkman & Wood 2001 (page 4) which references
Expand All @@ -76,8 +76,7 @@ def create_blackbody_packet_nus(T, no_of_packets, rng, l_samples=1000):
Parameters
----------
T : float
temperature
temperature : float
no_of_packets : int
l_samples : int
number of l_samples needed in the algorithm
Expand All @@ -100,7 +99,7 @@ def create_blackbody_packet_nus(T, no_of_packets, rng, l_samples=1000):
xis_prod = np.prod(xis[1:], 0)
x = ne.evaluate("-log(xis_prod)/l")

return x * (const.k_B.cgs.value * T) / const.h.cgs.value
return x * (const.k_B.cgs.value * temperature) / const.h.cgs.value


class BlackBodySimpleSource(BasePacketSource):
Expand All @@ -109,13 +108,12 @@ class BlackBodySimpleSource(BasePacketSource):
part.
"""

def create_packets(self, T, no_of_packets, rng, radius):
def create_packets(self, temperature, no_of_packets, rng, radius):
"""Generate black-body packet properties as arrays
Parameters
----------
T : float64
Temperature
temperature : float64
no_of_packets : int
Number of packets
rng : numpy random number generator
Expand All @@ -134,21 +132,20 @@ def create_packets(self, T, no_of_packets, rng, radius):
Packet energies
"""
radii = np.ones(no_of_packets) * radius
nus = self.create_blackbody_packet_nus(T, no_of_packets, rng)
nus = self.create_blackbody_packet_nus(temperature, no_of_packets, rng)
mus = self.create_zero_limb_darkening_packet_mus(no_of_packets, rng)
energies = self.create_uniform_packet_energies(no_of_packets, rng)

return radii, nus, mus, energies


class BlackBodySimpleSourceRelativistic(BlackBodySimpleSource):
def create_packets(self, T, no_of_packets, rng, radius, time_explosion):
def create_packets(self, temperature, no_of_packets, rng, radius, time_explosion):
"""Generate relativistic black-body packet properties as arrays
Parameters
----------
T : float64
Temperature
temperature : float64
no_of_packets : int
Number of packets
rng : numpy random number generator
Expand All @@ -169,7 +166,7 @@ def create_packets(self, T, no_of_packets, rng, radius, time_explosion):
Packet energies
"""
self.beta = ((radius / time_explosion) / const.c).to("")
return super().create_packets(T, no_of_packets, rng, radius)
return super().create_packets(temperature, no_of_packets, rng, radius)

def create_zero_limb_darkening_packet_mus(self, no_of_packets, rng):
"""
Expand Down
4 changes: 2 additions & 2 deletions tardis/plasma/properties/continuum_processes.py
Original file line number Diff line number Diff line change
Expand Up @@ -990,9 +990,9 @@ def calculate(self, t_electrons):
@njit(error_model="numpy", fastmath=True)
def nu_ff(shell):

T = t_electrons[shell]
temperature = t_electrons[shell]
zrand = np.random.random()
return -K_B * T / H * np.log(zrand)
return -K_B * temperature / H * np.log(zrand)

return nu_ff

Expand Down
8 changes: 4 additions & 4 deletions tardis/util/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -285,27 +285,27 @@ def create_synpp_yaml(radial1d_mdl, fname, shell_no=0, lines_db=None):
yaml.dump(yaml_reference, stream=f, explicit_start=True)


def intensity_black_body(nu, T):
def intensity_black_body(nu, temperature):
"""
Calculate the intensity of a black-body according to the following formula
.. math::
I(\\nu, T) = \\frac{2h\\nu^3}{c^2}\\frac{1}
I(\\nu, temperature) = \\frac{2h\\nu^3}{c^2}\\frac{1}
{e^{h\\nu \\beta_\\textrm{rad}} - 1}
Parameters
----------
nu : float
Frequency of light
T : float
temperature : float
Temperature in kelvin
Returns
-------
Intensity : float
Returns the intensity of the black body
"""
beta_rad = 1 / (k_B_cgs * T)
beta_rad = 1 / (k_B_cgs * temperature)
coefficient = 2 * h_cgs / c_cgs**2
intensity = ne.evaluate(
"coefficient * nu**3 / " "(exp(h_cgs * nu * beta_rad) -1 )"
Expand Down

0 comments on commit 97047ce

Please sign in to comment.