Reference documentation for deal.II version 9.1.1
\(\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  template <int dim, class DOFINFO>
38  class DoFInfoBox;
39 
40 
73  template <int dim, int spacedim = dim, typename number = double>
74  class DoFInfo : public LocalResults<number>
75  {
76  public:
79 
82 
89  unsigned int face_number;
90 
97  unsigned int sub_number;
98 
103  std::vector<types::global_dof_index> indices;
104 
109  std::vector<std::vector<types::global_dof_index>> indices_by_block;
110 
114  DoFInfo(const BlockInfo &block_info);
115 
120  DoFInfo(const DoFHandler<dim, spacedim> &dof_handler);
121 
125  template <class DHCellIterator>
126  void
127  reinit(const DHCellIterator &c);
128 
132  template <class DHCellIterator, class DHFaceIterator>
133  void
134  reinit(const DHCellIterator &c,
135  const DHFaceIterator &f,
136  const unsigned int face_no);
137 
141  template <class DHCellIterator, class DHFaceIterator>
142  void
143  reinit(const DHCellIterator &c,
144  const DHFaceIterator &f,
145  const unsigned int face_no,
146  const unsigned int subface_no);
147 
152  template <class DHFaceIterator>
153  void
154  set_face(const DHFaceIterator &f, const unsigned int face_no);
155 
160  template <class DHFaceIterator>
161  void
162  set_subface(const DHFaceIterator &f,
163  const unsigned int face_no,
164  const unsigned int subface_no);
165 
166  const BlockIndices &
167  local_indices() const;
168 
169 
172 
177 
178  private:
184  DoFInfo();
185 
187  void
189 
191  template <class DHCellIterator>
192  void
193  get_indices(const DHCellIterator &c);
194 
196  std::vector<types::global_dof_index> indices_org;
197 
204 
205  friend class DoFInfoBox<dim, DoFInfo<dim, spacedim, number>>;
206  };
207 
208 
220  template <int dim, class DOFINFO>
221  class DoFInfoBox
222  {
223  public:
227  DoFInfoBox(const DOFINFO &seed);
228 
234 
238  void
239  reset();
240 
246  template <class ASSEMBLER>
247  void
248  assemble(ASSEMBLER &ass) const;
249 
250 
254  DOFINFO cell;
263 
269 
275 
280  };
281 
282  //----------------------------------------------------------------------//
283 
284  template <int dim, int spacedim, typename number>
286  : face_number(numbers::invalid_unsigned_int)
287  , sub_number(numbers::invalid_unsigned_int)
288  , level_cell(false)
289  {}
290 
291 
292 
293  template <int dim, int spacedim, typename number>
295  const DoFHandler<dim, spacedim> &dof_handler)
296  : face_number(numbers::invalid_unsigned_int)
297  , sub_number(numbers::invalid_unsigned_int)
298  , level_cell(false)
299  {
300  std::vector<types::global_dof_index> aux(1);
301  aux[0] = dof_handler.get_fe().dofs_per_cell;
303  }
304 
305 
306  template <int dim, int spacedim, typename number>
307  template <class DHCellIterator>
308  inline void
310  {
311  indices.resize(c->get_fe().dofs_per_cell);
312  if (block_info == nullptr || block_info->local().size() == 0)
313  c->get_active_or_mg_dof_indices(indices);
314  else
315  {
316  indices_org.resize(c->get_fe().dofs_per_cell);
317  c->get_active_or_mg_dof_indices(indices_org);
318  set_block_indices();
319  }
320  }
321 
322 
323  template <int dim, int spacedim, typename number>
324  template <class DHCellIterator>
325  inline void
326  DoFInfo<dim, spacedim, number>::reinit(const DHCellIterator &c)
327  {
328  get_indices(c);
329  level_cell = c->is_level_cell();
330 
332  face_number = numbers::invalid_unsigned_int;
333  sub_number = numbers::invalid_unsigned_int;
334  if (block_info)
335  LocalResults<number>::reinit(block_info->local());
336  else
337  LocalResults<number>::reinit(aux_local_indices);
338  }
339 
340 
341  template <int dim, int spacedim, typename number>
342  template <class DHFaceIterator>
343  inline void
345  const unsigned int face_no)
346  {
347  face = static_cast<typename Triangulation<dim, spacedim>::face_iterator>(f);
348  face_number = face_no;
349  sub_number = numbers::invalid_unsigned_int;
350  }
351 
352 
353  template <int dim, int spacedim, typename number>
354  template <class DHCellIterator, class DHFaceIterator>
355  inline void
356  DoFInfo<dim, spacedim, number>::reinit(const DHCellIterator &c,
357  const DHFaceIterator &f,
358  const unsigned int face_no)
359  {
360  if ((cell.state() != IteratorState::valid) ||
361  cell != typename Triangulation<dim, spacedim>::cell_iterator(*c))
362  get_indices(c);
363  level_cell = c->is_level_cell();
364 
366  set_face(f, face_no);
367 
368  if (block_info)
369  LocalResults<number>::reinit(block_info->local());
370  else
371  LocalResults<number>::reinit(aux_local_indices);
372  }
373 
374 
375  template <int dim, int spacedim, typename number>
376  template <class DHFaceIterator>
377  inline void
379  const unsigned int face_no,
380  const unsigned int subface_no)
381  {
382  face = static_cast<typename Triangulation<dim, spacedim>::face_iterator>(f);
383  face_number = face_no;
384  sub_number = subface_no;
385  }
386 
387 
388  template <int dim, int spacedim, typename number>
389  template <class DHCellIterator, class DHFaceIterator>
390  inline void
391  DoFInfo<dim, spacedim, number>::reinit(const DHCellIterator &c,
392  const DHFaceIterator &f,
393  const unsigned int face_no,
394  const unsigned int subface_no)
395  {
396  if (cell.state() != IteratorState::valid ||
397  cell !=
398  static_cast<typename Triangulation<dim, spacedim>::cell_iterator>(c))
399  get_indices(c);
400  level_cell = c->is_level_cell();
401 
402  cell = static_cast<typename Triangulation<dim, spacedim>::cell_iterator>(c);
403  set_subface(f, face_no, subface_no);
404 
405  if (block_info)
406  LocalResults<number>::reinit(block_info->local());
407  else
408  LocalResults<number>::reinit(aux_local_indices);
409  }
410 
411 
412  template <int dim, int spacedim, typename number>
413  inline const BlockIndices &
415  {
416  if (block_info)
417  return block_info->local();
418  return aux_local_indices;
419  }
420 
421  //----------------------------------------------------------------------//
422 
423  template <int dim, class DOFINFO>
424  inline DoFInfoBox<dim, DOFINFO>::DoFInfoBox(const DOFINFO &seed)
425  : cell(seed)
426  , cell_valid(true)
427  {
428  for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
429  {
430  exterior[i] = seed;
431  interior[i] = seed;
432  interior_face_available[i] = false;
433  exterior_face_available[i] = false;
434  }
435  }
436 
437 
438  template <int dim, class DOFINFO>
440  const DoFInfoBox<dim, DOFINFO> &other)
441  : cell(other.cell)
442  , cell_valid(other.cell_valid)
443  {
444  for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
445  {
446  exterior[i] = other.exterior[i];
447  interior[i] = other.interior[i];
448  interior_face_available[i] = false;
449  exterior_face_available[i] = false;
450  }
451  }
452 
453 
454  template <int dim, class DOFINFO>
455  inline void
457  {
458  cell_valid = false;
459  for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
460  {
461  interior_face_available[i] = false;
462  exterior_face_available[i] = false;
463  }
464  }
465 
466 
467  template <int dim, class DOFINFO>
468  template <class ASSEMBLER>
469  inline void
470  DoFInfoBox<dim, DOFINFO>::assemble(ASSEMBLER &assembler) const
471  {
472  if (!cell_valid)
473  return;
474 
475  assembler.assemble(cell);
476  for (unsigned int i = 0; i < GeometryInfo<dim>::faces_per_cell; ++i)
477  {
478  // Only do something if data available
479  if (interior_face_available[i])
480  {
481  // If both data
482  // available, it is an
483  // interior face
484  if (exterior_face_available[i])
485  assembler.assemble(interior[i], exterior[i]);
486  else
487  assembler.assemble(interior[i]);
488  }
489  }
490  }
491 } // namespace MeshWorker
492 
493 DEAL_II_NAMESPACE_CLOSE
494 
495 #endif
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
Definition: dof_info.h:171
void get_indices(const DHCellIterator &c)
Fill index vector with active indices.
Definition: dof_info.h:309
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
Definition: dof_info.h:78
static const unsigned int invalid_unsigned_int
Definition: types.h:173
std::vector< types::global_dof_index > indices_org
Auxiliary vector.
Definition: dof_info.h:196
std::vector< types::global_dof_index > indices
Definition: dof_info.h:103
void set_block_indices()
Set up local block indices.
unsigned int sub_number
Definition: dof_info.h:97
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:274
BlockIndices aux_local_indices
Definition: dof_info.h:203
void set_subface(const DHFaceIterator &f, const unsigned int face_no, const unsigned int subface_no)
Definition: dof_info.h:378
bool interior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:268
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
DoFInfoBox(const DOFINFO &seed)
Definition: dof_info.h:424
void assemble(ASSEMBLER &ass) const
Definition: dof_info.h:470
std::vector< std::vector< types::global_dof_index > > indices_by_block
Definition: dof_info.h:109
void reinit(const BlockIndices &local_sizes)
Definition: mesh_worker.cc:28
DOFINFO interior[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:258
unsigned int face_number
Definition: dof_info.h:89
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition: block_info.h:96
Iterator points to a valid object.
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator
Definition: tria.h:1511
DOFINFO exterior[GeometryInfo< dim >::faces_per_cell]
Definition: dof_info.h:262
Triangulation< dim, spacedim >::face_iterator face
The current face.
Definition: dof_info.h:81
void reinit(const DHCellIterator &c)
Definition: dof_info.h:326
void set_face(const DHFaceIterator &f, const unsigned int face_no)
Definition: dof_info.h:344