Reference documentation for deal.II version Git 4abc4a1666 2020-07-04 19:58:34 +0200
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
dof_accessor.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2019 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.md at
12 // the top level directory of deal.II.
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 
24 
26 
27 #include <deal.II/hp/dof_handler.h>
28 
29 #include <vector>
30 
32 
33 // Forward declarations
34 #ifndef DOXYGEN
35 template <typename number>
36 class FullMatrix;
37 template <typename number>
38 class Vector;
39 template <typename number>
40 class AffineConstraints;
41 
42 template <typename Accessor>
43 class TriaRawIterator;
44 
45 template <int, int>
46 class FiniteElement;
47 
48 namespace internal
49 {
50  namespace DoFCellAccessorImplementation
51  {
52  struct Implementation;
53  }
54 
55  namespace DoFHandlerImplementation
56  {
57  struct Implementation;
58  namespace Policy
59  {
60  struct Implementation;
61  }
62  } // namespace DoFHandlerImplementation
63 
64  namespace hp
65  {
66  namespace DoFHandlerImplementation
67  {
68  struct Implementation;
69  }
70  } // namespace hp
71 } // namespace internal
72 #endif
73 
74 // note: the file dof_accessor.templates.h is included at the end of
75 // this file. this includes a lot of templates and thus makes
76 // compilation slower, but at the same time allows for more aggressive
77 // inlining and thus faster code.
78 
79 
80 namespace internal
81 {
82  namespace DoFAccessorImplementation
83  {
100  template <int structdim, int dim, int spacedim>
101  struct Inheritance
102  {
108  };
109 
110 
116  template <int dim, int spacedim>
117  struct Inheritance<dim, dim, spacedim>
118  {
124  };
125 
126  struct Implementation;
127  } // namespace DoFAccessorImplementation
128 } // namespace internal
129 
130 
131 /* -------------------------------------------------------------------------- */
132 
133 
134 
206 template <int structdim, int dim, int spacedim, bool level_dof_access>
207 class DoFAccessor : public ::internal::DoFAccessorImplementation::
208  Inheritance<structdim, dim, spacedim>::BaseClass
209 {
210 public:
215  static const unsigned int dimension = dim;
216 
221  static const unsigned int space_dimension = spacedim;
222 
227  using BaseClass = typename ::internal::DoFAccessorImplementation::
228  Inheritance<structdim, dimension, space_dimension>::BaseClass;
229 
234 
245  DoFAccessor();
246 
263  const int level,
264  const int index,
265  const DoFHandler<dim, spacedim> * dof_handler);
266 
279  template <int structdim2, int dim2, int spacedim2>
281 
286  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
287  DoFAccessor(
289 
293  template <bool level_dof_access2>
295 
307  delete;
308 
317  get_dof_handler() const;
318 
322  template <bool level_dof_access2>
323  void
325 
330  void
331  copy_from(const TriaAccessorBase<structdim, dim, spacedim> &da);
332 
337  static bool
338  is_level_cell();
339 
351  child(const unsigned int c) const;
352 
358  typename ::internal::DoFHandlerImplementation::
359  Iterators<dim, spacedim, level_dof_access>::line_iterator
360  line(const unsigned int i) const;
361 
367  typename ::internal::DoFHandlerImplementation::
368  Iterators<dim, spacedim, level_dof_access>::quad_iterator
369  quad(const unsigned int i) const;
370 
416  void
417  get_dof_indices(std::vector<types::global_dof_index> &dof_indices,
418  const unsigned int fe_index =
420 
427  void
428  get_mg_dof_indices(const int level,
429  std::vector<types::global_dof_index> &dof_indices,
430  const unsigned int fe_index =
432 
436  void
437  set_mg_dof_indices(
438  const int level,
439  const std::vector<types::global_dof_index> &dof_indices,
440  const unsigned int fe_index = DoFHandler<dim, spacedim>::invalid_fe_index);
441 
460  vertex_dof_index(const unsigned int vertex,
461  const unsigned int i,
462  const unsigned int fe_index =
464 
471  mg_vertex_dof_index(const int level,
472  const unsigned int vertex,
473  const unsigned int i,
474  const unsigned int fe_index =
476 
505  dof_index(const unsigned int i,
506  const unsigned int fe_index =
508 
513  mg_dof_index(const int level, const unsigned int i) const;
514 
536  unsigned int
537  n_active_fe_indices() const;
538 
546  unsigned int
547  nth_active_fe_index(const unsigned int n) const;
548 
555  std::set<unsigned int>
556  get_active_fe_indices() const;
557 
567  bool
568  fe_index_is_active(const unsigned int fe_index) const;
569 
576  get_fe(const unsigned int fe_index) const;
577 
587  DeclExceptionMsg(ExcInvalidObject,
588  "This accessor object has not been "
589  "associated with any DoFHandler object.");
595  DeclException0(ExcVectorNotEmpty);
601  DeclException0(ExcVectorDoesNotMatch);
607  DeclException0(ExcMatrixDoesNotMatch);
615  DeclException0(ExcNotActive);
622 
623 protected:
628 
629 public:
643  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
644  bool
645  operator==(
647 
651  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
652  bool
653  operator!=(
655 
656 protected:
660  void
661  set_dof_handler(DoFHandler<dim, spacedim> *dh);
662 
681  void
682  set_dof_index(const unsigned int i,
683  const types::global_dof_index index,
684  const unsigned int fe_index =
686 
687  void
688  set_mg_dof_index(const int level,
689  const unsigned int i,
690  const types::global_dof_index index) const;
691 
710  void
711  set_vertex_dof_index(const unsigned int vertex,
712  const unsigned int i,
713  const types::global_dof_index index,
714  const unsigned int fe_index =
716 
717  void
718  set_mg_vertex_dof_index(const int level,
719  const unsigned int vertex,
720  const unsigned int i,
721  const types::global_dof_index index,
722  const unsigned int fe_index =
724 
725  // Iterator classes need to be friends because they need to access
726  // operator== and operator!=.
727  template <typename>
728  friend class TriaRawIterator;
729  template <int, int, int, bool>
730  friend class DoFAccessor;
731 
732 private:
733  // Make the DoFHandler class a friend so that it can call the set_xxx()
734  // functions.
735  template <int, int>
736  friend class DoFHandler;
737  template <int, int>
738  friend class hp::DoFHandler;
739 
740  friend struct ::internal::DoFHandlerImplementation::Policy::
741  Implementation;
742  friend struct ::internal::DoFHandlerImplementation::Implementation;
743  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
744  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
745  friend struct ::internal::DoFAccessorImplementation::Implementation;
746 };
747 
748 
749 
757 template <int spacedim, bool level_dof_access>
758 class DoFAccessor<0, 1, spacedim, level_dof_access>
759  : public TriaAccessor<0, 1, spacedim>
760 {
761 public:
766  static const unsigned int dimension = 1;
767 
772  static const unsigned int space_dimension = spacedim;
773 
779 
784 
795  DoFAccessor();
796 
813  DoFAccessor(
814  const Triangulation<1, spacedim> * tria,
815  const typename TriaAccessor<0, 1, spacedim>::VertexKind vertex_kind,
816  const unsigned int vertex_index,
817  const DoFHandler<1, spacedim> * dof_handler);
818 
825  const int = 0,
826  const int = 0,
827  const DoFHandler<1, spacedim> *dof_handler = 0);
828 
841  template <int structdim2, int dim2, int spacedim2>
843 
848  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
849  DoFAccessor(
851 
862  operator=(const DoFAccessor<0, 1, spacedim, level_dof_access> &da) = delete;
863 
872  get_dof_handler() const;
873 
877  template <bool level_dof_access2>
878  void
880 
885  void
886  copy_from(const TriaAccessorBase<0, 1, spacedim> &da);
887 
901  child(const unsigned int c) const;
902 
909  typename ::internal::DoFHandlerImplementation::
910  Iterators<1, spacedim, level_dof_access>::line_iterator
911  line(const unsigned int i) const;
912 
919  typename ::internal::DoFHandlerImplementation::
920  Iterators<1, spacedim, level_dof_access>::quad_iterator
921  quad(const unsigned int i) const;
922 
966  void
967  get_dof_indices(
968  std::vector<types::global_dof_index> &dof_indices,
969  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
970 
977  void
978  get_mg_dof_indices(
979  const int level,
980  std::vector<types::global_dof_index> &dof_indices,
981  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
982 
1002  vertex_dof_index(
1003  const unsigned int vertex,
1004  const unsigned int i,
1005  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1006 
1025  dof_index(const unsigned int i,
1026  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1027 
1046  unsigned int
1047  n_active_fe_indices() const;
1048 
1056  unsigned int
1057  nth_active_fe_index(const unsigned int n) const;
1058 
1067  bool
1068  fe_index_is_active(const unsigned int fe_index) const;
1069 
1076  get_fe(const unsigned int fe_index) const;
1077 
1087  DeclExceptionMsg(ExcInvalidObject,
1088  "This accessor object has not been "
1089  "associated with any DoFHandler object.");
1095  DeclException0(ExcVectorNotEmpty);
1101  DeclException0(ExcVectorDoesNotMatch);
1107  DeclException0(ExcMatrixDoesNotMatch);
1115  DeclException0(ExcNotActive);
1122 
1123 protected:
1128 
1132  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1133  bool
1134  operator==(
1136 
1140  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1141  bool
1142  operator!=(
1144 
1148  void set_dof_handler(DoFHandler<1, spacedim> *dh);
1149 
1168  void
1169  set_dof_index(
1170  const unsigned int i,
1171  const types::global_dof_index index,
1172  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1173 
1192  void
1193  set_vertex_dof_index(
1194  const unsigned int vertex,
1195  const unsigned int i,
1196  const types::global_dof_index index,
1197  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1198 
1199  // Iterator classes need to be friends because they need to access
1200  // operator== and operator!=.
1201  template <typename>
1202  friend class TriaRawIterator;
1203 
1204 
1205  // Make the DoFHandler class a friend so that it can call the set_xxx()
1206  // functions.
1207  template <int, int>
1208  friend class DoFHandler;
1209  template <int, int>
1210  friend class hp::DoFHandler;
1211 
1212  friend struct ::internal::DoFHandlerImplementation::Policy::
1213  Implementation;
1214  friend struct ::internal::DoFHandlerImplementation::Implementation;
1215  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1216  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1217 };
1218 
1219 
1220 
1221 /* -------------------------------------------------------------------------- */
1222 
1223 
1244 template <int structdim, int dim, int spacedim = dim>
1245 class DoFInvalidAccessor : public InvalidAccessor<structdim, dim, spacedim>
1246 {
1247 public:
1251  using AccessorData =
1253 
1262  const int level = -1,
1263  const int index = -1,
1264  const AccessorData * local_data = 0);
1265 
1274 
1279  template <typename OtherAccessor>
1280  DoFInvalidAccessor(const OtherAccessor &);
1281 
1288  dof_index(const unsigned int i,
1289  const unsigned int fe_index =
1291 
1297  void
1298  set_dof_index(const unsigned int i,
1299  const types::global_dof_index index,
1300  const unsigned int fe_index =
1302 };
1303 
1304 
1305 
1306 /* -------------------------------------------------------------------------- */
1307 
1308 
1320 template <int dimension_, int space_dimension_, bool level_dof_access>
1321 class DoFCellAccessor : public DoFAccessor<dimension_,
1322  dimension_,
1323  space_dimension_,
1324  level_dof_access>
1325 {
1326 public:
1330  static const unsigned int dim = dimension_;
1331 
1335  static const unsigned int spacedim = space_dimension_;
1336 
1337 
1342 
1347  using BaseClass =
1349 
1354 
1359  using face_iterator = TriaIterator<DoFAccessor<dimension_ - 1,
1360  dimension_,
1361  space_dimension_,
1362  level_dof_access>>;
1363 
1375  const int level,
1376  const int index,
1377  const AccessorData *local_data);
1378 
1391  template <int structdim2, int dim2, int spacedim2>
1393 
1398  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1399  explicit DoFCellAccessor(
1401 
1412  operator=(
1414  delete;
1415 
1430  parent() const;
1431 
1445  neighbor(const unsigned int i) const;
1446 
1453  periodic_neighbor(const unsigned int i) const;
1454 
1461  neighbor_or_periodic_neighbor(const unsigned int i) const;
1462 
1469  child(const unsigned int i) const;
1470 
1478  face(const unsigned int i) const;
1479 
1483  inline std::array<face_iterator, GeometryInfo<dimension_>::faces_per_cell>
1484  face_iterators() const;
1485 
1493  neighbor_child_on_subface(const unsigned int face_no,
1494  const unsigned int subface_no) const;
1495 
1503  periodic_neighbor_child_on_subface(const unsigned int face_no,
1504  const unsigned int subface_no) const;
1505 
1531  template <class InputVector, typename number>
1532  void
1533  get_dof_values(const InputVector &values, Vector<number> &local_values) const;
1534 
1549  template <class InputVector, typename ForwardIterator>
1550  void
1551  get_dof_values(const InputVector &values,
1552  ForwardIterator local_values_begin,
1553  ForwardIterator local_values_end) const;
1554 
1572  template <class InputVector, typename ForwardIterator>
1573  void
1574  get_dof_values(
1576  const InputVector & values,
1577  ForwardIterator local_values_begin,
1578  ForwardIterator local_values_end) const;
1579 
1601  template <class OutputVector, typename number>
1602  void
1603  set_dof_values(const Vector<number> &local_values,
1604  OutputVector & values) const;
1605 
1637  template <class InputVector, typename number>
1638  void
1639  get_interpolated_dof_values(
1640  const InputVector &values,
1641  Vector<number> & interpolated_values,
1642  const unsigned int fe_index =
1644 
1697  template <class OutputVector, typename number>
1698  void
1699  set_dof_values_by_interpolation(
1700  const Vector<number> &local_values,
1701  OutputVector & values,
1702  const unsigned int fe_index =
1704 
1713  template <typename number, typename OutputVector>
1714  void
1715  distribute_local_to_global(const Vector<number> &local_source,
1716  OutputVector & global_destination) const;
1717 
1726  template <typename ForwardIterator, typename OutputVector>
1727  void
1728  distribute_local_to_global(ForwardIterator local_source_begin,
1729  ForwardIterator local_source_end,
1730  OutputVector & global_destination) const;
1731 
1742  template <typename ForwardIterator, typename OutputVector>
1743  void
1744  distribute_local_to_global(
1746  ForwardIterator local_source_begin,
1747  ForwardIterator local_source_end,
1748  OutputVector & global_destination) const;
1749 
1756  template <typename number, typename OutputMatrix>
1757  void
1758  distribute_local_to_global(const FullMatrix<number> &local_source,
1759  OutputMatrix &global_destination) const;
1760 
1765  template <typename number, typename OutputMatrix, typename OutputVector>
1766  void
1767  distribute_local_to_global(const FullMatrix<number> &local_matrix,
1768  const Vector<number> & local_vector,
1769  OutputMatrix & global_matrix,
1770  OutputVector & global_vector) const;
1771 
1796  void
1797  get_active_or_mg_dof_indices(
1798  std::vector<types::global_dof_index> &dof_indices) const;
1799 
1835  void
1836  get_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1837 
1842  void
1843  get_mg_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1844 
1870  get_fe() const;
1871 
1897  unsigned int
1898  active_fe_index() const;
1899 
1926  void
1927  set_active_fe_index(const unsigned int i) const;
1936  void
1937  set_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
1938 
1942  void
1943  set_mg_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
1944 
1948  void
1949  update_cell_dof_indices_cache() const;
1950 
1977  get_future_fe() const;
1978 
1997  unsigned int
1998  future_fe_index() const;
1999 
2008  void
2009  set_future_fe_index(const unsigned int i) const;
2010 
2017  bool
2018  future_fe_index_set() const;
2019 
2028  void
2029  clear_future_fe_index() const;
2034 private:
2035  // Make the DoFHandler class a friend so that it can call the
2036  // update_cell_dof_indices_cache() function
2037  template <int, int>
2038  friend class DoFHandler;
2039  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
2040 };
2041 
2042 
2043 template <int structdim, int dim, int spacedim, bool level_dof_access>
2044 inline bool
2046 {
2047  return level_dof_access;
2048 }
2049 
2050 
2051 
2052 template <int structdim, int dim, int spacedim>
2053 template <typename OtherAccessor>
2055  const OtherAccessor &)
2056 {
2057  Assert(false,
2058  ExcMessage("You are attempting an illegal conversion between "
2059  "iterator/accessor types. The constructor you call "
2060  "only exists to make certain template constructs "
2061  "easier to write as dimension independent code but "
2062  "the conversion is not valid in the current context."));
2063 }
2064 
2065 
2066 
2068 
2069 // include more templates
2070 #include "dof_accessor.templates.h"
2071 
2072 
2073 #endif
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
typename InvalidAccessor< structdim, dim, spacedim >::AccessorData AccessorData
typename ::internal::DoFAccessorImplementation::Inheritance< structdim, dimension, space_dimension >::BaseClass BaseClass
Definition: dof_accessor.h:228
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
static ::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1403
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
unsigned int level
Definition: grid_out.cc:4355
Definition: hp.h:117
static bool is_level_cell()
DoFHandler< dim, spacedim > * dof_handler
Definition: dof_accessor.h:627
DoFInvalidAccessor(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0)
Definition: dof_accessor.cc:46
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
void get_active_fe_indices(const DoFHandler< dim, spacedim > &dof_handler, std::vector< unsigned int > &active_fe_indices)
Definition: dof_tools.cc:1350
static ::ExceptionBase & ExcCantCompareIterators()