fe_enriched.h

Go to the documentation of this file.

// ---------------------------------------------------------------------
//
// Copyright (C) 2016 - 2020 by the deal.II authors
//
// This file is part of the deal.II library.
//
// The deal.II library is free software; you can use it, redistribute
// it, and/or modify it under the terms of the GNU Lesser General
// Public License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// The full text of the license can be found in the file LICENSE.md at
// the top level directory of deal.II.
//
// ---------------------------------------------------------------------
 
#ifndef dealii_fe_enriched_h
#define dealii_fe_enriched_h
 
#include <deal.II/base/config.h>
 
#include <deal.II/base/function.h>
#include <deal.II/base/quadrature.h>
#include <deal.II/base/symmetric_tensor.h>
 
#include <deal.II/fe/fe.h>
#include <deal.II/fe/fe_nothing.h>
#include <deal.II/fe/fe_q.h>
#include <deal.II/fe/fe_system.h>
#include <deal.II/fe/fe_update_flags.h>
 
#include <deal.II/hp/dof_handler.h>
#include <deal.II/hp/fe_collection.h>
 
#include <map>
#include <numeric>
#include <utility>
#include <vector>
 
DEAL_II_NAMESPACE_OPEN
 
template <int dim, int spacedim = dim>
class FE_Enriched : public FiniteElement<dim, spacedim>
{
public:
  FE_Enriched(const FiniteElement<dim, spacedim> &fe_base,
              const FiniteElement<dim, spacedim> &fe_enriched,
              const Function<spacedim> *          enrichment_function);
 
  FE_Enriched(const FiniteElement<dim, spacedim> &fe_base);
 
  FE_Enriched(
    const FiniteElement<dim, spacedim> *                     fe_base,
    const std::vector<const FiniteElement<dim, spacedim> *> &fe_enriched,
    const std::vector<std::vector<std::function<const Function<spacedim> *(
      const typename Triangulation<dim, spacedim>::cell_iterator &)>>>
      &functions);
 
private:
  FE_Enriched(
    const std::vector<const FiniteElement<dim, spacedim> *> &fes,
    const std::vector<unsigned int> &                        multiplicities,
    const std::vector<std::vector<std::function<const Function<spacedim> *(
      const typename Triangulation<dim, spacedim>::cell_iterator &)>>>
      &functions);
 
public:
  virtual std::unique_ptr<FiniteElement<dim, spacedim>>
  clone() const override;
 
  virtual UpdateFlags
  requires_update_flags(const UpdateFlags update_flags) const override;
 
  virtual std::string
  get_name() const override;
 
  virtual const FiniteElement<dim, spacedim> &
  base_element(const unsigned int index) const override;
 
  virtual double
  shape_value(const unsigned int i, const Point<dim> &p) const override;
 
  virtual const FullMatrix<double> &
  get_restriction_matrix(
    const unsigned int         child,
    const RefinementCase<dim> &refinement_case =
      RefinementCase<dim>::isotropic_refinement) const override;
 
  virtual const FullMatrix<double> &
  get_prolongation_matrix(
    const unsigned int         child,
    const RefinementCase<dim> &refinement_case =
      RefinementCase<dim>::isotropic_refinement) const override;
 
 
  virtual bool
  hp_constraints_are_implemented() const override;
 
  virtual void
  get_face_interpolation_matrix(const FiniteElement<dim, spacedim> &source,
                                FullMatrix<double> &matrix) const override;
 
  virtual void
  get_subface_interpolation_matrix(const FiniteElement<dim, spacedim> &source,
                                   const unsigned int                  subface,
                                   FullMatrix<double> &matrix) const override;
 
  virtual std::vector<std::pair<unsigned int, unsigned int>>
  hp_vertex_dof_identities(
    const FiniteElement<dim, spacedim> &fe_other) const override;
 
  virtual std::vector<std::pair<unsigned int, unsigned int>>
  hp_line_dof_identities(
    const FiniteElement<dim, spacedim> &fe_other) const override;
 
  virtual std::vector<std::pair<unsigned int, unsigned int>>
  hp_quad_dof_identities(
    const FiniteElement<dim, spacedim> &fe_other) const override;
 
  virtual FiniteElementDomination::Domination
  compare_for_domination(const FiniteElement<dim, spacedim> &fe_other,
                         const unsigned int codim = 0) const override final;
 
 
 
  const std::vector<std::vector<std::function<const Function<spacedim> *(
    const typename Triangulation<dim, spacedim>::cell_iterator &)>>>
  get_enrichments() const;
 
