Update openmc mesh tallies

include/userobject/FunctionUserObject.h

@@ -2,6 +2,7 @@
 
 // MOOSE includes
 #include "GeneralUserObject.h"
+#include "DisplacedProblem.h"
 
 // libMesh includes
 #include "libmesh/mesh_function.h"
@@ -34,9 +35,15 @@ class FunctionUserObject : public GeneralUserObject
 
 private:
 
+  // Reference to the mesh
+  MooseMesh& mesh();
+
   // Name of variable we want to turn into a function
   std::string _var_name;
 
+  // Tolerance to pass to point locator
+  double tolerance;
+
   // Pointer to the libMesh system containing our variable
   System* sysPtr;
 

openmc/include/executioners/OpenMCExecutioner.h

@@ -156,6 +156,8 @@ class OpenMCExecutioner : public Transient
   void updateMaterials();
   void updateMaterialDensities();
 
+  // Update mesh tallies
+  void updateMeshTallies();
 
   // Set up OpenMC tallies
   void setupTallies(openmc::Filter* filter_ptr);

openmc/src/executioners/OpenMCExecutioner.C

@@ -475,6 +475,8 @@ OpenMCExecutioner::setupTallies(openmc::Filter* filter_ptr)
   }
 
 }
+
+
 void
 OpenMCExecutioner::setupTally(int32_t& tally_id,
                               openmc::Filter* filter_ptr,
@@ -865,6 +867,8 @@ OpenMCExecutioner::resetOpenMC()
 
   updateMaterials();
 
+  updateMeshTallies();
+
   return true;
 }
 
@@ -1073,6 +1077,49 @@ OpenMCExecutioner::updateMaterialDensities()
   }
 }
 
+void
+OpenMCExecutioner::updateMeshTallies()
+{
+  // Retrieve the current mesh id
+  int32_t mesh_id = openmc::model::meshes.back()->id_;
+
+  // Update in place the mesh pointer
+  openmc::model::meshes.back().reset(new openmc::MOABMesh(moab().moabPtr));
+
+  // Set mesh ID to what is was before
+  openmc::model::meshes.back()->id_ = mesh_id;
+
+  // Get a pointer to current mesh filter
+  openmc::Filter* filter_ptr = openmc::model::tally_filters.back().get();
+
+  // Upcast pointer type
+  openmc::MeshFilter* mesh_filter = dynamic_cast<openmc::MeshFilter*>(filter_ptr);
+
+  int nBinsBefore = mesh_filter->n_bins();
+
+  // Update the mesh in the mesh_filter
+  int32_t mesh_idx = openmc::model::meshes.size() -1;
+  mesh_filter->set_mesh(mesh_idx);
+
+  int nBinsAfter = mesh_filter->n_bins();
+
+  // Update strides in tallies if number of bins has changed
+  if(nBinsAfter != nBinsBefore){
+    for(auto tally_pair : tally_ids_to_scores){
+      int32_t tally_id = tally_pair.first;
+
+      // Look up tally index
+      auto tally_it = openmc::model::tally_map.find(tally_id);
+      if(tally_it == openmc::model::tally_map.end()){
+        mooseError("Could not find tally ID in map.");
+      }
+      int32_t tally_index = tally_it->second;
+
+      // Re-calcuate strides
+      openmc::model::tallies.at(tally_index)->set_strides();
+    }
+  }
+}
 
 bool
 OpenMCExecutioner::setupCells()

openmc/src/userobject/MoabUserObject.C

@@ -115,7 +115,7 @@ MoabUserObject::MoabUserObject(const InputParameters & parameters) :
     powMin = int(floor(log10(var_min)));
     powMax = int(ceil(log10(var_max)));
     nPow = std::max(powMax-powMin, 1);
