Reference documentation for deal.II version 8.1.0
dof_accessor.h
1 // ---------------------------------------------------------------------
2 // @f$Id: dof_accessor.h 30040 2013-07-18 17:06:48Z maier @f$
3 //
4 // Copyright (C) 1998 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__dof_accessor_h
18 #define __deal2__dof_accessor_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/grid/tria_accessor.h>
23 #include <deal.II/dofs/dof_handler.h>
24 #include <deal.II/hp/dof_handler.h>
25 
26 #include <vector>
27 
29 
30 template <typename number> class FullMatrix;
31 template <typename number> class SparseMatrix;
32 template <typename number> class Vector;
33 class ConstraintMatrix;
34 
35 template <typename Accessor> class TriaRawIterator;
36 
37 template <int, int> class FiniteElement;
38 
39 
40 namespace internal
41 {
42  namespace DoFCellAccessor
43  {
44  struct Implementation;
45  }
46 
47  namespace DoFHandler
48  {
49  struct Implementation;
50  namespace Policy
51  {
52  struct Implementation;
53  }
54  }
55 
56  namespace hp
57  {
58  namespace DoFHandler
59  {
60  struct Implementation;
61  }
62  }
63 }
64 
65 // note: the file dof_accessor.templates.h is included at the end of
66 // this file. this includes a lot of templates and thus makes
67 // compilation slower, but at the same time allows for more aggressive
68 // inlining and thus faster code.
69 
70 
71 namespace internal
72 {
73  namespace DoFAccessor
74  {
92  template <int structdim, int dim, int spacedim>
93  struct Inheritance
94  {
100  typedef ::TriaAccessor<structdim,dim,spacedim> BaseClass;
101  };
102 
103 
109  template <int dim, int spacedim>
110  struct Inheritance<dim,dim,spacedim>
111  {
117  typedef ::CellAccessor<dim,spacedim> BaseClass;
118  };
119  }
120 }
121 
122 
123 /* -------------------------------------------------------------------------- */
124 
125 
126 
188 template <int structdim, class DH, bool level_dof_access>
189 class DoFAccessor : public ::internal::DoFAccessor::Inheritance<structdim, DH::dimension, DH::space_dimension>::BaseClass
190 {
191 public:
192 
198  static const unsigned int dimension=DH::dimension;
199 
205  static const unsigned int space_dimension=DH::space_dimension;
206 
213  typedef
214  typename ::internal::DoFAccessor::Inheritance<structdim, dimension, space_dimension>::BaseClass
216 
220  typedef DH AccessorData;
221 
234  DoFAccessor ();
235 
240  const int level,
241  const int index,
242  const DH *local_data);
243 
265  template <int structdim2, int dim2, int spacedim2>
267 
274  template <int dim2, class DH2, bool level_dof_access2>
276 
282  template <bool level_dof_access2>
293  const DH &
294  get_dof_handler () const;
295 
300  template <bool level_dof_access2>
302 
311 
317  static bool is_level_cell();
318 
324  parent () const;
325 
338  child (const unsigned int c) const;
339 
349  typename ::internal::DoFHandler::Iterators<DH, level_dof_access>::line_iterator
350  line (const unsigned int i) const;
351 
361  typename ::internal::DoFHandler::Iterators<DH, level_dof_access>::quad_iterator
362  quad (const unsigned int i) const;
363 
412  void get_dof_indices (std::vector<types::global_dof_index> &dof_indices,
413  const unsigned int fe_index = DH::default_fe_index) const;
414 
421  void get_mg_dof_indices (const int level,
422  std::vector<types::global_dof_index> &dof_indices,
423  const unsigned int fe_index = DH::default_fe_index) const;
424 
428  void set_mg_dof_indices (const int level,
429  const std::vector<types::global_dof_index> &dof_indices,
430  const unsigned int fe_index = DH::default_fe_index);
431 
468  types::global_dof_index vertex_dof_index (const unsigned int vertex,
469  const unsigned int i,
470  const unsigned int fe_index = DH::default_fe_index) const;
471 
478  const unsigned int vertex,
479  const unsigned int i,
480  const unsigned int fe_index = DH::default_fe_index) const;
481 
536  types::global_dof_index dof_index (const unsigned int i,
537  const unsigned int fe_index = DH::default_fe_index) const;
538 
542  types::global_dof_index mg_dof_index (const int level, const unsigned int i) const;
543 
575  unsigned int
576  n_active_fe_indices () const;
577 
592  unsigned int
593  nth_active_fe_index (const unsigned int n) const;
594 
612  bool
613  fe_index_is_active (const unsigned int fe_index) const;
614 
624  get_fe (const unsigned int fe_index) const;
625 
635  DeclException0 (ExcInvalidObject);
641  DeclException0 (ExcVectorNotEmpty);
647  DeclException0 (ExcVectorDoesNotMatch);
653  DeclException0 (ExcMatrixDoesNotMatch);
661  DeclException0 (ExcNotActive);
667  DeclException0 (ExcCantCompareIterators);
668 
669 protected:
670 
676 public:
701  template <int dim2, class DH2, bool level_dof_access2>
703 
708  template <int dim2, class DH2, bool level_dof_access2>
710 protected:
714  void set_dof_handler (DH *dh);
715 
752  void set_dof_index (const unsigned int i,
754  const unsigned int fe_index = DH::default_fe_index) const;
755 
756  void set_mg_dof_index (const int level, const unsigned int i, const types::global_dof_index index) const;
757 
794  void set_vertex_dof_index (const unsigned int vertex,
795  const unsigned int i,
797  const unsigned int fe_index = DH::default_fe_index) const;
798 
799  void set_mg_vertex_dof_index (const int level, const unsigned int vertex, const unsigned int i, const types::global_dof_index index, const unsigned int fe_index = DH::default_fe_index) const;
800 
806  template <typename> friend class TriaRawIterator;
807  template <int, class, bool> friend class DoFAccessor;
808 
809 private:
829 
835  template <int dim, int spacedim> friend class DoFHandler;
836  template <int dim, int spacedim> friend class hp::DoFHandler;
837 
838  friend struct ::internal::DoFHandler::Policy::Implementation;
839  friend struct ::internal::DoFHandler::Implementation;
840  friend struct ::internal::hp::DoFHandler::Implementation;
841  friend struct ::internal::DoFCellAccessor::Implementation;
842  friend struct ::internal::DoFAccessor::Implementation;
843 };
844 
845 
846 
857 template <template <int, int> class DH, int spacedim, bool level_dof_access>
858 class DoFAccessor<0,DH<1,spacedim>, level_dof_access> : public TriaAccessor<0,1,spacedim>
859 {
860 public:
861 
867  static const unsigned int dimension=1;
868 
874  static const unsigned int space_dimension=spacedim;
875 
883 
887  typedef DH<1,spacedim> AccessorData;
888 
901  DoFAccessor ();
902 
936  const typename TriaAccessor<0,1,spacedim>::VertexKind vertex_kind,
937  const unsigned int vertex_index,
938  const DH<1,spacedim> *dof_handler);
939 
950  const int = 0,
951  const int = 0,
952  const DH<1,spacedim> *dof_handler = 0);
953 
975  template <int structdim2, int dim2, int spacedim2>
977 
984  template <int dim2, class DH2, bool level_dof_access2>
986 
996  const DH<1,spacedim> &
997  get_dof_handler () const;
998 
1002  DoFAccessor<0,DH<1,spacedim>, level_dof_access> &
1003  operator = (const DoFAccessor<0,DH<1,spacedim>, level_dof_access> &da);
1004 
1009  template <bool level_dof_access2>
1010  void copy_from (const DoFAccessor<0, DH<1,spacedim>, level_dof_access2> &a);
1011 
1020 
1025  TriaIterator<DoFAccessor<0,DH<1,spacedim>, level_dof_access> >
1026  parent () const;
1027 
1039  TriaIterator<DoFAccessor<0,DH<1,spacedim>, level_dof_access > >
1040  child (const unsigned int c) const;
1041 
1051  typename ::internal::DoFHandler::Iterators<DH<1,spacedim>, level_dof_access>::line_iterator
1052  line (const unsigned int i) const;
1053 
1063  typename ::internal::DoFHandler::Iterators<DH<1,spacedim>, level_dof_access>::quad_iterator
1064  quad (const unsigned int i) const;
1065 
1135  void get_dof_indices (std::vector<types::global_dof_index> &dof_indices,
1136  const unsigned int fe_index = AccessorData::default_fe_index) const;
1137 
1174  types::global_dof_index vertex_dof_index (const unsigned int vertex,
1175  const unsigned int i,
1176  const unsigned int fe_index = AccessorData::default_fe_index) const;
1177 
1232  types::global_dof_index dof_index (const unsigned int i,
1233  const unsigned int fe_index = AccessorData::default_fe_index) const;
1234 
1266  unsigned int
1267  n_active_fe_indices () const;
1268 
1283  unsigned int
1284  nth_active_fe_index (const unsigned int n) const;
1285 
1303  bool
1304  fe_index_is_active (const unsigned int fe_index) const;
1305 
1314  const FiniteElement<DH<1,spacedim>::dimension,DH<1,spacedim>::space_dimension> &
1315  get_fe (const unsigned int fe_index) const;
1316 
1326  DeclException0 (ExcInvalidObject);
1332  DeclException0 (ExcVectorNotEmpty);
1338  DeclException0 (ExcVectorDoesNotMatch);
1344  DeclException0 (ExcMatrixDoesNotMatch);
1352  DeclException0 (ExcNotActive);
1358  DeclException0 (ExcCantCompareIterators);
1359 
1360 protected:
1361 
1366  DH<1,spacedim> *dof_handler;
1367 
1371  template <int dim2, class DH2, bool level_dof_access2>
1373 
1377  template <int dim2, class DH2, bool level_dof_access2>
1379 
1383  void set_dof_handler (DH<1,spacedim> *dh);
1384 
1421  void set_dof_index (const unsigned int i,
1423  const unsigned int fe_index = AccessorData::default_fe_index) const;
1424 
1461  void set_vertex_dof_index (const unsigned int vertex,
1462  const unsigned int i,
1464  const unsigned int fe_index = AccessorData::default_fe_index) const;
1465 
1471  template <typename> friend class TriaRawIterator;
1472 
1473 
1479  template <int, int> friend class DoFHandler;
1480  template <int, int> friend class hp::DoFHandler;
1481 
1482  friend struct ::internal::DoFHandler::Policy::Implementation;
1483  friend struct ::internal::DoFHandler::Implementation;
1484  friend struct ::internal::hp::DoFHandler::Implementation;
1485  friend struct ::internal::DoFCellAccessor::Implementation;
1486 };
1487 
1488 
1489 /* -------------------------------------------------------------------------- */
1490 
1491 
1504 template <class DH, bool level_dof_access>
1505 class DoFCellAccessor : public DoFAccessor<DH::dimension,DH, level_dof_access>
1506 {
1507 public:
1511  static const unsigned int dim = DH::dimension;
1512 
1516  static const unsigned int spacedim = DH::space_dimension;
1517 
1518 
1522  typedef DH AccessorData;
1523 
1531 
1536  typedef DH Container;
1537 
1549  const int level,
1550  const int index,
1551  const AccessorData *local_data);
1552 
1574  template <int structdim2, int dim2, int spacedim2>
1576 
1583  template <int dim2, class DH2, bool level_dof_access2>
1584  explicit
1586 
1601  parent () const;
1602 
1620  neighbor (const unsigned int) const;
1621 
1632  child (const unsigned int) const;
1633 
1642  TriaIterator<DoFAccessor<DH::dimension-1,DH, level_dof_access> >
1643  face (const unsigned int i) const;
1644 
1657  neighbor_child_on_subface (const unsigned int face_no,
1658  const unsigned int subface_no) const;
1659 
1699  template <class InputVector, typename number>
1700  void get_dof_values (const InputVector &values,
1701  Vector<number> &local_values) const;
1702 
1731  template <class InputVector, typename ForwardIterator>
1732  void get_dof_values (const InputVector &values,
1733  ForwardIterator local_values_begin,
1734  ForwardIterator local_values_end) const;
1735 
1767  template <class InputVector, typename ForwardIterator>
1768  void get_dof_values (const ConstraintMatrix &constraints,
1769  const InputVector &values,
1770  ForwardIterator local_values_begin,
1771  ForwardIterator local_values_end) const;
1772 
1807  template <class OutputVector, typename number>
1808  void set_dof_values (const Vector<number> &local_values,
1809  OutputVector &values) const;
1810 
1843  template <class InputVector, typename number>
1844  void get_interpolated_dof_values (const InputVector &values,
1845  Vector<number> &interpolated_values) const;
1846 
1936  template <class OutputVector, typename number>
1937  void set_dof_values_by_interpolation (const Vector<number> &local_values,
1938  OutputVector &values) const;
1939 
1956  template <typename number, typename OutputVector>
1957  void
1958  distribute_local_to_global (const Vector<number> &local_source,
1959  OutputVector &global_destination) const;
1960 
1974  template <typename ForwardIterator, typename OutputVector>
1975  void
1976  distribute_local_to_global (ForwardIterator local_source_begin,
1977  ForwardIterator local_source_end,
1978  OutputVector &global_destination) const;
1979 
1998  template <typename ForwardIterator, typename OutputVector>
1999  void
2000  distribute_local_to_global (const ConstraintMatrix &constraints,
2001  ForwardIterator local_source_begin,
2002  ForwardIterator local_source_end,
2003  OutputVector &global_destination) const;
2004 
2016  template <typename number, typename OutputMatrix>
2017  void
2018  distribute_local_to_global (const FullMatrix<number> &local_source,
2019  OutputMatrix &global_destination) const;
2020 
2027  template <typename number, typename OutputMatrix, typename OutputVector>
2028  void
2029  distribute_local_to_global (const FullMatrix<number> &local_matrix,
2030  const Vector<number> &local_vector,
2031  OutputMatrix &global_matrix,
2032  OutputVector &global_vector) const;
2033 
2062  void get_active_or_mg_dof_indices (std::vector<types::global_dof_index> &dof_indices) const;
2063 
2097  void get_dof_indices (std::vector<types::global_dof_index> &dof_indices) const;
2098 
2105  void get_mg_dof_indices (std::vector<types::global_dof_index> &dof_indices) const;
2106 
2130  get_fe () const;
2131 
2137  unsigned int active_fe_index () const;
2138 
2143  void set_active_fe_index (const unsigned int i);
2155  void set_dof_indices (const std::vector<types::global_dof_index> &dof_indices);
2156 
2161  void set_mg_dof_indices (const std::vector<types::global_dof_index> &dof_indices);
2162 
2168  void update_cell_dof_indices_cache () const;
2169 
2170 private:
2190 
2197  template <int dim, int spacedim> friend class DoFHandler;
2198  friend struct ::internal::DoFCellAccessor::Implementation;
2199 };
2200 
2201 
2202 template <int sd, class DH, bool lda>
2203 inline
2204 bool
2206 {
2207  return lda;
2208 }
2209 
2210 
2211 
2212 DEAL_II_NAMESPACE_CLOSE
2213 
2214 // include more templates
2215 #include "dof_accessor.templates.h"
2216 
2217 
2218 #endif
void set_dof_indices(const std::vector< types::global_dof_index > &dof_indices)
types::global_dof_index mg_dof_index(const int level, const unsigned int i) const
TriaIterator< DoFCellAccessor< DH, level_dof_access > > child(const unsigned int) const
::TriaAccessor< structdim, dim, spacedim > BaseClass
Definition: dof_accessor.h:100
TriaIterator< DoFCellAccessor< DH, level_dof_access > > neighbor(const unsigned int) const
void distribute_local_to_global(const Vector< number > &local_source, OutputVector &global_destination) const
const DH & get_dof_handler() const
void set_active_fe_index(const unsigned int i)
DoFAccessor< structdim, DH, level_dof_access > & operator=(const DoFAccessor< structdim, DH, level_dof_access > &da)
TriaIterator< DoFAccessor< structdim, DH, level_dof_access > > parent() const
unsigned int active_fe_index() const
bool operator==(const DoFAccessor< dim2, DH2, level_dof_access2 > &) const
void set_dof_handler(DH *dh)
void get_active_or_mg_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
typename::internal::DoFAccessor::Inheritance< structdim, dimension, space_dimension >::BaseClass BaseClass
Definition: dof_accessor.h:215
typename::internal::DoFHandler::Iterators< DH, level_dof_access >::line_iterator line(const unsigned int i) const
DeclException0(ExcInvalidObject)
DH * dof_handler
Definition: dof_accessor.h:675
types::global_dof_index mg_vertex_dof_index(const int level, const unsigned int vertex, const unsigned int i, const unsigned int fe_index=DH::default_fe_index) const
TriaIterator< DoFAccessor< structdim, DH, level_dof_access > > child(const unsigned int c) const
void set_mg_dof_indices(const int level, const std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DH::default_fe_index)
types::global_dof_index dof_index(const unsigned int i, const unsigned int fe_index=DH::default_fe_index) const
DoFCellAccessor< DH, level_dof_access > & operator=(const DoFCellAccessor< DH, level_dof_access > &da)
static const unsigned int dimension
Definition: dof_accessor.h:198
void set_dof_values_by_interpolation(const Vector< number > &local_values, OutputVector &values) const
unsigned int global_dof_index
Definition: types.h:100
void get_mg_dof_indices(const int level, std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DH::default_fe_index) const
static const unsigned int space_dimension
Definition: dof_accessor.h:205
void get_dof_values(const InputVector &values, Vector< number > &local_values) const
unsigned int vertex_index(const unsigned int i) const
const Triangulation< dim, spacedim > * tria
void set_mg_dof_indices(const std::vector< types::global_dof_index > &dof_indices)
void set_dof_values(const Vector< number > &local_values, OutputVector &values) const
DoFAccessor< DH::dimension, DH, level_dof_access > BaseClass
static const unsigned int dim
bool fe_index_is_active(const unsigned int fe_index) const
types::global_dof_index vertex_dof_index(const unsigned int vertex, const unsigned int i, const unsigned int fe_index=DH::default_fe_index) const
void get_interpolated_dof_values(const InputVector &values, Vector< number > &interpolated_values) const
Definition: hp.h:103
TriaIterator< DoFCellAccessor< DH, level_dof_access > > neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
const FiniteElement< DH::dimension, DH::space_dimension > & get_fe(const unsigned int fe_index) const
unsigned int nth_active_fe_index(const unsigned int n) const
void set_dof_index(const unsigned int i, const types::global_dof_index index, const unsigned int fe_index=DH::default_fe_index) const
void update_cell_dof_indices_cache() const
TriaIterator< DoFCellAccessor< DH, level_dof_access > > parent() const
static const unsigned int spacedim
TriaIterator< DoFAccessor< DH::dimension-1, DH, level_dof_access > > face(const unsigned int i) const
bool operator!=(const DoFAccessor< dim2, DH2, level_dof_access2 > &) const
static bool is_level_cell()
const FiniteElement< DH::dimension, DH::space_dimension > & get_fe() const
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DH::default_fe_index) const
void get_mg_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
DoFCellAccessor(const Triangulation< DH::dimension, DH::space_dimension > *tria, const int level, const int index, const AccessorData *local_data)
void copy_from(const DoFAccessor< structdim, DH, level_dof_access2 > &a)
void set_vertex_dof_index(const unsigned int vertex, const unsigned int i, const types::global_dof_index index, const unsigned int fe_index=DH::default_fe_index) const
typename::internal::DoFHandler::Iterators< DH, level_dof_access >::quad_iterator quad(const unsigned int i) const
Point< spacedim > & vertex(const unsigned int i) const
unsigned int n_active_fe_indices() const