Reference documentation for deal.II version Git 74ec1a6c4b 2021-01-19 10:11:04 -0500
\(\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 
740  friend struct ::internal::DoFHandlerImplementation::Policy::
741  Implementation;
742  friend struct ::internal::DoFHandlerImplementation::Implementation;
743  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
744  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
745  friend struct ::internal::DoFAccessorImplementation::Implementation;
746 };
747 
748 
749 
757 template <int spacedim, bool level_dof_access>
758 class DoFAccessor<0, 1, spacedim, level_dof_access>
759  : public TriaAccessor<0, 1, spacedim>
760 {
761 public:
766  static const unsigned int dimension = 1;
767 
772  static const unsigned int space_dimension = spacedim;
773 
779 
784 
795  DoFAccessor();
796 
813  DoFAccessor(
814  const Triangulation<1, spacedim> * tria,
815  const typename TriaAccessor<0, 1, spacedim>::VertexKind vertex_kind,
816  const unsigned int vertex_index,
817  const DoFHandler<1, spacedim> * dof_handler);
818 
825  const int = 0,
826  const int = 0,
827  const DoFHandler<1, spacedim> *dof_handler = 0);
828 
841  template <int structdim2, int dim2, int spacedim2>
843 
848  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
849  DoFAccessor(
851 
862  operator=(const DoFAccessor<0, 1, spacedim, level_dof_access> &da) = delete;
863 
872  get_dof_handler() const;
873 
877  template <bool level_dof_access2>
878  void
880 
885  void
886  copy_from(const TriaAccessorBase<0, 1, spacedim> &da);
887 
901  child(const unsigned int c) const;
902 
909  typename ::internal::DoFHandlerImplementation::
910  Iterators<1, spacedim, level_dof_access>::line_iterator
911  line(const unsigned int i) const;
912 
919  typename ::internal::DoFHandlerImplementation::
920  Iterators<1, spacedim, level_dof_access>::quad_iterator
921  quad(const unsigned int i) const;
922 
966  void
967  get_dof_indices(
968  std::vector<types::global_dof_index> &dof_indices,
969  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
970 
977  void
978  get_mg_dof_indices(
979  const int level,
980  std::vector<types::global_dof_index> &dof_indices,
981  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
982 
1002  vertex_dof_index(
1003  const unsigned int vertex,
1004  const unsigned int i,
1005  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1006 
1025  dof_index(const unsigned int i,
1026  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1027 
1046  unsigned int
1047  n_active_fe_indices() const;
1048 
1056  unsigned int
1057  nth_active_fe_index(const unsigned int n) const;
1058 
1067  bool
1068  fe_index_is_active(const unsigned int fe_index) const;
1069 
1076  get_fe(const unsigned int fe_index) const;
1077 
1087  DeclExceptionMsg(ExcInvalidObject,
1088  "This accessor object has not been "
1089  "associated with any DoFHandler object.");
1095  DeclException0(ExcVectorNotEmpty);
1101  DeclException0(ExcVectorDoesNotMatch);
1107  DeclException0(ExcMatrixDoesNotMatch);
1115  DeclException0(ExcNotActive);
1122 
1123 protected:
1128 
1132  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1133  bool
1134  operator==(
1136 
1140  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1141  bool
1142  operator!=(
1144 
1148  void set_dof_handler(DoFHandler<1, spacedim> *dh);
1149 
1168  void
1169  set_dof_index(
1170  const unsigned int i,
1171  const types::global_dof_index index,
1172  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1173 
1192  void
1193  set_vertex_dof_index(
1194  const unsigned int vertex,
1195  const unsigned int i,
1196  const types::global_dof_index index,
1197  const unsigned int fe_index = AccessorData::invalid_fe_index) const;
1198 
1199  // Iterator classes need to be friends because they need to access
1200  // operator== and operator!=.
1201  template <typename>
1202  friend class TriaRawIterator;
1203 
1204 
1205  // Make the DoFHandler class a friend so that it can call the set_xxx()
1206  // functions.
1207  template <int, int>
1208  friend class DoFHandler;
1209 
1210  friend struct ::internal::DoFHandlerImplementation::Policy::
1211  Implementation;
1212  friend struct ::internal::DoFHandlerImplementation::Implementation;
1213  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1214  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1215 };
1216 
1217 
1218 
1219 /* -------------------------------------------------------------------------- */
1220 
1221 
1242 template <int structdim, int dim, int spacedim = dim>
1243 class DoFInvalidAccessor : public InvalidAccessor<structdim, dim, spacedim>
1244 {
1245 public:
1249  using AccessorData =
1251 
1260  const int level = -1,
1261  const int index = -1,
1262  const AccessorData * local_data = 0);
1263 
1272 
1277  template <typename OtherAccessor>
1278  DoFInvalidAccessor(const OtherAccessor &);
1279 
1286  dof_index(const unsigned int i,
1287  const unsigned int fe_index =
1289 
1295  void
1296  set_dof_index(const unsigned int i,
1297  const types::global_dof_index index,
1298  const unsigned int fe_index =
1300 };
1301 
1302 
1303 
1304 /* -------------------------------------------------------------------------- */
1305 
1306 
1318 template <int dimension_, int space_dimension_, bool level_dof_access>
1319 class DoFCellAccessor : public DoFAccessor<dimension_,
1320  dimension_,
1321  space_dimension_,
1322  level_dof_access>
1323 {
1324 public:
1328  static const unsigned int dim = dimension_;
1329 
1333  static const unsigned int spacedim = space_dimension_;
1334 
1335 
1340 
1345  using BaseClass =
1347 
1352 
1357  using face_iterator = TriaIterator<DoFAccessor<dimension_ - 1,
1358  dimension_,
1359  space_dimension_,
1360  level_dof_access>>;
1361 
1373  const int level,
1374  const int index,
1375  const AccessorData *local_data);
1376 
1389  template <int structdim2, int dim2, int spacedim2>
1391 
1396  template <int structdim2, int dim2, int spacedim2, bool level_dof_access2>
1397  explicit DoFCellAccessor(
1399 
1410  operator=(
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 
1472  boost::container::small_vector<
1473  TriaIterator<
1476  child_iterators() const;
1477 
1485  face(const unsigned int i) const;
1486 
1490  boost::container::small_vector<face_iterator,
1492  face_iterators() const;
1493 
1501  neighbor_child_on_subface(const unsigned int face_no,
1502  const unsigned int subface_no) const;
1503 
1511  periodic_neighbor_child_on_subface(const unsigned int face_no,
1512  const unsigned int subface_no) const;
1513 
1542  template <class InputVector, typename number>
1543  void
1544  get_dof_values(const InputVector &values, Vector<number> &local_values) const;
1545 
1563  template <class InputVector, typename ForwardIterator>
1564  void
1565  get_dof_values(const InputVector &values,
1566  ForwardIterator local_values_begin,
1567  ForwardIterator local_values_end) const;
1568 
1589  template <class InputVector, typename ForwardIterator>
1590  void
1591  get_dof_values(
1593  const InputVector & values,
1594  ForwardIterator local_values_begin,
1595  ForwardIterator local_values_end) const;
1596 
1621  template <class OutputVector, typename number>
1622  void
1623  set_dof_values(const Vector<number> &local_values,
1624  OutputVector & values) const;
1625 
1657  template <class InputVector, typename number>
1658  void
1659  get_interpolated_dof_values(
1660  const InputVector &values,
1661  Vector<number> & interpolated_values,
1662  const unsigned int fe_index =
1664 
1717  template <class OutputVector, typename number>
1718  void
1719  set_dof_values_by_interpolation(
1720  const Vector<number> &local_values,
1721  OutputVector & values,
1722  const unsigned int fe_index =
1724 
1739  template <typename number, typename OutputVector>
1740  void
1741  distribute_local_to_global(const Vector<number> &local_source,
1742  OutputVector & global_destination) const;
1743 
1758  template <typename ForwardIterator, typename OutputVector>
1759  void
1760  distribute_local_to_global(ForwardIterator local_source_begin,
1761  ForwardIterator local_source_end,
1762  OutputVector & global_destination) const;
1763 
1777  template <typename ForwardIterator, typename OutputVector>
1778  void
1779  distribute_local_to_global(
1781  ForwardIterator local_source_begin,
1782  ForwardIterator local_source_end,
1783  OutputVector & global_destination) const;
1784 
1791  template <typename number, typename OutputMatrix>
1792  void
1793  distribute_local_to_global(const FullMatrix<number> &local_source,
1794  OutputMatrix &global_destination) const;
1795 
1800  template <typename number, typename OutputMatrix, typename OutputVector>
1801  void
1802  distribute_local_to_global(const FullMatrix<number> &local_matrix,
1803  const Vector<number> & local_vector,
1804  OutputMatrix & global_matrix,
1805  OutputVector & global_vector) const;
1806 
1831  void
1832  get_active_or_mg_dof_indices(
1833  std::vector<types::global_dof_index> &dof_indices) const;
1834 
1870  void
1871  get_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1872 
1877  void
1878  get_mg_dof_indices(std::vector<types::global_dof_index> &dof_indices) const;
1879 
1905  get_fe() const;
1906 
1932  unsigned int
1933  active_fe_index() const;
1934 
1961  void
1962  set_active_fe_index(const unsigned int i) const;
1971  void
1972  set_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
1973 
1977  void
1978  set_mg_dof_indices(const std::vector<types::global_dof_index> &dof_indices);
1979 
1983  void
1984  update_cell_dof_indices_cache() const;
1985 
2012  get_future_fe() const;
2013 
2032  unsigned int
2033  future_fe_index() const;
2034 
2043  void
2044  set_future_fe_index(const unsigned int i) const;
2045 
2052  bool
2053  future_fe_index_set() const;
2054 
2063  void
2064  clear_future_fe_index() const;
2065 
2084  unsigned int
2085  dominated_future_fe_on_children() const;
2103  ExcNoDominatedFiniteElementOnChildren,
2104  "No FiniteElement has been found in your FECollection that is "
2105  "dominated by all children of a cell you are trying to coarsen!");
2106 
2111 private:
2112  // Make the DoFHandler class a friend so that it can call the
2113  // update_cell_dof_indices_cache() function
2114  template <int, int>
2115  friend class DoFHandler;
2116  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
2117 };
2118 
2119 
2120 template <int structdim, int dim, int spacedim, bool level_dof_access>
2121 inline bool
2123 {
2124  return level_dof_access;
2125 }
2126 
2127 
2128 
2129 template <int structdim, int dim, int spacedim>
2130 template <typename OtherAccessor>
2132  const OtherAccessor &)
2133 {
2134  Assert(false,
2135  ExcMessage("You are attempting an illegal conversion between "
2136  "iterator/accessor types. The constructor you call "
2137  "only exists to make certain template constructs "
2138  "easier to write as dimension independent code but "
2139  "the conversion is not valid in the current context."));
2140 }
2141 
2142 
2143 
2145 
2146 // include more templates
2147 #include "dof_accessor.templates.h"
2148 
2149 
2150 #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:1466
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:380
unsigned int level
Definition: grid_out.cc:4362
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:379
void get_active_fe_indices(const DoFHandler< dim, spacedim > &dof_handler, std::vector< unsigned int > &active_fe_indices)
Definition: dof_tools.cc:1389
static ::ExceptionBase & ExcCantCompareIterators()