product_fe_space_factory.h

#ifndef PROJECTS_DPG_PRODUCT_FE_SPACE_FACTORY
#define PROJECTS_DPG_PRODUCT_FE_SPACE_FACTORY
 
#include <lf/base/base.h>
#include <lf/fe/fe.h>
#include <lf/mesh/mesh.h>
#include <lf/uscalfe/uscalfe.h>
 
#include "discontinuous_fe_constant.h"
#include "discontinuous_scalar_reference_finite_element.h"
#include "dpg.h"
#include "product_fe_space.h"
#include "trace_scalar_reference_finite_element.h"
 
namespace projects::dpg {
 
template <typename SCALAR>
class ProductUniformFESpaceFactory {
 public:
  ProductUniformFESpaceFactory(const ProductUniformFESpaceFactory&) = delete;
  ProductUniformFESpaceFactory(ProductUniformFESpaceFactory&&) noexcept =
      delete;
  ProductUniformFESpaceFactory& operator=(const ProductUniformFESpaceFactory&) =
      delete;
  ProductUniformFESpaceFactory& operator=(ProductUniformFESpaceFactory&&) =
      delete;
 
  explicit ProductUniformFESpaceFactory(
      std::shared_ptr<const lf::mesh::Mesh> mesh_p)
      : mesh_p_(std::move(mesh_p)) {}
 
  size_type AddH1Component(size_type degree);
 
  size_type AddL2Component(size_type degree);
 
  size_type AddFluxComponent(size_type degree);
 
  size_type AddTraceComponent(size_type degree);
 
  std::shared_ptr<ProductUniformFESpace<SCALAR>> Build();
 
  size_type NumComponents() { return num_components_; }
 
  virtual ~ProductUniformFESpaceFactory() = default;
 
 private:
  std::shared_ptr<const lf::fe::ScalarReferenceFiniteElement<SCALAR>>
  ReferenceFiniteElement(lf::base::RefEl ref_el, size_type degree);
 
  std::shared_ptr<const lf::mesh::Mesh> mesh_p_;
 
  std::vector<
      std::shared_ptr<const lf::fe::ScalarReferenceFiniteElement<SCALAR>>>
      rfs_tria_v_;
  std::vector<
      std::shared_ptr<const lf::fe::ScalarReferenceFiniteElement<SCALAR>>>
      rfs_quad_v_;
  std::vector<
      std::shared_ptr<const lf::fe::ScalarReferenceFiniteElement<SCALAR>>>
      rfs_segment_v_;
  std::vector<
      std::shared_ptr<const lf::fe::ScalarReferenceFiniteElement<SCALAR>>>
      rfs_point_v_;
 
