Keep the lifetime clock and delayed group of (n,xn) and split neutrons

[ondrejch]

Description

Fixes #4177

Problem

Two kinds of neutron are banked as SourceSites and started again by Particle::from_source(), although they continue the generation of the neutron that produced them:

• the extra neutrons of an (n,xn) reaction with an integral yield, banked by create_secondary() in inelastic_scatter();
• the neutrons made by weight-window splitting, banked by split() in apply_weight_window().

They lose two pieces of state that the iterated fission probability (IFP) method reads:

• The lifetime clock. from_source() sets lifetime() = 0.0 for every particle, while it restores time() from the site. ifp() stores p.lifetime() with each fission site, and ifp-time-numerator scores it. The IFP generation time is biased low.
• The delayed group. split() and create_secondary() do not copy delayed_group, so the copy starts with the site's default of 0. In continuous-energy mode a neutron's delayed_group() does not change during transport, so it is the neutron's birth tag. ifp() stores it with each fission site, and ifp-beta-numerator counts a fission when the tag of its ancestor is nonzero. The copies of a delayed neutron therefore count as prompt, and the IFP beta_eff is biased low.

The two (n,xn) paths of inelastic_scatter() are also inconsistent. An integral yield banks yield - 1 identical copies, which lose the clock and the tag. A non-integral yield multiplies the weight of the incident neutron instead, and that neutron keeps both.

Size of the bias:

• Weight windows. In a model where most fission neutrons are split, the generation time is about 14 % low and beta_eff 24 % low.
• (n,xn), generation time. About 1e-4 relative in a bare U-235 sphere (U-235 has (n,2n) and (n,3n) reactions), and 3.2e-3 in a U-235/Be-9 mixture.
• (n,xn), beta_eff. Negligible in most systems, because delayed neutrons are born below most (n,2n) thresholds (see "The (n,xn) case" below).

k-effective, ifp-denominator and every non-IFP tally are not affected.

Physics

The IFP kinetics parameters follow the adjoint-weighted formulation of MCNP's KOPTS option [1], with

Lambda   = <psi^dagger, v^-1 psi> / <psi^dagger, F psi>      (Eq. (1) of [1])
beta_eff = <psi^dagger, F_d psi>  / <psi^dagger, F psi>

where F is the fission production operator and F_d its delayed part. Eq. (8) of [1] scores the Lambda numerator by track length. It sums track length over speed over every track of the original generation, each weighted by the asymptotic population of the progeny that carry the track's progenitor index. Section 2 of [1] (p. 2) treats particle splitting and (n,2n) reactions as branches of a history within the original generation, which then score separately. So a split or an (n,xn) reaction does not start a new generation.

The clock. OpenMC scores the same numerator with a clock instead of track lengths. Each fission is scored with the lifetime of the neutron that caused it, weighted by that fission's progeny in the asymptotic generation. This equals the track-length form in expectation only if the clock runs from the birth of the generation.

• The time a neutron spends before a branch then enters the score of every branch, each weighted by that branch's own progeny.
• Summed over the branches, that time is weighted by the progeny of the neutron before the branch, as in the track-length form.
• A branch whose clock restarts at zero drops the time t_b from the birth of the generation to the branching event. IFP then loses t_b times that branch's progeny weight from the numerator.
• The error is always negative. It is small when few neutrons branch ((n,xn)), and large when most do (weight-window splitting).

The delayed group. The beta_eff numerator weights each neutron born delayed in a generation by its asymptotic progeny. That progeny includes the progeny of every branch of the neutron's history: the neutron that remains after a split, its split copies, and the extra neutrons of an (n,xn) reaction. The birth tag must therefore follow every branch, as the clock does. A branch whose tag is reset to 0 has its progeny counted as the progeny of a prompt neutron, and the numerator loses it.

The non-integral-yield path already behaves this way. Multiplying the weight of the incident neutron is the expected value of banking copies, and the neutron keeps its clock and its tag. With this change the integral-yield path agrees with it.

[1] B. C. Kiedrowski, "Theory, Interface, Verification, Validation, and Performance of the Adjoint-Weighted Point Reactor Kinetics Parameter Calculations in MCNP," LA-UR-10-01700, Los Alamos National Laboratory (2010), https://mcnp.lanl.gov/pdf_files/TechReport_2010_LANL_LA-UR-10-01700_Kiedrowski.pdf

Design

• Field. double lifetime {0.0}; is added at the end of SourceSite, after n_collision. The default of zero keeps the old behavior for every site that does not set it. delayed_group is already a member, with a default of 0.
• Producers.
  • create_secondary() sets bank.lifetime = lifetime() and bank.delayed_group = delayed_group() when a neutron banks a neutron (type.is_neutron() && this->type().is_neutron()). Its only neutron caller is the integral-yield branch of inelastic_scatter(). Multigroup mode has no create_secondary() call, because its scattering multiplicity is a weight factor.
  • split() sets both fields when the particle being split is a neutron. Its only caller is apply_weight_window().
• Neutrons only. Photons, electrons and positrons keep 0 for both fields, as secondaries and when split, and so would a neutron banked by another particle type. Such a particle is not part of a neutron generation, and nothing reads its clock or its tag.
• Consumer. from_source() sets lifetime() = src->lifetime instead of 0. It already set delayed_group() = src->delayed_group. Every path that starts a particle from a site goes through it:
  • source sites;
  • eigenvalue fission sites;
  • the local secondary bank (event_revive_from_secondary());
  • the shared secondary bank;
  • event-based revival.
