Reference documentation for deal.II version 9.5.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\}}\)
Loading...
Searching...
No Matches
dof_info.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2006 - 2023 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
23
25
27
30
31#include <memory>
32
34
35namespace MeshWorker
36{
37 // Forward declaration
38#ifndef DOXYGEN
39 template <int dim, class DOFINFO>
40 class DoFInfoBox;
41#endif
42
74 template <int dim, int spacedim = dim, typename number = double>
75 class DoFInfo : public LocalResults<number>
76 {
77 public:
80
83
90 unsigned int face_number;
91
98 unsigned int sub_number;
99
104 std::vector<types::global_dof_index> indices;
105
110 std::vector<std::vector<types::global_dof_index>> indices_by_block;
111
116
121 DoFInfo(const DoFHandler<dim, spacedim> &dof_handler);
122
126 template <class DHCellIterator>
127 void
128 reinit(const DHCellIterator &c);
129
133 template <class DHCellIterator, class DHFaceIterator>
134 void
135 reinit(const DHCellIterator &c,
136 const DHFaceIterator &f,
137 const unsigned int face_no);
138
142 template <class DHCellIterator, class DHFaceIterator>
143 void
144 reinit(const DHCellIterator &c,
145 const DHFaceIterator &f,
146 const unsigned int face_no,
147 const unsigned int subface_no);
148
153 template <class DHFaceIterator>
154 void
155 set_face(const DHFaceIterator &f, const unsigned int face_no);
156
161 template <class DHFaceIterator>
162 void
163 set_subface(const DHFaceIterator &f,
164 const unsigned int face_no,
165 const unsigned int subface_no);
166
167 const BlockIndices &
168 local_indices() const;
169
170
173
178
179 private:
185 DoFInfo();
186
188 void
190
192 template <class DHCellIterator>
193 void
194 get_indices(const DHCellIterator &c);
195
197 std::vector<types::global_dof_index> indices_org;
198
205
206 friend class DoFInfoBox<dim, DoFInfo<dim, spacedim, number>>;
207 };
208
209
220 template <int dim, class DOFINFO>
222 {
223 public:
227 DoFInfoBox(const DOFINFO &seed);
228
234
238 DoFInfoBox &
240
244 void
245 reset();
246
252 template <class ASSEMBLER>
253 void
254 assemble(ASSEMBLER &ass) const;
255
256
260 DOFINFO cell;
269
275
281
286 };
287
288 //----------------------------------------------------------------------//
289
290 template <int dim, int spacedim, typename number>
292 : face_number(numbers::invalid_unsigned_int)
293 , sub_number(numbers::invalid_unsigned_int)
294 , level_cell(false)
295 {}
296
297
298
299 template <int dim, int spacedim, typename number>
301 const DoFHandler<dim, spacedim> &dof_handler)
302 : face_number(numbers::invalid_unsigned_int)
303 , sub_number(numbers::invalid_unsigned_int)
304 , level_cell(false)
305 {
306 std::vector<types::global_dof_index> aux(1);
307 aux[0] = dof_handler.get_fe().n_dofs_per_cell();
309 }
310
311
312 template <int dim, int spacedim, typename number>
313 template <class DHCellIterator>
314 inline void
316 {
317 indices.resize(c->get_fe().n_dofs_per_cell());
318 if (block_info == nullptr || block_info->local().size() == 0)
319 c->get_active_or_mg_dof_indices(indices);
320 else
321 {
322 indices_org.resize(c->get_fe().n_dofs_per_cell());
323 c->get_active_or_mg_dof_indices(indices_org);
324 set_block_indices();
325 }
326 }
327
328
329 template <int dim, int spacedim, typename number>
330 template <class DHCellIterator>
331 inline void
333 {
334 get_indices(c);
335 level_cell = c->is_level_cell();
336
338 face_number = numbers::invalid_unsigned_int;
340 if (block_info)
341 LocalResults<number>::reinit(block_info->local());
342 else
343 LocalResults<number>::reinit(aux_local_indices);
344 }
345
346
347 template <int dim, int spacedim, typename number>
348 template <class DHFaceIterator>
349 inline void
351 const unsigned int face_no)
352 {
353 face = static_cast<typename Triangulation<dim, spacedim>::face_iterator>(f);
354 face_number = face_no;
356 }
357
358
359 template <int dim, int spacedim, typename number>
360 template <class DHCellIterator, class DHFaceIterator>
361 inline void
363 const DHFaceIterator &f,
364 const unsigned int face_no)
365 {
366 if ((cell.state() != IteratorState::valid) ||
368 get_indices(c);
369 level_cell = c->is_level_cell();
370
372 set_face(f, face_no);
373
374 if (block_info)
375 LocalResults<number>::reinit(block_info->local());
376 else
377 LocalResults<number>::reinit(aux_local_indices);
378 }
379
380
381 template <int dim, int spacedim, typename number>
382 template <class DHFaceIterator>
383 inline void
385 const unsigned int face_no,
386 const unsigned int subface_no)
387 {
388 face = static_cast<typename Triangulation<dim, spacedim>::face_iterator>(f);
389 face_number = face_no;
390 sub_number = subface_no;
391 }
392
393
394 template <int dim, int spacedim, typename number>
395 template <class DHCellIterator, class DHFaceIterator>
396 inline void
398 const DHFaceIterator &f,
399 const unsigned int face_no,
400 const unsigned int subface_no)
401 {
402 if (cell.state() != IteratorState::valid ||
403 cell !=
404 static_cast<typename Triangulation<dim, spacedim>::cell_iterator>(c))
405 get_indices(c);
406 level_cell = c->is_level_cell();
407
408 cell = static_cast<typename Triangulation<dim, spacedim>::cell_iterator>(c);
409 set_subface(f, face_no, subface_no);
410
411 if (block_info)
412 LocalResults<number>::reinit(block_info->local());
413 else
414 LocalResults<number>::reinit(aux_local_indices);
415 }
416
417
418 template <int dim, int spacedim, typename number>
419 inline const BlockIndices &
421 {
422 if (block_info)
423 return block_info->local();
424 return aux_local_indices;
425 }
426
427 //----------------------------------------------------------------------//
428
429 template <int dim, class DOFINFO>
430 inline DoFInfoBox<dim, DOFINFO>::DoFInfoBox(const DOFINFO &seed)
431 : cell(seed)
432 , cell_valid(true)
433 {
434 for (const unsigned int i : GeometryInfo<dim>::face_indices())
435 {
436 exterior[i] = seed;
437 interior[i] = seed;
438 interior_face_available[i] = false;
439 exterior_face_available[i] = false;
440 }
441 }
442
443
444 template <int dim, class DOFINFO>
446 const DoFInfoBox<dim, DOFINFO> &other)
447 : cell(other.cell)
448 , cell_valid(other.cell_valid)
449 {
450 for (const unsigned int i : GeometryInfo<dim>::face_indices())
451 {
452 exterior[i] = other.exterior[i];
453 interior[i] = other.interior[i];
454 interior_face_available[i] = false;
455 exterior_face_available[i] = false;
456 }
457 }
458
459
460 template <int dim, class DOFINFO>
463 {
464 cell = other.cell;
465 cell_valid = other.cell_valid;
466 for (const unsigned int i : GeometryInfo<dim>::face_indices())
467 {
468 exterior[i] = other.exterior[i];
469 interior[i] = other.interior[i];
470 interior_face_available[i] = false;
471 exterior_face_available[i] = false;
472 }
473 return *this;
474 }
475
476
477 template <int dim, class DOFINFO>
478 inline void
480 {
481 cell_valid = false;
482 for (const unsigned int i : GeometryInfo<dim>::face_indices())
483 {
484 interior_face_available[i] = false;
485 exterior_face_available[i] = false;
486 }
487 }
488
489
490 template <int dim, class DOFINFO>
491 template <class ASSEMBLER>
492 inline void
493 DoFInfoBox<dim, DOFINFO>::assemble(ASSEMBLER &assembler) const
494 {
495 if (!cell_valid)
496 return;
497
498 assembler.assemble(cell);
499 for (const unsigned int i : GeometryInfo<dim>::face_indices())
500 {
501 // Only do something if data available
502 if (interior_face_available[i])
503 {
504 // If both data
505 // available, it is an
506 // interior face
507 if (exterior_face_available[i])
508 assembler.assemble(interior[i], exterior[i]);
509 else
510 assembler.assemble(interior[i]);
511 }
512 }
513 }
514} // namespace MeshWorker
515
517
518#endif
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition block_info.h:96
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
unsigned int n_dofs_per_cell() const
DOFINFO interior[GeometryInfo< dim >::faces_per_cell]
Definition dof_info.h:264
void assemble(ASSEMBLER &ass) const
Definition dof_info.h:493
bool exterior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition dof_info.h:280
bool interior_face_available[GeometryInfo< dim >::faces_per_cell]
Definition dof_info.h:274
DOFINFO exterior[GeometryInfo< dim >::faces_per_cell]
Definition dof_info.h:268
DoFInfoBox & operator=(const DoFInfoBox< dim, DOFINFO > &)
Definition dof_info.h:462
DoFInfoBox(const DOFINFO &seed)
Definition dof_info.h:430
unsigned int sub_number
Definition dof_info.h:98
void set_face(const DHFaceIterator &f, const unsigned int face_no)
Definition dof_info.h:350
std::vector< types::global_dof_index > indices_org
Auxiliary vector.
Definition dof_info.h:197
BlockIndices aux_local_indices
Definition dof_info.h:204
DoFInfo(const BlockInfo &block_info)
std::vector< types::global_dof_index > indices
Definition dof_info.h:104
void set_subface(const DHFaceIterator &f, const unsigned int face_no, const unsigned int subface_no)
Definition dof_info.h:384
void set_block_indices()
Set up local block indices.
const BlockIndices & local_indices() const
Definition dof_info.h:420
void get_indices(const DHCellIterator &c)
Fill index vector with active indices.
Definition dof_info.h:315
Triangulation< dim, spacedim >::face_iterator face
The current face.
Definition dof_info.h:82
std::vector< std::vector< types::global_dof_index > > indices_by_block
Definition dof_info.h:110
SmartPointer< const BlockInfo, DoFInfo< dim, spacedim > > block_info
The block structure of the system.
Definition dof_info.h:172
unsigned int face_number
Definition dof_info.h:90
void reinit(const DHCellIterator &c)
Definition dof_info.h:332
Triangulation< dim, spacedim >::cell_iterator cell
The current cell.
Definition dof_info.h:79
void reinit(const BlockIndices &local_sizes)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator
Definition tria.h:1370
@ valid
Iterator points to a valid object.
static const unsigned int invalid_unsigned_int
Definition types.h:213
static std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices()