Reference documentation for deal.II version GIT 6da2e5d553 2022-07-01 18:55:02+00:00
\(\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.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2020 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 
19 
20 #include <deal.II/fe/fe.h>
21 
23 #include <deal.II/grid/tria_iterator.templates.h>
24 
25 #include <vector>
26 
28 
29 /* --------------------- Static variables: DoFAccessor --------------------- */
30 
31 template <int structdim, int dim, int spacedim, bool level_dof_access>
32 const unsigned int
34 
35 template <int structdim, int dim, int spacedim, bool level_dof_access>
36 const unsigned int
38 
39 
40 
41 /* ---------------------- Functions: DoFInvalidAccessor -------------------- */
42 
43 template <int structdim, int dim, int spacedim>
46  const int,
47  const int,
48  const AccessorData *)
49 {
50  Assert(false,
51  ExcMessage("You are attempting an illegal conversion between "
52  "iterator/accessor types. The constructor you call "
53  "only exists to make certain template constructs "
54  "easier to write as dimension independent code but "
55  "the conversion is not valid in the current context."));
56 }
57 
58 
59 
60 template <int structdim, int dim, int spacedim>
62  const DoFInvalidAccessor &i)
63  : InvalidAccessor<structdim, dim, spacedim>(
64  static_cast<const InvalidAccessor<structdim, dim, spacedim> &>(i))
65 {
66  Assert(false,
67  ExcMessage("You are attempting an illegal conversion between "
68  "iterator/accessor types. The constructor you call "
69  "only exists to make certain template constructs "
70  "easier to write as dimension independent code but "
71  "the conversion is not valid in the current context."));
72 }
73 
74 
75 
76 template <int structdim, int dim, int spacedim>
79  const unsigned int,
80  const unsigned int) const
81 {
82  Assert(false, ExcInternalError());
83  return 0;
84 }
85 
86 
87 
88 template <int structdim, int dim, int spacedim>
89 void
91  const unsigned int,
93  const unsigned int) const
94 {
95  Assert(false, ExcInternalError());
96 }
97 
98 
99 
100 /*------------------------- Functions: DoFCellAccessor -----------------------*/
101 
102 
103 
104 template <int dim, int spacedim, bool lda>
105 void
107 {
108  Assert(static_cast<unsigned int>(this->present_level) <
109  this->dof_handler->object_dof_indices.size(),
110  ExcMessage("DoFHandler not initialized"));
111 
112  Assert(this->dof_handler != nullptr, typename BaseClass::ExcInvalidObject());
113 
114  internal::DoFCellAccessorImplementation::Implementation::
115  update_cell_dof_indices_cache(*this);
116 }
117 
118 
119 
120 template <int dim, int spacedim, bool lda>
121 void
123  const std::vector<types::global_dof_index> &local_dof_indices)
124 {
125  Assert(static_cast<unsigned int>(this->present_level) <
126  this->dof_handler->object_dof_indices.size(),
127  ExcMessage("DoFHandler not initialized"));
128 
129  Assert(this->dof_handler != nullptr, typename BaseClass::ExcInvalidObject());
130 
131  internal::DoFAccessorImplementation::Implementation::
132  template set_dof_indices<dim, spacedim, lda, dim>(*this,
133  local_dof_indices,
134  this->active_fe_index());
135 }
136 
137 
138 
139 template <int dim, int spacedim, bool lda>
142  const unsigned int face,
143  const unsigned int subface) const
144 {
148  this->dof_handler);
149 }
150 
151 
152 
153 template <int dim, int spacedim, bool lda>
156  const unsigned int face,
157  const unsigned int subface) const
158 {
161  subface);
163  this->dof_handler);
164 }
165 
166 
167 
168 template <int dim, int spacedim, bool lda>
171  const unsigned int face) const
172 {
176  this->dof_handler);
177 }
178 
179 
180 
181 template <int dim, int spacedim, bool lda>
184  const unsigned int face) const
185 {
189  this->dof_handler);
190 }
191 
192 
193 
194 template <int structdim, int dim, int spacedim, bool level_dof_access>
195 inline void
197  std::vector<types::global_dof_index> &dof_indices,
198  const unsigned int fe_index_) const
199 {
200  Assert(this->dof_handler != nullptr, ExcInvalidObject());
201 
202  const auto fe_index =
203  internal::DoFAccessorImplementation::get_fe_index_or_default(*this,
204  fe_index_);
205 
206  Assert(static_cast<unsigned int>(this->level()) <
207  this->dof_handler->object_dof_indices.size(),
208  ExcMessage(
209  "The DoFHandler to which this accessor points has not "
210  "been initialized, i.e., it doesn't appear that DoF indices "
211  "have been distributed on it."));
212 
213  // this function really only makes sense if either a) there are degrees of
214  // freedom defined on the present object, or b) the object is non-active
215  // objects but all degrees of freedom are located on vertices, since
216  // otherwise there are degrees of freedom on sub-objects which are not
217  // allocated for this non-active thing
218  Assert(this->fe_index_is_active(fe_index) ||
219  (this->dof_handler->get_fe(fe_index).n_dofs_per_cell() ==
220  this->n_vertices() *
221  this->dof_handler->get_fe(fe_index).n_dofs_per_vertex()),
222  ExcInternalError());
223 
224  // now do the actual work
225  ::internal::DoFAccessorImplementation::Implementation::get_dof_indices(
226  *this, dof_indices, fe_index);
227 }
228 
229 
230 
231 template <int structdim, int dim, int spacedim, bool level_dof_access>
232 inline void
234  const int level,
235  std::vector<types::global_dof_index> &dof_indices,
236  const unsigned int fe_index_) const
237 {
238  Assert(this->dof_handler != nullptr, ExcInvalidObject());
239  Assert(this->dof_handler->mg_vertex_dofs.size() > 0,
240  ExcMessage("Multigrid DoF indices can only be accessed after "
241  "DoFHandler::distribute_mg_dofs() has been called!"));
242 
243  const auto fe_index =
244  internal::DoFAccessorImplementation::get_fe_index_or_default(*this,
245  fe_index_);
246 
247  internal::DoFAccessorImplementation::Implementation::get_mg_dof_indices(
248  *this, level, dof_indices, fe_index);
249 }
250 
251 
252 
253 template <int structdim, int dim, int spacedim, bool level_dof_access>
254 inline void
256  const int level,
257  const std::vector<types::global_dof_index> &dof_indices,
258  const unsigned int fe_index_)
259 {
260  Assert(this->dof_handler != nullptr, ExcInvalidObject());
261 
262  const auto fe_index =
263  internal::DoFAccessorImplementation::get_fe_index_or_default(*this,
264  fe_index_);
265 
266  internal::DoFAccessorImplementation::Implementation::set_mg_dof_indices(
267  *this, level, dof_indices, fe_index);
268 }
269 
270 
271 // --------------------------------------------------------------------------
272 // explicit instantiations
273 #include "dof_accessor.inst"
274 
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor(const unsigned int i) const
TriaIterator< CellAccessor< dim, spacedim > > neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
TriaIterator< CellAccessor< dim, spacedim > > neighbor_or_periodic_neighbor(const unsigned int i) const
void set_mg_dof_indices(const int level, const std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index)
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index) const
void get_mg_dof_indices(const int level, std::vector< types::global_dof_index > &dof_indices, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index) const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > periodic_neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
void update_cell_dof_indices_cache() const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > neighbor_or_periodic_neighbor(const unsigned int i) const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > periodic_neighbor(const unsigned int i) const
void set_dof_indices(const std::vector< types::global_dof_index > &dof_indices)
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
DoFInvalidAccessor(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0)
Definition: dof_accessor.cc:44
types::global_dof_index dof_index(const unsigned int i, const unsigned int fe_index=DoFHandler< dim, spacedim >::default_fe_index) const
Definition: dof_accessor.cc:78
typename InvalidAccessor< structdim, dim, spacedim >::AccessorData AccessorData
void set_dof_index(const unsigned int i, const types::global_dof_index index, const unsigned int fe_index=DoFHandler< dim, spacedim >::invalid_fe_index) const
Definition: dof_accessor.cc:90
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
unsigned int level
Definition: grid_out.cc:4607
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int global_dof_index
Definition: types.h:76