• Fission neutrons start a new generation. create_fission_sites() is not changed. It builds each site from a default SourceSite, so the clock starts at zero, and the site gets the delayed group sampled for its own emission. This holds in eigenvalue mode, in fixed-source mode with create_fission_neutrons (where the site goes to the secondary bank), and in continuous-energy and multigroup mode alike.
• External sources are unchanged. IndependentSource, file sources and MCPL sources produce sites with a zero clock and a delayed group of 0. So do compiled sources, unless they set the fields.
• MPI. initialize_mpi() lists every field of SourceSite in mpi::source_site, and the new lifetime field is its 16th entry. delayed_group was already in it.
  • A field missing from the datatype would make its extent shorter than sizeof(SourceSite), and arrays of sites sent over MPI would be misread.
  • Sites cross ranks in the eigenvalue bank exchange, the source bank write and MCPL output.
  • In fixed-source mode they also cross ranks in the shared secondary bank, which is on by default with weight windows. That bank carries split neutrons, and with them their clock and tag.
• Python. openmc.lib.core._SourceSite gains the lifetime field, so the structured arrays of openmc.lib.source_bank() and openmc.lib.sample_external_source(as_array=True) keep matching the C++ layout.
• No file format change. The source and surface-source file type is a fixed list of 8 fields (h5banktype()), and only its memory type uses sizeof(SourceSite). The clock is in-memory state. delayed_group is one of the 8 fields already (see "File contents" below).
• Restart is unaffected.
  • Particle restart files carry neither the clock nor the delayed group, before or after this change. Particle restart clears all tallies (particle_restart.cpp), so the clock is not observable there.
  • In eigenvalue mode, the sites that reach a statepoint or source file are fission or external-source sites. Their clock is zero by construction.
• Documentation. Two sentences in docs/source/usersguide/kinetics.rst say when the IFP lifetime clock starts and that such a neutron also keeps the delayed group of the neutron that produced it.

No IFP code is changed (ifp.cpp, ifp.h, bank.cpp, eigenvalue.cpp, tally_scoring.cpp).

The (n,xn) case

split() must copy the tag: a split only divides weight. For create_secondary(), this PR copies it as well, for three reasons.

1. The tag is a generation label. The IFP tag records how the neutron that started this line of the generation was born, not how this particular neutron was emitted. An (n,2n) neutron produced by a delayed neutron is emitted promptly. It still continues the delayed neutron's generation, and its progeny are part of that delayed neutron's importance. The same reasoning keeps the clock running.
2. The two yield representations agree. A non-integral yield keeps the tag on the weight-multiplied neutron. Banking integral-yield copies with a tag of 0 would give two representations of one reaction different beta_eff.
3. One rule for both fields. Every neutron that continues a generation keeps both the clock and the tag. That is simpler to state and to review.

The alternative is to leave (n,xn) neutrons at 0, on the ground that an (n,xn) neutron is a prompt neutron. We did not take it, for the reasons above.

In practice the choice matters little. In the local ENDF/B-VII.1 data, the delayed-neutron spectra of U-235 and U-238 end at 3.0 MeV. That is below the (n,2n) thresholds of U-235 (5.32 MeV) and U-238 (6.18 MeV), so a delayed neutron cannot cause an (n,2n) reaction in either. Be-9, with a threshold of about 1.75 MeV, is the exception that a delayed spectrum's tail can reach. The C++ test pins the copy. Leaving (n,xn) neutrons at 0 would be a one-line change in create_secondary() and in the test.

Readers of the particle's delayed group

This PR changes the value of p.delayed_group() only for split and (n,xn) neutrons. These are all the code paths that read it:

reader	effect for a split or (n,xn) neutron
ifp() (src/ifp.cpp), scored by ifp-beta-numerator with or without a DelayedGroupFilter	the producing neutron's group instead of 0: the fix
add_surf_source_to_bank() (src/particle.cpp), surface-source sites	a site written by such a neutron records the producing neutron's group instead of 0
collision tracking (src/collision_track.cpp)	the same, for collision-track sites
print_particle() (src/output.cpp)	debug output

Not affected:

• delayed-nu-fission scores and the DelayedGroupFilter on them read the group of the sites created at the collision (p.nu_bank()), not the particle's tag;
• fission sites set their own group;
• the time cutoff and the delayed-neutron counts in create_fission_sites() use the created site;
• external sources set 0;
• sort_bank() orders by parent and progeny IDs only;
• the multigroup create_fission_sites() reads only the value it has just set (see "Multigroup mode" below);
• the particle restart file does not store the tag;
• the MPI datatype and _SourceSite already carry delayed_group, and the HDF5 source and surface-source types are unchanged.

File contents. Surface-source files and collision-track files keep their formats. Their content changes only for sites written by split or (n,xn) neutrons, which now carry the delayed group of the neutron that produced them instead of 0. In the surface_source_write, surface_source and collision_track regression tests, the generated files are unchanged (see "Regression references" below).

Multigroup mode (behavior before this PR, not changed). The multigroup create_fission_sites() (src/physics_mg.cpp) also sets the colliding particle's delayed_group() to the group of each site it creates. After a multigroup fission collision, a particle's delayed_group() is therefore no longer its birth tag, and a later split copies that value, the same as the parent's. IFP does not read it there: the multigroup create_fission_sites() does not call ifp(). Only the surface-source, collision-track and print readers above see it.

Tests

C++ unit test, tests/cpp_unit_tests/test_particle.cpp (new)

3 test cases and 55 assertions, all exact. The neutron has lifetime() = 1.25e-6 s and time() = 3.0e-6 s:

