dof_accessor.h

// ---------------------------------------------------------------------
//
// Copyright (C) 1998 - 2018 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_dof_accessor_h
#define dealii_dof_accessor_h


#include <deal.II/base/config.h>
#include <deal.II/grid/tria_accessor.h>
#include <deal.II/dofs/dof_handler.h>
#include <deal.II/hp/dof_handler.h>

#include <vector>

DEAL_II_NAMESPACE_OPEN

template <typename number> class FullMatrix;
template <typename number> class Vector;
class ConstraintMatrix;

template <typename Accessor> class TriaRawIterator;

template <int, int> class FiniteElement;


namespace internal
{
  namespace DoFCellAccessorImplementation
  {
    struct Implementation;
  }

  namespace DoFHandlerImplementation
  {
    struct Implementation;
    namespace Policy
    {
      struct Implementation;
    }
  }

  namespace hp
  {
    namespace DoFHandlerImplementation
    {
      struct Implementation;
    }
  }
}

// note: the file dof_accessor.templates.h is included at the end of
// this file.  this includes a lot of templates and thus makes
// compilation slower, but at the same time allows for more aggressive
// inlining and thus faster code.


namespace internal
{
  namespace DoFAccessorImplementation
  {
    template <int structdim, int dim, int spacedim>
    struct Inheritance
    {
      typedef ::TriaAccessor<structdim,dim,spacedim> BaseClass;
    };


    template <int dim, int spacedim>
    struct Inheritance<dim,dim,spacedim>
    {
      typedef ::CellAccessor<dim,spacedim> BaseClass;
    };
  }
}


/* -------------------------------------------------------------------------- */



template <int structdim, typename DoFHandlerType, bool level_dof_access>
class DoFAccessor : public ::internal::DoFAccessorImplementation::Inheritance<structdim, DoFHandlerType::dimension, DoFHandlerType::space_dimension>::BaseClass
{
public:

  static const unsigned int dimension=DoFHandlerType::dimension;

  static const unsigned int space_dimension=DoFHandlerType::space_dimension;

  typedef
  typename ::internal::DoFAccessorImplementation::Inheritance<structdim, dimension, space_dimension>::BaseClass
  BaseClass;

  typedef DoFHandlerType AccessorData;

  DoFAccessor ();

  DoFAccessor (const Triangulation<DoFHandlerType::dimension,DoFHandlerType::space_dimension> *tria,
               const int             level,
               const int             index,
               const DoFHandlerType *dof_handler);

  template <int structdim2, int dim2, int spacedim2>
  DoFAccessor (const InvalidAccessor<structdim2,dim2,spacedim2> &);

  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
  DoFAccessor (const DoFAccessor<dim2, DoFHandlerType2, level_dof_access2> &);

  template <bool level_dof_access2>
  DoFAccessor(const DoFAccessor<structdim, DoFHandlerType, level_dof_access2> &);

  DoFAccessor<structdim,DoFHandlerType, level_dof_access> &
  operator = (const DoFAccessor<structdim,DoFHandlerType, level_dof_access> &da) = delete;

  const DoFHandlerType &
  get_dof_handler () const;

  template <bool level_dof_access2>
  void copy_from (const DoFAccessor<structdim, DoFHandlerType, level_dof_access2> &a);

  void copy_from (const TriaAccessorBase<structdim, DoFHandlerType::dimension, DoFHandlerType::space_dimension> &da);

  static bool is_level_cell();

  TriaIterator<DoFAccessor<structdim,DoFHandlerType, level_dof_access> >
  child (const unsigned int c) const;

  typename ::internal::DoFHandlerImplementation::Iterators<DoFHandlerType, level_dof_access>::line_iterator
  line (const unsigned int i) const;

  typename ::internal::DoFHandlerImplementation::Iterators<DoFHandlerType, level_dof_access>::quad_iterator
  quad (const unsigned int i) const;

  void get_dof_indices (std::vector<types::global_dof_index> &dof_indices,
                        const unsigned int fe_index = DoFHandlerType::default_fe_index) const;

  void get_mg_dof_indices (const int level,
                           std::vector<types::global_dof_index> &dof_indices,
                           const unsigned int fe_index = DoFHandlerType::default_fe_index) const;

  void set_mg_dof_indices (const int level,
                           const std::vector<types::global_dof_index> &dof_indices,
                           const unsigned int fe_index = DoFHandlerType::default_fe_index);

  types::global_dof_index vertex_dof_index
  (const unsigned int vertex,
   const unsigned int i,
   const unsigned int fe_index = DoFHandlerType::default_fe_index) const;