  size_type num_components_ = 0;
};
 
template <typename SCALAR>
std::shared_ptr<const lf::fe::ScalarReferenceFiniteElement<SCALAR>>
ProductUniformFESpaceFactory<SCALAR>::ReferenceFiniteElement(
    lf::base::RefEl ref_el, size_type degree) {
  switch (ref_el) {
    case (lf::base::RefEl::kTria()): {
      switch (degree) {
        case (0):
          return std::make_shared<FeDiscontinuousO0Tria<SCALAR>>();
        case (1):
          return std::make_shared<lf::uscalfe::FeLagrangeO1Tria<SCALAR>>();
        case (2):
          return std::make_shared<lf::uscalfe::FeLagrangeO2Tria<SCALAR>>();
        case (3):
          return std::make_shared<lf::uscalfe::FeLagrangeO3Tria<SCALAR>>();
        default:
          LF_ASSERT_MSG(false, "no fe of degree " << degree
                                                  << " available on ref_el "
                                                  << ref_el);
      }
    }
    case (lf::base::RefEl::kQuad()): {
      switch (degree) {
        case (0):
          return std::make_shared<FeDiscontinuousO0Quad<SCALAR>>();
        case (1):
          return std::make_shared<lf::uscalfe::FeLagrangeO1Quad<SCALAR>>();
        case (2):
          return std::make_shared<lf::uscalfe::FeLagrangeO2Quad<SCALAR>>();
        case (3):
          return std::make_shared<lf::uscalfe::FeLagrangeO3Quad<SCALAR>>();
        default:
          LF_ASSERT_MSG(false, "no fe of degree " << degree
                                                  << "available on ref_el "
                                                  << ref_el);
      }
    }
    case (lf::base::RefEl::kSegment()): {
      switch (degree) {
        case (0):
          return std::make_shared<FeDiscontinuousO0Segment<SCALAR>>();
        case (1):
          return std::make_shared<lf::uscalfe::FeLagrangeO1Segment<SCALAR>>();
        case (2):
          return std::make_shared<lf::uscalfe::FeLagrangeO2Segment<SCALAR>>();
        case (3):
          return std::make_shared<lf::uscalfe::FeLagrangeO3Segment<SCALAR>>();
        default:
          LF_ASSERT_MSG(false, "no fe of degree " << degree
                                                  << "available on ref_el "
                                                  << ref_el);
      }
    }
    case (lf::base::RefEl::kPoint()): {
      return std::make_shared<lf::fe::FePoint<SCALAR>>(degree);
    }
    default: {
      LF_ASSERT_MSG(false, "unsupported reference element " << ref_el);
      return nullptr;
    }
  }
}
 
template <typename SCALAR>
size_type ProductUniformFESpaceFactory<SCALAR>::AddH1Component(
    size_type degree) {
  LF_ASSERT_MSG(
      degree > 0,
      "Can only represent H1 conforming finite element with degree > 0");
  // add continuous reference shape functions on triangles, quads and segments.
  rfs_tria_v_.push_back(
      ReferenceFiniteElement(lf::base::RefEl::kTria(), degree));
  rfs_quad_v_.push_back(
      ReferenceFiniteElement(lf::base::RefEl::kQuad(), degree));
  rfs_segment_v_.push_back(
      ReferenceFiniteElement(lf::base::RefEl::kSegment(), degree));
  rfs_point_v_.push_back(nullptr);
  return num_components_++;
}
 
template <typename SCALAR>
size_type ProductUniformFESpaceFactory<SCALAR>::AddL2Component(
    size_type degree) {
  // add discontinuous shape functions to triangles and quads.
  rfs_tria_v_.push_back(
      std::make_shared<DiscontinuousScalarReferenceFiniteElement<SCALAR>>(
          ReferenceFiniteElement(lf::base::RefEl::kTria(), degree)));
  rfs_quad_v_.push_back(
      std::make_shared<DiscontinuousScalarReferenceFiniteElement<SCALAR>>(
          ReferenceFiniteElement(lf::base::RefEl::kQuad(), degree)));
  rfs_segment_v_.push_back(nullptr);
  rfs_point_v_.push_back(nullptr);
  return num_components_++;
}
 
template <typename SCALAR>
size_type ProductUniformFESpaceFactory<SCALAR>::AddFluxComponent(
    size_type degree) {
  // add discontinuous shape functions to segments
  rfs_tria_v_.push_back(nullptr);
  rfs_quad_v_.push_back(nullptr);
  rfs_segment_v_.push_back(
      std::make_shared<DiscontinuousScalarReferenceFiniteElement<SCALAR>>(
          ReferenceFiniteElement(lf::base::RefEl::kSegment(), degree)));
  rfs_point_v_.push_back(nullptr);
  return num_components_++;
}
 
template <typename SCALAR>
size_type ProductUniformFESpaceFactory<SCALAR>::AddTraceComponent(
    size_type degree) {
  // add trace representing shape functions to triangles and quadrilaterals.
  // add full shape function spec. to segments.
  rfs_tria_v_.push_back(
      std::make_shared<TraceScalarReferenceFiniteElement<SCALAR>>(
          ReferenceFiniteElement(lf::base::RefEl::kTria(), degree)));
  rfs_quad_v_.push_back(
      std::make_shared<TraceScalarReferenceFiniteElement<SCALAR>>(
          ReferenceFiniteElement(lf::base::RefEl::kQuad(), degree)));
  rfs_segment_v_.push_back(
      ReferenceFiniteElement(lf::base::RefEl::kSegment(), degree));
  rfs_point_v_.push_back(nullptr);
  return num_components_++;
}
 
template <typename SCALAR>
std::shared_ptr<ProductUniformFESpace<SCALAR>>
ProductUniformFESpaceFactory<SCALAR>::Build() {
  // construct the fe space based on the provided local reference shape function
  // layouts.
  auto fe_space = std::make_shared<ProductUniformFESpace<SCALAR>>(
      mesh_p_, rfs_tria_v_, rfs_quad_v_, rfs_segment_v_, rfs_point_v_);
 
  // clear supplied information:
  rfs_tria_v_ = {};
  rfs_quad_v_ = {};
  rfs_segment_v_ = {};
  rfs_point_v_ = {};
  return fe_space;
}
 
}  // namespace projects::dpg
#endif  // PROJECTS_DPG_PRODUCT_FE_SPACE_FACTORY