  const FESystem<dim, spacedim> &
  get_fe_system() const;
 
protected:
  class InternalData : public FiniteElement<dim, spacedim>::InternalDataBase
  {
  public:
    struct EnrichmentValues
    {
      std::vector<double>                       values;
      std::vector<Tensor<1, spacedim>>          gradients;
      std::vector<SymmetricTensor<2, spacedim>> hessians;
    };
 
    InternalData(std::unique_ptr<typename FESystem<dim, spacedim>::InternalData>
                   fesystem_data);
 
    typename FiniteElement<dim, spacedim>::InternalDataBase &
    get_fe_data(const unsigned int base_no) const;
 
    internal::FEValuesImplementation::FiniteElementRelatedData<dim, spacedim> &
    get_fe_output_object(const unsigned int base_no) const;
 
    std::unique_ptr<typename FESystem<dim, spacedim>::InternalData>
      fesystem_data;
 
    mutable std::vector<std::vector<EnrichmentValues>> enrichment;
  };
 
  std::vector<std::vector<std::vector<unsigned int>>>
    base_no_mult_local_enriched_dofs;
 
  const std::vector<std::vector<std::function<const Function<spacedim> *(
    const typename Triangulation<dim, spacedim>::cell_iterator &)>>>
    enrichments;
 
  const bool is_enriched;
 
  template <int dim_1>
  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
  setup_data(
    std::unique_ptr<typename FESystem<dim, spacedim>::InternalData> fes_data,
    const UpdateFlags                                               flags,
    const Quadrature<dim_1> &quadrature) const;
 
  virtual std::unique_ptr<
    typename FiniteElement<dim, spacedim>::InternalDataBase>
  get_data(
    const UpdateFlags             flags,
    const Mapping<dim, spacedim> &mapping,
    const Quadrature<dim> &       quadrature,
    ::internal::FEValuesImplementation::FiniteElementRelatedData<dim,
                                                                       spacedim>
      &output_data) const override;
 
  virtual std::unique_ptr<
    typename FiniteElement<dim, spacedim>::InternalDataBase>
  get_face_data(
    const UpdateFlags             update_flags,
    const Mapping<dim, spacedim> &mapping,
    const Quadrature<dim - 1> &   quadrature,
    ::internal::FEValuesImplementation::FiniteElementRelatedData<dim,
                                                                       spacedim>
      &output_data) const override;
 
  virtual std::unique_ptr<
    typename FiniteElement<dim, spacedim>::InternalDataBase>
  get_subface_data(
    const UpdateFlags             update_flags,
    const Mapping<dim, spacedim> &mapping,
    const Quadrature<dim - 1> &   quadrature,
    ::internal::FEValuesImplementation::FiniteElementRelatedData<dim,
                                                                       spacedim>
      &output_data) const override;
 
  virtual void
  fill_fe_values(
    const typename Triangulation<dim, spacedim>::cell_iterator &cell,
    const CellSimilarity::Similarity                            cell_similarity,
    const Quadrature<dim> &                                     quadrature,
    const Mapping<dim, spacedim> &                              mapping,
    const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
    const ::internal::FEValuesImplementation::MappingRelatedData<dim,
                                                                       spacedim>
      &                                                            mapping_data,
    const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
    ::internal::FEValuesImplementation::FiniteElementRelatedData<dim,
                                                                       spacedim>
      &output_data) const override;
 
  virtual void
  fill_fe_face_values(
    const typename Triangulation<dim, spacedim>::cell_iterator &cell,
    const unsigned int                                          face_no,
    const Quadrature<dim - 1> &                                 quadrature,
    const Mapping<dim, spacedim> &                              mapping,
    const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
    const ::internal::FEValuesImplementation::MappingRelatedData<dim,
                                                                       spacedim>
      &                                                            mapping_data,
    const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
    ::internal::FEValuesImplementation::FiniteElementRelatedData<dim,
                                                                       spacedim>
      &output_data) const override;
 
  virtual void
  fill_fe_subface_values(
    const typename Triangulation<dim, spacedim>::cell_iterator &cell,
    const unsigned int                                          face_no,
    const unsigned int                                          sub_no,
    const Quadrature<dim - 1> &                                 quadrature,
    const Mapping<dim, spacedim> &                              mapping,
    const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
    const ::internal::FEValuesImplementation::MappingRelatedData<dim,
                                                                       spacedim>
      &                                                            mapping_data,
    const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
    ::internal::FEValuesImplementation::FiniteElementRelatedData<dim,
                                                                       spacedim>
      &output_data) const override;
 
private:
  void
  initialize(const std::vector<const FiniteElement<dim, spacedim> *> &fes,
             const std::vector<unsigned int> &multiplicities);
 
  const std::unique_ptr<const FESystem<dim, spacedim>> fe_system;
 