  types::global_dof_index mg_vertex_dof_index
  (const int level,
   const unsigned int vertex,
   const unsigned int i,
   const unsigned int fe_index = DoFHandlerType::default_fe_index) const;

  types::global_dof_index dof_index
  (const unsigned int i,
   const unsigned int fe_index = DoFHandlerType::default_fe_index) const;

  types::global_dof_index mg_dof_index (const int level, const unsigned int i) const;

  unsigned int
  n_active_fe_indices () const;

  unsigned int
  nth_active_fe_index (const unsigned int n) const;

  bool
  fe_index_is_active (const unsigned int fe_index) const;

  const FiniteElement<DoFHandlerType::dimension,DoFHandlerType::space_dimension> &
  get_fe (const unsigned int fe_index) const;

  DeclExceptionMsg (ExcInvalidObject, "This accessor object has not been "
                    "associated with any DoFHandler object.");
  DeclException0 (ExcVectorNotEmpty);
  DeclException0 (ExcVectorDoesNotMatch);
  DeclException0 (ExcMatrixDoesNotMatch);
  DeclException0 (ExcNotActive);
  DeclException0 (ExcCantCompareIterators);

protected:

  DoFHandlerType *dof_handler;
public:
  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
  bool operator == (const DoFAccessor<dim2,DoFHandlerType2,level_dof_access2> &) const;

  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
  bool operator != (const DoFAccessor<dim2,DoFHandlerType2,level_dof_access2> &) const;
protected:
  void set_dof_handler (DoFHandlerType *dh);

  void set_dof_index (const unsigned int i,
                      const types::global_dof_index index,
                      const unsigned int fe_index = DoFHandlerType::default_fe_index) const;

  void set_mg_dof_index (const int level, const unsigned int i, const types::global_dof_index index) const;

  void set_vertex_dof_index
  (const unsigned int            vertex,
   const unsigned int            i,
   const types::global_dof_index index,
   const unsigned int            fe_index = DoFHandlerType::default_fe_index) const;

  void set_mg_vertex_dof_index
  (const int level,
   const unsigned int vertex,
   const unsigned int i,
   const types::global_dof_index index,
   const unsigned int fe_index = DoFHandlerType::default_fe_index) const;

  template <typename> friend class TriaRawIterator;
  template <int, class, bool> friend class DoFAccessor;

private:

  template <int dim, int spacedim> friend class DoFHandler;
  template <int dim, int spacedim> friend class hp::DoFHandler;

  friend struct ::internal::DoFHandlerImplementation::Policy::Implementation;
  friend struct ::internal::DoFHandlerImplementation::Implementation;
  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
  friend struct ::internal::DoFAccessorImplementation::Implementation;
};



template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access>
class DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access> : public TriaAccessor<0,1,spacedim>
{
public:

  static const unsigned int dimension=1;

  static const unsigned int space_dimension=spacedim;

  typedef TriaAccessor<0,1,spacedim> BaseClass;

  typedef DoFHandlerType<1,spacedim> AccessorData;

  DoFAccessor ();

  DoFAccessor (const Triangulation<1,spacedim>                       *tria,
               const typename TriaAccessor<0,1,spacedim>::VertexKind  vertex_kind,
               const unsigned int                                     vertex_index,
               const DoFHandlerType<1,spacedim>                      *dof_handler);

  DoFAccessor (const Triangulation<1,spacedim> *,
               const                             int         = 0,
               const                             int         = 0,
               const DoFHandlerType<1,spacedim> *dof_handler = 0);

  template <int structdim2, int dim2, int spacedim2>
  DoFAccessor (const InvalidAccessor<structdim2,dim2,spacedim2> &);

  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
  DoFAccessor (const DoFAccessor<dim2, DoFHandlerType2, level_dof_access2> &);

  DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access> &
  operator = (const DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access> &da) = delete;

  const DoFHandlerType<1,spacedim> &
  get_dof_handler () const;

  template <bool level_dof_access2>
  void copy_from (const DoFAccessor<0, DoFHandlerType<1,spacedim>, level_dof_access2> &a);

  void copy_from (const TriaAccessorBase<0, 1, spacedim> &da);

  TriaIterator<DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access > >
  child (const unsigned int c) const;

  typename ::internal::DoFHandlerImplementation::Iterators<DoFHandlerType<1,spacedim>, level_dof_access>::line_iterator
  line (const unsigned int i) const;

  typename ::internal::DoFHandlerImplementation::Iterators<DoFHandlerType<1,spacedim>, level_dof_access>::quad_iterator
  quad (const unsigned int i) const;

