Reference documentation for deal.II version Git 8f0f58a98a 2021-10-21 17:29:01 -0600
\(\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 - 2021 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 
275  default;
276 
280  DoFAccessor( // NOLINT
282  default; // NOLINT
283 
287  ~DoFAccessor() = default;
288 
301  template <int structdim2, int dim2, int spacedim2>
303 
308  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
309  DoFAccessor(
311 
315  template <bool level_dof_access2>
317 
329  delete;
330 
335  operator=( // NOLINT
337  default; // NOLINT
338 
347  get_dof_handler() const;
348 
352  template <bool level_dof_access2>
353  void
355 
360  void
361  copy_from(const TriaAccessorBase<structdim, dim, spacedim> &da);
362 
367  static bool
368  is_level_cell();
369 
381  child(const unsigned int c) const;
382 
388  typename ::internal::DoFHandlerImplementation::
389  Iterators<dim, spacedim, level_dof_access>::line_iterator
390  line(const unsigned int i) const;
391 
397  typename ::internal::DoFHandlerImplementation::
398  Iterators<dim, spacedim, level_dof_access>::quad_iterator
399  quad(const unsigned int i) const;
400 
446  void
447  get_dof_indices(std::vector<types::global_dof_index> &dof_indices,
448  const unsigned int fe_index =
450 
457  void
458  get_mg_dof_indices(const int level,
459  std::vector<types::global_dof_index> &dof_indices,
460  const unsigned int fe_index =
462 
466  void
467  set_mg_dof_indices(
468  const int level,
469  const std::vector<types::global_dof_index> &dof_indices,
470  const unsigned int fe_index = DoFHandler<dim, spacedim>::invalid_fe_index);
471 
494  vertex_dof_index(const unsigned int vertex,
495  const unsigned int i,
496  const unsigned int fe_index =
498 
505  mg_vertex_dof_index(const int level,
506  const unsigned int vertex,
507  const unsigned int i,
508  const unsigned int fe_index =
510 
539  dof_index(const unsigned int i,
540  const unsigned int fe_index =
542 
547  mg_dof_index(const int level, const unsigned int i) const;
548 
570  unsigned int
571  n_active_fe_indices() const;
572 
580  unsigned int
581  nth_active_fe_index(const unsigned int n) const;
582 
589  std::set<unsigned int>
590  get_active_fe_indices() const;
591 
601  bool
602  fe_index_is_active(const unsigned int fe_index) const;
603 
610  get_fe(const unsigned int fe_index) const;
611 
621  DeclExceptionMsg(ExcInvalidObject,
622  "This accessor object has not been "
623  "associated with any DoFHandler object.");
629  DeclException0(ExcVectorNotEmpty);
635  DeclException0(ExcVectorDoesNotMatch);
641  DeclException0(ExcMatrixDoesNotMatch);
649  DeclException0(ExcNotActive);
656 
657 protected:
662 
663 public:
677  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
678  bool
679  operator==(
681 
685  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
686  bool
687  operator!=(
689 
690 protected:
694  void
695  set_dof_handler(DoFHandler<dim, spacedim> *dh);
696 
715  void
716  set_dof_index(const unsigned int i,
717  const types::global_dof_index index,
718  const unsigned int fe_index =
720 
721  void
722  set_mg_dof_index(const int level,
723  const unsigned int i,
724  const types::global_dof_index index) const;
725 
744  void
745  set_vertex_dof_index(const unsigned int vertex,
746  const unsigned int i,
747  const types::global_dof_index index,
748  const unsigned int fe_index =
750 
751  void
752  set_mg_vertex_dof_index(const int level,
753  const unsigned int vertex,
754  const unsigned int i,
755  const types::global_dof_index index,
756  const unsigned int fe_index =
758 
759  // Iterator classes need to be friends because they need to access
760  // operator== and operator!=.
761  template <typename>
762  friend class TriaRawIterator;
763  template <int, int, int, bool>
764  friend class DoFAccessor;
765 
766 private:
767  // Make the DoFHandler class a friend so that it can call the set_xxx()
768  // functions.
769  template <int, int>
770  friend class DoFHandler;
771 
772  friend struct ::internal::DoFHandlerImplementation::Policy::
773  Implementation;
774  friend struct ::internal::DoFHandlerImplementation::Implementation;
775  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
776  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
777  friend struct ::internal::DoFAccessorImplementation::Implementation;
778 };
779 
780 
781 
789 template <int spacedim, bool level_dof_access>
790 class DoFAccessor<0, 1, spacedim, level_dof_access>
791  : public TriaAccessor<0, 1, spacedim>
792 {
793 public:
798  static const unsigned int dimension = 1;
799 
804  static const unsigned int space_dimension = spacedim;
805 
811 
816 
827  DoFAccessor();
828 
845  DoFAccessor(
847  const typename TriaAccessor<0, 1, spacedim>::VertexKind vertex_kind,
848  const unsigned int vertex_index,
849  const DoFHandler<1, spacedim> * dof_handler);
850 
858  const int = 0,
859  const int = 0,
860  const DoFHandler<1, spacedim> *dof_handler = 0);
861 
874  template <int structdim2, int dim2, int spacedim2>
876 
881  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
882  DoFAccessor(
884 
889 
893  // NOLINTNEXTLINE OSX does not compile with noexcept
895 
899  ~DoFAccessor() = default;
900 
911  operator=(const DoFAccessor<0, 1, spacedim, level_dof_access> &da) = delete;
912 
917  operator=(DoFAccessor<0, 1, spacedim, level_dof_access> &&) noexcept =
918  default;
919 
927  const DoFHandler<1, spacedim> &
928  get_dof_handler() const;
929 
933  template <bool level_dof_access2>
934  void
935  copy_from(const DoFAccessor<0, 1, spacedim, level_dof_access2> &a);
936 
941  void
942  copy_from(const TriaAccessorBase<0, 1, spacedim> &da);
943 
956  TriaIterator<DoFAccessor<0, 1, spacedim, level_dof_access>>
957  child(const unsigned int c) const;
958 
965  typename ::internal::DoFHandlerImplementation::
966  Iterators<1, spacedim, level_dof_access>::line_iterator
967  line(const unsigned int i) const;
968 
975  typename ::internal::DoFHandlerImplementation::
976  Iterators<1, spacedim, level_dof_access>::quad_iterator
977  quad(const unsigned int i) const;
978 
1022  void
1023  get_dof_indices(
1024  std::vector<types::global_dof_index> &dof_indices,
1025  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1026 
1033  void
1034  get_mg_dof_indices(
1035  const int level,
1036  std::vector<types::global_dof_index> &dof_indices,
1037  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1038 
1058  vertex_dof_index(
1059  const unsigned int vertex,
1060  const unsigned int i,
1061  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1062 
1081  dof_index(const unsigned int i,
1082  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1083 
1102  unsigned int
1103  n_active_fe_indices() const;
1104 
1112  unsigned int
1113  nth_active_fe_index(const unsigned int n) const;
1114 
1123  bool
1124  fe_index_is_active(const unsigned int fe_index) const;
1125 
1131  const FiniteElement<1, spacedim> &
1132  get_fe(const unsigned int fe_index) const;
1133 
1143  DeclExceptionMsg(ExcInvalidObject,
1144  "This accessor object has not been "
1145  "associated with any DoFHandler object.");
1151  DeclException0(ExcVectorNotEmpty);
1157  DeclException0(ExcVectorDoesNotMatch);
1163  DeclException0(ExcMatrixDoesNotMatch);
1171  DeclException0(ExcNotActive);
1178 
1179 protected:
1183  DoFHandler<1, spacedim> *dof_handler;
1184 
1188  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1189  bool
1190  operator==(
1191  const DoFAccessor<structdim2, dim2, spacedim2, level_dof_access2> &) const;
1192 
1196  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1197  bool
1198  operator!=(
1199  const DoFAccessor<structdim2, dim2, spacedim2, level_dof_access2> &) const;
1200 
1204  void
1205  set_dof_handler(DoFHandler<1, spacedim> *dh);
1206 
1225  void
1226  set_dof_index(
1227  const unsigned int i,
1228  const types::global_dof_index index,
1229  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1230 
1249  void
1250  set_vertex_dof_index(
1251  const unsigned int vertex,
1252  const unsigned int i,
1253  const types::global_dof_index index,
1254  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1255 
1256  // Iterator classes need to be friends because they need to access
1257  // operator== and operator!=.
1258  template <typename>
1259  friend class TriaRawIterator;
1260 
1261 
1262  // Make the DoFHandler class a friend so that it can call the set_xxx()
1263  // functions.
1264  template <int, int>
1265  friend class DoFHandler;
1266 
1267  friend struct ::internal::DoFHandlerImplementation::Policy::
1268  Implementation;
1269  friend struct ::internal::DoFHandlerImplementation::Implementation;
1270  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1271  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1272 };
1273 
1274 
1275 
1276 /* -------------------------------------------------------------------------- */
1277 
1278 
1299 template <int structdim, int dim, int spacedim = dim>
1300 class DoFInvalidAccessor : public InvalidAccessor<structdim, dim, spacedim>
1301 {
1302 public:
1306  using AccessorData =
1308 
1317  const int level = -1,
1318  const int index = -1,
1319  const AccessorData * local_data = 0);
1320 
1329 
1334  template <typename OtherAccessor>
1335  DoFInvalidAccessor(const OtherAccessor &);
1336 
1343  dof_index(const unsigned int i,
1344  const unsigned int fe_index =
1346 
1352  void
1353  set_dof_index(const unsigned int i,
1354  const types::global_dof_index index,
1355  const unsigned int fe_index =
1357 };
1358 
1359 
1360 
1361 /* -------------------------------------------------------------------------- */
1362 
1363 
1375 template <int dimension_, int space_dimension_, bool level_dof_access>
1376 class DoFCellAccessor : public DoFAccessor<dimension_,
1377  dimension_,
1378  space_dimension_,
1379  level_dof_access>
1380 {
1381 public:
1385  static const unsigned int dim = dimension_;
1386 
1390  static const unsigned int spacedim = space_dimension_;
1391 
1392 
1397 
1402  using BaseClass =
1404 
1409 
1414  using face_iterator = TriaIterator<DoFAccessor<dimension_ - 1,
1415  dimension_,
1416  space_dimension_,
1417  level_dof_access>>;
1418 
1430  const int level,
1431  const int index,
1432  const AccessorData *local_data);
1433 
1446  template <int structdim2, int dim2, int spacedim2>
1448 
1453  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1454  explicit DoFCellAccessor(
1456 
1462  default;
1463 
1467  DoFCellAccessor( // NOLINT
1469  &&) = default; // NOLINT
1470 
1474  ~DoFCellAccessor() = default;
1475 
1486  operator=(
1488  delete;
1489 
1494  operator=( // NOLINT
1496  &&) = default; // NOLINT
1497 
1512  parent() const;
1513 
1527  neighbor(const unsigned int i) const;
1528 
1535  periodic_neighbor(const unsigned int i) const;
1536 
1543  neighbor_or_periodic_neighbor(const unsigned int i) const;
1544 
1551  child(const unsigned int i) const;
1552 
1556  boost::container::small_vector<
1557  TriaIterator<
1560  child_iterators() const;
1561 
1569  face(const unsigned int i) const;
1570 
1574  boost::container::small_vector<face_iterator,
1576  face_iterators() const;
1577 
1585  neighbor_child_on_subface(const unsigned int face_no,
1586  const unsigned int subface_no) const;
1587 
1595  periodic_neighbor_child_on_subface(const unsigned int face_no,
1596  const unsigned int subface_no) const;
1597 
1626  template <class InputVector, typename number>
1627  void
1628  get_dof_values(const InputVector &values, Vector<number> &local_values) const;
1629 
1647  template <class InputVector, typename ForwardIterator>
1648  void
1649  get_dof_values(const InputVector &values,
1650  ForwardIterator local_values_begin,
1651  ForwardIterator local_values_end) const;
1652 
1673  template <class InputVector, typename ForwardIterator>
1674  void
1675  get_dof_values(
1677  const InputVector & values,
1678  ForwardIterator local_values_begin,
1679  ForwardIterator local_values_end) const;
1680 
1705  template <class OutputVector, typename number>
1706  void
1707  set_dof_values(const Vector<number> &local_values,
1708  OutputVector & values) const;
1709 
1741  template <class InputVector, typename number>
1742  void
1743  get_interpolated_dof_values(
1744  const InputVector &values,
1745  Vector<number> & interpolated_values,
1746  const unsigned int fe_index =
1748 
1801  template <class OutputVector, typename number>
1802  void
1803  set_dof_values_by_interpolation(
1804  const Vector<number> &local_values,
1805  OutputVector & values,
1806  const unsigned int fe_index =
1808 
1823  template <typename number, typename OutputVector>
1824  void
1825  distribute_local_to_global(const Vector<number> &local_source,
1826  OutputVector & global_destination) const;
1827 
1842  template <typename ForwardIterator, typename OutputVector>
1843  void
1844  distribute_local_to_global(ForwardIterator local_source_begin,
1845  ForwardIterator local_source_end,
1846  OutputVector & global_destination) const;
1847 
1861  template <typename ForwardIterator, typename OutputVector>
1862  void
1863  distribute_local_to_global(
1865  ForwardIterator local_source_begin,
1866  ForwardIterator local_source_end,
1867  OutputVector & global_destination) const;
1868 
1875  template <typename number, typename OutputMatrix>
1876  void
1877  distribute_local_to_global(const FullMatrix<number> &local_source,
1878  OutputMatrix &global_destination) const;
1879 
1884  template <typename number, typename OutputMatrix, typename OutputVector>
1885  void
1886  distribute_local_to_global(const FullMatrix<number> &local_matrix,
1887  const Vector<number> & local_vector,
1888  OutputMatrix & global_matrix,
1889  OutputVector & global_vector) const;
1890 
1915  void
1916  get_active_or_mg_dof_indices(
1917  std::vector<types::global_dof_index> &dof_indices) const;
1918 
1950  void
1951  get_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1952 
1957  void
1958  get_mg_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1959 
1985  get_fe() const;
1986 
2012  unsigned int
2013  active_fe_index() const;
2014 
2041  void
2042  set_active_fe_index(const unsigned int i) const;
2051  void
2052  set_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
2053 
2057  void
2058  set_mg_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
2059 
2063  void
2064  update_cell_dof_indices_cache() const;
2065 
2092  get_future_fe() const;
2093 
2112  unsigned int
2113  future_fe_index() const;
2114 
2123  void
2124  set_future_fe_index(const unsigned int i) const;
2125 
2132  bool
2133  future_fe_index_set() const;
2134 
2143  void
2144  clear_future_fe_index() const;
2149 private:
2150  // Make the DoFHandler class a friend so that it can call the
2151  // update_cell_dof_indices_cache() function
2152  template <int, int>
2153  friend class DoFHandler;
2154  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
2155 };
2156 
2157 
2158 template <int structdim, int dim, int spacedim, bool level_dof_access>
2159 inline bool
2161 {
2162  return level_dof_access;
2163 }
2164 
2165 
2166 
2167 template <int structdim, int dim, int spacedim>
2168 template <typename OtherAccessor>
2170  const OtherAccessor &)
2171 {
2172  Assert(false,
2173  ExcMessage("You are attempting an illegal conversion between "
2174  "iterator/accessor types. The constructor you call "
2175  "only exists to make certain template constructs "
2176  "easier to write as dimension independent code but "
2177  "the conversion is not valid in the current context."));
2178 }
2179 
2180 
2181 
2183 
2184 // include more templates
2185 #include "dof_accessor.templates.h"
2186 
2187 
2188 #endif
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
const ::Triangulation< dim, spacedim > & tria
STL namespace.
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:414
static ::ExceptionBase & ExcMessage(std::string arg1)
Definition: types.h:31
#define Assert(cond, exc)
Definition: exceptions.h:1461
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
#define DeclException0(Exception0)
Definition: exceptions.h:464
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
unsigned int level
Definition: grid_out.cc:4589
Definition: hp.h:117
unsigned int global_dof_index
Definition: types.h:76
static bool is_level_cell()
DoFHandler< dim, spacedim > * dof_handler
Definition: dof_accessor.h:661
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:452
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
void get_active_fe_indices(const DoFHandler< dim, spacedim > &dof_handler, std::vector< unsigned int > &active_fe_indices)
Definition: dof_tools.cc:1322
static ::ExceptionBase & ExcCantCompareIterators()