From fd474019699b420c533d87fb4c04114db78f1941 Mon Sep 17 00:00:00 2001 From: GuySten Date: Tue, 26 Aug 2025 03:44:15 -0700 Subject: [PATCH 1/9] initial implementation wip --- include/openmc/particle.h | 3 +- include/openmc/particle_data.h | 3 ++ include/openmc/tallies/filter.h | 1 + include/openmc/tallies/filter_particle.h | 40 ++++++++++++++ include/openmc/tallies/tally.h | 1 + include/openmc/tallies/tally_scoring.h | 6 ++- openmc/filter.py | 29 ++++++++-- openmc/lib/filter.py | 10 ++-- src/mesh.cpp | 2 +- src/particle.cpp | 21 ++++++-- src/particle_restart.cpp | 1 + src/random_ray/random_ray.cpp | 2 +- src/simulation.cpp | 5 +- src/tallies/filter.cpp | 2 + src/tallies/filter_particle.cpp | 69 +++++++++++++++++++++++- src/tallies/tally.cpp | 11 +++- src/tallies/tally_scoring.cpp | 8 +-- 17 files changed, 190 insertions(+), 24 deletions(-) diff --git a/include/openmc/particle.h b/include/openmc/particle.h index 0f37719b94f..cdfbfbb6b82 100644 --- a/include/openmc/particle.h +++ b/include/openmc/particle.h @@ -62,7 +62,8 @@ class Particle : public ParticleData { //! site may have been produced from an external source, from fission, or //! simply as a secondary particle. //! \param src Source site data - void from_source(const SourceSite* src); + //! \param particle ParticleType incident particle type + void from_source(const SourceSite* src, ParticleType particle); // Coarse-grained particle events void event_calculate_xs(); diff --git a/include/openmc/particle_data.h b/include/openmc/particle_data.h index ec393845f81..235fbf59ac7 100644 --- a/include/openmc/particle_data.h +++ b/include/openmc/particle_data.h @@ -469,6 +469,7 @@ class ParticleData : public GeometryState { CacheDataMG mg_xs_cache_; ParticleType type_ {ParticleType::neutron}; + ParticleType type_last_ {ParticleType::neutron}; double E_; double E_last_; @@ -566,6 +567,8 @@ class ParticleData : public GeometryState { // Particle type (n, p, e, gamma, etc) ParticleType& type() { return type_; } const ParticleType& type() const { return type_; } + ParticleType& type_last() { return type_last_; } + const ParticleType& type_last() const { return type_last_; } // Current particle energy, energy before collision, // and corresponding multigroup group indices. Energy diff --git a/include/openmc/tallies/filter.h b/include/openmc/tallies/filter.h index 65098597a5d..6fd9d3bf8d5 100644 --- a/include/openmc/tallies/filter.h +++ b/include/openmc/tallies/filter.h @@ -39,6 +39,7 @@ enum class FilterType { MUSURFACE, PARENT_NUCLIDE, PARTICLE, + PARTICLE_OUT, POLAR, SPHERICAL_HARMONICS, SPATIAL_LEGENDRE, diff --git a/include/openmc/tallies/filter_particle.h b/include/openmc/tallies/filter_particle.h index 863a6d282fe..fff955ee5f2 100644 --- a/include/openmc/tallies/filter_particle.h +++ b/include/openmc/tallies/filter_particle.h @@ -48,5 +48,45 @@ class ParticleFilter : public Filter { vector particles_; }; +//============================================================================== +//! Bins by type of outgoing particle (e.g. neutron, photon). +//============================================================================== + +class ParticleOutFilter : public Filter { +public: + //---------------------------------------------------------------------------- + // Constructors, destructors + + ~ParticleOutFilter() = default; + + //---------------------------------------------------------------------------- + // Methods + + std::string type_str() const override { return "particleout"; } + FilterType type() const override { return FilterType::PARTICLE_OUT; } + + void from_xml(pugi::xml_node node) override; + + void get_all_bins(const Particle& p, TallyEstimator estimator, + FilterMatch& match) const override; + + void to_statepoint(hid_t filter_group) const override; + + std::string text_label(int bin) const override; + + //---------------------------------------------------------------------------- + // Accessors + + const vector& particles() const { return particles_; } + + void set_particles(span particles); + +private: + //---------------------------------------------------------------------------- + // Data members + + vector particles_; +}; + } // namespace openmc #endif // OPENMC_TALLIES_FILTER_PARTICLE_H diff --git a/include/openmc/tallies/tally.h b/include/openmc/tallies/tally.h index 8c088c46089..4d4ed4f4735 100644 --- a/include/openmc/tallies/tally.h +++ b/include/openmc/tallies/tally.h @@ -202,6 +202,7 @@ extern std::unordered_map tally_map; extern vector> tallies; extern vector active_tallies; extern vector active_analog_tallies; +extern vector active_particleout_analog_tallies; extern vector active_tracklength_tallies; extern vector active_collision_tallies; extern vector active_meshsurf_tallies; diff --git a/include/openmc/tallies/tally_scoring.h b/include/openmc/tallies/tally_scoring.h index 28f1f16223d..043ea6d03d6 100644 --- a/include/openmc/tallies/tally_scoring.h +++ b/include/openmc/tallies/tally_scoring.h @@ -72,14 +72,16 @@ void score_collision_tally(Particle& p); //! Analog tallies are triggered at every collision, not every event. // //! \param p The particle being tracked -void score_analog_tally_ce(Particle& p); +//! \param tallies A vector of the indices of the tallies to score to +void score_analog_tally_ce(Particle& p, const vector& tallies); //! Score tallies based on a simple count of events (for multigroup). // //! Analog tallies are triggered at every collision, not every event. // //! \param p The particle being tracked -void score_analog_tally_mg(Particle& p); +//! \param tallies A vector of the indices of the tallies to score to +void score_analog_tally_mg(Particle& p, const vector& tallies); //! Score tallies using a tracklength estimate of the flux. // diff --git a/openmc/filter.py b/openmc/filter.py index 6a666d2a02c..ae039cbcef0 100644 --- a/openmc/filter.py +++ b/openmc/filter.py @@ -24,9 +24,9 @@ 'universe', 'material', 'cell', 'cellborn', 'surface', 'mesh', 'energy', 'energyout', 'mu', 'musurface', 'polar', 'azimuthal', 'distribcell', 'delayedgroup', 'energyfunction', 'cellfrom', 'materialfrom', 'legendre', 'spatiallegendre', - 'sphericalharmonics', 'zernike', 'zernikeradial', 'particle', 'cellinstance', - 'collision', 'time', 'parentnuclide', 'weight', 'meshborn', 'meshsurface', - 'meshmaterial', + 'sphericalharmonics', 'zernike', 'zernikeradial', 'particle', 'particleout', + 'cellinstance', 'collision', 'time', 'parentnuclide', 'weight', 'meshborn', + 'meshsurface', 'meshmaterial', ) _CURRENT_NAMES = ( @@ -785,6 +785,29 @@ def from_xml_element(cls, elem, **kwargs): filter_id = int(get_text(elem, "id")) bins = get_elem_list(elem, "bins", str) or [] return cls(bins, filter_id=filter_id) + + +class ParticleoutFilter(ParticleFilter): + """Bins tally events based on the outgoing particle type. + + Parameters + ---------- + bins : str, or sequence of str + The particles to tally represented as strings ('neutron', 'photon', + 'electron', 'positron'). + filter_id : int + Unique identifier for the filter + + Attributes + ---------- + bins : sequence of str + The particles to tally + id : int + Unique identifier for the filter + num_bins : Integral + The number of filter bins + + """ class ParentNuclideFilter(ParticleFilter): diff --git a/openmc/lib/filter.py b/openmc/lib/filter.py index 55b681d89ad..3ed03ec590a 100644 --- a/openmc/lib/filter.py +++ b/openmc/lib/filter.py @@ -22,9 +22,9 @@ 'EnergyFilter', 'EnergyoutFilter', 'EnergyFunctionFilter', 'LegendreFilter', 'MaterialFilter', 'MaterialFromFilter', 'MeshFilter', 'MeshBornFilter', 'MeshMaterialFilter', 'MeshSurfaceFilter', 'MuFilter', 'MuSurfaceFilter', - 'ParentNuclideFilter', 'ParticleFilter', 'PolarFilter', 'SphericalHarmonicsFilter', - 'SpatialLegendreFilter', 'SurfaceFilter', 'TimeFilter', 'UniverseFilter', - 'WeightFilter', 'ZernikeFilter', 'ZernikeRadialFilter', 'filters' + 'ParentNuclideFilter', 'ParticleFilter', 'ParticleoutFilter', 'PolarFilter', + 'SphericalHarmonicsFilter', 'SpatialLegendreFilter', 'SurfaceFilter', 'TimeFilter', + 'UniverseFilter', 'WeightFilter', 'ZernikeFilter', 'ZernikeRadialFilter', 'filters' ] # Tally functions @@ -564,6 +564,10 @@ def bins(self): return [ParticleType(i) for i in particle_i] +class ParticleoutFilter(ParticleFilter): + filter_type = 'particleout' + + class PolarFilter(Filter): filter_type = 'polar' diff --git a/src/mesh.cpp b/src/mesh.cpp index 3e4ff1a3eca..97243e0d7d4 100644 --- a/src/mesh.cpp +++ b/src/mesh.cpp @@ -358,7 +358,7 @@ void Mesh::material_volumes(int nx, int ny, int nz, int table_size, site.r[ax1] = min1 + (i1 + 0.5) * d1; site.r[ax2] = min2 + (i2 + 0.5) * d2; - p.from_source(&site); + p.from_source(&site, site.particle); // Determine particle's location if (!exhaustive_find_cell(p)) { diff --git a/src/particle.cpp b/src/particle.cpp index fb41d82e950..6dfb6c5ff95 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -110,7 +110,7 @@ void Particle::split(double wgt) } } -void Particle::from_source(const SourceSite* src) +void Particle::from_source(const SourceSite* src, ParticleType particle) { // Reset some attributes clear(); @@ -141,6 +141,7 @@ void Particle::from_source(const SourceSite* src) E() = data::mg.energy_bin_avg_[g()]; } E_last() = E(); + type_last() = particle; time() = src->time; time_last() = src->time; parent_nuclide() = src->parent_nuclide; @@ -160,6 +161,7 @@ void Particle::event_calculate_xs() // Store pre-collision particle properties wgt_last() = wgt(); E_last() = E(); + type_last() = type(); u_last() = u(); r_last() = r(); time_last() = time(); @@ -348,9 +350,9 @@ void Particle::event_collide() score_collision_tally(*this); if (!model::active_analog_tallies.empty()) { if (settings::run_CE) { - score_analog_tally_ce(*this); + score_analog_tally_ce(*this, model::active_analog_tallies); } else { - score_analog_tally_mg(*this); + score_analog_tally_mg(*this, model::active_analog_tallies); } } @@ -408,6 +410,8 @@ void Particle::event_revive_from_secondary() wgt() = 0.0; } + auto type = this->type(); + // Check for secondary particles if this particle is dead if (!alive()) { // Write final position for this particle @@ -419,11 +423,20 @@ void Particle::event_revive_from_secondary() if (secondary_bank().empty()) return; - from_source(&secondary_bank().back()); + from_source(&secondary_bank().back(), type); secondary_bank().pop_back(); n_event() = 0; bank_second_E() = 0.0; + // Score tallies affected by secondary particles + if (!model::active_particleout_analog_tallies.empty()) { + if (settings::run_CE) { + score_analog_tally_ce(*this, model::active_particleout_analog_tallies); + } else { + score_analog_tally_mg(*this, model::active_particleout_analog_tallies); + } + } + // Subtract secondary particle energy from interim pulse-height results if (!model::active_pulse_height_tallies.empty() && this->type() == ParticleType::photon) { diff --git a/src/particle_restart.cpp b/src/particle_restart.cpp index 6b7778211af..ea4070b951c 100644 --- a/src/particle_restart.cpp +++ b/src/particle_restart.cpp @@ -64,6 +64,7 @@ void read_particle_restart(Particle& p, RunMode& previous_run_mode) p.r_last() = p.r(); p.u_last() = p.u(); p.E_last() = p.E(); + p.type_last() = p.type(); p.g_last() = p.g(); p.time_last() = p.time(); diff --git a/src/random_ray/random_ray.cpp b/src/random_ray/random_ray.cpp index 27f674c4278..5e59acb3691 100644 --- a/src/random_ray/random_ray.cpp +++ b/src/random_ray/random_ray.cpp @@ -796,7 +796,7 @@ void RandomRay::initialize_ray(uint64_t ray_id, FlatSourceDomain* domain) } site.E = 0.0; - this->from_source(&site); + this->from_source(&site, site.particle); // Locate ray if (lowest_coord().cell() == C_NONE) { diff --git a/src/simulation.cpp b/src/simulation.cpp index b9986f47375..a0adbaac9f9 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -567,7 +567,8 @@ void initialize_history(Particle& p, int64_t index_source) // set defaults if (settings::run_mode == RunMode::EIGENVALUE) { // set defaults for eigenvalue simulations from primary bank - p.from_source(&simulation::source_bank[index_source - 1]); + p.from_source( + &simulation::source_bank[index_source - 1], ParticleType::neutron); } else if (settings::run_mode == RunMode::FIXED_SOURCE) { // initialize random number seed int64_t id = (simulation::total_gen + overall_generation() - 1) * @@ -576,7 +577,7 @@ void initialize_history(Particle& p, int64_t index_source) uint64_t seed = init_seed(id, STREAM_SOURCE); // sample from external source distribution or custom library then set auto site = sample_external_source(&seed); - p.from_source(&site); + p.from_source(&site, site.particle); } p.current_work() = index_source; diff --git a/src/tallies/filter.cpp b/src/tallies/filter.cpp index 79817981db0..af23718b24b 100644 --- a/src/tallies/filter.cpp +++ b/src/tallies/filter.cpp @@ -146,6 +146,8 @@ Filter* Filter::create(const std::string& type, int32_t id) return Filter::create(id); } else if (type == "particle") { return Filter::create(id); + } else if (type == "particleout") { + return Filter::create(id); } else if (type == "polar") { return Filter::create(id); } else if (type == "surface") { diff --git a/src/tallies/filter_particle.cpp b/src/tallies/filter_particle.cpp index eef1d1e63c5..80295cded8c 100644 --- a/src/tallies/filter_particle.cpp +++ b/src/tallies/filter_particle.cpp @@ -6,6 +6,10 @@ namespace openmc { +//============================================================================== +// ParticleFilter implementation +//============================================================================== + void ParticleFilter::from_xml(pugi::xml_node node) { auto particles = get_node_array(node, "bins"); @@ -34,8 +38,11 @@ void ParticleFilter::set_particles(span particles) void ParticleFilter::get_all_bins( const Particle& p, TallyEstimator estimator, FilterMatch& match) const { + // Get the pre-collision type of the particle. + auto type = p.type_last(); + for (auto i = 0; i < particles_.size(); i++) { - if (particles_[i] == p.type()) { + if (particles_[i] == type) { match.bins_.push_back(i); match.weights_.push_back(1.0); } @@ -58,6 +65,66 @@ std::string ParticleFilter::text_label(int bin) const return fmt::format("Particle: {}", particle_type_to_str(p)); } +//============================================================================== +// ParticleOutFilter implementation +//============================================================================== + +void ParticleOutFilter::from_xml(pugi::xml_node node) +{ + auto particles = get_node_array(node, "bins"); + + // Convert to vector of ParticleType + vector types; + for (auto& p : particles) { + types.push_back(str_to_particle_type(p)); + } + this->set_particles(types); +} + +void ParticleOutFilter::set_particles(span particles) +{ + // Clear existing particles + particles_.clear(); + particles_.reserve(particles.size()); + + // Set particles and number of bins + for (auto p : particles) { + particles_.push_back(p); + } + n_bins_ = particles_.size(); +} + +void ParticleOutFilter::get_all_bins( + const Particle& p, TallyEstimator estimator, FilterMatch& match) const +{ + for (auto i = 0; i < particles_.size(); i++) { + if (particles_[i] == p.type()) { + match.bins_.push_back(i); + match.weights_.push_back(1.0); + } + } +} + +void ParticleOutFilter::to_statepoint(hid_t filter_group) const +{ + Filter::to_statepoint(filter_group); + vector particles; + for (auto p : particles_) { + particles.push_back(particle_type_to_str(p)); + } + write_dataset(filter_group, "bins", particles); +} + +std::string ParticleOutFilter::text_label(int bin) const +{ + const auto& p = particles_.at(bin); + return fmt::format("ParticleOut: {}", particle_type_to_str(p)); +} + +//============================================================================== +// C-API functions +//============================================================================== + extern "C" int openmc_particle_filter_get_bins(int32_t idx, int bins[]) { if (int err = verify_filter(idx)) diff --git a/src/tallies/tally.cpp b/src/tallies/tally.cpp index ba183860269..cf6760bda19 100644 --- a/src/tallies/tally.cpp +++ b/src/tallies/tally.cpp @@ -55,6 +55,7 @@ std::unordered_map tally_map; vector> tallies; vector active_tallies; vector active_analog_tallies; +vector active_particleout_analog_tallies; vector active_tracklength_tallies; vector active_collision_tallies; vector active_meshsurf_tallies; @@ -152,7 +153,8 @@ Tally::Tally(pugi::xml_node node) // Change the tally estimator if a filter demands it FilterType filt_type = f->type(); if (filt_type == FilterType::ENERGY_OUT || - filt_type == FilterType::LEGENDRE) { + filt_type == FilterType::LEGENDRE || + filt_type == FilterType::PARTICLE_OUT) { estimator_ = TallyEstimator::ANALOG; } else if (filt_type == FilterType::SPHERICAL_HARMONICS) { auto sf = dynamic_cast(f); @@ -569,7 +571,6 @@ void Tally::set_scores(const vector& scores) fatal_error("Cannot tally flux with an outgoing energy filter."); break; - case SCORE_TOTAL: case SCORE_ABSORPTION: case SCORE_FISSION: if (energyout_present) @@ -578,6 +579,7 @@ void Tally::set_scores(const vector& scores) "outgoing energy filter"); break; + case SCORE_TOTAL: case SCORE_SCATTER: if (legendre_present) estimator_ = TallyEstimator::ANALOG; @@ -1068,6 +1070,7 @@ void setup_active_tallies() { model::active_tallies.clear(); model::active_analog_tallies.clear(); + model::active_particleout_analog_tallies.clear(); model::active_tracklength_tallies.clear(); model::active_collision_tallies.clear(); model::active_meshsurf_tallies.clear(); @@ -1079,12 +1082,15 @@ void setup_active_tallies() if (tally.active_) { model::active_tallies.push_back(i); + bool particleout_present = tally.has_filter(FilterType::PARTICLE_OUT); switch (tally.type_) { case TallyType::VOLUME: switch (tally.estimator_) { case TallyEstimator::ANALOG: model::active_analog_tallies.push_back(i); + if (particleout_present) + model::active_particleout_analog_tallies.push_back(i); break; case TallyEstimator::TRACKLENGTH: model::active_tracklength_tallies.push_back(i); @@ -1122,6 +1128,7 @@ void free_memory_tally() model::active_tallies.clear(); model::active_analog_tallies.clear(); + model::active_particleout_analog_tallies.clear(); model::active_tracklength_tallies.clear(); model::active_collision_tallies.clear(); model::active_meshsurf_tallies.clear(); diff --git a/src/tallies/tally_scoring.cpp b/src/tallies/tally_scoring.cpp index e73fb90f311..a1220f2139a 100644 --- a/src/tallies/tally_scoring.cpp +++ b/src/tallies/tally_scoring.cpp @@ -2303,7 +2303,7 @@ void score_general_mg(Particle& p, int i_tally, int start_index, } } -void score_analog_tally_ce(Particle& p) +void score_analog_tally_ce(Particle& p, const vector& tallies) { // Since electrons/positrons are not transported, we assign a flux of zero. // Note that the heating score does NOT use the flux and will be non-zero for @@ -2313,7 +2313,7 @@ void score_analog_tally_ce(Particle& p) ? 1.0 : 0.0; - for (auto i_tally : model::active_analog_tallies) { + for (auto i_tally : tallies) { const Tally& tally {*model::tallies[i_tally]}; // Initialize an iterator over valid filter bin combinations. If there are @@ -2355,9 +2355,9 @@ void score_analog_tally_ce(Particle& p) match.bins_present_ = false; } -void score_analog_tally_mg(Particle& p) +void score_analog_tally_mg(Particle& p, const vector& tallies) { - for (auto i_tally : model::active_analog_tallies) { + for (auto i_tally : tallies) { const Tally& tally {*model::tallies[i_tally]}; // Initialize an iterator over valid filter bin combinations. If there are From e60d78dd92e1cc6fc4bec981f102541f070b477f Mon Sep 17 00:00:00 2001 From: GuySten Date: Tue, 26 Aug 2025 21:21:51 +0300 Subject: [PATCH 2/9] store progeny energy in SourceSite E_last --- include/openmc/particle_data.h | 1 + openmc/lib/core.py | 1 + src/initialize.cpp | 28 +++++++++++++++------------- src/mesh.cpp | 1 + src/particle.cpp | 2 +- src/physics.cpp | 1 + src/physics_mg.cpp | 1 + src/random_ray/random_ray.cpp | 1 + 8 files changed, 22 insertions(+), 14 deletions(-) diff --git a/include/openmc/particle_data.h b/include/openmc/particle_data.h index 235fbf59ac7..109631a7b5a 100644 --- a/include/openmc/particle_data.h +++ b/include/openmc/particle_data.h @@ -44,6 +44,7 @@ struct SourceSite { Position r; Direction u; double E; + double E_last; double time {0.0}; double wgt {1.0}; int delayed_group {0}; diff --git a/openmc/lib/core.py b/openmc/lib/core.py index 9f8db69d577..f65742d8d63 100644 --- a/openmc/lib/core.py +++ b/openmc/lib/core.py @@ -22,6 +22,7 @@ class _SourceSite(Structure): _fields_ = [('r', c_double*3), ('u', c_double*3), ('E', c_double), + ('E_last', c_double), ('time', c_double), ('wgt', c_double), ('delayed_group', c_int), diff --git a/src/initialize.cpp b/src/initialize.cpp index 2da78137d10..bab50aa8c09 100644 --- a/src/initialize.cpp +++ b/src/initialize.cpp @@ -157,26 +157,28 @@ void initialize_mpi(MPI_Comm intracomm) // Create bank datatype SourceSite b; - MPI_Aint disp[11]; + MPI_Aint disp[12]; MPI_Get_address(&b.r, &disp[0]); MPI_Get_address(&b.u, &disp[1]); MPI_Get_address(&b.E, &disp[2]); - MPI_Get_address(&b.time, &disp[3]); - MPI_Get_address(&b.wgt, &disp[4]); - MPI_Get_address(&b.delayed_group, &disp[5]); - MPI_Get_address(&b.surf_id, &disp[6]); - MPI_Get_address(&b.particle, &disp[7]); - MPI_Get_address(&b.parent_nuclide, &disp[8]); - MPI_Get_address(&b.parent_id, &disp[9]); - MPI_Get_address(&b.progeny_id, &disp[10]); - for (int i = 10; i >= 0; --i) { + MPI_Get_address(&b.E_last, &disp[3]); + MPI_Get_address(&b.time, &disp[4]); + MPI_Get_address(&b.wgt, &disp[5]); + MPI_Get_address(&b.delayed_group, &disp[6]); + MPI_Get_address(&b.surf_id, &disp[7]); + MPI_Get_address(&b.particle, &disp[8]); + MPI_Get_address(&b.parent_nuclide, &disp[9]); + MPI_Get_address(&b.parent_id, &disp[10]); + MPI_Get_address(&b.progeny_id, &disp[11]); + for (int i = 11; i >= 0; --i) { disp[i] -= disp[0]; } - int blocks[] {3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1}; + int blocks[] {3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}; MPI_Datatype types[] {MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, - MPI_DOUBLE, MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_LONG, MPI_LONG}; - MPI_Type_create_struct(11, blocks, disp, types, &mpi::source_site); + MPI_DOUBLE, MPI_DOUBLE, MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_LONG, + MPI_LONG}; + MPI_Type_create_struct(12, blocks, disp, types, &mpi::source_site); MPI_Type_commit(&mpi::source_site); } #endif // OPENMC_MPI diff --git a/src/mesh.cpp b/src/mesh.cpp index 97243e0d7d4..4586177f739 100644 --- a/src/mesh.cpp +++ b/src/mesh.cpp @@ -321,6 +321,7 @@ void Mesh::material_volumes(int nx, int ny, int nz, int table_size, SourceSite site; site.E = 1.0; + site.E_last = 1.0; site.particle = ParticleType::neutron; for (int axis = 0; axis < 3; ++axis) { diff --git a/src/particle.cpp b/src/particle.cpp index 6dfb6c5ff95..1ebfbf2fc40 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -132,6 +132,7 @@ void Particle::from_source(const SourceSite* src, ParticleType particle) r_last_current() = src->r; r_last() = src->r; u_last() = src->u; + E_last() = src->E_last; if (settings::run_CE) { E() = src->E; g() = 0; @@ -140,7 +141,6 @@ void Particle::from_source(const SourceSite* src, ParticleType particle) g_last() = static_cast(src->E); E() = data::mg.energy_bin_avg_[g()]; } - E_last() = E(); type_last() = particle; time() = src->time; time_last() = src->time; diff --git a/src/physics.cpp b/src/physics.cpp index 3c06e543de1..1d463174988 100644 --- a/src/physics.cpp +++ b/src/physics.cpp @@ -209,6 +209,7 @@ void create_fission_sites(Particle& p, int i_nuclide, const Reaction& rx) site.r = p.r(); site.particle = ParticleType::neutron; site.time = p.time(); + site.E_last = p.E(); site.wgt = 1. / weight; site.parent_id = p.id(); site.progeny_id = p.n_progeny()++; diff --git a/src/physics_mg.cpp b/src/physics_mg.cpp index 5ebef9c1417..6445ccea356 100644 --- a/src/physics_mg.cpp +++ b/src/physics_mg.cpp @@ -138,6 +138,7 @@ void create_fission_sites(Particle& p) site.r = p.r(); site.particle = ParticleType::neutron; site.time = p.time(); + site.E_last = p.E(); site.wgt = 1. / weight; site.parent_id = p.id(); site.progeny_id = p.n_progeny()++; diff --git a/src/random_ray/random_ray.cpp b/src/random_ray/random_ray.cpp index 5e59acb3691..a3f88345622 100644 --- a/src/random_ray/random_ray.cpp +++ b/src/random_ray/random_ray.cpp @@ -796,6 +796,7 @@ void RandomRay::initialize_ray(uint64_t ray_id, FlatSourceDomain* domain) } site.E = 0.0; + site.E_last = 0.0; this->from_source(&site, site.particle); // Locate ray From fac1c90d84c27eecceefb0e09f3d721d7019012b Mon Sep 17 00:00:00 2001 From: GuySten Date: Tue, 26 Aug 2025 22:06:13 +0300 Subject: [PATCH 3/9] move E_last to be later in SourceSite structs --- include/openmc/particle_data.h | 2 +- openmc/lib/core.py | 2 +- src/initialize.cpp | 14 +++++++------- 3 files changed, 9 insertions(+), 9 deletions(-) diff --git a/include/openmc/particle_data.h b/include/openmc/particle_data.h index 109631a7b5a..c45edadf700 100644 --- a/include/openmc/particle_data.h +++ b/include/openmc/particle_data.h @@ -44,7 +44,6 @@ struct SourceSite { Position r; Direction u; double E; - double E_last; double time {0.0}; double wgt {1.0}; int delayed_group {0}; @@ -52,6 +51,7 @@ struct SourceSite { ParticleType particle; // Extra attributes that don't show up in source written to file + double E_last; int parent_nuclide {-1}; int64_t parent_id; int64_t progeny_id; diff --git a/openmc/lib/core.py b/openmc/lib/core.py index f65742d8d63..68ab289aa1e 100644 --- a/openmc/lib/core.py +++ b/openmc/lib/core.py @@ -22,12 +22,12 @@ class _SourceSite(Structure): _fields_ = [('r', c_double*3), ('u', c_double*3), ('E', c_double), - ('E_last', c_double), ('time', c_double), ('wgt', c_double), ('delayed_group', c_int), ('surf_id', c_int), ('particle', c_int), + ('E_last', c_double), ('parent_nuclide', c_int), ('parent_id', c_int64), ('progeny_id', c_int64)] diff --git a/src/initialize.cpp b/src/initialize.cpp index bab50aa8c09..6d28f0bf84b 100644 --- a/src/initialize.cpp +++ b/src/initialize.cpp @@ -161,12 +161,12 @@ void initialize_mpi(MPI_Comm intracomm) MPI_Get_address(&b.r, &disp[0]); MPI_Get_address(&b.u, &disp[1]); MPI_Get_address(&b.E, &disp[2]); - MPI_Get_address(&b.E_last, &disp[3]); - MPI_Get_address(&b.time, &disp[4]); - MPI_Get_address(&b.wgt, &disp[5]); - MPI_Get_address(&b.delayed_group, &disp[6]); - MPI_Get_address(&b.surf_id, &disp[7]); - MPI_Get_address(&b.particle, &disp[8]); + MPI_Get_address(&b.time, &disp[3]); + MPI_Get_address(&b.wgt, &disp[4]); + MPI_Get_address(&b.delayed_group, &disp[5]); + MPI_Get_address(&b.surf_id, &disp[6]); + MPI_Get_address(&b.particle, &disp[7]); + MPI_Get_address(&b.E_last, &disp[8]); MPI_Get_address(&b.parent_nuclide, &disp[9]); MPI_Get_address(&b.parent_id, &disp[10]); MPI_Get_address(&b.progeny_id, &disp[11]); @@ -176,7 +176,7 @@ void initialize_mpi(MPI_Comm intracomm) int blocks[] {3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}; MPI_Datatype types[] {MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, - MPI_DOUBLE, MPI_DOUBLE, MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_LONG, + MPI_DOUBLE, MPI_INT, MPI_INT, MPI_INT, MPI_DOUBLE, MPI_INT, MPI_LONG, MPI_LONG}; MPI_Type_create_struct(12, blocks, disp, types, &mpi::source_site); MPI_Type_commit(&mpi::source_site); From b3020d1a30692fe97d0da8aec498af260b0d69dd Mon Sep 17 00:00:00 2001 From: GuySten Date: Tue, 26 Aug 2025 22:46:09 +0300 Subject: [PATCH 4/9] try new approach --- include/openmc/particle.h | 3 +-- include/openmc/particle_data.h | 1 - openmc/lib/core.py | 1 - src/initialize.cpp | 18 ++++++++--------- src/mesh.cpp | 3 +-- src/particle.cpp | 36 ++++++++++++++++++++-------------- src/physics.cpp | 1 - src/physics_mg.cpp | 1 - src/random_ray/random_ray.cpp | 3 +-- src/simulation.cpp | 5 ++--- 10 files changed, 34 insertions(+), 38 deletions(-) diff --git a/include/openmc/particle.h b/include/openmc/particle.h index cdfbfbb6b82..0f37719b94f 100644 --- a/include/openmc/particle.h +++ b/include/openmc/particle.h @@ -62,8 +62,7 @@ class Particle : public ParticleData { //! site may have been produced from an external source, from fission, or //! simply as a secondary particle. //! \param src Source site data - //! \param particle ParticleType incident particle type - void from_source(const SourceSite* src, ParticleType particle); + void from_source(const SourceSite* src); // Coarse-grained particle events void event_calculate_xs(); diff --git a/include/openmc/particle_data.h b/include/openmc/particle_data.h index c45edadf700..235fbf59ac7 100644 --- a/include/openmc/particle_data.h +++ b/include/openmc/particle_data.h @@ -51,7 +51,6 @@ struct SourceSite { ParticleType particle; // Extra attributes that don't show up in source written to file - double E_last; int parent_nuclide {-1}; int64_t parent_id; int64_t progeny_id; diff --git a/openmc/lib/core.py b/openmc/lib/core.py index 68ab289aa1e..9f8db69d577 100644 --- a/openmc/lib/core.py +++ b/openmc/lib/core.py @@ -27,7 +27,6 @@ class _SourceSite(Structure): ('delayed_group', c_int), ('surf_id', c_int), ('particle', c_int), - ('E_last', c_double), ('parent_nuclide', c_int), ('parent_id', c_int64), ('progeny_id', c_int64)] diff --git a/src/initialize.cpp b/src/initialize.cpp index 6d28f0bf84b..2da78137d10 100644 --- a/src/initialize.cpp +++ b/src/initialize.cpp @@ -157,7 +157,7 @@ void initialize_mpi(MPI_Comm intracomm) // Create bank datatype SourceSite b; - MPI_Aint disp[12]; + MPI_Aint disp[11]; MPI_Get_address(&b.r, &disp[0]); MPI_Get_address(&b.u, &disp[1]); MPI_Get_address(&b.E, &disp[2]); @@ -166,19 +166,17 @@ void initialize_mpi(MPI_Comm intracomm) MPI_Get_address(&b.delayed_group, &disp[5]); MPI_Get_address(&b.surf_id, &disp[6]); MPI_Get_address(&b.particle, &disp[7]); - MPI_Get_address(&b.E_last, &disp[8]); - MPI_Get_address(&b.parent_nuclide, &disp[9]); - MPI_Get_address(&b.parent_id, &disp[10]); - MPI_Get_address(&b.progeny_id, &disp[11]); - for (int i = 11; i >= 0; --i) { + MPI_Get_address(&b.parent_nuclide, &disp[8]); + MPI_Get_address(&b.parent_id, &disp[9]); + MPI_Get_address(&b.progeny_id, &disp[10]); + for (int i = 10; i >= 0; --i) { disp[i] -= disp[0]; } - int blocks[] {3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1}; + int blocks[] {3, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1}; MPI_Datatype types[] {MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, MPI_DOUBLE, - MPI_DOUBLE, MPI_INT, MPI_INT, MPI_INT, MPI_DOUBLE, MPI_INT, MPI_LONG, - MPI_LONG}; - MPI_Type_create_struct(12, blocks, disp, types, &mpi::source_site); + MPI_DOUBLE, MPI_INT, MPI_INT, MPI_INT, MPI_INT, MPI_LONG, MPI_LONG}; + MPI_Type_create_struct(11, blocks, disp, types, &mpi::source_site); MPI_Type_commit(&mpi::source_site); } #endif // OPENMC_MPI diff --git a/src/mesh.cpp b/src/mesh.cpp index 4586177f739..3e4ff1a3eca 100644 --- a/src/mesh.cpp +++ b/src/mesh.cpp @@ -321,7 +321,6 @@ void Mesh::material_volumes(int nx, int ny, int nz, int table_size, SourceSite site; site.E = 1.0; - site.E_last = 1.0; site.particle = ParticleType::neutron; for (int axis = 0; axis < 3; ++axis) { @@ -359,7 +358,7 @@ void Mesh::material_volumes(int nx, int ny, int nz, int table_size, site.r[ax1] = min1 + (i1 + 0.5) * d1; site.r[ax2] = min2 + (i2 + 0.5) * d2; - p.from_source(&site, site.particle); + p.from_source(&site); // Determine particle's location if (!exhaustive_find_cell(p)) { diff --git a/src/particle.cpp b/src/particle.cpp index 1ebfbf2fc40..f6002bfd23f 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -88,6 +88,23 @@ bool Particle::create_secondary( bank.E = settings::run_CE ? E : g(); bank.time = time(); bank_second_E() += bank.E; + + // Score tallies affected by secondary particles + if (!model::active_particleout_analog_tallies.empty()) { + // Create secondary particle for tallying purposes only + Particle p; + p.from_source(&bank); + p.u_last() = this->u(); + p.r_last() = this->r(); + p.E_last() = this->E(); + p.type_last() = this->type(); + + if (settings::run_CE) { + score_analog_tally_ce(p, model::active_particleout_analog_tallies); + } else { + score_analog_tally_mg(p, model::active_particleout_analog_tallies); + } + } return true; } @@ -110,7 +127,7 @@ void Particle::split(double wgt) } } -void Particle::from_source(const SourceSite* src, ParticleType particle) +void Particle::from_source(const SourceSite* src) { // Reset some attributes clear(); @@ -124,6 +141,7 @@ void Particle::from_source(const SourceSite* src, ParticleType particle) // Copy attributes from source bank site type() = src->particle; + type_last() = src->particle; wgt() = src->wgt; wgt_last() = src->wgt; r() = src->r; @@ -132,7 +150,6 @@ void Particle::from_source(const SourceSite* src, ParticleType particle) r_last_current() = src->r; r_last() = src->r; u_last() = src->u; - E_last() = src->E_last; if (settings::run_CE) { E() = src->E; g() = 0; @@ -141,7 +158,7 @@ void Particle::from_source(const SourceSite* src, ParticleType particle) g_last() = static_cast(src->E); E() = data::mg.energy_bin_avg_[g()]; } - type_last() = particle; + E_last() = E(); time() = src->time; time_last() = src->time; parent_nuclide() = src->parent_nuclide; @@ -410,8 +427,6 @@ void Particle::event_revive_from_secondary() wgt() = 0.0; } - auto type = this->type(); - // Check for secondary particles if this particle is dead if (!alive()) { // Write final position for this particle @@ -423,20 +438,11 @@ void Particle::event_revive_from_secondary() if (secondary_bank().empty()) return; - from_source(&secondary_bank().back(), type); + from_source(&secondary_bank().back()); secondary_bank().pop_back(); n_event() = 0; bank_second_E() = 0.0; - // Score tallies affected by secondary particles - if (!model::active_particleout_analog_tallies.empty()) { - if (settings::run_CE) { - score_analog_tally_ce(*this, model::active_particleout_analog_tallies); - } else { - score_analog_tally_mg(*this, model::active_particleout_analog_tallies); - } - } - // Subtract secondary particle energy from interim pulse-height results if (!model::active_pulse_height_tallies.empty() && this->type() == ParticleType::photon) { diff --git a/src/physics.cpp b/src/physics.cpp index 1d463174988..3c06e543de1 100644 --- a/src/physics.cpp +++ b/src/physics.cpp @@ -209,7 +209,6 @@ void create_fission_sites(Particle& p, int i_nuclide, const Reaction& rx) site.r = p.r(); site.particle = ParticleType::neutron; site.time = p.time(); - site.E_last = p.E(); site.wgt = 1. / weight; site.parent_id = p.id(); site.progeny_id = p.n_progeny()++; diff --git a/src/physics_mg.cpp b/src/physics_mg.cpp index 6445ccea356..5ebef9c1417 100644 --- a/src/physics_mg.cpp +++ b/src/physics_mg.cpp @@ -138,7 +138,6 @@ void create_fission_sites(Particle& p) site.r = p.r(); site.particle = ParticleType::neutron; site.time = p.time(); - site.E_last = p.E(); site.wgt = 1. / weight; site.parent_id = p.id(); site.progeny_id = p.n_progeny()++; diff --git a/src/random_ray/random_ray.cpp b/src/random_ray/random_ray.cpp index a3f88345622..27f674c4278 100644 --- a/src/random_ray/random_ray.cpp +++ b/src/random_ray/random_ray.cpp @@ -796,8 +796,7 @@ void RandomRay::initialize_ray(uint64_t ray_id, FlatSourceDomain* domain) } site.E = 0.0; - site.E_last = 0.0; - this->from_source(&site, site.particle); + this->from_source(&site); // Locate ray if (lowest_coord().cell() == C_NONE) { diff --git a/src/simulation.cpp b/src/simulation.cpp index a0adbaac9f9..b9986f47375 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -567,8 +567,7 @@ void initialize_history(Particle& p, int64_t index_source) // set defaults if (settings::run_mode == RunMode::EIGENVALUE) { // set defaults for eigenvalue simulations from primary bank - p.from_source( - &simulation::source_bank[index_source - 1], ParticleType::neutron); + p.from_source(&simulation::source_bank[index_source - 1]); } else if (settings::run_mode == RunMode::FIXED_SOURCE) { // initialize random number seed int64_t id = (simulation::total_gen + overall_generation() - 1) * @@ -577,7 +576,7 @@ void initialize_history(Particle& p, int64_t index_source) uint64_t seed = init_seed(id, STREAM_SOURCE); // sample from external source distribution or custom library then set auto site = sample_external_source(&seed); - p.from_source(&site, site.particle); + p.from_source(&site); } p.current_work() = index_source; From 4e47f041b7968f572418a0121e7613af7f499049 Mon Sep 17 00:00:00 2001 From: GuySten Date: Wed, 27 Aug 2025 18:30:10 +0300 Subject: [PATCH 5/9] pre allocate tmp parricle --- src/particle.cpp | 21 +++++++++++++-------- 1 file changed, 13 insertions(+), 8 deletions(-) diff --git a/src/particle.cpp b/src/particle.cpp index f6002bfd23f..9d723bba858 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -38,6 +38,10 @@ namespace openmc { +namespace simulation { +thread_local Particle tmp_particle; +} + //============================================================================== // Particle implementation //============================================================================== @@ -92,17 +96,18 @@ bool Particle::create_secondary( // Score tallies affected by secondary particles if (!model::active_particleout_analog_tallies.empty()) { // Create secondary particle for tallying purposes only - Particle p; - p.from_source(&bank); - p.u_last() = this->u(); - p.r_last() = this->r(); - p.E_last() = this->E(); - p.type_last() = this->type(); + simulation::tmp_particle.from_source(&bank); + simulation::tmp_particle.u_last() = this->u(); + simulation::tmp_particle.r_last() = this->r(); + simulation::tmp_particle.E_last() = this->E(); + simulation::tmp_particle.type_last() = this->type(); if (settings::run_CE) { - score_analog_tally_ce(p, model::active_particleout_analog_tallies); + score_analog_tally_ce( + simulation::tmp_particle, model::active_particleout_analog_tallies); } else { - score_analog_tally_mg(p, model::active_particleout_analog_tallies); + score_analog_tally_mg( + simulation::tmp_particle, model::active_particleout_analog_tallies); } } return true; From 32be041fc74be6e451e4325010438b37214294cb Mon Sep 17 00:00:00 2001 From: GuySten Date: Sat, 31 Jan 2026 21:59:06 +0200 Subject: [PATCH 6/9] commit review suggestions --- src/particle.cpp | 25 ++++++++++++------------- 1 file changed, 12 insertions(+), 13 deletions(-) diff --git a/src/particle.cpp b/src/particle.cpp index 2a9dea19d19..e3b4405f7a2 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -39,10 +39,6 @@ namespace openmc { -namespace simulation { -thread_local Particle tmp_particle; -} - //============================================================================== // Particle implementation //============================================================================== @@ -97,18 +93,21 @@ bool Particle::create_secondary( // Score tallies affected by secondary particles if (!model::active_particleout_analog_tallies.empty()) { // Create secondary particle for tallying purposes only - simulation::tmp_particle.from_source(&bank); - simulation::tmp_particle.u_last() = this->u(); - simulation::tmp_particle.r_last() = this->r(); - simulation::tmp_particle.E_last() = this->E(); - simulation::tmp_particle.type_last() = this->type(); + Particle tmp; + tmp.from_source(&bank); + tmp.u_last() = this->u(); + tmp.r_last() = this->r(); + tmp.E_last() = this->E(); + tmp.type_last() = this->type(); + + // Load geometry info + if (!exhaustive_find_cell(tmp)) + fatal_error("Cannot find temporary particle in model."); if (settings::run_CE) { - score_analog_tally_ce( - simulation::tmp_particle, model::active_particleout_analog_tallies); + score_analog_tally_ce(tmp, model::active_particleout_analog_tallies); } else { - score_analog_tally_mg( - simulation::tmp_particle, model::active_particleout_analog_tallies); + score_analog_tally_mg(tmp, model::active_particleout_analog_tallies); } } return true; From f9a5b4a6feb0d0dcd4c76bfd141f3ffa42c41c6a Mon Sep 17 00:00:00 2001 From: GuySten Date: Sat, 31 Jan 2026 22:29:07 +0200 Subject: [PATCH 7/9] add a regression test --- .../filter_particleout/__init__.py | 0 .../filter_particleout/inputs_true.dat | 43 ++++ .../filter_particleout/results_true.dat | 199 ++++++++++++++++++ .../filter_particleout/test.py | 66 ++++++ 4 files changed, 308 insertions(+) create mode 100644 tests/regression_tests/filter_particleout/__init__.py create mode 100644 tests/regression_tests/filter_particleout/inputs_true.dat create mode 100644 tests/regression_tests/filter_particleout/results_true.dat create mode 100644 tests/regression_tests/filter_particleout/test.py diff --git a/tests/regression_tests/filter_particleout/__init__.py b/tests/regression_tests/filter_particleout/__init__.py new file mode 100644 index 00000000000..e69de29bb2d diff --git a/tests/regression_tests/filter_particleout/inputs_true.dat b/tests/regression_tests/filter_particleout/inputs_true.dat new file mode 100644 index 00000000000..6e87cc1f22d --- /dev/null +++ b/tests/regression_tests/filter_particleout/inputs_true.dat @@ -0,0 +1,43 @@ + + + + + + + + + + + + + + fixed source + 10000 + 10 + + + 0.0 0.0 0.0 + + + 2000000.0 1.0 + + + true + + + + neutron + + + photon + + + 10000.0 10974.987654930568 12045.035402587811 13219.411484660288 14508.287784959402 15922.82793341094 17475.28400007683 19179.10261672489 21049.04144512022 23101.29700083158 25353.644939701113 27825.59402207126 30538.555088334124 33516.02650938841 36783.79771828634 40370.17258596558 44306.21457583878 48626.015800653535 53366.99231206313 58570.20818056673 64280.73117284319 70548.02310718645 77426.36826811278 84975.34359086439 93260.334688322 102353.10218990268 112332.40329780265 123284.67394420659 135304.77745798076 148496.82622544636 162975.08346206436 178864.9529057435 196304.06500402724 215443.46900318866 236448.9412645407 259502.42113997374 284803.5868435805 312571.58496882353 343046.9286314919 376493.58067924716 413201.24001153343 453487.8508128582 497702.35643321136 546227.7217684337 599484.2503189408 657933.2246575683 722080.9018385471 792482.8983539186 869749.0026177834 954548.4566618347 1047615.7527896662 1149756.9953977356 1261856.8830660211 1384886.3713938745 1519911.0829529332 1668100.537200059 1830738.2802953697 2009233.0025650458 2205130.739903046 2420128.264794383 2656087.782946684 2915053.0628251815 3199267.1377973845 3511191.7342151348 3853528.593710535 4229242.8743894985 4641588.833612782 5094138.014816386 5590810.182512223 6135907.273413176 6734150.657750828 7390722.033525775 8111308.307896872 8902150.854450393 9770099.572992247 10722672.220103253 11768119.524349991 12915496.650148828 14174741.629268076 15556761.439304722 17073526.47470692 18738174.22860387 20565123.083486516 22570197.196339216 24770763.55991714 27185882.427329402 29836472.402833402 32745491.628777318 35938136.63804626 39442060.59437664 43287612.810830615 47508101.62102793 52140082.879996955 57223676.5935022 62802914.41834259 68926121.04349709 75646332.7554629 83021756.81319752 91116275.61154905 100000000.0 + + + 1 2 3 + total + analog + + + diff --git a/tests/regression_tests/filter_particleout/results_true.dat b/tests/regression_tests/filter_particleout/results_true.dat new file mode 100644 index 00000000000..27610573cd3 --- /dev/null +++ b/tests/regression_tests/filter_particleout/results_true.dat @@ -0,0 +1,199 @@ +tally 1: +2.803000E-01 +7.899090E-03 +3.483000E-01 +1.214713E-02 +3.907000E-01 +1.532231E-02 +4.527000E-01 +2.054189E-02 +5.254000E-01 +2.765018E-02 +5.925000E-01 +3.513577E-02 +7.129000E-01 +5.086499E-02 +8.202000E-01 +6.733918E-02 +1.017600E+00 +1.036527E-01 +1.247900E+00 +1.559066E-01 +1.520700E+00 +2.313731E-01 +1.823900E+00 +3.328213E-01 +2.282400E+00 +5.211727E-01 +2.852200E+00 +8.138458E-01 +3.568000E+00 +1.273357E+00 +4.329900E+00 +1.875591E+00 +5.602350E+01 +3.138822E+02 +1.124200E+00 +1.264838E-01 +1.167100E+00 +1.364348E-01 +1.299500E+00 +1.690742E-01 +1.525800E+00 +2.330320E-01 +1.913600E+00 +3.663497E-01 +2.032500E+00 +4.133832E-01 +2.613100E+00 +6.830790E-01 +3.794300E+00 +1.439925E+00 +6.613200E+00 +4.374276E+00 +3.068100E+00 +9.415793E-01 +2.519200E+00 +6.348918E-01 +2.344800E+00 +5.499640E-01 +2.259700E+00 +5.106904E-01 +2.151800E+00 +4.632270E-01 +2.317600E+00 +5.374687E-01 +2.557800E+00 +6.544854E-01 +2.704600E+00 +7.319093E-01 +2.852900E+00 +8.141188E-01 +3.125500E+00 +9.772850E-01 +3.300000E+00 +1.089410E+00 +3.476100E+00 +1.208627E+00 +3.596300E+00 +1.293913E+00 +3.766100E+00 +1.418948E+00 +3.907500E+00 +1.527112E+00 +4.048700E+00 +1.639931E+00 +4.825700E+00 +2.329247E+00 +6.058500E+00 +3.671578E+00 +6.091100E+00 +3.711118E+00 +5.054100E+00 +2.554896E+00 +4.319200E+00 +1.865836E+00 +4.292900E+00 +1.843096E+00 +4.404600E+00 +1.940307E+00 +3.878800E+00 +1.504872E+00 +3.778500E+00 +1.428066E+00 +3.875500E+00 +1.502577E+00 +3.373500E+00 +1.138456E+00 +3.549000E+00 +1.259773E+00 +3.510600E+00 +1.232725E+00 +3.731000E+00 +1.392373E+00 +3.907600E+00 +1.527663E+00 +3.640600E+00 +1.325788E+00 +3.090100E+00 +9.549517E-01 +2.517400E+00 +6.341261E-01 +2.075400E+00 +4.308890E-01 +1.555100E+00 +2.420523E-01 +1.005800E+00 +1.012240E-01 +8.566000E-01 +7.343774E-02 +1.405700E+00 +1.978421E-01 +2.084000E-01 +4.349900E-03 +1.401000E-01 +1.984950E-03 +6.610000E-02 +4.420700E-04 +2.790000E-02 +8.125000E-05 +1.790000E-02 +3.465000E-05 +1.670000E-02 +3.023000E-05 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 +0.000000E+00 diff --git a/tests/regression_tests/filter_particleout/test.py b/tests/regression_tests/filter_particleout/test.py new file mode 100644 index 00000000000..6bf437fac0e --- /dev/null +++ b/tests/regression_tests/filter_particleout/test.py @@ -0,0 +1,66 @@ +import numpy as np +import openmc +import pytest + +from tests.testing_harness import PyAPITestHarness + +@pytest.fixture +def model(): + openmc.reset_auto_ids() + model = openmc.Model() + + # ============================================================================= + # Materials + # ============================================================================= + + u238 = openmc.Material(name='U-238') + u238.add_nuclide('U238', 1.0) + u238.set_density('g/cm3', 19.1) + model.materials = openmc.Materials([u238]) + + # ============================================================================= + # Geometry + # ============================================================================= + + sphere = openmc.Sphere(r=1e90, boundary_type='vacuum') + cell = openmc.Cell(fill=u238, region=-sphere) + universe = openmc.Universe(cells=[cell]) + model.geometry = openmc.Geometry(universe) + + # ============================================================================= + # Settings + # ============================================================================= + + source = openmc.Source() + source.space = openmc.stats.Point() + source.energy = openmc.stats.Discrete([2.0e6], [1.0]) + model.settings = openmc.Settings() + model.settings.run_mode = 'fixed source' + model.settings.particles = 10_000 + model.settings.batches = 10 + model.settings.source = source + model.settings.photon_transport = True + + # ============================================================================= + # Tallies + # ============================================================================= + + energy_bins = np.logspace(4, 8, 100) # 10 keV to 100 MeV + particle_filter = openmc.ParticleFilter(['neutron']) + particleout_filter = openmc.ParticleoutFilter(['photon']) + energyout_filter = openmc.EnergyoutFilter(energy_bins) + cell_filter = openmc.CellFilter([cell]) + material_filter = openmc.MaterialFilter([u238]) + + # Define tally + tally = openmc.Tally(name='photon spectrum') + tally.filters = [particle_filter, particleout_filter, energyout_filter] + tally.scores = ['total'] + tally.estimator = 'analog' + model.tallies = openmc.Tallies([tally]) + return model + + +def test_filter_particleout(model): + harness = PyAPITestHarness("statepoint.10.h5", model=model) + harness.main() From a20d1123db4b98a53b890fee5449c705e48c269b Mon Sep 17 00:00:00 2001 From: GuySten Date: Sat, 31 Jan 2026 22:37:55 +0200 Subject: [PATCH 8/9] ran clang format --- src/particle.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/particle.cpp b/src/particle.cpp index e3b4405f7a2..5536feca46a 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -99,7 +99,7 @@ bool Particle::create_secondary( tmp.r_last() = this->r(); tmp.E_last() = this->E(); tmp.type_last() = this->type(); - + // Load geometry info if (!exhaustive_find_cell(tmp)) fatal_error("Cannot find temporary particle in model."); From 6e8322eaa68f69ee8a7f32cdf90d0cdc3a9453c3 Mon Sep 17 00:00:00 2001 From: GuySten Date: Sun, 1 Feb 2026 01:54:21 +0200 Subject: [PATCH 9/9] explicitly store energy bounds --- .../filter_particleout/inputs_true.dat | 2 +- .../filter_particleout/test.py | 28 ++++++++++++++++++- 2 files changed, 28 insertions(+), 2 deletions(-) diff --git a/tests/regression_tests/filter_particleout/inputs_true.dat b/tests/regression_tests/filter_particleout/inputs_true.dat index 6e87cc1f22d..42619630960 100644 --- a/tests/regression_tests/filter_particleout/inputs_true.dat +++ b/tests/regression_tests/filter_particleout/inputs_true.dat @@ -32,7 +32,7 @@ photon - 10000.0 10974.987654930568 12045.035402587811 13219.411484660288 14508.287784959402 15922.82793341094 17475.28400007683 19179.10261672489 21049.04144512022 23101.29700083158 25353.644939701113 27825.59402207126 30538.555088334124 33516.02650938841 36783.79771828634 40370.17258596558 44306.21457583878 48626.015800653535 53366.99231206313 58570.20818056673 64280.73117284319 70548.02310718645 77426.36826811278 84975.34359086439 93260.334688322 102353.10218990268 112332.40329780265 123284.67394420659 135304.77745798076 148496.82622544636 162975.08346206436 178864.9529057435 196304.06500402724 215443.46900318866 236448.9412645407 259502.42113997374 284803.5868435805 312571.58496882353 343046.9286314919 376493.58067924716 413201.24001153343 453487.8508128582 497702.35643321136 546227.7217684337 599484.2503189408 657933.2246575683 722080.9018385471 792482.8983539186 869749.0026177834 954548.4566618347 1047615.7527896662 1149756.9953977356 1261856.8830660211 1384886.3713938745 1519911.0829529332 1668100.537200059 1830738.2802953697 2009233.0025650458 2205130.739903046 2420128.264794383 2656087.782946684 2915053.0628251815 3199267.1377973845 3511191.7342151348 3853528.593710535 4229242.8743894985 4641588.833612782 5094138.014816386 5590810.182512223 6135907.273413176 6734150.657750828 7390722.033525775 8111308.307896872 8902150.854450393 9770099.572992247 10722672.220103253 11768119.524349991 12915496.650148828 14174741.629268076 15556761.439304722 17073526.47470692 18738174.22860387 20565123.083486516 22570197.196339216 24770763.55991714 27185882.427329402 29836472.402833402 32745491.628777318 35938136.63804626 39442060.59437664 43287612.810830615 47508101.62102793 52140082.879996955 57223676.5935022 62802914.41834259 68926121.04349709 75646332.7554629 83021756.81319752 91116275.61154905 100000000.0 + 10000.0 10974.9877 12045.0354 13219.4115 14508.2878 15922.8279 17475.284 19179.1026 21049.0414 23101.297 25353.6449 27825.594 30538.5551 33516.0265 36783.7977 40370.1726 44306.2146 48626.0158 53366.9923 58570.2082 64280.7312 70548.0231 77426.3683 84975.3436 93260.3347 102353.102 112332.403 123284.674 135304.777 148496.826 162975.083 178864.953 196304.065 215443.469 236448.941 259502.421 284803.587 312571.585 343046.929 376493.581 413201.24 453487.851 497702.356 546227.722 599484.25 657933.225 722080.902 792482.898 869749.003 954548.457 1047615.75 1149757.0 1261856.88 1384886.37 1519911.08 1668100.54 1830738.28 2009233.0 2205130.74 2420128.26 2656087.78 2915053.06 3199267.14 3511191.73 3853528.59 4229242.87 4641588.83 5094138.01 5590810.18 6135907.27 6734150.66 7390722.03 8111308.31 8902150.85 9770099.57 10722672.2 11768119.5 12915496.7 14174741.6 15556761.4 17073526.5 18738174.2 20565123.1 22570197.2 24770763.6 27185882.4 29836472.4 32745491.6 35938136.6 39442060.6 43287612.8 47508101.6 52140082.9 57223676.6 62802914.4 68926121.0 75646332.8 83021756.8 91116275.6 100000000.0 1 2 3 diff --git a/tests/regression_tests/filter_particleout/test.py b/tests/regression_tests/filter_particleout/test.py index 6bf437fac0e..45f4d43e262 100644 --- a/tests/regression_tests/filter_particleout/test.py +++ b/tests/regression_tests/filter_particleout/test.py @@ -45,7 +45,33 @@ def model(): # Tallies # ============================================================================= - energy_bins = np.logspace(4, 8, 100) # 10 keV to 100 MeV + energy_bins = np.array( + [1.00000000e+04, 1.09749877e+04, 1.20450354e+04, 1.32194115e+04, + 1.45082878e+04, 1.59228279e+04, 1.74752840e+04, 1.91791026e+04, + 2.10490414e+04, 2.31012970e+04, 2.53536449e+04, 2.78255940e+04, + 3.05385551e+04, 3.35160265e+04, 3.67837977e+04, 4.03701726e+04, + 4.43062146e+04, 4.86260158e+04, 5.33669923e+04, 5.85702082e+04, + 6.42807312e+04, 7.05480231e+04, 7.74263683e+04, 8.49753436e+04, + 9.32603347e+04, 1.02353102e+05, 1.12332403e+05, 1.23284674e+05, + 1.35304777e+05, 1.48496826e+05, 1.62975083e+05, 1.78864953e+05, + 1.96304065e+05, 2.15443469e+05, 2.36448941e+05, 2.59502421e+05, + 2.84803587e+05, 3.12571585e+05, 3.43046929e+05, 3.76493581e+05, + 4.13201240e+05, 4.53487851e+05, 4.97702356e+05, 5.46227722e+05, + 5.99484250e+05, 6.57933225e+05, 7.22080902e+05, 7.92482898e+05, + 8.69749003e+05, 9.54548457e+05, 1.04761575e+06, 1.14975700e+06, + 1.26185688e+06, 1.38488637e+06, 1.51991108e+06, 1.66810054e+06, + 1.83073828e+06, 2.00923300e+06, 2.20513074e+06, 2.42012826e+06, + 2.65608778e+06, 2.91505306e+06, 3.19926714e+06, 3.51119173e+06, + 3.85352859e+06, 4.22924287e+06, 4.64158883e+06, 5.09413801e+06, + 5.59081018e+06, 6.13590727e+06, 6.73415066e+06, 7.39072203e+06, + 8.11130831e+06, 8.90215085e+06, 9.77009957e+06, 1.07226722e+07, + 1.17681195e+07, 1.29154967e+07, 1.41747416e+07, 1.55567614e+07, + 1.70735265e+07, 1.87381742e+07, 2.05651231e+07, 2.25701972e+07, + 2.47707636e+07, 2.71858824e+07, 2.98364724e+07, 3.27454916e+07, + 3.59381366e+07, 3.94420606e+07, 4.32876128e+07, 4.75081016e+07, + 5.21400829e+07, 5.72236766e+07, 6.28029144e+07, 6.89261210e+07, + 7.56463328e+07, 8.30217568e+07, 9.11162756e+07, 1.00000000e+08]) + particle_filter = openmc.ParticleFilter(['neutron']) particleout_filter = openmc.ParticleoutFilter(['photon']) energyout_filter = openmc.EnergyoutFilter(energy_bins)