• A default SourceSite has lifetime == 0.0 and delayed_group == 0.
• Lifetime:
  • create_secondary(..., neutron) and split() bank sites with lifetime == 1.25e-6 and time == 3.0e-6;
  • a photon secondary, split() of a photon, and a neutron banked by a photon bank lifetime == 0.0 and keep time;
  • from_source() of each banked site gives p.lifetime() == site.lifetime.
• Delayed group, for a neutron with delayed_group() = 4:
  • the (n,xn) and split copies bank delayed_group == 4;
  • a photon secondary, split() of a photon, and a neutron banked by a photon bank 0;
  • from_source() and event_revive_from_secondary() of each banked site restore p.delayed_group() == site.delayed_group.

The lifetime test cases do not compile on develop, which has no SourceSite::lifetime. The delayed-group test case fails on the lifetime-only part of this PR, with bank[i].delayed_group == 4 evaluating to 0 == 4.

C++ unit test, tests/cpp_unit_tests/test_secondary_bank.cpp (new)

The test runs with the settings of a fixed-source run with a shared secondary bank.

• It banks split and (n,2n) neutrons, photons, electrons and positrons, plus a site that is not a secondary, such as a fission site.
• It redistributes the bank with synchronize_global_secondary_bank() and revives every site with event_revive_from_secondary().
• A revived neutron must resume the lifetime of the neutron that banked it. The other particles and the non-secondary site must start at zero.

With MPI, the test also checks that the extent of the SourceSite MPI datatype equals sizeof(SourceSite). The CTest entry test_secondary_bank_mpi runs it on two ranks, so that sites banked on the first rank are sent to the second. This is the direct test of the MPI layout.

The CTest entry passes --order decl, because since version 3.9 Catch2 runs test cases in random order with a seed drawn by each process, and the two ranks must run the two test cases that call synchronize_global_secondary_bank() in the same order. With the Catch2 3.16 that CI fetches, 10 of 20 local runs of test_secondary_bank_mpi failed without the option, as the MPI CI job of this PR did, and 20 of 20 passed with it, with Catch2 3.16 and with 3.4.

A second test case follows the delayed group through the same path.

• Neutrons born from precursors of delayed groups 1 to 4 each bank a split neutron, an (n,2n) neutron, a photon, an electron and a positron.
• A photon that was given a nonzero delayed group, which transport never gives a photon, banks a split photon and a neutron.
• After synchronize_global_secondary_bank() (rank 0 to rank 1 under test_secondary_bank_mpi) and event_revive_from_secondary():
  • a neutron banked by a neutron must resume that neutron's delayed group;
  • every other banked particle must have none;
  • a site that is not a secondary, such as a fission site, keeps its own group.

On the lifetime-only part of this PR this case fails, serially and on 2 ranks, while the lifetime cases pass.

Python unit test, tests/unit_tests/test_secondary_delayed_group.py (new)

This test checks the delayed group in transport.

• Model. Every history starts from a source-file site of delayed group 3 in a fixed-source run with photon transport. The geometry is a Be-9 sphere (r = 9 cm) with a Fe-56 shell out to 10 cm. Nothing fissions, so every neutron continues the history of a source neutron.
• Cases. Four, each with a local and with a shared secondary bank:
  • n2n: 14 MeV source neutrons without weight windows, so the only neutrons banked during transport are the extra neutrons of (n,2n) reactions;
  • split: 1 MeV source neutrons, below the (n,2n) thresholds, with weight windows that split the neutrons colliding outside the central mesh cell into 10.
• Check. It reads the particles that leave the sphere from the surface source file. Every neutron must have delayed group 3 and every photon none.
• Modes. The test follows the --event and --mpi options of the test suite.

The test passes with NNDC and ENDF/B-VII.1 data in history mode, with --event, with --mpi (2 ranks) and with --mpi --event. On the lifetime-only part of this PR, all four cases fail in every mode. With NNDC data in history mode, these neutrons had delayed group 0:

case	local bank	shared bank
n2n	640 of 1581	589 of 1524
split	8116 of 9099	8187 of 9172

Python unit test, tests/unit_tests/test_ifp_weight_windows.py (new)

The test runs the U-235 sphere of the ifp regression tests twice, without and with weight windows. It asserts |Lambda_ww / Lambda_analog - 1| <= 0.025. Settings: 50000 particles, 30 batches (10 inactive), ifp_n_generation = 1.

Weight windows.

• OpenMC scales a weight window by the ratio of the particle's birth weight to the centre of its birth window (apply_weight_window()). So a window that is uniform in space never splits a weight-1 fission neutron, whatever its bounds.
• A neutron born where no window applies takes 1 as its reference weight instead (apply_weight_windows()), and the absolute bounds then apply.
• The test uses this. On a 20 x 20 x 20 mesh over the sphere, every fourth cell, (i + j + k) % 4 == 0, has the window [0.05, 0.1], and the other cells have none (-1).
• Most fission neutrons are born where no window applies. At their first collision in a cell with a window, each is split into 10.

Single-seed results, NNDC data, 2 threads, measured on develop 3cded0f before the rebase:

Build	Lambda analog (s)	Lambda with weight windows (s)	Ratio	Test
develop 3cded0f	6.3008e-09	5.4012e-09	0.857	fails
this PR	6.3033e-09	6.2716e-09	0.995	passes

With ENDF/B-VII.1 data the ratios are 0.869 (develop, fails) and 1.0065 (this PR, passes). At OMP_NUM_THREADS=2 the test takes 8 s to 10 s in history mode and 17 s to 23 s with --event. It passes on the final branch with NNDC data, in both modes.

