Reference documentation for deal.II version Git 7026f387cc 2019-10-15 14:19:01 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
dof_info.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2006 - 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 
17 #ifndef dealii_mesh_worker_dof_info_h
18 #define dealii_mesh_worker_dof_info_h
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/quadrature_lib.h>
23 
24 #include <deal.II/dofs/block_info.h>
25 
26 #include <deal.II/fe/fe_values.h>
27 
28 #include <deal.II/meshworker/local_results.h>
29 #include <deal.II/meshworker/vector_selector.h>
30 
31 #include <memory>
32 
33 DEAL_II_NAMESPACE_OPEN
34 
35 namespace MeshWorker
36 {
37  // Forward declaration
38 #ifndef DOXYGEN
39  template <int dim, class DOFINFO>
40  class DoFInfoBox;
41 #endif
42 
75  template <int dim, int spacedim = dim, typename number = double>
76  class DoFInfo : public LocalResults<number>
77  {
78  public:
81 
84 
91  unsigned int face_number;
92 
99  unsigned int sub_number;
100 
105  std::vector<types::global_dof_index> indices;
106 
111  std::vector<std::vector<types::global_dof_index>> indices_by_block;
112 
116  DoFInfo(const BlockInfo &block_info);
117 
122  DoFInfo(const DoFHandler<dim, spacedim> &dof_handler);
123 
127  template <class DHCellIterator>
128  void
129  reinit(const DHCellIterator &c);
130 
134  template <class DHCellIterator, class DHFaceIterator>
135  void
136  reinit(const DHCellIterator &c,
137  const DHFaceIterator &f,
138  const unsigned int face_no);
139 
143  template <class DHCellIterator, class DHFaceIterator>
144  void
145  reinit(const DHCellIterator &c,
146  const DHFaceIterator &f,
147  const unsigned int face_no,
148  const unsigned int subface_no);
149 
154  template <class DHFaceIterator>
155  void
156  set_face(const DHFaceIterator &f, const unsigned int face_no);
157 
162  template <class DHFaceIterator>
163  void
164  set_subface(const DHFaceIterator &f,
165  const unsigned int face_no,
166  const unsigned int subface_no);
167 
168  const BlockIndices &
169  local_indices() const;
170 
171 
174 
179 
180  private:
186  DoFInfo();
187 
189  void
191 
193  template <class DHCellIterator>
194  void
195  get_indices(const DHCellIterator &c);
196 
198  std::vector<types::global_dof_index> indices_org;
199 
206 
207  friend class DoFInfoBox<dim, DoFInfo<dim, spacedim, number>>;
208  };
209 
210 
222  template <int dim, class DOFINFO>
224  {
225  public:
229  DoFInfoBox(const DOFINFO &seed);
230 
236 
240  void
241  reset();
242 
248  template <class ASSEMBLER>
249  void
250  assemble(ASSEMBLER &ass) const;
251 
252 
256  DOFINFO cell;
265 
270  bool interior_face_available[GeometryInfo<dim>::faces_per_cell];
271 
276  bool exterior_face_available[GeometryInfo<dim>::faces_per_cell];
277 
282  };
283 
284  //----------------------------------------------------------------------//
285 
286  template <int dim, int spacedim, typename number>
288  : face_number(numbers::invalid_unsigned_int)
289  , sub_number(numbers::invalid_unsigned_int)
290  , level_cell(false)
291  {}
292 
293 
294 
295  template <int dim, int spacedim, typename number>
297  const DoFHandler<dim, spacedim> &dof_handler)
298  : face_number(numbers::invalid_unsigned_int)
299  , sub_number(numbers::invalid_unsigned_int)
300  , level_cell(false)
301  {
302  std::vector<types::global_dof_index> aux(1);
303  aux[0] = dof_handler.get_fe().dofs_per_cell;
305  }
306 
307 
308  template <int dim, int spacedim, typename number>
309  template <class DHCellIterator>
310  inline void
312  {
313  indices.resize(c->get_fe().dofs_per_cell);
314  if (block_info == nullptr || block_info->local().size() == 0)
315  c->get_active_or_mg_dof_indices(indices);
316  else
317  {
318  indices_org.resize(c->get_fe().dofs_per_cell);
319  c->get_active_or_mg_dof_indices(indices_org);
321  }
322  }
323 
324 
325  template <int dim, int spacedim, typename number>
326  template <class DHCellIterator>
327  inline void
328  DoFInfo<dim, spacedim, number>::reinit(const DHCellIterator &c)
329  {
330  get_indices(c);
331  level_cell = c->is_level_cell();
332 
336  if (block_info)
338  else
340  }
341 
342 
343  template <int dim, int spacedim, typename number>
344  template <class DHFaceIterator>
345  inline void
347  const unsigned int face_no)
348  {
349  face = static_cast<typename Triangulation<dim, spacedim>::face_iterator>(f);
350  face_number = face_no;
352  }
353 
354 
355  template <int dim, int spacedim, typename number>
356  template <class DHCellIterator, class DHFaceIterator>
357  inline void
358  DoFInfo<dim, spacedim, number>::reinit(const DHCellIterator &c,
359  const DHFaceIterator &f,
360  const unsigned int face_no)
361  {
362  if ((cell.state() != IteratorState::valid) ||
364  get_indices(c);
365  level_cell = c->is_level_cell();
366 
368  set_face(f, face_no);
369 
370  if (block_info)
372  else
374  }
375 
376 
377  template <int dim, int spacedim, typename number>
378  template <class DHFaceIterator>
379  inline void
381  const unsigned int face_no,
382  const unsigned int subface_no)
383  {
384  face = static_cast<typename Triangulation<dim, spacedim>::face_iterator>(f);
385  face_number = face_no;
386  sub_number = subface_no;
387  }
388 
389 
390  template <int dim, int spacedim, typename number>
391  template <class DHCellIterator, class DHFaceIterator>
392  inline void
393  DoFInfo<dim, spacedim, number>::reinit(const DHCellIterator &c,
394  const DHFaceIterator &f,
395  const unsigned int face_no,
396  const unsigned int subface_no)
397  {
398  if (cell.state() != IteratorState::valid ||
399  cell !=
400  static_cast<typename Triangulation<dim, spacedim>::cell_iterator>(c))
401  get_indices(c);
402  level_cell = c->is_level_cell();
403 
404  cell = static_cast<typename Triangulation<dim, spacedim>::cell_iterator>(c);
405  set_subface(f, face_no, subface_no);
406 
407  if (block_info)
409  else
411  }
412 
413 
414  template <int dim, int spacedim, typename number>
415  inline const BlockIndices &
417  {
418  if (block_info)
419  return block_info->local();
420  return aux_local_indices;
421  }
422 
423  //----------------------------------------------------------------------//
424 
425  template <int dim, class DOFINFO>
426  inline DoFInfoBox<dim, DOFINFO>::DoFInfoBox(const DOFINFO &seed)
427  : cell(seed)
428  , cell_valid(true)
429  {
430  for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
431  {
432  exterior[i] = seed;
433  interior[i] = seed;
434  interior_face_available[i] = false;
435  exterior_face_available[i] = false;
436  }
437  }
438 
439 
440  template <int dim, class DOFINFO>
442  const DoFInfoBox<dim, DOFINFO> &other)
443  : cell(other.cell)
444  , cell_valid(other.cell_valid)
445  {
446  for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
447  {
448  exterior[i] = other.exterior[i];
449  interior[i] = other.interior[i];
450  interior_face_available[i] = false;
451  exterior_face_available[i] = false;
452  }
453  }
454 
455 
456  template <int dim, class DOFINFO>
457  inline void
459  {
460  cell_valid = false;
461  for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
462  {
463  interior_face_available[i] = false;
464  exterior_face_available[i] = false;
465  }
466  }
467 
468 
469  template <int dim, class DOFINFO>
470  template <class ASSEMBLER>
471  inline void
472  DoFInfoBox<dim, DOFINFO>::assemble(ASSEMBLER &assembler) const
473  {
474  if (!cell_valid)
475  return;
476 
477  assembler.assemble(cell);
478  for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
479  {
480  // Only do something if data available
482  {
483  // If both data
484  // available, it is an
485  // interior face
487  assembler.assemble(interior[i], exterior[i]);
488  else
489  assembler.assemble(interior[i]);
490  }
491  }
492  }
493 } // namespace MeshWorker
494 
495 DEAL_II_NAMESPACE_CLOSE
496 
497 #endif
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
Definition: dof_info.h:173
void get_indices(const DHCellIterator &c)
Fill index vector with active indices.
Definition: dof_info.h:311
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
Definition: dof_info.h:80
static const unsigned int invalid_unsigned_int
Definition: types.h:187
std::vector< types::global_dof_index > indices_org
Auxiliary vector.
Definition: dof_info.h:198
std::vector< types::global_dof_index > indices
Definition: dof_info.h:105
void set_block_indices()
Set up local block indices.
unsigned int sub_number
Definition: dof_info.h:99
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
bool exterior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:276
BlockIndices aux_local_indices
Definition: dof_info.h:205
void set_subface(const DHFaceIterator &f, const unsigned int face_no, const unsigned int subface_no)
Definition: dof_info.h:380
bool interior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:270
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
DoFInfoBox(const DOFINFO &seed)
Definition: dof_info.h:426
void assemble(ASSEMBLER &ass) const
Definition: dof_info.h:472
std::vector< std::vector< types::global_dof_index > > indices_by_block
Definition: dof_info.h:111
void reinit(const BlockIndices &local_sizes)
Definition: mesh_worker.cc:28
DOFINFO interior[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:260
unsigned int face_number
Definition: dof_info.h:91
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition: block_info.h:97
Iterator points to a valid object.
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator
Definition: tria.h:1529
DOFINFO exterior[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:264
Triangulation< dim, spacedim >::face_iterator face
The current face.
Definition: dof_info.h:83
void reinit(const DHCellIterator &c)
Definition: dof_info.h:328
void set_face(const DHFaceIterator &f, const unsigned int face_no)
Definition: dof_info.h:346