diff --git a/include/openmc/particle_data.h b/include/openmc/particle_data.h index a07a99ffb..82b101c22 100644 --- a/include/openmc/particle_data.h +++ b/include/openmc/particle_data.h @@ -435,6 +435,7 @@ private: double wgt_ {1.0}; double wgt_born_ {1.0}; + double wgt_ww_born_ {-1.0}; double mu_; double time_ {0.0}; double time_last_ {0.0}; @@ -545,6 +546,10 @@ public: double& wgt_born() { return wgt_born_; } double wgt_born() const { return wgt_born_; } + // Weight window value at birth + double& wgt_ww_born() { return wgt_ww_born_; } + const double& wgt_ww_born() const { return wgt_ww_born_; } + // Statistic weight of particle at last collision double& wgt_last() { return wgt_last_; } const double& wgt_last() const { return wgt_last_; } diff --git a/src/simulation.cpp b/src/simulation.cpp index e4521c06d..f8d4e5470 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -614,6 +614,10 @@ void initialize_history(Particle& p, int64_t index_source) // Set particle track. p.write_track() = check_track_criteria(p); + // Set the particle's initial weight window value. + p.wgt_ww_born() = -1.0; + apply_weight_windows(p); + // Display message if high verbosity or trace is on if (settings::verbosity >= 9 || p.trace()) { write_message("Simulating Particle {}", p.id()); diff --git a/src/weight_windows.cpp b/src/weight_windows.cpp index 0efa5a6ae..50f3fb2ed 100644 --- a/src/weight_windows.cpp +++ b/src/weight_windows.cpp @@ -73,10 +73,26 @@ void apply_weight_windows(Particle& p) if (weight_window.is_valid()) break; } + + // If particle has not yet had its birth weight window value set, set it to + // the current weight window (or 1.0 if not born in a weight window). + if (p.wgt_ww_born() == -1.0) { + if (weight_window.is_valid()) { + p.wgt_ww_born() = + (weight_window.lower_weight + weight_window.upper_weight) / 2; + } else { + p.wgt_ww_born() = 1.0; + } + } + // particle is not in any of the ww domains, do nothing if (!weight_window.is_valid()) return; + // Normalize weight windows based on particle's starting weight + // and the value of the weight window the particle was born in. + weight_window.scale(p.wgt_born() / p.wgt_ww_born()); + // get the paramters double weight = p.wgt();