Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix bug in relativistic packet source #2453

Merged
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions tardis/montecarlo/packet_source.py
Original file line number Diff line number Diff line change
Expand Up @@ -313,7 +313,7 @@ def create_packets(self, no_of_packets):
self.beta = ((self.radius / self.time_explosion) / const.c).to("")
return super().create_packets(no_of_packets)

def create_packet_nus(self, no_of_packets):
def create_packet_mus(self, no_of_packets):
"""
Create zero-limb-darkening packet :math:`\mu^\prime` distributed
according to :math:`\\mu^\\prime=2 \\frac{\\mu^\\prime + \\beta}{2 \\beta + 1}`.
Expand All @@ -327,7 +327,7 @@ def create_packet_nus(self, no_of_packets):

Returns
-------
array of frequencies
Directions for packets
numpy.ndarray
"""
z = self.rng.random(no_of_packets)
Expand Down
68 changes: 65 additions & 3 deletions tardis/montecarlo/tests/test_packet_source.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,8 +3,12 @@
import numpy as np
import pandas as pd
import pytest
from numpy.testing import assert_allclose

from tardis.montecarlo.packet_source import BlackBodySimpleSource
from tardis.montecarlo.packet_source import (
BlackBodySimpleSource,
BlackBodySimpleSourceRelativistic,
)
from tardis.montecarlo import (
montecarlo_configuration as montecarlo_configuration,
)
Expand Down Expand Up @@ -43,6 +47,22 @@ def blackbodysimplesource(self, request):
yield cls.bb
montecarlo_configuration.LEGACY_MODE_ENABLED = False

@pytest.fixture(scope="class")
def blackbody_simplesource_relativistic(self, request):
"""
Create BlackBodySimpleSourceRelativistic instance.

Yields
-------
tardis.montecarlo.packet_source.BlackBodySimpleSourceRelativistic
"""
montecarlo_configuration.LEGACY_MODE_ENABLED = True
bb_rel = BlackBodySimpleSourceRelativistic(
base_seed=1963, legacy_second_seed=2508
)
yield bb_rel
montecarlo_configuration.LEGACY_MODE_ENABLED = False

def test_bb_packet_sampling(
self,
request,
Expand All @@ -51,8 +71,6 @@ def test_bb_packet_sampling(
blackbodysimplesource,
):
"""
Test generate_plot_mpl method.

Parameters
----------
request : _pytest.fixtures.SubRequest
Expand All @@ -75,3 +93,47 @@ def test_bb_packet_sampling(
assert np.all(np.isclose(nus, ref_df["nus"]))
assert np.all(np.isclose(mus, ref_df["mus"]))
assert np.all(np.isclose(unif_energies, ref_df["energies"]))

def test_bb_packet_sampling_relativistic(
self,
tardis_ref_data,
blackbody_simplesource_relativistic,
):
"""
Parameters
----------
tardis_ref_data : pd.HDFStore
blackbody_simplesource_relativistic : tardis.montecarlo.packet_source.BlackBodySimpleSourceRelativistic
"""
blackbody_simplesource_relativistic.temperature = 10000
blackbody_simplesource_relativistic.beta = 0.25

nus = blackbody_simplesource_relativistic.create_packet_nus(100)
unif_energies = (
blackbody_simplesource_relativistic.create_packet_energies(100)
)
blackbody_simplesource_relativistic._reseed(2508)
mus = blackbody_simplesource_relativistic.create_packet_mus(10)

gamma = np.sqrt(1 - blackbody_simplesource_relativistic.beta**2) ** -1
ref_df = tardis_ref_data["/packet_unittest/blackbody"]
expected_nus = ref_df["nus"]
expected_unif_energies = ref_df["energies"] * 1.6 / gamma
chvogl marked this conversation as resolved.
Show resolved Hide resolved
expected_mus = np.array(
[
0.60420546,
0.49899691,
0.69583288,
0.96812652,
0.01544154,
0.93562304,
0.44306545,
0.77010037,
0.896973,
0.67876489,
]
)

assert_allclose(nus, expected_nus)
assert_allclose(unif_energies, expected_unif_energies)
assert_allclose(mus, expected_mus, rtol=1e-6)
Loading