Reference documentation for deal.II version Git f81eda9982 2020-03-28 21:30:57 -0400
\(\newcommand{\vcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\vcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
dof_accessor.h
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 
22 #include <deal.II/dofs/dof_handler.h>
23 
24 #include <deal.II/grid/tria_accessor.h>
25 
26 #include <deal.II/hp/dof_handler.h>
27 
28 #include <vector>
29 
30 DEAL_II_NAMESPACE_OPEN
31 
32 // Forward declarations
33 #ifndef DOXYGEN
34 template <typename number>
35 class FullMatrix;
36 template <typename number>
37 class Vector;
38 template <typename number>
39 class AffineConstraints;
40 
41 template <typename Accessor>
42 class TriaRawIterator;
43 
44 template <int, int>
45 class FiniteElement;
46 
47 namespace internal
48 {
49  namespace DoFCellAccessorImplementation
50  {
51  struct Implementation;
52  }
53 
54  namespace DoFHandlerImplementation
55  {
56  struct Implementation;
57  namespace Policy
58  {
59  struct Implementation;
60  }
61  } // namespace DoFHandlerImplementation
62 
63  namespace hp
64  {
65  namespace DoFHandlerImplementation
66  {
67  struct Implementation;
68  }
69  } // namespace hp
70 } // namespace internal
71 #endif
72 
73 // note: the file dof_accessor.templates.h is included at the end of
74 // this file. this includes a lot of templates and thus makes
75 // compilation slower, but at the same time allows for more aggressive
76 // inlining and thus faster code.
77 
78 
79 namespace internal
80 {
81  namespace DoFAccessorImplementation
82  {
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  } // namespace DoFAccessorImplementation
126 } // namespace internal
127 
128 
129 /* -------------------------------------------------------------------------- */
130 
131 
132 
210 template <int structdim, typename DoFHandlerType, bool level_dof_access>
213  structdim,
214  DoFHandlerType::dimension,
215  DoFHandlerType::space_dimension>::BaseClass
216 {
217 public:
222  static const unsigned int dimension = DoFHandlerType::dimension;
223 
228  static const unsigned int space_dimension = DoFHandlerType::space_dimension;
229 
234  using BaseClass = typename ::internal::DoFAccessorImplementation::
235  Inheritance<structdim, dimension, space_dimension>::BaseClass;
236 
240  using AccessorData = DoFHandlerType;
241 
252  DoFAccessor();
253 
270  DoFAccessor(const Triangulation<DoFHandlerType::dimension,
271  DoFHandlerType::space_dimension> *tria,
272  const int level,
273  const int index,
274  const DoFHandlerType *dof_handler);
275 
288  template <int structdim2, int dim2, int spacedim2>
290 
295  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
297 
301  template <bool level_dof_access2>
302  DoFAccessor(
304 
316  &da) = delete;
317 
325  const DoFHandlerType &
326  get_dof_handler() const;
327 
331  template <bool level_dof_access2>
332  void
334 
339  void
340  copy_from(const TriaAccessorBase<structdim,
341  DoFHandlerType::dimension,
342  DoFHandlerType::space_dimension> &da);
343 
348  static bool
349  is_level_cell();
350 
362  child(const unsigned int c) const;
363 
369  typename ::internal::DoFHandlerImplementation::
370  Iterators<DoFHandlerType, level_dof_access>::line_iterator
371  line(const unsigned int i) const;
372 
378  typename ::internal::DoFHandlerImplementation::
379  Iterators<DoFHandlerType, level_dof_access>::quad_iterator
380  quad(const unsigned int i) const;
381 
427  void
428  get_dof_indices(
429  std::vector<types::global_dof_index> &dof_indices,
430  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
431 
438  void
439  get_mg_dof_indices(
440  const int level,
441  std::vector<types::global_dof_index> &dof_indices,
442  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
443 
447  void
448  set_mg_dof_indices(
449  const int level,
450  const std::vector<types::global_dof_index> &dof_indices,
451  const unsigned int fe_index = DoFHandlerType::default_fe_index);
452 
471  vertex_dof_index(
472  const unsigned int vertex,
473  const unsigned int i,
474  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
475 
482  mg_vertex_dof_index(
483  const int level,
484  const unsigned int vertex,
485  const unsigned int i,
486  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
487 
516  dof_index(
517  const unsigned int i,
518  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
519 
524  mg_dof_index(const int level, const unsigned int i) const;
525 
547  unsigned int
548  n_active_fe_indices() const;
549 
557  unsigned int
558  nth_active_fe_index(const unsigned int n) const;
559 
566  std::set<unsigned int>
567  get_active_fe_indices() const;
568 
577  bool
578  fe_index_is_active(const unsigned int fe_index) const;
579 
585  const FiniteElement<DoFHandlerType::dimension,
586  DoFHandlerType::space_dimension> &
587  get_fe(const unsigned int fe_index) const;
588 
598  DeclExceptionMsg(ExcInvalidObject,
599  "This accessor object has not been "
600  "associated with any DoFHandler object.");
606  DeclException0(ExcVectorNotEmpty);
612  DeclException0(ExcVectorDoesNotMatch);
618  DeclException0(ExcMatrixDoesNotMatch);
626  DeclException0(ExcNotActive);
633 
634 protected:
638  DoFHandlerType *dof_handler;
639 
640 public:
654  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
655  bool
656  operator==(
658 
662  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
663  bool
664  operator!=(
666 
667 protected:
671  void
672  set_dof_handler(DoFHandlerType *dh);
673 
691  void
692  set_dof_index(
693  const unsigned int i,
694  const types::global_dof_index index,
695  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
696 
697  void
698  set_mg_dof_index(const int level,
699  const unsigned int i,
700  const types::global_dof_index index) const;
701 
719  void
720  set_vertex_dof_index(
721  const unsigned int vertex,
722  const unsigned int i,
723  const types::global_dof_index index,
724  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
725 
726  void
727  set_mg_vertex_dof_index(
728  const int level,
729  const unsigned int vertex,
730  const unsigned int i,
731  const types::global_dof_index index,
732  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
733 
734  // Iterator classes need to be friends because they need to access
735  // operator== and operator!=.
736  template <typename>
737  friend class TriaRawIterator;
738  template <int, class, bool>
739  friend class DoFAccessor;
740 
741 private:
742  // Make the DoFHandler class a friend so that it can call the set_xxx()
743  // functions.
744  template <int dim, int spacedim>
745  friend class DoFHandler;
746  template <int dim, int spacedim>
747  friend class hp::DoFHandler;
748 
749  friend struct ::internal::DoFHandlerImplementation::Policy::
750  Implementation;
751  friend struct ::internal::DoFHandlerImplementation::Implementation;
752  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
753  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
754  friend struct ::internal::DoFAccessorImplementation::Implementation;
755 };
756 
757 
758 
768 template <template <int, int> class DoFHandlerType,
769  int spacedim,
770  bool level_dof_access>
771 class DoFAccessor<0, DoFHandlerType<1, spacedim>, level_dof_access>
772  : public TriaAccessor<0, 1, spacedim>
773 {
774 public:
779  static const unsigned int dimension = 1;
780 
785  static const unsigned int space_dimension = spacedim;
786 
792 
796  using AccessorData = DoFHandlerType<1, spacedim>;
797 
808  DoFAccessor();
809 
826  DoFAccessor(
827  const Triangulation<1, spacedim> * tria,
828  const typename TriaAccessor<0, 1, spacedim>::VertexKind vertex_kind,
829  const unsigned int vertex_index,
830  const DoFHandlerType<1, spacedim> * dof_handler);
831 
838  const int = 0,
839  const int = 0,
840  const DoFHandlerType<1, spacedim> *dof_handler = 0);
841 
854  template <int structdim2, int dim2, int spacedim2>
856 
861  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
863 
873  DoFAccessor<0, DoFHandlerType<1, spacedim>, level_dof_access> &
874  operator=(const DoFAccessor<0, DoFHandlerType<1, spacedim>, level_dof_access>
875  &da) = delete;
876 
884  const DoFHandlerType<1, spacedim> &
885  get_dof_handler() const;
886 
890  template <bool level_dof_access2>
891  void
892  copy_from(
893  const DoFAccessor<0, DoFHandlerType<1, spacedim>, level_dof_access2> &a);
894 
899  void
900  copy_from(const TriaAccessorBase<0, 1, spacedim> &da);
901 
915  child(const unsigned int c) const;
916 
923  typename ::internal::DoFHandlerImplementation::
924  Iterators<DoFHandlerType<1, spacedim>, level_dof_access>::line_iterator
925  line(const unsigned int i) const;
926 
933  typename ::internal::DoFHandlerImplementation::
934  Iterators<DoFHandlerType<1, spacedim>, level_dof_access>::quad_iterator
935  quad(const unsigned int i) const;
936 
979  void
980  get_dof_indices(
981  std::vector<types::global_dof_index> &dof_indices,
982  const unsigned int fe_index = AccessorData::default_fe_index) const;
983 
990  void
991  get_mg_dof_indices(
992  const int level,
993  std::vector<types::global_dof_index> &dof_indices,
994  const unsigned int fe_index = AccessorData::default_fe_index) const;
995 
1014  vertex_dof_index(
1015  const unsigned int vertex,
1016  const unsigned int i,
1017  const unsigned int fe_index = AccessorData::default_fe_index) const;
1018 
1034  dof_index(const unsigned int i,
1035  const unsigned int fe_index = AccessorData::default_fe_index) const;
1036 
1055  unsigned int
1056  n_active_fe_indices() const;
1057 
1065  unsigned int
1066  nth_active_fe_index(const unsigned int n) const;
1067 
1076  bool
1077  fe_index_is_active(const unsigned int fe_index) const;
1078 
1085  DoFHandlerType<1, spacedim>::space_dimension> &
1086  get_fe(const unsigned int fe_index) const;
1087 
1097  DeclExceptionMsg(ExcInvalidObject,
1098  "This accessor object has not been "
1099  "associated with any DoFHandler object.");
1105  DeclException0(ExcVectorNotEmpty);
1111  DeclException0(ExcVectorDoesNotMatch);
1117  DeclException0(ExcMatrixDoesNotMatch);
1125  DeclException0(ExcNotActive);
1132 
1133 protected:
1137  DoFHandlerType<1, spacedim> *dof_handler;
1138 
1142  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
1143  bool
1144  operator==(
1146 
1150  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
1151  bool
1152  operator!=(
1154 
1158  void set_dof_handler(DoFHandlerType<1, spacedim> *dh);
1159 
1177  void
1178  set_dof_index(
1179  const unsigned int i,
1180  const types::global_dof_index index,
1181  const unsigned int fe_index = AccessorData::default_fe_index) const;
1182 
1200  void
1201  set_vertex_dof_index(
1202  const unsigned int vertex,
1203  const unsigned int i,
1204  const types::global_dof_index index,
1205  const unsigned int fe_index = AccessorData::default_fe_index) const;
1206 
1207  // Iterator classes need to be friends because they need to access
1208  // operator== and operator!=.
1209  template <typename>
1210  friend class TriaRawIterator;
1211 
1212 
1213  // Make the DoFHandler class a friend so that it can call the set_xxx()
1214  // functions.
1215  template <int, int>
1216  friend class DoFHandler;
1217  template <int, int>
1218  friend class hp::DoFHandler;
1219 
1220  friend struct ::internal::DoFHandlerImplementation::Policy::
1221  Implementation;
1222  friend struct ::internal::DoFHandlerImplementation::Implementation;
1223  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1224  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1225 };
1226 
1227 
1228 
1229 /* -------------------------------------------------------------------------- */
1230 
1231 
1253 template <int structdim, int dim, int spacedim = dim>
1254 class DoFInvalidAccessor : public InvalidAccessor<structdim, dim, spacedim>
1255 {
1256 public:
1260  using AccessorData =
1262 
1271  const int level = -1,
1272  const int index = -1,
1273  const AccessorData * local_data = 0);
1274 
1283 
1288  template <typename OtherAccessor>
1289  DoFInvalidAccessor(const OtherAccessor &);
1290 
1296  void
1297  set_dof_index(const unsigned int i,
1298  const types::global_dof_index index,
1299  const unsigned int fe_index =
1301 };
1302 
1303 
1304 
1305 /* -------------------------------------------------------------------------- */
1306 
1307 
1320 template <typename DoFHandlerType, bool level_dof_access>
1321 class DoFCellAccessor : public DoFAccessor<DoFHandlerType::dimension,
1322  DoFHandlerType,
1323  level_dof_access>
1324 {
1325 public:
1329  static const unsigned int dim = DoFHandlerType::dimension;
1330 
1334  static const unsigned int spacedim = DoFHandlerType::space_dimension;
1335 
1336 
1340  using AccessorData = DoFHandlerType;
1341 
1346  using BaseClass =
1348 
1352  using Container = DoFHandlerType;
1353 
1358  using face_iterator = TriaIterator<DoFAccessor<DoFHandlerType::dimension - 1,
1359  DoFHandlerType,
1360  level_dof_access>>;
1361 
1372  DoFCellAccessor(const Triangulation<DoFHandlerType::dimension,
1373  DoFHandlerType::space_dimension> *tria,
1374  const int level,
1375  const int index,
1376  const AccessorData *local_data);
1377 
1390  template <int structdim2, int dim2, int spacedim2>
1392 
1397  template <int dim2, class DoFHandlerType2, bool level_dof_access2>
1398  explicit DoFCellAccessor(
1400 
1412  delete;
1413 
1428  parent() const;
1429 
1443  neighbor(const unsigned int i) const;
1444 
1451  periodic_neighbor(const unsigned int i) const;
1452 
1459  neighbor_or_periodic_neighbor(const unsigned int i) const;
1460 
1467  child(const unsigned int i) const;
1468 
1476  face(const unsigned int i) const;
1477 
1481  inline std::array<face_iterator,
1483  face_iterators() const;
1484 
1492  neighbor_child_on_subface(const unsigned int face_no,
1493  const unsigned int subface_no) const;
1494 
1502  periodic_neighbor_child_on_subface(const unsigned int face_no,
1503  const unsigned int subface_no) const;
1504 
1530  template <class InputVector, typename number>
1531  void
1532  get_dof_values(const InputVector &values, Vector<number> &local_values) const;
1533 
1548  template <class InputVector, typename ForwardIterator>
1549  void
1550  get_dof_values(const InputVector &values,
1551  ForwardIterator local_values_begin,
1552  ForwardIterator local_values_end) const;
1553 
1571  template <class InputVector, typename ForwardIterator>
1572  void
1573  get_dof_values(
1575  const InputVector & values,
1576  ForwardIterator local_values_begin,
1577  ForwardIterator local_values_end) const;
1578 
1600  template <class OutputVector, typename number>
1601  void
1602  set_dof_values(const Vector<number> &local_values,
1603  OutputVector & values) const;
1604 
1636  template <class InputVector, typename number>
1637  void
1638  get_interpolated_dof_values(
1639  const InputVector &values,
1640  Vector<number> & interpolated_values,
1641  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
1642 
1695  template <class OutputVector, typename number>
1696  void
1697  set_dof_values_by_interpolation(
1698  const Vector<number> &local_values,
1699  OutputVector & values,
1700  const unsigned int fe_index = DoFHandlerType::default_fe_index) const;
1701 
1710  template <typename number, typename OutputVector>
1711  void
1712  distribute_local_to_global(const Vector<number> &local_source,
1713  OutputVector & global_destination) const;
1714 
1723  template <typename ForwardIterator, typename OutputVector>
1724  void
1725  distribute_local_to_global(ForwardIterator local_source_begin,
1726  ForwardIterator local_source_end,
1727  OutputVector & global_destination) const;
1728 
1739  template <typename ForwardIterator, typename OutputVector>
1740  void
1741  distribute_local_to_global(
1743  ForwardIterator local_source_begin,
1744  ForwardIterator local_source_end,
1745  OutputVector & global_destination) const;
1746 
1753  template <typename number, typename OutputMatrix>
1754  void
1755  distribute_local_to_global(const FullMatrix<number> &local_source,
1756  OutputMatrix &global_destination) const;
1757 
1762  template <typename number, typename OutputMatrix, typename OutputVector>
1763  void
1764  distribute_local_to_global(const FullMatrix<number> &local_matrix,
1765  const Vector<number> & local_vector,
1766  OutputMatrix & global_matrix,
1767  OutputVector & global_vector) const;
1768 
1793  void
1794  get_active_or_mg_dof_indices(
1795  std::vector<types::global_dof_index> &dof_indices) const;
1796 
1832  void
1833  get_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1834 
1839  void
1840  get_mg_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1841 
1866  const FiniteElement<DoFHandlerType::dimension,
1867  DoFHandlerType::space_dimension> &
1868  get_fe() const;
1869 
1895  unsigned int
1896  active_fe_index() const;
1897 
1924  void
1925  set_active_fe_index(const unsigned int i) const;
1934  void
1935  set_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
1936 
1940  void
1941  set_mg_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
1942 
1946  void
1947  update_cell_dof_indices_cache() const;
1948 
1974  const FiniteElement<DoFHandlerType::dimension,
1975  DoFHandlerType::space_dimension> &
1976  get_future_fe() const;
1977 
1996  unsigned int
1997  future_fe_index() const;
1998 
2007  void
2008  set_future_fe_index(const unsigned int i) const;
2009 
2016  bool
2017  future_fe_index_set() const;
2018 
2027  void
2028  clear_future_fe_index() const;
2033 private:
2034  // Make the DoFHandler class a friend so that it can call the
2035  // update_cell_dof_indices_cache() function
2036  template <int dim, int spacedim>
2037  friend class DoFHandler;
2038  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
2039 };
2040 
2041 
2042 template <int sd, typename DoFHandlerType, bool level_dof_access>
2043 inline bool
2045 {
2046  return level_dof_access;
2047 }
2048 
2049 
2050 
2051 template <int structdim, int dim, int spacedim>
2052 template <typename OtherAccessor>
2054  const OtherAccessor &)
2055 {
2056  Assert(false,
2057  ExcMessage("You are attempting an illegal conversion between "
2058  "iterator/accessor types. The constructor you call "
2059  "only exists to make certain template constructs "
2060  "easier to write as dimension independent code but "
2061  "the conversion is not valid in the current context."));
2062 }
2063 
2064 
2065 
2066 DEAL_II_NAMESPACE_CLOSE
2067 
2068 // include more templates
2069 #include "dof_accessor.templates.h"
2070 
2071 
2072 #endif
DoFHandlerType * dof_handler
Definition: dof_accessor.h:638
static bool is_level_cell()
typename ::internal::DoFAccessorImplementation::Inheritance< structdim, dimension, space_dimension >::BaseClass BaseClass
Definition: dof_accessor.h:235
typename InvalidAccessor< structdim, dim, spacedim >::AccessorData AccessorData
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
LinearAlgebra::distributed::Vector< Number > Vector
static ::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1419
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
#define DeclException0(Exception0)
Definition: exceptions.h:473
Definition: hp.h:117
unsigned int global_dof_index
Definition: types.h:77
DoFHandlerType Container
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
DoFHandlerType AccessorData
static ::ExceptionBase & ExcCantCompareIterators()