  void get_dof_indices (std::vector<types::global_dof_index> &dof_indices,
                        const unsigned int fe_index = AccessorData::default_fe_index) const;

  void get_mg_dof_indices (const int level,
                           std::vector<types::global_dof_index> &dof_indices,
                           const unsigned int fe_index = AccessorData::default_fe_index) const;

  types::global_dof_index vertex_dof_index (const unsigned int vertex,
                                            const unsigned int i,
                                            const unsigned int fe_index = AccessorData::default_fe_index) const;

  types::global_dof_index dof_index (const unsigned int i,
                                     const unsigned int fe_index = AccessorData::default_fe_index) const;

  unsigned int
  n_active_fe_indices () const;

  unsigned int
  nth_active_fe_index (const unsigned int n) const;

  bool
  fe_index_is_active (const unsigned int fe_index) const;

  const FiniteElement<DoFHandlerType<1,spacedim>::dimension,DoFHandlerType<1,spacedim>::space_dimension> &
  get_fe (const unsigned int fe_index) const;

  DeclExceptionMsg (ExcInvalidObject, "This accessor object has not been "
                    "associated with any DoFHandler object.");
  DeclException0 (ExcVectorNotEmpty);
  DeclException0 (ExcVectorDoesNotMatch);
  DeclException0 (ExcMatrixDoesNotMatch);
  DeclException0 (ExcNotActive);
  DeclException0 (ExcCantCompareIterators);

protected:

  DoFHandlerType<1,spacedim> *dof_handler;

  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
  bool operator == (const DoFAccessor<dim2,DoFHandlerType2,level_dof_access2> &) const;

  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
  bool operator != (const DoFAccessor<dim2,DoFHandlerType2,level_dof_access2> &) const;

  void set_dof_handler (DoFHandlerType<1,spacedim> *dh);

  void set_dof_index (const unsigned int i,
                      const types::global_dof_index index,
                      const unsigned int fe_index = AccessorData::default_fe_index) const;

  void set_vertex_dof_index (const unsigned int vertex,
                             const unsigned int i,
                             const types::global_dof_index index,
                             const unsigned int fe_index = AccessorData::default_fe_index) const;

  template <typename> friend class TriaRawIterator;


  template <int, int> friend class DoFHandler;
  template <int, int> friend class hp::DoFHandler;

  friend struct ::internal::DoFHandlerImplementation::Policy::Implementation;
  friend struct ::internal::DoFHandlerImplementation::Implementation;
  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
};



/* -------------------------------------------------------------------------- */


template <int structdim, int dim, int spacedim=dim>
class DoFInvalidAccessor : public InvalidAccessor<structdim,dim,spacedim>
{
public:
  typedef typename InvalidAccessor<structdim,dim,spacedim>::AccessorData AccessorData;

  DoFInvalidAccessor (const Triangulation<dim,spacedim> *parent     =  0,
                      const int                 level      = -1,
                      const int                 index      = -1,
                      const AccessorData       *local_data =  0);

  DoFInvalidAccessor (const DoFInvalidAccessor<structdim,dim,spacedim> &);

  template <typename OtherAccessor>
  DoFInvalidAccessor (const OtherAccessor &);

  void set_dof_index (const unsigned int i,
                      const types::global_dof_index index,
                      const unsigned int fe_index = DoFHandler<dim,spacedim>::default_fe_index) const;
};




/* -------------------------------------------------------------------------- */


template <typename DoFHandlerType, bool level_dof_access>
class DoFCellAccessor :  public DoFAccessor<DoFHandlerType::dimension,DoFHandlerType, level_dof_access>
{
public:
  static const unsigned int dim = DoFHandlerType::dimension;

  static const unsigned int spacedim = DoFHandlerType::space_dimension;


  typedef DoFHandlerType AccessorData;

  typedef DoFAccessor<DoFHandlerType::dimension,DoFHandlerType, level_dof_access> BaseClass;

  typedef DoFHandlerType Container;

  typedef
  TriaIterator<DoFAccessor<DoFHandlerType::dimension-1, DoFHandlerType, level_dof_access> >
  face_iterator;

  DoFCellAccessor (const Triangulation<DoFHandlerType::dimension,DoFHandlerType::space_dimension> *tria,
                   const int           level,
                   const int           index,
                   const AccessorData *local_data);

  template <int structdim2, int dim2, int spacedim2>
  DoFCellAccessor (const InvalidAccessor<structdim2,dim2,spacedim2> &);

  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
  explicit
  DoFCellAccessor (const DoFAccessor<dim2, DoFHandlerType2, level_dof_access2> &);