-    if(nPow > nVarBins){
+    if(nPow > nVarBins && logscale){
       mooseError("Please ensure number of powers for variable is less than the number of bins");
     }
     nMinor = nVarBins/nPow;
@@ -419,8 +419,8 @@ MoabUserObject::createElems(std::map<dof_id_type,moab::EntityHandle>& node_id_to
   moab::Range all_elems;
 
   // Iterate over elements in the mesh
-  auto itelem = mesh().elements_begin();
-  auto endelem = mesh().elements_end();
+  auto itelem = mesh().active_elements_begin();
+  auto endelem = mesh().active_elements_end();
   for( ; itelem!=endelem; ++itelem){
 
     // Get a reference to current elem
@@ -722,8 +722,8 @@ MoabUserObject::setSolution(unsigned int iSysNow,  unsigned int iVarNow, std::ve
   bool procHasNonZeroResult=false;
 
   // When we set the solution, we only want to set dofs that belong to this process
-  auto itelem  = mesh().local_elements_begin();
-  auto endelem = mesh().local_elements_end();
+  auto itelem  = mesh().active_local_elements_begin();
+  auto endelem = mesh().active_local_elements_end();
   for( ; itelem!=endelem; ++itelem){
 
     Elem& elem = **itelem;
@@ -1056,7 +1056,8 @@ MoabUserObject::sortElemsByResults()
   for(const auto & elemSet : sortedElems){
     elemCountCheck += elemSet.size();
   }
-  if(elemCountCheck != mesh().n_elem()){
+
+  if(elemCountCheck != mesh().n_active_elem()){
     mooseError("Disparity in number of sorted elements.");
   }
 
@@ -1285,6 +1286,15 @@ MoabUserObject::resetContainers()
   // Update the serial solutions
   for(const auto& sol :  serial_solutions){
     System & sys = 	systems().get_system(sol.first);
+
+    // Check if solution vector size has changed, e.g. due to mesh refinement
+    if(sys.n_dofs() != sol.second->size()){
+      // clear
+      sol.second->init(0,false,SERIAL);
+      // resize
+      sol.second->init(sys.n_dofs(),false,SERIAL);
+    }
+
     sys.solution->localize(*sol.second);
   }
 }

openmc/unit/include/BasicTest.h

@@ -150,9 +150,16 @@ class InputFileTest : public BasicTest {
     ASSERT_TRUE(fileExists(newName));
   }
 
-  void fetchInput(const std::vector<std::string>& inputFiles){
-    for(const auto file : inputFiles){
-      fetchInputFile(file);
+  void fetchInput(const std::vector<std::string>& inputFilesSrc,
+                  const std::vector<std::string>& inputFilesDest){
+
+    ASSERT_EQ(inputFilesSrc.size(),inputFilesDest.size());
+
+    size_t nFiles = inputFilesSrc.size();
+    for(size_t iFile=0; iFile<nFiles; iFile++){
+      std::string fileSrc = inputFilesSrc.at(iFile);
+      std::string fileDest = inputFilesDest.at(iFile);
+      fetchInputFile(fileSrc,fileDest);
     }
   }
 
@@ -194,6 +201,9 @@ class OpenMCRunTest : public InputFileTest {
     openmcInputXMLFiles.push_back("settings.xml");
     openmcInputXMLFiles.push_back("materials.xml");
 
+    openmcInputXMLFilesSrc.push_back("settings.xml");
+    openmcInputXMLFilesSrc.push_back("materials.xml");
+
     openmcOutputFiles.push_back("summary.h5");
     openmcOutputFiles.push_back("tallies.out");
     openmcOutputFiles.push_back("statepoint.2.h5");
@@ -202,7 +212,8 @@ class OpenMCRunTest : public InputFileTest {
 
   virtual void SetUp() override {
     // Get the input xml files from dir input to current working dir
-    fetchInput(openmcInputXMLFiles);
+    fetchInput(openmcInputXMLFilesSrc,
+               openmcInputXMLFiles);
   }
 
   virtual void TearDown() override {
@@ -213,6 +224,7 @@ class OpenMCRunTest : public InputFileTest {
   }
 
   std::vector<std::string> openmcInputXMLFiles;
+  std::vector<std::string> openmcInputXMLFilesSrc;
   std::vector<std::string> openmcOutputFiles;
   std::string dagmcFilename;
 

src/userobject/FunctionUserObject.C

@@ -11,20 +11,24 @@ validParams<FunctionUserObject>()
 
   params.addRequiredParam<std::string>(
       "variable", "Variable name to evaluate the mesh function on.");
+  params.addParam<double>("tolerance", 1e-9,"Set point locator tolerance.");
+
   return params;
 }
 
 FunctionUserObject::FunctionUserObject(const InputParameters & parameters) :
   GeneralUserObject(parameters),
-  _var_name(getParam<std::string>("variable"))
+  _var_name(getParam<std::string>("variable")),
+  tolerance(getParam<double>("tolerance"))
 {}
 
 void
 FunctionUserObject::execute()
 {
 
   // Fetch a pointer to the point locator object
-  point_locator_ptr = _fe_problem.mesh().getPointLocator();
+  point_locator_ptr = mesh().getPointLocator();
+  point_locator_ptr->set_close_to_point_tol(tolerance);
 
   // Fetch a pointer to the system containing named variable
   sysPtr = &_fe_problem.getSystem(_var_name);
@@ -69,3 +73,12 @@ FunctionUserObject::value (const Point &p) const
     throw std::logic_error("Point locator is null in FunctionUserObject.");
   }
 }
+
+MooseMesh&
+FunctionUserObject::mesh()
+{
+  if(_fe_problem.haveDisplaced()){
+    return _fe_problem.getDisplacedProblem()->mesh();
+  }
+  return _fe_problem.mesh();
+}

unit/Makefile

@@ -37,7 +37,7 @@ RDG                       := no
 RICHARDS                  := no
 SOLID_MECHANICS           := no
 STOCHASTIC_TOOLS          := no
-TENSOR_MECHANICS          := no
+TENSOR_MECHANICS          := yes
 XFEM                      := no
 POROUS_FLOW               := no
 LEVEL_SET                 := no
@@ -68,7 +68,7 @@ OPENMC_APP_LIB = -Wl,-rpath,$(OPENMC_APP_DIR)/lib -L$(OPENMC_APP_DIR)/lib -l$(OP
 
 include $(OPENMC_APP_DIR)/config.inc
 
-ADDITIONAL_INCLUDES += $(OPENMC_APP_INC) $(OPENMC_INC) $(MOAB_INC)
+ADDITIONAL_INCLUDES += $(OPENMC_APP_INC) $(OPENMC_INC) $(MOAB_INC) $(HDF5_INC)
 EXTERNAL_FLAGS += $(OPENMC_APP_LIB) $(OPENMC_LIB) $(MOAB_LIB)