From 895249b6ca2ed992d979c472c824f42b68f9bd22 Mon Sep 17 00:00:00 2001 From: Patrick Shriwise Date: Wed, 25 Mar 2020 03:21:39 -0500 Subject: [PATCH] Addressing changes from PR review. --- docs/source/usersguide/beginners.rst | 4 +- .../jupyter/unstructured-mesh-part-i.ipynb | 32 +--- .../jupyter/unstructured-mesh-part-ii.ipynb | 2 +- include/openmc/mesh.h | 17 +- openmc/mesh.py | 35 +++- openmc/statepoint.py | 2 +- src/mesh.cpp | 57 ++---- src/state_point.cpp | 99 +++++----- src/tallies/tally.cpp | 170 +----------------- .../unstructured_mesh/test.py | 38 ++-- 10 files changed, 129 insertions(+), 327 deletions(-) diff --git a/docs/source/usersguide/beginners.rst b/docs/source/usersguide/beginners.rst index a87078cb7..0fbb4bd33 100644 --- a/docs/source/usersguide/beginners.rst +++ b/docs/source/usersguide/beginners.rst @@ -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. diff --git a/examples/jupyter/unstructured-mesh-part-i.ipynb b/examples/jupyter/unstructured-mesh-part-i.ipynb index e28bded76..dab06db99 100644 --- a/examples/jupyter/unstructured-mesh-part-i.ipynb +++ b/examples/jupyter/unstructured-mesh-part-i.ipynb @@ -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!" ] }, { diff --git a/examples/jupyter/unstructured-mesh-part-ii.ipynb b/examples/jupyter/unstructured-mesh-part-ii.ipynb index 60fe335df..8373656db 100644 --- a/examples/jupyter/unstructured-mesh-part-ii.ipynb +++ b/examples/jupyter/unstructured-mesh-part-ii.ipynb @@ -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", diff --git a/include/openmc/mesh.h b/include/openmc/mesh.h index ba68086ce..fe5a79583 100644 --- a/include/openmc/mesh.h +++ b/include/openmc/mesh.h @@ -258,18 +258,10 @@ public: UnstructuredMesh(pugi::xml_node); ~UnstructuredMesh() = default; - void bins_crossed(const Particle* p, std::vector& bins, + void bins_crossed(const Particle* p, + std::vector& bins, std::vector& 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_; // mbi_; //!< MOAB instance std::unique_ptr kdtree_; //!< MOAB KDTree instance diff --git a/openmc/mesh.py b/openmc/mesh.py index d11b60dec..4bde4a52a 100644 --- a/openmc/mesh.py +++ b/openmc/mesh.py @@ -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 diff --git a/openmc/statepoint.py b/openmc/statepoint.py index 84512d70e..5802a74ca 100644 --- a/openmc/statepoint.py +++ b/openmc/statepoint.py @@ -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): diff --git a/src/mesh.cpp b/src/mesh.cpp index e8facb421..6be85685e 100644 --- a/src/mesh.cpp +++ b/src/mesh.cpp @@ -3,7 +3,6 @@ #include // for copy, equal, min, min_element #include // for size_t #include // for ceil -#include // for fmt #include // for allocator #include @@ -16,6 +15,7 @@ #include "xtensor/xsort.hpp" #include "xtensor/xtensor.hpp" #include "xtensor/xview.hpp" +#include // 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(); - // 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& 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); } } diff --git a/src/state_point.cpp b/src/state_point.cpp index e20b508b2..a74b1d685 100644 --- a/src/state_point.cpp +++ b/src/state_point.cpp @@ -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(filter.get()); - auto mesh_idx = mesh_filter->mesh(); - auto umesh = dynamic_cast(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(filter.get()); + auto mesh_idx = mesh_filter->mesh(); + auto umesh = dynamic_cast(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 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 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); } } } diff --git a/src/tallies/tally.cpp b/src/tallies/tally.cpp index bbc671d21..6b89136a0 100644 --- a/src/tallies/tally.cpp +++ b/src/tallies/tally.cpp @@ -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]); } //============================================================================== diff --git a/tests/regression_tests/unstructured_mesh/test.py b/tests/regression_tests/unstructured_mesh/test.py index 1b1b98a46..e07a0f1de 100644 --- a/tests/regression_tests/unstructured_mesh/test.py +++ b/tests/regression_tests/unstructured_mesh/test.py @@ -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()