Tolerance. The tolerance comes from a 5-seed study of this exact configuration (ENDF/B-VII.1 data, 5 seeds per arm, Welch t on the per-seed Lambda):

	develop	this PR
mean ratio, weight windows / analog	0.866 (bias -13.4 %)	1.0041
Welch t	-104	3.19 (nu = 7.2, t_0.995 = 3.46)
standard deviation of the ratio (5 seeds)	0.18 %	0.19 %

• 15 further seeds were run with this PR as a check only. Over the 20 seeds, the mean ratio is 0.9996 (Welch t = -0.31), and the standard deviation of the ratio is 0.53 %. The tolerance rests on this estimate.
• The tolerance of 0.025 is 4.7 standard deviations, and just under one fifth of the bias on develop (0.0267).
• With ifp_n_generation = 5, the same study gives -13.6 % on develop (t = -35.1) and +0.17 % with this PR (t = 0.41, t_0.995 = 3.57).
• The test uses one generation because the scatter of the ratio is several times smaller there. With five generations, no test that fits in 60 s keeps the tolerance both above 4 standard deviations and below one fifth of the bias.

No beta_eff assertion. The same two runs give beta_eff, but at this size its per-seed scatter is too large for a test. Over 5 seeds the standard deviation of the beta_eff ratio with this PR is 1.8 %. A tolerance of 4 standard deviations (0.072) would exceed one fifth of the bias that it must detect (0.048), so the test asserts only Lambda. The delayed group is covered by the exact C++ tests and by the paired runs below.

beta_eff with weight windows

The model, weight windows and settings are those of the Python test, with ENDF/B-VII.1 data, 5 seeds per arm, 2 threads and a Welch t on the per-seed beta_eff. "Before" is the lifetime-only part of this PR, the first two commits. The lifetime part does not change beta_eff: in the paired runs of the lifetime study, ifp-beta-numerator agrees with develop to 4.5e-16 relative.

ifp_n_generation		before	this PR
1	beta_eff with weight windows / analog, means	-23.9 %	-0.58 %
1	Welch t (nu, t_0.995)	-15.4 (7.9, 3.37)	-0.39 (7.7, 3.40)
1	standard deviation of the per-seed ratio	0.88 %	1.8 %
5	beta_eff with weight windows / analog, means	-25.0 %	-2.9 %
5	Welch t (nu, t_0.995)	-7.17 (8.0, 3.36)	-0.79 (8.0, 3.36)

With one generation and a DelayedGroupFilter on groups 1 to 6, every group is 19 % to 27 % low before (Welch t from -2.5 to -13.1), and within |t| <= 1.1 with this PR. With five generations every group is within |t| <= 1.91 with this PR.

Python unit test, tests/unit_tests/test_ifp_n2n_lifetime.py (new)

This test checks the lifetime clock of (n,2n) neutrons in transport. The effect of (n,xn) reactions on the generation time is too small for a statistical comparison (see "Size of the bias" above), so the test checks an exact lower bound instead.

• Model. An eigenvalue run with 5000 particles, 2 batches (1 inactive) and ifp_n_generation = 1. Every source neutron starts at the origin at t = 1 s with 14 MeV, moving along +z. It crosses 100 m of void before it reaches a 5 cm Be-9 slab (1.85 g/cm3), followed by a 10 cm U-235 slab (16.0 g/cm3). All cells are 30 cm wide in x and y, and every outer surface is a vacuum boundary. One tally without filters scores ifp-time-numerator (N) and ifp-denominator (D).
• Why the bound is exact.
  • Nothing is split (there are no weight windows), and fission neutrons start the next generation. Every neutron of the first generation is therefore a source neutron or an extra neutron of an (n,xn) reaction with an integral yield, which continues the history of the neutron that produced it.
  • The clock of a source neutron reaches the flight time L / v, with L = 100 m and v its speed in transport, before its first collision. If the extra neutrons keep the clock of the neutron that produced them, every fission site of the first generation records at least that time.
  • With one IFP generation, N/D in the second generation is the mean of those recorded times over positive weights, so it is at least L / v too, whatever the random numbers.
  • The bound T_min = L / (c sqrt(2 E / (m c^2))) uses the non-relativistic speed, with the speed of light and neutron mass of constants.h. For every E > 0 that speed exceeds the relativistic speed that OpenMC uses in transport, c sqrt(E (E + 2 m c^2)) / (E + m c^2). T_min = 1.9323e-6 s is therefore below OpenMC's flight time of 1.9538e-6 s, which is 1.1 % longer.
• Check. D > 0, N/D >= T_min and N/D <= 1 ms. The upper bound fails for a clock that copies the absolute time time() instead of the lifetime, because the source neutrons start at t = 1 s.
• Modes. The test follows the --event and --mpi options of the test suite.

If an extra neutron restarts its clock at the reaction, in the target, its fission sites record only the time since that reaction, and N/D falls below the bound. The void is convex, so a neutron that leaves the target cannot come back and add flight time.

Results, NNDC data, 2 threads, N/D / T_min. Seed 1 is the test's own run; seeds 2 to 5 run the same model with another seed.

Build	seed 1	seed 2	seed 3	seed 4	seed 5
this PR	1.0195	1.0160	1.0189	1.0168	1.0168
develop a5bc348 (same code as the base, 4e69968)	0.8213	0.8172	0.8259	0.8332	0.8319
this PR with bank.lifetime = lifetime(); removed from create_secondary()	0.8213	0.8172	0.8259	0.8332	0.8319

• With this PR the three checks hold for every seed, and the test passes in history mode, with --event and with --mpi (2 ranks).
• On develop and with the clock copy removed, the test fails in history mode, with --event and with --mpi (2 ranks), with IFP mean lifetime (1.5870e-06 s) is below the flight time of the source neutrons (1.9323e-06 s).
• The random walk is the same with all three builds. At 1 thread (seed 1 and six other seeds), k of both generations, ifp-denominator and the statepoint source banks are bitwise equal for all three, and N is bitwise equal between develop and the build without the clock copy: only the time numerator of this PR differs. At 2 threads, in history mode, ifp-denominator and the source banks are still bitwise equal, while k, and N between develop and the build without the clock copy, differ only in the last bits, by up to about 2e-15 between builds as between two runs of the same build, from the order of the OpenMP reductions.
• On the lifetime-only part of this PR, the first two commits, create_secondary() already copies the clock, and the test passes.

At OMP_NUM_THREADS=2 the test takes 0.11 s in history mode, 0.13 s with --event and 0.5 s with --mpi. The --mpi runs used Open MPI 4.1.6 with process binding off, through a small wrapper around mpiexec that unsets all OMPI_*, PMIX_*, OPAL_* and ORTE_* environment variables but the binding policy (OMPI_MCA_hwloc_base_binding_policy) before it calls mpiexec.

On the sixth commit, ctest gives the counts listed under "Tests run on the final branch", serially and with MPI, and ifp/total, ifp/groupwise, test_ifp_weight_windows.py and test_secondary_delayed_group.py pass in history mode.

Regression references, ifp/total and ifp/groupwise

Both models are a bare U-235 sphere. Only the ifp-time-numerator rows of both results_true.dat files change. They were regenerated with the NNDC library and a strict-FP build. The develop build on the same machine reproduces the old references.

ifp-time-numerator	before	after	relative change
sum	9.109384E-08	9.110594E-08	+1.33e-4
sum of squares	5.667165E-16	5.668788E-16	+2.86e-4

All other values in both files are unchanged, including k-effective, ifp-beta-numerator and ifp-denominator, and inputs_true.dat does not change.

The delayed-group change leaves every reference unchanged. These models have no weight windows, so nothing is split, and no U-235 delayed neutron reaches the (n,2n) threshold. The lifetime-only part and the full PR were run with --update on ENDF/B-VII.1 data, at OMP_NUM_THREADS=2 except where a test sets its own thread count. The generated files are identical for the tests below; in some HDF5 files only the version, commit and date attributes differ.

• ifp;
• weightwindows, weightwindows_fw_cadis, weightwindows_fw_cadis_local, weightwindows_fw_cadis_mesh and weightwindows_pulse_height;
• surface_source_write and mg_fixed_source_ww_fission_shared_secondary;
• surface_source and collision_track, run at 1 thread to remove thread-scheduling differences.

The one exception is surface_source_write/case-a01, which runs at 2 threads by design. Its 24 sites differ in order only. They are equal after sorting, which is how the test compares them.

Paired checks before and after

Each pair uses the same seed at 2 threads, and both builds are RelWithDebInfo with strict FP. The random walk is the same in both runs of a pair, so a paired difference is not statistical. The lifetime checks were run with develop 3cded0f as "before" and the lifetime commit as "after", before the rebase. The delayed-group checks were run on the rebased branch, with the first two commits as "before" and the full PR as "after".

• Nothing else changes (lifetime).
  • Every regression test that this build can run was run with --update, with both builds, in two copies of tests/ (ENDF/B-VII.1 data). Only the two ifp-time-numerator rows above differ for a reason other than run-to-run variation.
  • Every other file that differed also differed between two runs of develop at 2 threads, because it depends on thread scheduling: collision_track, the depletion tests, microxs, mgxs_library_hdf5, surface_source, surface_source_write, sourcepoint_latest and random_ray_diagonal_stabilization.
  • Nine of those test directories were rerun at one thread. They are identical between the builds, apart from the wall-clock depletion times and a spreadsheet timestamp.
  • With NNDC data, on the first machine, each build ran the full regression suite against its own committed references. 303 tests pass on both. The same 7 fail or error on both:
    • cpp_driver, source_dlopen and source_parameterized_dlopen need an installed OpenMC;
    • weightwindows_fw_cadis, weightwindows_fw_cadis_local and weightwindows_fw_cadis_mesh (flat, linear) fail on that machine only. On the machine of the full test suite under "Tests run on the final branch", they pass with both builds, and the regression part there has 307 passes.
• ifp/total statepoints (ENDF/B-VII.1).
  • k-combined is equal to 2.8e-15 relative.
  • ifp-denominator and ifp-beta-numerator are bitwise equal, and so is the source bank.
  • ifp-time-numerator increases by 1.89e-4.
• (n,xn)-rich model. U-235 and Be-9 (atom ratio 1:2, 10 g/cm3, r = 10 cm, k = 0.770), ENDF/B-VII.1 data, 1e5 particles, 30 batches (10 inactive), 5 IFP generations:
  • there are 0.0295 (n,2n) and 1.4e-6 (n,3n) reactions per source neutron;
  • k-combined is equal to 7.9e-14 relative;
  • ifp-denominator, ifp-beta-numerator and the source bank are bitwise equal;
  • the generation time increases by 3.23e-3 relative. That is a paired difference; the statistical standard deviation of Lambda in one run is 0.45 %.
