Reference documentation for deal.II version Git 78a8940608 2021-04-17 09:24:19 -0400
\(\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 
30 #include <boost/container/small_vector.hpp>
32 
33 #include <vector>
34 
36 
37 // Forward declarations
38 #ifndef DOXYGEN
39 template <typename number>
40 class FullMatrix;
41 template <typename number>
42 class Vector;
43 template <typename number>
44 class AffineConstraints;
45 
46 template <typename Accessor>
47 class TriaRawIterator;
48 
49 template <int, int>
50 class FiniteElement;
51 
52 namespace internal
53 {
54  namespace DoFCellAccessorImplementation
55  {
56  struct Implementation;
57  }
58 
59  namespace DoFHandlerImplementation
60  {
61  struct Implementation;
62  namespace Policy
63  {
64  struct Implementation;
65  }
66  } // namespace DoFHandlerImplementation
67 
68  namespace hp
69  {
70  namespace DoFHandlerImplementation
71  {
72  struct Implementation;
73  }
74  } // namespace hp
75 } // namespace internal
76 #endif
77 
78 // note: the file dof_accessor.templates.h is included at the end of
79 // this file. this includes a lot of templates and thus makes
80 // compilation slower, but at the same time allows for more aggressive
81 // inlining and thus faster code.
82 
83 
84 namespace internal
85 {
86  namespace DoFAccessorImplementation
87  {
104  template <int structdim, int dim, int spacedim>
105  struct Inheritance
106  {
112  };
113 
114 
120  template <int dim, int spacedim>
121  struct Inheritance<dim, dim, spacedim>
122  {
128  };
129 
130  struct Implementation;
131  } // namespace DoFAccessorImplementation
132 } // namespace internal
133 
134 
135 /* -------------------------------------------------------------------------- */
136 
137 
138 
210 template <int structdim, int dim, int spacedim, bool level_dof_access>
211 class DoFAccessor : public ::internal::DoFAccessorImplementation::
212  Inheritance<structdim, dim, spacedim>::BaseClass
213 {
214 public:
219  static const unsigned int dimension = dim;
220 
225  static const unsigned int space_dimension = spacedim;
226 
231  using BaseClass = typename ::internal::DoFAccessorImplementation::
232  Inheritance<structdim, dimension, space_dimension>::BaseClass;
233 
238 
249  DoFAccessor();
250 
267  const int level,
268  const int index,
269  const DoFHandler<dim, spacedim> * dof_handler);
270 
283  template <int structdim2, int dim2, int spacedim2>
285 
290  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
291  DoFAccessor(
293 
297  template <bool level_dof_access2>
299 
311  delete;
312 
321  get_dof_handler() const;
322 
326  template <bool level_dof_access2>
327  void
329 
334  void
335  copy_from(const TriaAccessorBase<structdim, dim, spacedim> &da);
336 
341  static bool
342  is_level_cell();
343 
355  child(const unsigned int c) const;
356 
362  typename ::internal::DoFHandlerImplementation::
363  Iterators<dim, spacedim, level_dof_access>::line_iterator
364  line(const unsigned int i) const;
365 
371  typename ::internal::DoFHandlerImplementation::
372  Iterators<dim, spacedim, level_dof_access>::quad_iterator
373  quad(const unsigned int i) const;
374 
420  void
421  get_dof_indices(std::vector<types::global_dof_index> &dof_indices,
422  const unsigned int fe_index =
424 
431  void
432  get_mg_dof_indices(const int level,
433  std::vector<types::global_dof_index> &dof_indices,
434  const unsigned int fe_index =
436 
440  void
441  set_mg_dof_indices(
442  const int level,
443  const std::vector<types::global_dof_index> &dof_indices,
444  const unsigned int fe_index = DoFHandler<dim, spacedim>::invalid_fe_index);
445 
464  vertex_dof_index(const unsigned int vertex,
465  const unsigned int i,
466  const unsigned int fe_index =
468 
475  mg_vertex_dof_index(const int level,
476  const unsigned int vertex,
477  const unsigned int i,
478  const unsigned int fe_index =
480 
509  dof_index(const unsigned int i,
510  const unsigned int fe_index =
512 
517  mg_dof_index(const int level, const unsigned int i) const;
518 
540  unsigned int
541  n_active_fe_indices() const;
542 
550  unsigned int
551  nth_active_fe_index(const unsigned int n) const;
552 
559  std::set<unsigned int>
560  get_active_fe_indices() const;
561 
571  bool
572  fe_index_is_active(const unsigned int fe_index) const;
573 
580  get_fe(const unsigned int fe_index) const;
581 
591  DeclExceptionMsg(ExcInvalidObject,
592  "This accessor object has not been "
593  "associated with any DoFHandler object.");
599  DeclException0(ExcVectorNotEmpty);
605  DeclException0(ExcVectorDoesNotMatch);
611  DeclException0(ExcMatrixDoesNotMatch);
619  DeclException0(ExcNotActive);
626 
627 protected:
632 
633 public:
647  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
648  bool
649  operator==(
651 
655  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
656  bool
657  operator!=(
659 
660 protected:
664  void
665  set_dof_handler(DoFHandler<dim, spacedim> *dh);
666 
685  void
686  set_dof_index(const unsigned int i,
687  const types::global_dof_index index,
688  const unsigned int fe_index =
690 
691  void
692  set_mg_dof_index(const int level,
693  const unsigned int i,
694  const types::global_dof_index index) const;
695 
714  void
715  set_vertex_dof_index(const unsigned int vertex,
716  const unsigned int i,
717  const types::global_dof_index index,
718  const unsigned int fe_index =
720 
721  void
722  set_mg_vertex_dof_index(const int level,
723  const unsigned int vertex,
724  const unsigned int i,
725  const types::global_dof_index index,
726  const unsigned int fe_index =
728 
729  // Iterator classes need to be friends because they need to access
730  // operator== and operator!=.
731  template <typename>
732  friend class TriaRawIterator;
733  template <int, int, int, bool>
734  friend class DoFAccessor;
735 
736 private:
737  // Make the DoFHandler class a friend so that it can call the set_xxx()
738  // functions.
739  template <int, int>
740  friend class DoFHandler;
741 
742  friend struct ::internal::DoFHandlerImplementation::Policy::
743  Implementation;
744  friend struct ::internal::DoFHandlerImplementation::Implementation;
745  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
746  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
747  friend struct ::internal::DoFAccessorImplementation::Implementation;
748 };
749 
750 
751 
759 template <int spacedim, bool level_dof_access>
760 class DoFAccessor<0, 1, spacedim, level_dof_access>
761  : public TriaAccessor<0, 1, spacedim>
762 {
763 public:
768  static const unsigned int dimension = 1;
769 
774  static const unsigned int space_dimension = spacedim;
775 
781 
786 
797  DoFAccessor();
798 
815  DoFAccessor(
816  const Triangulation<1, spacedim> * tria,
817  const typename TriaAccessor<0, 1, spacedim>::VertexKind vertex_kind,
818  const unsigned int vertex_index,
819  const DoFHandler<1, spacedim> * dof_handler);
820 
827  const int = 0,
828  const int = 0,
829  const DoFHandler<1, spacedim> *dof_handler = 0);
830 
843  template <int structdim2, int dim2, int spacedim2>
845 
850  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
851  DoFAccessor(
853 
864  operator=(const DoFAccessor<0, 1, spacedim, level_dof_access> &da) = delete;
865 
874  get_dof_handler() const;
875 
879  template <bool level_dof_access2>
880  void
882 
887  void
888  copy_from(const TriaAccessorBase<0, 1, spacedim> &da);
889 
903  child(const unsigned int c) const;
904 
911  typename ::internal::DoFHandlerImplementation::
912  Iterators<1, spacedim, level_dof_access>::line_iterator
913  line(const unsigned int i) const;
914 
921  typename ::internal::DoFHandlerImplementation::
922  Iterators<1, spacedim, level_dof_access>::quad_iterator
923  quad(const unsigned int i) const;
924 
968  void
969  get_dof_indices(
970  std::vector<types::global_dof_index> &dof_indices,
971  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
972 
979  void
980  get_mg_dof_indices(
981  const int level,
982  std::vector<types::global_dof_index> &dof_indices,
983  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
984 
1004  vertex_dof_index(
1005  const unsigned int vertex,
1006  const unsigned int i,
1007  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1008 
1027  dof_index(const unsigned int i,
1028  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1029 
1048  unsigned int
1049  n_active_fe_indices() const;
1050 
1058  unsigned int
1059  nth_active_fe_index(const unsigned int n) const;
1060 
1069  bool
1070  fe_index_is_active(const unsigned int fe_index) const;
1071 
1078  get_fe(const unsigned int fe_index) const;
1079 
1089  DeclExceptionMsg(ExcInvalidObject,
1090  "This accessor object has not been "
1091  "associated with any DoFHandler object.");
1097  DeclException0(ExcVectorNotEmpty);
1103  DeclException0(ExcVectorDoesNotMatch);
1109  DeclException0(ExcMatrixDoesNotMatch);
1117  DeclException0(ExcNotActive);
1124 
1125 protected:
1130 
1134  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1135  bool
1136  operator==(
1138 
1142  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1143  bool
1144  operator!=(
1146 
1150  void set_dof_handler(DoFHandler<1, spacedim> *dh);
1151 
1170  void
1171  set_dof_index(
1172  const unsigned int i,
1173  const types::global_dof_index index,
1174  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1175 
1194  void
1195  set_vertex_dof_index(
1196  const unsigned int vertex,
1197  const unsigned int i,
1198  const types::global_dof_index index,
1199  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1200 
1201  // Iterator classes need to be friends because they need to access
1202  // operator== and operator!=.
1203  template <typename>
1204  friend class TriaRawIterator;
1205 
1206 
1207  // Make the DoFHandler class a friend so that it can call the set_xxx()
1208  // functions.
1209  template <int, int>
1210  friend class 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 
1474  boost::container::small_vector<
1475  TriaIterator<
1478  child_iterators() const;
1479 
1487  face(const unsigned int i) const;
1488 
1492  boost::container::small_vector<face_iterator,
1494  face_iterators() const;
1495 
1503  neighbor_child_on_subface(const unsigned int face_no,
1504  const unsigned int subface_no) const;
1505 
1513  periodic_neighbor_child_on_subface(const unsigned int face_no,
1514  const unsigned int subface_no) const;
1515 
1544  template <class InputVector, typename number>
1545  void
1546  get_dof_values(const InputVector &values, Vector<number> &local_values) const;
1547 
1565  template <class InputVector, typename ForwardIterator>
1566  void
1567  get_dof_values(const InputVector &values,
1568  ForwardIterator local_values_begin,
1569  ForwardIterator local_values_end) const;
1570 
1591  template <class InputVector, typename ForwardIterator>
1592  void
1593  get_dof_values(
1595  const InputVector & values,
1596  ForwardIterator local_values_begin,
1597  ForwardIterator local_values_end) const;
1598 
1623  template <class OutputVector, typename number>
1624  void
1625  set_dof_values(const Vector<number> &local_values,
1626  OutputVector & values) const;
1627 
1659  template <class InputVector, typename number>
1660  void
1661  get_interpolated_dof_values(
1662  const InputVector &values,
1663  Vector<number> & interpolated_values,
1664  const unsigned int fe_index =
1666 
1719  template <class OutputVector, typename number>
1720  void
1721  set_dof_values_by_interpolation(
1722  const Vector<number> &local_values,
1723  OutputVector & values,
1724  const unsigned int fe_index =
1726 
1741  template <typename number, typename OutputVector>
1742  void
1743  distribute_local_to_global(const Vector<number> &local_source,
1744  OutputVector & global_destination) const;
1745 
1760  template <typename ForwardIterator, typename OutputVector>
1761  void
1762  distribute_local_to_global(ForwardIterator local_source_begin,
1763  ForwardIterator local_source_end,
1764  OutputVector & global_destination) const;
1765 
1779  template <typename ForwardIterator, typename OutputVector>
1780  void
1781  distribute_local_to_global(
1783  ForwardIterator local_source_begin,
1784  ForwardIterator local_source_end,
1785  OutputVector & global_destination) const;
1786 
1793  template <typename number, typename OutputMatrix>
1794  void
1795  distribute_local_to_global(const FullMatrix<number> &local_source,
1796  OutputMatrix &global_destination) const;
1797 
1802  template <typename number, typename OutputMatrix, typename OutputVector>
1803  void
1804  distribute_local_to_global(const FullMatrix<number> &local_matrix,
1805  const Vector<number> & local_vector,
1806  OutputMatrix & global_matrix,
1807  OutputVector & global_vector) const;
1808 
1833  void
1834  get_active_or_mg_dof_indices(
1835  std::vector<types::global_dof_index> &dof_indices) const;
1836 
1868  void
1869  get_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1870 
1875  void
1876  get_mg_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1877 
1903  get_fe() const;
1904 
1930  unsigned int
1931  active_fe_index() const;
1932 
1959  void
1960  set_active_fe_index(const unsigned int i) const;
1969  void
1970  set_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
1971 
1975  void
1976  set_mg_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
1977 
1981  void
1982  update_cell_dof_indices_cache() const;
1983 
2010  get_future_fe() const;
2011 
2030  unsigned int
2031  future_fe_index() const;
2032 
2041  void
2042  set_future_fe_index(const unsigned int i) const;
2043 
2050  bool
2051  future_fe_index_set() const;
2052 
2061  void
2062  clear_future_fe_index() const;
2067 private:
2068  // Make the DoFHandler class a friend so that it can call the
2069  // update_cell_dof_indices_cache() function
2070  template <int, int>
2071  friend class DoFHandler;
2072  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
2073 };
2074 
2075 
2076 template <int structdim, int dim, int spacedim, bool level_dof_access>
2077 inline bool
2079 {
2080  return level_dof_access;
2081 }
2082 
2083 
2084 
2085 template <int structdim, int dim, int spacedim>
2086 template <typename OtherAccessor>
2088  const OtherAccessor &)
2089 {
2090  Assert(false,
2091  ExcMessage("You are attempting an illegal conversion between "
2092  "iterator/accessor types. The constructor you call "
2093  "only exists to make certain template constructs "
2094  "easier to write as dimension independent code but "
2095  "the conversion is not valid in the current context."));
2096 }
2097 
2098 
2099 
2101 
2102 // include more templates
2103 #include "dof_accessor.templates.h"
2104 
2105 
2106 #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:232
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
Definition: config.h:409
static ::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1466
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:396
unsigned int level
Definition: grid_out.cc:4590
Definition: hp.h:117
static bool is_level_cell()
DoFHandler< dim, spacedim > * dof_handler
Definition: dof_accessor.h:631
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_ENABLE_EXTRA_DIAGNOSTICS
Definition: config.h:446
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:395
void get_active_fe_indices(const DoFHandler< dim, spacedim > &dof_handler, std::vector< unsigned int > &active_fe_indices)
Definition: dof_tools.cc:1399
static ::ExceptionBase & ExcCantCompareIterators()