mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-27 05:35:49 -04:00
Addressing changes from PR review.
This commit is contained in:
Patrick Shriwise
parent
ca7cba57df
commit
895249b6ca
10 changed files with 129 additions and 327 deletions
|
|
@ -59,8 +59,8 @@ deterministic methods:
|
|||
Now let's look at the pros and cons of Monte Carlo methods:
|
||||
|
||||
- **Pro**: No mesh generation is required to build geometry. By using
|
||||
`constructive solid geometry`_, it's possible to build arbitrarily complex
|
||||
reactor models with curved surfaces.
|
||||
`constructive solid geometry`_, it's possible to build complex
|
||||
xmodels with curved surfaces.
|
||||
|
||||
- **Pro**: Monte Carlo methods can be used with either continuous-energy or
|
||||
multi-group cross sections.
|
||||
|
|
|
|||
|
|
@ -33,34 +33,6 @@
|
|||
"We'll need to download the unstructured mesh file used in this notebook. We'll be retrieving those using the function and URLs below."
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 2,
|
||||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"%matplotlib inline\n",
|
||||
"from matplotlib import pyplot as plt\n",
|
||||
"plt.rcParams[\"figure.figsize\"] = (30,10)\n",
|
||||
"\n",
|
||||
"import urllib.request\n",
|
||||
"\n",
|
||||
"pin_mesh_url = 'https://tinyurl.com/u9ce9d7' # 1.2 MB\n",
|
||||
"\n",
|
||||
"def download(url, filename='dagmc.h5m'):\n",
|
||||
" \"\"\"\n",
|
||||
" Helper function for retrieving dagmc models\n",
|
||||
" \"\"\"\n",
|
||||
" u = urllib.request.urlopen(url)\n",
|
||||
" \n",
|
||||
" if u.status != 200:\n",
|
||||
" raise RuntimeError(\"Failed to download file.\")\n",
|
||||
" \n",
|
||||
" # save file as dagmc.h5m\n",
|
||||
" with open(filename, 'wb') as f:\n",
|
||||
" f.write(u.read())"
|
||||
]
|
||||
},
|
||||
{
|
||||
"cell_type": "code",
|
||||
"execution_count": 3,
|
||||
|
|
@ -467,7 +439,7 @@
|
|||
"outputs": [],
|
||||
"source": [
|
||||
"download(pin_mesh_url, \"pins1-4.h5m\")\n",
|
||||
"umesh = openmc.UnstructuredMesh(filename=\"pins1-4.h5m\")\n",
|
||||
"umesh = openmc.UnstructuredMesh(\"pins1-4.h5m\")\n",
|
||||
"mesh_filter = openmc.MeshFilter(umesh)"
|
||||
]
|
||||
},
|
||||
|
|
@ -893,7 +865,7 @@
|
|||
"cell_type": "markdown",
|
||||
"metadata": {},
|
||||
"source": [
|
||||
"We hope you've found this example notebook useful. More unstructured mesh features are under development and will be included in additional examples soon."
|
||||
"We hope you've found this example notebook useful!"
|
||||
]
|
||||
},
|
||||
{
|
||||
|
|
|
|||
|
|
@ -351,7 +351,7 @@
|
|||
"metadata": {},
|
||||
"outputs": [],
|
||||
"source": [
|
||||
"unstructured_mesh = openmc.UnstructuredMesh(filename=\"manifold.h5m\")\n",
|
||||
"unstructured_mesh = openmc.UnstructuredMesh(\"manifold.h5m\")\n",
|
||||
"\n",
|
||||
"mesh_filter = openmc.MeshFilter(unstructured_mesh)\n",
|
||||
"\n",
|
||||
|
|
|
|||
|
|
@ -258,18 +258,10 @@ public:
|
|||
UnstructuredMesh(pugi::xml_node);
|
||||
~UnstructuredMesh() = default;
|
||||
|
||||
void bins_crossed(const Particle* p, std::vector<int>& bins,
|
||||
void bins_crossed(const Particle* p,
|
||||
std::vector<int>& bins,
|
||||
std::vector<double>& lengths) const override;
|
||||
|
||||
//! Check where a line segment intersects the mesh and if it intersects at all
|
||||
//
|
||||
//! \param[in,out] r0 In: starting position, out: intersection point
|
||||
//! \param[in] r1 Ending position
|
||||
//! \param[out] ijk Indices of the mesh bin containing the intersection point
|
||||
//! \return Whether the line segment connecting r0 and r1 intersects mesh
|
||||
bool intersects(Position& r0, Position r1, int* ijk);
|
||||
|
||||
|
||||
private:
|
||||
|
||||
//! Find all intersections with faces of the mesh.
|
||||
|
|
@ -410,12 +402,11 @@ public:
|
|||
//! Write the mesh with any current tally data
|
||||
void write(std::string base_filename) const;
|
||||
|
||||
std::string filename_; //<! Path to unstructured mesh file
|
||||
std::string filename_; //!< Path to unstructured mesh file
|
||||
|
||||
private:
|
||||
moab::Range ehs_; //!< Range of tetrahedra EntityHandle's in the mesh
|
||||
moab::EntityHandle meshset_; //!< EntitySet containing all elements
|
||||
moab::EntityHandle tet_set_; //! < EntitySet containing all tetrahedra
|
||||
moab::EntityHandle tetset_; //!< EntitySet containing all tetrahedra
|
||||
moab::EntityHandle kdtree_root_; //!< Root of the MOAB KDTree
|
||||
std::unique_ptr<moab::Interface> mbi_; //!< MOAB instance
|
||||
std::unique_ptr<moab::AdaptiveKDTree> kdtree_; //!< MOAB KDTree instance
|
||||
|
|
|
|||
|
|
@ -619,7 +619,7 @@ class UnstructuredMesh(MeshBase):
|
|||
(1.0, 1.0, 1.0), ...]
|
||||
"""
|
||||
|
||||
def __init__(self, mesh_id=None, name='', filename=''):
|
||||
def __init__(self, filename, mesh_id=None, name=''):
|
||||
super().__init__(mesh_id, name)
|
||||
self._filename = filename
|
||||
self._volumes = []
|
||||
|
|
@ -631,11 +631,9 @@ class UnstructuredMesh(MeshBase):
|
|||
|
||||
@filename.setter
|
||||
def filename(self, filename):
|
||||
if filename is not None:
|
||||
cv.check_type('Unstructured Mesh filename', filename, str)
|
||||
self._filename = filename
|
||||
else:
|
||||
self.filename = ''
|
||||
cv.check_type('Unstructured Mesh filename: {}'.format(filename),
|
||||
filename, str)
|
||||
self._filename = filename
|
||||
|
||||
@property
|
||||
def volumes(self):
|
||||
|
|
@ -667,9 +665,9 @@ class UnstructuredMesh(MeshBase):
|
|||
@classmethod
|
||||
def from_hdf5(cls, group):
|
||||
mesh_id = int(group.name.split('/')[-1].lstrip('mesh '))
|
||||
filename = group['filename'][()].decode()
|
||||
|
||||
mesh = cls(mesh_id)
|
||||
mesh.filename = group['filename'][()].decode()
|
||||
mesh = cls(filename, mesh_id=mesh_id)
|
||||
vol_data = group['volumes'][()]
|
||||
centroids = group['centroids'][()]
|
||||
mesh.volumes = np.reshape(vol_data, (vol_data.shape[0],))
|
||||
|
|
@ -695,3 +693,24 @@ class UnstructuredMesh(MeshBase):
|
|||
subelement.text = self.filename
|
||||
|
||||
return element
|
||||
|
||||
@classmethod
|
||||
def from_xml_element(cls, elem):
|
||||
"""Generate unstructured mesh object from XML element
|
||||
|
||||
Parameters
|
||||
----------
|
||||
elem : xml.etree.ElementTree.Element
|
||||
XML element
|
||||
|
||||
Returns
|
||||
-------
|
||||
openmc.UnstructuredMesh
|
||||
UnstructuredMesh generated from an XML element
|
||||
"""
|
||||
mesh_id = int(get_text(elem, 'id'))
|
||||
filename = get_text(elem, 'mesh_file')
|
||||
|
||||
mesh = cls(filename, mesh_id)
|
||||
|
||||
return mesh
|
||||
|
|
|
|||
|
|
@ -310,7 +310,7 @@ class StatePoint(object):
|
|||
if self.run_mode == 'eigenvalue':
|
||||
return self._f['n_inactive'][()]
|
||||
else:
|
||||
return None
|
||||
return None
|
||||
|
||||
@property
|
||||
def n_particles(self):
|
||||
|
|
|
|||
57
src/mesh.cpp
57
src/mesh.cpp
|
|
@ -3,7 +3,6 @@
|
|||
#include <algorithm> // for copy, equal, min, min_element
|
||||
#include <cstddef> // for size_t
|
||||
#include <cmath> // for ceil
|
||||
#include <fmt/core.h> // for fmt
|
||||
#include <memory> // for allocator
|
||||
#include <string>
|
||||
|
||||
|
|
@ -16,6 +15,7 @@
|
|||
#include "xtensor/xsort.hpp"
|
||||
#include "xtensor/xtensor.hpp"
|
||||
#include "xtensor/xview.hpp"
|
||||
#include <fmt/core.h> // for fmt
|
||||
|
||||
#include "openmc/capi.h"
|
||||
#include "openmc/constants.h"
|
||||
|
|
@ -1456,10 +1456,6 @@ openmc_mesh_set_id(int32_t index, int32_t id)
|
|||
extern "C" int
|
||||
openmc_mesh_get_dimension(int32_t index, int** dims, int* n)
|
||||
{
|
||||
if (index < 0 || index >= model::meshes.size()) {
|
||||
set_errmsg("Index in meshes array is out of bounds.");
|
||||
return OPENMC_E_OUT_OF_BOUNDS;
|
||||
}
|
||||
RegularMesh* mesh;
|
||||
if (int err = check_regular_mesh(index, &mesh)) return err;
|
||||
*dims = mesh->shape_.data();
|
||||
|
|
@ -1471,11 +1467,6 @@ openmc_mesh_get_dimension(int32_t index, int** dims, int* n)
|
|||
extern "C" int
|
||||
openmc_mesh_set_dimension(int32_t index, int n, const int* dims)
|
||||
{
|
||||
if (index < 0 || index >= model::meshes.size()) {
|
||||
set_errmsg("Index in meshes array is out of bounds.");
|
||||
return OPENMC_E_OUT_OF_BOUNDS;
|
||||
}
|
||||
|
||||
RegularMesh* mesh;
|
||||
if (int err = check_regular_mesh(index, &mesh)) return err;
|
||||
|
||||
|
|
@ -1559,19 +1550,14 @@ UnstructuredMesh::UnstructuredMesh(pugi::xml_node node) : Mesh(node)
|
|||
|
||||
// create MOAB instance
|
||||
mbi_ = std::make_unique<moab::Core>();
|
||||
// create meshset to load mesh into
|
||||
moab::ErrorCode rval = mbi_->create_meshset(moab::MESHSET_SET, meshset_);
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
fatal_error("Failed to create fileset for umesh: " + filename_);
|
||||
}
|
||||
// load unstructured mesh file
|
||||
rval = mbi_->load_file(filename_.c_str(), &meshset_);
|
||||
moab::ErrorCode rval = mbi_->load_file(filename_.c_str());
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
fatal_error("Failed to load the unstructured mesh file: " + filename_);
|
||||
}
|
||||
|
||||
// set member range of tetrahedral entities
|
||||
rval = mbi_->get_entities_by_dimension(meshset_, n_dimension_, ehs_);
|
||||
rval = mbi_->get_entities_by_dimension(0, n_dimension_, ehs_);
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
fatal_error("Failed to get all tetrahedral elements");
|
||||
}
|
||||
|
|
@ -1583,12 +1569,12 @@ UnstructuredMesh::UnstructuredMesh(pugi::xml_node node) : Mesh(node)
|
|||
|
||||
// make an entity set for all tetrahedra
|
||||
// this is used for convenience later in output
|
||||
rval = mbi_->create_meshset(moab::MESHSET_SET, tet_set_);
|
||||
rval = mbi_->create_meshset(moab::MESHSET_SET, tetset_);
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
fatal_error("Failed to create an entity set for the tetrahedral elements");
|
||||
}
|
||||
|
||||
rval = mbi_->add_entities(tet_set_, ehs_);
|
||||
rval = mbi_->add_entities(tetset_, ehs_);
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
fatal_error("Failed to add tetrahedra to an entity set.");
|
||||
}
|
||||
|
|
@ -1672,7 +1658,7 @@ UnstructuredMesh::bins_crossed(const Particle* p,
|
|||
moab::ErrorCode rval;
|
||||
Position last_r{p->r_last_};
|
||||
Position r{p->r()};
|
||||
Position u{p->u()};
|
||||
Direction u{p->u()};
|
||||
u /= u.norm();
|
||||
moab::CartVect r0(last_r.x, last_r.y, last_r.z);
|
||||
moab::CartVect r1(r.x, r.y, r.z);
|
||||
|
|
@ -1680,8 +1666,8 @@ UnstructuredMesh::bins_crossed(const Particle* p,
|
|||
|
||||
double track_len = (r1 - r0).length();
|
||||
|
||||
r0 += TINY_BIT*dir;
|
||||
r1 -= TINY_BIT*dir;
|
||||
r0 += TINY_BIT * dir;
|
||||
r1 -= TINY_BIT * dir;
|
||||
|
||||
UnstructuredMeshHits hits;
|
||||
intersect_track(r0, dir, track_len, hits);
|
||||
|
|
@ -1729,7 +1715,6 @@ UnstructuredMesh::bins_crossed(const Particle* p,
|
|||
bins.push_back(get_bin_from_ent_handle(tet));
|
||||
lengths.push_back((hit.first - last_dist) / track_len);
|
||||
} else {
|
||||
// if in the loop, we should always find a tet
|
||||
warning("No tet found for location between triangle hits");
|
||||
}
|
||||
|
||||
|
|
@ -1805,7 +1790,7 @@ double UnstructuredMesh::tet_volume(moab::EntityHandle tet) const {
|
|||
}
|
||||
|
||||
moab::CartVect p[4];
|
||||
rval = mbi_->get_coords(&(conn[0]), (int)conn.size(), p[0].array());
|
||||
rval = mbi_->get_coords(conn.data(), conn.size(), p[0].array());
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
fatal_error("Failed to get tet coords");
|
||||
}
|
||||
|
|
@ -1815,7 +1800,7 @@ double UnstructuredMesh::tet_volume(moab::EntityHandle tet) const {
|
|||
|
||||
void UnstructuredMesh::surface_bins_crossed(const Particle* p, std::vector<int>& bins) const {
|
||||
// TODO: Implement triangle crossings here
|
||||
return;
|
||||
throw std::runtime_error{"Unstructured mesh surface tallies are not implemented."};
|
||||
}
|
||||
|
||||
int
|
||||
|
|
@ -1845,7 +1830,7 @@ UnstructuredMesh::compute_barycentric_data(const moab::Range& tets) {
|
|||
}
|
||||
|
||||
moab::CartVect p[4];
|
||||
rval = mbi_->get_coords(&(verts[0]), (int)verts.size(), p[0].array());
|
||||
rval = mbi_->get_coords(verts.data(), verts.size(), p[0].array());
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
fatal_error("Failed to get coordinates of a tet in umesh: " + filename_);
|
||||
}
|
||||
|
|
@ -2045,9 +2030,8 @@ UnstructuredMesh::get_score_tags(std::string score) const {
|
|||
moab::MB_TAG_DENSE|moab::MB_TAG_CREAT,
|
||||
&default_val);
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
std::stringstream msg;
|
||||
msg << "Could not create or retrieve the error tag for the score " << score
|
||||
<< " on unstructured mesh " << id_;
|
||||
auto msg = fmt::format("Could not create or retrieve the error tag for the score {}"
|
||||
" on unstructured mesh {}", score, id_);
|
||||
fatal_error(msg);
|
||||
}
|
||||
|
||||
|
|
@ -2078,18 +2062,16 @@ UnstructuredMesh::set_score_data(const std::string& score,
|
|||
// set the score value
|
||||
rval = mbi_->tag_set_data(score_tags.first, ehs_, values.data());
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
std::stringstream msg;
|
||||
msg << "Failed to set the tally value for score '" << score << "' "
|
||||
<< " on unstructured mesh " << id_;
|
||||
auto msg = fmt::format("Failed to set the tally value for score '{}' "
|
||||
"on unstructured mesh {}", score, id_);
|
||||
warning(msg);
|
||||
}
|
||||
|
||||
// set the error value
|
||||
rval = mbi_->tag_set_data(score_tags.second, ehs_, std_dev.data());
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
std::stringstream msg;
|
||||
msg << "Failed to set the tally value for score '" << score << "' "
|
||||
<< " on unstructured mesh " << id_;
|
||||
auto msg = fmt::format("Failed to set the tally error for score '{}' "
|
||||
"on unstructured mesh {}", score, id_);
|
||||
warning(msg);
|
||||
}
|
||||
}
|
||||
|
|
@ -2105,10 +2087,9 @@ UnstructuredMesh::write(std::string base_filename) const {
|
|||
// to avoid clutter from zero-value data on other
|
||||
// elements during visualization
|
||||
moab::ErrorCode rval;
|
||||
rval = mbi_->write_mesh(filename.c_str(), &tet_set_, 1);
|
||||
rval = mbi_->write_mesh(filename.c_str(), &tetset_, 1);
|
||||
if (rval != moab::MB_SUCCESS) {
|
||||
std::stringstream msg;
|
||||
msg << "Failed to write unstructured mesh " << id_;
|
||||
auto msg = fmt::format("Failed to write unstructured mesh {}", id_);
|
||||
warning(msg);
|
||||
}
|
||||
}
|
||||
|
|
|
|||
|
|
@ -689,64 +689,71 @@ void write_unstructured_mesh_results() {
|
|||
for (auto& tally : model::tallies) {
|
||||
for (auto filter_idx : tally->filters()) {
|
||||
auto& filter = model::tally_filters[filter_idx];
|
||||
if (filter->type() == "mesh") {
|
||||
// check if the filter uses an unstructured mesh
|
||||
auto mesh_filter = dynamic_cast<MeshFilter*>(filter.get());
|
||||
auto mesh_idx = mesh_filter->mesh();
|
||||
auto umesh = dynamic_cast<UnstructuredMesh*>(model::meshes[mesh_idx].get());
|
||||
if (umesh) {
|
||||
// if this tally has more than one filter, print
|
||||
// warning and skip writing the mesh
|
||||
if (tally->filters().size() > 1) {
|
||||
warning(fmt::format("Skipping unstructured mesh writing for tally "
|
||||
"{}. More than one filter is present on the tally.",
|
||||
tally->id_));
|
||||
break;
|
||||
}
|
||||
if (filter->type() != "mesh") continue;
|
||||
|
||||
int n_realizations = tally->n_realizations_;
|
||||
// check if the filter uses an unstructured mesh
|
||||
auto mesh_filter = dynamic_cast<MeshFilter*>(filter.get());
|
||||
auto mesh_idx = mesh_filter->mesh();
|
||||
auto umesh = dynamic_cast<UnstructuredMesh*>(model::meshes[mesh_idx].get());
|
||||
|
||||
// write each score/nuclide combination for this tally
|
||||
for (int i_score = 0; i_score < tally->scores_.size(); i_score++) {
|
||||
for (int i_nuc = 0; i_nuc < tally->nuclides_.size(); i_nuc++) {
|
||||
if (!umesh) continue;
|
||||
|
||||
// index for this nuclide and score
|
||||
int nuc_score_idx = i_score + i_nuc*tally->scores_.size();
|
||||
// if this tally has more than one filter, print
|
||||
// warning and skip writing the mesh
|
||||
if (tally->filters().size() > 1) {
|
||||
warning(fmt::format("Skipping unstructured mesh writing for tally "
|
||||
"{0}. More than one filter is present on the tally.",
|
||||
tally->id_));
|
||||
break;
|
||||
}
|
||||
|
||||
// construct result vectors
|
||||
std::vector<double> mean_vec, std_dev_vec;
|
||||
for (int j = 0; j < tally->results_.shape()[0]; j++) {
|
||||
double mean = tally->results_(j, nuc_score_idx, TallyResult::SUM) / n_realizations;
|
||||
double sum_sq = tally->results_(j , nuc_score_idx, TallyResult::SUM_SQ);
|
||||
double std_dev = sum_sq / n_realizations - std::pow(mean, 2) / (n_realizations - 1);
|
||||
std_dev_vec.push_back(std_dev);
|
||||
mean_vec.push_back(mean);
|
||||
}
|
||||
int n_realizations = tally->n_realizations_;
|
||||
|
||||
// generate a name for the value
|
||||
std::string nuclide_name = "total"; // start with total by default
|
||||
if (tally->nuclides_[i_nuc] > -1) {
|
||||
nuclide_name = data::nuclides[tally->nuclides_[i_nuc]]->name_;
|
||||
}
|
||||
// write each score/nuclide combination for this tally
|
||||
for (int i_score = 0; i_score < tally->scores_.size(); i_score++) {
|
||||
for (int i_nuc = 0; i_nuc < tally->nuclides_.size(); i_nuc++) {
|
||||
|
||||
std::string score_name = tally->score_name(i_score);
|
||||
// index for this nuclide and score
|
||||
int nuc_score_idx = i_score + i_nuc * tally->scores_.size();
|
||||
|
||||
auto score_str = fmt::format("{}_{}", score_name, nuclide_name);
|
||||
umesh->set_score_data(score_str, mean_vec, std_dev_vec);
|
||||
// construct result vectors
|
||||
std::vector<double> mean_vec, std_dev_vec;
|
||||
for (int j = 0; j < tally->results_.shape()[0]; j++) {
|
||||
// mean
|
||||
double mean = tally->results_(j, nuc_score_idx, TallyResult::SUM) / n_realizations;
|
||||
mean_vec.push_back(mean);
|
||||
// std. dev.
|
||||
double sum_sq = tally->results_(j , nuc_score_idx, TallyResult::SUM_SQ);
|
||||
if (n_realizations > 1) {
|
||||
double std_dev = sum_sq / n_realizations - (mean* mean);
|
||||
std_dev = std::sqrt(std_dev / (n_realizations - 1));
|
||||
std_dev_vec.push_back(std_dev);
|
||||
} else {
|
||||
std_dev_vec.push_back(0.0);
|
||||
}
|
||||
}
|
||||
|
||||
// Generate a file name based on the tally id
|
||||
// and the current batch number
|
||||
int w = std::to_string(settings::n_max_batches).size();
|
||||
std::string filename = fmt::format("tally_{0}.{1:0{2}}",
|
||||
tally->id_,
|
||||
simulation::current_batch,
|
||||
w);
|
||||
// Write the unstructured mesh and data to file
|
||||
umesh->write(filename);
|
||||
// generate a name for the value
|
||||
std::string nuclide_name = "total"; // start with total by default
|
||||
if (tally->nuclides_[i_nuc] > -1) {
|
||||
nuclide_name = data::nuclides[tally->nuclides_[i_nuc]]->name_;
|
||||
}
|
||||
|
||||
std::string score_name = tally->score_name(i_score);
|
||||
auto score_str = fmt::format("{}_{}", score_name, nuclide_name);
|
||||
umesh->set_score_data(score_str, mean_vec, std_dev_vec);
|
||||
}
|
||||
}
|
||||
|
||||
// Generate a file name based on the tally id
|
||||
// and the current batch number
|
||||
int w = std::to_string(settings::n_max_batches).size();
|
||||
std::string filename = fmt::format("tally_{0}.{1:0{2}}",
|
||||
tally->id_,
|
||||
simulation::current_batch,
|
||||
w);
|
||||
// Write the unstructured mesh and data to file
|
||||
umesh->write(filename);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
|
|
|||
|
|
@ -65,171 +65,6 @@ double global_tally_collision;
|
|||
double global_tally_tracklength;
|
||||
double global_tally_leakage;
|
||||
|
||||
std::string
|
||||
score_int_to_str(int score_int) {
|
||||
if (score_int == SCORE_FLUX)
|
||||
return "flux";
|
||||
|
||||
if (score_int == SCORE_SCATTER)
|
||||
return "scatter";
|
||||
|
||||
if (score_int == SCORE_TOTAL)
|
||||
return "total";
|
||||
|
||||
if (score_int == SCORE_SCATTER)
|
||||
return "scatter";
|
||||
|
||||
if (score_int == SCORE_NU_SCATTER)
|
||||
return "nu-scatter";
|
||||
|
||||
if (score_int == SCORE_ABSORPTION)
|
||||
return "absorption";
|
||||
|
||||
if (score_int == SCORE_FISSION)
|
||||
return "fission";
|
||||
|
||||
if (score_int == SCORE_NU_FISSION)
|
||||
return "nu-fission";
|
||||
|
||||
if (score_int == SCORE_DECAY_RATE)
|
||||
return "decay-rate";
|
||||
|
||||
if (score_int == SCORE_DELAYED_NU_FISSION)
|
||||
return "delayed-nu-fission";
|
||||
|
||||
if (score_int == SCORE_PROMPT_NU_FISSION)
|
||||
return "prompt-nu-fission";
|
||||
|
||||
if (score_int == SCORE_KAPPA_FISSION)
|
||||
return "kappa-fission";
|
||||
|
||||
if (score_int == SCORE_INVERSE_VELOCITY)
|
||||
return "inverse-velocity";
|
||||
|
||||
if (score_int == SCORE_FISS_Q_PROMPT)
|
||||
return "fission-q-prompt";
|
||||
|
||||
if (score_int == SCORE_FISS_Q_RECOV)
|
||||
return "fission-q-recoverable";
|
||||
|
||||
if (score_int == HEATING)
|
||||
return "heating";
|
||||
|
||||
if (score_int == HEATING_LOCAL)
|
||||
return "heating-local";
|
||||
|
||||
if (score_int == SCORE_CURRENT)
|
||||
return "current";
|
||||
|
||||
if (score_int == SCORE_EVENTS)
|
||||
return "events";
|
||||
|
||||
if (score_int == ELASTIC)
|
||||
return "(n,elastic)";
|
||||
|
||||
if (score_int == N_2N)
|
||||
return "(n,2n)";
|
||||
if (score_int == N_3N)
|
||||
return "(n,3n)";
|
||||
if (score_int == N_4N)
|
||||
return "(n,4n)";
|
||||
if (score_int == N_2ND)
|
||||
return "(n,2nd)";
|
||||
if (score_int == N_2NA)
|
||||
return "(n,na)";
|
||||
if (score_int == N_N3A)
|
||||
return "(n,n3a)";
|
||||
if (score_int == N_2NA)
|
||||
return "(n,2na)";
|
||||
if (score_int == N_3NA)
|
||||
return "(n,3na)";
|
||||
if (score_int == N_NP)
|
||||
return "(n,np)";
|
||||
if (score_int == N_N2A)
|
||||
return "(n,n2a)";
|
||||
if (score_int == N_2N2A)
|
||||
return "(n,2n2a)";
|
||||
if (score_int == N_ND)
|
||||
return "(n,nd)";
|
||||
if (score_int == N_NT)
|
||||
return "(n,nt)";
|
||||
if (score_int == N_N3HE)
|
||||
return "(n,nHe-3)";
|
||||
if (score_int == N_ND2A)
|
||||
return "(n,nd2a)";
|
||||
if (score_int == N_NT2A)
|
||||
return "(n,nt2a)";
|
||||
if (score_int == N_3NF)
|
||||
return "(n,3nf)";
|
||||
if (score_int == N_2NP)
|
||||
return "(n,2np)";
|
||||
if (score_int == N_3NP)
|
||||
return "(n,3np)";
|
||||
if (score_int == N_N2P)
|
||||
return "(n,n2p)";
|
||||
if (score_int == N_NPA)
|
||||
return "(n,npa)";
|
||||
if (score_int == N_N1)
|
||||
return "(n,n1)";
|
||||
if (score_int == N_NC)
|
||||
return "(n,nc)";
|
||||
if (score_int == N_GAMMA)
|
||||
return "(n,gamma)";
|
||||
if (score_int == N_P)
|
||||
return "(n,p)";
|
||||
if (score_int == N_D)
|
||||
return "(n,d)";
|
||||
if (score_int == N_T)
|
||||
return "(n,t)";
|
||||
if (score_int == N_3HE)
|
||||
return "(n,3He)";
|
||||
if (score_int == N_A)
|
||||
return "(n,a)";
|
||||
if (score_int == N_2A)
|
||||
return "(n,2a)";
|
||||
if (score_int == N_3A)
|
||||
return "(n,3a)";
|
||||
if (score_int == N_2P)
|
||||
return "(n,2p)";
|
||||
if (score_int == N_PA)
|
||||
return "(n,pa)";
|
||||
if (score_int == N_T2A)
|
||||
return "(n,t2a)";
|
||||
if (score_int == N_D2A)
|
||||
return "(n,d2a)";
|
||||
if (score_int == N_PD)
|
||||
return "(n,pd)";
|
||||
if (score_int == N_PT)
|
||||
return "(n,pt)";
|
||||
if (score_int == N_DA)
|
||||
return "(n,da)";
|
||||
if (score_int == N_XP)
|
||||
return "H1-production";
|
||||
if (score_int == N_XD)
|
||||
return "H2-production";
|
||||
if (score_int == N_XT)
|
||||
return "H3-production";
|
||||
if (score_int == N_X3HE)
|
||||
return "He3-production";
|
||||
if (score_int == N_XA)
|
||||
return "He4-production";
|
||||
if (score_int == DAMAGE_ENERGY)
|
||||
return "damage-energy";
|
||||
|
||||
// Assume the given int is a reaction MT number. Make sure it's a natural
|
||||
// number then return.
|
||||
std::string score_as_str = std::to_string(score_int);
|
||||
int MT;
|
||||
try {
|
||||
MT = std::stoi(score_as_str);
|
||||
} catch (const std::invalid_argument& ex) {
|
||||
throw std::invalid_argument("Invalid tally score \"" + score_as_str + "\"");
|
||||
}
|
||||
if (MT < 1)
|
||||
throw std::invalid_argument("Invalid tally score \"" + score_as_str + "\"");
|
||||
return "MT" + score_as_str;
|
||||
}
|
||||
|
||||
int
|
||||
score_str_to_int(std::string score_str)
|
||||
{
|
||||
|
|
@ -1001,10 +836,9 @@ void Tally::accumulate()
|
|||
std::string
|
||||
Tally::score_name(int score_idx) const {
|
||||
if (score_idx < 0 || score_idx >= scores_.size()) {
|
||||
warning("Index in scores array is out of bounds.");
|
||||
return "";
|
||||
fatal_error("Index in scores array is out of bounds.");
|
||||
}
|
||||
return score_int_to_str(scores_[score_idx]);
|
||||
return reaction_name(scores_[score_idx]);
|
||||
}
|
||||
|
||||
//==============================================================================
|
||||
|
|
|
|||
|
|
@ -17,22 +17,23 @@ TETS_PER_VOXEL = 12
|
|||
|
||||
|
||||
class UnstructuredMeshTest(PyAPITestHarness):
|
||||
def __init__(self, statepoint_name, **kwargs):
|
||||
|
||||
def __init__(self,
|
||||
statepoint_name,
|
||||
estimator='collision',
|
||||
external_geom=False,
|
||||
holes=None):
|
||||
|
||||
super().__init__(statepoint_name)
|
||||
|
||||
# defaults
|
||||
self.estimator = "collision" # tally estimator type
|
||||
self.external_geom = False # geometry size matches mesh
|
||||
self.mesh_has_holes = False # holes in the mesh
|
||||
self.holes = ()
|
||||
self.mesh_filename = "test_mesh_tets.h5m" # mesh file to use
|
||||
|
||||
# set parameters for the test
|
||||
self.estimator = kwargs.get('estimator', self.estimator)
|
||||
self.external_geom = kwargs.get('external_geom', self.external_geom)
|
||||
self.holes = kwargs.get('holes', self.holes)
|
||||
self.estimator = estimator # tally estimator type
|
||||
self.external_geom = external_geom # geometry size matches mesh
|
||||
self.holes = holes # holes in the test mesh
|
||||
if self.holes:
|
||||
self.mesh_filename = "test_mesh_tets_w_holes.h5m"
|
||||
else:
|
||||
self.mesh_filename = "test_mesh_tets.h5m" # mesh file to use
|
||||
print(self.estimator, self.external_geom, self.holes, self.mesh_filename)
|
||||
|
||||
def _build_inputs(self):
|
||||
### Materials ###
|
||||
|
|
@ -140,8 +141,7 @@ class UnstructuredMeshTest(PyAPITestHarness):
|
|||
|
||||
regular_mesh_filter = openmc.MeshFilter(mesh=regular_mesh)
|
||||
|
||||
uscd_mesh = openmc.UnstructuredMesh()
|
||||
uscd_mesh.filename = self.mesh_filename
|
||||
uscd_mesh = openmc.UnstructuredMesh(self.mesh_filename)
|
||||
uscd_mesh.mesh_lib = 'moab'
|
||||
uscd_filter = openmc.MeshFilter(mesh=uscd_mesh)
|
||||
|
||||
|
|
@ -194,7 +194,6 @@ class UnstructuredMeshTest(PyAPITestHarness):
|
|||
|
||||
def _compare_results(self):
|
||||
with openmc.StatePoint(self._sp_name) as sp:
|
||||
|
||||
# loop over the tallies and get data
|
||||
for tally in sp.tallies.values():
|
||||
# find the regular and unstructured meshes
|
||||
|
|
@ -204,8 +203,8 @@ class UnstructuredMeshTest(PyAPITestHarness):
|
|||
if isinstance(flt.mesh, openmc.RegularMesh):
|
||||
reg_mesh_data, reg_mesh_std_dev = self.get_mesh_tally_data(tally)
|
||||
if self.holes:
|
||||
reg_mesh_data = np.delete(reg_mesh_data, holes)
|
||||
reg_mesh_std_dev = np.delete(reg_mesh_std_dev, holes)
|
||||
reg_mesh_data = np.delete(reg_mesh_data, self.holes)
|
||||
reg_mesh_std_dev = np.delete(reg_mesh_std_dev, self.holes)
|
||||
else:
|
||||
unstructured_data, unstructured_std_dev = self.get_mesh_tally_data(tally, True)
|
||||
|
||||
|
|
@ -218,7 +217,6 @@ class UnstructuredMeshTest(PyAPITestHarness):
|
|||
decimals)
|
||||
except AssertionError as ae:
|
||||
print(ae)
|
||||
print()
|
||||
break
|
||||
# increment decimals
|
||||
decimals += 1
|
||||
|
|
@ -251,7 +249,7 @@ param_values = (['collision', 'tracklength'], # estimators
|
|||
[True, False], # geometry outside of the mesh
|
||||
[(333, 90, 77), (,)]) # location of holes in the mesh
|
||||
test_cases = []
|
||||
for estimator, holes, ext_geom in product(*param_values):
|
||||
for estimator, ext_geom, holes in product(*param_values):
|
||||
test_cases.append({'estimator' : estimator,
|
||||
'external_geom' : ext_geom,
|
||||
'holes' : holes})
|
||||
|
|
@ -259,5 +257,5 @@ for estimator, holes, ext_geom in product(*param_values):
|
|||
|
||||
@pytest.mark.parametrize("opts", test_cases)
|
||||
def test_unstructured_mesh(opts):
|
||||
harness = UnstructuredMeshTest('statepoint.10.h5', kwargs=opts)
|
||||
harness = UnstructuredMeshTest('statepoint.10.h5', **opts)
|
||||
harness.main()
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue