Neml2 czm compare clean

framework/include/neml2/userobjects/MOOSEMaterialPropertyToNEML2.h

@@ -27,10 +27,21 @@ class MOOSEMaterialPropertyToNEML2 : public MOOSEToNEML2Batched<T>
 
 #ifdef NEML2_ENABLED
 protected:
-  const MooseArray<T> & elemMOOSEData() const override { return _mat_prop.get(); }
+  const MooseArray<T> & elemMOOSEData() const override { return _volume_mat_prop.get(); }
+  const MooseArray<T> & elemSideMOOSEData() const override { return _face_mat_prop.get(); }
+  const MooseArray<T> & elemNeighborSideMOOSEData() const override
+  {
+    return _neighbor_face_mat_prop.get();
+  }
 
   /// MOOSE material property to read data from
-  const MaterialProperty<T> & _mat_prop;
+  const MaterialProperty<T> & _volume_mat_prop;
+
+  /// MOOSE material property to read data from (face)
+  const MaterialProperty<T> & _face_mat_prop;
+
+  /// MOOSE material property to read data from (neighbor face)
+  const MaterialProperty<T> & _neighbor_face_mat_prop;
 #endif
 };
 

framework/include/neml2/userobjects/MOOSEToNEML2Batched.h

@@ -10,7 +10,10 @@
 #pragma once
 
 #include "MOOSEToNEML2.h"
-#include "ElementUserObject.h"
+#include "DomainUserObject.h"
+
+#include <set>
+#include <tuple>
 
 /**
  * @brief Generic gatherer for collecting "batched" MOOSE data for NEML2
@@ -19,29 +22,37 @@
  * MooseArray<T>.
  *
  * It is not so generic in the sense that the collected data is always a std::vector of
- * MooseArray<T>, where the vector size is generally the number of elements this ElementUserObject
+ * MooseArray<T>, where the vector size is generally the number of elements this DomainUserObject
  * operates on, and the MooseArray<T> size is generally the number of quadrature points in each
  * element.
  *
  * @tparam T Type of the underlying MOOSE data, e.g., Real, SymmetricRankTwoTensor, etc.
  */
 template <typename T>
-class MOOSEToNEML2Batched : public MOOSEToNEML2, public ElementUserObject
+class MOOSEToNEML2Batched : public MOOSEToNEML2, public DomainUserObject
 {
 public:
   static InputParameters validParams();
 
   MOOSEToNEML2Batched(const InputParameters & params);
 
+  void finalize() override {}
+
 #ifndef NEML2_ENABLED
   void initialize() override {}
-  void execute() override {}
-  void finalize() override {}
+  void executeOnElement() override {}
+  void executeOnBoundary() override {}
+  void executeOnInternalSide() override {}
+  void executeOnExternalSide(const Elem * elem, unsigned int side) override {}
+  void executeOnInterface() override {}
   void threadJoin(const UserObject &) override {}
 #else
   void initialize() override;
-  void execute() override;
-  void finalize() override {}
+  void executeOnElement() override;
+  void executeOnBoundary() override;
+  void executeOnInternalSide() override;
+  void executeOnExternalSide(const Elem * elem, unsigned int side) override;
+  void executeOnInterface() override;
   void threadJoin(const UserObject &) override;
 
   neml2::Tensor gatheredData() const override;
@@ -50,11 +61,22 @@ class MOOSEToNEML2Batched : public MOOSEToNEML2, public ElementUserObject
   std::size_t size() const { return _buffer.size(); }
 
 protected:
+  using ElemSide = std::tuple<dof_id_type, unsigned int>;
+
   /// MOOSE data for the current element
   virtual const MooseArray<T> & elemMOOSEData() const = 0;
 
+  /// MOOSE data for the current element side
+  virtual const MooseArray<T> & elemSideMOOSEData() const = 0;
+
+  /// MOOSE data for the neighboring element side on internal/interface sides
+  virtual const MooseArray<T> & elemNeighborSideMOOSEData() const = 0;
+
   /// Intermediate data buffer, filled during the element loop
   std::vector<T> _buffer;
+
+  /// Element-side keys already gathered in this iteration
+  std::set<ElemSide> _visited_elem_sides;
 #endif
 };
 
