Reference documentation for deal.II version Git 50c3491829 2021-08-01 13:40:40 +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 - 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(
846  const Triangulation<1, spacedim> * tria,
847  const typename TriaAccessor<0, 1, spacedim>::VertexKind vertex_kind,
848  const unsigned int vertex_index,
849  const DoFHandler<1, spacedim> * dof_handler);
850 
857  const int = 0,
858  const int = 0,
859  const DoFHandler<1, spacedim> *dof_handler = 0);
860 
873  template <int structdim2, int dim2, int spacedim2>
875 
880  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
881  DoFAccessor(
883 
888 
892  // NOLINTNEXTLINE OSX does not compile with noexcept
894 
898  ~DoFAccessor() = default;
899 
910  operator=(const DoFAccessor<0, 1, spacedim, level_dof_access> &da) = delete;
911 
916  operator=(DoFAccessor<0, 1, spacedim, level_dof_access> &&) noexcept =
917  default;
918 
926  const DoFHandler<1, spacedim> &
927  get_dof_handler() const;
928 
932  template <bool level_dof_access2>
933  void
934  copy_from(const DoFAccessor<0, 1, spacedim, level_dof_access2> &a);
935 
940  void
941  copy_from(const TriaAccessorBase<0, 1, spacedim> &da);
942 
955  TriaIterator<DoFAccessor<0, 1, spacedim, level_dof_access>>
956  child(const unsigned int c) const;
957 
964  typename ::internal::DoFHandlerImplementation::
965  Iterators<1, spacedim, level_dof_access>::line_iterator
966  line(const unsigned int i) const;
967 
974  typename ::internal::DoFHandlerImplementation::
975  Iterators<1, spacedim, level_dof_access>::quad_iterator
976  quad(const unsigned int i) const;
977 
1021  void
1022  get_dof_indices(
1023  std::vector<types::global_dof_index> &dof_indices,
1024  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1025 
1032  void
1033  get_mg_dof_indices(
1034  const int level,
1035  std::vector<types::global_dof_index> &dof_indices,
1036  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1037 
1057  vertex_dof_index(
1058  const unsigned int vertex,
1059  const unsigned int i,
1060  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1061 
1080  dof_index(const unsigned int i,
1081  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1082 
1101  unsigned int
1102  n_active_fe_indices() const;
1103 
1111  unsigned int
1112  nth_active_fe_index(const unsigned int n) const;
1113 
1122  bool
1123  fe_index_is_active(const unsigned int fe_index) const;
1124 
1130  const FiniteElement<1, spacedim> &
1131  get_fe(const unsigned int fe_index) const;
1132 
1142  DeclExceptionMsg(ExcInvalidObject,
1143  "This accessor object has not been "
1144  "associated with any DoFHandler object.");
1150  DeclException0(ExcVectorNotEmpty);
1156  DeclException0(ExcVectorDoesNotMatch);
1162  DeclException0(ExcMatrixDoesNotMatch);
1170  DeclException0(ExcNotActive);
1177 
1178 protected:
1182  DoFHandler<1, spacedim> *dof_handler;
1183 
1187  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1188  bool
1189  operator==(
1190  const DoFAccessor<structdim2, dim2, spacedim2, level_dof_access2> &) const;
1191 
1195  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1196  bool
1197  operator!=(
1198  const DoFAccessor<structdim2, dim2, spacedim2, level_dof_access2> &) const;
1199 
1203  void
1204  set_dof_handler(DoFHandler<1, spacedim> *dh);
1205 
1224  void
1225  set_dof_index(
1226  const unsigned int i,
1227  const types::global_dof_index index,
1228  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1229 
1248  void
1249  set_vertex_dof_index(
1250  const unsigned int vertex,
1251  const unsigned int i,
1252  const types::global_dof_index index,
1253  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1254 
1255  // Iterator classes need to be friends because they need to access
1256  // operator== and operator!=.
1257  template <typename>
1258  friend class TriaRawIterator;
1259 
1260 
1261  // Make the DoFHandler class a friend so that it can call the set_xxx()
1262  // functions.
1263  template <int, int>
1264  friend class DoFHandler;
1265 
1266  friend struct ::internal::DoFHandlerImplementation::Policy::
1267  Implementation;
1268  friend struct ::internal::DoFHandlerImplementation::Implementation;
1269  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1270  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1271 };
1272 
1273 
1274 
1275 /* -------------------------------------------------------------------------- */
1276 
1277 
1298 template <int structdim, int dim, int spacedim = dim>
1299 class DoFInvalidAccessor : public InvalidAccessor<structdim, dim, spacedim>
1300 {
1301 public:
1305  using AccessorData =
1307 
1316  const int level = -1,
1317  const int index = -1,
1318  const AccessorData * local_data = 0);
1319 
1328 
1333  template <typename OtherAccessor>
1334  DoFInvalidAccessor(const OtherAccessor &);
1335 
1342  dof_index(const unsigned int i,
1343  const unsigned int fe_index =
1345 
1351  void
1352  set_dof_index(const unsigned int i,
1353  const types::global_dof_index index,
1354  const unsigned int fe_index =
1356 };
1357 
1358 
1359 
1360 /* -------------------------------------------------------------------------- */
1361 
1362 
1374 template <int dimension_, int space_dimension_, bool level_dof_access>
1375 class DoFCellAccessor : public DoFAccessor<dimension_,
1376  dimension_,
1377  space_dimension_,
1378  level_dof_access>
1379 {
1380 public:
1384  static const unsigned int dim = dimension_;
1385 
1389  static const unsigned int spacedim = space_dimension_;
1390 
1391 
1396 
1401  using BaseClass =
1403 
1408 
1413  using face_iterator = TriaIterator<DoFAccessor<dimension_ - 1,
1414  dimension_,
1415  space_dimension_,
1416  level_dof_access>>;
1417 
1429  const int level,
1430  const int index,
1431  const AccessorData *local_data);
1432 
1445  template <int structdim2, int dim2, int spacedim2>
1447 
1452  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1453  explicit DoFCellAccessor(
1455 
1461  default;
1462 
1466  DoFCellAccessor( // NOLINT
1468  &&) = default; // NOLINT
1469 
1473  ~DoFCellAccessor() = default;
1474 
1485  operator=(
1487  delete;
1488 
1493  operator=( // NOLINT
1495  &&) = default; // NOLINT
1496 
1511  parent() const;
1512 
1526  neighbor(const unsigned int i) const;
1527 
1534  periodic_neighbor(const unsigned int i) const;
1535 
1542  neighbor_or_periodic_neighbor(const unsigned int i) const;
1543 
1550  child(const unsigned int i) const;
1551 
1555  boost::container::small_vector<
1556  TriaIterator<
1559  child_iterators() const;
1560 
1568  face(const unsigned int i) const;
1569 
1573  boost::container::small_vector<face_iterator,
1575  face_iterators() const;
1576 
1584  neighbor_child_on_subface(const unsigned int face_no,
1585  const unsigned int subface_no) const;
1586 
1594  periodic_neighbor_child_on_subface(const unsigned int face_no,
1595  const unsigned int subface_no) const;
1596 
1625  template <class InputVector, typename number>
1626  void
1627  get_dof_values(const InputVector &values, Vector<number> &local_values) const;
1628 
1646  template <class InputVector, typename ForwardIterator>
1647  void
1648  get_dof_values(const InputVector &values,
1649  ForwardIterator local_values_begin,
1650  ForwardIterator local_values_end) const;
1651 
1672  template <class InputVector, typename ForwardIterator>
1673  void
1674  get_dof_values(
1676  const InputVector & values,
1677  ForwardIterator local_values_begin,
1678  ForwardIterator local_values_end) const;
1679 
1704  template <class OutputVector, typename number>
1705  void
1706  set_dof_values(const Vector<number> &local_values,
1707  OutputVector & values) const;
1708 
1740  template <class InputVector, typename number>
1741  void
1742  get_interpolated_dof_values(
1743  const InputVector &values,
1744  Vector<number> & interpolated_values,
1745  const unsigned int fe_index =
1747 
1800  template <class OutputVector, typename number>
1801  void
1802  set_dof_values_by_interpolation(
1803  const Vector<number> &local_values,
1804  OutputVector & values,
1805  const unsigned int fe_index =
1807 
1822  template <typename number, typename OutputVector>
1823  void
1824  distribute_local_to_global(const Vector<number> &local_source,
1825  OutputVector & global_destination) const;
1826 
1841  template <typename ForwardIterator, typename OutputVector>
1842  void
1843  distribute_local_to_global(ForwardIterator local_source_begin,
1844  ForwardIterator local_source_end,
1845  OutputVector & global_destination) const;
1846 
1860  template <typename ForwardIterator, typename OutputVector>
1861  void
1862  distribute_local_to_global(
1864  ForwardIterator local_source_begin,
1865  ForwardIterator local_source_end,
1866  OutputVector & global_destination) const;
1867 
1874  template <typename number, typename OutputMatrix>
1875  void
1876  distribute_local_to_global(const FullMatrix<number> &local_source,
1877  OutputMatrix &global_destination) const;
1878 
1883  template <typename number, typename OutputMatrix, typename OutputVector>
1884  void
1885  distribute_local_to_global(const FullMatrix<number> &local_matrix,
1886  const Vector<number> & local_vector,
1887  OutputMatrix & global_matrix,
1888  OutputVector & global_vector) const;
1889 
1914  void
1915  get_active_or_mg_dof_indices(
1916  std::vector<types::global_dof_index> &dof_indices) const;
1917 
1949  void
1950  get_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1951 
1956  void
1957  get_mg_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1958 
1984  get_fe() const;
1985 
2011  unsigned int
2012  active_fe_index() const;
2013 
2040  void
2041  set_active_fe_index(const unsigned int i) const;
2050  void
2051  set_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
2052 
2056  void
2057  set_mg_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
2058 
2062  void
2063  update_cell_dof_indices_cache() const;
2064 
2091  get_future_fe() const;
2092 
2111  unsigned int
2112  future_fe_index() const;
2113 
2122  void
2123  set_future_fe_index(const unsigned int i) const;
2124 
2131  bool
2132  future_fe_index_set() const;
2133 
2142  void
2143  clear_future_fe_index() const;
2148 private:
2149  // Make the DoFHandler class a friend so that it can call the
2150  // update_cell_dof_indices_cache() function
2151  template <int, int>
2152  friend class DoFHandler;
2153  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
2154 };
2155 
2156 
2157 template <int structdim, int dim, int spacedim, bool level_dof_access>
2158 inline bool
2160 {
2161  return level_dof_access;
2162 }
2163 
2164 
2165 
2166 template <int structdim, int dim, int spacedim>
2167 template <typename OtherAccessor>
2169  const OtherAccessor &)
2170 {
2171  Assert(false,
2172  ExcMessage("You are attempting an illegal conversion between "
2173  "iterator/accessor types. The constructor you call "
2174  "only exists to make certain template constructs "
2175  "easier to write as dimension independent code but "
2176  "the conversion is not valid in the current context."));
2177 }
2178 
2179 
2180 
2182 
2183 // include more templates
2184 #include "dof_accessor.templates.h"
2185 
2186 
2187 #endif
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
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:1467
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
#define DeclException0(Exception0)
Definition: exceptions.h:470
#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()