Reference documentation for deal.II version 9.3.0
\(\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_objects.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2006 - 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_objects_h
17 #define dealii_dof_objects_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 
24 #include <vector>
25 
27 
28 // Forward declarations
29 #ifndef DOXYGEN
30 template <int, int>
31 class DoFHandler;
32 #endif
33 
34 namespace internal
35 {
36  namespace DoFHandlerImplementation
37  {
38 #ifndef DOXYGEN
39  template <int>
40  class DoFLevel;
41  template <int>
42  class DoFFaces;
43 #endif
44 
68  template <int dim>
69  class DoFObjects
70  {
71  public:
75  std::vector<types::global_dof_index> dofs;
76 
77  public:
91  template <int dh_dim, int spacedim>
92  void
93  set_dof_index(const ::DoFHandler<dh_dim, spacedim> &dof_handler,
94  const unsigned int obj_index,
95  const unsigned int fe_index,
96  const unsigned int local_index,
98 
111  template <int dh_dim, int spacedim>
113  get_dof_index(const ::DoFHandler<dh_dim, spacedim> &dof_handler,
114  const unsigned int obj_index,
115  const unsigned int fe_index,
116  const unsigned int local_index) const;
117 
123  template <int dh_dim, int spacedim>
124  unsigned int
126  const ::DoFHandler<dh_dim, spacedim> &dof_handler,
127  const types::global_dof_index index) const;
128 
133  template <int dh_dim, int spacedim>
134  bool
136  const ::DoFHandler<dh_dim, spacedim> &dof_handler,
137  const types::global_dof_index index,
138  const unsigned int fe_index) const;
139 
144  std::size_t
145  memory_consumption() const;
146 
152  template <class Archive>
153  void
154  serialize(Archive &ar, const unsigned int version);
155 
156  // Declare the classes that store levels and faces of DoFs friends so
157  // that they can resize arrays.
158  template <int>
159  friend class DoFLevel;
160  template <int>
161  friend class DoFFaces;
162  };
163 
164 
165  // --------------------- template and inline functions ------------------
166 
167  template <int dim>
168  template <int dh_dim, int spacedim>
169  inline unsigned int
171  const ::DoFHandler<dh_dim, spacedim> &,
172  const types::global_dof_index) const
173  {
174  return 1;
175  }
176 
177 
178 
179  template <int dim>
180  template <int dh_dim, int spacedim>
181  inline bool
183  const ::DoFHandler<dh_dim, spacedim> &,
185  const unsigned int fe_index) const
186  {
187  (void)fe_index;
188  Assert((fe_index ==
190  ExcMessage("Only zero fe_index values are allowed for "
191  "non-hp-DoFHandlers."));
192  return true;
193  }
194 
195 
196 
197  template <int dim>
198  template <int dh_dim, int spacedim>
201  const ::DoFHandler<dh_dim, spacedim> &dof_handler,
202  const unsigned int obj_index,
203  const unsigned int fe_index,
204  const unsigned int local_index) const
205  {
206  (void)fe_index;
207  Assert(
209  ExcMessage(
210  "Only the default FE index is allowed for non-hp-DoFHandler objects"));
211  Assert(
212  local_index < dof_handler.get_fe().template n_dofs_per_object<dim>(),
213  ExcIndexRange(local_index,
214  0,
215  dof_handler.get_fe().template n_dofs_per_object<dim>()));
216  Assert(obj_index *
217  dof_handler.get_fe().template n_dofs_per_object<dim>() +
218  local_index <
219  dofs.size(),
220  ExcInternalError());
221 
222  return dofs[obj_index *
223  dof_handler.get_fe().template n_dofs_per_object<dim>() +
224  local_index];
225  }
226 
227 
228  template <int dim>
229  template <class Archive>
230  void
231  DoFObjects<dim>::serialize(Archive &ar, const unsigned int)
232  {
233  ar &dofs;
234  }
235 
236  } // namespace DoFHandlerImplementation
237 } // namespace internal
238 
240 
241 #endif
void serialize(Archive &ar, const unsigned int version)
Definition: dof_objects.h:231
types::global_dof_index get_dof_index(const ::DoFHandler< dh_dim, spacedim > &dof_handler, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index) const
Definition: dof_objects.h:200
void set_dof_index(const ::DoFHandler< dh_dim, spacedim > &dof_handler, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index)
Definition: dof_objects.cc:42
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
static ::ExceptionBase & ExcMessage(std::string arg1)
bool fe_index_is_active(const ::DoFHandler< dh_dim, spacedim > &dof_handler, const types::global_dof_index index, const unsigned int fe_index) const
Definition: dof_objects.h:182
#define Assert(cond, exc)
Definition: exceptions.h:1465
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:395
unsigned int n_active_fe_indices(const ::DoFHandler< dh_dim, spacedim > &dof_handler, const types::global_dof_index index) const
Definition: dof_objects.h:170
std::vector< types::global_dof_index > dofs
Definition: dof_objects.h:75
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:394
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
static ::ExceptionBase & ExcInternalError()