• Integral versus non-integral yield.
  • A copy of the Be-9 data has its (n,2n) yield changed from 2 to 2(1 + 1e-9), so that the reaction takes the weight branch.
  • With this PR the generation times from the two libraries agree: the weight branch is +0.09 % from the integral branch (Welch t = 0.25 over 5 seeds, t_0.995 = 3.39).
  • On develop the integral branch is low. 5 seeds cannot resolve the difference (+0.41 %, t = 1.11), but the paired run above measures it: 0.32 %.
• Weight-window study, lifetime, paired. In all 20 paired runs of the 5-seed study (with and without weight windows, 1 and 5 generations):
  • k agrees between the builds to 6.8e-14 relative;
  • ifp-denominator and ifp-beta-numerator agree to 4.5e-16;
  • the source banks are bitwise equal;
  • the generation time is higher with this PR.
• Weight-window study, delayed group, paired. The 20 paired runs of the beta_eff study above:
  • with weight windows, k agrees to 1.4e-13 relative, ifp-denominator to 5.6e-16 and ifp-time-numerator to 2.9e-15;
  • ifp-beta-numerator is higher with this PR in every active batch of every run, by 27 % to 33 % in total;
  • without weight windows, ifp-beta-numerator and ifp-denominator are bitwise equal, k agrees to 8.2e-14 and ifp-time-numerator to 7.9e-15, so no delayed neutron caused an (n,xn) reaction in these runs;
  • the final source banks are bitwise equal in every pair.

Tests run on the final branch

These were run locally (RelWithDebInfo, strict FP, OMP_NUM_THREADS=2). The last four commits change only tests and a comment, so the product code is that of the third commit. The fourth commit compiles to the same code as the fifth, sixth and seventh, and these ran on it: ctest, the new Python test, the ifp regression tests, the unit modules named below and the full test suite. The sixth commit adds test_ifp_n2n_lifetime.py, whose runs are given in its own subsection above, and the seventh passes --order decl to test_secondary_bank_mpi (see the test_secondary_bank.cpp subsection). The MPI regression tests and the full set of unit modules ran on the third commit:

• ctest. 16 tests serially, of which 15 pass and test_mcpl_stat_sum is skipped (no MCPL). With MPI, 17 tests, including test_secondary_bank_mpi; 16 pass and the same one is skipped. test_secondary_bank passes on 1 and 2 ranks, with both of its delayed-group and lifetime cases.
• IFP regression tests, NNDC data. ifp/total and ifp/groupwise pass against the committed references, together with test_ifp_weight_windows.py, in history mode and with --event.
• New delayed-group test. test_secondary_delayed_group.py passes with NNDC and ENDF/B-VII.1 data, in history mode, with --event, with --mpi (2 ranks) and with --mpi --event. test_ifp_weight_windows.py, test_surface_source_write.py and test_source_file.py pass with it (25 passed, 2 skipped).
• MPI regression tests. A build with -DOPENMC_USE_MPI=on passes mg_fixed_source_ww_fission_shared_secondary, mg_basic and mg_basic_delayed with --mpi at 2 ranks, against their committed references. The first is a fixed-source run with weight windows, fission neutrons and the shared secondary bank, so split neutrons cross ranks there.
• Weight-window, surface-source and collision-track regression tests. These pass against the committed NNDC references in the full test suite below. They also ran as the paired --update comparison on ENDF/B-VII.1 data described under "Regression references".
• Unit tests, NNDC data, on the third commit. 200 passed and 24 skipped (DAGMC, libMesh and XDG are not built). The modules run were:
  • test_ifp.py and test_ifp_weight_windows.py;
  • test_collision_filter_ww.py and weightwindows/;
  • test_lib.py;
  • the test_source*.py modules and test_surface_source_write.py;
  • test_uniform_source_sampling.py and test_deplete_external_source_rates.py.
• Struct layout. sizeof(SourceSite) is 144 bytes, and ctypes.sizeof(openmc.lib.core._SourceSite) is also 144, with lifetime at offset 136. On develop both are 136.
• Style. clang-format 18 reports no changes for the 5 changed C++ files. pyflakes and pycodestyle report nothing for the three new Python tests or for the changed lines of openmc/lib/core.py.
• Full test suite, NNDC data.
  • What ran. CI's targets: tests/test_matplotlib_import.py, tests/unit_tests and tests/regression_tests. History mode, OMP_NUM_THREADS=2, in 4 parallel shards, on Ubuntu 24.04 with an i7-8700K.
  • This branch: 1532 passed, 141 skipped, 64 errors and no failures.
  • develop 4e69968, the same way: 1527 passed, 141 skipped, the same 64 errors and no failures. The 5 extra passes are the new tests.
  • The errors are the same on both:
    • cpp_driver, source_dlopen and source_parameterized_dlopen need an installed OpenMC;
    • 61 data tests need OPENMC_ENDF_DATA (ENDF-6 files), which is not set here.
  • The skips are features not built or not installed: DAGMC, libMesh, XDG, MCPL, NCrystal, UWUW and NJOY. Three more checks are skipped because they apply to event-based mode.
  • Subsets. The IFP, weight-window, surface-source, collision-track and shared-secondary regression tests and the related unit tests, the new tests included:
    • with --event: 119 passed (develop: 114);
    • with --mpi on 2 ranks (OpenMPI): 134 passed (develop: 133).
  • One caveat under OpenMPI. In that same session the 4 cases of test_secondary_delayed_group.py failed, because they ran after tests that initialize openmc.lib in-process. Run alone with --mpi, they pass. Upstream's test_tracks.py fails the same way on develop when it runs after test_ww_gen.py.

Not run:

