mesh_function_fe.h

 
#ifndef INCG4ee2d6e8004446558bc6d2186596e392
#define INCG4ee2d6e8004446558bc6d2186596e392
 
#include <memory>
 
#include "scalar_fe_space.h"
 
namespace lf::fe {

template <class SCALAR_FE, class SCALAR_COEFF>
class MeshFunctionFE {
 public:
  // Why this? Because we can use real-valued finite element spaces to represent
  // complex-valued functions by using complex degrees of freedom!
  using Scalar = decltype(SCALAR_FE(0) * SCALAR_COEFF(0));  // NOLINT

  MeshFunctionFE(std::shared_ptr<const ScalarFESpace<SCALAR_FE>> fe_space,
                 Eigen::Matrix<SCALAR_COEFF, Eigen::Dynamic, 1> coeff_vector)
      : fe_space_(std::move(fe_space)), dof_vector_(std::move(coeff_vector)) {}

  std::vector<Scalar> operator()(const lf::mesh::Entity& e,
                                 const Eigen::MatrixXd& local) const {
    // Obtain values of all shape functions in the evaluation points
    auto sf_eval =
        fe_space_->ShapeFunctionLayout(e)->EvalReferenceShapeFunctions(local);
    // Extract the d.o.f.s for the current entity from the vector of global
    // d.o.f. values
    Eigen::Matrix<SCALAR_COEFF, 1, Eigen::Dynamic> local_dofs(1,
                                                              sf_eval.rows());
    auto global_dofs = fe_space_->LocGlobMap().GlobalDofIndices(e);
    for (Eigen::Index i = 0; i < sf_eval.rows(); ++i) {
      local_dofs(i) = dof_vector_(global_dofs[i]);
    }
    // Trick to combine Eigen data types with STL containers
    std::vector<Scalar> result(local.cols());
    Eigen::Map<Eigen::Matrix<Scalar, 1, Eigen::Dynamic>> temp(result.data(), 1,
                                                              local.cols());
    temp = local_dofs * sf_eval;
    return result;
  }

  [[nodiscard]] std::shared_ptr<const mesh::Mesh> getMesh() const {
    return fe_space_->Mesh();
  }

  [[nodiscard]] std::shared_ptr<const ScalarFESpace<SCALAR_FE>> getFESpace()
      const {
    return fe_space_;
  }
 
 private:
  std::shared_ptr<const ScalarFESpace<SCALAR_FE>> fe_space_;
  Eigen::Matrix<SCALAR_COEFF, Eigen::Dynamic, 1> dof_vector_;
};
 
// deduction guide
template <class T, class SCALAR_COEFF>
MeshFunctionFE(std::shared_ptr<T>,
               const Eigen::Matrix<SCALAR_COEFF, Eigen::Dynamic, 1>&)
    -> MeshFunctionFE<typename T::Scalar, SCALAR_COEFF>;
 
}  // namespace lf::fe
 
#endif  // INCG4ee2d6e8004446558bc6d2186596e392