trace_scalar_reference_finite_element.h

#ifndef PROJECTS_DPG_TRACE_SCALAR_REFERENCE_FINITE_ELEMENT
#define PROJECTS_DPG_TRACE_SCALAR_REFERENCE_FINITE_ELEMENT
 
#include <lf/fe/fe.h>
#include <lf/uscalfe/uscalfe.h>
 
#include <iostream>
#include <typeinfo>
 
#include "dpg.h"
 
namespace projects::dpg {
 
template <typename SCALAR>
class TraceScalarReferenceFiniteElement
    : public lf::fe::ScalarReferenceFiniteElement<SCALAR> {
 public:
  TraceScalarReferenceFiniteElement() = default;
 
  TraceScalarReferenceFiniteElement(const TraceScalarReferenceFiniteElement&) =
      delete;
  TraceScalarReferenceFiniteElement(
      TraceScalarReferenceFiniteElement&&) noexcept = default;
  TraceScalarReferenceFiniteElement& operator=(
      const TraceScalarReferenceFiniteElement&) = delete;
  TraceScalarReferenceFiniteElement& operator=(
      TraceScalarReferenceFiniteElement&&) noexcept = default;
 
  explicit TraceScalarReferenceFiniteElement(
      std::shared_ptr<const lf::fe::ScalarReferenceFiniteElement<SCALAR>> cfe)
      : lf::fe::ScalarReferenceFiniteElement<SCALAR>(), cfe_(std::move(cfe)) {}
 
  [[nodiscard]] inline bool isInitialized() const { return (cfe_ != nullptr); }
 
  [[nodiscard]] lf::base::RefEl RefEl() const override {
    LF_ASSERT_MSG(isInitialized(), "Not initialized!");
    return cfe_->RefEl();
  }
 
  [[nodiscard]] unsigned Degree() const override {
    LF_ASSERT_MSG(isInitialized(), "Not Initialized");
    return cfe_->Degree();
  }
 
  [[nodiscard]] size_type NumRefShapeFunctions() const override {
    LF_ASSERT_MSG(isInitialized(), "Not initialized");
    return cfe_->NumRefShapeFunctions() - cfe_->NumRefShapeFunctions(0);
  }
 
  [[nodiscard]] size_type NumRefShapeFunctions(dim_t codim) const override {
    LF_ASSERT_MSG(isInitialized(), "Not initialized");
    return (codim == 0) ? 0 : cfe_->NumRefShapeFunctions(codim);
  }
 
  [[nodiscard]] size_type NumRefShapeFunctions(
      dim_t codim, sub_idx_t subidx) const override {
    LF_ASSERT_MSG(isInitialized(), "Not initialized");
    return (codim == 0) ? 0 : cfe_->NumRefShapeFunctions(codim, subidx);
  }
 
  [[nodiscard]] Eigen::Matrix<SCALAR, Eigen::Dynamic, Eigen::Dynamic>
  EvalReferenceShapeFunctions(const Eigen::MatrixXd& local) const override {
    LF_ASSERT_MSG(isInitialized(), "Not initialized");
    return cfe_->EvalReferenceShapeFunctions(local).topRows(
        NumRefShapeFunctions());
  }
 
  [[nodiscard]] Eigen::Matrix<SCALAR, Eigen::Dynamic, Eigen::Dynamic>
  GradientsReferenceShapeFunctions(
      const Eigen::MatrixXd& local) const override {
    LF_ASSERT_MSG(isInitialized(), "Not initialized");
    return cfe_->GradientsReferenceShapeFunctions(local).topRows(
        NumRefShapeFunctions());
  }
 
  [[nodiscard]] Eigen::MatrixXd EvaluationNodes() const override {
    LF_ASSERT_MSG(isInitialized(), "Not initialized");
    return cfe_->EvaluationNodes().leftCols(NumRefShapeFunctions());
  }
 
  [[nodiscard]] lf::fe::size_type NumEvaluationNodes() const override {
    LF_ASSERT_MSG(isInitialized(), "Not initialized");
    return cfe_->NumEvaluationNodes() - cfe_->NumRefShapeFunctions(0);
  }
 
  [[nodiscard]] Eigen::Matrix<SCALAR, 1, Eigen::Dynamic> NodalValuesToDofs(
      const Eigen::Matrix<SCALAR, 1, Eigen::Dynamic>& nodvals) const override {
    LF_ASSERT_MSG(isInitialized(), "Not initialized");
    LF_ASSERT_MSG(nodvals.cols() == NumEvaluationNodes(),
                  "nodvals = " << nodvals << " <-> " << NumEvaluationNodes());
 
    // the inner element expects more nodal values, than the trace element
    // since the internal shape functions were dropped, so
    // the nodal values are extended by zeroes
    Eigen::Matrix<SCALAR, 1, Eigen::Dynamic> extendedNodvals(
        cfe_->NumEvaluationNodes());
    Eigen::Matrix<SCALAR, 1, Eigen::Dynamic> zeros(cfe_->NumEvaluationNodes() -
                                                   NumEvaluationNodes());
    zeros.setZero();
    extendedNodvals << nodvals, zeros;
    return cfe_->NodalValuesToDofs(extendedNodvals)
        .leftCols(NumRefShapeFunctions());
  }
 
  ~TraceScalarReferenceFiniteElement() override = default;
 
 private:
  std::shared_ptr<const lf::fe::ScalarReferenceFiniteElement<SCALAR>> cfe_;
};
 
}  // namespace projects::dpg
 
#endif  // PROJECTS_DPG_TRACE_SCALAR_REFERENCE_FINITE_ELEMENT