Merge pull request #25769 from lindsayad/fe-face-arg-eval-19420

Implement FaceArg evaluation for finite element variables

framework/include/variables/MooseVariableFE.h

@@ -698,16 +698,22 @@ class MooseVariableFE : public MooseVariableField<OutputType>
   using NodeArg = Moose::NodeArg;
   using ElemPointArg = Moose::ElemPointArg;
 
+  /**
+   * A common method that both evaluate(FaceArg) and evaluateDot(FaceArg) can call. A value
+   * evaluation vs dot evaluation is delineated via the passed-in \p cache_data, e.g. if the
+   * passed-in cache data is the sln data member then this will return a value evaluation and if the
+   * cache data is the dot data member then this will return a dot evaluation
+   */
+  ValueType
+  faceEvaluate(const FaceArg &, const StateArg &, const std::vector<ValueType> & cache_data) const;
+
   ValueType evaluate(const ElemQpArg & elem_qp, const StateArg & state) const override final;
   ValueType evaluate(const ElemSideQpArg & elem_side_qp,
                      const StateArg & state) const override final;
   ValueType evaluate(const ElemArg &, const StateArg &) const override final;
   ValueType evaluate(const ElemPointArg &, const StateArg &) const override final;
   ValueType evaluate(const NodeArg & node_arg, const StateArg & state) const override final;
-  ValueType evaluate(const FaceArg &, const StateArg &) const override final
-  {
-    mooseError("Face info functor overload not yet implemented for finite element variables");
-  }
+  ValueType evaluate(const FaceArg &, const StateArg &) const override final;
 
   GradientType evaluateGradient(const ElemQpArg & elem_qp, const StateArg & state) const override;
   GradientType evaluateGradient(const ElemSideQpArg & elem_side_qp,

framework/src/variables/MooseVariableFE.C

@@ -1031,6 +1031,52 @@ MooseVariableFE<OutputType>::evaluate(const ElemArg & elem_arg, const StateArg &
   return _current_elem_qp_functor_sln[0];
 }
 
+template <typename OutputType>
+typename MooseVariableFE<OutputType>::ValueType
+MooseVariableFE<OutputType>::faceEvaluate(const FaceArg & face_arg,
+                                          const StateArg & state,
+                                          const std::vector<ValueType> & cache_data) const
+{
+  const QMonomial qrule(face_arg.fi->elem().dim() - 1, CONSTANT);
+  auto side_evaluate =
+      [this, &qrule, &state, &cache_data](const Elem * const elem, const unsigned int side)
+  {
+    // We can use whatever we want for the point argument since it won't be used
+    const ElemSideQpArg elem_side_qp_arg{elem, side, /*qp=*/0, &qrule, Point(0, 0, 0)};
+    evaluateOnElementSide(elem_side_qp_arg, state, /*cache_eligible=*/false);
+    return cache_data[0];
+  };
+
+  const auto continuity = this->getContinuity();
+  const bool on_elem = !face_arg.face_side || (face_arg.face_side == face_arg.fi->elemPtr());
+  const bool on_neighbor =
+      !face_arg.face_side || (face_arg.face_side == face_arg.fi->neighborPtr());
+  if (on_neighbor)
+    mooseAssert(
+        face_arg.fi->neighborPtr(),
+        "If we are signaling we should evaluate on the neighbor, we better have a neighbor");
+
+  // Only do multiple evaluations if we are not continuous and we are on an internal face
+  if ((continuity != C_ZERO && continuity != C_ONE) && on_elem && on_neighbor)
+    return (side_evaluate(face_arg.fi->elemPtr(), face_arg.fi->elemSideID()) +
+            side_evaluate(face_arg.fi->neighborPtr(), face_arg.fi->neighborSideID())) /
+           2;
+  else if (on_elem)
+    return side_evaluate(face_arg.fi->elemPtr(), face_arg.fi->elemSideID());
+  else if (on_neighbor)
+    return side_evaluate(face_arg.fi->neighborPtr(), face_arg.fi->neighborSideID());
+  else
+    mooseError(
+        "Attempted to evaluate a moose finite element variable on a face where it is not defined");
+}
+
+template <typename OutputType>
+typename MooseVariableFE<OutputType>::ValueType
+MooseVariableFE<OutputType>::evaluate(const FaceArg & face_arg, const StateArg & state) const
+{
+  return faceEvaluate(face_arg, state, _current_elem_side_qp_functor_sln);
+}
+
 template <typename OutputType>
 typename MooseVariableFE<OutputType>::ValueType
 MooseVariableFE<OutputType>::evaluate(const ElemPointArg & elem_point_arg,
@@ -1214,12 +1260,7 @@ MooseVariableFE<OutputType>::evaluateDot(const FaceArg & face_arg, const StateAr
   mooseAssert(_time_integrator && _time_integrator->dt(),
               "A time derivative is being requested but we do not have a time integrator so we'll "
               "have no idea how to compute it");
-  const QMonomial qrule(face_arg.fi->elem().dim() - 1, CONSTANT);
-  // We can use whatever we want for the point argument since it won't be used
-  const ElemSideQpArg elem_side_qp_arg{
-      face_arg.fi->elemPtr(), face_arg.fi->elemSideID(), /*qp=*/0, &qrule, Point(0, 0, 0)};
-  evaluateOnElementSide(elem_side_qp_arg, state, /*cache_eligible=*/false);
-  return _current_elem_side_qp_functor_dot[0];
+  return faceEvaluate(face_arg, state, _current_elem_side_qp_functor_dot);
 }
 
 template <>

test/tests/functors/fe-var-for-fv-neumann/test.i

@@ -0,0 +1,68 @@
+[Mesh]
+  [gen]
+    type = GeneratedMeshGenerator
+    dim = 1
+    nx = 20
+  []
+[]
+
+[Variables]
+  [fe][]
+  [fv]
+    type = MooseVariableFVReal
+  []
+[]
+
+[Kernels]
+  [diff]
+    type = Diffusion
+    variable = fe
+  []
+[]
+
+[FVKernels]
+  [diff]
+    type = FVDiffusion
+    variable = fv
+    coeff = 1
+  []
+[]
+
+[BCs]
+  [left]
+    type = DirichletBC
+    variable = fe
+    value = 0
+    boundary = left
+   []
+   [right]
+     type = DirichletBC
+     variable = fe
+     value = 1
+     boundary = right
+   []
+ []
+
+[FVBCs]
+  [left]
+    type = FVDirichletBC
+    variable = fv
+    value = 0
+    boundary = left
+  []
+  [right]
+    type = FVFunctorNeumannBC
+    variable = fv
+    functor = fe
+    boundary = right
+  []
+[]
+
+[Executioner]
+  type = Steady
+  solve_type = NEWTON
+[]
+
+[Outputs]
+  exodus = true
+[]

test/tests/functors/fe-var-for-fv-neumann/tests

@@ -0,0 +1,18 @@
+[Tests]
+  design = 'Functors/index.md'
+  issues = '#19420'
+  [exo]
+    type = Exodiff
+    input = test.i
+    exodiff = test_out.e
+    requirement = 'The system shall be able to accurately evaluate a finite element variable through the functor system in a finite volume Dirichlet boundary condition.'
+  []
+  [jac]
+    type = PetscJacobianTester
+    ratio_tol = 1e-7
+    difference_tol = 1e-5
+    input = test.i
+    run_sim = True
+    requirement = 'The system shall compute the correct Jacobian when evaluating a finite element variable through the functor system in a finite volume Dirichlet boundary condition.'
+  []
+[]