@@ -63,7 +85,7 @@ InputParameters
 MOOSEToNEML2Batched<T>::validParams()
 {
   auto params = MOOSEToNEML2::validParams();
-  params += ElementUserObject::validParams();
+  params += DomainUserObject::validParams();
 
   // Since we use the NEML2 model to evaluate the residual AND the Jacobian at the same time, we
   // want to execute this user object only at execute_on = LINEAR (i.e. during residual evaluation).
@@ -77,7 +99,7 @@ MOOSEToNEML2Batched<T>::validParams()
 
 template <typename T>
 MOOSEToNEML2Batched<T>::MOOSEToNEML2Batched(const InputParameters & params)
-  : MOOSEToNEML2(params), ElementUserObject(params)
+  : MOOSEToNEML2(params), DomainUserObject(params)
 {
 }
 
@@ -87,24 +109,95 @@ void
 MOOSEToNEML2Batched<T>::initialize()
 {
   _buffer.clear();
+  _visited_elem_sides.clear();
 }
 
 template <typename T>
 void
-MOOSEToNEML2Batched<T>::execute()
+MOOSEToNEML2Batched<T>::executeOnElement()
 {
   const auto & elem_data = this->elemMOOSEData();
   for (auto i : index_range(elem_data))
     _buffer.push_back(elem_data[i]);
 }
 
+template <typename T>
+void
+MOOSEToNEML2Batched<T>::executeOnBoundary()
+{
+  // Keep boundary gathering consistent with the side-index generator:
+  // sides that have a neighbor are gathered in internal/interface callbacks.
+  if (_neighbor_elem)
+    return;
+
+  const auto elem_side = ElemSide(_current_elem->id(), _current_side);
+  if (_visited_elem_sides.insert(elem_side).second)
+  {
+    const auto & elem_data = this->elemSideMOOSEData();
+    for (auto i : index_range(elem_data))
+      _buffer.push_back(elem_data[i]);
+  }
+}
+
+template <typename T>
+void
+MOOSEToNEML2Batched<T>::executeOnInternalSide()
+{
+  const auto elem_side = ElemSide(_current_elem->id(), _current_side);
+  if (_visited_elem_sides.insert(elem_side).second)
+  {
+    const auto & elem_data = this->elemSideMOOSEData();
+    for (auto i : index_range(elem_data))
+      _buffer.push_back(elem_data[i]);
+  }
+
+  if (_neighbor_elem)
+  {
+    const auto neighbor_side = _neighbor_elem->which_neighbor_am_i(_current_elem);
+    const auto neighbor_elem_side = ElemSide(_neighbor_elem->id(), neighbor_side);
+    if (_visited_elem_sides.insert(neighbor_elem_side).second)
+    {
+      const auto & neighbor_elem_data = this->elemNeighborSideMOOSEData();
+      for (auto i : index_range(neighbor_elem_data))
+        _buffer.push_back(neighbor_elem_data[i]);
+    }
+  }
+}
+
+template <typename T>
+void
+MOOSEToNEML2Batched<T>::executeOnExternalSide(const Elem * /*elem*/, unsigned int /*side*/)
+{
+  // DomainUserObject does not reinit side/qp data for external-side callbacks.
+}
+
+template <typename T>
+void
+MOOSEToNEML2Batched<T>::executeOnInterface()
+{
+  // We do not want to double count the internal side data in the loop, so we only gather data from
+  // the "elem" side of the interface, and skip the "neighbor elem" side
+  const auto elem_side = ElemSide(_current_elem->id(), _current_side);
+  if (_visited_elem_sides.insert(elem_side).second)
+  {
+    const auto & elem_data = this->elemSideMOOSEData();
+    for (auto i : index_range(elem_data))
+      _buffer.push_back(elem_data[i]);
+  }
+
+  // const auto & neighbor_elem_data = this->elemNeighborSideMOOSEData();
+  // for (auto i : index_range(neighbor_elem_data))
+  //   _buffer.push_back(neighbor_elem_data[i]);
+}
+
 template <typename T>
 void
 MOOSEToNEML2Batched<T>::threadJoin(const UserObject & uo)
 {
   // append vectors
   const auto & m2n = static_cast<const MOOSEToNEML2Batched<T> &>(uo);
   _buffer.insert(_buffer.end(), m2n._buffer.begin(), m2n._buffer.end());
+  _visited_elem_sides.insert(m2n._visited_elem_sides.begin(), m2n._visited_elem_sides.end());
 }
 
 template <typename T>

framework/include/neml2/userobjects/MOOSEVariableToNEML2.h

@@ -25,9 +25,17 @@ class MOOSEVariableToNEML2Templ : public MOOSEToNEML2Batched<Real>
 #ifdef NEML2_ENABLED
 protected:
   const MooseArray<Real> & elemMOOSEData() const override { return _moose_variable; }
+  const MooseArray<Real> & elemSideMOOSEData() const override { return _moose_variable; }
+  const MooseArray<Real> & elemNeighborSideMOOSEData() const override
+  {
+    return _moose_variable_neighbor;
+  }
 
   /// Coupled MOOSE variable to read data from
   const VariableValue & _moose_variable;
+
+  /// Coupled neighboring MOOSE variable to read data from
+  const VariableValue & _moose_variable_neighbor;
 #endif
 };
 