  template <int dim_1>
  void
  multiply_by_enrichment(
    const Quadrature<dim_1> &quadrature,
    const InternalData &     fe_data,
    const internal::FEValuesImplementation::MappingRelatedData<dim, spacedim>
      &                                                         mapping_data,
    const typename Triangulation<dim, spacedim>::cell_iterator &cell,
    internal::FEValuesImplementation::FiniteElementRelatedData<dim, spacedim>
      &output_data) const;
};
 
 
 
namespace ColorEnriched
{
  template <int dim, int spacedim = dim>
  using predicate_function = std::function<bool(
    const typename Triangulation<dim, spacedim>::cell_iterator &)>;
 
#ifndef DOXYGEN
  namespace internal
  {
    template <int dim, int spacedim>
    bool
    find_connection_between_subdomains(
      const hp::DoFHandler<dim, spacedim> &    dof_handler,
      const predicate_function<dim, spacedim> &predicate_1,
      const predicate_function<dim, spacedim> &predicate_2);
 
    template <int dim, int spacedim>
    unsigned int
    color_predicates(
      const hp::DoFHandler<dim, spacedim> &                 dof_handler,
      const std::vector<predicate_function<dim, spacedim>> &predicates,
      std::vector<unsigned int> &                           predicate_colors);
 
    template <int dim, int spacedim>
    void
    set_cellwise_color_set_and_fe_index(
      hp::DoFHandler<dim, spacedim> &                       dof_handler,
      const std::vector<predicate_function<dim, spacedim>> &predicates,
      const std::vector<unsigned int> &                     predicate_colors,
      std::map<unsigned int, std::map<unsigned int, unsigned int>>
        &                                  cellwise_color_predicate_map,
      std::vector<std::set<unsigned int>> &fe_sets);
 
    template <int dim, int spacedim>
    void
    make_colorwise_enrichment_functions(
      const unsigned int                                      n_colors,
      const std::vector<std::shared_ptr<Function<spacedim>>> &enrichments,
      const std::map<unsigned int, std::map<unsigned int, unsigned int>>
        &cellwise_color_predicate_map,
      std::vector<std::function<const Function<spacedim> *(
        const typename Triangulation<dim, spacedim>::cell_iterator &)>>
        &color_enrichments);
 
 
    template <int dim, int spacedim>
    void
    make_fe_collection_from_colored_enrichments(
      const unsigned int n_colors,
      const std::vector<std::set<unsigned int>>
        &fe_sets, // total list of color sets possible
      const std::vector<std::function<const Function<spacedim> *(
        const typename Triangulation<dim, spacedim>::cell_iterator &)>>
        &color_enrichments, // color wise enrichment functions
      const FiniteElement<dim, spacedim> &fe_base, // basic FE element
      const FiniteElement<dim, spacedim>
        &fe_enriched, // FE multiplied by enrichment function
      const FE_Nothing<dim, spacedim> &fe_nothing,
      hp::FECollection<dim, spacedim> &fe_collection);
  }    // namespace internal
#endif // DOXYGEN
 
 
 
  template <int dim, int spacedim = dim>
  struct Helper
  {
    Helper(const FiniteElement<dim, spacedim> &                    fe_base,
           const FiniteElement<dim, spacedim> &                    fe_enriched,
           const std::vector<predicate_function<dim, spacedim>> &  predicates,
           const std::vector<std::shared_ptr<Function<spacedim>>> &enrichments);
 
    const hp::FECollection<dim, spacedim> &
    build_fe_collection(hp::DoFHandler<dim, spacedim> &dof_handler);
 
  private:
    hp::FECollection<dim, spacedim> fe_collection;
 
    const FiniteElement<dim, spacedim> &fe_base;
 
    const FiniteElement<dim, spacedim> &fe_enriched;
 
    const FE_Nothing<dim, spacedim> fe_nothing;
 
    const std::vector<predicate_function<dim, spacedim>> predicates;
 
    const std::vector<std::shared_ptr<Function<spacedim>>> enrichments;
 
    using cell_iterator_function = std::function<const Function<spacedim> *(
      const typename Triangulation<dim, spacedim>::cell_iterator &)>;
 
    std::vector<cell_iterator_function> color_enrichments;
 
    std::vector<unsigned int> predicate_colors;
 
    unsigned int n_colors;
 
    std::map<unsigned int, std::map<unsigned int, unsigned int>>
      cellwise_color_predicate_map;
 
    std::vector<std::set<unsigned int>> fe_sets;
  };
} // namespace ColorEnriched
 
DEAL_II_NAMESPACE_CLOSE
 
#endif // dealii_fe_enriched_h