Reference documentation for deal.II version 8.4.1
dof_accessor.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2016 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii__dof_accessor_h
17 #define dealii__dof_accessor_h
18 
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/grid/tria_accessor.h>
22 #include <deal.II/dofs/dof_handler.h>
23 #include <deal.II/hp/dof_handler.h>
24 
25 #include <vector>
26 
27 DEAL_II_NAMESPACE_OPEN
28 
29 template <typename number> class FullMatrix;
30 template <typename number> class SparseMatrix;
31 template <typename number> class Vector;
32 class ConstraintMatrix;
33 
34 template <typename Accessor> class TriaRawIterator;
35 
36 template <int, int> class FiniteElement;
37 
38 
39 namespace internal
40 {
41  namespace DoFCellAccessor
42  {
43  struct Implementation;
44  }
45 
46  namespace DoFHandler
47  {
48  struct Implementation;
49  namespace Policy
50  {
51  struct Implementation;
52  }
53  }
54 
55  namespace hp
56  {
57  namespace DoFHandler
58  {
59  struct Implementation;
60  }
61  }
62 }
63 
64 // note: the file dof_accessor.templates.h is included at the end of
65 // this file. this includes a lot of templates and thus makes
66 // compilation slower, but at the same time allows for more aggressive
67 // inlining and thus faster code.
68 
69 
70 namespace internal
71 {
72  namespace DoFAccessor
73  {
91  template <int structdim, int dim, int spacedim>
92  struct Inheritance
93  {
98  typedef ::TriaAccessor<structdim,dim,spacedim> BaseClass;
99  };
100 
101 
107  template <int dim, int spacedim>
108  struct Inheritance<dim,dim,spacedim>
109  {
114  typedef ::CellAccessor<dim,spacedim> BaseClass;
115  };
116  }
117 }
118 
119 
120 /* -------------------------------------------------------------------------- */
121 
122 
123 
184 template <int structdim, typename DoFHandlerType, bool level_dof_access>
185 class DoFAccessor : public ::internal::DoFAccessor::Inheritance<structdim, DoFHandlerType::dimension, DoFHandlerType::space_dimension>::BaseClass
186 {
187 public:
188 
193  static const unsigned int dimension=DoFHandlerType::dimension;
194 
199  static const unsigned int space_dimension=DoFHandlerType::space_dimension;
200 
205  typedef
206  typename ::internal::DoFAccessor::Inheritance<structdim, dimension, space_dimension>::BaseClass
208 
212  typedef DoFHandlerType AccessorData;
213 
224  DoFAccessor ();
225 
230  const int level,
231  const int index,
232  const DoFHandlerType *local_data);
233 
246  template <int structdim2, int dim2, int spacedim2>
248 
253  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
255 
259  template <bool level_dof_access2>
268  const DoFHandlerType &
269  get_dof_handler () const;
270 
274  template <bool level_dof_access2>
276 
282 
287  static bool is_level_cell();
288 
300  child (const unsigned int c) const;
301 
307  typename ::internal::DoFHandler::Iterators<DoFHandlerType, level_dof_access>::line_iterator
308  line (const unsigned int i) const;
309 
315  typename ::internal::DoFHandler::Iterators<DoFHandlerType, level_dof_access>::quad_iterator
316  quad (const unsigned int i) const;
317 
363  void get_dof_indices (std::vector<types::global_dof_index> &dof_indices,
364  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
365 
372  void get_mg_dof_indices (const int level,
373  std::vector<types::global_dof_index> &dof_indices,
374  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
375 
379  void set_mg_dof_indices (const int level,
380  const std::vector<types::global_dof_index> &dof_indices,
381  const unsigned int fe_index = DoFHandlerType::default_fe_index);
382 
401  (const unsigned int vertex,
402  const unsigned int i,
403  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
404 
411  (const int level,
412  const unsigned int vertex,
413  const unsigned int i,
414  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
415 
444  (const unsigned int i,
445  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
446 
450  types::global_dof_index mg_dof_index (const int level, const unsigned int i) const;
451 
473  unsigned int
474  n_active_fe_indices () const;
475 
483  unsigned int
484  nth_active_fe_index (const unsigned int n) const;
485 
494  bool
495  fe_index_is_active (const unsigned int fe_index) const;
496 
503  get_fe (const unsigned int fe_index) const;
504 
514  DeclException0 (ExcInvalidObject);
520  DeclException0 (ExcVectorNotEmpty);
526  DeclException0 (ExcVectorDoesNotMatch);
532  DeclException0 (ExcMatrixDoesNotMatch);
540  DeclException0 (ExcNotActive);
546  DeclException0 (ExcCantCompareIterators);
547 
548 protected:
549 
553  DoFHandlerType *dof_handler;
554 public:
568  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
570 
574  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
576 protected:
580  void set_dof_handler (DoFHandlerType *dh);
581 
599  void set_dof_index
600  (const unsigned int i,
602  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
603 
604  void set_mg_dof_index (const int level, const unsigned int i, const types::global_dof_index index) const;
605 
624  (const unsigned int vertex,
625  const unsigned int i,
627  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
628 
629  void set_mg_vertex_dof_index
630  (const int level,
631  const unsigned int vertex,
632  const unsigned int i,
634  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
635 
640  template <typename> friend class TriaRawIterator;
641  template <int, class, bool> friend class DoFAccessor;
642 
643 private:
656 
661  template <int dim, int spacedim> friend class DoFHandler;
662  template <int dim, int spacedim> friend class hp::DoFHandler;
663 
664  friend struct ::internal::DoFHandler::Policy::Implementation;
665  friend struct ::internal::DoFHandler::Implementation;
666  friend struct ::internal::hp::DoFHandler::Implementation;
667  friend struct ::internal::DoFCellAccessor::Implementation;
668  friend struct ::internal::DoFAccessor::Implementation;
669 };
670 
671 
672 
682 template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access>
683 class DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access> : public TriaAccessor<0,1,spacedim>
684 {
685 public:
686 
691  static const unsigned int dimension=1;
692 
697  static const unsigned int space_dimension=spacedim;
698 
704 
708  typedef DoFHandlerType<1,spacedim> AccessorData;
709 
720  DoFAccessor ();
721 
739  const typename TriaAccessor<0,1,spacedim>::VertexKind vertex_kind,
740  const unsigned int vertex_index,
741  const DoFHandlerType<1,spacedim> *dof_handler);
742 
749  const int = 0,
750  const int = 0,
751  const DoFHandlerType<1,spacedim> *dof_handler = 0);
752 
765  template <int structdim2, int dim2, int spacedim2>
767 
772  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
774 
782  const DoFHandlerType<1,spacedim> &
783  get_dof_handler () const;
784 
788  DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access> &
789  operator = (const DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access> &da);
790 
794  template <bool level_dof_access2>
795  void copy_from (const DoFAccessor<0, DoFHandlerType<1,spacedim>, level_dof_access2> &a);
796 
802 
816  child (const unsigned int c) const;
817 
823  typename ::internal::DoFHandler::Iterators<DoFHandlerType<1,spacedim>, level_dof_access>::line_iterator
824  line (const unsigned int i) const;
825 
831  typename ::internal::DoFHandler::Iterators<DoFHandlerType<1,spacedim>, level_dof_access>::quad_iterator
832  quad (const unsigned int i) const;
833 
878  void get_dof_indices (std::vector<types::global_dof_index> &dof_indices,
879  const unsigned int fe_index = AccessorData::default_fe_index) const;
880 
898  types::global_dof_index vertex_dof_index (const unsigned int vertex,
899  const unsigned int i,
900  const unsigned int fe_index = AccessorData::default_fe_index) const;
901 
929  types::global_dof_index dof_index (const unsigned int i,
930  const unsigned int fe_index = AccessorData::default_fe_index) const;
931 
953  unsigned int
954  n_active_fe_indices () const;
955 
963  unsigned int
964  nth_active_fe_index (const unsigned int n) const;
965 
974  bool
975  fe_index_is_active (const unsigned int fe_index) const;
976 
982  const FiniteElement<DoFHandlerType<1,spacedim>::dimension,DoFHandlerType<1,spacedim>::space_dimension> &
983  get_fe (const unsigned int fe_index) const;
984 
994  DeclException0 (ExcInvalidObject);
1000  DeclException0 (ExcVectorNotEmpty);
1006  DeclException0 (ExcVectorDoesNotMatch);
1012  DeclException0 (ExcMatrixDoesNotMatch);
1020  DeclException0 (ExcNotActive);
1026  DeclException0 (ExcCantCompareIterators);
1027 
1028 protected:
1029 
1033  DoFHandlerType<1,spacedim> *dof_handler;
1034 
1038  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
1040 
1044  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
1046 
1050  void set_dof_handler (DoFHandlerType<1,spacedim> *dh);
1051 
1069  void set_dof_index (const unsigned int i,
1071  const unsigned int fe_index = AccessorData::default_fe_index) const;
1072 
1090  void set_vertex_dof_index (const unsigned int vertex,
1091  const unsigned int i,
1093  const unsigned int fe_index = AccessorData::default_fe_index) const;
1094 
1099  template <typename> friend class TriaRawIterator;
1100 
1101 
1106  template <int, int> friend class DoFHandler;
1107  template <int, int> friend class hp::DoFHandler;
1108 
1109  friend struct ::internal::DoFHandler::Policy::Implementation;
1110  friend struct ::internal::DoFHandler::Implementation;
1111  friend struct ::internal::hp::DoFHandler::Implementation;
1112  friend struct ::internal::DoFCellAccessor::Implementation;
1113 };
1114 
1115 
1116 /* -------------------------------------------------------------------------- */
1117 
1118 
1131 template <typename DoFHandlerType, bool level_dof_access>
1132 class DoFCellAccessor : public DoFAccessor<DoFHandlerType::dimension,DoFHandlerType, level_dof_access>
1133 {
1134 public:
1138  static const unsigned int dim = DoFHandlerType::dimension;
1139 
1143  static const unsigned int spacedim = DoFHandlerType::space_dimension;
1144 
1145 
1149  typedef DoFHandlerType AccessorData;
1150 
1156 
1160  typedef DoFHandlerType Container;
1161 
1166  typedef
1167  TriaIterator<DoFAccessor<DoFHandlerType::dimension-1, DoFHandlerType, level_dof_access> >
1169 
1181  const int level,
1182  const int index,
1183  const AccessorData *local_data);
1184 
1197  template <int structdim2, int dim2, int spacedim2>
1199 
1204  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
1205  explicit
1207 
1222  parent () const;
1223 
1237  neighbor (const unsigned int) const;
1238 
1245  child (const unsigned int) const;
1246 
1254  face (const unsigned int i) const;
1255 
1263  neighbor_child_on_subface (const unsigned int face_no,
1264  const unsigned int subface_no) const;
1265 
1291  template <class InputVector, typename number>
1292  void get_dof_values (const InputVector &values,
1293  Vector<number> &local_values) const;
1294 
1309  template <class InputVector, typename ForwardIterator>
1310  void get_dof_values (const InputVector &values,
1311  ForwardIterator local_values_begin,
1312  ForwardIterator local_values_end) const;
1313 
1330  template <class InputVector, typename ForwardIterator>
1331  void get_dof_values (const ConstraintMatrix &constraints,
1332  const InputVector &values,
1333  ForwardIterator local_values_begin,
1334  ForwardIterator local_values_end) const;
1335 
1357  template <class OutputVector, typename number>
1358  void set_dof_values (const Vector<number> &local_values,
1359  OutputVector &values) const;
1360 
1392  template <class InputVector, typename number>
1393  void get_interpolated_dof_values (const InputVector &values,
1394  Vector<number> &interpolated_values,
1395  const unsigned int fe_index
1396  = DoFHandlerType::default_fe_index) const;
1397 
1450  template <class OutputVector, typename number>
1451  void set_dof_values_by_interpolation (const Vector<number> &local_values,
1452  OutputVector &values,
1453  const unsigned int fe_index
1454  = DoFHandlerType::default_fe_index) const;
1455 
1464  template <typename number, typename OutputVector>
1465  void
1466  distribute_local_to_global (const Vector<number> &local_source,
1467  OutputVector &global_destination) const;
1468 
1477  template <typename ForwardIterator, typename OutputVector>
1478  void
1479  distribute_local_to_global (ForwardIterator local_source_begin,
1480  ForwardIterator local_source_end,
1481  OutputVector &global_destination) const;
1482 
1493  template <typename ForwardIterator, typename OutputVector>
1494  void
1495  distribute_local_to_global (const ConstraintMatrix &constraints,
1496  ForwardIterator local_source_begin,
1497  ForwardIterator local_source_end,
1498  OutputVector &global_destination) const;
1499 
1506  template <typename number, typename OutputMatrix>
1507  void
1508  distribute_local_to_global (const FullMatrix<number> &local_source,
1509  OutputMatrix &global_destination) const;
1510 
1515  template <typename number, typename OutputMatrix, typename OutputVector>
1516  void
1517  distribute_local_to_global (const FullMatrix<number> &local_matrix,
1518  const Vector<number> &local_vector,
1519  OutputMatrix &global_matrix,
1520  OutputVector &global_vector) const;
1521 
1546  void get_active_or_mg_dof_indices (std::vector<types::global_dof_index> &dof_indices) const;
1547 
1583  void get_dof_indices (std::vector<types::global_dof_index> &dof_indices) const;
1584 
1592  void get_mg_dof_indices (std::vector<types::global_dof_index> &dof_indices) const;
1593 
1619  get_fe () const;
1620 
1633  unsigned int active_fe_index () const;
1634 
1648  void set_active_fe_index (const unsigned int i);
1657  void set_dof_indices (const std::vector<types::global_dof_index> &dof_indices);
1658 
1662  void set_mg_dof_indices (const std::vector<types::global_dof_index> &dof_indices);
1663 
1667  void update_cell_dof_indices_cache () const;
1668 
1669 private:
1682 
1687  template <int dim, int spacedim> friend class DoFHandler;
1688  friend struct ::internal::DoFCellAccessor::Implementation;
1689 };
1690 
1691 
1692 template <int sd, typename DoFHandlerType, bool level_dof_access>
1693 inline
1694 bool
1696 {
1697  return level_dof_access;
1698 }
1699 
1700 
1701 
1702 DEAL_II_NAMESPACE_CLOSE
1703 
1704 // include more templates
1705 #include "dof_accessor.templates.h"
1706 
1707 
1708 #endif
face_iterator face(const unsigned int i) const
DoFHandlerType * dof_handler
Definition: dof_accessor.h:553
unsigned int n_active_fe_indices() const
types::global_dof_index mg_dof_index(const int level, const unsigned int i) const
::TriaAccessor< structdim, dim, spacedim > BaseClass
Definition: dof_accessor.h:98
unsigned int nth_active_fe_index(const unsigned int n) const
static bool is_level_cell()
DeclException0(ExcInvalidObject)
static const unsigned int spacedim
TriaIterator< DoFAccessor< structdim, DoFHandlerType, level_dof_access > > child(const unsigned int c) const
void get_interpolated_dof_values(const InputVector &values, Vector< number > &interpolated_values, const unsigned int fe_index=DoFHandlerType::default_fe_index) const
TriaIterator< DoFCellAccessor< DoFHandlerType, level_dof_access > > neighbor(const unsigned int) 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
DoFHandlerType Container
void copy_from(const DoFAccessor< structdim, DoFHandlerType, level_dof_access2 > &a)
void set_dof_values(const Vector< number > &local_values, OutputVector &values) const
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
TriaIterator< DoFCellAccessor< DoFHandlerType, level_dof_access > > parent() const
TriaIterator< DoFAccessor< DoFHandlerType::dimension-1, DoFHandlerType, level_dof_access > > face_iterator
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)
TriaIterator< DoFCellAccessor< DoFHandlerType, level_dof_access > > neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
Definition: dof_accessor.cc:77
void set_dof_handler(DoFHandlerType *dh)
TriaIterator< DoFCellAccessor< DoFHandlerType, level_dof_access > > child(const unsigned int) const
const FiniteElement< DoFHandlerType::dimension, DoFHandlerType::space_dimension > & get_fe() const
DoFHandlerType AccessorData
Definition: dof_accessor.h:212
unsigned int active_fe_index() const
void set_active_fe_index(const unsigned int i)
typename::internal::DoFAccessor::Inheritance< structdim, dimension, space_dimension >::BaseClass BaseClass
Definition: dof_accessor.h:207
unsigned int global_dof_index
Definition: types.h:88
void get_dof_values(const InputVector &values, Vector< number > &local_values) const
types::global_dof_index dof_index(const unsigned int i, const unsigned int fe_index=DoFHandlerType::default_fe_index) const
void set_mg_dof_indices(const std::vector< types::global_dof_index > &dof_indices)
static const unsigned int dimension
Definition: dof_accessor.h:193
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
unsigned int vertex_index(const unsigned int i) const
const Triangulation< dim, spacedim > * tria
void get_mg_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
typename::internal::DoFHandler::Iterators< DoFHandlerType, level_dof_access >::quad_iterator quad(const unsigned int i) const
void set_dof_indices(const std::vector< types::global_dof_index > &dof_indices)
Definition: dof_accessor.cc:61
DoFCellAccessor< DoFHandlerType, level_dof_access > & operator=(const DoFCellAccessor< DoFHandlerType, level_dof_access > &da)
Definition: hp.h:102
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
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 get_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
typename::internal::DoFHandler::Iterators< DoFHandlerType, level_dof_access >::line_iterator line(const unsigned int i) const
static const unsigned int dim
const FiniteElement< DoFHandlerType::dimension, DoFHandlerType::space_dimension > & get_fe(const unsigned int fe_index) const
DoFHandlerType AccessorData
void get_active_or_mg_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
void update_cell_dof_indices_cache() const
Definition: dof_accessor.cc:46
Point< spacedim > & vertex(const unsigned int i) const
bool operator!=(const DoFAccessor< dim2, DoFHandlerType2, level_dof_access2 > &) const
int index() const
DoFCellAccessor(const Triangulation< DoFHandlerType::dimension, DoFHandlerType::space_dimension > *tria, const int level, const int index, const AccessorData *local_data)
static const unsigned int space_dimension
Definition: dof_accessor.h:199
bool operator==(const DoFAccessor< dim2, DoFHandlerType2, level_dof_access2 > &) const
bool fe_index_is_active(const unsigned int fe_index) const
DoFAccessor< structdim, DoFHandlerType, level_dof_access > & operator=(const DoFAccessor< structdim, DoFHandlerType, level_dof_access > &da)
void distribute_local_to_global(const Vector< number > &local_source, OutputVector &global_destination) const
void set_dof_values_by_interpolation(const Vector< number > &local_values, OutputVector &values, const unsigned int fe_index=DoFHandlerType::default_fe_index) const
const DoFHandlerType & get_dof_handler() const
int level() const
DoFAccessor< DoFHandlerType::dimension, DoFHandlerType, level_dof_access > BaseClass
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DoFHandlerType::default_fe_index) const