fe_system.h

// ---------------------------------------------------------------------
//
// Copyright (C) 1999 - 2016 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 at
// the top level of the deal.II distribution.
//
// ---------------------------------------------------------------------

#ifndef dealii__fe_system_h
#define dealii__fe_system_h


/*----------------------------   fe_system.h     ---------------------------*/


#include <deal.II/base/config.h>
#include <deal.II/base/thread_management.h>
#include <deal.II/fe/fe.h>
#include <vector>
#include <utility>

DEAL_II_NAMESPACE_OPEN

template <int dim, int spacedim> class FE_Enriched;


template <int dim, int spacedim=dim>
class FESystem : public FiniteElement<dim,spacedim>
{
public:

  FESystem (const FiniteElement<dim,spacedim> &fe,
            const unsigned int n_elements);

  FESystem (const FiniteElement<dim,spacedim> &fe1, const unsigned int n1,
            const FiniteElement<dim,spacedim> &fe2, const unsigned int n2);

  FESystem (const FiniteElement<dim,spacedim> &fe1, const unsigned int n1,
            const FiniteElement<dim,spacedim> &fe2, const unsigned int n2,
            const FiniteElement<dim,spacedim> &fe3, const unsigned int n3);

  FESystem (const FiniteElement<dim,spacedim> &fe1, const unsigned int n1,
            const FiniteElement<dim,spacedim> &fe2, const unsigned int n2,
            const FiniteElement<dim,spacedim> &fe3, const unsigned int n3,
            const FiniteElement<dim,spacedim> &fe4, const unsigned int n4);

  FESystem (const FiniteElement<dim,spacedim> &fe1, const unsigned int n1,
            const FiniteElement<dim,spacedim> &fe2, const unsigned int n2,
            const FiniteElement<dim,spacedim> &fe3, const unsigned int n3,
            const FiniteElement<dim,spacedim> &fe4, const unsigned int n4,
            const FiniteElement<dim,spacedim> &fe5, const unsigned int n5);

  FESystem (const std::vector<const FiniteElement<dim,spacedim>*> &fes,
            const std::vector<unsigned int>                   &multiplicities);

  virtual ~FESystem ();

  virtual std::string get_name () const;

  // for documentation, see the FiniteElement base class
  virtual
  UpdateFlags
  requires_update_flags (const UpdateFlags update_flags) const;

  virtual double shape_value (const unsigned int i,
                              const Point<dim> &p) const;

  virtual double shape_value_component (const unsigned int i,
                                        const Point<dim> &p,
                                        const unsigned int component) const;

  virtual Tensor<1,dim> shape_grad (const unsigned int  i,
                                    const Point<dim>   &p) const;

  virtual Tensor<1,dim> shape_grad_component (const unsigned int i,
                                              const Point<dim> &p,
                                              const unsigned int component) const;

  virtual Tensor<2,dim> shape_grad_grad (const unsigned int  i,
                                         const Point<dim> &p) const;

  virtual
  Tensor<2,dim>
  shape_grad_grad_component (const unsigned int i,
                             const Point<dim> &p,
                             const unsigned int component) const;

  virtual Tensor<3,dim> shape_3rd_derivative (const unsigned int  i,
                                              const Point<dim>   &p) const;

  virtual Tensor<3,dim> shape_3rd_derivative_component (const unsigned int i,
                                                        const Point<dim>   &p,
                                                        const unsigned int component) const;

  virtual Tensor<4,dim> shape_4th_derivative (const unsigned int  i,
                                              const Point<dim>   &p) const;

  virtual Tensor<4,dim> shape_4th_derivative_component (const unsigned int i,
                                                        const Point<dim>   &p,
                                                        const unsigned int component) const;

  virtual void
  get_interpolation_matrix (const FiniteElement<dim,spacedim> &source,
                            FullMatrix<double>           &matrix) const;

  virtual const FiniteElement<dim,spacedim> &
  base_element (const unsigned int index) const;

  virtual bool has_support_on_face (const unsigned int shape_index,
                                    const unsigned int face_index) const;

  virtual const FullMatrix<double> &
  get_restriction_matrix (const unsigned int child,
                          const RefinementCase<dim> &refinement_case=RefinementCase<dim>::isotropic_refinement) const;

  virtual const FullMatrix<double> &
  get_prolongation_matrix (const unsigned int child,
                           const RefinementCase<dim> &refinement_case=RefinementCase<dim>::isotropic_refinement) const;

  virtual
  unsigned int face_to_cell_index (const unsigned int face_dof_index,
                                   const unsigned int face,
                                   const bool face_orientation = true,
                                   const bool face_flip        = false,
                                   const bool face_rotation    = false) const;

  virtual
  Point<dim>
  unit_support_point (const unsigned int index) const;

  virtual
  Point<dim-1>
  unit_face_support_point (const unsigned int index) const;

  virtual std::pair<Table<2,bool>, std::vector<unsigned int> >
  get_constant_modes () const;

  virtual bool hp_constraints_are_implemented () const;

  virtual void
  get_face_interpolation_matrix (const FiniteElement<dim,spacedim> &source,
                                 FullMatrix<double>       &matrix) const;