framework/include/neml2/userobjects/NEML2BatchIndexGenerator.h

@@ -9,34 +9,46 @@
 
 #pragma once
 
-#include "ElementUserObject.h"
+#include "DomainUserObject.h"
 
 #include <map>
 
 /**
  * NEML2BatchIndexGenerator iterates over the mesh and generates a map from element ID to batch
  * index which is used by NEML2ModelExecutor for transfer data between MOOSE and NEML2.
  */
-class NEML2BatchIndexGenerator : public ElementUserObject
+class NEML2BatchIndexGenerator : public DomainUserObject
 {
 public:
   static InputParameters validParams();
 
   NEML2BatchIndexGenerator(const InputParameters & params);
 
   void initialize() override;
-  void execute() override;
+  void executeOnElement() override;
+  void executeOnBoundary() override;
+  void executeOnInternalSide() override;
+  void executeOnExternalSide(const Elem * elem, unsigned int side) override;
+  void executeOnInterface() override;
   void threadJoin(const UserObject &) override;
   void finalize() override;
 
   void meshChanged() override;
 
+  using ElemSide = std::tuple<dof_id_type, unsigned int>;
+
   /// Get the current batch index (in almost all cases this is the total batch size)
   std::size_t getBatchIndex() const { return _batch_index; }
 
   /// Get the batch index for the given element ID
   std::size_t getBatchIndex(dof_id_type elem_id) const;
 
+  /// Get the batch index for the given element side
+  std::size_t getSideBatchIndex(const ElemSide & elem_side) const;
+
+  /// Check if a batch index exists for the given element side
+  bool isSideBatchIndexExist(const ElemSide & elem_side) const;
+
   /// Whether the batch is empty
   bool isEmpty() const { return _batch_index == 0; }
 
@@ -50,6 +62,12 @@ class NEML2BatchIndexGenerator : public ElementUserObject
   /// Map from element IDs to batch indices
   std::map<dof_id_type, std::size_t> _elem_to_batch_index;
 
+  /// Map from element sides to batch indices
+  std::map<ElemSide, std::size_t> _elemside_to_batch_index;
+
   /// cache the index for the current element
   mutable std::pair<dof_id_type, std::size_t> _elem_to_batch_index_cache;
+
+  /// cache the index for the current element side
+  mutable std::pair<ElemSide, std::size_t> _elemside_to_batch_index_cache;
 };

framework/include/neml2/userobjects/NEML2ModelExecutor.h

@@ -45,6 +45,12 @@ class NEML2ModelExecutor : public NEML2ModelInterface<GeneralUserObject>
   /// Get the batch index for the given element ID
   std::size_t getBatchIndex(dof_id_type elem_id) const;
 
+  /// Get the batch index for the given element side
+  std::size_t getSideBatchIndex(const NEML2BatchIndexGenerator::ElemSide & elem_side) const;
+
+  /// Check if a batch index exists for the given element side
+  bool isSideBatchIndexExist(const NEML2BatchIndexGenerator::ElemSide & elem_side) const;
+
   /// Get a reference(!) to the requested output view
   const neml2::Tensor & getOutput(const neml2::VariableName & output_name) const;
 
