Recover boundary elements (idaholab#32198)

Author: dschwen

framework/include/utils/MooseUtils.h

@@ -648,14 +648,10 @@ getUnion(const std::vector<T> & vector1, const std::vector<T> & vector2, std::ve
 {
   std::unordered_set<T> unique_elements;
   unique_elements.reserve(vector1.size() + vector2.size());
-
-  for (const T & entry : vector1)
-    unique_elements.insert(entry);
-  for (const T & entry : vector2)
-    unique_elements.insert(entry);
+  unique_elements.insert(vector1.begin(), vector1.end());
+  unique_elements.insert(vector2.begin(), vector2.end());
 
   // Now populate the common vector with the union
-  common.clear();
   common.assign(unique_elements.begin(), unique_elements.end());
 }
 

modules/solid_mechanics/include/timeintegrators/NEML2CentralDifference.h

@@ -36,7 +36,7 @@ class NEML2CentralDifference : public ExplicitMixedOrder
 
 private:
   /// Empty element vector to help zero out the algebraic range
-  std::vector<const Elem *> _no_elem = {};
+  std::vector<const Elem *> _boundary_elements = {};
 
   /// Empty node vector to help zero out the algebraic range
   std::vector<const Node *> _no_node = {};

modules/solid_mechanics/src/timeintegrators/NEML2CentralDifference.C

@@ -7,6 +7,9 @@
 //* Licensed under LGPL 2.1, please see LICENSE for details
 //* https://www.gnu.org/licenses/lgpl-2.1.html
 
+#include "MooseTypes.h"
+#include "IntegratedBCBase.h"
+
 #ifdef NEML2_ENABLED
 
 // MOOSE includes
@@ -35,6 +38,27 @@ NEML2CentralDifference::NEML2CentralDifference(const InputParameters & parameter
 void
 NEML2CentralDifference::initialSetup()
 {
+  // build boundary element list by iterating over all integrated BCs
+  const auto & ibcs = _nl->getIntegratedBCWarehouse();
+  std::unordered_set<BoundaryID> bnds;
+  for (const auto & ibc : ibcs.getObjects())
+    bnds.insert(ibc->boundaryIDs().begin(), ibc->boundaryIDs().end());
+
+  if (!bnds.empty())
+  {
+    mooseInfo("Dectected BCs on ", bnds.size(), " boundaries.");
+
+    // deduplicate elements that have multiple boundaries
+    std::unordered_set<const Elem *> unique_elems;
+    const auto end = _fe_problem.mesh().bndElemsEnd();
+    for (auto it = _fe_problem.mesh().bndElemsBegin(); it != end; ++it)
+      if (bnds.find((*it)->_bnd_id) != bnds.end())
+        unique_elems.insert((*it)->_elem);
+
+    _boundary_elems.assign(unique_elems.begin(), unique_elems.end());
+    mooseInfo("Adding ", _boundary_elems.size(), " elements to the algebraic range.");
+  }
+
   ExplicitMixedOrder::initialSetup();
   _neml2_assembly = &_fe_problem.getUserObject<NEML2Assembly>("assembly", /*tid=*/0);
   _fe = &_fe_problem.getUserObject<NEML2FEInterpolation>("fe", /*tid=*/0);
@@ -52,8 +76,8 @@ NEML2CentralDifference::evaluateRHSResidual()
 {
   if (_fe->contextUpToDate() && _neml2_assembly->upToDate())
   {
-    libMesh::ConstElemRange null_elem_range(&_no_elem);
-    _fe_problem.setCurrentAlgebraicElementRange(&null_elem_range);
+    libMesh::ConstElemRange boundary_elem_range(&_boundary_elems);
+    _fe_problem.setCurrentAlgebraicElementRange(&boundary_elem_range);
 
     libMesh::ConstNodeRange null_node_range(&_no_node);
     _fe_problem.setCurrentAlgebraicNodeRange(&null_node_range);