Apply suggestions from code review (idaholab#26053)

Co-authored-by: Alex Lindsay <alexander.lindsay@inl.gov>

address_more_comments

framework/include/actions/ProjectedStatefulMaterialStorageAction.h

@@ -18,9 +18,9 @@
 #include "libmesh/fe_type.h"
 
 // created types
-#include <InterpolatedStatefulMaterial.h>
-#include <ProjectedStatefulMaterialAux.h>
-#include <ProjectedStatefulMaterialNodalPatchRecovery.h>
+#include "InterpolatedStatefulMaterial.h"
+#include "ProjectedStatefulMaterialAux.h"
+#include "ProjectedStatefulMaterialNodalPatchRecovery.h"
 
 /**
  * Set up AuxKernels and AuxVariables for projected material property storage (PoMPS).
@@ -37,7 +37,7 @@ class ProjectedStatefulMaterialStorageAction : public Action
   using Action::addRelationshipManagers;
   virtual void addRelationshipManagers(Moose::RelationshipManagerType input_rm_type) override;
 
-  /// List of supported types
+  /// List of supported raw types (equivalent AD types are also supported)
   typedef Moose::TypeList<Real, RealVectorValue, RankTwoTensor, RankFourTensor> SupportedTypes;
 
   /// get an enum with all supported types
@@ -61,12 +61,17 @@ class ProjectedStatefulMaterialStorageAction : public Action
   template <typename... Ts>
   static MooseEnum getTypeEnum(typename Moose::TypeList<Ts...>);
 
+  /// Property names for which we will do stateful material property projection
   const std::vector<MaterialPropertyName> & _prop_names;
 
+  /// Variable order to project the properties into
   const MooseEnum _order;
+
+  /// FEType for the variables to project the properties into
   FEType _fe_type;
+
+  /// MooseObject name for the underlying variable type
   const std::string _var_type;
-  const std::string _pomps_prefix;
 };
 
 template <typename... Ts>
@@ -81,8 +86,10 @@ void
 ProjectedStatefulMaterialStorageAction::processProperty(const MaterialPropertyName & prop_name)
 {
   // check if a property of type T exists
-  const auto & data = _problem->getMaterialData(Moose::BLOCK_MATERIAL_DATA);
-  if (!data.haveGenericProperty<T, is_ad>(prop_name))
+  const auto & block_data = _problem->getMaterialData(Moose::BLOCK_MATERIAL_DATA);
+  const auto & boundary_data = _problem->getMaterialData(Moose::BOUNDARY_MATERIAL_DATA);
+  if (!block_data.haveGenericProperty<T, is_ad>(prop_name) &&
+      !boundary_data.haveGenericProperty<T, is_ad>(prop_name))
     return;
 
   // number of scalar components

framework/include/auxkernels/NodalPatchRecoveryAuxBase.h

@@ -29,8 +29,9 @@ class NodalPatchRecoveryAuxBase : public AuxKernel
 protected:
   virtual Real computeValue() override;
 
-  /// Override this to get the fitted value form a Nodal Patch Recovery User Object
+  /// Override this to get the fitted value from a Nodal Patch Recovery User Object
   virtual Real nodalPatchRecovery() = 0;
 
+  /// local patch of elements used for recovery
   std::vector<dof_id_type> _elem_ids;
 };

framework/include/auxkernels/ProjectedMaterialPropertyNodalPatchRecoveryAux.h

@@ -25,6 +25,9 @@ class ProjectedMaterialPropertyNodalPatchRecoveryAux : public NodalPatchRecovery
   virtual Real nodalPatchRecovery() override;
 
 private:
+  /// User object holding the data needed for patch recovery
   const ProjectedStatefulMaterialNodalPatchRecoveryBase & _npr;
+
+  /// Property component (index into a serialized representation of the property)
   const unsigned int _component;
 };

framework/include/auxkernels/ProjectedStatefulMaterialAux.h

@@ -10,9 +10,8 @@
 #pragma once
 
 #include "AuxKernel.h"
-#include "IndexableProperty.h"
 
-#include <deque>
+#include <unordered_map>
 
 class MaterialBase;
 
@@ -29,12 +28,16 @@ class ProjectedStatefulMaterialAuxTempl : public AuxKernel
 protected:
   virtual Real computeValue() override;
 
+  /// ID of the subdomain currently being iterated over
   const SubdomainID & _current_subdomain_id;
 
-  std::map<SubdomainID, std::vector<MaterialBase *>> _required_materials;
+  /// A sorted list of all material objects the processed property depends on on each subdomain
+  std::unordered_map<SubdomainID, std::vector<MaterialBase *>> _required_materials;
 
+  /// The property a component of which is being projected
   const GenericMaterialProperty<T, is_ad> & _prop;
 
+  /// Property component (index into a serialized representation of the property)
   const unsigned int _component;
 };
 

framework/include/materials/InterpolatedStatefulMaterial.h

@@ -14,8 +14,8 @@
 #include "RankFourTensorForward.h"
 
 /**
- * Reconstitute a materal property from the old and older states of projected AuxVariables. Use
- * though the ProjectedStatefulMaterialStorageAction.
+ * Reconstitute a materal property from the old and older states of projected AuxVariables. Used
+ * through the ProjectedStatefulMaterialStorageAction.
  */
 template <typename T>
 class InterpolatedStatefulMaterialTempl : public Material
@@ -34,7 +34,7 @@ class InterpolatedStatefulMaterialTempl : public Material
   /// Older projected state
   const std::vector<const VariableValue *> _older_state;
 
-  /// total number of components
+  /// Total number of components (e.g. 1 for Real, LIBMESH_DIM for RealVectorValue, LIBMESH_DIM^2 for RankTwoTensor)
   const std::size_t _size;
 
   /// emitted property name

framework/include/materials/Material.h

@@ -277,9 +277,9 @@ Material::getGenericMaterialPropertyByName(const std::string & prop_name_in,
   if (_use_interpolated_state)
   {
     if (state == 1)
-      return getGenericMaterialPropertyByName<T, is_ad>(prop_name_in + "_interpolated_old", 0);
+      return getGenericMaterialPropertyByName<T, is_ad>(prop_name_in + _interpolated_old, 0);
     if (state == 2)
-      return getGenericMaterialPropertyByName<T, is_ad>(prop_name_in + "_interpolated_older", 0);
+      return getGenericMaterialPropertyByName<T, is_ad>(prop_name_in + _interpolated_older, 0);
   }
 
   MaterialBase::checkExecutionStage();

framework/include/materials/MaterialPropertyInterface.h

@@ -18,7 +18,7 @@
 #include "InputParameters.h"
 #include "SubProblem.h"
 
-#include <deque>
+#include <unordered_map>
 
 // Forward declarations
 class MooseObject;
@@ -265,7 +265,7 @@ class MaterialPropertyInterface
   ///@}
 
   /// get a map of MaterialBase pointers for all material objects that this object depends on for each block
-  std::map<SubdomainID, std::vector<MaterialBase *>>
+  std::unordered_map<SubdomainID, std::vector<MaterialBase *>>
   buildRequiredMaterials(bool allow_stateful = true);
 
   ///@{
@@ -538,6 +538,11 @@ class MaterialPropertyInterface
   /// Use the interpolated state set up through the ProjectedStatefulMaterialStorageAction
   const bool _use_interpolated_state;
 
+  ///@{ name suffixes for interpolated old and older properties
+  static const std::string _interpolated_old;
+  static const std::string _interpolated_older;
+  ///@}
+
 private:
   /**
    * @returns The MaterialDataType given the interface's parameters
@@ -775,10 +780,10 @@ MaterialPropertyInterface::getGenericMaterialPropertyByName(const MaterialProper
   {
     if (state == 1)
       return getGenericMaterialPropertyByName<T, is_ad>(
-          name_in + "_interpolated_old", material_data, 0);
+          name_in + _interpolated_old, material_data, 0);
     if (state == 2)
       return getGenericMaterialPropertyByName<T, is_ad>(
-          name_in + "_interpolated_older", material_data, 0);
+          name_in + _interpolated_older, material_data, 0);
   }
 
   const auto name = _get_suffix.empty()

framework/include/userobjects/ProjectedStatefulMaterialNodalPatchRecovery.h

@@ -13,6 +13,8 @@
 #include "ElementUserObject.h"
 #include "SerialAccess.h"
 
+#include <unordered_map>
+
 class ProjectedStatefulMaterialNodalPatchRecoveryBase : public ElementUserObject
 {
 public:
@@ -38,7 +40,7 @@ class ProjectedStatefulMaterialNodalPatchRecoveryTempl
   ProjectedStatefulMaterialNodalPatchRecoveryTempl(const InputParameters & parameters);
 
   /**
-   * Solve the least-squares problem. Use the fitted coefficients to calculate  the value at the
+   * Solve the least-squares problem. Use the fitted coefficients to calculate the value at the
    * requested point.
    *
    * @param p           Point at which to compute the fitted value
@@ -94,13 +96,13 @@ class ProjectedStatefulMaterialNodalPatchRecoveryTempl
   const GenericMaterialProperty<T, is_ad> & _prop;
 
   /// The element-level A matrix and the element-level b vectors for each component
-  std::map<dof_id_type, ElementData> _abs;
+  std::unordered_map<dof_id_type, ElementData> _abs;
 
   /// Current subdomain
   const SubdomainID & _current_subdomain_id;
 
-  /// list of require materials that need to be explicityly initialized at step zero
-  std::map<SubdomainID, std::vector<MaterialBase *>> _required_materials;
+  /// list of required materials that need to be explicitly initialized at step zero
+  std::unordered_map<SubdomainID, std::vector<MaterialBase *>> _required_materials;
 };
 
 typedef ProjectedStatefulMaterialNodalPatchRecoveryTempl<Real, false>

framework/include/utils/SerialAccess.h

@@ -80,28 +80,30 @@ SERIAL_ACCESS_DYNAMIC_SIZE(DenseVector, &obj(0u), obj.size());
  * member or where value_type doesn't have a suitable meaning (ADReal)).
  */
 template <typename T>
-struct SerialAccessVlaueTypeHelper
+struct SerialAccessValueTypeHelper
 {
   typedef typename T::value_type value_type;
 };
 template <>
-struct SerialAccessVlaueTypeHelper<ADReal>
+struct SerialAccessValueTypeHelper<ADReal>
 {
   typedef ADReal value_type;
 };
 template <>
-struct SerialAccessVlaueTypeHelper<Real>
+struct SerialAccessValueTypeHelper<Real>
 {
   typedef Real value_type;
 };
 
 template <typename T>
 class SerialAccessRange
 {
-  typedef typename SerialAccessVlaueTypeHelper<typename std::remove_const<T>::type>::value_type R;
+public:
+  /// Value type of the components of T
+  typedef typename SerialAccessValueTypeHelper<typename std::remove_const<T>::type>::value_type R;
+  /// Value type with the correct constness
   typedef typename std::conditional<std::is_const_v<T>, const R, R>::type V;
 
-public:
   class iterator
   {
   public:

framework/src/actions/ProjectedStatefulMaterialStorageAction.C

@@ -34,7 +34,7 @@ ProjectedStatefulMaterialStorageAction::validParams()
 
   params.addParam<MooseEnum>(
       "family",
-      MooseEnum("LAGRANGE MONOMIAL", "LAGRANGE"),
+      MooseEnum("LAGRANGE MONOMIAL L2_LAGRANGE L2_HIERARCHIC", "LAGRANGE"),
       "Finite element variable family to project the material properties onto");
   params.addParam<MooseEnum>(
       "order",
@@ -56,8 +56,7 @@ ProjectedStatefulMaterialStorageAction::ProjectedStatefulMaterialStorageAction(
     _order(params.get<MooseEnum>("order")),
     _fe_type({Utility::string_to_enum<Order>(_order),
               Utility::string_to_enum<FEFamily>(params.get<MooseEnum>("family"))}),
-    _var_type(AddVariableAction::variableType(_fe_type)),
-    _pomps_prefix("_pomps_")
+    _var_type(AddVariableAction::variableType(_fe_type))
 {
 }
 

framework/src/auxkernels/NodalPatchRecoveryAuxBase.C

@@ -32,7 +32,7 @@ Real
 NodalPatchRecoveryAuxBase::computeValue()
 {
   // get node-to-conneted-elem map
-  const std::map<dof_id_type, std::vector<dof_id_type>> & node_to_elem_map = _mesh.nodeToElemMap();
+  const auto & node_to_elem_map = _mesh.nodeToElemMap();
   auto node_to_elem_pair = node_to_elem_map.find(_current_node->id());
   mooseAssert(node_to_elem_pair != node_to_elem_map.end(), "Missing entry in node to elem map");
 
@@ -42,7 +42,7 @@ NodalPatchRecoveryAuxBase::computeValue()
   // consider the case for corner node
   if (_elem_ids.size() == 1)
   {
-    dof_id_type elem_id = _elem_ids[0];
+    const dof_id_type elem_id = _elem_ids[0];
     for (auto & n : _mesh.elemPtr(elem_id)->node_ref_range())
     {
       node_to_elem_pair = node_to_elem_map.find(n.id());

framework/src/auxkernels/ProjectedMaterialPropertyNodalPatchRecoveryAux.C

@@ -33,6 +33,7 @@ ProjectedMaterialPropertyNodalPatchRecoveryAux::ProjectedMaterialPropertyNodalPa
                "Nodal patch recovery auxiliary kernel is not defined in a subset of blocks of the "
                "associated user object. Revise your input file.");
 }
+
 Real
 ProjectedMaterialPropertyNodalPatchRecoveryAux::nodalPatchRecovery()
 {

framework/src/materials/InterpolatedStatefulMaterial.C

@@ -34,8 +34,8 @@ InterpolatedStatefulMaterialTempl<T>::InterpolatedStatefulMaterialTempl(
     _older_state(coupledValuesOlder("old_state")),
     _size(Moose::SerialAccess<T>::size()),
     _prop_name(getParam<MaterialPropertyName>("prop_name")),
-    _prop_old(declareProperty<T>(_prop_name + "_interpolated_old")),
-    _prop_older(declareProperty<T>(_prop_name + "_interpolated_older"))
+    _prop_old(declareProperty<T>(_prop_name + _interpolated_old)),
+    _prop_older(declareProperty<T>(_prop_name + _interpolated_older))
 {
   if (_old_state.size() != _size)
     paramError("old_state", "Wrong number of component AuxVariables passed in.");

framework/src/materials/MaterialPropertyInterface.C

@@ -12,6 +12,9 @@
 #include "MooseApp.h"
 #include "MaterialBase.h"
 
+const std::string MaterialPropertyInterface::_interpolated_old = "_interpolated_old";
+const std::string MaterialPropertyInterface::_interpolated_older = "_interpolated_older";
+
 InputParameters
 MaterialPropertyInterface::validParams()
 {
@@ -198,10 +201,10 @@ MaterialPropertyInterface::getMaterialByName(const std::string & name, bool no_w
   return *discrete;
 }
 
-std::map<SubdomainID, std::vector<MaterialBase *>>
+std::unordered_map<SubdomainID, std::vector<MaterialBase *>>
 MaterialPropertyInterface::buildRequiredMaterials(bool allow_stateful)
 {
-  std::map<SubdomainID, std::vector<MaterialBase *>> required_mats;
+  std::unordered_map<SubdomainID, std::vector<MaterialBase *>> required_mats;
   const auto & mwh = _mi_feproblem.getMaterialWarehouse();
   for (const auto id : _mi_block_ids)
   {

framework/src/userobjects/ProjectedStatefulMaterialNodalPatchRecovery.C

@@ -10,6 +10,7 @@
 #include "ProjectedStatefulMaterialNodalPatchRecovery.h"
 #include "ElementUserObject.h"
 #include "MaterialBase.h"
+#include "MathUtils.h"
 #include "Assembly.h"
 
 // TIMPI includes
@@ -101,18 +102,18 @@ ProjectedStatefulMaterialNodalPatchRecoveryTempl<T, is_ad>::nodalPatchRecovery(
   // Assemble the least squares problem over the patch
   RealEigenMatrix A = RealEigenMatrix::Zero(_q, _q);
   RealEigenVector b = RealEigenVector::Zero(_q);
-  for (auto elem_id : elem_ids)
+  for (const auto elem_id : elem_ids)
   {
     const auto abs = libmesh_map_find(_abs, elem_id);
     A += abs.first;
     b += abs.second[component];
   }
 
   // Solve the least squares fitting
-  RealEigenVector coef = A.completeOrthogonalDecomposition().solve(b);
+  const RealEigenVector coef = A.completeOrthogonalDecomposition().solve(b);
 
   // Compute the fitted nodal value
-  RealEigenVector p = evaluateBasisFunctions(x);
+  const RealEigenVector p = evaluateBasisFunctions(x);
   return p.dot(coef);
 }
 
@@ -123,13 +124,12 @@ ProjectedStatefulMaterialNodalPatchRecoveryTempl<T, is_ad>::evaluateBasisFunctio
 {
   RealEigenVector p(_q);
   Real polynomial;
-  for (unsigned int r = 0; r < _multi_index.size(); r++)
+  for (const auto r : index_range(_multi_index))
   {
     polynomial = 1.0;
     mooseAssert(_multi_index[r].size() == _mesh.dimension(), "Wrong multi-index size.");
-    for (unsigned int c = 0; c < _multi_index[r].size(); c++)
-      for (unsigned int p = 0; p < _multi_index[r][c]; p++)
-        polynomial *= q_point(c);
+    for (const auto c : index_range(_multi_index[r]))
+      polynomial *= MathUtils::pow(q_point(c), _multi_index[r][c]);
     p(r) = polynomial;
   }
   return p;
@@ -163,7 +163,7 @@ ProjectedStatefulMaterialNodalPatchRecoveryTempl<T, is_ad>::execute()
       bs[index++] += MetaPhysicL::raw_value(v) * p;
   }
 
-  dof_id_type elem_id = _current_elem->id();
+  const dof_id_type elem_id = _current_elem->id();
   _abs[elem_id] = {Ae, bs};
 }
 
@@ -187,7 +187,7 @@ ProjectedStatefulMaterialNodalPatchRecoveryTempl<T, is_ad>::finalize()
   // Populate algebraically ghosted elements to query
   std::unordered_map<processor_id_type, std::vector<dof_id_type>> query_ids;
   const ConstElemRange evaluable_elem_range = _fe_problem.getEvaluableElementRange();
-  for (auto elem : evaluable_elem_range)
+  for (const auto * const elem : evaluable_elem_range)
     if (elem->processor_id() != processor_id())
       query_ids[elem->processor_id()].push_back(elem->id());