@@ -59,6 +65,9 @@ class NEML2ModelExecutor : public NEML2ModelInterface<GeneralUserObject>
   /// check if the output is fully computed and ready to be fetched
   bool outputReady() const { return _output_ready; }
 
+  /// Whether to only execute the NEML2 model at volume Gauss points
+  bool volumeGPOnly() const { return _volume_gp_only; }
+
 protected:
   /// Register a NEML2 input variable gathered by a gatherer
   virtual void addGatheredVariable(const UserObjectName &, const neml2::VariableName &);
@@ -150,5 +159,7 @@ class NEML2ModelExecutor : public NEML2ModelInterface<GeneralUserObject>
   bool _error;
   /// Error message
   std::string _error_message;
+
+  bool _volume_gp_only;
 #endif
 };

framework/src/neml2/actions/NEML2ActionCommon.C

@@ -19,6 +19,8 @@ NEML2ActionCommon::commonParams()
 {
   auto params = NEML2ModelInterface<Action>::validParams();
   params += NEML2ModelExecutor::actionParams();
+  params.addParam<bool>(
+      "volume_gp_only", false, "Whether to only execute the NEML2 model at volume Gauss points.");
 
   MultiMooseEnum moose_types("MATERIAL VARIABLE POSTPROCESSOR");
 

framework/src/neml2/materials/NEML2ToMOOSEMaterialProperty.C

@@ -94,7 +94,10 @@ NEML2ToMOOSEMaterialProperty<T>::computeProperties()
     return;
 
   // look up start index for current element
-  const auto i = _execute_neml2_model.getBatchIndex(_current_elem->id());
+  const auto i = (_bnd && !_execute_neml2_model.volumeGPOnly())
+                     ? _execute_neml2_model.getSideBatchIndex(
+                           NEML2BatchIndexGenerator::ElemSide(_current_elem->id(), _current_side))
+                     : _execute_neml2_model.getBatchIndex(_current_elem->id());
   for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
     NEML2Utils::copyTensorToMOOSEData(_value.batch_index({neml2::Size(i + _qp)}), _prop[_qp]);
 }

framework/src/neml2/userobjects/MOOSEMaterialPropertyToNEML2.C

@@ -36,7 +36,20 @@ MOOSEMaterialPropertyToNEML2<T, state>::MOOSEMaterialPropertyToNEML2(const Input
   : MOOSEToNEML2Batched<T>(params)
 #ifdef NEML2_ENABLED
     ,
-    _mat_prop(this->template getGenericMaterialProperty<T, false>("from_moose", state))
+    _volume_mat_prop(this->template getGenericMaterialProperty<T, false>("from_moose", state)),
+    _face_mat_prop([this]() -> const MaterialProperty<T> & {
+      if constexpr (state == 0)
+        return this->template getFaceMaterialProperty<T>("from_moose");
+      else if constexpr (state == 1)
+        return this->template getFaceMaterialPropertyOld<T>("from_moose");
+      else
+      {
+        static_assert(state == 2, "Unsupported face material property state");
+        return this->template getFaceMaterialPropertyOlder<T>("from_moose");
+      }
+    }()),
+    _neighbor_face_mat_prop(
+        this->template getGenericNeighborMaterialProperty<T, false>("from_moose", state))
 #endif
 {
 }

framework/src/neml2/userobjects/MOOSEVariableToNEML2.C

@@ -28,7 +28,8 @@ MOOSEVariableToNEML2Templ<0>::MOOSEVariableToNEML2Templ(const InputParameters &
   : MOOSEToNEML2Batched(params)
 #ifdef NEML2_ENABLED
     ,
-    _moose_variable(coupledValue("from_moose"))
+    _moose_variable(coupledValue("from_moose")),
+    _moose_variable_neighbor(coupledNeighborValue("from_moose"))
 #endif
 {
 }
@@ -38,7 +39,8 @@ MOOSEVariableToNEML2Templ<1>::MOOSEVariableToNEML2Templ(const InputParameters &
   : MOOSEToNEML2Batched(params)
 #ifdef NEML2_ENABLED
     ,
-    _moose_variable(coupledValueOld("from_moose"))
+    _moose_variable(coupledValueOld("from_moose")),
+    _moose_variable_neighbor(coupledNeighborValueOld("from_moose"))
 #endif
 {
 }