[PL] support FunctionParameter in Neumann bc

ProcessLib/BoundaryCondition/GenericNaturalBoundaryConditionLocalAssembler.h

@@ -9,6 +9,7 @@
 
 #pragma once
 
+#include "MeshLib/Elements/Element.h"
 #include "NumLib/DOF/LocalToGlobalIndexMap.h"
 #include "NumLib/Fem/ShapeMatrixPolicy.h"
 #include "ProcessLib/Utils/InitShapeMatrices.h"
@@ -44,7 +45,8 @@ public:
         : _integration_method(integration_order),
           _shape_matrices(initShapeMatrices<ShapeFunction, ShapeMatricesType,
                                             IntegrationMethod, GlobalDim>(
-              e, is_axially_symmetric, _integration_method))
+              e, is_axially_symmetric, _integration_method)),
+          _element(e)
     {
     }
 
@@ -54,6 +56,7 @@ protected:
                 Eigen::aligned_allocator<
                     typename ShapeMatricesType::ShapeMatrices>> const
         _shape_matrices;
+    MeshLib::Element const& _element;
 };
 
 }  // namespace ProcessLib

ProcessLib/BoundaryCondition/NeumannBoundaryConditionLocalAssembler.h

@@ -11,7 +11,9 @@
 
 #include "GenericNaturalBoundaryConditionLocalAssembler.h"
 #include "NumLib/DOF/DOFTableUtil.h"
+#include "NumLib/Fem/ShapeMatrixPolicy.h"
 #include "ProcessLib/Parameter/Parameter.h"
+#include "ProcessLib/Utils/InitShapeMatrices.h"
 
 namespace ProcessLib
 {
@@ -23,6 +25,7 @@ class NeumannBoundaryConditionLocalAssembler final
 {
     using Base = GenericNaturalBoundaryConditionLocalAssembler<
         ShapeFunction, IntegrationMethod, GlobalDim>;
+    using ShapeMatricesType = ShapeMatrixPolicyType<ShapeFunction, GlobalDim>;
 
 public:
     /// The neumann_bc_value factor is directly integrated into the local
@@ -53,11 +56,20 @@ public:
         SpatialPosition pos;
         pos.setElementID(id);
 
+        using FemType =
+            NumLib::TemplateIsoparametric<ShapeFunction, ShapeMatricesType>;
+        FemType fe(
+             static_cast<const typename ShapeFunction::MeshElement&>(Base::_element));
+
         for (unsigned ip = 0; ip < n_integration_points; ip++)
         {
             pos.setIntegrationPoint(ip);
             auto const& sm = Base::_shape_matrices[ip];
             auto const& wp = Base::_integration_method.getWeightedPoint(ip);
+
+            MathLib::TemplatePoint<double, 3> coords_ip(fe.interpolateCoordinates(sm.N));
+            pos.setCoordinates(coords_ip);
+
             _local_rhs.noalias() += sm.N * _neumann_bc_parameter(t, pos)[0] *
                                     sm.detJ * wp.getWeight() *
                                     sm.integralMeasure;