  DoFCellAccessor<DoFHandlerType, level_dof_access> &
  operator = (const DoFCellAccessor<DoFHandlerType, level_dof_access> &da) = delete;

  TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
  parent () const;

  TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
  neighbor (const unsigned int i) const;

  TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
  periodic_neighbor (const unsigned int i) const;

  TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
  neighbor_or_periodic_neighbor (const unsigned int i) const;

  TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
  child (const unsigned int i) const;

  face_iterator
  face (const unsigned int i) const;

  TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
  neighbor_child_on_subface (const unsigned int face_no,
                             const unsigned int subface_no) const;

  TriaIterator<DoFCellAccessor<DoFHandlerType, level_dof_access> >
  periodic_neighbor_child_on_subface (const unsigned int face_no,
                                      const unsigned int subface_no) const;

  template <class InputVector, typename number>
  void get_dof_values (const InputVector &values,
                       Vector<number>    &local_values) const;

  template <class InputVector, typename ForwardIterator>
  void get_dof_values (const InputVector &values,
                       ForwardIterator    local_values_begin,
                       ForwardIterator    local_values_end) const;

  template <class InputVector, typename ForwardIterator>
  void get_dof_values (const ConstraintMatrix &constraints,
                       const InputVector      &values,
                       ForwardIterator         local_values_begin,
                       ForwardIterator         local_values_end) const;

  template <class OutputVector, typename number>
  void set_dof_values (const Vector<number> &local_values,
                       OutputVector         &values) const;

  template <class InputVector, typename number>
  void get_interpolated_dof_values (const InputVector &values,
                                    Vector<number>    &interpolated_values,
                                    const unsigned int fe_index
                                    = DoFHandlerType::default_fe_index) const;

  template <class OutputVector, typename number>
  void set_dof_values_by_interpolation (const Vector<number> &local_values,
                                        OutputVector         &values,
                                        const unsigned int    fe_index
                                        = DoFHandlerType::default_fe_index) const;

  template <typename number, typename OutputVector>
  void
  distribute_local_to_global (const Vector<number> &local_source,
                              OutputVector         &global_destination) const;

  template <typename ForwardIterator, typename OutputVector>
  void
  distribute_local_to_global (ForwardIterator   local_source_begin,
                              ForwardIterator   local_source_end,
                              OutputVector     &global_destination) const;

  template <typename ForwardIterator, typename OutputVector>
  void
  distribute_local_to_global (const ConstraintMatrix &constraints,
                              ForwardIterator         local_source_begin,
                              ForwardIterator         local_source_end,
                              OutputVector           &global_destination) const;

  template <typename number, typename OutputMatrix>
  void
  distribute_local_to_global (const FullMatrix<number> &local_source,
                              OutputMatrix             &global_destination) const;

  template <typename number, typename OutputMatrix, typename OutputVector>
  void
  distribute_local_to_global (const FullMatrix<number> &local_matrix,
                              const Vector<number>     &local_vector,
                              OutputMatrix             &global_matrix,
                              OutputVector             &global_vector) const;

  void get_active_or_mg_dof_indices (std::vector<types::global_dof_index> &dof_indices) const;

  void get_dof_indices (std::vector<types::global_dof_index> &dof_indices) const;

  void get_mg_dof_indices (std::vector<types::global_dof_index> &dof_indices) const;

  const FiniteElement<DoFHandlerType::dimension,DoFHandlerType::space_dimension> &
  get_fe () const;

  unsigned int active_fe_index () const;

  void set_active_fe_index (const unsigned int i) const;
  void set_dof_indices (const std::vector<types::global_dof_index> &dof_indices);

  void set_mg_dof_indices (const std::vector<types::global_dof_index> &dof_indices);

  void update_cell_dof_indices_cache () const;

private:

  template <int dim, int spacedim> friend class DoFHandler;
  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
};


template <int sd, typename DoFHandlerType, bool level_dof_access>
inline
bool
DoFAccessor<sd, DoFHandlerType, level_dof_access>::is_level_cell()
{
  return level_dof_access;
}



template <int structdim, int dim, int spacedim>
template <typename OtherAccessor>
DoFInvalidAccessor<structdim, dim, spacedim>::
DoFInvalidAccessor (const OtherAccessor &)
{
  Assert (false,
          ExcMessage ("You are attempting an illegal conversion between "
                      "iterator/accessor types. The constructor you call "
                      "only exists to make certain template constructs "
                      "easier to write as dimension independent code but "
                      "the conversion is not valid in the current context."));
}



DEAL_II_NAMESPACE_CLOSE

// include more templates
#include "dof_accessor.templates.h"


#endif