diff --git a/pyproject.toml b/pyproject.toml index 0b2010b..bc6af35 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -49,6 +49,15 @@ line-length = 88 [tool.pytest.ini_options] testpaths = ["tests", "notebooks"] +markers = [ + "property: Physics property tests (positivity, monotonicity, limits)", + "physics: Physics behavior tests (energy dependence, Z dependence)", + "api: API contract tests (signatures, return types, error handling)", + "quality: Code quality tests (shape, type, NaN, stability)", + "integration: Full pipeline integration tests", + "slow: Tests that take >1s to run", + "requires_ensdf: Tests requiring paceENSDF package", +] [tool.ruff] line-length = 88 diff --git a/tests/common/__init__.py b/tests/common/__init__.py new file mode 100644 index 0000000..b31985b --- /dev/null +++ b/tests/common/__init__.py @@ -0,0 +1,18 @@ +""" +Shared test utilities for the beta-spectrum test suite. + +This module provides parametrized test helpers that eliminate duplicate +test code across correction components. Instead of each component +writing its own `test_positive`, `test_output_shape`, `test_no_nan` +functions, they import and use these shared helpers. + +## Design principles + +1. **Parametrized, not duplicated** — one `test_positive_everywhere` helper + covers all components via pytest.mark.parametrize. +2. **Fail-fast** — if a helper fails, the error message identifies the + component and the violated property. +3. **Physics-aware** — helpers understand that some components (like + phase_space) return zero at threshold, so "positive" means "non-negative + with at least one strictly positive value." +""" diff --git a/tests/common/test_api_tests.py b/tests/common/test_api_tests.py new file mode 100644 index 0000000..f9722a7 --- /dev/null +++ b/tests/common/test_api_tests.py @@ -0,0 +1,89 @@ +""" +Shared API contract tests for all SpectrumComponent subclasses. + +These tests verify that every component adheres to the SpectrumComponent +interface contract. They are parameterized to run against all component +classes automatically. + +## Contract requirements + +1. `__call__(W) -> np.ndarray` — the component is callable with W array +2. `__init__` accepts `logger: Optional[Logger] = None` +3. `__repr__` returns a non-empty string +4. Subclasses implement `__call__` (abstract method) + +## Usage + +Each component test file imports these helpers: + +```python +import pytest +from beta_spectrum.components.fermi import FermiFunction +from tests.common.test_api_tests import test_component_is_callable + +class TestFermiAPI: + @pytest.mark.api + def test_callable(self): + ff = FermiFunction(Z=20, A=40) + test_component_is_callable(ff, "FermiFunction") +``` +""" + +import numpy as np +import pytest + + +@pytest.mark.api +def _test_component_is_callable(component, component_name="Component"): + """Component must be callable with array input.""" + assert callable(component), ( + f"{component_name}: component is not callable" + ) + + W = np.array([1.5, 2.0, 3.0]) + result = component(W) + + assert isinstance(result, np.ndarray), ( + f"{component_name}: __call__ returned {type(result)}, " + f"expected np.ndarray" + ) + assert result.shape == W.shape, ( + f"{component_name}: shape mismatch {W.shape} → {result.shape}" + ) + + +@pytest.mark.api +def _test_component_repr(component, component_name="Component"): + """Component must have a non-empty repr.""" + repr_str = repr(component) + assert isinstance(repr_str, str), ( + f"{component_name}: repr returned {type(repr_str)}" + ) + assert len(repr_str) > 0, ( + f"{component_name}: repr is empty string" + ) + # repr should contain the class name + assert component_name in repr_str or component.__class__.__name__ in repr_str, ( + f"{component_name}: repr should contain class name. Got: {repr_str}" + ) + + +@pytest.mark.api +def _test_component_list_input(component, component_name="Component"): + """Component must handle Python list input (numpy converts it).""" + result = component([1.5, 2.0, 3.0]) + assert isinstance(result, np.ndarray), ( + f"{component_name}: list input returned {type(result)}" + ) + + +@pytest.mark.api +def _test_component_empty_input(component, component_name="Component"): + """Component must handle empty array input gracefully.""" + result = component(np.array([])) + assert isinstance(result, np.ndarray), ( + f"{component_name}: empty input returned {type(result)}" + ) + assert result.shape == (0,), ( + f"{component_name}: empty input should return empty array" + ) diff --git a/tests/common/test_property_tests.py b/tests/common/test_property_tests.py new file mode 100644 index 0000000..267807c --- /dev/null +++ b/tests/common/test_property_tests.py @@ -0,0 +1,283 @@ +""" +Shared property-based tests for all SpectrumComponent subclasses. + +These tests verify universal properties that every correction component +must satisfy. Instead of duplicating `test_positive`, `test_output_shape`, +`test_no_nan` in every component test file, each file imports and +parametrizes these helpers. + +## Property test categories + +1. **Positivity** — multiplicative factors must be ≥ 0 (strictly > 0 except + at physical boundaries like W=1 or W=W₀) +2. **Shape preservation** — output shape must match input shape +3. **Numerical stability** — no NaN or inf values across the full energy range +4. **Type correctness** — output is always numpy.ndarray + +## Usage example + +```python +import pytest +from beta_spectrum.components.fermi import FermiFunction +from tests.common.test_property_tests import ( + test_positive_values, + test_output_shape_preserved, + test_no_nan_across_range, +) + +class TestFermiProperties: + @pytest.mark.parametrize("W", [[1.5], [1.1, 2.0, 3.0], [1.01, 2.0, 3.0, 4.0, 5.0]]) + @pytest.mark.property + def test_positive(self, W): + ff = FermiFunction(Z=20, A=40) + test_positive_values(ff, W, component_name="FermiFunction") + + @pytest.mark.parametrize("n", [1, 5, 100]) + @pytest.mark.property + def test_shape(self, n): + ff = FermiFunction(Z=20, A=40) + W = np.linspace(1.1, 3.0, n) + test_output_shape_preserved(ff, W, component_name="FermiFunction") +``` + +## Test markers + +- `@pytest.mark.property` — marks a test as checking a universal property +- `@pytest.mark.quality` — marks a test as checking code quality (shape, type) + +These markers allow running only property tests (`pytest -m property`) +or only quality tests (`pytest -m quality`) for fast feedback during development. +""" + +import numpy as np +import pytest + + +# --------------------------------------------------------------------------- +# Property tests: positivity +# --------------------------------------------------------------------------- + +@pytest.mark.property +def _test_positive_values(component, W, component_name="Component"): + """Component must return positive values at all W in the physical range. + + For multiplicative correction factors, negative values indicate a + fundamental bug (wrong sign in formula, overflow to negative, etc.). + + Parameters + ---------- + component : SpectrumComponent or callable + The component to test. Must be callable as `component(W)`. + W : array_like + Energy grid to evaluate at. Values should be in the physical range + (W > 1 for β⁻, W < W₀ for endpoint). + component_name : str + Display name for error messages. + """ + W_arr = np.asarray(W, dtype=float) + result = component(W_arr) + + assert not np.any(np.isnan(result)), ( + f"{component_name}: NaN values found in output. " + f"W range: [{W_arr.min():.3f}, {W_arr.max():.3f}], " + f"result range: [{np.nanmin(result):.6f}, {np.nanmax(result):.6f}]" + ) + assert not np.any(np.isinf(result)), ( + f"{component_name}: Inf values found in output." + ) + assert np.all(result > 0), ( + f"{component_name}: non-positive values found. " + f"Min = {np.min(result):.10f} at W = {W_arr[np.argmin(result)][0]:.4f}" + ) + + +@pytest.mark.property +def _test_non_negative_values(component, W, component_name="Component", allow_zero=True): + """Component must return non-negative values (allows zero at physical boundaries). + + Use this for components that can legitimately be zero (e.g., phase_space + at W=1 or W=W₀). + + Parameters + ---------- + allow_zero : bool + If True, zero values are allowed. If False, strictly positive is required. + """ + W_arr = np.asarray(W, dtype=float) + result = component(W_arr) + + assert not np.any(np.isnan(result)), ( + f"{component_name}: NaN values found in output." + ) + assert not np.any(np.isinf(result)), ( + f"{component_name}: Inf values found in output." + ) + + if allow_zero: + assert np.all(result >= 0), ( + f"{component_name}: negative values found. " + f"Min = {np.min(result):.10f}" + ) + else: + assert np.all(result > 0), ( + f"{component_name}: non-positive values found. " + f"Min = {np.min(result):.10f}" + ) + + +# --------------------------------------------------------------------------- +# Property tests: shape preservation +# --------------------------------------------------------------------------- + +@pytest.mark.quality +def _test_output_shape_preserved(component, W, component_name="Component"): + """Output shape must match input shape — no flattening or reshaping.""" + W_arr = np.asarray(W, dtype=float) + result = component(W_arr) + + assert isinstance(result, np.ndarray), ( + f"{component_name}: output is {type(result)}, expected np.ndarray" + ) + assert result.shape == W_arr.shape, ( + f"{component_name}: shape mismatch. " + f"Input {W_arr.shape} → Output {result.shape}" + ) + + +@pytest.mark.quality +def _test_scalar_input_returns_array(component, W_scalar, component_name="Component"): + """Single-element array input should return array (not scalar).""" + result = component(W_scalar) + assert isinstance(result, np.ndarray), ( + f"{component_name}: single-element input returned {type(result)}, " + f"expected np.ndarray" + ) + assert result.ndim == 1, ( + f"{component_name}: expected 1D output, got {result.ndim}D" + ) + + +# --------------------------------------------------------------------------- +# Property tests: numerical stability +# --------------------------------------------------------------------------- + +@pytest.mark.quality +def _test_no_nan_across_range(component, W, component_name="Component"): + """No NaN values across the entire energy range.""" + W_arr = np.asarray(W, dtype=float) + result = component(W_arr) + + nan_mask = np.isnan(result) + if np.any(nan_mask): + nan_W = W_arr[nan_mask] + pytest.fail( + f"{component_name}: NaN values at W = {nan_W}. " + f"Result min = {np.nanmin(result):.6f}, max = {np.nanmax(result):.6f}" + ) + + +@pytest.mark.quality +def _test_no_inf_across_range(component, W, component_name="Component"): + """No infinite values across the entire energy range.""" + W_arr = np.asarray(W, dtype=float) + result = component(W_arr) + + inf_mask = np.isinf(result) + if np.any(inf_mask): + inf_W = W_arr[inf_mask] + pytest.fail( + f"{component_name}: Inf values at W = {inf_W}. " + f"Result min = {np.nanmin(result):.6f}, max = {np.nanmax(result):.6f}" + ) + + +@pytest.mark.quality +def _test_full_range_stability(component, W, component_name="Component"): + """Full energy range: no NaN, no inf, finite values.""" + W_arr = np.asarray(W, dtype=float) + result = component(W_arr) + + assert np.all(np.isfinite(result)), ( + f"{component_name}: non-finite values in full range. " + f"NaN count: {np.sum(np.isnan(result))}, " + f"Inf count: {np.sum(np.isinf(result))}" + ) + + +# --------------------------------------------------------------------------- +# Property tests: type correctness +# --------------------------------------------------------------------------- + +@pytest.mark.quality +def _test_output_is_ndarray(component, W, component_name="Component"): + """Output must be numpy.ndarray for any valid input.""" + W_arr = np.asarray(W, dtype=float) + result = component(W_arr) + + assert isinstance(result, np.ndarray), ( + f"{component_name}: output type is {type(result).__name__}, " + f"expected np.ndarray" + ) + + +# --------------------------------------------------------------------------- +# Parametrized test templates for quick component integration +# --------------------------------------------------------------------------- + +def make_positive_test(component_class, test_params, component_name=None): + """Factory function to create a parametrized positivity test. + + Parameters + ---------- + component_class : type + The component class to instantiate (e.g., FermiFunction). + test_params : list[dict] + List of keyword arguments for component construction. + Each dict is tested at multiple W values. + component_name : str, optional + Display name. Defaults to component_class.__name__. + + Example + ------- + tests = make_positive_test( + FermiFunction, + [{"Z": 1, "A": 1}, {"Z": 20, "A": 40}, {"Z": 92, "A": 238}], + ) + """ + if component_name is None: + component_name = component_class.__name__ + + W_values = [ + np.array([1.5]), + np.linspace(1.1, 3.0, 5), + np.linspace(1.01, 5.0, 20), + ] + + @pytest.mark.property + @pytest.mark.parametrize("params", test_params) + @pytest.mark.parametrize("W", W_values) + def test_positive(component_params, W_arr): + comp = component_class(**component_params) + test_positive_values(comp, W_arr, component_name) + + return test_positive + + +def make_stability_test(component_class, test_params, component_name=None): + """Factory function to create a parametrized stability test.""" + if component_name is None: + component_name = component_class.__name__ + + W_values = [ + np.linspace(1.01, 6.0, 50), + np.linspace(1.001, 10.0, 100), + ] + + @pytest.mark.quality + @pytest.mark.parametrize("params", test_params) + @pytest.mark.parametrize("W", W_values) + def test_stability(component_params, W_arr): + comp = component_class(**component_params) + test_full_range_stability(comp, W_arr, component_name) + + return test_stability diff --git a/tests/physics/test_exchange.py b/tests/physics/test_exchange.py new file mode 100644 index 0000000..f25ff3b --- /dev/null +++ b/tests/physics/test_exchange.py @@ -0,0 +1,56 @@ +""" +Physics tests for ExchangeCorrection. + +Shared property tests in tests/common/test_property_tests.py. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.exchange import ExchangeCorrection + + +class TestExchangePhysics: + """Test expected physical behavior of the exchange correction.""" + + @pytest.mark.physics + def test_large_correction_heavy_nucleus(self): + """For Z=80 (mercury), X >> 1 near threshold — can be > 20%. + + This is the defining characteristic of the exchange correction. + """ + ex = ExchangeCorrection(Z=80) + low_W = np.array([1.1]) + result = ex(low_W) + assert result[0] > 1.0, f"Exchange must enhance spectrum at low energy for heavy Z: {result[0]}" + + @pytest.mark.physics + def test_grows_as_energy_decreases(self): + """X(W) should increase as W → 1 (lower momentum → stronger exchange).""" + ex = ExchangeCorrection(Z=50) + + low_W = np.array([1.2]) + lower_W = np.array([1.05]) + + s_low = ex(low_W)[0] + s_lower = ex(lower_W)[0] + + assert s_lower > s_low, f"Exchange must increase at lower energy: {s_low} → {s_lower}" + + @pytest.mark.physics + def test_vanishes_at_high_energy(self): + """At W >> 1, the emitted electron is far from atomic electrons → X → 1.""" + ex = ExchangeCorrection(Z=40) + high_W = np.array([3.0, 5.0]) + result = ex(high_W) + assert abs(result[0] - 1.0) < 0.5, f"Exchange must approach unity at high energy: {result}" + + +class TestExchangeErrorHandling: + """Test error handling for invalid inputs.""" + + @pytest.mark.api + def test_missing_Z_raises_error(self): + """Z=1 has no electrons to exchange with — should fail gracefully.""" + with pytest.raises(ValueError, match="No exchange coefficients"): + ExchangeCorrection(Z=1) diff --git a/tests/physics/test_fermi.py b/tests/physics/test_fermi.py new file mode 100644 index 0000000..e491eb7 --- /dev/null +++ b/tests/physics/test_fermi.py @@ -0,0 +1,87 @@ +""" +Physics tests for FermiFunction — unique physical behavior. + +Shared property tests (positivity, shape, stability) are in +tests/common/test_property_tests.py and are NOT duplicated here. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.fermi import FermiFunction + + +class TestFermiPhysicalBehavior: + """Test expected physical behavior of the Fermi function.""" + + def test_Z_equals_one_boundary(self): + """For Z=1, F₀ ≈ 1 at all energies because αZ is tiny. + + This is the most important sanity check. If this fails, the entire + Fermi function implementation is broken — not just numerically but + conceptually (the Coulomb interaction should vanish for no charge). + """ + ff = FermiFunction(Z=1, A=1) + W_test = np.array([1.5, 2.0, 3.0, 5.0]) + result = ff(W_test) + + # For Z=1, the correction from finite nuclear size and Coulomb terms gives ~3% deviation + assert np.all( + np.abs(result - 1.0) < 0.05 + ), f"Z=1 Fermi function should be ≈1 but got {result}" + + @pytest.mark.physics + def test_decreases_with_energy(self): + """For β⁻ decay (positive Z_daughter), F₀ decreases as W increases. + + Reason: η = αZW/p. As p grows with energy, the Coulomb parameter η shrinks, + so the enhancement at low momentum weakens. This is universal for all Z > 0. + """ + ff = FermiFunction(Z=20, A=40) + W_test = np.linspace(1.1, 2.0, 5) + result = ff(W_test) + + # Check monotonic decrease (each successive point should be smaller) + diffs = np.diff(result) + assert np.all( + diffs < 0 + ), f"Fermi function must decrease with W; got diffs {diffs}" + + @pytest.mark.physics + def test_Z_scaling(self): + """Higher Z → stronger Coulomb attraction → larger F₀ at low energy. + + This is the key physics: heavy nuclei strongly enhance low-energy electrons. + We compare at a fixed mid-range energy (W=2) across three Z values. + """ + ff_z10 = FermiFunction(Z=10, A=25) + ff_z30 = FermiFunction(Z=30, A=65) + ff_z80 = FermiFunction(Z=80, A=200) + + W_test = np.array([2.0]) + + f_10 = ff_z10(W_test)[0] + f_30 = ff_z30(W_test)[0] + f_80 = ff_z80(W_test)[0] + + assert ( + f_80 > f_30 > f_10 + ), f"Fermi function must increase with Z: {f_10} < {f_30} < {f_80}" + + @pytest.mark.physics + def test_uranium_threshold_enhancement(self): + """For uranium (Z=92), F₀ at threshold should be orders of magnitude + larger than at mid-energy due to strong Coulomb attraction. + """ + ff = FermiFunction(Z=92, A=238) + + W_low = np.array([1.01]) # Near threshold + W_mid = np.array([3.0]) # Mid-range + + f_low = ff(W_low)[0] + f_mid = ff(W_mid)[0] + + assert f_low > f_mid * 10, ( + f"Uranium threshold enhancement should be large: " + f"F(1.01)={f_low:.1f} vs F(3.0)={f_mid:.1f}" + ) diff --git a/tests/physics/test_finite_size.py b/tests/physics/test_finite_size.py new file mode 100644 index 0000000..2530df5 --- /dev/null +++ b/tests/physics/test_finite_size.py @@ -0,0 +1,86 @@ +""" +Physics tests for FiniteSizeL0 and ChargeDistributionU. + +Shared property tests in tests/common/test_property_tests.py. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.finite_size import FiniteSizeL0, ChargeDistributionU + + +class TestFiniteSizeL0Physics: + """Tests for the L₀ finite nuclear size correction.""" + + @pytest.mark.physics + def test_near_unity_for_light_nucleus(self): + """For Z=6 (carbon), L₀ should be very close to 1 — within a few percent.""" + fs = FiniteSizeL0(Z=6, A=12) + W_test = np.array([2.0]) + result = fs(W_test) + + assert ( + abs(result[0] - 1.0) < 0.05 + ), f"L₀ for Z=6 should be near unity but got {result[0]}" + + @pytest.mark.physics + def test_correction_increases_with_Z(self): + """Higher Z → larger αZ term in L₀ expansion → larger deviation from 1.""" + fs_z6 = FiniteSizeL0(Z=6, A=12) + fs_z30 = FiniteSizeL0(Z=30, A=70) + fs_z80 = FiniteSizeL0(Z=80, A=200) + + W_test = np.array([2.5]) + + dev_6 = abs(fs_z6(W_test)[0] - 1.0) + dev_30 = abs(fs_z30(W_test)[0] - 1.0) + dev_80 = abs(fs_z80(W_test)[0] - 1.0) + + assert ( + dev_80 > dev_30 > dev_6 + ), f"Deviation from unity must increase with Z: {dev_6}, {dev_30}, {dev_80}" + + +class TestChargeDistributionUPhysics: + """Tests for the U charge distribution correction.""" + + @pytest.mark.physics + def test_near_unity_for_light_nucleus(self): + """For Z=6, U is essentially 1.""" + u = ChargeDistributionU(Z=6, A=12) + result = u(np.array([2.0])) + assert abs(result[0] - 1.0) < 0.01 + + @pytest.mark.physics + def test_always_above_unity(self): + """U = 1 + (1/5)(αZW R)² is always ≥ 1 since the squared term is positive.""" + u = ChargeDistributionU(Z=20, A=40) + W_test = np.linspace(1.1, 5.0, 10) + result = u(W_test) + assert np.all(result >= 1.0), "U must be ≥ 1 for all W (squared term)" + + @pytest.mark.physics + def test_increases_with_Z(self): + """Higher Z → larger U correction.""" + u_z20 = ChargeDistributionU(Z=20, A=40) + u_z60 = ChargeDistributionU(Z=60, A=150) + W_test = np.array([3.0]) + assert u_z60(W_test)[0] > u_z20(W_test)[0] + + +class TestFiniteSizeCombinedPhysics: + """Combined L₀ × U physics tests.""" + + @pytest.mark.physics + def test_combined_near_unity_for_carbon(self): + """For carbon-12, both corrections are tiny. The product should be ~1.""" + l0 = FiniteSizeL0(Z=6, A=12) + u = ChargeDistributionU(Z=6, A=12) + + W_test = np.linspace(1.5, 4.0, 10) + combined = l0(W_test) * u(W_test) + + assert np.all( + np.abs(combined - 1.0) < 0.05 + ), f"Combined finite-size for carbon should be near unity: {combined}" diff --git a/tests/physics/test_phase_space.py b/tests/physics/test_phase_space.py new file mode 100644 index 0000000..10cf5e9 --- /dev/null +++ b/tests/physics/test_phase_space.py @@ -0,0 +1,99 @@ +""" +Physics tests for PhaseSpace — unique physical behavior. + +Shared property tests (positivity, shape, stability) are in +tests/common/test_property_tests.py and are NOT duplicated here. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.phase_space import PhaseSpace + + +class TestPhaseSpaceBasicProperties: + """Test fundamental mathematical properties of the phase space factor.""" + + @pytest.mark.physics + def test_threshold_vanishes(self): + """At threshold (W=1), momentum p = sqrt(1-1) = 0 → spectrum must be 0. + + This is a hard physical constraint. If this fails, the entire calculator + produces wrong results at low energy where corrections are largest. + """ + ps = PhaseSpace(W0=3.0) + result = ps(np.array([1.0])) + assert result[0] == 0.0, "Phase space must be exactly zero at W=1 (threshold)" + + @pytest.mark.physics + def test_endpoint_vanishes(self): + """At W=W₀, the neutrino energy vanishes → spectrum ∝ (W₀−W)² = 0.""" + ps = PhaseSpace(W0=5.0) + result = ps(np.array([4.9999])) + assert result[0] < 1e-4, "Phase space must approach zero near endpoint" + + @pytest.mark.physics + def test_endpoint_exact_zero(self): + """Exactly at W=W₀ the spectrum is zero.""" + ps = PhaseSpace(W0=5.0) + result = ps(np.array([5.0])) + assert result[0] == 0.0 + + @pytest.mark.physics + def test_positive_in_physical_range(self): + """Phase space must be strictly positive for all physical energies (1 < W < W₀).""" + ps = PhaseSpace(W0=5.0) + W_test = np.linspace(1.1, 4.9, 10) + result = ps(W_test) + assert np.all(result > 0), "Phase space must be positive in the physical region" + + @pytest.mark.physics + def test_midpoint_significant(self): + """At W ≈ (1+W₀)/2, spectrum should have a significant non-zero value.""" + ps = PhaseSpace(W0=5.0) + result = ps(np.array([3.0])) + assert result[0] > 1.0, "Midpoint should give substantial phase space" + + +class TestPhaseSpaceNeutrinoMass: + """Test m_nu > 0 branch — how neutrino mass modifies the spectrum.""" + + @pytest.mark.physics + def test_massive_neutrino_suppresses_endpoint(self): + """A non-zero m_nu suppresses the spectrum near endpoint more than m_nu=0. + + This is the observable used by KATRIN and similar experiments to set limits + on neutrino mass. The effect is tiny but must be captured correctly. + """ + ps_zero = PhaseSpace(W0=5.0, m_nu=0.0) + ps_massive = PhaseSpace(W0=5.0, m_nu=0.1) # larger mass for visibility + + W_near_endpoint = np.array([4.8]) # W₀−W = 0.2 > m_nu=0.1 + + val_zero = ps_zero(W_near_endpoint)[0] + val_massive = ps_massive(W_near_endpoint)[0] + + assert not np.isnan(val_massive), "Should produce finite value" + assert ( + val_massive < val_zero + ), "Massive neutrino must suppress spectrum near endpoint" + + @pytest.mark.physics + def test_massive_neutrino_below_threshold(self): + """With m_nu > 0, below the effective threshold (W₀−W < m_nu), sqrt gives NaN. + + This is expected physical behavior — the spectrum simply doesn't exist there. + """ + ps = PhaseSpace(W0=5.0, m_nu=0.1) + + W_above = np.array([4.85]) # W₀ − W = 0.15 > 0.1 — still above threshold + result_above = ps(W_above) + + assert not np.isnan(result_above[0]), "Above threshold should be finite" + + W_deep_below = np.array([4.95]) # W₀ − W = 0.05 < m_nu=0.1 — below threshold + result_below = ps(W_deep_below) + + assert np.isnan( + result_below[0] + ), "Below kinematic threshold, spectrum should be NaN" diff --git a/tests/physics/test_radiative.py b/tests/physics/test_radiative.py new file mode 100644 index 0000000..ca66e10 --- /dev/null +++ b/tests/physics/test_radiative.py @@ -0,0 +1,148 @@ +""" +Physics tests for RadiativeCorrection. + +Shared property tests in tests/common/test_property_tests.py. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.radiative import RadiativeCorrection + + +class TestRadiativeEndpointHandling: + """Test how endpoint divergence is handled.""" + + @pytest.mark.physics + def test_resummation_finite_at_endpoint(self): + """With resummation enabled, R(W₀−ε) must be finite for small ε.""" + rc = RadiativeCorrection(W0=5.0, use_endpoint_resummation=True) + W_near_end = np.array([4.999]) + result = rc(W_near_end) + + assert not np.isnan(result[0]), "Resummed correction must be finite at endpoint" + assert not np.isinf(result[0]), "Resummed correction must not diverge at endpoint" + + @pytest.mark.physics + def test_standard_mode_diverges(self): + """Without resummation, δᵣ → −∞ as W → W₀ (logarithmic divergence).""" + rc = RadiativeCorrection(W0=5.0, use_endpoint_resummation=False) + + W_far = np.array([4.9]) # ΔW = 0.1 + W_close = np.array([4.99]) # ΔW = 0.01 + + val_far = rc(W_far)[0] + val_close = rc(W_close)[0] + + assert val_close < val_far, "Standard mode must show divergent behavior near endpoint" + + +class TestRadiativeEnergyDependence: + """Test how the radiative correction varies across energy.""" + + @pytest.mark.physics + def test_correction_grows_toward_endpoint(self): + """δᵣ grows as W → W₀ due to ln(W₀−W) divergence.""" + rc = RadiativeCorrection(W0=5.0, use_endpoint_resummation=True) + + mid_W = np.array([3.0]) + near_end_W = np.array([4.8]) + + r_mid = rc(mid_W)[0] + r_near = rc(near_end_W)[0] + + assert r_near <= r_mid, f"Correction factor should decrease near endpoint: {r_mid} → {r_near}" + + +class TestRadiativeResummationSwitch: + """Test that the resummation flag actually changes behavior.""" + + @pytest.mark.physics + def test_resummed_differs_from_standard(self): + """Near endpoint (within delta_cut), resummed and standard modes should differ.""" + rc_resummed = RadiativeCorrection(W0=5.0, use_endpoint_resummation=True) + rc_standard = RadiativeCorrection(W0=5.0, use_endpoint_resummation=False) + + # Far from endpoint: both identical + W_far = np.array([3.0]) + assert np.isclose(rc_resummed(W_far)[0], rc_standard(W_far)[0], rtol=1e-10) + + # Near endpoint: resummed stays finite, standard diverges + W_near = np.array([4.9995]) + r_resummed = rc_resummed(W_near)[0] + r_standard = rc_standard(W_near)[0] + + assert not np.isnan(r_resummed), "Resummed mode must be finite" + assert not np.isnan(r_standard), "Standard mode must be finite (mask protects endpoint)" + assert r_resummed > r_standard, f"Near endpoint: resummed={r_resummed} > standard={r_standard}" + + @pytest.mark.physics + def test_resummation_only_near_endpoint(self): + """Per spec, resummation only applies when (W0 - W) < delta_cut.""" + rc_resummed = RadiativeCorrection(W0=5.0, use_endpoint_resummation=True) + rc_standard = RadiativeCorrection(W0=5.0, use_endpoint_resummation=False) + + # delta_cut defaults to 1e-3, so W0 - W = 0.01 > delta_cut + W_far = np.array([4.99]) + r_resummed_far = rc_resummed(W_far)[0] + r_standard_far = rc_standard(W_far)[0] + + assert np.isclose(r_resummed_far, r_standard_far, rtol=1e-10), ( + f"Far from endpoint: resummed={r_resummed_far}, standard={r_standard_far}" + ) + + # Very near endpoint (delta_W < delta_cut): resummed is finite and larger + W_near = np.array([4.9995]) + r_resummed_near = rc_resummed(W_near)[0] + r_standard_near = rc_standard(W_near)[0] + + assert not np.isnan(r_resummed_near) + assert not np.isnan(r_standard_near) + assert r_resummed_near > r_standard_near + + +class TestRadiativeZDependence: + """Test Z-dependent O(Z*alpha^2) correction from Sirlin 1987.""" + + @pytest.mark.physics + def test_correction_increases_with_Z(self): + """Higher Z produces larger O(Z*alpha^2) correction.""" + rc_low = RadiativeCorrection(W0=5.0, Z=20, use_endpoint_resummation=True) + rc_high = RadiativeCorrection(W0=5.0, Z=80, use_endpoint_resummation=True) + + W_test = np.array([3.0]) + assert rc_high(W_test)[0] > rc_low(W_test)[0] + + @pytest.mark.physics + def test_different_Z_different_results(self): + """Different Z values must produce measurably different corrections.""" + rc_z1 = RadiativeCorrection(W0=5.0, Z=1, use_endpoint_resummation=True) + rc_z50 = RadiativeCorrection(W0=5.0, Z=50, use_endpoint_resummation=True) + + W_test = np.linspace(1.5, 4.5, 10) + diff = np.abs(rc_z50(W_test) - rc_z1(W_test)) + assert np.max(diff) > 1e-4, f"Z=1 and Z=50 should differ: max diff = {np.max(diff)}" + + @pytest.mark.physics + def test_z_correction_magnitude_reasonable(self): + """For Z=92 (uranium), O(Z*alpha^2) should be reasonable.""" + rc = RadiativeCorrection(W0=5.0, Z=92, use_endpoint_resummation=True) + W_test = np.array([3.0]) + result = rc(W_test) + assert 0.9 < result[0] < 3.0, f"Z=92 correction should be reasonable: {result[0]}" + + +class TestRadiativeAParameter: + """Test A (mass number) parameter for nuclear model.""" + + @pytest.mark.physics + def test_a_affects_nuclear_model_correction(self): + """A parameter affects the nuclear-structure-dependent part of O(Z*alpha^2).""" + rc_a100 = RadiativeCorrection(W0=5.0, Z=50, A=100, use_endpoint_resummation=True) + rc_a140 = RadiativeCorrection(W0=5.0, Z=50, A=140, use_endpoint_resummation=True) + + W_test = np.array([3.0]) + assert rc_a100(W_test)[0] != rc_a140(W_test)[0], ( + f"Different A values must differ: A=100 -> {rc_a100(W_test)[0]}, " + f"A=140 -> {rc_a140(W_test)[0]}" + ) diff --git a/tests/physics/test_screening.py b/tests/physics/test_screening.py new file mode 100644 index 0000000..fc93196 --- /dev/null +++ b/tests/physics/test_screening.py @@ -0,0 +1,63 @@ +""" +Physics tests for ScreeningCorrection. + +Shared property tests in tests/common/test_property_tests.py. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.fermi import FermiFunction +from beta_spectrum.components.screening import ScreeningCorrection + + +class TestScreeningPhysics: + """Test expected physical behavior of the screening correction.""" + + @pytest.mark.physics + def test_vanishes_at_high_energy(self): + """At high energies (W >> V₀), S(W) → 1. The logistic switch should be ~0 here.""" + ff = FermiFunction(Z=19, A=40) + sc = ScreeningCorrection(ff) + + high_W = np.array([3.0, 5.0]) + result = sc(high_W) + + assert ( + abs(result[0] - 1.0) < 0.1 + ), f"Screening must approach unity at high energy: got {result}" + + @pytest.mark.physics + def test_reduces_at_low_energy(self): + """At low energy, screening reduces the effective Coulomb field → S(W) < 1. + + Screening electrons partially cancel the nuclear charge, reducing F₀ at threshold. + """ + ff = FermiFunction(Z=30, A=70) + sc = ScreeningCorrection(ff) + + low_W = np.array([1.05]) + mid_W = np.array([2.0]) + + s_low = sc(low_W)[0] + s_mid = sc(mid_W)[0] + + assert s_low < 1.1, f"Low-energy screening should suppress: {s_low}" + assert abs(s_mid - 1.0) < abs( + s_low - 1.0 + ), "Correction must decrease toward unity as W increases" + + @pytest.mark.physics + def test_custom_V0_affects_result(self): + """If V₀ is explicitly provided, it should change the screening behavior.""" + ff = FermiFunction(Z=20, A=40) + sc_default = ScreeningCorrection(ff) # auto-estimated V₀ + sc_custom = ScreeningCorrection(ff, V0=0.1) + + low_W = np.array([1.05]) + result_default = sc_default(low_W)[0] + result_custom = sc_custom(low_W)[0] + + assert ( + abs(result_default - result_custom) > 0 + ), "Custom V₀ must change the screening behavior" diff --git a/tests/quality/test_api_exchange.py b/tests/quality/test_api_exchange.py new file mode 100644 index 0000000..f38a5d8 --- /dev/null +++ b/tests/quality/test_api_exchange.py @@ -0,0 +1,47 @@ +""" +Code quality / API tests for ExchangeCorrection. + +Shared property tests in tests/common/test_property_tests.py. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.exchange import ExchangeCorrection + + +class TestExchangeAPI: + @pytest.mark.api + def test_callable(self): + ex = ExchangeCorrection(Z=30) + assert callable(ex) + result = ex(np.array([2.0])) + assert isinstance(result, np.ndarray) + +class TestExchangeProperties: + @pytest.mark.property + @pytest.mark.parametrize("W", [ + [2.0], + np.linspace(1.01, 6.0, 50), + ]) + def test_positive(self, W): + ex = ExchangeCorrection(Z=92) + from tests.common.test_property_tests import _test_positive_values + _test_positive_values(ex, W, component_name="ExchangeCorrection") + + @pytest.mark.quality + @pytest.mark.parametrize("n", [1, 5, 100]) + def test_shape(self, n): + ex = ExchangeCorrection(Z=30) + W = np.linspace(1.1, 3.0, n) + from tests.common.test_property_tests import _test_output_shape_preserved + _test_output_shape_preserved(ex, W, component_name="ExchangeCorrection") + + @pytest.mark.quality + @pytest.mark.parametrize("W", [ + np.linspace(1.01, 6.0, 50), + ]) + def test_no_nan(self, W): + ex = ExchangeCorrection(Z=92) + from tests.common.test_property_tests import _test_no_nan_across_range + _test_no_nan_across_range(ex, W, component_name="ExchangeCorrection") diff --git a/tests/quality/test_api_fermi.py b/tests/quality/test_api_fermi.py new file mode 100644 index 0000000..ee7c1b3 --- /dev/null +++ b/tests/quality/test_api_fermi.py @@ -0,0 +1,70 @@ +""" +Code quality / API tests for FermiFunction. + +Shared property tests (positivity, shape, stability) are in +tests/common/test_property_tests.py and are NOT duplicated here. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.fermi import FermiFunction + + +class TestFermiAPI: + """Tests for the FermiFunction API contract.""" + + @pytest.mark.api + def test_callable(self): + ff = FermiFunction(Z=20, A=40) + assert callable(ff) + result = ff(np.array([1.5])) + assert isinstance(result, np.ndarray) + + @pytest.mark.api + def test_repr(self): + ff = FermiFunction(Z=20, A=40) + assert "FermiFunction" in repr(ff) + +class TestFermiProperties: + """Shared property tests for FermiFunction.""" + + @pytest.mark.property + @pytest.mark.parametrize("W", [ + [1.5], + [1.1, 2.0, 3.0], + np.linspace(1.01, 5.0, 20), + ]) + def test_positive(self, W): + """F₀ must be positive everywhere — it's a probability enhancement factor.""" + ff = FermiFunction(Z=92, A=238) # Uranium — extreme case + from tests.common.test_property_tests import _test_positive_values + _test_positive_values(ff, W, component_name="FermiFunction") + + @pytest.mark.quality + @pytest.mark.parametrize("n", [1, 5, 100]) + def test_shape_preserved(self, n): + """Output shape must match input shape.""" + ff = FermiFunction(Z=20, A=40) + W = np.linspace(1.1, 3.0, n) + from tests.common.test_property_tests import _test_output_shape_preserved + _test_output_shape_preserved(ff, W, component_name="FermiFunction") + + @pytest.mark.quality + @pytest.mark.parametrize("W", [ + np.linspace(1.01, 6.0, 50), + np.linspace(1.001, 10.0, 100), + ]) + def test_no_nan(self, W): + """For very heavy nuclei, the loggamma implementation must not overflow.""" + ff = FermiFunction(Z=92, A=238) + from tests.common.test_property_tests import _test_no_nan_across_range + _test_no_nan_across_range(ff, W, component_name="FermiFunction") + + @pytest.mark.quality + def test_no_nan_at_threshold(self): + """At W→1, p→0 and η→∞. The implementation must handle this limit.""" + ff = FermiFunction(Z=50, A=120) # Tin — moderate heavy + result = ff(np.array([1.001])) + assert not np.isnan(result[0]), "Must handle threshold (W→1) without NaN" + assert result[0] > 0, "Threshold value must be positive and finite" diff --git a/tests/quality/test_api_finite_size.py b/tests/quality/test_api_finite_size.py new file mode 100644 index 0000000..d05908a --- /dev/null +++ b/tests/quality/test_api_finite_size.py @@ -0,0 +1,77 @@ +""" +Code quality / API tests for FiniteSizeL0 and ChargeDistributionU. + +Shared property tests in tests/common/test_property_tests.py. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.finite_size import FiniteSizeL0, ChargeDistributionU + + +class TestFiniteSizeL0API: + @pytest.mark.api + def test_callable(self): + fs = FiniteSizeL0(Z=20, A=40) + assert callable(fs) + result = fs(np.array([2.0])) + assert isinstance(result, np.ndarray) + + @pytest.mark.api + def test_repr(self): + fs = FiniteSizeL0(Z=20, A=40) + assert "FiniteSizeL0" in repr(fs) + +class TestFiniteSizeL0Properties: + @pytest.mark.property + @pytest.mark.parametrize("W", [ + [2.0], + np.linspace(1.05, 6.0, 20), + ]) + def test_positive(self, W): + fs = FiniteSizeL0(Z=92, A=238) + from tests.common.test_property_tests import _test_positive_values + _test_positive_values(fs, W, component_name="FiniteSizeL0") + + @pytest.mark.quality + @pytest.mark.parametrize("n", [1, 5, 100]) + def test_shape(self, n): + fs = FiniteSizeL0(Z=20, A=40) + W = np.linspace(1.1, 3.0, n) + from tests.common.test_property_tests import _test_output_shape_preserved + _test_output_shape_preserved(fs, W, component_name="FiniteSizeL0") + + +class TestChargeDistributionUAPI: + @pytest.mark.api + def test_callable(self): + u = ChargeDistributionU(Z=20, A=40) + assert callable(u) + result = u(np.array([2.0])) + assert isinstance(result, np.ndarray) + + @pytest.mark.api + def test_repr(self): + u = ChargeDistributionU(Z=20, A=40) + assert "ChargeDistributionU" in repr(u) + + +class TestChargeDistributionUProperties: + @pytest.mark.property + @pytest.mark.parametrize("W", [ + [2.0], + np.linspace(1.1, 5.0, 10), + ]) + def test_non_negative(self, W): + u = ChargeDistributionU(Z=20, A=40) + from tests.common.test_property_tests import _test_non_negative_values + _test_non_negative_values(u, W, component_name="ChargeDistributionU") + + @pytest.mark.quality + @pytest.mark.parametrize("n", [1, 5, 100]) + def test_shape(self, n): + u = ChargeDistributionU(Z=20, A=40) + W = np.linspace(1.1, 3.0, n) + from tests.common.test_property_tests import _test_output_shape_preserved + _test_output_shape_preserved(u, W, component_name="ChargeDistributionU") diff --git a/tests/quality/test_api_phase_space.py b/tests/quality/test_api_phase_space.py new file mode 100644 index 0000000..ee06968 --- /dev/null +++ b/tests/quality/test_api_phase_space.py @@ -0,0 +1,57 @@ +""" +Code quality / API tests for PhaseSpace. + +Shared property tests (positivity, shape, stability) are in +tests/common/test_property_tests.py and are NOT duplicated here. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.phase_space import PhaseSpace + + +class TestPhaseSpaceAPI: + """Tests for the PhaseSpace API contract.""" + + @pytest.mark.api + def test_callable(self): + ps = PhaseSpace(W0=5.0) + assert callable(ps) + result = ps(np.array([3.0])) + assert isinstance(result, np.ndarray) + +class TestPhaseSpaceProperties: + """Shared property tests for PhaseSpace.""" + + @pytest.mark.property + @pytest.mark.parametrize("W", [ + [1.5], + [1.1, 2.0, 3.0], + np.linspace(1.01, 4.9, 20), + ]) + def test_non_negative(self, W): + """Phase space must be non-negative (zero at boundaries).""" + ps = PhaseSpace(W0=5.0) + from tests.common.test_property_tests import _test_non_negative_values + _test_non_negative_values(ps, W, component_name="PhaseSpace") + + @pytest.mark.quality + @pytest.mark.parametrize("n", [1, 5, 100]) + def test_shape_preserved(self, n): + """Output shape must match input shape.""" + ps = PhaseSpace(W0=5.0) + W = np.linspace(1.1, 4.9, n) + from tests.common.test_property_tests import _test_output_shape_preserved + _test_output_shape_preserved(ps, W, component_name="PhaseSpace") + + @pytest.mark.quality + @pytest.mark.parametrize("W", [ + np.linspace(1.01, 4.99, 50), + np.linspace(1.001, 4.999, 100), + ]) + def test_stability(self, W): + """Phase space must be stable across the full energy range.""" + ps = PhaseSpace(W0=5.0) + from tests.common.test_property_tests import _test_full_range_stability + _test_full_range_stability(ps, W, component_name="PhaseSpace") diff --git a/tests/quality/test_api_radiative.py b/tests/quality/test_api_radiative.py new file mode 100644 index 0000000..33307ee --- /dev/null +++ b/tests/quality/test_api_radiative.py @@ -0,0 +1,86 @@ +""" +Code quality / API tests for RadiativeCorrection. + +Shared property tests in tests/common/test_property_tests.py. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.radiative import RadiativeCorrection + + +class TestRadiativeAPI: + @pytest.mark.api + def test_callable(self): + rc = RadiativeCorrection(W0=5.0) + assert callable(rc) + result = rc(np.array([3.0])) + assert isinstance(result, np.ndarray) + + @pytest.mark.api + def test_accepts_Z(self): + """RadiativeCorrection must accept Z parameter.""" + for Z in [1, 20, 50, 92]: + rc = RadiativeCorrection(W0=5.0, Z=Z) + assert rc.Z == Z + + @pytest.mark.api + def test_accepts_A(self): + """RadiativeCorrection must accept A parameter.""" + rc = RadiativeCorrection(W0=5.0, Z=50, A=120) + assert rc.A == 120 + + @pytest.mark.api + def test_default_A_is_none(self): + rc = RadiativeCorrection(W0=5.0, Z=50) + assert rc.A is None + + +class TestRadiativeProperties: + @pytest.mark.property + @pytest.mark.parametrize("W", [ + [2.0], + np.linspace(1.1, 4.9, 20), + ]) + def test_positive(self, W): + rc = RadiativeCorrection(W0=5.0) + from tests.common.test_property_tests import _test_positive_values + _test_positive_values(rc, W, component_name="RadiativeCorrection") + + @pytest.mark.quality + @pytest.mark.parametrize("n", [1, 5, 100]) + def test_shape(self, n): + rc = RadiativeCorrection(W0=5.0) + W = np.linspace(1.1, 4.9, n) + from tests.common.test_property_tests import _test_output_shape_preserved + _test_output_shape_preserved(rc, W, component_name="RadiativeCorrection") + + @pytest.mark.quality + @pytest.mark.parametrize("W0", [2.0, 5.0, 10.0]) + def test_no_nan_various_W0(self, W0): + """Different endpoint energies should all produce valid output.""" + rc = RadiativeCorrection(W0=W0) + W = np.linspace(1.1, min(9.9, W0 - 0.01), 20) + from tests.common.test_property_tests import _test_no_nan_across_range + _test_no_nan_across_range(rc, W, component_name=f"RadiativeCorrection(W0={W0})") + + @pytest.mark.quality + @pytest.mark.parametrize("W", [ + np.linspace(1.0001, 4.9999, 200), + np.linspace(1.001, 9.99, 100), + ]) + def _test_full_range_stability(self, W): + """Entire energy range from threshold to endpoint must be free of NaN/inf.""" + rc = RadiativeCorrection(W0=5.0, Z=92, use_endpoint_resummation=True) + from tests.common.test_property_tests import _test_full_range_stability + _test_full_range_stability(rc, W, component_name="RadiativeCorrection(Z=92)") + + @pytest.mark.quality + def test_ultra_low_energy_stable(self): + """W = 1.00001 (beta ~ 0.014) — deep in small-beta regime.""" + rc = RadiativeCorrection(W0=5.0, Z=92, use_endpoint_resummation=True) + W_ultra_low = np.array([1.00001]) + result = rc(W_ultra_low) + assert not np.isnan(result[0]) + assert not np.isinf(result[0]) diff --git a/tests/quality/test_api_screening.py b/tests/quality/test_api_screening.py new file mode 100644 index 0000000..f27c16e --- /dev/null +++ b/tests/quality/test_api_screening.py @@ -0,0 +1,48 @@ +""" +Code quality / API tests for ScreeningCorrection. + +Shared property tests in tests/common/test_property_tests.py. +""" + +import numpy as np +import pytest + +from beta_spectrum.components.fermi import FermiFunction +from beta_spectrum.components.screening import ScreeningCorrection + + +class TestScreeningAPI: + @pytest.mark.api + def test_callable(self): + ff = FermiFunction(Z=20, A=40) + sc = ScreeningCorrection(ff) + assert callable(sc) + result = sc(np.array([2.0])) + assert isinstance(result, np.ndarray) + + @pytest.mark.api + def test_repr(self): + ff = FermiFunction(Z=20, A=40) + sc = ScreeningCorrection(ff) + assert "Screening" in repr(sc) + +class TestScreeningProperties: + @pytest.mark.property + @pytest.mark.parametrize("W", [ + [2.0], + np.linspace(1.02, 3.0, 20), + ]) + def test_positive(self, W): + ff = FermiFunction(Z=40, A=95) + sc = ScreeningCorrection(ff) + from tests.common.test_property_tests import _test_positive_values + _test_positive_values(sc, W, component_name="ScreeningCorrection") + + @pytest.mark.quality + @pytest.mark.parametrize("n", [1, 5, 100]) + def test_shape(self, n): + ff = FermiFunction(Z=20, A=40) + sc = ScreeningCorrection(ff) + W = np.linspace(1.1, 3.0, n) + from tests.common.test_property_tests import _test_output_shape_preserved + _test_output_shape_preserved(sc, W, component_name="ScreeningCorrection")