Reference documentation for deal.II version Git 0096380e24 2020-08-06 21:03:09 -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 
29 #include <boost/container/small_vector.hpp>
30 
31 #include <vector>
32 
34 
35 // Forward declarations
36 #ifndef DOXYGEN
37 template <typename number>
38 class FullMatrix;
39 template <typename number>
40 class Vector;
41 template <typename number>
42 class AffineConstraints;
43 
44 template <typename Accessor>
45 class TriaRawIterator;
46 
47 template <int, int>
48 class FiniteElement;
49 
50 namespace internal
51 {
52  namespace DoFCellAccessorImplementation
53  {
54  struct Implementation;
55  }
56 
57  namespace DoFHandlerImplementation
58  {
59  struct Implementation;
60  namespace Policy
61  {
62  struct Implementation;
63  }
64  } // namespace DoFHandlerImplementation
65 
66  namespace hp
67  {
68  namespace DoFHandlerImplementation
69  {
70  struct Implementation;
71  }
72  } // namespace hp
73 } // namespace internal
74 #endif
75 
76 // note: the file dof_accessor.templates.h is included at the end of
77 // this file. this includes a lot of templates and thus makes
78 // compilation slower, but at the same time allows for more aggressive
79 // inlining and thus faster code.
80 
81 
82 namespace internal
83 {
84  namespace DoFAccessorImplementation
85  {
102  template <int structdim, int dim, int spacedim>
103  struct Inheritance
104  {
110  };
111 
112 
118  template <int dim, int spacedim>
119  struct Inheritance<dim, dim, spacedim>
120  {
126  };
127 
128  struct Implementation;
129  } // namespace DoFAccessorImplementation
130 } // namespace internal
131 
132 
133 /* -------------------------------------------------------------------------- */
134 
135 
136 
208 template <int structdim, int dim, int spacedim, bool level_dof_access>
209 class DoFAccessor : public ::internal::DoFAccessorImplementation::
210  Inheritance<structdim, dim, spacedim>::BaseClass
211 {
212 public:
217  static const unsigned int dimension = dim;
218 
223  static const unsigned int space_dimension = spacedim;
224 
229  using BaseClass = typename ::internal::DoFAccessorImplementation::
230  Inheritance<structdim, dimension, space_dimension>::BaseClass;
231 
236 
247  DoFAccessor();
248 
265  const int level,
266  const int index,
267  const DoFHandler<dim, spacedim> * dof_handler);
268 
281  template <int structdim2, int dim2, int spacedim2>
283 
288  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
289  DoFAccessor(
291 
295  template <bool level_dof_access2>
297 
309  delete;
310 
319  get_dof_handler() const;
320 
324  template <bool level_dof_access2>
325  void
327 
332  void
333  copy_from(const TriaAccessorBase<structdim, dim, spacedim> &da);
334 
339  static bool
340  is_level_cell();
341 
353  child(const unsigned int c) const;
354 
360  typename ::internal::DoFHandlerImplementation::
361  Iterators<dim, spacedim, level_dof_access>::line_iterator
362  line(const unsigned int i) const;
363 
369  typename ::internal::DoFHandlerImplementation::
370  Iterators<dim, spacedim, level_dof_access>::quad_iterator
371  quad(const unsigned int i) const;
372 
418  void
419  get_dof_indices(std::vector<types::global_dof_index> &dof_indices,
420  const unsigned int fe_index =
422 
429  void
430  get_mg_dof_indices(const int level,
431  std::vector<types::global_dof_index> &dof_indices,
432  const unsigned int fe_index =
434 
438  void
439  set_mg_dof_indices(
440  const int level,
441  const std::vector<types::global_dof_index> &dof_indices,
442  const unsigned int fe_index = DoFHandler<dim, spacedim>::invalid_fe_index);
443 
462  vertex_dof_index(const unsigned int vertex,
463  const unsigned int i,
464  const unsigned int fe_index =
466 
473  mg_vertex_dof_index(const int level,
474  const unsigned int vertex,
475  const unsigned int i,
476  const unsigned int fe_index =
478 
507  dof_index(const unsigned int i,
508  const unsigned int fe_index =
510 
515  mg_dof_index(const int level, const unsigned int i) const;
516 
538  unsigned int
539  n_active_fe_indices() const;
540 
548  unsigned int
549  nth_active_fe_index(const unsigned int n) const;
550 
557  std::set<unsigned int>
558  get_active_fe_indices() const;
559 
569  bool
570  fe_index_is_active(const unsigned int fe_index) const;
571 
578  get_fe(const unsigned int fe_index) const;
579 
589  DeclExceptionMsg(ExcInvalidObject,
590  "This accessor object has not been "
591  "associated with any DoFHandler object.");
597  DeclException0(ExcVectorNotEmpty);
603  DeclException0(ExcVectorDoesNotMatch);
609  DeclException0(ExcMatrixDoesNotMatch);
617  DeclException0(ExcNotActive);
624 
625 protected:
630 
631 public:
645  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
646  bool
647  operator==(
649 
653  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
654  bool
655  operator!=(
657 
658 protected:
662  void
663  set_dof_handler(DoFHandler<dim, spacedim> *dh);
664 
683  void
684  set_dof_index(const unsigned int i,
685  const types::global_dof_index index,
686  const unsigned int fe_index =
688 
689  void
690  set_mg_dof_index(const int level,
691  const unsigned int i,
692  const types::global_dof_index index) const;
693 
712  void
713  set_vertex_dof_index(const unsigned int vertex,
714  const unsigned int i,
715  const types::global_dof_index index,
716  const unsigned int fe_index =
718 
719  void
720  set_mg_vertex_dof_index(const int level,
721  const unsigned int vertex,
722  const unsigned int i,
723  const types::global_dof_index index,
724  const unsigned int fe_index =
726 
727  // Iterator classes need to be friends because they need to access
728  // operator== and operator!=.
729  template <typename>
730  friend class TriaRawIterator;
731  template <int, int, int, bool>
732  friend class DoFAccessor;
733 
734 private:
735  // Make the DoFHandler class a friend so that it can call the set_xxx()
736  // functions.
737  template <int, int>
738  friend class DoFHandler;
739  template <int, int>
740  friend class hp::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  template <int, int>
1212  friend class hp::DoFHandler;
1213 
1214  friend struct ::internal::DoFHandlerImplementation::Policy::
1215  Implementation;
1216  friend struct ::internal::DoFHandlerImplementation::Implementation;
1217  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1218  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1219 };
1220 
1221 
1222 
1223 /* -------------------------------------------------------------------------- */
1224 
1225 
1246 template <int structdim, int dim, int spacedim = dim>
1247 class DoFInvalidAccessor : public InvalidAccessor<structdim, dim, spacedim>
1248 {
1249 public:
1253  using AccessorData =
1255 
1264  const int level = -1,
1265  const int index = -1,
1266  const AccessorData * local_data = 0);
1267 
1276 
1281  template <typename OtherAccessor>
1282  DoFInvalidAccessor(const OtherAccessor &);
1283 
1290  dof_index(const unsigned int i,
1291  const unsigned int fe_index =
1293 
1299  void
1300  set_dof_index(const unsigned int i,
1301  const types::global_dof_index index,
1302  const unsigned int fe_index =
1304 };
1305 
1306 
1307 
1308 /* -------------------------------------------------------------------------- */
1309 
1310 
1322 template <int dimension_, int space_dimension_, bool level_dof_access>
1323 class DoFCellAccessor : public DoFAccessor<dimension_,
1324  dimension_,
1325  space_dimension_,
1326  level_dof_access>
1327 {
1328 public:
1332  static const unsigned int dim = dimension_;
1333 
1337  static const unsigned int spacedim = space_dimension_;
1338 
1339 
1344 
1349  using BaseClass =
1351 
1356 
1361  using face_iterator = TriaIterator<DoFAccessor<dimension_ - 1,
1362  dimension_,
1363  space_dimension_,
1364  level_dof_access>>;
1365 
1377  const int level,
1378  const int index,
1379  const AccessorData *local_data);
1380 
1393  template <int structdim2, int dim2, int spacedim2>
1395 
1400  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1401  explicit DoFCellAccessor(
1403 
1414  operator=(
1416  delete;
1417 
1432  parent() const;
1433 
1447  neighbor(const unsigned int i) const;
1448 
1455  periodic_neighbor(const unsigned int i) const;
1456 
1463  neighbor_or_periodic_neighbor(const unsigned int i) const;
1464 
1471  child(const unsigned int i) const;
1472 
1480  face(const unsigned int i) const;
1481 
1485  inline boost::container::
1486  small_vector<face_iterator, GeometryInfo<dimension_>::faces_per_cell>
1487  face_iterators() const;
1488 
1496  neighbor_child_on_subface(const unsigned int face_no,
1497  const unsigned int subface_no) const;
1498 
1506  periodic_neighbor_child_on_subface(const unsigned int face_no,
1507  const unsigned int subface_no) const;
1508 
1534  template <class InputVector, typename number>
1535  void
1536  get_dof_values(const InputVector &values, Vector<number> &local_values) const;
1537 
1552  template <class InputVector, typename ForwardIterator>
1553  void
1554  get_dof_values(const InputVector &values,
1555  ForwardIterator local_values_begin,
1556  ForwardIterator local_values_end) const;
1557 
1575  template <class InputVector, typename ForwardIterator>
1576  void
1577  get_dof_values(
1579  const InputVector & values,
1580  ForwardIterator local_values_begin,
1581  ForwardIterator local_values_end) const;
1582 
1604  template <class OutputVector, typename number>
1605  void
1606  set_dof_values(const Vector<number> &local_values,
1607  OutputVector & values) const;
1608 
1640  template <class InputVector, typename number>
1641  void
1642  get_interpolated_dof_values(
1643  const InputVector &values,
1644  Vector<number> & interpolated_values,
1645  const unsigned int fe_index =
1647 
1700  template <class OutputVector, typename number>
1701  void
1702  set_dof_values_by_interpolation(
1703  const Vector<number> &local_values,
1704  OutputVector & values,
1705  const unsigned int fe_index =
1707 
1716  template <typename number, typename OutputVector>
1717  void
1718  distribute_local_to_global(const Vector<number> &local_source,
1719  OutputVector & global_destination) const;
1720 
1729  template <typename ForwardIterator, typename OutputVector>
1730  void
1731  distribute_local_to_global(ForwardIterator local_source_begin,
1732  ForwardIterator local_source_end,
1733  OutputVector & global_destination) const;
1734 
1745  template <typename ForwardIterator, typename OutputVector>
1746  void
1747  distribute_local_to_global(
1749  ForwardIterator local_source_begin,
1750  ForwardIterator local_source_end,
1751  OutputVector & global_destination) const;
1752 
1759  template <typename number, typename OutputMatrix>
1760  void
1761  distribute_local_to_global(const FullMatrix<number> &local_source,
1762  OutputMatrix &global_destination) const;
1763 
1768  template <typename number, typename OutputMatrix, typename OutputVector>
1769  void
1770  distribute_local_to_global(const FullMatrix<number> &local_matrix,
1771  const Vector<number> & local_vector,
1772  OutputMatrix & global_matrix,
1773  OutputVector & global_vector) const;
1774 
1799  void
1800  get_active_or_mg_dof_indices(
1801  std::vector<types::global_dof_index> &dof_indices) const;
1802 
1838  void
1839  get_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1840 
1845  void
1846  get_mg_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1847 
1873  get_fe() const;
1874 
1900  unsigned int
1901  active_fe_index() const;
1902 
1929  void
1930  set_active_fe_index(const unsigned int i) const;
1939  void
1940  set_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
1941 
1945  void
1946  set_mg_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
1947 
1951  void
1952  update_cell_dof_indices_cache() const;
1953 
1980  get_future_fe() const;
1981 
2000  unsigned int
2001  future_fe_index() const;
2002 
2011  void
2012  set_future_fe_index(const unsigned int i) const;
2013 
2020  bool
2021  future_fe_index_set() const;
2022 
2031  void
2032  clear_future_fe_index() const;
2037 private:
2038  // Make the DoFHandler class a friend so that it can call the
2039  // update_cell_dof_indices_cache() function
2040  template <int, int>
2041  friend class DoFHandler;
2042  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
2043 };
2044 
2045 
2046 template <int structdim, int dim, int spacedim, bool level_dof_access>
2047 inline bool
2049 {
2050  return level_dof_access;
2051 }
2052 
2053 
2054 
2055 template <int structdim, int dim, int spacedim>
2056 template <typename OtherAccessor>
2058  const OtherAccessor &)
2059 {
2060  Assert(false,
2061  ExcMessage("You are attempting an illegal conversion between "
2062  "iterator/accessor types. The constructor you call "
2063  "only exists to make certain template constructs "
2064  "easier to write as dimension independent code but "
2065  "the conversion is not valid in the current context."));
2066 }
2067 
2068 
2069 
2071 
2072 // include more templates
2073 #include "dof_accessor.templates.h"
2074 
2075 
2076 #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:230
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:1411
#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:4341
Definition: hp.h:117
static bool is_level_cell()
DoFHandler< dim, spacedim > * dof_handler
Definition: dof_accessor.h:629
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:1352
static ::ExceptionBase & ExcCantCompareIterators()