  virtual void
  get_subface_interpolation_matrix (const FiniteElement<dim,spacedim> &source,
                                    const unsigned int        subface,
                                    FullMatrix<double>       &matrix) const;

  virtual
  std::vector<std::pair<unsigned int, unsigned int> >
  hp_vertex_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;

  virtual
  std::vector<std::pair<unsigned int, unsigned int> >
  hp_line_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;

  virtual
  std::vector<std::pair<unsigned int, unsigned int> >
  hp_quad_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;

  virtual
  FiniteElementDomination::Domination
  compare_for_face_domination (const FiniteElement<dim,spacedim> &fe_other) const;

  virtual std::size_t memory_consumption () const;

protected:

  virtual FiniteElement<dim,spacedim> *clone() const;


  virtual typename FiniteElement<dim,spacedim>::InternalDataBase *
  get_data (const UpdateFlags                                                    update_flags,
            const Mapping<dim,spacedim>                                         &mapping,
            const Quadrature<dim>                                               &quadrature,
            ::internal::FEValues::FiniteElementRelatedData<dim, spacedim> &output_data) const;

  virtual
  typename FiniteElement<dim,spacedim>::InternalDataBase *
  get_face_data (const UpdateFlags                                                    update_flags,
                 const Mapping<dim,spacedim>                                         &mapping,
                 const Quadrature<dim-1>                                             &quadrature,
                 ::internal::FEValues::FiniteElementRelatedData<dim, spacedim> &output_data) const;

  virtual
  typename FiniteElement<dim,spacedim>::InternalDataBase *
  get_subface_data (const UpdateFlags                                                    update_flags,
                    const Mapping<dim,spacedim>                                         &mapping,
                    const Quadrature<dim-1>                                             &quadrature,
                    ::internal::FEValues::FiniteElementRelatedData<dim, spacedim> &output_data) const;

  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::FEValues::MappingRelatedData<dim, spacedim> &mapping_data,
                  const typename FiniteElement<dim,spacedim>::InternalDataBase        &fe_internal,
                  ::internal::FEValues::FiniteElementRelatedData<dim, spacedim> &output_data) const;

  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::FEValues::MappingRelatedData<dim, spacedim> &mapping_data,
                       const typename FiniteElement<dim,spacedim>::InternalDataBase        &fe_internal,
                       ::internal::FEValues::FiniteElementRelatedData<dim, spacedim> &output_data) const;

  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::FEValues::MappingRelatedData<dim, spacedim> &mapping_data,
                          const typename FiniteElement<dim,spacedim>::InternalDataBase        &fe_internal,
                          ::internal::FEValues::FiniteElementRelatedData<dim, spacedim> &output_data) const;

  template <int dim_1>
  void compute_fill (const Mapping<dim,spacedim>                      &mapping,
                     const typename Triangulation<dim,spacedim>::cell_iterator &cell,
                     const unsigned int                                face_no,
                     const unsigned int                                sub_no,
                     const Quadrature<dim_1>                          &quadrature,
                     const CellSimilarity::Similarity                   cell_similarity,
                     const typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
                     const typename FiniteElement<dim,spacedim>::InternalDataBase &fe_data,
                     const internal::FEValues::MappingRelatedData<dim,spacedim> &mapping_data,
                     internal::FEValues::FiniteElementRelatedData<dim,spacedim> &output_data) const;

private:

  static const unsigned int invalid_face_number = numbers::invalid_unsigned_int;

  std::vector<std::pair<std_cxx11::shared_ptr<const FiniteElement<dim,spacedim> >,
      unsigned int> >
      base_elements;

  void initialize_unit_support_points ();

  void initialize_unit_face_support_points ();

  void initialize_quad_dof_index_permutation ();
  void initialize (const std::vector<const FiniteElement<dim,spacedim>*> &fes,
                   const std::vector<unsigned int> &multiplicities);

  void build_interface_constraints ();

  template <int structdim>
  std::vector<std::pair<unsigned int, unsigned int> >
  hp_object_dof_identities (const FiniteElement<dim,spacedim> &fe_other) const;

  class InternalData : public FiniteElement<dim,spacedim>::InternalDataBase
  {
  public:
    InternalData (const unsigned int n_base_elements);

    ~InternalData();

    void set_fe_data(const unsigned int                        base_no,
                     typename FiniteElement<dim,spacedim>::InternalDataBase *);

    typename FiniteElement<dim,spacedim>::InternalDataBase &
    get_fe_data (const unsigned int base_no) const;

    internal::FEValues::FiniteElementRelatedData<dim,spacedim> &
    get_fe_output_object (const unsigned int base_no) const;

  private:

    typename std::vector<typename FiniteElement<dim,spacedim>::InternalDataBase *> base_fe_datas;

    mutable std::vector<internal::FEValues::FiniteElementRelatedData<dim,spacedim> > base_fe_output_objects;
  };

  /*
   * Mutex for protecting initialization of restriction and embedding matrix.
   */
  mutable Threads::Mutex mutex;

  friend class FE_Enriched<dim,spacedim>;
};


DEAL_II_NAMESPACE_CLOSE

/*----------------------------  fe_system.h  ---------------------------*/
#endif
/*----------------------------  fe_system.h  ---------------------------*/