• tests that need DAGMC, libMesh, XDG, NCrystal, MCPL, UWUW or NJOY;
• the 61 tests that need ENDF-6 data (OPENMC_ENDF_DATA);
• the three tests that compile against an installed OpenMC (cpp_driver, source_dlopen, source_parameterized_dlopen), because nothing was installed.

Compatibility

• Results. Only the IFP kinetics parameters change, and only upward.
  • The generation time (ifp-time-numerator) changes in eigenvalue models with integral-yield (n,xn) reactions or with weight-window splitting: by about 1e-4 relative for bare U-235, and by much more with splitting.
  • beta_eff (ifp-beta-numerator) changes in eigenvalue models with weight-window splitting. A change from (n,xn) reactions needs delayed neutrons above an (n,2n) threshold, as in beryllium.
  • The random walk depends on neither field, and k-effective does not change.
• File contents. Surface-source and collision-track sites written by split or (n,xn) neutrons carry the producing neutron's delayed group instead of 0.
• SourceSite layout. The struct grows from 136 to 144 bytes.
  • The MPI datatype and the ctypes structure in openmc.lib are updated with it.
  • Compiled custom sources and C++ drivers that use SourceSite must be rebuilt against the new headers, as after Preserve collision count when splitting particles #3958.
• Files. The source, surface-source, collision-track, statepoint and particle restart file formats do not change.
• Python API. The only change is the new lifetime field, in _SourceSite and in the structured arrays returned by openmc.lib.source_bank() and openmc.lib.sample_external_source(as_array=True).

Related issues and PRs

• Preserve collision count when splitting particles #3958 (merged, "Preserve collision count when splitting particles") is the precedent.
  • It added n_collision to SourceSite, its MPI datatype and _SourceSite in the same way, so that split particles keep their collision count.
  • It added test_collision_filter_ww.py, an analog-versus-weight-window comparison like the new IFP test.
• Introducing a functionality to reset the time for produced particles (neutrons) #2355 (closed in 2023, "Introducing a functionality to reset the time for produced particles (neutrons)") proposed a setting that resets time() for neutrons produced by (n,xn) reactions and fission. This PR does not change time(). It keeps the IFP clock running through (n,xn) reactions and splits and restarts it at fission, as the adjoint-weighted definition requires.
• Replace the IFP parameter enum with independent flags #4094 (merged, in the base, "Replace the IFP parameter enum with independent flags"). This change does not interact with it.
• Bundle each IFP data stream with the flag and banks that govern it #4097 (open, "Bundle each IFP data stream with the flag and banks that govern it") refactors how IFP stores its data.
  • It shares no file with this PR.
  • ifp() still stores p.lifetime() and p.delayed_group() there.
  • It changes no reference file, so the references regenerated here should stand after it.
• Support pulse-height tallies in shared secondary mode #4089 (open, "Support pulse-height tallies in shared secondary mode") changes code next to the lines changed here.
  • It renames SourceSite::parent_id, adds members after n_collision, and edits create_secondary(), split(), initialize_mpi() and _SourceSite.
  • The two changes are independent. Whichever is merged second needs a small textual merge in SourceSite, create_secondary() and split().
• Delayed neutron fraction by isotope with IFP method #3659 (open, "Delayed neutron fraction by isotope with IFP method") adds a regression test, ifp/nuclidewise.
  • The test is a U-235/Pu-239 sphere whose reference includes an ifp-time-numerator tally (tally 1).
  • Whichever of the two PRs is merged second must regenerate those two rows.
  • Both PRs also edit kinetics.rst and particle_data.h, in different places.
• Split storage optimization #3656 (draft, "Split storage optimization") replaces split() with replicate().
  • If it goes ahead, the replicated sites must carry the clock and the delayed group too.
  • The new C++ tests would need the new call, and the new Python test would catch a clock lost there.

AI Assistance

1. Code, tests, analysis and this text:
   • Harness: Claude Code
   • Model: Claude Opus 5.5 (claude-opus-5-5)
   • Reasoning effort: xhigh
2. Review of the change:
   • Harness: OpenCode
   • Model: grok-4.7
   • Reasoning effort: not recorded
3. Review of the PR formalities:
   • Harness: OpenCode
   • Model: muse-spark
   • Reasoning effort: not recorded

Checklist

• I have performed a self-review of my own code
• I have run clang-format (version 18) on any C++ source files (if applicable)
• I have followed the style guidelines for Python source files (if applicable)
• I have made corresponding changes to the documentation (if applicable)
• I have added tests that prove my fix is effective or that my feature works (if applicable)

Summary by CodeRabbit

• New Features
  • IFP lifetime scoring now preserves neutron travel time through (n,xn) reactions and weight-window splitting, improving generation-time estimates. Fission neutrons start with a new lifetime clock.
  • Delayed-group data is retained by secondary neutrons from neutron parents, including during splitting; non-neutron secondaries receive group zero.
  • Documentation now describes how lifetime and delayed-group data carry across generations.

[coderabbitai]

No actionable comments were generated in the recent review. 🎉

ℹ️ Recent review info

⚙️ Run configuration

• Configuration used: defaults
• Review profile: CHILL
• Plan: Advanced
• Run ID: dbc79306-ca93-4743-9101-ded814244b3e

📥 Commits

Reviewing files that changed from the base of the PR and between ba4b30a and b9931a7.

⛔ Files ignored due to path filters (2)

• tests/regression_tests/ifp/groupwise/results_true.dat is excluded by !**/*.dat
• tests/regression_tests/ifp/total/results_true.dat is excluded by !**/*.dat

📒 Files selected for processing (11)

