From a96aef19c79be302d2bc6963df801d1fe4817f5c Mon Sep 17 00:00:00 2001 From: "Micah D. Gale" Date: Tue, 4 Nov 2025 09:27:14 -0600 Subject: [PATCH 1/4] Added test example that keeps showing up. --- tests/test_nuclide.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/tests/test_nuclide.py b/tests/test_nuclide.py index b3b761d7..88978881 100644 --- a/tests/test_nuclide.py +++ b/tests/test_nuclide.py @@ -1,6 +1,6 @@ # Copyright 2024, Battelle Energy Alliance, LLC All Rights Reserved. import pytest -from hypothesis import assume, given, note, strategies as st, settings +from hypothesis import assume, example, given, note, strategies as st, settings import montepy @@ -383,6 +383,7 @@ def test_particle_str(_): class TestNucleus: + @example(Z=97, A=185, meta=2) @given(Z=st.integers(1, 99), A=st.integers(0, 300), meta=st.integers(0, 4)) def test_nucleus_init_eq_hash(_, Z, A, meta): # avoid metastable elemental From e5ca4d868c4bd92004de4b9bac2fffb61cca92c7 Mon Sep 17 00:00:00 2001 From: "Micah D. Gale" Date: Tue, 4 Nov 2025 10:08:27 -0600 Subject: [PATCH 2/4] Implemented semi-emperical mass formula. --- tests/test_nuclide.py | 42 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 42 insertions(+) diff --git a/tests/test_nuclide.py b/tests/test_nuclide.py index 88978881..54253094 100644 --- a/tests/test_nuclide.py +++ b/tests/test_nuclide.py @@ -1,6 +1,7 @@ # Copyright 2024, Battelle Energy Alliance, LLC All Rights Reserved. import pytest from hypothesis import assume, example, given, note, strategies as st, settings +import numpy as np import montepy @@ -383,11 +384,52 @@ def test_particle_str(_): class TestNucleus: + @staticmethod + def _semi_emper_binding_energy_p_nucleon(Z, A): + """ + Calculates the base nuclide's binding energy using the semi-empirical mass-formula. + + Based on Wikipedia: + """ + CONSTANTS = { + "A_V": 15.8, + "A_S": 18.3, + "A_C": 0.714, + "A_A": 23.2, + "A_P": 12, + "K_P": -1 / 2, + } + """ + Units MeV + based on "least-squares (1)" from Alonso, Finn. + """ + FUNCTIONS = { + "A_V": lambda z, a: a, + "A_S": lambda z, a: np.pow(a, 2 / 3), + "A_C": lambda z, a: (z * (z - 1)) / (np.pow(a, 1 / 3)), + "A_A": lambda z, a: (a - 2 * z) ** 2 / a, + } + energy = 0.0 + for constant, func in zip(CONSTANTS.values(), FUNCTIONS.values()): + energy += constant * func(Z, A) + # handle even odd stuff + N = A - Z + delta = 34 / (np.pow(A, 3 / 4)) + even_n = N % 2 == 0 + even_z = Z % 2 == 0 + parity_match = even_z == even_z + # 0 for parity mismatch + if not parity_match: + return energy / A + energy += delta * 1 if even_n else -1 + return energy / A + @example(Z=97, A=185, meta=2) @given(Z=st.integers(1, 99), A=st.integers(0, 300), meta=st.integers(0, 4)) def test_nucleus_init_eq_hash(_, Z, A, meta): # avoid metastable elemental assume((A == 0) == (meta == 0)) + print(_._semi_emper_binding_energy_p_nucleon(Z, A)) nucleus = Nucleus(Element(Z), A, meta) assert nucleus.Z == Z assert nucleus.A == A From 54382af34a0e69633c7b05de27839189bdc6446c Mon Sep 17 00:00:00 2001 From: "Micah D. Gale" Date: Tue, 4 Nov 2025 10:35:28 -0600 Subject: [PATCH 3/4] Debugged SEMF. --- tests/test_nuclide.py | 21 +++++++++++---------- 1 file changed, 11 insertions(+), 10 deletions(-) diff --git a/tests/test_nuclide.py b/tests/test_nuclide.py index 54253094..1e5bcb84 100644 --- a/tests/test_nuclide.py +++ b/tests/test_nuclide.py @@ -392,29 +392,30 @@ def _semi_emper_binding_energy_p_nucleon(Z, A): Based on Wikipedia: """ CONSTANTS = { - "A_V": 15.8, - "A_S": 18.3, - "A_C": 0.714, - "A_A": 23.2, - "A_P": 12, - "K_P": -1 / 2, + "A_V": 15.75, + "A_S": 17.8, + "A_C": 0.711, + "A_A": 23.7, + "A_P": 11.18, } """ Units MeV - based on "least-squares (1)" from Alonso, Finn. + based on Rohlf """ FUNCTIONS = { "A_V": lambda z, a: a, - "A_S": lambda z, a: np.pow(a, 2 / 3), - "A_C": lambda z, a: (z * (z - 1)) / (np.pow(a, 1 / 3)), + "A_S": lambda z, a: np.pow(a, 2.0 / 3), + "A_C": lambda z, a: (z * (z - 1)) / (np.pow(a, 1.0 / 3)), "A_A": lambda z, a: (a - 2 * z) ** 2 / a, } energy = 0.0 for constant, func in zip(CONSTANTS.values(), FUNCTIONS.values()): + print(func(Z, A), constant * func(Z, A)) energy += constant * func(Z, A) + print(energy, energy / A) # handle even odd stuff N = A - Z - delta = 34 / (np.pow(A, 3 / 4)) + delta = CONSTANTS["A_P"] / (np.pow(A, 1 / 2)) even_n = N % 2 == 0 even_z = Z % 2 == 0 parity_match = even_z == even_z From ff0772d901d444640e516092d5d8abab3a2564e1 Mon Sep 17 00:00:00 2001 From: "Micah D. Gale" Date: Tue, 4 Nov 2025 10:55:38 -0600 Subject: [PATCH 4/4] Fixed SEMF, I missed the minus signs. --- tests/test_nuclide.py | 17 ++++++++++------- 1 file changed, 10 insertions(+), 7 deletions(-) diff --git a/tests/test_nuclide.py b/tests/test_nuclide.py index 1e5bcb84..8f3124c5 100644 --- a/tests/test_nuclide.py +++ b/tests/test_nuclide.py @@ -393,9 +393,9 @@ def _semi_emper_binding_energy_p_nucleon(Z, A): """ CONSTANTS = { "A_V": 15.75, - "A_S": 17.8, - "A_C": 0.711, - "A_A": 23.7, + "A_S": -17.8, + "A_C": -0.711, + "A_A": -23.7, "A_P": 11.18, } """ @@ -409,10 +409,8 @@ def _semi_emper_binding_energy_p_nucleon(Z, A): "A_A": lambda z, a: (a - 2 * z) ** 2 / a, } energy = 0.0 - for constant, func in zip(CONSTANTS.values(), FUNCTIONS.values()): - print(func(Z, A), constant * func(Z, A)) + for (name, constant), func in zip(CONSTANTS.items(), FUNCTIONS.values()): energy += constant * func(Z, A) - print(energy, energy / A) # handle even odd stuff N = A - Z delta = CONSTANTS["A_P"] / (np.pow(A, 1 / 2)) @@ -430,7 +428,12 @@ def _semi_emper_binding_energy_p_nucleon(Z, A): def test_nucleus_init_eq_hash(_, Z, A, meta): # avoid metastable elemental assume((A == 0) == (meta == 0)) - print(_._semi_emper_binding_energy_p_nucleon(Z, A)) + # use SEMF to rule out non-physical nuclides + # start caring after iron + if Z > 26 and A > 0: + assume( + _._semi_emper_binding_energy_p_nucleon(Z, A) > 7.5 + ) # 7.5 MeV per nucleon cutoff nucleus = Nucleus(Element(Z), A, meta) assert nucleus.Z == Z assert nucleus.A == A