mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-21 14:35:27 -04:00
Updates to object lifecycle for WeightWindows and WeightWindowGenerators (#2582)
This commit is contained in:
Patrick Shriwise
GitHub
parent
36f229eb01
commit
382bcb2e8e
6 changed files with 275 additions and 211 deletions
|
|
@ -15,211 +15,215 @@
|
|||
#include "openmc/tallies/tally.h"
|
||||
#include "openmc/vector.h"
|
||||
|
||||
namespace openmc
|
||||
{
|
||||
namespace openmc {
|
||||
|
||||
enum class WeightWindowUpdateMethod {
|
||||
MAGIC,
|
||||
};
|
||||
enum class WeightWindowUpdateMethod {
|
||||
MAGIC,
|
||||
};
|
||||
|
||||
//==============================================================================
|
||||
// Constants
|
||||
//==============================================================================
|
||||
//==============================================================================
|
||||
// Constants
|
||||
//==============================================================================
|
||||
|
||||
constexpr double DEFAULT_WEIGHT_CUTOFF {1.0e-38}; // default low weight cutoff
|
||||
constexpr double DEFAULT_WEIGHT_CUTOFF {1.0e-38}; // default low weight cutoff
|
||||
|
||||
//==============================================================================
|
||||
// Non-member functions
|
||||
//==============================================================================
|
||||
//==============================================================================
|
||||
// Non-member functions
|
||||
//==============================================================================
|
||||
|
||||
//! Apply weight windows to a particle
|
||||
//! \param[in] p Particle to apply weight windows to
|
||||
void apply_weight_windows(Particle & p);
|
||||
//! Apply weight windows to a particle
|
||||
//! \param[in] p Particle to apply weight windows to
|
||||
void apply_weight_windows(Particle& p);
|
||||
|
||||
//! Free memory associated with weight windows
|
||||
void free_memory_weight_windows();
|
||||
//! Free memory associated with weight windows
|
||||
void free_memory_weight_windows();
|
||||
|
||||
//==============================================================================
|
||||
// Global variables
|
||||
//==============================================================================
|
||||
//==============================================================================
|
||||
// Global variables
|
||||
//==============================================================================
|
||||
|
||||
class WeightWindows;
|
||||
class WeightWindowsGenerator;
|
||||
class WeightWindows;
|
||||
class WeightWindowsGenerator;
|
||||
|
||||
namespace variance_reduction {
|
||||
namespace variance_reduction {
|
||||
|
||||
extern std::unordered_map<int32_t, int32_t> ww_map;
|
||||
extern vector<unique_ptr<WeightWindows>> weight_windows;
|
||||
extern vector<unique_ptr<WeightWindowsGenerator>> weight_windows_generators;
|
||||
extern std::unordered_map<int32_t, int32_t> ww_map;
|
||||
extern vector<unique_ptr<WeightWindows>> weight_windows;
|
||||
extern vector<unique_ptr<WeightWindowsGenerator>> weight_windows_generators;
|
||||
|
||||
} // namespace variance_reduction
|
||||
} // namespace variance_reduction
|
||||
|
||||
//==============================================================================
|
||||
//! Individual weight window information
|
||||
//==============================================================================
|
||||
//==============================================================================
|
||||
//! Individual weight window information
|
||||
//==============================================================================
|
||||
|
||||
struct WeightWindow {
|
||||
double lower_weight {-1}; // -1 indicates invalid state
|
||||
double upper_weight {1};
|
||||
double max_lb_ratio {1};
|
||||
double survival_weight {0.5};
|
||||
double weight_cutoff {DEFAULT_WEIGHT_CUTOFF};
|
||||
int max_split {1};
|
||||
struct WeightWindow {
|
||||
double lower_weight {-1}; // -1 indicates invalid state
|
||||
double upper_weight {1};
|
||||
double max_lb_ratio {1};
|
||||
double survival_weight {0.5};
|
||||
double weight_cutoff {DEFAULT_WEIGHT_CUTOFF};
|
||||
int max_split {10};
|
||||
|
||||
//! Whether the weight window is in a valid state
|
||||
bool is_valid() const { return lower_weight >= 0.0; }
|
||||
//! Whether the weight window is in a valid state
|
||||
bool is_valid() const { return lower_weight >= 0.0; }
|
||||
|
||||
//! Adjust the weight window by a constant factor
|
||||
void scale(double factor)
|
||||
{
|
||||
lower_weight *= factor;
|
||||
upper_weight *= factor;
|
||||
}
|
||||
};
|
||||
//! Adjust the weight window by a constant factor
|
||||
void scale(double factor)
|
||||
{
|
||||
lower_weight *= factor;
|
||||
upper_weight *= factor;
|
||||
}
|
||||
};
|
||||
|
||||
//==============================================================================
|
||||
//! Weight window settings
|
||||
//==============================================================================
|
||||
//==============================================================================
|
||||
//! Weight window settings
|
||||
//==============================================================================
|
||||
|
||||
class WeightWindows {
|
||||
public:
|
||||
//----------------------------------------------------------------------------
|
||||
// Constructors
|
||||
WeightWindows(int32_t id = -1);
|
||||
WeightWindows(pugi::xml_node node);
|
||||
~WeightWindows();
|
||||
static WeightWindows* create(int32_t id = -1);
|
||||
static WeightWindows* from_hdf5(
|
||||
hid_t wws_group, const std::string& group_name);
|
||||
class WeightWindows {
|
||||
public:
|
||||
//----------------------------------------------------------------------------
|
||||
// Constructors
|
||||
WeightWindows(int32_t id = -1);
|
||||
WeightWindows(pugi::xml_node node);
|
||||
~WeightWindows();
|
||||
static WeightWindows* create(int32_t id = -1);
|
||||
static WeightWindows* from_hdf5(
|
||||
hid_t wws_group, const std::string& group_name);
|
||||
|
||||
//----------------------------------------------------------------------------
|
||||
// Methods
|
||||
private:
|
||||
template<class T>
|
||||
void check_bounds(const T& lower, const T& upper) const;
|
||||
//----------------------------------------------------------------------------
|
||||
// Methods
|
||||
private:
|
||||
template<class T>
|
||||
void check_bounds(const T& lower, const T& upper) const;
|
||||
|
||||
template<class T>
|
||||
void check_bounds(const T& lower) const;
|
||||
template<class T>
|
||||
void check_bounds(const T& lower) const;
|
||||
|
||||
void check_tally_update_compatibility(const Tally* tally);
|
||||
void check_tally_update_compatibility(const Tally* tally);
|
||||
|
||||
public:
|
||||
//! Set the weight window ID
|
||||
void set_id(int32_t id = -1);
|
||||
public:
|
||||
//! Set the weight window ID
|
||||
void set_id(int32_t id = -1);
|
||||
|
||||
void set_energy_bounds(gsl::span<const double> bounds);
|
||||
void set_energy_bounds(gsl::span<const double> bounds);
|
||||
|
||||
void set_mesh(const std::unique_ptr<Mesh>& mesh);
|
||||
void set_mesh(const std::unique_ptr<Mesh>& mesh);
|
||||
|
||||
void set_mesh(const Mesh* mesh);
|
||||
void set_mesh(const Mesh* mesh);
|
||||
|
||||
void set_mesh(int32_t mesh_idx);
|
||||
void set_mesh(int32_t mesh_idx);
|
||||
|
||||
//! Ready the weight window class for use
|
||||
void set_defaults();
|
||||
//! Ready the weight window class for use
|
||||
void set_defaults();
|
||||
|
||||
//! Update weight window boundaries using tally results
|
||||
//! \param[in] tally Pointer to the tally whose results will be used to
|
||||
//! update weight windows \param[in] value String representing the type of
|
||||
//! value to use for weight window generation (one of "mean" or "rel_err")
|
||||
//! \param[in] threshold Relative error threshold. Results over this
|
||||
//! threshold will be ignored \param[in] ratio Ratio of upper to lower
|
||||
//! weight window bounds
|
||||
void update_magic(const Tally* tally, const std::string& value = "mean",
|
||||
double threshold = 1.0, double ratio = 5.0);
|
||||
//! Ensure the weight window lower bounds are properly allocated
|
||||
void allocate_ww_bounds();
|
||||
|
||||
// NOTE: This is unused for now but may be used in the future
|
||||
//! Write weight window settings to an HDF5 file
|
||||
//! \param[in] group HDF5 group to write to
|
||||
void to_hdf5(hid_t group) const;
|
||||
//! Update weight window boundaries using tally results
|
||||
//! \param[in] tally Pointer to the tally whose results will be used to
|
||||
//! update weight windows \param[in] value String representing the type of
|
||||
//! value to use for weight window generation (one of "mean" or "rel_err")
|
||||
//! \param[in] threshold Relative error threshold. Results over this
|
||||
//! threshold will be ignored \param[in] ratio Ratio of upper to lower
|
||||
//! weight window bounds
|
||||
void update_magic(const Tally* tally, const std::string& value = "mean",
|
||||
double threshold = 1.0, double ratio = 5.0);
|
||||
|
||||
//! Retrieve the weight window for a particle
|
||||
//! \param[in] p Particle to get weight window for
|
||||
WeightWindow get_weight_window(const Particle& p) const;
|
||||
// NOTE: This is unused for now but may be used in the future
|
||||
//! Write weight window settings to an HDF5 file
|
||||
//! \param[in] group HDF5 group to write to
|
||||
void to_hdf5(hid_t group) const;
|
||||
|
||||
std::array<int, 2> bounds_size() const;
|
||||
//! Retrieve the weight window for a particle
|
||||
//! \param[in] p Particle to get weight window for
|
||||
WeightWindow get_weight_window(const Particle& p) const;
|
||||
|
||||
const vector<double>& energy_bounds() const { return energy_bounds_; }
|
||||
std::array<int, 2> bounds_size() const;
|
||||
|
||||
void set_bounds(const xt::xtensor<double, 2>& lower_ww_bounds,
|
||||
const xt::xtensor<double, 2>& upper_bounds);
|
||||
const vector<double>& energy_bounds() const { return energy_bounds_; }
|
||||
|
||||
void set_bounds(const xt::xtensor<double, 2>& lower_bounds, double ratio);
|
||||
void set_bounds(const xt::xtensor<double, 2>& lower_ww_bounds,
|
||||
const xt::xtensor<double, 2>& upper_bounds);
|
||||
|
||||
void set_bounds(gsl::span<const double> lower_bounds,
|
||||
gsl::span<const double> upper_bounds);
|
||||
void set_bounds(const xt::xtensor<double, 2>& lower_bounds, double ratio);
|
||||
|
||||
void set_bounds(gsl::span<const double> lower_bounds, double ratio);
|
||||
void set_bounds(
|
||||
gsl::span<const double> lower_bounds, gsl::span<const double> upper_bounds);
|
||||
|
||||
void set_particle_type(ParticleType p_type);
|
||||
void set_bounds(gsl::span<const double> lower_bounds, double ratio);
|
||||
|
||||
//----------------------------------------------------------------------------
|
||||
// Accessors
|
||||
int32_t id() const { return id_; }
|
||||
int32_t& id() { return id_; }
|
||||
void set_particle_type(ParticleType p_type);
|
||||
|
||||
int32_t index() const { return index_; }
|
||||
//----------------------------------------------------------------------------
|
||||
// Accessors
|
||||
int32_t id() const { return id_; }
|
||||
int32_t& id() { return id_; }
|
||||
|
||||
vector<double>& energy_bounds() { return energy_bounds_; }
|
||||
int32_t index() const { return index_; }
|
||||
|
||||
const std::unique_ptr<Mesh>& mesh() const
|
||||
{
|
||||
return model::meshes[mesh_idx_];
|
||||
}
|
||||
vector<double>& energy_bounds() { return energy_bounds_; }
|
||||
|
||||
const xt::xtensor<double, 2>& lower_ww_bounds() const { return lower_ww_; }
|
||||
xt::xtensor<double, 2>& lower_ww_bounds() { return lower_ww_; }
|
||||
const std::unique_ptr<Mesh>& mesh() const { return model::meshes[mesh_idx_]; }
|
||||
|
||||
const xt::xtensor<double, 2>& upper_ww_bounds() const { return upper_ww_; }
|
||||
xt::xtensor<double, 2>& upper_ww_bounds() { return upper_ww_; }
|
||||
const xt::xtensor<double, 2>& lower_ww_bounds() const { return lower_ww_; }
|
||||
xt::xtensor<double, 2>& lower_ww_bounds() { return lower_ww_; }
|
||||
|
||||
ParticleType particle_type() const { return particle_type_; }
|
||||
const xt::xtensor<double, 2>& upper_ww_bounds() const { return upper_ww_; }
|
||||
xt::xtensor<double, 2>& upper_ww_bounds() { return upper_ww_; }
|
||||
|
||||
private:
|
||||
//----------------------------------------------------------------------------
|
||||
// Data members
|
||||
int32_t id_; //!< Unique ID
|
||||
gsl::index index_; //!< Index into weight windows vector
|
||||
ParticleType particle_type_ {
|
||||
ParticleType::neutron}; //!< Particle type to apply weight windows to
|
||||
vector<double> energy_bounds_; //!< Energy boundaries [eV]
|
||||
xt::xtensor<double, 2> lower_ww_; //!< Lower weight window bounds (shape:
|
||||
//!< energy_bins, mesh_bins (k, j, i))
|
||||
xt::xtensor<double, 2>
|
||||
upper_ww_; //!< Upper weight window bounds (shape: energy_bins, mesh_bins)
|
||||
double survival_ratio_ {3.0}; //!< Survival weight ratio
|
||||
double max_lb_ratio_ {
|
||||
1.0}; //!< Maximum lower bound to particle weight ratio
|
||||
double weight_cutoff_ {DEFAULT_WEIGHT_CUTOFF}; //!< Weight cutoff
|
||||
int max_split_ {10}; //!< Maximum value for particle splitting
|
||||
int32_t mesh_idx_; //!< Index in meshes vector
|
||||
};
|
||||
ParticleType particle_type() const { return particle_type_; }
|
||||
|
||||
class WeightWindowsGenerator {
|
||||
public:
|
||||
// Constructors
|
||||
WeightWindowsGenerator(pugi::xml_node node);
|
||||
private:
|
||||
//----------------------------------------------------------------------------
|
||||
// Data members
|
||||
int32_t id_; //!< Unique ID
|
||||
gsl::index index_; //!< Index into weight windows vector
|
||||
ParticleType particle_type_ {
|
||||
ParticleType::neutron}; //!< Particle type to apply weight windows to
|
||||
vector<double> energy_bounds_; //!< Energy boundaries [eV]
|
||||
xt::xtensor<double, 2> lower_ww_; //!< Lower weight window bounds (shape:
|
||||
//!< energy_bins, mesh_bins (k, j, i))
|
||||
xt::xtensor<double, 2>
|
||||
upper_ww_; //!< Upper weight window bounds (shape: energy_bins, mesh_bins)
|
||||
double survival_ratio_ {3.0}; //!< Survival weight ratio
|
||||
double max_lb_ratio_ {1.0}; //!< Maximum lower bound to particle weight ratio
|
||||
double weight_cutoff_ {DEFAULT_WEIGHT_CUTOFF}; //!< Weight cutoff
|
||||
int max_split_ {10}; //!< Maximum value for particle splitting
|
||||
int32_t mesh_idx_ {-1}; //!< Index in meshes vector
|
||||
};
|
||||
|
||||
// Methods
|
||||
void update() const;
|
||||
class WeightWindowsGenerator {
|
||||
public:
|
||||
// Constructors
|
||||
WeightWindowsGenerator(pugi::xml_node node);
|
||||
|
||||
// Data members
|
||||
int32_t
|
||||
tally_idx_; //!< Index of the tally used to update the weight windows
|
||||
int32_t ww_idx_; //!< Index of the weight windows object being generated
|
||||
std::string method_; //!< Method used to update weight window. Only "magic"
|
||||
//!< is valid for now.
|
||||
int32_t max_realizations_; //!< Maximum number of tally realizations
|
||||
int32_t update_interval_; //!< Determines how often updates occur
|
||||
bool on_the_fly_; //!< Whether or not weight windows
|
||||
// Methods
|
||||
void update() const;
|
||||
|
||||
// MAGIC update parameters
|
||||
std::string tally_value_ {
|
||||
"mean"}; //<! Tally value to use (one of {"mean", "rel_err"})
|
||||
double threshold_ {1.0}; //<! Relative error threshold for values used to
|
||||
//update weight windows
|
||||
double ratio_ {5.0}; //<! ratio of lower to upper weight window bounds
|
||||
};
|
||||
//! Create the tally used for weight window generation
|
||||
void create_tally();
|
||||
|
||||
// Data members
|
||||
int32_t tally_idx_; //!< Index of the tally used to update the weight windows
|
||||
int32_t ww_idx_; //!< Index of the weight windows object being generated
|
||||
std::string method_; //!< Method used to update weight window. Only "magic"
|
||||
//!< is valid for now.
|
||||
int32_t max_realizations_; //!< Maximum number of tally realizations
|
||||
int32_t update_interval_; //!< Determines how often updates occur
|
||||
bool on_the_fly_; //!< Whether or not to keep tally results between batches or
|
||||
//!< realizations
|
||||
|
||||
// MAGIC update parameters
|
||||
std::string tally_value_ {
|
||||
"mean"}; //<! Tally value to use (one of {"mean", "rel_err"})
|
||||
double threshold_ {1.0}; //<! Relative error threshold for values used to
|
||||
// update weight windows
|
||||
double ratio_ {5.0}; //<! ratio of lower to upper weight window bounds
|
||||
};
|
||||
|
||||
//! Finalize variance reduction objects after all inputs have been read
|
||||
void finalize_variance_reduction();
|
||||
|
||||
} // namespace openmc
|
||||
#endif // OPENMC_WEIGHT_WINDOWS_H
|
||||
|
|
|
|||
|
|
@ -49,7 +49,8 @@ class WeightWindows(IDManagerMixin):
|
|||
Ratio of the lower to upper weight window bounds
|
||||
energy_bounds : Iterable of Real
|
||||
A list of values for which each successive pair constitutes a range of
|
||||
energies in [eV] for a single bin
|
||||
energies in [eV] for a single bin. If no energy bins are provided, the
|
||||
maximum and minimum energy for the data available at runtime.
|
||||
particle_type : {'neutron', 'photon'}
|
||||
Particle type the weight windows apply to
|
||||
survival_ratio : float
|
||||
|
|
@ -125,7 +126,9 @@ class WeightWindows(IDManagerMixin):
|
|||
self.mesh = mesh
|
||||
self.id = id
|
||||
self.particle_type = particle_type
|
||||
self.energy_bounds = energy_bounds
|
||||
self._energy_bounds = None
|
||||
if energy_bounds is not None:
|
||||
self.energy_bounds = energy_bounds
|
||||
self.lower_ww_bounds = lower_ww_bounds
|
||||
|
||||
if upper_ww_bounds is not None and upper_bound_ratio:
|
||||
|
|
@ -229,7 +232,7 @@ class WeightWindows(IDManagerMixin):
|
|||
@property
|
||||
def num_energy_bins(self) -> int:
|
||||
if self.energy_bounds is None:
|
||||
raise ValueError('Energy bounds are not set')
|
||||
return 1
|
||||
return self.energy_bounds.size - 1
|
||||
|
||||
@property
|
||||
|
|
@ -650,7 +653,8 @@ class WeightWindowGenerator:
|
|||
Mesh used to represent the weight windows spatially
|
||||
energy_bounds : Iterable of Real
|
||||
A list of values for which each successive pair constitutes a range of
|
||||
energies in [eV] for a single bin
|
||||
energies in [eV] for a single bin. If no energy bins are provided, the
|
||||
maximum and minimum energy for the data available at runtime.
|
||||
particle_type : {'neutron', 'photon'}
|
||||
Particle type the weight windows apply to
|
||||
|
||||
|
|
@ -681,6 +685,7 @@ class WeightWindowGenerator:
|
|||
|
||||
def __init__(self, mesh, energy_bounds=None, particle_type='neutron'):
|
||||
self.mesh = mesh
|
||||
self._energy_bounds = None
|
||||
if energy_bounds is not None:
|
||||
self.energy_bounds = energy_bounds
|
||||
self.particle_type = particle_type
|
||||
|
|
|
|||
|
|
@ -33,6 +33,7 @@
|
|||
#include "openmc/thermal.h"
|
||||
#include "openmc/timer.h"
|
||||
#include "openmc/vector.h"
|
||||
#include "openmc/weight_windows.h"
|
||||
|
||||
#ifdef LIBMESH
|
||||
#include "libmesh/libmesh.h"
|
||||
|
|
@ -404,6 +405,8 @@ bool read_model_xml()
|
|||
}
|
||||
}
|
||||
|
||||
finalize_variance_reduction();
|
||||
|
||||
return true;
|
||||
}
|
||||
|
||||
|
|
@ -427,6 +430,8 @@ void read_separate_xml_files()
|
|||
// Read the plots.xml regardless of plot mode in case plots are requested
|
||||
// via the API
|
||||
read_plots_xml();
|
||||
|
||||
finalize_variance_reduction();
|
||||
}
|
||||
|
||||
void initial_output()
|
||||
|
|
|
|||
|
|
@ -257,21 +257,26 @@ WeightWindows* WeightWindows::from_hdf5(
|
|||
|
||||
void WeightWindows::set_defaults()
|
||||
{
|
||||
// ensure default values are set
|
||||
// set energy bounds to the min/max energy supported by the data
|
||||
if (energy_bounds_.size() == 0) {
|
||||
int p_type = static_cast<int>(particle_type_);
|
||||
energy_bounds_.push_back(data::energy_min[p_type]);
|
||||
energy_bounds_.push_back(data::energy_max[p_type]);
|
||||
}
|
||||
}
|
||||
|
||||
// some constructors won't allocate space for the bounds
|
||||
// do that here so the object is valid
|
||||
if (lower_ww_.size() == 0 || upper_ww_.size() == 0) {
|
||||
lower_ww_ = xt::empty<double>(bounds_size());
|
||||
lower_ww_.fill(-1);
|
||||
upper_ww_ = xt::empty<double>(bounds_size());
|
||||
upper_ww_.fill(-1);
|
||||
void WeightWindows::allocate_ww_bounds()
|
||||
{
|
||||
auto shape = bounds_size();
|
||||
if (shape[0] * shape[1] == 0) {
|
||||
auto msg = fmt::format(
|
||||
"Size of weight window bounds is zero for WeightWindows {}", id());
|
||||
warning(msg);
|
||||
}
|
||||
lower_ww_ = xt::empty<double>(shape);
|
||||
lower_ww_.fill(-1);
|
||||
upper_ww_ = xt::empty<double>(shape);
|
||||
upper_ww_.fill(-1);
|
||||
}
|
||||
|
||||
void WeightWindows::set_id(int32_t id)
|
||||
|
|
@ -308,6 +313,9 @@ void WeightWindows::set_energy_bounds(gsl::span<const double> bounds)
|
|||
{
|
||||
energy_bounds_.clear();
|
||||
energy_bounds_.insert(energy_bounds_.begin(), bounds.begin(), bounds.end());
|
||||
// if the mesh is set, allocate space for weight window bounds
|
||||
if (mesh_idx_ != C_NONE)
|
||||
allocate_ww_bounds();
|
||||
}
|
||||
|
||||
void WeightWindows::set_particle_type(ParticleType p_type)
|
||||
|
|
@ -325,6 +333,7 @@ void WeightWindows::set_mesh(int32_t mesh_idx)
|
|||
fatal_error(fmt::format("Could not find a mesh for index {}", mesh_idx));
|
||||
|
||||
mesh_idx_ = mesh_idx;
|
||||
allocate_ww_bounds();
|
||||
}
|
||||
|
||||
void WeightWindows::set_mesh(const std::unique_ptr<Mesh>& mesh)
|
||||
|
|
@ -743,43 +752,6 @@ WeightWindowsGenerator::WeightWindowsGenerator(pugi::xml_node node)
|
|||
e_bounds.push_back(data::energy_max[p_type]);
|
||||
}
|
||||
|
||||
// create a tally based on the WWG information
|
||||
Tally* ww_tally = Tally::create();
|
||||
tally_idx_ = model::tally_map[ww_tally->id()];
|
||||
ww_tally->set_scores({"flux"});
|
||||
|
||||
// see if there's already a mesh filter using this mesh
|
||||
bool found_mesh_filter = false;
|
||||
for (const auto& f : model::tally_filters) {
|
||||
if (f->type() == FilterType::MESH) {
|
||||
const auto* mesh_filter = dynamic_cast<MeshFilter*>(f.get());
|
||||
if (mesh_filter->mesh() == mesh_idx && !mesh_filter->translated()) {
|
||||
ww_tally->add_filter(f.get());
|
||||
found_mesh_filter = true;
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if (!found_mesh_filter) {
|
||||
auto mesh_filter = Filter::create("mesh");
|
||||
openmc_mesh_filter_set_mesh(mesh_filter->index(), model::mesh_map[mesh_id]);
|
||||
ww_tally->add_filter(mesh_filter);
|
||||
}
|
||||
|
||||
if (e_bounds.size() > 0) {
|
||||
auto energy_filter = Filter::create("energy");
|
||||
openmc_energy_filter_set_bins(
|
||||
energy_filter->index(), e_bounds.size(), e_bounds.data());
|
||||
ww_tally->add_filter(energy_filter);
|
||||
}
|
||||
|
||||
// add a particle filter
|
||||
auto particle_filter = Filter::create("particle");
|
||||
auto pf = dynamic_cast<ParticleFilter*>(particle_filter);
|
||||
pf->set_particles({&particle_type, 1});
|
||||
ww_tally->add_filter(particle_filter);
|
||||
|
||||
// set method and parameters for updates
|
||||
method_ = get_node_value(node, "method");
|
||||
if (method_ == "magic") {
|
||||
|
|
@ -815,13 +787,59 @@ WeightWindowsGenerator::WeightWindowsGenerator(pugi::xml_node node)
|
|||
// create a matching weight windows object
|
||||
auto wws = WeightWindows::create();
|
||||
ww_idx_ = wws->index();
|
||||
wws->set_mesh(mesh_idx);
|
||||
if (e_bounds.size() > 0)
|
||||
wws->set_energy_bounds(e_bounds);
|
||||
wws->set_mesh(model::mesh_map[mesh_id]);
|
||||
wws->set_particle_type(particle_type);
|
||||
wws->set_defaults();
|
||||
}
|
||||
|
||||
void WeightWindowsGenerator::create_tally()
|
||||
{
|
||||
const auto& wws = variance_reduction::weight_windows[ww_idx_];
|
||||
|
||||
// create a tally based on the WWG information
|
||||
Tally* ww_tally = Tally::create();
|
||||
tally_idx_ = model::tally_map[ww_tally->id()];
|
||||
ww_tally->set_scores({"flux"});
|
||||
|
||||
int32_t mesh_id = wws->mesh()->id();
|
||||
int32_t mesh_idx = model::mesh_map.at(mesh_id);
|
||||
// see if there's already a mesh filter using this mesh
|
||||
bool found_mesh_filter = false;
|
||||
for (const auto& f : model::tally_filters) {
|
||||
if (f->type() == FilterType::MESH) {
|
||||
const auto* mesh_filter = dynamic_cast<MeshFilter*>(f.get());
|
||||
if (mesh_filter->mesh() == mesh_idx && !mesh_filter->translated()) {
|
||||
ww_tally->add_filter(f.get());
|
||||
found_mesh_filter = true;
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if (!found_mesh_filter) {
|
||||
auto mesh_filter = Filter::create("mesh");
|
||||
openmc_mesh_filter_set_mesh(mesh_filter->index(), model::mesh_map[mesh_id]);
|
||||
ww_tally->add_filter(mesh_filter);
|
||||
}
|
||||
|
||||
const auto& e_bounds = wws->energy_bounds();
|
||||
if (e_bounds.size() > 0) {
|
||||
auto energy_filter = Filter::create("energy");
|
||||
openmc_energy_filter_set_bins(
|
||||
energy_filter->index(), e_bounds.size(), e_bounds.data());
|
||||
ww_tally->add_filter(energy_filter);
|
||||
}
|
||||
|
||||
// add a particle filter
|
||||
auto particle_type = wws->particle_type();
|
||||
auto particle_filter = Filter::create("particle");
|
||||
auto pf = dynamic_cast<ParticleFilter*>(particle_filter);
|
||||
pf->set_particles({&particle_type, 1});
|
||||
ww_tally->add_filter(particle_filter);
|
||||
}
|
||||
|
||||
void WeightWindowsGenerator::update() const
|
||||
{
|
||||
const auto& wws = variance_reduction::weight_windows[ww_idx_];
|
||||
|
|
@ -845,6 +863,17 @@ void WeightWindowsGenerator::update() const
|
|||
// complete
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
// Non-member functions
|
||||
//==============================================================================
|
||||
|
||||
void finalize_variance_reduction()
|
||||
{
|
||||
for (const auto& wwg : variance_reduction::weight_windows_generators) {
|
||||
wwg->create_tally();
|
||||
}
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
// C API
|
||||
//==============================================================================
|
||||
|
|
|
|||
|
|
@ -28,13 +28,20 @@ def test_ww_generator(run_in_tmpdir):
|
|||
energy_bounds = np.linspace(0.0, 1e6, 70)
|
||||
particle = 'neutron'
|
||||
|
||||
# include another tally to make sure user-specified tallies and those automaticaly
|
||||
# created by weight window generators can coexist
|
||||
tally = openmc.Tally()
|
||||
ef = openmc.EnergyFilter(energy_bounds)
|
||||
tally.filters = [ef]
|
||||
tally.scores = ['flux']
|
||||
model.tallies = [tally]
|
||||
|
||||
wwg = openmc.WeightWindowGenerator(mesh, energy_bounds, particle)
|
||||
wwg.update_parameters = {'ratio' : 5.0, 'threshold': 0.8, 'value' : 'mean'}
|
||||
wwg.update_parameters = {'ratio': 5.0, 'threshold': 0.8, 'value': 'mean'}
|
||||
|
||||
model.settings.weight_window_generators = wwg
|
||||
model.export_to_xml()
|
||||
|
||||
model.run()
|
||||
|
||||
# we test the effectiveness of the update method elsewhere, so
|
||||
# just test that the generation happens successfully here
|
||||
assert os.path.exists('weight_windows.h5')
|
||||
|
|
|
|||
|
|
@ -241,7 +241,21 @@ def test_roundtrip(run_in_tmpdir, model, wws):
|
|||
assert(ww_out == ww_in)
|
||||
|
||||
|
||||
def test_ww_attrs(run_in_tmpdir, model):
|
||||
def test_ww_attrs_python(model):
|
||||
mesh = openmc.RegularMesh.from_domain(model.geometry)
|
||||
lower_bounds = np.ones(mesh.dimension)
|
||||
|
||||
# ensure that creation of weight window objects with default arg values
|
||||
# is successful
|
||||
wws = openmc.WeightWindows(mesh, lower_bounds, upper_bound_ratio=10.0)
|
||||
|
||||
assert wws.energy_bounds == None
|
||||
|
||||
wwg = openmc.WeightWindowGenerator(mesh)
|
||||
|
||||
assert wwg.energy_bounds == None
|
||||
|
||||
def test_ww_attrs_capi(run_in_tmpdir, model):
|
||||
model.export_to_xml()
|
||||
|
||||
openmc.lib.init()
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue