update strides method

framework/include/problems/DisplacedProblem.h

@@ -374,6 +374,8 @@ class DisplacedProblem : public SubProblem
 
   virtual bool hasNonlocalCoupling() const override;
 
+  void setStrides();
+
 protected:
   FEProblemBase & _mproblem;
   MooseMesh & _mesh;

framework/include/problems/FEProblemBase.h

@@ -2750,6 +2750,8 @@ class FEProblemBase : public SubProblem, public Restartable
    */
   virtual void resetState();
 
+  void setStrides();
+
   bool _error_on_jacobian_nonzero_reallocation;
   bool _ignore_zeros_in_jacobian;
   const bool _force_restart;

framework/include/systems/SystemBase.h

@@ -940,6 +940,11 @@ class SystemBase : public libMesh::ParallelObject, public ConsoleStreamInterface
    */
   virtual void compute(ExecFlagType type) = 0;
 
+  /**
+   * Initialize the system
+   */
+  virtual void setStrides();
+
 protected:
   /**
    * Internal getter for solution owned by libMesh.

framework/include/variables/MooseVariableBase.h

@@ -182,14 +182,14 @@ class MooseVariableBase : public MooseObject,
    */
   bool isLowerD() const { return _is_lower_d; }
 
+  void setStride() const;
+
 protected:
   /**
    * @returns whether we should insert derivatives
    */
   bool doDerivatives() const;
 
-  void setStride() const;
-
   /// System this variable is part of
   SystemBase & _sys;
 

framework/src/problems/DisplacedProblem.C

@@ -1376,3 +1376,10 @@ DisplacedProblem::solverSysNum(const SolverSystemName & sys_name) const
 {
   return _mproblem.solverSysNum(sys_name);
 }
+void
+DisplacedProblem::setStrides()
+{
+  for (auto & sys : _displaced_solver_systems)
+    sys->setStrides();
+  _displaced_aux->setStrides();
+}

framework/src/problems/FEProblemBase.C

@@ -1336,6 +1336,7 @@ FEProblemBase::timestepSetup()
     // u4) Now that all the geometric searches have been done (both undisplaced and displaced),
     //     we're ready to update the sparsity pattern
     es().reinit_systems();
+    setStrides();
   }
 
   if (_line_search)
@@ -4709,6 +4710,8 @@ FEProblemBase::reinitBecauseOfGhostingOrNewGeomObjects(const bool mortar_changed
 
     if (_displaced_mesh)
       _displaced_problem->es().reinit();
+
+    setStrides();
   }
 }
 
@@ -5844,6 +5847,7 @@ FEProblemBase::init()
   if (_displaced_problem)
     _displaced_problem->init();
 
+  setStrides();
   _initialized = true;
 }
 
@@ -7472,6 +7476,7 @@ FEProblemBase::adaptMesh()
   if (mesh_changed)
     es().reinit_systems();
 
+  setStrides();
   return mesh_changed;
 }
 #endif // LIBMESH_ENABLE_AMR
@@ -7560,6 +7565,8 @@ FEProblemBase::meshChangedHelper(bool intermediate_change)
   else
     es().reinit();
 
+  setStrides();
+
   // Updating MooseMesh first breaks other adaptivity code, unless we
   // then *again* update the MooseMesh caches.  E.g. the definition of
   // "active" and "local" may have been *changed* by refinement and
@@ -8776,3 +8783,12 @@ FEProblemBase::setCurrentAlgebraicBndNodeRange(ConstBndNodeRange * range)
 
   _current_algebraic_bnd_node_range = std::make_unique<ConstBndNodeRange>(*range);
 }
+void
+FEProblemBase::setStrides()
+{
+  for (auto & sys : _solver_systems)
+    sys->setStrides();
+  _aux->setStrides();
+  if (_displaced_problem)
+    _displaced_problem->setStrides();
+}

framework/src/systems/SystemBase.C

@@ -1214,6 +1214,17 @@ SystemBase::copyVars(ExodusII_IO & io)
     solution().close();
 }
 
+void
+SystemBase::setStrides()
+{
+  for (MooseIndex(_vars) thread = 0; thread < _vars.size(); ++thread)
+  {
+    const std::vector<MooseVariableFieldBase *> & var_vec = _vars[thread].fieldVariables();
+    for (auto & var : var_vec)
+      var->setStride();
+  }
+}
+
 void
 SystemBase::update(const bool update_libmesh_system)
 {

framework/src/variables/MooseVariableBase.C

@@ -180,40 +180,35 @@ MooseVariableBase::componentDofIndices(const std::vector<dof_id_type> & dof_indi
 
   // Only perform mapping if the component is not 0
   if (component != 0)
-  {
-    // Set once. Need to have redo stride when mesh or dof map changes.
-    setStride();
-
     for (auto & id : new_dof_indices)
       id += component * _stride[id];
-  }
 
   return new_dof_indices;
 }
 
 void
 MooseVariableBase::setStride() const
 {
-  if (_stride.empty())
-  {
-    std::vector<dof_id_type> comp_0_dofs;
-    _dof_map.local_variable_indices(comp_0_dofs, _mesh, _var_num);
-    // Determine the variable number for the specified component
-    auto comp_num = _var_num + 1;
-
-    // Grab all DOFs for the specified component
-    std::vector<dof_id_type> comp_1_dofs;
-    _dof_map.local_variable_indices(comp_1_dofs, _mesh, comp_num);
-
-    // This hangs when a processor doesn't have the variable.
-    _communicator.allgather(comp_0_dofs);
-    _communicator.allgather(comp_1_dofs);
-
-    _stride.resize(*std::max_element(comp_0_dofs.begin(), comp_0_dofs.end()) + 1,
-                   std::numeric_limits<dof_id_type>::max());
-    for (const auto i : index_range(comp_0_dofs))
-      _stride[comp_0_dofs[i]] = comp_1_dofs[i] - comp_0_dofs[i];
-  }
+  if (!_is_array || _count == 1)
+    return;
+
+  std::vector<dof_id_type> comp_0_dofs;
+  _dof_map.local_variable_indices(comp_0_dofs, _mesh, _var_num);
+  // Determine the variable number for the specified component
+  auto comp_num = _var_num + 1;
+
+  // Grab all DOFs for the specified component
+  std::vector<dof_id_type> comp_1_dofs;
+  _dof_map.local_variable_indices(comp_1_dofs, _mesh, comp_num);
+
+  // This hangs when a processor doesn't have the variable.
+  _communicator.allgather(comp_0_dofs);
+  _communicator.allgather(comp_1_dofs);
+
+  _stride.resize(*std::max_element(comp_0_dofs.begin(), comp_0_dofs.end()) + 1,
+                 std::numeric_limits<dof_id_type>::max());
+  for (const auto i : index_range(comp_0_dofs))
+    _stride[comp_0_dofs[i]] = comp_1_dofs[i] - comp_0_dofs[i];
 }
 
 void

framework/src/variables/MooseVariableDataBase.C

@@ -405,22 +405,14 @@ MooseVariableDataBase<RealEigenVector>::insert(NumericVector<Number> & residual)
   if (_has_dof_values)
   {
     auto & dof_values = _vector_tags_dof_u[_solution_tag];
-    if (isNodal())
-    {
-      for (unsigned int i = 0; i < _dof_indices.size(); ++i)
-        for (unsigned int j = 0; j < _count; ++j)
-          residual.set(_dof_indices[i] + j, dof_values[i](j));
-    }
-    else
-    {
-      unsigned int n = 0;
+    std::vector<std::vector<dof_id_type>> all_dofs_indices(_count);
+    all_dofs_indices[0] = _dof_indices;
+    for (unsigned int comp = 1; comp < _count; ++comp)
+      all_dofs_indices[comp] = _var.componentDofIndices(_dof_indices, comp);
+
+    for (unsigned int i = 0; i < _dof_indices.size(); ++i)
       for (unsigned int j = 0; j < _count; ++j)
-      {
-        for (unsigned int i = 0; i < _dof_indices.size(); ++i)
-          residual.set(_dof_indices[i] + n, dof_values[i](j));
-        n += _dof_indices.size();
-      }
-    }
+        residual.set(all_dofs_indices[j][i], dof_values[i](j));
   }
 }
 
@@ -442,20 +434,17 @@ MooseVariableDataBase<RealEigenVector>::add(NumericVector<Number> & residual)
   if (_has_dof_values)
   {
     auto & dof_values = _vector_tags_dof_u[_solution_tag];
-    if (isNodal())
-    {
-      for (unsigned int i = 0; i < _dof_indices.size(); ++i)
-        for (unsigned int j = 0; j < _count; ++j)
-          residual.add(_dof_indices[i] + j, dof_values[i](j));
-    }
-    else
+
+    std::vector<std::vector<dof_id_type>> all_dofs_indices(_count);
+    all_dofs_indices[0] = _dof_indices;
+    for (unsigned int comp = 1; comp < _count; ++comp)
+      all_dofs_indices[comp] = _var.componentDofIndices(_dof_indices, comp);
+
+    for (unsigned int i = 0; i < _dof_indices.size(); ++i)
     {
-      unsigned int n = 0;
       for (unsigned int j = 0; j < _count; ++j)
       {
-        for (unsigned int i = 0; i < _dof_indices.size(); ++i)
-          residual.add(_dof_indices[i] + n, dof_values[i](j));
-        n += _dof_indices.size();
+        residual.add(all_dofs_indices[j][i], dof_values[i](j));
       }
     }
   }

test/tests/variables/side_hierarchic/side_hierarchic_array.i

@@ -0,0 +1,58 @@
+[Problem]
+  solve = false
+[]
+
+[Mesh]
+  type = GeneratedMesh
+  elem_type = QUAD9
+  dim = 2
+  nx = 2
+  ny = 2
+[]
+
+[Variables]
+  [side_var]
+    order = ELEVENTH
+    family = SIDE_HIERARCHIC
+    components = 2
+  []
+[]
+
+[AuxVariables]
+  [aux_side_var]
+    order = FIRST
+    family = SIDE_HIERARCHIC
+  []
+[]
+
+[Functions]
+  [nl_var]
+    type = ParsedFunction
+    expression = 'x+y+1'
+  []
+  [aux_var]
+    type = ParsedFunction
+    expression = 'x-y+10'
+  []
+[]
+
+[ICs]
+  [side_nl]
+    type = ArrayFunctionIC
+    variable = side_var
+    function = 'nl_var nl_var'
+  []
+  [side_aux]
+    type = FunctionIC
+    variable = aux_side_var
+    function = aux_var
+  []
+[]
+
+[Outputs]
+  exodus = true
+[]
+
+[Executioner]
+  type = Steady
+[]

test/tests/variables/side_hierarchic/tests

@@ -8,4 +8,13 @@
     design = 'syntax/Variables/index.md'
     requirement = 'The system shall initialize SIDE_HIERARCHIC variables as instructed'
   [../]
+  [side_hierarchic_array]
+      type = Exodiff
+      input = 'side_hierarchic_array.i'
+      exodiff = 'side_hierarchic_array_out.e'
+
+      issues = '#27639'
+      design = 'syntax/Variables/index.md'
+      requirement = 'The system shall initialize SIDE_HIERARCHIC array variables as instructed'
+  []
 []