• docs/source/usersguide/kinetics.rst
• include/openmc/particle_data.h
• openmc/lib/core.py
• src/initialize.cpp
• src/particle.cpp
• tests/cpp_unit_tests/CMakeLists.txt
• tests/cpp_unit_tests/test_particle.cpp
• tests/cpp_unit_tests/test_secondary_bank.cpp
• tests/unit_tests/test_ifp_n2n_lifetime.py
• tests/unit_tests/test_ifp_weight_windows.py
• tests/unit_tests/test_secondary_delayed_group.py

Included review availability: This review used your included allowance. Your plan provides up to 10 included reviews per hour; 9 remain after this review.

---

📝 Walkthrough

Walkthrough

Source sites now carry neutron lifetime. Neutron secondary creation and weight-window splitting preserve lifetime and delayed-group state through banking and restoration. The changes also update the MPI datatype and add C++ and transport tests for these paths and IFP scoring.

Changes

Neutron state propagation

Layer / File(s)	Summary
Source-site lifetime layout include/openmc/particle_data.h, openmc/lib/core.py, src/initialize.cpp	SourceSite and its ctypes representation include lifetime. The MPI datatype registers the field as its sixteenth entry.
Particle state propagation src/particle.cpp, docs/source/usersguide/kinetics.rst	Neutron secondaries and splits copy lifetime and delayed group. from_source restores lifetime from the source site. The guide describes clock continuation for (n,xn) reactions and weight-window splitting, and states that fission neutrons start at zero.
C++ particle and bank validation tests/cpp_unit_tests/test_particle.cpp, tests/cpp_unit_tests/test_secondary_bank.cpp, tests/cpp_unit_tests/CMakeLists.txt	C++ tests check lifetime and delayed-group values through particle creation, splitting, source restoration, and secondary-bank revival. The test configuration registers both test programs and adds an MPI run for test_secondary_bank.
Transport and IFP regression tests tests/unit_tests/test_secondary_delayed_group.py, tests/unit_tests/test_ifp_n2n_lifetime.py, tests/unit_tests/test_ifp_weight_windows.py	Transport tests check delayed groups on recorded neutrons and photons. IFP tests check lifetime after (n,2n) reactions and compare generation time with and without neutron weight windows.

Priority: ➖ Normal

Estimated code review effort: 3 (Moderate) | ~25 minutes

Change: Bug fix · Severity of issue fixed: Medium

Merge Risk: ⚪ Minimal · up to b9931

The change is mergeable after normal checks; no concrete issue requiring a pre-merge fix was identified.

Security Architecture Review

Security architecture risk: 🔵 Low · up to b9931

The updated representations agree, and the inspected banking paths preserve generation state while fresh fission particles reset it. Compatibility with independently compiled source-record consumers remains uncertain because the binary layout grows. No introduced privilege-boundary bypass was established.

Retained concerns
No architecture-level concerns identified.

Security review details

Security Blast Radius

• inferred — The demonstrated exposure is the simulation process, its in-process source consumers and participating MPI ranks. An incompatible native buffer could affect process memory, but no tenant, credential or service-authority expansion was established.

Trust Boundaries and Controls

• observed — The existing C source-bank accessor exposes mutable native records, and sampling writes native SourceSite elements into caller-provided storage. The PR updates the known Python consumer to the matching layout; these interfaces are documented as experimental rather than as a stable binary contract.

Resilience and Maintainability Implications

• inferred — A caller allocating records with the pre-change layout would provide insufficient storage for sampling writes made by the new library. The inspected API cannot detect that mismatch because it accepts a void pointer without element-size negotiation. Whether such stale callers are supported or deployed remains unknown; this is a conditional compatibility hazard, not a demonstrated attack.

Hardening Proposals

• proposed — If independently compiled consumers must remain supported, document coordinated rebuild and rollback requirements, or introduce a record-layout compatibility check before accepting their buffers.

🚥 Pre-merge checks | ✅ 4 | ❌ 1

❌ Failed checks (1 warning)

Check name	Status	Explanation	Resolution
Docstring Coverage		Docstring coverage is 28.57% which is insufficient. The required threshold is 80.00%. Docstring coverage is scoped to functions touched by this diff. Analyzed 21 functions across 9 files. (2 skipped: …	Write docstrings for the functions missing them to satisfy the coverage threshold.

✅ Passed checks (4 passed)

Check name	Status	Explanation
Linked Issues check		The PR meets the coding requirements in [#4177]. SourceSite stores lifetime, create_secondary() and split() copy lifetime and delayed_group for neutron-to-neutron branching, and `from_sour…
Out of Scope Changes check		The changes stay within [#4177]. The documentation, SourceSite layout updates, MPI and Python bindings, regression tests, and CTest ordering support state propagation, compatibility, or verification o…
Title check		The title clearly and concisely identifies the main change: preserving the lifetime clock and delayed group for (n,xn) and split neutrons.
Description check		The description follows the repository template. It explains the problem, motivation, design, compatibility impact, tests, AI assistance, linked issue, and completed checklist items.

Full details: Docstring Coverage

Explanation

Docstring coverage is 28.57% which is insufficient. The required threshold is 80.00%. Docstring coverage is scoped to functions touched by this diff. Analyzed 21 functions across 9 files. (2 skipped: 2 unsupported.)

• Fix all pre-merge checks with AI

✨ Finishing Touches

🧪 Generate unit tests (beta)

• Create a new PR

• Autopilot · Keep fixing CodeRabbit findings and required CI, and resolving merge conflicts

---

_{Comment @coderabbitai help to get the list of available commands.}