Reference documentation for deal.II version 9.3.3
\(\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_handler.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 1998 - 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#include <deal.II/base/config.h>
17
18#include <deal.II/base/geometry_info.h>
19#include <deal.II/base/memory_consumption.h>
20#include <deal.II/base/mpi.templates.h>
21
22#include <deal.II/distributed/cell_data_transfer.templates.h>
23#include <deal.II/distributed/fully_distributed_tria.h>
24#include <deal.II/distributed/shared_tria.h>
25#include <deal.II/distributed/tria.h>
26
27#include <deal.II/dofs/dof_handler.h>
28#include <deal.II/dofs/dof_handler_policy.h>
29
30#include <deal.II/grid/grid_tools.h>
31#include <deal.II/grid/tria.h>
32#include <deal.II/grid/tria_accessor.h>
33#include <deal.II/grid/tria_iterator.h>
34#include <deal.II/grid/tria_levels.h>
35
36#include <algorithm>
37#include <memory>
38#include <set>
39#include <unordered_set>
40
41DEAL_II_NAMESPACE_OPEN
42
43template <int dim, int spacedim>
44const typename DoFHandler<dim, spacedim>::active_fe_index_type
45 DoFHandler<dim, spacedim>::invalid_active_fe_index;
46
47namespace internal
48{
49 template <int dim, int spacedim>
50 std::string
51 policy_to_string(const ::internal::DoFHandlerImplementation::Policy::
52 PolicyBase<dim, spacedim> &policy)
53 {
54 std::string policy_name;
55 if (dynamic_cast<const typename ::internal::DoFHandlerImplementation::
56 Policy::Sequential<dim, spacedim> *>(&policy) ||
57 dynamic_cast<const typename ::internal::DoFHandlerImplementation::
58 Policy::Sequential<dim, spacedim> *>(&policy))
59 policy_name = "Policy::Sequential<";
60 else if (dynamic_cast<
61 const typename ::internal::DoFHandlerImplementation::
62 Policy::ParallelDistributed<dim, spacedim> *>(&policy) ||
63 dynamic_cast<
64 const typename ::internal::DoFHandlerImplementation::
65 Policy::ParallelDistributed<dim, spacedim> *>(&policy))
66 policy_name = "Policy::ParallelDistributed<";
67 else if (dynamic_cast<
68 const typename ::internal::DoFHandlerImplementation::
69 Policy::ParallelShared<dim, spacedim> *>(&policy) ||
70 dynamic_cast<
71 const typename ::internal::DoFHandlerImplementation::
72 Policy::ParallelShared<dim, spacedim> *>(&policy))
73 policy_name = "Policy::ParallelShared<";
74 else
75 AssertThrow(false, ExcNotImplemented());
76 policy_name += Utilities::int_to_string(dim) + "," +
77 Utilities::int_to_string(spacedim) + ">";
78 return policy_name;
79 }
80
81
82 namespace DoFHandlerImplementation
83 {
88 struct Implementation
89 {
94 template <int spacedim>
95 static unsigned int
96 max_couplings_between_dofs(const DoFHandler<1, spacedim> &dof_handler)
97 {
98 return std::min(static_cast<types::global_dof_index>(
99 3 * dof_handler.fe_collection.max_dofs_per_vertex() +
100 2 * dof_handler.fe_collection.max_dofs_per_line()),
101 dof_handler.n_dofs());
102 }
103
104 template <int spacedim>
105 static unsigned int
106 max_couplings_between_dofs(const DoFHandler<2, spacedim> &dof_handler)
107 {
108 // get these numbers by drawing pictures
109 // and counting...
110 // example:
111 // | | |
112 // --x-----x--x--X--
113 // | | | |
114 // | x--x--x
115 // | | | |
116 // --x--x--*--x--x--
117 // | | | |
118 // x--x--x |
119 // | | | |
120 // --X--x--x-----x--
121 // | | |
122 // x = vertices connected with center vertex *;
123 // = total of 19
124 // (the X vertices are connected with * if
125 // the vertices adjacent to X are hanging
126 // nodes)
127 // count lines -> 28 (don't forget to count
128 // mother and children separately!)
129 types::global_dof_index max_couplings;
130 switch (dof_handler.tria->max_adjacent_cells())
131 {
132 case 4:
133 max_couplings =
134 19 * dof_handler.fe_collection.max_dofs_per_vertex() +
135 28 * dof_handler.fe_collection.max_dofs_per_line() +
136 8 * dof_handler.fe_collection.max_dofs_per_quad();
137 break;
138 case 5:
139 max_couplings =
140 21 * dof_handler.fe_collection.max_dofs_per_vertex() +
141 31 * dof_handler.fe_collection.max_dofs_per_line() +
142 9 * dof_handler.fe_collection.max_dofs_per_quad();
143 break;
144 case 6:
145 max_couplings =
146 28 * dof_handler.fe_collection.max_dofs_per_vertex() +
147 42 * dof_handler.fe_collection.max_dofs_per_line() +
148 12 * dof_handler.fe_collection.max_dofs_per_quad();
149 break;
150 case 7:
151 max_couplings =
152 30 * dof_handler.fe_collection.max_dofs_per_vertex() +
153 45 * dof_handler.fe_collection.max_dofs_per_line() +
154 13 * dof_handler.fe_collection.max_dofs_per_quad();
155 break;
156 case 8:
157 max_couplings =
158 37 * dof_handler.fe_collection.max_dofs_per_vertex() +
159 56 * dof_handler.fe_collection.max_dofs_per_line() +
160 16 * dof_handler.fe_collection.max_dofs_per_quad();
161 break;
162
163 // the following numbers are not based on actual counting but by
164 // extrapolating the number sequences from the previous ones (for
165 // example, for n_dofs_per_vertex(), the sequence above is 19, 21,
166 // 28, 30, 37, and is continued as follows):
167 case 9:
168 max_couplings =
169 39 * dof_handler.fe_collection.max_dofs_per_vertex() +
170 59 * dof_handler.fe_collection.max_dofs_per_line() +
171 17 * dof_handler.fe_collection.max_dofs_per_quad();
172 break;
173 case 10:
174 max_couplings =
175 46 * dof_handler.fe_collection.max_dofs_per_vertex() +
176 70 * dof_handler.fe_collection.max_dofs_per_line() +
177 20 * dof_handler.fe_collection.max_dofs_per_quad();
178 break;
179 case 11:
180 max_couplings =
181 48 * dof_handler.fe_collection.max_dofs_per_vertex() +
182 73 * dof_handler.fe_collection.max_dofs_per_line() +
183 21 * dof_handler.fe_collection.max_dofs_per_quad();
184 break;
185 case 12:
186 max_couplings =
187 55 * dof_handler.fe_collection.max_dofs_per_vertex() +
188 84 * dof_handler.fe_collection.max_dofs_per_line() +
189 24 * dof_handler.fe_collection.max_dofs_per_quad();
190 break;
191 case 13:
192 max_couplings =
193 57 * dof_handler.fe_collection.max_dofs_per_vertex() +
194 87 * dof_handler.fe_collection.max_dofs_per_line() +
195 25 * dof_handler.fe_collection.max_dofs_per_quad();
196 break;
197 case 14:
198 max_couplings =
199 63 * dof_handler.fe_collection.max_dofs_per_vertex() +
200 98 * dof_handler.fe_collection.max_dofs_per_line() +
201 28 * dof_handler.fe_collection.max_dofs_per_quad();
202 break;
203 case 15:
204 max_couplings =
205 65 * dof_handler.fe_collection.max_dofs_per_vertex() +
206 103 * dof_handler.fe_collection.max_dofs_per_line() +
207 29 * dof_handler.fe_collection.max_dofs_per_quad();
208 break;
209 case 16:
210 max_couplings =
211 72 * dof_handler.fe_collection.max_dofs_per_vertex() +
212 114 * dof_handler.fe_collection.max_dofs_per_line() +
213 32 * dof_handler.fe_collection.max_dofs_per_quad();
214 break;
215
216 default:
217 Assert(false, ExcNotImplemented());
218 max_couplings = 0;
219 }
220 return std::min(max_couplings, dof_handler.n_dofs());
221 }
222
223 template <int spacedim>
224 static unsigned int
225 max_couplings_between_dofs(const DoFHandler<3, spacedim> &dof_handler)
226 {
227 // TODO:[?] Invent significantly better estimates than the ones in this
228 // function
229
230 // doing the same thing here is a rather complicated thing, compared
231 // to the 2d case, since it is hard to draw pictures with several
232 // refined hexahedra :-) so I presently only give a coarse
233 // estimate for the case that at most 8 hexes meet at each vertex
234 //
235 // can anyone give better estimate here?
236 const unsigned int max_adjacent_cells =
237 dof_handler.tria->max_adjacent_cells();
238
239 types::global_dof_index max_couplings;
240 if (max_adjacent_cells <= 8)
241 max_couplings =
242 7 * 7 * 7 * dof_handler.fe_collection.max_dofs_per_vertex() +
243 7 * 6 * 7 * 3 * dof_handler.fe_collection.max_dofs_per_line() +
244 9 * 4 * 7 * 3 * dof_handler.fe_collection.max_dofs_per_quad() +
245 27 * dof_handler.fe_collection.max_dofs_per_hex();
246 else
247 {
248 Assert(false, ExcNotImplemented());
249 max_couplings = 0;
250 }
251
252 return std::min(max_couplings, dof_handler.n_dofs());
253 }
254
259 template <int dim, int spacedim>
260 static void
261 reset_to_empty_objects(DoFHandler<dim, spacedim> &dof_handler)
262 {
263 dof_handler.object_dof_indices.clear();
264 dof_handler.object_dof_indices.resize(dof_handler.tria->n_levels());
265 dof_handler.object_dof_indices.shrink_to_fit();
266
267 dof_handler.object_dof_ptr.clear();
268 dof_handler.object_dof_ptr.resize(dof_handler.tria->n_levels());
269 dof_handler.object_dof_ptr.shrink_to_fit();
270
271 dof_handler.cell_dof_cache_indices.clear();
272 dof_handler.cell_dof_cache_indices.resize(dof_handler.tria->n_levels());
273 dof_handler.cell_dof_cache_indices.shrink_to_fit();
274
275 dof_handler.cell_dof_cache_ptr.clear();
276 dof_handler.cell_dof_cache_ptr.resize(dof_handler.tria->n_levels());
277 dof_handler.cell_dof_cache_ptr.shrink_to_fit();
278 }
279
283 template <int dim, int spacedim>
284 static void
285 reserve_cells(DoFHandler<dim, spacedim> &dof_handler,
286 const unsigned int n_inner_dofs_per_cell)
287 {
288 for (unsigned int i = 0; i < dof_handler.tria->n_levels(); ++i)
289 {
290 // 1) object_dof_indices
291 dof_handler.object_dof_ptr[i][dim].assign(
292 dof_handler.tria->n_raw_cells(i) + 1, 0);
293
294 for (const auto &cell :
295 dof_handler.tria->cell_iterators_on_level(i))
296 if (cell->is_active() && !cell->is_artificial())
297 dof_handler.object_dof_ptr[i][dim][cell->index() + 1] =
298 n_inner_dofs_per_cell;
299
300 for (unsigned int j = 0; j < dof_handler.tria->n_raw_cells(i); ++j)
301 dof_handler.object_dof_ptr[i][dim][j + 1] +=
302 dof_handler.object_dof_ptr[i][dim][j];
303
304 dof_handler.object_dof_indices[i][dim].resize(
305 dof_handler.object_dof_ptr[i][dim].back(),
306 numbers::invalid_dof_index);
307
308 // 2) cell_dof_cache_indices
309 dof_handler.cell_dof_cache_ptr[i].assign(
310 dof_handler.tria->n_raw_cells(i) + 1, 0);
311
312 for (const auto &cell :
313 dof_handler.tria->cell_iterators_on_level(i))
314 if (cell->is_active() && !cell->is_artificial())
315 dof_handler.cell_dof_cache_ptr[i][cell->index() + 1] =
316 dof_handler.get_fe().n_dofs_per_cell();
317
318 for (unsigned int j = 0; j < dof_handler.tria->n_raw_cells(i); ++j)
319 dof_handler.cell_dof_cache_ptr[i][j + 1] +=
320 dof_handler.cell_dof_cache_ptr[i][j];
321
322 dof_handler.cell_dof_cache_indices[i].resize(
323 dof_handler.cell_dof_cache_ptr[i].back(),
324 numbers::invalid_dof_index);
325 }
326 }
327
332 template <int dim, int spacedim, typename T>
333 static void
334 reserve_subentities(DoFHandler<dim, spacedim> &dof_handler,
335 const unsigned int structdim,
336 const unsigned int n_raw_entities,
337 const T & cell_process)
338 {
339 if (dof_handler.tria->n_cells() == 0)
340 return;
341
342 dof_handler.object_dof_ptr[0][structdim].assign(n_raw_entities + 1, -1);
343 // determine for each entity the number of dofs
344 for (const auto &cell : dof_handler.tria->cell_iterators())
345 if (cell->is_active() && !cell->is_artificial())
346 cell_process(
347 cell,
348 [&](const unsigned int n_dofs_per_entity,
349 const unsigned int index) {
350 auto &n_dofs_per_entity_target =
351 dof_handler.object_dof_ptr[0][structdim][index + 1];
352
353 // make sure that either the entity has not been visited or
354 // the entity has the same number of dofs assigned
355 Assert((n_dofs_per_entity_target ==
356 static_cast<
357 typename DoFHandler<dim, spacedim>::offset_type>(
358 -1) ||
359 n_dofs_per_entity_target == n_dofs_per_entity),
360 ExcNotImplemented());
361
362 n_dofs_per_entity_target = n_dofs_per_entity;
363 });
364
365 // convert the absolute numbers to CRS
366 dof_handler.object_dof_ptr[0][structdim][0] = 0;
367 for (unsigned int i = 1; i < n_raw_entities + 1; ++i)
368 {
369 if (dof_handler.object_dof_ptr[0][structdim][i] ==
370 static_cast<typename DoFHandler<dim, spacedim>::offset_type>(
371 -1))
372 dof_handler.object_dof_ptr[0][structdim][i] =
373 dof_handler.object_dof_ptr[0][structdim][i - 1];
374 else
375 dof_handler.object_dof_ptr[0][structdim][i] +=
376 dof_handler.object_dof_ptr[0][structdim][i - 1];
377 }
378
379 // allocate memory for indices
380 dof_handler.object_dof_indices[0][structdim].resize(
381 dof_handler.object_dof_ptr[0][structdim].back(),
382 numbers::invalid_dof_index);
383 }
384
391 template <int dim, int spacedim>
392 static void
393 reserve_space(DoFHandler<dim, spacedim> &dof_handler)
394 {
395 reset_to_empty_objects(dof_handler);
396
397 const auto &fe = dof_handler.get_fe();
398
399 // cell
400 reserve_cells(dof_handler,
401 dim == 1 ? fe.n_dofs_per_line() :
402 (dim == 2 ? fe.n_dofs_per_quad(0) :
403 fe.n_dofs_per_hex()));
404
405 // vertices
406 reserve_subentities(dof_handler,
407 0,
408 dof_handler.tria->n_vertices(),
409 [&](const auto &cell, const auto &process) {
410 for (const auto vertex_index :
411 cell->vertex_indices())
412 process(fe.n_dofs_per_vertex(),
413 cell->vertex_index(vertex_index));
414 });
415
416 // lines
417 if (dim == 2 || dim == 3)
418 reserve_subentities(dof_handler,
419 1,
420 dof_handler.tria->n_raw_lines(),
421 [&](const auto &cell, const auto &process) {
422 for (const auto line_index :
423 cell->line_indices())
424 process(fe.n_dofs_per_line(),
425 cell->line(line_index)->index());
426 });
427
428 // quads
429 if (dim == 3)
430 reserve_subentities(dof_handler,
431 2,
432 dof_handler.tria->n_raw_quads(),
433 [&](const auto &cell, const auto &process) {
434 for (const auto face_index :
435 cell->face_indices())
436 process(fe.n_dofs_per_quad(face_index),
437 cell->face(face_index)->index());
438 });
439 }
440
441 template <int spacedim>
442 static void reserve_space_mg(DoFHandler<1, spacedim> &dof_handler)
443 {
444 Assert(dof_handler.get_triangulation().n_levels() > 0,
445 ExcMessage("Invalid triangulation"));
446 dof_handler.clear_mg_space();
447
448 const ::Triangulation<1, spacedim> &tria =
449 dof_handler.get_triangulation();
450 const unsigned int dofs_per_line =
451 dof_handler.get_fe().n_dofs_per_line();
452 const unsigned int n_levels = tria.n_levels();
453
454 for (unsigned int i = 0; i < n_levels; ++i)
455 {
456 dof_handler.mg_levels.emplace_back(
457 new internal::DoFHandlerImplementation::DoFLevel<1>);
458 dof_handler.mg_levels.back()->dof_object.dofs =
459 std::vector<types::global_dof_index>(tria.n_raw_lines(i) *
460 dofs_per_line,
461 numbers::invalid_dof_index);
462 }
463
464 const unsigned int n_vertices = tria.n_vertices();
465
466 dof_handler.mg_vertex_dofs.resize(n_vertices);
467
468 std::vector<unsigned int> max_level(n_vertices, 0);
469 std::vector<unsigned int> min_level(n_vertices, n_levels);
470
471 for (typename ::Triangulation<1, spacedim>::cell_iterator cell =
472 tria.begin();
473 cell != tria.end();
474 ++cell)
475 {
476 const unsigned int level = cell->level();
477
478 for (const auto vertex : cell->vertex_indices())
479 {
480 const unsigned int vertex_index = cell->vertex_index(vertex);
481
482 if (min_level[vertex_index] > level)
483 min_level[vertex_index] = level;
484
485 if (max_level[vertex_index] < level)
486 max_level[vertex_index] = level;
487 }
488 }
489
490 for (unsigned int vertex = 0; vertex < n_vertices; ++vertex)
491 if (tria.vertex_used(vertex))
492 {
493 Assert(min_level[vertex] < n_levels, ExcInternalError());
494 Assert(max_level[vertex] >= min_level[vertex],
495 ExcInternalError());
496 dof_handler.mg_vertex_dofs[vertex].init(
497 min_level[vertex],
498 max_level[vertex],
499 dof_handler.get_fe().n_dofs_per_vertex());
500 }
501
502 else
503 {
504 Assert(min_level[vertex] == n_levels, ExcInternalError());
505 Assert(max_level[vertex] == 0, ExcInternalError());
506 dof_handler.mg_vertex_dofs[vertex].init(1, 0, 0);
507 }
508 }
509
510 template <int spacedim>
511 static void reserve_space_mg(DoFHandler<2, spacedim> &dof_handler)
512 {
513 Assert(dof_handler.get_triangulation().n_levels() > 0,
514 ExcMessage("Invalid triangulation"));
515 dof_handler.clear_mg_space();
516
517 const ::FiniteElement<2, spacedim> &fe = dof_handler.get_fe();
518 const ::Triangulation<2, spacedim> &tria =
519 dof_handler.get_triangulation();
520 const unsigned int n_levels = tria.n_levels();
521
522 for (unsigned int i = 0; i < n_levels; ++i)
523 {
524 dof_handler.mg_levels.emplace_back(
525 std::make_unique<
526 internal::DoFHandlerImplementation::DoFLevel<2>>());
527 dof_handler.mg_levels.back()->dof_object.dofs =
528 std::vector<types::global_dof_index>(
529 tria.n_raw_quads(i) *
530 fe.n_dofs_per_quad(0 /*note: in 2D there is only one quad*/),
531 numbers::invalid_dof_index);
532 }
533
534 dof_handler.mg_faces =
535 std::make_unique<internal::DoFHandlerImplementation::DoFFaces<2>>();
536 dof_handler.mg_faces->lines.dofs =
537 std::vector<types::global_dof_index>(tria.n_raw_lines() *
538 fe.n_dofs_per_line(),
539 numbers::invalid_dof_index);
540
541 const unsigned int n_vertices = tria.n_vertices();
542
543 dof_handler.mg_vertex_dofs.resize(n_vertices);
544
545 std::vector<unsigned int> max_level(n_vertices, 0);
546 std::vector<unsigned int> min_level(n_vertices, n_levels);
547
548 for (typename ::Triangulation<2, spacedim>::cell_iterator cell =
549 tria.begin();
550 cell != tria.end();
551 ++cell)
552 {
553 const unsigned int level = cell->level();
554
555 for (const auto vertex : cell->vertex_indices())
556 {
557 const unsigned int vertex_index = cell->vertex_index(vertex);
558
559 if (min_level[vertex_index] > level)
560 min_level[vertex_index] = level;
561
562 if (max_level[vertex_index] < level)
563 max_level[vertex_index] = level;
564 }
565 }
566
567 for (unsigned int vertex = 0; vertex < n_vertices; ++vertex)
568 if (tria.vertex_used(vertex))
569 {
570 Assert(min_level[vertex] < n_levels, ExcInternalError());
571 Assert(max_level[vertex] >= min_level[vertex],
572 ExcInternalError());
573 dof_handler.mg_vertex_dofs[vertex].init(min_level[vertex],
574 max_level[vertex],
575 fe.n_dofs_per_vertex());
576 }
577
578 else
579 {
580 Assert(min_level[vertex] == n_levels, ExcInternalError());
581 Assert(max_level[vertex] == 0, ExcInternalError());
582 dof_handler.mg_vertex_dofs[vertex].init(1, 0, 0);
583 }
584 }
585
586 template <int spacedim>
587 static void reserve_space_mg(DoFHandler<3, spacedim> &dof_handler)
588 {
589 Assert(dof_handler.get_triangulation().n_levels() > 0,
590 ExcMessage("Invalid triangulation"));
591 dof_handler.clear_mg_space();
592
593 const ::FiniteElement<3, spacedim> &fe = dof_handler.get_fe();
594 const ::Triangulation<3, spacedim> &tria =
595 dof_handler.get_triangulation();
596 const unsigned int n_levels = tria.n_levels();
597
598 for (unsigned int i = 0; i < n_levels; ++i)
599 {
600 dof_handler.mg_levels.emplace_back(
601 std::make_unique<
602 internal::DoFHandlerImplementation::DoFLevel<3>>());
603 dof_handler.mg_levels.back()->dof_object.dofs =
604 std::vector<types::global_dof_index>(tria.n_raw_hexs(i) *
605 fe.n_dofs_per_hex(),
606 numbers::invalid_dof_index);
607 }
608
609 dof_handler.mg_faces =
610 std::make_unique<internal::DoFHandlerImplementation::DoFFaces<3>>();
611 dof_handler.mg_faces->lines.dofs =
612 std::vector<types::global_dof_index>(tria.n_raw_lines() *
613 fe.n_dofs_per_line(),
614 numbers::invalid_dof_index);
615
616 // TODO: the implementation makes the assumption that all faces have the
617 // same number of dofs
618 AssertDimension(fe.n_unique_faces(), 1);
619 dof_handler.mg_faces->quads.dofs = std::vector<types::global_dof_index>(
620 tria.n_raw_quads() * fe.n_dofs_per_quad(0 /*=face_no*/),
621 numbers::invalid_dof_index);
622
623 const unsigned int n_vertices = tria.n_vertices();
624
625 dof_handler.mg_vertex_dofs.resize(n_vertices);
626
627 std::vector<unsigned int> max_level(n_vertices, 0);
628 std::vector<unsigned int> min_level(n_vertices, n_levels);
629
630 for (typename ::Triangulation<3, spacedim>::cell_iterator cell =
631 tria.begin();
632 cell != tria.end();
633 ++cell)
634 {
635 const unsigned int level = cell->level();
636
637 for (const auto vertex : cell->vertex_indices())
638 {
639 const unsigned int vertex_index = cell->vertex_index(vertex);
640
641 if (min_level[vertex_index] > level)
642 min_level[vertex_index] = level;
643
644 if (max_level[vertex_index] < level)
645 max_level[vertex_index] = level;
646 }
647 }
648
649 for (unsigned int vertex = 0; vertex < n_vertices; ++vertex)
650 if (tria.vertex_used(vertex))
651 {
652 Assert(min_level[vertex] < n_levels, ExcInternalError());
653 Assert(max_level[vertex] >= min_level[vertex],
654 ExcInternalError());
655 dof_handler.mg_vertex_dofs[vertex].init(min_level[vertex],
656 max_level[vertex],
657 fe.n_dofs_per_vertex());
658 }
659
660 else
661 {
662 Assert(min_level[vertex] == n_levels, ExcInternalError());
663 Assert(max_level[vertex] == 0, ExcInternalError());
664 dof_handler.mg_vertex_dofs[vertex].init(1, 0, 0);
665 }
666 }
667
668 template <int spacedim>
669 static types::global_dof_index
670 get_dof_index(
671 const DoFHandler<1, spacedim> &dof_handler,
672 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<1>>
673 &mg_level,
674 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<1>>
675 &,
676 const unsigned int obj_index,
677 const unsigned int fe_index,
678 const unsigned int local_index,
679 const std::integral_constant<int, 1>)
680 {
681 Assert(dof_handler.hp_capability_enabled == false,
682 (typename DoFHandler<1, spacedim>::ExcNotImplementedWithHP()));
683
684 return mg_level->dof_object.get_dof_index(
685 static_cast<const DoFHandler<1, spacedim> &>(dof_handler),
686 obj_index,
687 fe_index,
688 local_index);
689 }
690
691 template <int spacedim>
692 static types::global_dof_index
693 get_dof_index(
694 const DoFHandler<2, spacedim> &dof_handler,
695 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<2>>
696 &,
697 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<2>>
698 & mg_faces,
699 const unsigned int obj_index,
700 const unsigned int fe_index,
701 const unsigned int local_index,
702 const std::integral_constant<int, 1>)
703 {
704 return mg_faces->lines.get_dof_index(
705 static_cast<const DoFHandler<2, spacedim> &>(dof_handler),
706 obj_index,
707 fe_index,
708 local_index);
709 }
710
711 template <int spacedim>
712 static types::global_dof_index
713 get_dof_index(
714 const DoFHandler<2, spacedim> &dof_handler,
715 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<2>>
716 &mg_level,
717 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<2>>
718 &,
719 const unsigned int obj_index,
720 const unsigned int fe_index,
721 const unsigned int local_index,
722 const std::integral_constant<int, 2>)
723 {
724 Assert(dof_handler.hp_capability_enabled == false,
725 (typename DoFHandler<2, spacedim>::ExcNotImplementedWithHP()));
726 return mg_level->dof_object.get_dof_index(
727 static_cast<const DoFHandler<2, spacedim> &>(dof_handler),
728 obj_index,
729 fe_index,
730 local_index);
731 }
732
733 template <int spacedim>
734 static types::global_dof_index
735 get_dof_index(
736 const DoFHandler<3, spacedim> &dof_handler,
737 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<3>>
738 &,
739 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<3>>
740 & mg_faces,
741 const unsigned int obj_index,
742 const unsigned int fe_index,
743 const unsigned int local_index,
744 const std::integral_constant<int, 1>)
745 {
746 Assert(dof_handler.hp_capability_enabled == false,
747 (typename DoFHandler<3, spacedim>::ExcNotImplementedWithHP()));
748 return mg_faces->lines.get_dof_index(
749 static_cast<const DoFHandler<3, spacedim> &>(dof_handler),
750 obj_index,
751 fe_index,
752 local_index);
753 }
754
755 template <int spacedim>
756 static types::global_dof_index
757 get_dof_index(
758 const DoFHandler<3, spacedim> &dof_handler,
759 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<3>>
760 &,
761 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<3>>
762 & mg_faces,
763 const unsigned int obj_index,
764 const unsigned int fe_index,
765 const unsigned int local_index,
766 const std::integral_constant<int, 2>)
767 {
768 Assert(dof_handler.hp_capability_enabled == false,
769 (typename DoFHandler<3, spacedim>::ExcNotImplementedWithHP()));
770 return mg_faces->quads.get_dof_index(
771 static_cast<const DoFHandler<3, spacedim> &>(dof_handler),
772 obj_index,
773 fe_index,
774 local_index);
775 }
776
777 template <int spacedim>
778 static types::global_dof_index
779 get_dof_index(
780 const DoFHandler<3, spacedim> &dof_handler,
781 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<3>>
782 &mg_level,
783 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<3>>
784 &,
785 const unsigned int obj_index,
786 const unsigned int fe_index,
787 const unsigned int local_index,
788 const std::integral_constant<int, 3>)
789 {
790 Assert(dof_handler.hp_capability_enabled == false,
791 (typename DoFHandler<3, spacedim>::ExcNotImplementedWithHP()));
792 return mg_level->dof_object.get_dof_index(
793 static_cast<const DoFHandler<3, spacedim> &>(dof_handler),
794 obj_index,
795 fe_index,
796 local_index);
797 }
798
799 template <int spacedim>
800 static void
801 set_dof_index(
802 const DoFHandler<1, spacedim> &dof_handler,
803 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<1>>
804 &mg_level,
805 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<1>>
806 &,
807 const unsigned int obj_index,
808 const unsigned int fe_index,
809 const unsigned int local_index,
810 const types::global_dof_index global_index,
811 const std::integral_constant<int, 1>)
812 {
813 Assert(dof_handler.hp_capability_enabled == false,
814 (typename DoFHandler<1, spacedim>::ExcNotImplementedWithHP()));
815 mg_level->dof_object.set_dof_index(
816 static_cast<const DoFHandler<1, spacedim> &>(dof_handler),
817 obj_index,
818 fe_index,
819 local_index,
820 global_index);
821 }
822
823 template <int spacedim>
824 static void
825 set_dof_index(
826 const DoFHandler<2, spacedim> &dof_handler,
827 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<2>>
828 &,
829 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<2>>
830 & mg_faces,
831 const unsigned int obj_index,
832 const unsigned int fe_index,
833 const unsigned int local_index,
834 const types::global_dof_index global_index,
835 const std::integral_constant<int, 1>)
836 {
837 Assert(dof_handler.hp_capability_enabled == false,
838 (typename DoFHandler<2, spacedim>::ExcNotImplementedWithHP()));
839 mg_faces->lines.set_dof_index(
840 static_cast<const DoFHandler<2, spacedim> &>(dof_handler),
841 obj_index,
842 fe_index,
843 local_index,
844 global_index);
845 }
846
847 template <int spacedim>
848 static void
849 set_dof_index(
850 const DoFHandler<2, spacedim> &dof_handler,
851 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<2>>
852 &mg_level,
853 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<2>>
854 &,
855 const unsigned int obj_index,
856 const unsigned int fe_index,
857 const unsigned int local_index,
858 const types::global_dof_index global_index,
859 const std::integral_constant<int, 2>)
860 {
861 Assert(dof_handler.hp_capability_enabled == false,
862 (typename DoFHandler<2, spacedim>::ExcNotImplementedWithHP()));
863 mg_level->dof_object.set_dof_index(
864 static_cast<const DoFHandler<2, spacedim> &>(dof_handler),
865 obj_index,
866 fe_index,
867 local_index,
868 global_index);
869 }
870
871 template <int spacedim>
872 static void
873 set_dof_index(
874 const DoFHandler<3, spacedim> &dof_handler,
875 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<3>>
876 &,
877 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<3>>
878 & mg_faces,
879 const unsigned int obj_index,
880 const unsigned int fe_index,
881 const unsigned int local_index,
882 const types::global_dof_index global_index,
883 const std::integral_constant<int, 1>)
884 {
885 Assert(dof_handler.hp_capability_enabled == false,
886 (typename DoFHandler<3, spacedim>::ExcNotImplementedWithHP()));
887 mg_faces->lines.set_dof_index(
888 static_cast<const DoFHandler<3, spacedim> &>(dof_handler),
889 obj_index,
890 fe_index,
891 local_index,
892 global_index);
893 }
894
895 template <int spacedim>
896 static void
897 set_dof_index(
898 const DoFHandler<3, spacedim> &dof_handler,
899 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<3>>
900 &,
901 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<3>>
902 & mg_faces,
903 const unsigned int obj_index,
904 const unsigned int fe_index,
905 const unsigned int local_index,
906 const types::global_dof_index global_index,
907 const std::integral_constant<int, 2>)
908 {
909 Assert(dof_handler.hp_capability_enabled == false,
910 (typename DoFHandler<3, spacedim>::ExcNotImplementedWithHP()));
911 mg_faces->quads.set_dof_index(
912 static_cast<const DoFHandler<3, spacedim> &>(dof_handler),
913 obj_index,
914 fe_index,
915 local_index,
916 global_index);
917 }
918
919 template <int spacedim>
920 static void
921 set_dof_index(
922 const DoFHandler<3, spacedim> &dof_handler,
923 const std::unique_ptr<internal::DoFHandlerImplementation::DoFLevel<3>>
924 &mg_level,
925 const std::unique_ptr<internal::DoFHandlerImplementation::DoFFaces<3>>
926 &,
927 const unsigned int obj_index,
928 const unsigned int fe_index,
929 const unsigned int local_index,
930 const types::global_dof_index global_index,
931 const std::integral_constant<int, 3>)
932 {
933 Assert(dof_handler.hp_capability_enabled == false,
934 (typename DoFHandler<3, spacedim>::ExcNotImplementedWithHP()));
935 mg_level->dof_object.set_dof_index(
936 static_cast<const DoFHandler<3, spacedim> &>(dof_handler),
937 obj_index,
938 fe_index,
939 local_index,
940 global_index);
941 }
942 };
943 } // namespace DoFHandlerImplementation
944
945
946
947 namespace hp
948 {
949 namespace DoFHandlerImplementation
950 {
955 struct Implementation
956 {
962 template <int dim, int spacedim>
963 static void
964 ensure_absence_of_future_fe_indices(
965 DoFHandler<dim, spacedim> &dof_handler)
966 {
967 (void)dof_handler;
968 for (const auto &cell : dof_handler.active_cell_iterators())
969 if (cell->is_locally_owned())
970 Assert(
971 !cell->future_fe_index_set(),
972 ExcMessage(
973 "There shouldn't be any cells flagged for p-adaptation when partitioning."));
974 }
975
976
977
982 template <int dim, int spacedim>
983 static void
984 reserve_space_vertices(DoFHandler<dim, spacedim> &dof_handler)
985 {
986 // The final step in all of the reserve_space() functions is to set
987 // up vertex dof information. since vertices are sequentially
988 // numbered, what we do first is to set up an array in which
989 // we record whether a vertex is associated with any of the
990 // given fe's, by setting a bit. in a later step, we then
991 // actually allocate memory for the required dofs
992 //
993 // in the following, we only need to consider vertices that are
994 // adjacent to either a locally owned or a ghost cell; we never
995 // store anything on vertices that are only surrounded by
996 // artificial cells. so figure out that subset of vertices
997 // first
998 std::vector<bool> locally_used_vertices(
999 dof_handler.tria->n_vertices(), false);
1000 for (const auto &cell : dof_handler.active_cell_iterators())
1001 if (!cell->is_artificial())
1002 for (const auto v : cell->vertex_indices())
1003 locally_used_vertices[cell->vertex_index(v)] = true;
1004
1005 std::vector<std::vector<bool>> vertex_fe_association(
1006 dof_handler.fe_collection.size(),
1007 std::vector<bool>(dof_handler.tria->n_vertices(), false));
1008
1009 for (const auto &cell : dof_handler.active_cell_iterators())
1010 if (!cell->is_artificial())
1011 for (const auto v : cell->vertex_indices())
1012 vertex_fe_association[cell->active_fe_index()]
1013 [cell->vertex_index(v)] = true;
1014
1015 // in debug mode, make sure that each vertex is associated
1016 // with at least one FE (note that except for unused
1017 // vertices, all vertices are actually active). this is of
1018 // course only true for vertices that are part of either
1019 // ghost or locally owned cells
1020#ifdef DEBUG
1021 for (unsigned int v = 0; v < dof_handler.tria->n_vertices(); ++v)
1022 if (locally_used_vertices[v] == true)
1023 if (dof_handler.tria->vertex_used(v) == true)
1024 {
1025 unsigned int fe = 0;
1026 for (; fe < dof_handler.fe_collection.size(); ++fe)
1027 if (vertex_fe_association[fe][v] == true)
1028 break;
1029 Assert(fe != dof_handler.fe_collection.size(),
1030 ExcInternalError());
1031 }
1032#endif
1033
1034 const unsigned int d = 0;
1035 const unsigned int l = 0;
1036
1037 dof_handler.hp_object_fe_ptr[d].clear();
1038 dof_handler.hp_object_fe_indices[d].clear();
1039 dof_handler.object_dof_ptr[l][d].clear();
1040 dof_handler.object_dof_indices[l][d].clear();
1041
1042 dof_handler.hp_object_fe_ptr[d].reserve(
1043 dof_handler.tria->n_vertices() + 1);
1044
1045 unsigned int vertex_slots_needed = 0;
1046 unsigned int fe_slots_needed = 0;
1047
1048 for (unsigned int v = 0; v < dof_handler.tria->n_vertices(); ++v)
1049 {
1050 dof_handler.hp_object_fe_ptr[d].push_back(fe_slots_needed);
1051
1052 if (dof_handler.tria->vertex_used(v) && locally_used_vertices[v])
1053 {
1054 for (unsigned int fe = 0;
1055 fe < dof_handler.fe_collection.size();
1056 ++fe)
1057 if (vertex_fe_association[fe][v] == true)
1058 {
1059 fe_slots_needed++;
1060 vertex_slots_needed +=
1061 dof_handler.get_fe(fe).n_dofs_per_vertex();
1062 }
1063 }
1064 }
1065
1066 dof_handler.hp_object_fe_ptr[d].push_back(fe_slots_needed);
1067
1068 dof_handler.hp_object_fe_indices[d].reserve(fe_slots_needed);
1069 dof_handler.object_dof_ptr[l][d].reserve(fe_slots_needed + 1);
1070
1071 dof_handler.object_dof_indices[l][d].reserve(vertex_slots_needed);
1072
1073 for (unsigned int v = 0; v < dof_handler.tria->n_vertices(); ++v)
1074 if (dof_handler.tria->vertex_used(v) && locally_used_vertices[v])
1075 {
1076 for (unsigned int fe = 0; fe < dof_handler.fe_collection.size();
1077 ++fe)
1078 if (vertex_fe_association[fe][v] == true)
1079 {
1080 dof_handler.hp_object_fe_indices[d].push_back(fe);
1081 dof_handler.object_dof_ptr[l][d].push_back(
1082 dof_handler.object_dof_indices[l][d].size());
1083
1084 for (unsigned int i = 0;
1085 i < dof_handler.get_fe(fe).n_dofs_per_vertex();
1086 i++)
1087 dof_handler.object_dof_indices[l][d].push_back(
1088 numbers::invalid_dof_index);
1089 }
1090 }
1091
1092
1093 dof_handler.object_dof_ptr[l][d].push_back(
1094 dof_handler.object_dof_indices[l][d].size());
1095
1096 AssertDimension(vertex_slots_needed,
1097 dof_handler.object_dof_indices[l][d].size());
1098 AssertDimension(fe_slots_needed,
1099 dof_handler.hp_object_fe_indices[d].size());
1100 AssertDimension(fe_slots_needed + 1,
1101 dof_handler.object_dof_ptr[l][d].size());
1102 AssertDimension(dof_handler.tria->n_vertices() + 1,
1103 dof_handler.hp_object_fe_ptr[d].size());
1104
1105 dof_handler.object_dof_indices[l][d].assign(
1106 vertex_slots_needed, numbers::invalid_dof_index);
1107 }
1108
1109
1110
1115 template <int dim, int spacedim>
1116 static void
1117 reserve_space_cells(DoFHandler<dim, spacedim> &dof_handler)
1118 {
1119 (void)dof_handler;
1120 // count how much space we need on each level for the cell
1121 // dofs and set the dof_*_offsets data. initially set the
1122 // latter to an invalid index, and only later set it to
1123 // something reasonable for active dof_handler.cells
1124 //
1125 // note that for dof_handler.cells, the situation is simpler
1126 // than for other (lower dimensional) objects since exactly
1127 // one finite element is used for it
1128 for (unsigned int level = 0; level < dof_handler.tria->n_levels();
1129 ++level)
1130 {
1131 dof_handler.object_dof_ptr[level][dim] =
1132 std::vector<typename DoFHandler<dim, spacedim>::offset_type>(
1133 dof_handler.tria->n_raw_cells(level),
1134 static_cast<typename DoFHandler<dim, spacedim>::offset_type>(
1135 -1));
1136 dof_handler.cell_dof_cache_ptr[level] =
1137 std::vector<typename DoFHandler<dim, spacedim>::offset_type>(
1138 dof_handler.tria->n_raw_cells(level),
1139 static_cast<typename DoFHandler<dim, spacedim>::offset_type>(
1140 -1));
1141
1142 types::global_dof_index next_free_dof = 0;
1143 types::global_dof_index cache_size = 0;
1144
1145 for (auto cell :
1146 dof_handler.active_cell_iterators_on_level(level))
1147 if (cell->is_active() && !cell->is_artificial())
1148 {
1149 dof_handler.object_dof_ptr[level][dim][cell->index()] =
1150 next_free_dof;
1151 next_free_dof +=
1152 cell->get_fe().template n_dofs_per_object<dim>();
1153
1154 dof_handler.cell_dof_cache_ptr[level][cell->index()] =
1155 cache_size;
1156 cache_size += cell->get_fe().n_dofs_per_cell();
1157 }
1158
1159 dof_handler.object_dof_indices[level][dim] =
1160 std::vector<types::global_dof_index>(
1161 next_free_dof, numbers::invalid_dof_index);
1162 dof_handler.cell_dof_cache_indices[level] =
1163 std::vector<types::global_dof_index>(
1164 cache_size, numbers::invalid_dof_index);
1165 }
1166 }
1167
1168
1169
1174 template <int dim, int spacedim>
1175 static void
1176 reserve_space_faces(DoFHandler<dim, spacedim> &dof_handler)
1177 {
1178 // FACE DOFS
1179 //
1180 // Count face dofs, then allocate as much space
1181 // as we need and prime the linked list for faces (see the
1182 // description in hp::DoFLevel) with the indices we will
1183 // need. Note that our task is more complicated than for the
1184 // cell case above since two adjacent cells may have different
1185 // active FE indices, in which case we need to allocate
1186 // *two* sets of face dofs for the same face. But they don't
1187 // *have* to be different, and so we need to prepare for this
1188 // as well.
1189 //
1190 // The way we do things is that we loop over all active
1191 // cells (these are the only ones that have DoFs
1192 // anyway) and all their faces. We note in the
1193 // user flags whether we have previously visited a face and
1194 // if so skip it (consequently, we have to save and later
1195 // restore the face flags)
1196 {
1197 std::vector<bool> saved_face_user_flags;
1198 switch (dim)
1199 {
1200 case 2:
1201 {
1202 const_cast<::Triangulation<dim, spacedim> &>(
1203 *dof_handler.tria)
1204 .save_user_flags_line(saved_face_user_flags);
1205 const_cast<::Triangulation<dim, spacedim> &>(
1206 *dof_handler.tria)
1207 .clear_user_flags_line();
1208
1209 break;
1210 }
1211
1212 case 3:
1213 {
1214 const_cast<::Triangulation<dim, spacedim> &>(
1215 *dof_handler.tria)
1216 .save_user_flags_quad(saved_face_user_flags);
1217 const_cast<::Triangulation<dim, spacedim> &>(
1218 *dof_handler.tria)
1219 .clear_user_flags_quad();
1220
1221 break;
1222 }
1223
1224 default:
1225 Assert(false, ExcNotImplemented());
1226 }
1227
1228 const unsigned int d = dim - 1;
1229 const unsigned int l = 0;
1230
1231 dof_handler.hp_object_fe_ptr[d].clear();
1232 dof_handler.hp_object_fe_indices[d].clear();
1233 dof_handler.object_dof_ptr[l][d].clear();
1234 dof_handler.object_dof_indices[l][d].clear();
1235
1236 dof_handler.hp_object_fe_ptr[d].resize(
1237 dof_handler.tria->n_raw_faces() + 1);
1238
1239 // An array to hold how many slots (see the hp::DoFLevel
1240 // class) we will have to store on each level
1241 unsigned int n_face_slots = 0;
1242
1243 for (const auto &cell : dof_handler.active_cell_iterators())
1244 if (!cell->is_artificial())
1245 for (const auto face : cell->face_indices())
1246 if (cell->face(face)->user_flag_set() == false)
1247 {
1248 unsigned int fe_slots_needed = 0;
1249
1250 if (cell->at_boundary(face) ||
1251 cell->face(face)->has_children() ||
1252 cell->neighbor_is_coarser(face) ||
1253 (!cell->at_boundary(face) &&
1254 cell->neighbor(face)->is_artificial()) ||
1255 (!cell->at_boundary(face) &&
1256 !cell->neighbor(face)->is_artificial() &&
1257 (cell->active_fe_index() ==
1258 cell->neighbor(face)->active_fe_index())))
1259 {
1260 fe_slots_needed = 1;
1261 n_face_slots +=
1262 dof_handler.get_fe(cell->active_fe_index())
1263 .template n_dofs_per_object<dim - 1>(face);
1264 }
1265 else
1266 {
1267 fe_slots_needed = 2;
1268 n_face_slots +=
1269 dof_handler.get_fe(cell->active_fe_index())
1270 .template n_dofs_per_object<dim - 1>(face) +
1271 dof_handler
1272 .get_fe(cell->neighbor(face)->active_fe_index())
1273 .template n_dofs_per_object<dim - 1>(
1274 cell->neighbor_face_no(face));
1275 }
1276
1277 // mark this face as visited
1278 cell->face(face)->set_user_flag();
1279
1280 dof_handler
1281 .hp_object_fe_ptr[d][cell->face(face)->index() + 1] =
1282 fe_slots_needed;
1283 }
1284
1285 for (unsigned int i = 1; i < dof_handler.hp_object_fe_ptr[d].size();
1286 i++)
1287 dof_handler.hp_object_fe_ptr[d][i] +=
1288 dof_handler.hp_object_fe_ptr[d][i - 1];
1289
1290
1291 dof_handler.hp_object_fe_indices[d].resize(
1292 dof_handler.hp_object_fe_ptr[d].back());
1293 dof_handler.object_dof_ptr[l][d].resize(
1294 dof_handler.hp_object_fe_ptr[d].back() + 1);
1295
1296 dof_handler.object_dof_indices[l][d].reserve(n_face_slots);
1297
1298
1299 // With the memory now allocated, loop over the
1300 // dof_handler cells again and prime the _offset values as
1301 // well as the fe_index fields
1302 switch (dim)
1303 {
1304 case 2:
1305 {
1306 const_cast<::Triangulation<dim, spacedim> &>(
1307 *dof_handler.tria)
1308 .clear_user_flags_line();
1309
1310 break;
1311 }
1312
1313 case 3:
1314 {
1315 const_cast<::Triangulation<dim, spacedim> &>(
1316 *dof_handler.tria)
1317 .clear_user_flags_quad();
1318
1319 break;
1320 }
1321
1322 default:
1323 Assert(false, ExcNotImplemented());
1324 }
1325
1326 for (const auto &cell : dof_handler.active_cell_iterators())
1327 if (!cell->is_artificial())
1328 for (const auto face : cell->face_indices())
1329 if (!cell->face(face)->user_flag_set())
1330 {
1331 // Same decision tree as before
1332 if (cell->at_boundary(face) ||
1333 cell->face(face)->has_children() ||
1334 cell->neighbor_is_coarser(face) ||
1335 (!cell->at_boundary(face) &&
1336 cell->neighbor(face)->is_artificial()) ||
1337 (!cell->at_boundary(face) &&
1338 !cell->neighbor(face)->is_artificial() &&
1339 (cell->active_fe_index() ==
1340 cell->neighbor(face)->active_fe_index())))
1341 {
1342 const unsigned int fe = cell->active_fe_index();
1343 const unsigned int n_dofs =
1344 dof_handler.get_fe(fe)
1345 .template n_dofs_per_object<dim - 1>(face);
1346 const unsigned int offset =
1347 dof_handler
1348 .hp_object_fe_ptr[d][cell->face(face)->index()];
1349
1350 dof_handler.hp_object_fe_indices[d][offset] = fe;
1351 dof_handler.object_dof_ptr[l][d][offset + 1] = n_dofs;
1352
1353 for (unsigned int i = 0; i < n_dofs; i++)
1354 dof_handler.object_dof_indices[l][d].push_back(
1355 numbers::invalid_dof_index);
1356 }
1357 else
1358 {
1359 unsigned int fe_1 = cell->active_fe_index();
1360 unsigned int face_no_1 = face;
1361 unsigned int fe_2 =
1362 cell->neighbor(face)->active_fe_index();
1363 unsigned int face_no_2 = cell->neighbor_face_no(face);
1364
1365 if (fe_2 < fe_1)
1366 {
1367 std::swap(fe_1, fe_2);
1368 std::swap(face_no_1, face_no_2);
1369 }
1370
1371 const unsigned int n_dofs_1 =
1372 dof_handler.get_fe(fe_1)
1373 .template n_dofs_per_object<dim - 1>(face_no_1);
1374
1375 const unsigned int n_dofs_2 =
1376 dof_handler.get_fe(fe_2)
1377 .template n_dofs_per_object<dim - 1>(face_no_2);
1378
1379 const unsigned int offset =
1380 dof_handler
1381 .hp_object_fe_ptr[d][cell->face(face)->index()];
1382
1383 dof_handler.hp_object_fe_indices[d].push_back(
1384 cell->active_fe_index());
1385 dof_handler.object_dof_ptr[l][d].push_back(
1386 dof_handler.object_dof_indices[l][d].size());
1387
1388 dof_handler.hp_object_fe_indices[d][offset + 0] =
1389 fe_1;
1390 dof_handler.hp_object_fe_indices[d][offset + 1] =
1391 fe_2;
1392 dof_handler.object_dof_ptr[l][d][offset + 1] =
1393 n_dofs_1;
1394 dof_handler.object_dof_ptr[l][d][offset + 2] =
1395 n_dofs_2;
1396
1397
1398 for (unsigned int i = 0; i < n_dofs_1 + n_dofs_2; i++)
1399 dof_handler.object_dof_indices[l][d].push_back(
1400 numbers::invalid_dof_index);
1401 }
1402
1403 // mark this face as visited
1404 cell->face(face)->set_user_flag();
1405 }
1406
1407 for (unsigned int i = 1;
1408 i < dof_handler.object_dof_ptr[l][d].size();
1409 i++)
1410 dof_handler.object_dof_ptr[l][d][i] +=
1411 dof_handler.object_dof_ptr[l][d][i - 1];
1412
1413 // at the end, restore the user flags for the faces
1414 switch (dim)
1415 {
1416 case 2:
1417 {
1418 const_cast<::Triangulation<dim, spacedim> &>(
1419 *dof_handler.tria)
1420 .load_user_flags_line(saved_face_user_flags);
1421
1422 break;
1423 }
1424
1425 case 3:
1426 {
1427 const_cast<::Triangulation<dim, spacedim> &>(
1428 *dof_handler.tria)
1429 .load_user_flags_quad(saved_face_user_flags);
1430
1431 break;
1432 }
1433
1434 default:
1435 Assert(false, ExcNotImplemented());
1436 }
1437 }
1438 }
1439
1440
1441
1448 template <int spacedim>
1449 static void reserve_space(::DoFHandler<1, spacedim> &dof_handler)
1450 {
1451 Assert(dof_handler.fe_collection.size() > 0,
1452 (typename ::DoFHandler<1, spacedim>::ExcNoFESelected()));
1453 Assert(dof_handler.tria->n_levels() > 0,
1454 ExcMessage("The current Triangulation must not be empty."));
1455 Assert(dof_handler.tria->n_levels() ==
1456 dof_handler.hp_cell_future_fe_indices.size(),
1457 ExcInternalError());
1458
1459 internal::DoFHandlerImplementation::Implementation::
1460 reset_to_empty_objects(dof_handler);
1461
1462 Threads::TaskGroup<> tasks;
1463 tasks +=
1464 Threads::new_task(&reserve_space_cells<1, spacedim>, dof_handler);
1465 tasks += Threads::new_task(&reserve_space_vertices<1, spacedim>,
1466 dof_handler);
1467 tasks.join_all();
1468 }
1469
1470
1471
1472 template <int spacedim>
1473 static void reserve_space(::DoFHandler<2, spacedim> &dof_handler)
1474 {
1475 Assert(dof_handler.fe_collection.size() > 0,
1476 (typename ::DoFHandler<1, spacedim>::ExcNoFESelected()));
1477 Assert(dof_handler.tria->n_levels() > 0,
1478 ExcMessage("The current Triangulation must not be empty."));
1479 Assert(dof_handler.tria->n_levels() ==
1480 dof_handler.hp_cell_future_fe_indices.size(),
1481 ExcInternalError());
1482
1483 internal::DoFHandlerImplementation::Implementation::
1484 reset_to_empty_objects(dof_handler);
1485
1486 Threads::TaskGroup<> tasks;
1487 tasks +=
1488 Threads::new_task(&reserve_space_cells<2, spacedim>, dof_handler);
1489 tasks +=
1490 Threads::new_task(&reserve_space_faces<2, spacedim>, dof_handler);
1491 tasks += Threads::new_task(&reserve_space_vertices<2, spacedim>,
1492 dof_handler);
1493 tasks.join_all();
1494 }
1495
1496
1497
1498 template <int spacedim>
1499 static void reserve_space(::DoFHandler<3, spacedim> &dof_handler)
1500 {
1501 Assert(dof_handler.fe_collection.size() > 0,
1502 (typename ::DoFHandler<1, spacedim>::ExcNoFESelected()));
1503 Assert(dof_handler.tria->n_levels() > 0,
1504 ExcMessage("The current Triangulation must not be empty."));
1505 Assert(dof_handler.tria->n_levels() ==
1506 dof_handler.hp_cell_future_fe_indices.size(),
1507 ExcInternalError());
1508
1509 internal::DoFHandlerImplementation::Implementation::
1510 reset_to_empty_objects(dof_handler);
1511
1512 Threads::TaskGroup<> tasks;
1513 tasks +=
1514 Threads::new_task(&reserve_space_cells<3, spacedim>, dof_handler);
1515 tasks +=
1516 Threads::new_task(&reserve_space_faces<3, spacedim>, dof_handler);
1517 tasks += Threads::new_task(&reserve_space_vertices<3, spacedim>,
1518 dof_handler);
1519
1520 // While the tasks above are running, we can turn to line dofs
1521
1522 // the situation here is pretty much like with vertices:
1523 // there can be an arbitrary number of finite elements
1524 // associated with each line.
1525 //
1526 // the algorithm we use is somewhat similar to what we do in
1527 // reserve_space_vertices()
1528 {
1529 // what we do first is to set up an array in which we
1530 // record whether a line is associated with any of the
1531 // given fe's, by setting a bit. in a later step, we
1532 // then actually allocate memory for the required dofs
1533 std::vector<std::vector<bool>> line_fe_association(
1534 dof_handler.fe_collection.size(),
1535 std::vector<bool>(dof_handler.tria->n_raw_lines(), false));
1536
1537 for (const auto &cell : dof_handler.active_cell_iterators())
1538 if (!cell->is_artificial())
1539 for (const auto l : cell->line_indices())
1540 line_fe_association[cell->active_fe_index()]
1541 [cell->line_index(l)] = true;
1542
1543 // first check which of the lines is used at all,
1544 // i.e. is associated with a finite element. we do this
1545 // since not all lines may actually be used, in which
1546 // case we do not have to allocate any memory at all
1547 std::vector<bool> line_is_used(dof_handler.tria->n_raw_lines(),
1548 false);
1549 for (unsigned int line = 0; line < dof_handler.tria->n_raw_lines();
1550 ++line)
1551 for (unsigned int fe = 0; fe < dof_handler.fe_collection.size();
1552 ++fe)
1553 if (line_fe_association[fe][line] == true)
1554 {
1555 line_is_used[line] = true;
1556 break;
1557 }
1558
1559
1560
1561 const unsigned int d = 1;
1562 const unsigned int l = 0;
1563
1564 dof_handler.hp_object_fe_ptr[d].clear();
1565 dof_handler.hp_object_fe_indices[d].clear();
1566 dof_handler.object_dof_ptr[l][d].clear();
1567 dof_handler.object_dof_indices[l][d].clear();
1568
1569 dof_handler.hp_object_fe_ptr[d].reserve(
1570 dof_handler.tria->n_raw_lines() + 1);
1571
1572 unsigned int line_slots_needed = 0;
1573 unsigned int fe_slots_needed = 0;
1574
1575 for (unsigned int line = 0; line < dof_handler.tria->n_raw_lines();
1576 ++line)
1577 {
1578 dof_handler.hp_object_fe_ptr[d].push_back(fe_slots_needed);
1579
1580 if (line_is_used[line] == true)
1581 {
1582 for (unsigned int fe = 0;
1583 fe < dof_handler.fe_collection.size();
1584 ++fe)
1585 if (line_fe_association[fe][line] == true)
1586 {
1587 fe_slots_needed++;
1588 line_slots_needed +=
1589 dof_handler.get_fe(fe).n_dofs_per_line();
1590 }
1591 }
1592 }
1593
1594 dof_handler.hp_object_fe_ptr[d].push_back(fe_slots_needed);
1595
1596 // make sure that all entries have been set
1597 AssertDimension(dof_handler.hp_object_fe_ptr[d].size(),
1598 dof_handler.tria->n_raw_lines() + 1);
1599
1600 dof_handler.hp_object_fe_indices[d].reserve(fe_slots_needed);
1601 dof_handler.object_dof_ptr[l][d].reserve(fe_slots_needed + 1);
1602
1603 dof_handler.object_dof_indices[l][d].reserve(line_slots_needed);
1604
1605 for (unsigned int line = 0; line < dof_handler.tria->n_raw_lines();
1606 ++line)
1607 if (line_is_used[line] == true)
1608 {
1609 for (unsigned int fe = 0;
1610 fe < dof_handler.fe_collection.size();
1611 ++fe)
1612 if (line_fe_association[fe][line] == true)
1613 {
1614 dof_handler.hp_object_fe_indices[d].push_back(fe);
1615 dof_handler.object_dof_ptr[l][d].push_back(
1616 dof_handler.object_dof_indices[l][d].size());
1617
1618 for (unsigned int i = 0;
1619 i < dof_handler.get_fe(fe).n_dofs_per_line();
1620 i++)
1621 dof_handler.object_dof_indices[l][d].push_back(
1622 numbers::invalid_dof_index);
1623 }
1624 }
1625
1626 dof_handler.object_dof_ptr[l][d].push_back(
1627 dof_handler.object_dof_indices[l][d].size());
1628
1629 // make sure that all entries have been set
1630 AssertDimension(dof_handler.hp_object_fe_indices[d].size(),
1631 fe_slots_needed);
1632 AssertDimension(dof_handler.object_dof_ptr[l][d].size(),
1633 fe_slots_needed + 1);
1634 AssertDimension(dof_handler.object_dof_indices[l][d].size(),
1635 line_slots_needed);
1636 }
1637
1638 // Ensure that everything is done at this point.
1639 tasks.join_all();
1640 }
1641
1642
1643
1655 template <int dim, int spacedim>
1656 static void
1657 communicate_active_fe_indices(DoFHandler<dim, spacedim> &dof_handler)
1658 {
1659 Assert(
1660 dof_handler.hp_capability_enabled == true,
1661 (typename DoFHandler<dim, spacedim>::ExcOnlyAvailableWithHP()));
1662
1663 using active_fe_index_type =
1664 typename ::DoFHandler<dim, spacedim>::active_fe_index_type;
1665
1666 if (const ::parallel::shared::Triangulation<dim, spacedim> *tr =
1667 dynamic_cast<
1668 const ::parallel::shared::Triangulation<dim, spacedim>
1669 *>(&dof_handler.get_triangulation()))
1670 {
1671 // we have a shared triangulation. in this case, every processor
1672 // knows about all cells, but every processor only has knowledge
1673 // about the active FE index on the cells it owns.
1674 //
1675 // we can create a complete set of active FE indices by letting
1676 // every processor create a vector of indices for all cells,
1677 // filling only those on the cells it owns and setting the indices
1678 // on the other cells to zero. then we add all of these vectors
1679 // up, and because every vector entry has exactly one processor
1680 // that owns it, the sum is correct
1681 std::vector<active_fe_index_type> active_fe_indices(
1682 tr->n_active_cells(), 0u);
1683 for (const auto &cell : dof_handler.active_cell_iterators())
1684 if (cell->is_locally_owned())
1685 active_fe_indices[cell->active_cell_index()] =
1686 cell->active_fe_index();
1687
1688 Utilities::MPI::sum(active_fe_indices,
1689 tr->get_communicator(),
1690 active_fe_indices);
1691
1692 // now go back and fill the active FE index on all other
1693 // cells. we would like to call cell->set_active_fe_index(),
1694 // but that function does not allow setting these indices on
1695 // non-locally_owned cells. so we have to work around the
1696 // issue a little bit by accessing the underlying data
1697 // structures directly
1698 for (const auto &cell : dof_handler.active_cell_iterators())
1699 if (!cell->is_locally_owned())
1700 dof_handler
1701 .hp_cell_active_fe_indices[cell->level()][cell->index()] =
1702 active_fe_indices[cell->active_cell_index()];
1703 }
1704 else if (const ::parallel::
1705 DistributedTriangulationBase<dim, spacedim> *tr =
1706 dynamic_cast<
1707 const ::parallel::
1708 DistributedTriangulationBase<dim, spacedim> *>(
1709 &dof_handler.get_triangulation()))
1710 {
1711 // For completely distributed meshes, use the function that is
1712 // able to move data from locally owned cells on one processor to
1713 // the corresponding ghost cells on others. To this end, we need
1714 // to have functions that can pack and unpack the data we want to
1715 // transport -- namely, the single unsigned int active_fe_index
1716 // objects
1717 auto pack =
1718 [](const typename ::DoFHandler<dim, spacedim>::
1719 active_cell_iterator &cell) -> active_fe_index_type {
1720 return cell->active_fe_index();
1721 };
1722
1723 auto unpack =
1724 [&dof_handler](
1725 const typename ::DoFHandler<dim, spacedim>::
1726 active_cell_iterator & cell,
1727 const active_fe_index_type active_fe_index) -> void {
1728 // we would like to say
1729 // cell->set_active_fe_index(active_fe_index);
1730 // but this is not allowed on cells that are not
1731 // locally owned, and we are on a ghost cell
1732 dof_handler
1733 .hp_cell_active_fe_indices[cell->level()][cell->index()] =
1734 active_fe_index;
1735 };
1736
1737 GridTools::exchange_cell_data_to_ghosts<
1738 active_fe_index_type,
1739 ::DoFHandler<dim, spacedim>>(dof_handler, pack, unpack);
1740 }
1741 else
1742 {
1743 // a sequential triangulation. there is nothing we need to do here
1744 Assert(
1745 (dynamic_cast<
1746 const ::parallel::TriangulationBase<dim, spacedim> *>(
1747 &dof_handler.get_triangulation()) == nullptr),
1748 ExcInternalError());
1749 }
1750 }
1751
1752
1753
1767 template <int dim, int spacedim>
1768 static void
1769 communicate_future_fe_indices(DoFHandler<dim, spacedim> &dof_handler)
1770 {
1771 Assert(
1772 dof_handler.hp_capability_enabled == true,
1773 (typename DoFHandler<dim, spacedim>::ExcOnlyAvailableWithHP()));
1774
1775 using active_fe_index_type =
1776 typename ::DoFHandler<dim, spacedim>::active_fe_index_type;
1777
1778 if (const ::parallel::shared::Triangulation<dim, spacedim> *tr =
1779 dynamic_cast<
1780 const ::parallel::shared::Triangulation<dim, spacedim>
1781 *>(&dof_handler.get_triangulation()))
1782 {
1783 std::vector<active_fe_index_type> future_fe_indices(
1784 tr->n_active_cells(), 0u);
1785 for (const auto &cell : dof_handler.active_cell_iterators())
1786 if (cell->is_locally_owned())
1787 future_fe_indices[cell->active_cell_index()] =
1788 dof_handler
1789 .hp_cell_future_fe_indices[cell->level()][cell->index()];
1790
1791 Utilities::MPI::sum(future_fe_indices,
1792 tr->get_communicator(),
1793 future_fe_indices);
1794
1795 for (const auto &cell : dof_handler.active_cell_iterators())
1796 if (!cell->is_locally_owned())
1797 dof_handler
1798 .hp_cell_future_fe_indices[cell->level()][cell->index()] =
1799 future_fe_indices[cell->active_cell_index()];
1800 }
1801 else if (const ::parallel::
1802 DistributedTriangulationBase<dim, spacedim> *tr =
1803 dynamic_cast<
1804 const ::parallel::
1805 DistributedTriangulationBase<dim, spacedim> *>(
1806 &dof_handler.get_triangulation()))
1807 {
1808 auto pack =
1809 [&dof_handler](
1810 const typename ::DoFHandler<dim, spacedim>::
1811 active_cell_iterator &cell) -> active_fe_index_type {
1812 return dof_handler
1813 .hp_cell_future_fe_indices[cell->level()][cell->index()];
1814 };
1815
1816 auto unpack =
1817 [&dof_handler](
1818 const typename ::DoFHandler<dim, spacedim>::
1819 active_cell_iterator & cell,
1820 const active_fe_index_type future_fe_index) -> void {
1821 dof_handler
1822 .hp_cell_future_fe_indices[cell->level()][cell->index()] =
1823 future_fe_index;
1824 };
1825
1826 GridTools::exchange_cell_data_to_ghosts<
1827 active_fe_index_type,
1828 ::DoFHandler<dim, spacedim>>(dof_handler, pack, unpack);
1829 }
1830 else
1831 {
1832 Assert(
1833 (dynamic_cast<
1834 const ::parallel::TriangulationBase<dim, spacedim> *>(
1835 &dof_handler.get_triangulation()) == nullptr),
1836 ExcInternalError());
1837 }
1838 }
1839
1840
1841
1862 template <int dim, int spacedim>
1863 static void
1864 collect_fe_indices_on_cells_to_be_refined(
1865 DoFHandler<dim, spacedim> &dof_handler)
1866 {
1867 const auto &fe_transfer = dof_handler.active_fe_index_transfer;
1868
1869 for (const auto &cell : dof_handler.active_cell_iterators())
1870 if (cell->is_locally_owned())
1871 {
1872 if (cell->refine_flag_set())
1873 {
1874 // Store the active FE index of each cell that will be
1875 // refined to and distribute it later on its children.
1876 // Pick their future index if flagged for p-refinement.
1877 fe_transfer->refined_cells_fe_index.insert(
1878 {cell, cell->future_fe_index()});
1879 }
1880 else if (cell->coarsen_flag_set())
1881 {
1882 // From all cells that will be coarsened, determine their
1883 // parent and calculate its proper active FE index, so that
1884 // it can be set after refinement. But first, check if that
1885 // particular cell has a parent at all.
1886 Assert(cell->level() > 0, ExcInternalError());
1887 const auto &parent = cell->parent();
1888
1889 // Check if the active FE index for the current cell has
1890 // been determined already.
1891 if (fe_transfer->coarsened_cells_fe_index.find(parent) ==
1892 fe_transfer->coarsened_cells_fe_index.end())
1893 {
1894 // Find a suitable active FE index for the parent cell
1895 // based on the 'least dominant finite element' of its
1896 // children. Consider the childrens' hypothetical future
1897 // index when they have been flagged for p-refinement.
1898#ifdef DEBUG
1899 for (const auto &child : parent->child_iterators())
1900 Assert(child->is_active() &&
1901 child->coarsen_flag_set(),
1902 typename ::Triangulation<
1903 dim>::ExcInconsistentCoarseningFlags());
1904#endif
1905
1906 const unsigned int fe_index = ::internal::hp::
1907 DoFHandlerImplementation::Implementation::
1908 dominated_future_fe_on_children<dim, spacedim>(
1909 parent);
1910
1911 fe_transfer->coarsened_cells_fe_index.insert(
1912 {parent, fe_index});
1913 }
1914 }
1915 else
1916 {
1917 // No h-refinement is scheduled for this cell.
1918 // However, it may have p-refinement indicators, so we
1919 // choose a new active FE index based on its flags.
1920 if (cell->future_fe_index_set() == true)
1921 fe_transfer->persisting_cells_fe_index.insert(
1922 {cell, cell->future_fe_index()});
1923 }
1924 }
1925 }
1926
1927
1928
1933 template <int dim, int spacedim>
1934 static void
1935 distribute_fe_indices_on_refined_cells(
1936 DoFHandler<dim, spacedim> &dof_handler)
1937 {
1938 const auto &fe_transfer = dof_handler.active_fe_index_transfer;
1939
1940 // Set active FE indices on persisting cells.
1941 for (const auto &persist : fe_transfer->persisting_cells_fe_index)
1942 {
1943 const auto &cell = persist.first;
1944
1945 if (cell->is_locally_owned())
1946 {
1947 Assert(cell->is_active(), ExcInternalError());
1948 cell->set_active_fe_index(persist.second);
1949 }
1950 }
1951
1952 // Distribute active FE indices from all refined cells on their
1953 // respective children.
1954 for (const auto &refine : fe_transfer->refined_cells_fe_index)
1955 {
1956 const auto &parent = refine.first;
1957
1958 for (const auto &child : parent->child_iterators())
1959 if (child->is_locally_owned())
1960 {
1961 Assert(child->is_active(), ExcInternalError());
1962 child->set_active_fe_index(refine.second);
1963 }
1964 }
1965
1966 // Set active FE indices on coarsened cells that have been determined
1967 // before the actual coarsening happened.
1968 for (const auto &coarsen : fe_transfer->coarsened_cells_fe_index)
1969 {
1970 const auto &cell = coarsen.first;
1971
1972 if (cell->is_locally_owned())
1973 {
1974 Assert(cell->is_active(), ExcInternalError());
1975 cell->set_active_fe_index(coarsen.second);
1976 }
1977 }
1978 }
1979
1980
1991 template <int dim, int spacedim>
1992 static unsigned int
1993 determine_fe_from_children(
1994 const typename Triangulation<dim, spacedim>::cell_iterator &,
1995 const std::vector<unsigned int> & children_fe_indices,
1996 const ::hp::FECollection<dim, spacedim> &fe_collection)
1997 {
1998 Assert(!children_fe_indices.empty(), ExcInternalError());
1999
2000 // convert vector to set
2001 const std::set<unsigned int> children_fe_indices_set(
2002 children_fe_indices.begin(), children_fe_indices.end());
2003
2004 const unsigned int dominated_fe_index =
2005 fe_collection.find_dominated_fe_extended(children_fe_indices_set,
2006 /*codim=*/0);
2007
2008 Assert(dominated_fe_index != numbers::invalid_unsigned_int,
2009 ExcNoDominatedFiniteElementOnChildren());
2010
2011 return dominated_fe_index;
2012 }
2013
2014
2022 template <int dim, int spacedim>
2023 static unsigned int
2024 dominated_future_fe_on_children(
2025 const typename DoFHandler<dim, spacedim>::cell_iterator &parent)
2026 {
2027 Assert(
2028 !parent->is_active(),
2029 ExcMessage(
2030 "You ask for information on children of this cell which is only "
2031 "available for active cells. This cell has no children."));
2032
2033 const auto &dof_handler = parent->get_dof_handler();
2034 Assert(
2035 dof_handler.has_hp_capabilities(),
2036 (typename DoFHandler<dim, spacedim>::ExcOnlyAvailableWithHP()));
2037
2038 std::set<unsigned int> future_fe_indices_children;
2039 for (const auto &child : parent->child_iterators())
2040 {
2041 Assert(
2042 child->is_active(),
2043 ExcMessage(
2044 "You ask for information on children of this cell which is only "
2045 "available for active cells. One of its children is not active."));
2046
2047 // Ghost siblings might occur on parallel::shared::Triangulation
2048 // objects. The public interface does not allow to access future
2049 // FE indices on ghost cells. However, we need this information
2050 // here and thus call the internal function that does not check
2051 // for cell ownership. This requires that future FE indices have
2052 // been communicated prior to calling this function.
2053 const unsigned int future_fe_index_child =
2054 ::internal::DoFCellAccessorImplementation::
2055 Implementation::future_fe_index<dim, spacedim, false>(*child);
2056
2057 future_fe_indices_children.insert(future_fe_index_child);
2058 }
2059 Assert(!future_fe_indices_children.empty(), ExcInternalError());
2060
2061 const unsigned int future_fe_index =
2062 dof_handler.fe_collection.find_dominated_fe_extended(
2063 future_fe_indices_children,
2064 /*codim=*/0);
2065
2066 Assert(future_fe_index != numbers::invalid_unsigned_int,
2067 ExcNoDominatedFiniteElementOnChildren());
2068
2069 return future_fe_index;
2070 }
2071 };
2072
2073
2074
2078 template <int dim, int spacedim>
2079 void
2080 communicate_future_fe_indices(DoFHandler<dim, spacedim> &dof_handler)
2081 {
2082 Implementation::communicate_future_fe_indices<dim, spacedim>(
2083 dof_handler);
2084 }
2085
2086
2087
2091 template <int dim, int spacedim>
2092 unsigned int
2093 dominated_future_fe_on_children(
2094 const typename DoFHandler<dim, spacedim>::cell_iterator &parent)
2095 {
2096 return Implementation::dominated_future_fe_on_children<dim, spacedim>(
2097 parent);
2098 }
2099 } // namespace DoFHandlerImplementation
2100 } // namespace hp
2101} // namespace internal
2102
2103
2104
2105template <int dim, int spacedim>
2106DoFHandler<dim, spacedim>::DoFHandler()
2107 : hp_capability_enabled(true)
2108 , tria(nullptr, typeid(*this).name())
2109 , mg_faces(nullptr)
2110{}
2111
2112
2113
2114template <int dim, int spacedim>
2115DoFHandler<dim, spacedim>::DoFHandler(const Triangulation<dim, spacedim> &tria)
2116 : DoFHandler()
2117{
2118 reinit(tria);
2119}
2120
2121
2122
2123template <int dim, int spacedim>
2124DoFHandler<dim, spacedim>::~DoFHandler()
2125{
2126 // unsubscribe all attachments to signals of the underlying triangulation
2127 for (auto &connection : this->tria_listeners)
2128 connection.disconnect();
2129 this->tria_listeners.clear();
2130
2131 for (auto &connection : this->tria_listeners_for_transfer)
2132 connection.disconnect();
2133 this->tria_listeners_for_transfer.clear();
2134
2135 // release allocated memory
2136 // virtual functions called in constructors and destructors never use the
2137 // override in a derived class
2138 // for clarity be explicit on which function is called
2139 DoFHandler<dim, spacedim>::clear();
2140
2141 // also release the policy. this needs to happen before the
2142 // current object disappears because the policy objects
2143 // store references to the DoFhandler object they work on
2144 this->policy.reset();
2145}
2146
2147
2148
2149template <int dim, int spacedim>
2150void
2151DoFHandler<dim, spacedim>::initialize(const Triangulation<dim, spacedim> &tria,
2152 const FiniteElement<dim, spacedim> &fe)
2153{
2154 this->initialize(tria, hp::FECollection<dim, spacedim>(fe));
2155}
2156
2157
2158
2159template <int dim, int spacedim>
2160void
2161DoFHandler<dim, spacedim>::initialize(const Triangulation<dim, spacedim> &tria,
2162 const hp::FECollection<dim, spacedim> &fe)
2163{
2164 this->reinit(tria);
2165 this->distribute_dofs(fe);
2166}
2167
2168
2169
2170template <int dim, int spacedim>
2171void
2172DoFHandler<dim, spacedim>::reinit(const Triangulation<dim, spacedim> &tria)
2173{
2174 //
2175 // call destructor
2176 //
2177 // remove association with old triangulation
2178 for (auto &connection : this->tria_listeners)
2179 connection.disconnect();
2180 this->tria_listeners.clear();
2181
2182 for (auto &connection : this->tria_listeners_for_transfer)
2183 connection.disconnect();
2184 this->tria_listeners_for_transfer.clear();
2185
2186 // release allocated memory and policy
2187 DoFHandler<dim, spacedim>::clear();
2188 this->policy.reset();
2189
2190 // reset the finite element collection
2191 this->fe_collection = hp::FECollection<dim, spacedim>();
2192
2193 //
2194 // call constructor
2195 //
2196 // establish connection to new triangulation
2197 this->tria = &tria;
2198 this->setup_policy();
2199
2200 // start in hp-mode and let distribute_dofs toggle it if necessary
2201 hp_capability_enabled = true;
2202 this->connect_to_triangulation_signals();
2203 this->create_active_fe_table();
2204}
2205
2206
2207
2208/*------------------------ Cell iterator functions ------------------------*/
2209
2210template <int dim, int spacedim>
2211typename DoFHandler<dim, spacedim>::cell_iterator
2212DoFHandler<dim, spacedim>::begin(const unsigned int level) const
2213{
2214 typename Triangulation<dim, spacedim>::cell_iterator cell =
2215 this->get_triangulation().begin(level);
2216 if (cell == this->get_triangulation().end(level))
2217 return end(level);
2218 return cell_iterator(*cell, this);
2219}
2220
2221
2222
2223template <int dim, int spacedim>
2224typename DoFHandler<dim, spacedim>::active_cell_iterator
2225DoFHandler<dim, spacedim>::begin_active(const unsigned int level) const
2226{
2227 // level is checked in begin
2228 cell_iterator i = begin(level);
2229 if (i.state() != IteratorState::valid)
2230 return i;
2231 while (i->has_children())
2232 if ((++i).state() != IteratorState::valid)
2233 return i;
2234 return i;
2235}
2236
2237
2238
2239template <int dim, int spacedim>
2240typename DoFHandler<dim, spacedim>::cell_iterator
2241DoFHandler<dim, spacedim>::end() const
2242{
2243 return cell_iterator(&this->get_triangulation(), -1, -1, this);
2244}
2245
2246
2247
2248template <int dim, int spacedim>
2249typename DoFHandler<dim, spacedim>::cell_iterator
2250DoFHandler<dim, spacedim>::end(const unsigned int level) const
2251{
2252 typename Triangulation<dim, spacedim>::cell_iterator cell =
2253 this->get_triangulation().end(level);
2254 if (cell.state() != IteratorState::valid)
2255 return end();
2256 return cell_iterator(*cell, this);
2257}
2258
2259
2260
2261template <int dim, int spacedim>
2262typename DoFHandler<dim, spacedim>::active_cell_iterator
2263DoFHandler<dim, spacedim>::end_active(const unsigned int level) const
2264{
2265 typename Triangulation<dim, spacedim>::cell_iterator cell =
2266 this->get_triangulation().end_active(level);
2267 if (cell.state() != IteratorState::valid)
2268 return active_cell_iterator(end());
2269 return active_cell_iterator(*cell, this);
2270}
2271
2272
2273
2274template <int dim, int spacedim>
2275typename DoFHandler<dim, spacedim>::level_cell_iterator
2276DoFHandler<dim, spacedim>::begin_mg(const unsigned int level) const
2277{
2278 Assert(this->has_level_dofs(),
2279 ExcMessage("You can only iterate over mg "
2280 "levels if mg dofs got distributed."));
2281 typename Triangulation<dim, spacedim>::cell_iterator cell =
2282 this->get_triangulation().begin(level);
2283 if (cell == this->get_triangulation().end(level))
2284 return end_mg(level);
2285 return level_cell_iterator(*cell, this);
2286}
2287
2288
2289
2290template <int dim, int spacedim>
2291typename DoFHandler<dim, spacedim>::level_cell_iterator
2292DoFHandler<dim, spacedim>::end_mg(const unsigned int level) const
2293{
2294 Assert(this->has_level_dofs(),
2295 ExcMessage("You can only iterate over mg "
2296 "levels if mg dofs got distributed."));
2297 typename Triangulation<dim, spacedim>::cell_iterator cell =
2298 this->get_triangulation().end(level);
2299 if (cell.state() != IteratorState::valid)
2300 return end();
2301 return level_cell_iterator(*cell, this);
2302}
2303
2304
2305
2306template <int dim, int spacedim>
2307typename DoFHandler<dim, spacedim>::level_cell_iterator
2308DoFHandler<dim, spacedim>::end_mg() const
2309{
2310 return level_cell_iterator(&this->get_triangulation(), -1, -1, this);
2311}
2312
2313
2314
2315template <int dim, int spacedim>
2316IteratorRange<typename DoFHandler<dim, spacedim>::cell_iterator>
2317DoFHandler<dim, spacedim>::cell_iterators() const
2318{
2319 return IteratorRange<typename DoFHandler<dim, spacedim>::cell_iterator>(
2320 begin(), end());
2321}
2322
2323
2324
2325template <int dim, int spacedim>
2326IteratorRange<typename DoFHandler<dim, spacedim>::active_cell_iterator>
2327DoFHandler<dim, spacedim>::active_cell_iterators() const
2328{
2329 return IteratorRange<
2330 typename DoFHandler<dim, spacedim>::active_cell_iterator>(begin_active(),
2331 end());
2332}
2333
2334
2335
2336template <int dim, int spacedim>
2337IteratorRange<typename DoFHandler<dim, spacedim>::level_cell_iterator>
2338DoFHandler<dim, spacedim>::mg_cell_iterators() const
2339{
2340 return IteratorRange<typename DoFHandler<dim, spacedim>::level_cell_iterator>(
2341 begin_mg(), end_mg());
2342}
2343
2344
2345
2346template <int dim, int spacedim>
2347IteratorRange<typename DoFHandler<dim, spacedim>::cell_iterator>
2348DoFHandler<dim, spacedim>::cell_iterators_on_level(
2349 const unsigned int level) const
2350{
2351 return IteratorRange<typename DoFHandler<dim, spacedim>::cell_iterator>(
2352 begin(level), end(level));
2353}
2354
2355
2356
2357template <int dim, int spacedim>
2358IteratorRange<typename DoFHandler<dim, spacedim>::active_cell_iterator>
2359DoFHandler<dim, spacedim>::active_cell_iterators_on_level(
2360 const unsigned int level) const
2361{
2362 return IteratorRange<
2363 typename DoFHandler<dim, spacedim>::active_cell_iterator>(
2364 begin_active(level), end_active(level));
2365}
2366
2367
2368
2369template <int dim, int spacedim>
2370IteratorRange<typename DoFHandler<dim, spacedim>::level_cell_iterator>
2371DoFHandler<dim, spacedim>::mg_cell_iterators_on_level(
2372 const unsigned int level) const
2373{
2374 return IteratorRange<typename DoFHandler<dim, spacedim>::level_cell_iterator>(
2375 begin_mg(level), end_mg(level));
2376}
2377
2378
2379
2380//---------------------------------------------------------------------------
2381
2382
2383
2384template <int dim, int spacedim>
2385types::global_dof_index
2386DoFHandler<dim, spacedim>::n_boundary_dofs() const
2387{
2388 Assert(!(dim == 2 && spacedim == 3) || hp_capability_enabled == false,
2389 ExcNotImplementedWithHP());
2390
2391 Assert(this->fe_collection.size() > 0, ExcNoFESelected());
2392
2393 std::unordered_set<types::global_dof_index> boundary_dofs;
2394 std::vector<types::global_dof_index> dofs_on_face;
2395 dofs_on_face.reserve(this->get_fe_collection().max_dofs_per_face());
2396
2397 const IndexSet &owned_dofs = locally_owned_dofs();
2398
2399 // loop over all faces to check whether they are at a
2400 // boundary. note that we need not take special care of single
2401 // lines in 3d (using @p{cell->has_boundary_lines}), since we do
2402 // not support boundaries of dimension dim-2, and so every
2403 // boundary line is also part of a boundary face.
2404 for (const auto &cell : this->active_cell_iterators())
2405 if (cell->is_locally_owned() && cell->at_boundary())
2406 {
2407 for (const auto iface : cell->face_indices())
2408 {
2409 const auto face = cell->face(iface);
2410 if (face->at_boundary())
2411 {
2412 const unsigned int dofs_per_face =
2413 cell->get_fe().n_dofs_per_face(iface);
2414 dofs_on_face.resize(dofs_per_face);
2415
2416 face->get_dof_indices(dofs_on_face, cell->active_fe_index());
2417 for (unsigned int i = 0; i < dofs_per_face; ++i)
2418 {
2419 const unsigned int global_idof_index = dofs_on_face[i];
2420 if (owned_dofs.is_element(global_idof_index))
2421 {
2422 boundary_dofs.insert(global_idof_index);
2423 }
2424 }
2425 }
2426 }
2427 }
2428 return boundary_dofs.size();
2429}
2430
2431
2432
2433template <int dim, int spacedim>
2434types::global_dof_index
2435DoFHandler<dim, spacedim>::n_boundary_dofs(
2436 const std::set<types::boundary_id> &boundary_ids) const
2437{
2438 Assert(!(dim == 2 && spacedim == 3) || hp_capability_enabled == false,
2439 ExcNotImplementedWithHP());
2440
2441 Assert(this->fe_collection.size() > 0, ExcNoFESelected());
2442 Assert(boundary_ids.find(numbers::internal_face_boundary_id) ==
2443 boundary_ids.end(),
2444 ExcInvalidBoundaryIndicator());
2445
2446 // same as above, but with additional checks for set of boundary
2447 // indicators
2448 std::unordered_set<types::global_dof_index> boundary_dofs;
2449 std::vector<types::global_dof_index> dofs_on_face;
2450 dofs_on_face.reserve(this->get_fe_collection().max_dofs_per_face());
2451
2452 const IndexSet &owned_dofs = locally_owned_dofs();
2453
2454 for (const auto &cell : this->active_cell_iterators())
2455 if (cell->is_locally_owned() && cell->at_boundary())
2456 {
2457 for (const auto iface : cell->face_indices())
2458 {
2459 const auto face = cell->face(iface);
2460 const unsigned int boundary_id = face->boundary_id();
2461 if (face->at_boundary() &&
2462 (boundary_ids.find(boundary_id) != boundary_ids.end()))
2463 {
2464 const unsigned int dofs_per_face =
2465 cell->get_fe().n_dofs_per_face(iface);
2466 dofs_on_face.resize(dofs_per_face);
2467
2468 face->get_dof_indices(dofs_on_face, cell->active_fe_index());
2469 for (unsigned int i = 0; i < dofs_per_face; ++i)
2470 {
2471 const unsigned int global_idof_index = dofs_on_face[i];
2472 if (owned_dofs.is_element(global_idof_index))
2473 {
2474 boundary_dofs.insert(global_idof_index);
2475 }
2476 }
2477 }
2478 }
2479 }
2480 return boundary_dofs.size();
2481}
2482
2483
2484
2485template <int dim, int spacedim>
2486std::size_t
2487DoFHandler<dim, spacedim>::memory_consumption() const
2488{
2489 std::size_t mem = MemoryConsumption::memory_consumption(this->tria) +
2490 MemoryConsumption::memory_consumption(this->fe_collection) +
2491 MemoryConsumption::memory_consumption(this->number_cache);
2492
2493 mem += MemoryConsumption::memory_consumption(cell_dof_cache_indices) +
2494 MemoryConsumption::memory_consumption(cell_dof_cache_ptr) +
2495 MemoryConsumption::memory_consumption(object_dof_indices) +
2496 MemoryConsumption::memory_consumption(object_dof_ptr) +
2497 MemoryConsumption::memory_consumption(hp_object_fe_indices) +
2498 MemoryConsumption::memory_consumption(hp_object_fe_ptr) +
2499 MemoryConsumption::memory_consumption(hp_cell_active_fe_indices) +
2500 MemoryConsumption::memory_consumption(hp_cell_future_fe_indices);
2501
2502
2503 if (hp_capability_enabled)
2504 {
2505 // nothing to add
2506 }
2507 else
2508 {
2509 // collect size of multigrid data structures
2510
2511 mem += MemoryConsumption::memory_consumption(this->block_info_object);
2512
2513 for (unsigned int level = 0; level < this->mg_levels.size(); ++level)
2514 mem += this->mg_levels[level]->memory_consumption();
2515
2516 if (this->mg_faces != nullptr)
2517 mem += MemoryConsumption::memory_consumption(*this->mg_faces);
2518
2519 for (unsigned int i = 0; i < this->mg_vertex_dofs.size(); ++i)
2520 mem += sizeof(MGVertexDoFs) +
2521 (1 + this->mg_vertex_dofs[i].get_finest_level() -
2522 this->mg_vertex_dofs[i].get_coarsest_level()) *
2523 sizeof(types::global_dof_index);
2524 }
2525
2526 return mem;
2527}
2528
2529
2530
2531template <int dim, int spacedim>
2532void
2533DoFHandler<dim, spacedim>::set_fe(const FiniteElement<dim, spacedim> &fe)
2534{
2535 this->set_fe(hp::FECollection<dim, spacedim>(fe));
2536}
2537
2538
2539
2540template <int dim, int spacedim>
2541void
2542DoFHandler<dim, spacedim>::set_fe(const hp::FECollection<dim, spacedim> &ff)
2543{
2544 Assert(
2545 this->tria != nullptr,
2546 ExcMessage(
2547 "You need to set the Triangulation in the DoFHandler using reinit() or "
2548 "in the constructor before you can distribute DoFs."));
2549 Assert(this->tria->n_levels() > 0,
2550 ExcMessage("The Triangulation you are using is empty!"));
2551 Assert(ff.size() > 0, ExcMessage("The hp::FECollection given is empty!"));
2552
2553 // don't create a new object if the one we have is already appropriate
2554 if (this->fe_collection != ff)
2555 {
2556 this->fe_collection = hp::FECollection<dim, spacedim>(ff);
2557
2558 const bool contains_multiple_fes = (this->fe_collection.size() > 1);
2559
2560 // disable hp-mode if only a single finite element has been registered
2561 if (hp_capability_enabled && !contains_multiple_fes)
2562 {
2563 hp_capability_enabled = false;
2564
2565 // unsubscribe connections to signals that are only relevant for
2566 // hp-mode, since we only have a single element here
2567 for (auto &connection : this->tria_listeners_for_transfer)
2568 connection.disconnect();
2569 this->tria_listeners_for_transfer.clear();
2570
2571 // release active and future finite element tables
2572 this->hp_cell_active_fe_indices.clear();
2573 this->hp_cell_active_fe_indices.shrink_to_fit();
2574 this->hp_cell_future_fe_indices.clear();
2575 this->hp_cell_future_fe_indices.shrink_to_fit();
2576 }
2577
2578 // re-enabling hp-mode is not permitted since the active and future FE
2579 // tables are no longer available
2580 AssertThrow(
2581 hp_capability_enabled || !contains_multiple_fes,
2582 ExcMessage(
2583 "You cannot re-enable hp-capabilities after you registered a single "
2584 "finite element. Please create a new DoFHandler object instead."));
2585 }
2586
2587 if (hp_capability_enabled)
2588 {
2589 // make sure every processor knows the active FE indices
2590 // on both its own cells and all ghost cells
2591 ::internal::hp::DoFHandlerImplementation::Implementation::
2592 communicate_active_fe_indices(*this);
2593
2594 // make sure that the FE collection is large enough to
2595 // cover all FE indices presently in use on the mesh
2596 for (const auto &cell : this->active_cell_iterators())
2597 if (!cell->is_artificial())
2598 Assert(cell->active_fe_index() < this->fe_collection.size(),
2599 ExcInvalidFEIndex(cell->active_fe_index(),
2600 this->fe_collection.size()));
2601 }
2602}
2603
2604
2605
2606template <int dim, int spacedim>
2607void
2608DoFHandler<dim, spacedim>::distribute_dofs(
2609 const FiniteElement<dim, spacedim> &fe)
2610{
2611 this->distribute_dofs(hp::FECollection<dim, spacedim>(fe));
2612}
2613
2614
2615
2616template <int dim, int spacedim>
2617void
2618DoFHandler<dim, spacedim>::distribute_dofs(
2619 const hp::FECollection<dim, spacedim> &ff)
2620{
2621 Assert(
2622 this->tria != nullptr,
2623 ExcMessage(
2624 "You need to set the Triangulation in the DoFHandler using reinit() or "
2625 "in the constructor before you can distribute DoFs."));
2626 Assert(this->tria->n_levels() > 0,
2627 ExcMessage("The Triangulation you are using is empty!"));
2628 Assert(ff.size() > 0, ExcMessage("The hp::FECollection given is empty!"));
2629
2630 //
2631 // register the new finite element collection
2632 //
2633 // don't create a new object if the one we have is identical
2634 if (this->fe_collection != ff)
2635 {
2636 this->fe_collection = hp::FECollection<dim, spacedim>(ff);
2637
2638 const bool contains_multiple_fes = (this->fe_collection.size() > 1);
2639
2640 // disable hp-mode if only a single finite element has been registered
2641 if (hp_capability_enabled && !contains_multiple_fes)
2642 {
2643 hp_capability_enabled = false;
2644
2645 // unsubscribe connections to signals that are only relevant for
2646 // hp-mode, since we only have a single element here
2647 for (auto &connection : this->tria_listeners_for_transfer)
2648 connection.disconnect();
2649 this->tria_listeners_for_transfer.clear();
2650
2651 // release active and future finite element tables
2652 this->hp_cell_active_fe_indices.clear();
2653 this->hp_cell_active_fe_indices.shrink_to_fit();
2654 this->hp_cell_future_fe_indices.clear();
2655 this->hp_cell_future_fe_indices.shrink_to_fit();
2656 }
2657
2658 // re-enabling hp-mode is not permitted since the active and future FE
2659 // tables are no longer available
2660 AssertThrow(
2661 hp_capability_enabled || !contains_multiple_fes,
2662 ExcMessage(
2663 "You cannot re-enable hp-capabilities after you registered a single "
2664 "finite element. Please call reinit() or create a new DoFHandler "
2665 "object instead."));
2666 }
2667
2668 //
2669 // enumerate all degrees of freedom
2670 //
2671 if (hp_capability_enabled)
2672 {
2673 // make sure every processor knows the active FE indices
2674 // on both its own cells and all ghost cells
2675 internal::hp::DoFHandlerImplementation::Implementation::
2676 communicate_active_fe_indices(*this);
2677
2678#ifdef DEBUG
2679 // make sure that the FE collection is large enough to
2680 // cover all FE indices presently in use on the mesh
2681 for (const auto &cell : this->active_cell_iterators())
2682 {
2683 if (!cell->is_artificial())
2684 Assert(cell->active_fe_index() < this->fe_collection.size(),
2685 ExcInvalidFEIndex(cell->active_fe_index(),
2686 this->fe_collection.size()));
2687 if (cell->is_locally_owned())
2688 Assert(cell->future_fe_index() < this->fe_collection.size(),
2689 ExcInvalidFEIndex(cell->future_fe_index(),
2690 this->fe_collection.size()));
2691 }
2692#endif
2693 }
2694
2695 {
2696 // We would like to enumerate all dofs for shared::Triangulations. If an
2697 // underlying shared::Tria allows artificial cells, we need to restore the
2698 // true cell owners temporarily.
2699 // We use the TemporarilyRestoreSubdomainIds class for this purpose: we save
2700 // the current set of subdomain ids, set subdomain ids to the "true" owner
2701 // of each cell upon construction of the TemporarilyRestoreSubdomainIds
2702 // object, and later restore these flags when it is destroyed.
2703 const internal::parallel::shared::TemporarilyRestoreSubdomainIds<dim,
2704 spacedim>
2705 subdomain_modifier(this->get_triangulation());
2706
2707 // Adjust size of levels to the triangulation. Note that we still have to
2708 // allocate space for all degrees of freedom on this mesh (including ghost
2709 // and cells that are entirely stored on different processors), though we
2710 // may not assign numbers to some of them (i.e. they will remain at
2711 // invalid_dof_index). We need to allocate the space because we will want
2712 // to be able to query the dof_indices on each cell, and simply be told
2713 // that we don't know them on some cell (i.e. get back invalid_dof_index)
2714 if (hp_capability_enabled)
2715 internal::hp::DoFHandlerImplementation::Implementation::reserve_space(
2716 *this);
2717 else
2718 internal::DoFHandlerImplementation::Implementation::reserve_space(*this);
2719 }
2720
2721 // hand the actual work over to the policy
2722 this->number_cache = this->policy->distribute_dofs();
2723
2724 // do some housekeeping: compress indices
2725 // if(hp_capability_enabled)
2726 // {
2727 // Threads::TaskGroup<> tg;
2728 // for (int level = this->levels_hp.size() - 1; level >= 0; --level)
2729 // tg += Threads::new_task(
2730 // &::internal::hp::DoFLevel::compress_data<dim, spacedim>,
2731 // *this->levels_hp[level],
2732 // this->fe_collection);
2733 // tg.join_all();
2734 // }
2735
2736 // Initialize the block info object only if this is a sequential
2737 // triangulation. It doesn't work correctly yet if it is parallel and has not
2738 // yet been implemented for hp-mode.
2739 if (!hp_capability_enabled &&
2740 dynamic_cast<const parallel::DistributedTriangulationBase<dim, spacedim>
2741 *>(&*this->tria) == nullptr)
2742 this->block_info_object.initialize(*this, false, true);
2743}
2744
2745
2746
2747template <int dim, int spacedim>
2748void
2749DoFHandler<dim, spacedim>::distribute_mg_dofs()
2750{
2751 AssertThrow(hp_capability_enabled == false, ExcNotImplementedWithHP());
2752
2753 Assert(
2754 this->object_dof_indices.size() > 0,
2755 ExcMessage(
2756 "Distribute active DoFs using distribute_dofs() before calling distribute_mg_dofs()."));
2757
2758 Assert(
2759 ((this->tria->get_mesh_smoothing() &
2760 Triangulation<dim, spacedim>::limit_level_difference_at_vertices) !=
2761 Triangulation<dim, spacedim>::none),
2762 ExcMessage(
2763 "The mesh smoothing requirement 'limit_level_difference_at_vertices' has to be set for using multigrid!"));
2764
2765 this->clear_mg_space();
2766
2767 internal::DoFHandlerImplementation::Implementation::reserve_space_mg(*this);
2768 this->mg_number_cache = this->policy->distribute_mg_dofs();
2769
2770 // initialize the block info object only if this is a sequential
2771 // triangulation. it doesn't work correctly yet if it is parallel
2772 if (dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
2773 &*this->tria) == nullptr)
2774 this->block_info_object.initialize(*this, true, false);
2775}
2776
2777
2778
2779template <int dim, int spacedim>
2780void
2781DoFHandler<dim, spacedim>::initialize_local_block_info()
2782{
2783 AssertThrow(hp_capability_enabled == false, ExcNotImplementedWithHP());
2784
2785 this->block_info_object.initialize_local(*this);
2786}
2787
2788
2789
2790template <int dim, int spacedim>
2791void
2792DoFHandler<dim, spacedim>::setup_policy()
2793{
2794 // decide whether we need a sequential or a parallel distributed policy
2795 if (dynamic_cast<const ::parallel::shared::Triangulation<dim, spacedim>
2796 *>(&this->get_triangulation()) != nullptr)
2797 this->policy = std::make_unique<internal::DoFHandlerImplementation::Policy::
2798 ParallelShared<dim, spacedim>>(*this);
2799 else if (dynamic_cast<
2800 const ::parallel::DistributedTriangulationBase<dim, spacedim>
2801 *>(&this->get_triangulation()) == nullptr)
2802 this->policy = std::make_unique<
2803 internal::DoFHandlerImplementation::Policy::Sequential<dim, spacedim>>(
2804 *this);
2805 else
2806 this->policy =
2807 std::make_unique<internal::DoFHandlerImplementation::Policy::
2808 ParallelDistributed<dim, spacedim>>(*this);
2809}
2810
2811
2812
2813template <int dim, int spacedim>
2814void
2815DoFHandler<dim, spacedim>::clear()
2816{
2817 // release memory
2818 this->clear_space();
2819 this->clear_mg_space();
2820}
2821
2822
2823
2824template <int dim, int spacedim>
2825void
2826DoFHandler<dim, spacedim>::clear_space()
2827{
2828 cell_dof_cache_indices.clear();
2829
2830 cell_dof_cache_ptr.clear();
2831
2832 object_dof_indices.clear();
2833
2834 object_dof_ptr.clear();
2835
2836 this->number_cache.clear();
2837
2838 this->hp_cell_active_fe_indices.clear();
2839 this->hp_cell_future_fe_indices.clear();
2840}
2841
2842
2843
2844template <int dim, int spacedim>
2845void
2846DoFHandler<dim, spacedim>::clear_mg_space()
2847{
2848 this->mg_levels.clear();
2849 this->mg_faces.reset();
2850
2851 std::vector<MGVertexDoFs> tmp;
2852
2853 std::swap(this->mg_vertex_dofs, tmp);
2854
2855 this->mg_number_cache.clear();
2856}
2857
2858
2859
2860template <int dim, int spacedim>
2861void
2862DoFHandler<dim, spacedim>::renumber_dofs(
2863 const std::vector<types::global_dof_index> &new_numbers)
2864{
2865 if (hp_capability_enabled)
2866 {
2867 Assert(this->hp_cell_future_fe_indices.size() > 0,
2868 ExcMessage(
2869 "You need to distribute DoFs before you can renumber them."));
2870
2871 AssertDimension(new_numbers.size(), this->n_locally_owned_dofs());
2872
2873#ifdef DEBUG
2874 // assert that the new indices are consecutively numbered if we are
2875 // working on a single processor. this doesn't need to
2876 // hold in the case of a parallel mesh since we map the interval
2877 // [0...n_dofs()) into itself but only globally, not on each processor
2878 if (this->n_locally_owned_dofs() == this->n_dofs())
2879 {
2880 std::vector<types::global_dof_index> tmp(new_numbers);
2881 std::sort(tmp.begin(), tmp.end());
2882 std::vector<types::global_dof_index>::const_iterator p = tmp.begin();
2883 types::global_dof_index i = 0;
2884 for (; p != tmp.end(); ++p, ++i)
2885 Assert(*p == i, ExcNewNumbersNotConsecutive(i));
2886 }
2887 else
2888 for (const auto new_number : new_numbers)
2889 Assert(new_number < this->n_dofs(),
2890 ExcMessage(
2891 "New DoF index is not less than the total number of dofs."));
2892#endif
2893
2894 // uncompress the internal storage scheme of dofs on cells so that
2895 // we can access dofs in turns. uncompress in parallel, starting
2896 // with the most expensive levels (the highest ones)
2897 //{
2898 // Threads::TaskGroup<> tg;
2899 // for (int level = this->levels_hp.size() - 1; level >= 0; --level)
2900 // tg += Threads::new_task(
2901 // &::internal::hp::DoFLevel::uncompress_data<dim, spacedim>,
2902 // *this->levels_hp[level],
2903 // this->fe_collection);
2904 // tg.join_all();
2905 //}
2906
2907 // do the renumbering
2908 this->number_cache = this->policy->renumber_dofs(new_numbers);
2909
2910 // now re-compress the dof indices
2911 //{
2912 // Threads::TaskGroup<> tg;
2913 // for (int level = this->levels_hp.size() - 1; level >= 0; --level)
2914 // tg += Threads::new_task(
2915 // &::internal::hp::DoFLevel::compress_data<dim, spacedim>,
2916 // *this->levels_hp[level],
2917 // this->fe_collection);
2918 // tg.join_all();
2919 //}
2920 }
2921 else
2922 {
2923 Assert(this->object_dof_indices.size() > 0,
2924 ExcMessage(
2925 "You need to distribute DoFs before you can renumber them."));
2926
2927#ifdef DEBUG
2928 if (dynamic_cast<const parallel::shared::Triangulation<dim, spacedim> *>(
2929 &*this->tria) != nullptr)
2930 {
2931 Assert(new_numbers.size() == this->n_dofs() ||
2932 new_numbers.size() == this->n_locally_owned_dofs(),
2933 ExcMessage("Incorrect size of the input array."));
2934 }
2935 else if (dynamic_cast<
2936 const parallel::DistributedTriangulationBase<dim, spacedim> *>(
2937 &*this->tria) != nullptr)
2938 {
2939 AssertDimension(new_numbers.size(), this->n_locally_owned_dofs());
2940 }
2941 else
2942 {
2943 AssertDimension(new_numbers.size(), this->n_dofs());
2944 }
2945
2946 // assert that the new indices are consecutively numbered if we are
2947 // working on a single processor. this doesn't need to
2948 // hold in the case of a parallel mesh since we map the interval
2949 // [0...n_dofs()) into itself but only globally, not on each processor
2950 if (this->n_locally_owned_dofs() == this->n_dofs())
2951 {
2952 std::vector<types::global_dof_index> tmp(new_numbers);
2953 std::sort(tmp.begin(), tmp.end());
2954 std::vector<types::global_dof_index>::const_iterator p = tmp.begin();
2955 types::global_dof_index i = 0;
2956 for (; p != tmp.end(); ++p, ++i)
2957 Assert(*p == i, ExcNewNumbersNotConsecutive(i));
2958 }
2959 else
2960 for (const auto new_number : new_numbers)
2961 Assert(new_number < this->n_dofs(),
2962 ExcMessage(
2963 "New DoF index is not less than the total number of dofs."));
2964#endif
2965
2966 this->number_cache = this->policy->renumber_dofs(new_numbers);
2967 }
2968}
2969
2970
2971
2972template <int dim, int spacedim>
2973void
2974DoFHandler<dim, spacedim>::renumber_dofs(
2975 const unsigned int level,
2976 const std::vector<types::global_dof_index> &new_numbers)
2977{
2978 AssertThrow(hp_capability_enabled == false, ExcNotImplementedWithHP());
2979
2980 Assert(
2981 this->mg_levels.size() > 0 && this->object_dof_indices.size() > 0,
2982 ExcMessage(
2983 "You need to distribute active and level DoFs before you can renumber level DoFs."));
2984 AssertIndexRange(level, this->get_triangulation().n_global_levels());
2985 AssertDimension(new_numbers.size(),
2986 this->locally_owned_mg_dofs(level).n_elements());
2987
2988#ifdef DEBUG
2989 // assert that the new indices are consecutively numbered if we are working
2990 // on a single processor. this doesn't need to hold in the case of a
2991 // parallel mesh since we map the interval [0...n_dofs(level)) into itself
2992 // but only globally, not on each processor
2993 if (this->n_locally_owned_dofs() == this->n_dofs())
2994 {
2995 std::vector<types::global_dof_index> tmp(new_numbers);
2996 std::sort(tmp.begin(), tmp.end());
2997 std::vector<types::global_dof_index>::const_iterator p = tmp.begin();
2998 types::global_dof_index i = 0;
2999 for (; p != tmp.end(); ++p, ++i)
3000 Assert(*p == i, ExcNewNumbersNotConsecutive(i));
3001 }
3002 else
3003 for (const auto new_number : new_numbers)
3004 Assert(new_number < this->n_dofs(level),
3005 ExcMessage(
3006 "New DoF index is not less than the total number of dofs."));
3007#endif
3008
3009 this->mg_number_cache[level] =
3010 this->policy->renumber_mg_dofs(level, new_numbers);
3011}
3012
3013
3014
3015template <int dim, int spacedim>
3016unsigned int
3017DoFHandler<dim, spacedim>::max_couplings_between_boundary_dofs() const
3018{
3019 Assert(this->fe_collection.size() > 0, ExcNoFESelected());
3020
3021 switch (dim)
3022 {
3023 case 1:
3024 return this->fe_collection.max_dofs_per_vertex();
3025 case 2:
3026 return (3 * this->fe_collection.max_dofs_per_vertex() +
3027 2 * this->fe_collection.max_dofs_per_line());
3028 case 3:
3029 // we need to take refinement of one boundary face into
3030 // consideration here; in fact, this function returns what
3031 // #max_coupling_between_dofs<2> returns
3032 //
3033 // we assume here, that only four faces meet at the boundary;
3034 // this assumption is not justified and needs to be fixed some
3035 // time. fortunately, omitting it for now does no harm since
3036 // the matrix will cry foul if its requirements are not
3037 // satisfied
3038 return (19 * this->fe_collection.max_dofs_per_vertex() +
3039 28 * this->fe_collection.max_dofs_per_line() +
3040 8 * this->fe_collection.max_dofs_per_quad());
3041 default:
3042 Assert(false, ExcNotImplemented());
3043 return 0;
3044 }
3045}
3046
3047
3048
3049template <int dim, int spacedim>
3050unsigned int
3051DoFHandler<dim, spacedim>::max_couplings_between_dofs() const
3052{
3053 Assert(this->fe_collection.size() > 0, ExcNoFESelected());
3054 return internal::DoFHandlerImplementation::Implementation::
3055 max_couplings_between_dofs(*this);
3056}
3057
3058
3059
3060template <int dim, int spacedim>
3061template <int structdim>
3062types::global_dof_index
3063DoFHandler<dim, spacedim>::get_dof_index(const unsigned int obj_level,
3064 const unsigned int obj_index,
3065 const unsigned int fe_index,
3066 const unsigned int local_index) const
3067{
3068 if (hp_capability_enabled)
3069 {
3070 Assert(false, ExcNotImplemented());
3071 return numbers::invalid_dof_index;
3072 }
3073 else
3074 {
3075 return internal::DoFHandlerImplementation::Implementation::get_dof_index(
3076 *this,
3077 this->mg_levels[obj_level],
3078 this->mg_faces,
3079 obj_index,
3080 fe_index,
3081 local_index,
3082 std::integral_constant<int, structdim>());
3083 }
3084}
3085
3086
3087
3088template <int dim, int spacedim>
3089template <int structdim>
3090void
3091DoFHandler<dim, spacedim>::set_dof_index(
3092 const unsigned int obj_level,
3093 const unsigned int obj_index,
3094 const unsigned int fe_index,
3095 const unsigned int local_index,
3096 const types::global_dof_index global_index) const
3097{
3098 if (hp_capability_enabled)
3099 {
3100 Assert(false, ExcNotImplemented());
3101 return;
3102 }
3103 else
3104 {
3105 internal::DoFHandlerImplementation::Implementation::set_dof_index(
3106 *this,
3107 this->mg_levels[obj_level],
3108 this->mg_faces,
3109 obj_index,
3110 fe_index,
3111 local_index,
3112 global_index,
3113 std::integral_constant<int, structdim>());
3114 }
3115}
3116
3117
3118
3119template <int dim, int spacedim>
3120void
3121DoFHandler<dim, spacedim>::set_active_fe_indices(
3122 const std::vector<unsigned int> &active_fe_indices)
3123{
3124 Assert(active_fe_indices.size() == this->get_triangulation().n_active_cells(),
3125 ExcDimensionMismatch(active_fe_indices.size(),
3126 this->get_triangulation().n_active_cells()));
3127
3128 this->create_active_fe_table();
3129 // we could set the values directly, since they are stored as
3130 // protected data of this object, but for simplicity we use the
3131 // cell-wise access. this way we also have to pass some debug-mode
3132 // tests which we would have to duplicate ourselves otherwise
3133 for (const auto &cell : this->active_cell_iterators())
3134 if (cell->is_locally_owned())
3135 cell->set_active_fe_index(active_fe_indices[cell->active_cell_index()]);
3136}
3137
3138
3139
3140template <int dim, int spacedim>
3141void
3142DoFHandler<dim, spacedim>::get_active_fe_indices(
3143 std::vector<unsigned int> &active_fe_indices) const
3144{
3145 active_fe_indices.resize(this->get_triangulation().n_active_cells());
3146
3147 // we could try to extract the values directly, since they are
3148 // stored as protected data of this object, but for simplicity we
3149 // use the cell-wise access.
3150 for (const auto &cell : this->active_cell_iterators())
3151 if (!cell->is_artificial())
3152 active_fe_indices[cell->active_cell_index()] = cell->active_fe_index();
3153}
3154
3155
3156
3157template <int dim, int spacedim>
3158void
3159DoFHandler<dim, spacedim>::connect_to_triangulation_signals()
3160{
3161 // make sure this is called during initialization in hp-mode
3162 Assert(hp_capability_enabled, ExcOnlyAvailableWithHP());
3163
3164 // connect functions to signals of the underlying triangulation
3165 this->tria_listeners.push_back(this->tria->signals.create.connect(
3166 [this]() { this->reinit(*(this->tria)); }));
3167 this->tria_listeners.push_back(
3168 this->tria->signals.clear.connect([this]() { this->clear(); }));
3169
3170 // attach corresponding callback functions dealing with the transfer of
3171 // active FE indices depending on the type of triangulation
3172 if (dynamic_cast<
3173 const ::parallel::fullydistributed::Triangulation<dim, spacedim>
3174 *>(&this->get_triangulation()))
3175 {
3176 // no transfer of active FE indices for this class
3177 }
3178 else if (dynamic_cast<
3179 const ::parallel::distributed::Triangulation<dim, spacedim>
3180 *>(&this->get_triangulation()))
3181 {
3182 // repartitioning signals
3183 this->tria_listeners_for_transfer.push_back(
3184 this->tria->signals.pre_distributed_repartition.connect([this]() {
3185 internal::hp::DoFHandlerImplementation::Implementation::
3186 ensure_absence_of_future_fe_indices<dim, spacedim>(*this);
3187 }));
3188 this->tria_listeners_for_transfer.push_back(
3189 this->tria->signals.pre_distributed_repartition.connect(
3190 [this]() { this->pre_distributed_transfer_action(); }));
3191 this->tria_listeners_for_transfer.push_back(
3192 this->tria->signals.post_distributed_repartition.connect(
3193 [this] { this->post_distributed_transfer_action(); }));
3194
3195 // refinement signals
3196 this->tria_listeners_for_transfer.push_back(
3197 this->tria->signals.post_p4est_refinement.connect(
3198 [this]() { this->pre_distributed_transfer_action(); }));
3199 this->tria_listeners_for_transfer.push_back(
3200 this->tria->signals.post_distributed_refinement.connect(
3201 [this]() { this->post_distributed_transfer_action(); }));
3202
3203 // serialization signals
3204 this->tria_listeners_for_transfer.push_back(
3205 this->tria->signals.post_distributed_save.connect(
3206 [this]() { this->active_fe_index_transfer.reset(); }));
3207 this->tria_listeners_for_transfer.push_back(
3208 this->tria->signals.post_distributed_load.connect(
3209 [this]() { this->update_active_fe_table(); }));
3210 }
3211 else if (dynamic_cast<
3212 const ::parallel::shared::Triangulation<dim, spacedim> *>(
3213 &this->get_triangulation()) != nullptr)
3214 {
3215 // partitioning signals
3216 this->tria_listeners_for_transfer.push_back(
3217 this->tria->signals.pre_partition.connect([this]() {
3218 internal::hp::DoFHandlerImplementation::Implementation::
3219 ensure_absence_of_future_fe_indices(*this);
3220 }));
3221
3222 // refinement signals
3223 this->tria_listeners_for_transfer.push_back(
3224 this->tria->signals.pre_refinement.connect([this]() {
3225 internal::hp::DoFHandlerImplementation::Implementation::
3226 communicate_future_fe_indices(*this);
3227 }));
3228 this->tria_listeners_for_transfer.push_back(
3229 this->tria->signals.pre_refinement.connect(
3230 [this] { this->pre_transfer_action(); }));
3231 this->tria_listeners_for_transfer.push_back(
3232 this->tria->signals.post_refinement.connect(
3233 [this] { this->post_transfer_action(); }));
3234 }
3235 else
3236 {
3237 // refinement signals
3238 this->tria_listeners_for_transfer.push_back(
3239 this->tria->signals.pre_refinement.connect(
3240 [this] { this->pre_transfer_action(); }));
3241 this->tria_listeners_for_transfer.push_back(
3242 this->tria->signals.post_refinement.connect(
3243 [this] { this->post_transfer_action(); }));
3244 }
3245}
3246
3247
3248
3249template <int dim, int spacedim>
3250void
3251DoFHandler<dim, spacedim>::create_active_fe_table()
3252{
3253 AssertThrow(hp_capability_enabled == true, ExcOnlyAvailableWithHP());
3254
3255
3256 // Create sufficiently many hp::DoFLevels.
3257 // while (this->levels_hp.size() < this->tria->n_levels())
3258 // this->levels_hp.emplace_back(new ::internal::hp::DoFLevel);
3259
3260 this->hp_cell_active_fe_indices.resize(this->tria->n_levels());
3261 this->hp_cell_future_fe_indices.resize(this->tria->n_levels());
3262
3263 // then make sure that on each level we have the appropriate size
3264 // of active FE indices; preset them to zero, i.e. the default FE
3265 for (unsigned int level = 0; level < this->hp_cell_future_fe_indices.size();
3266 ++level)
3267 {
3268 if (this->hp_cell_active_fe_indices[level].size() == 0 &&
3269 this->hp_cell_future_fe_indices[level].size() == 0)
3270 {
3271 this->hp_cell_active_fe_indices[level].resize(
3272 this->tria->n_raw_cells(level), 0);
3273 this->hp_cell_future_fe_indices[level].resize(
3274 this->tria->n_raw_cells(level), invalid_active_fe_index);
3275 }
3276 else
3277 {
3278 // Either the active FE indices have size zero because
3279 // they were just created, or the correct size. Other
3280 // sizes indicate that something went wrong.
3281 Assert(this->hp_cell_active_fe_indices[level].size() ==
3282 this->tria->n_raw_cells(level) &&
3283 this->hp_cell_future_fe_indices[level].size() ==
3284 this->tria->n_raw_cells(level),
3285 ExcInternalError());
3286 }
3287
3288 // it may be that the previous table was compressed; in that
3289 // case, restore the correct active FE index. the fact that
3290 // this no longer matches the indices in the table is of no
3291 // importance because the current function is called at a
3292 // point where we have to recreate the dof_indices tables in
3293 // the levels anyway
3294 // this->levels_hp[level]->normalize_active_fe_indices();
3295 }
3296}
3297
3298
3299
3300template <int dim, int spacedim>
3301void
3302DoFHandler<dim, spacedim>::update_active_fe_table()
3303{
3304 // // Normally only one level is added, but if this Triangulation
3305 // // is created by copy_triangulation, it can be more than one level.
3306 // while (this->levels_hp.size() < this->tria->n_levels())
3307 // this->levels_hp.emplace_back(new ::internal::hp::DoFLevel);
3308 //
3309 // // Coarsening can lead to the loss of levels. Hence remove them.
3310 // while (this->levels_hp.size() > this->tria->n_levels())
3311 // {
3312 // // drop the last element. that also releases the memory pointed to
3313 // this->levels_hp.pop_back();
3314 // }
3315
3316 this->hp_cell_active_fe_indices.resize(this->tria->n_levels());
3317 this->hp_cell_active_fe_indices.shrink_to_fit();
3318
3319 this->hp_cell_future_fe_indices.resize(this->tria->n_levels());
3320 this->hp_cell_future_fe_indices.shrink_to_fit();
3321
3322 for (unsigned int i = 0; i < this->hp_cell_future_fe_indices.size(); ++i)
3323 {
3324 // Resize active FE indices vectors. Use zero indicator to extend.
3325 this->hp_cell_active_fe_indices[i].resize(this->tria->n_raw_cells(i), 0);
3326
3327 // Resize future FE indices vectors. Make sure that all
3328 // future FE indices have been cleared after refinement happened.
3329 //
3330 // We have used future FE indices to update all active FE indices
3331 // before refinement happened, thus we are safe to clear them now.
3332 this->hp_cell_future_fe_indices[i].assign(this->tria->n_raw_cells(i),
3333 invalid_active_fe_index);
3334 }
3335}
3336
3337
3338template <int dim, int spacedim>
3339void
3340DoFHandler<dim, spacedim>::pre_transfer_action()
3341{
3342 Assert(this->active_fe_index_transfer == nullptr, ExcInternalError());
3343
3344 this->active_fe_index_transfer = std::make_unique<ActiveFEIndexTransfer>();
3345
3346 ::internal::hp::DoFHandlerImplementation::Implementation::
3347 collect_fe_indices_on_cells_to_be_refined(*this);
3348}
3349
3350
3351
3352template <int dim, int spacedim>
3353void
3354DoFHandler<dim, spacedim>::pre_distributed_transfer_action()
3355{
3356#ifndef DEAL_II_WITH_P4EST
3357 Assert(false,
3358 ExcMessage(
3359 "You are attempting to use a functionality that is only available "
3360 "if deal.II was configured to use p4est, but cmake did not find a "
3361 "valid p4est library."));
3362#else
3363 // the implementation below requires a p:d:T currently
3364 Assert(
3365 (dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim> *>(
3366 &this->get_triangulation()) != nullptr),
3367 ExcNotImplemented());
3368
3369 Assert(active_fe_index_transfer == nullptr, ExcInternalError());
3370
3371 active_fe_index_transfer = std::make_unique<ActiveFEIndexTransfer>();
3372
3373 // If we work on a p::d::Triangulation, we have to transfer all
3374 // active FE indices since ownership of cells may change. We will
3375 // use our p::d::CellDataTransfer member to achieve this. Further,
3376 // we prepare the values in such a way that they will correspond to
3377 // the active FE indices on the new mesh.
3378
3379 // Gather all current future FE indices.
3380 active_fe_index_transfer->active_fe_indices.resize(
3381 get_triangulation().n_active_cells(), numbers::invalid_unsigned_int);
3382
3383 for (const auto &cell : active_cell_iterators())
3384 if (cell->is_locally_owned())
3385 active_fe_index_transfer->active_fe_indices[cell->active_cell_index()] =
3386 cell->future_fe_index();
3387
3388 // Create transfer object and attach to it.
3389 const auto *distributed_tria =
3390 dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim> *>(
3391 &this->get_triangulation());
3392
3393 active_fe_index_transfer->cell_data_transfer = std::make_unique<
3394 parallel::distributed::
3395 CellDataTransfer<dim, spacedim, std::vector<unsigned int>>>(
3396 *distributed_tria,
3397 /*transfer_variable_size_data=*/false,
3398 /*refinement_strategy=*/
3399 &::AdaptationStrategies::Refinement::
3400 preserve<dim, spacedim, unsigned int>,
3401 /*coarsening_strategy=*/
3402 [this](const typename Triangulation<dim, spacedim>::cell_iterator &parent,
3403 const std::vector<unsigned int> &children_fe_indices)
3404 -> unsigned int {
3405 return ::internal::hp::DoFHandlerImplementation::Implementation::
3406 determine_fe_from_children<dim, spacedim>(parent,
3407 children_fe_indices,
3408 fe_collection);
3409 });
3410
3411 active_fe_index_transfer->cell_data_transfer
3412 ->prepare_for_coarsening_and_refinement(
3413 active_fe_index_transfer->active_fe_indices);
3414#endif
3415}
3416
3417
3418
3419template <int dim, int spacedim>
3420void
3421DoFHandler<dim, spacedim>::post_transfer_action()
3422{
3423 update_active_fe_table();
3424
3425 Assert(this->active_fe_index_transfer != nullptr, ExcInternalError());
3426
3427 ::internal::hp::DoFHandlerImplementation::Implementation::
3428 distribute_fe_indices_on_refined_cells(*this);
3429
3430 // We have to distribute the information about active FE indices
3431 // of all cells (including the artificial ones) on all processors,
3432 // if a parallel::shared::Triangulation has been used.
3433 ::internal::hp::DoFHandlerImplementation::Implementation::
3434 communicate_active_fe_indices(*this);
3435
3436 // Free memory.
3437 this->active_fe_index_transfer.reset();
3438}
3439
3440
3441
3442template <int dim, int spacedim>
3443void
3444DoFHandler<dim, spacedim>::post_distributed_transfer_action()
3445{
3446#ifndef DEAL_II_WITH_P4EST
3447 Assert(false, ExcInternalError());
3448#else
3449 update_active_fe_table();
3450
3451 Assert(this->active_fe_index_transfer != nullptr, ExcInternalError());
3452
3453 // Unpack active FE indices.
3454 this->active_fe_index_transfer->active_fe_indices.resize(
3455 this->get_triangulation().n_active_cells(), numbers::invalid_unsigned_int);
3456 this->active_fe_index_transfer->cell_data_transfer->unpack(
3457 this->active_fe_index_transfer->active_fe_indices);
3458
3459 // Update all locally owned active FE indices.
3460 this->set_active_fe_indices(
3461 this->active_fe_index_transfer->active_fe_indices);
3462
3463 // Update active FE indices on ghost cells.
3464 ::internal::hp::DoFHandlerImplementation::Implementation::
3465 communicate_active_fe_indices(*this);
3466
3467 // Free memory.
3468 this->active_fe_index_transfer.reset();
3469#endif
3470}
3471
3472
3473
3474template <int dim, int spacedim>
3475void
3476DoFHandler<dim, spacedim>::prepare_for_serialization_of_active_fe_indices()
3477{
3478#ifndef DEAL_II_WITH_P4EST
3479 Assert(false,
3480 ExcMessage(
3481 "You are attempting to use a functionality that is only available "
3482 "if deal.II was configured to use p4est, but cmake did not find a "
3483 "valid p4est library."));
3484#else
3485 // the implementation below requires a p:d:T currently
3486 Assert(
3487 (dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim> *>(
3488 &this->get_triangulation()) != nullptr),
3489 ExcNotImplemented());
3490
3491 Assert(active_fe_index_transfer == nullptr, ExcInternalError());
3492
3493 active_fe_index_transfer = std::make_unique<ActiveFEIndexTransfer>();
3494
3495 // Create transfer object and attach to it.
3496 const auto *distributed_tria =
3497 dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim> *>(
3498 &this->get_triangulation());
3499
3500 active_fe_index_transfer->cell_data_transfer = std::make_unique<
3501 parallel::distributed::
3502 CellDataTransfer<dim, spacedim, std::vector<unsigned int>>>(
3503 *distributed_tria,
3504 /*transfer_variable_size_data=*/false,
3505 /*refinement_strategy=*/
3506 &::AdaptationStrategies::Refinement::
3507 preserve<dim, spacedim, unsigned int>,
3508 /*coarsening_strategy=*/
3509 [this](const typename Triangulation<dim, spacedim>::cell_iterator &parent,
3510 const std::vector<unsigned int> &children_fe_indices)
3511 -> unsigned int {
3512 return ::internal::hp::DoFHandlerImplementation::Implementation::
3513 determine_fe_from_children<dim, spacedim>(parent,
3514 children_fe_indices,
3515 fe_collection);
3516 });
3517
3518 // If we work on a p::d::Triangulation, we have to transfer all
3519 // active FE indices since ownership of cells may change.
3520
3521 // Gather all current active FE indices
3522 get_active_fe_indices(active_fe_index_transfer->active_fe_indices);
3523
3524 // Attach to transfer object
3525 active_fe_index_transfer->cell_data_transfer->prepare_for_serialization(
3526 active_fe_index_transfer->active_fe_indices);
3527#endif
3528}
3529
3530
3531
3532template <int dim, int spacedim>
3533void
3534DoFHandler<dim, spacedim>::deserialize_active_fe_indices()
3535{
3536#ifndef DEAL_II_WITH_P4EST
3537 Assert(false,
3538 ExcMessage(
3539 "You are attempting to use a functionality that is only available "
3540 "if deal.II was configured to use p4est, but cmake did not find a "
3541 "valid p4est library."));
3542#else
3543 // the implementation below requires a p:d:T currently
3544 Assert(
3545 (dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim> *>(
3546 &this->get_triangulation()) != nullptr),
3547 ExcNotImplemented());
3548
3549 Assert(active_fe_index_transfer == nullptr, ExcInternalError());
3550
3551 active_fe_index_transfer = std::make_unique<ActiveFEIndexTransfer>();
3552
3553 // Create transfer object and attach to it.
3554 const auto *distributed_tria =
3555 dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim> *>(
3556 &this->get_triangulation());
3557
3558 active_fe_index_transfer->cell_data_transfer = std::make_unique<
3559 parallel::distributed::
3560 CellDataTransfer<dim, spacedim, std::vector<unsigned int>>>(
3561 *distributed_tria,
3562 /*transfer_variable_size_data=*/false,
3563 /*refinement_strategy=*/
3564 &::AdaptationStrategies::Refinement::
3565 preserve<dim, spacedim, unsigned int>,
3566 /*coarsening_strategy=*/
3567 [this](const typename Triangulation<dim, spacedim>::cell_iterator &parent,
3568 const std::vector<unsigned int> &children_fe_indices)
3569 -> unsigned int {
3570 return ::internal::hp::DoFHandlerImplementation::Implementation::
3571 determine_fe_from_children<dim, spacedim>(parent,
3572 children_fe_indices,
3573 fe_collection);
3574 });
3575
3576 // Unpack active FE indices.
3577 active_fe_index_transfer->active_fe_indices.resize(
3578 get_triangulation().n_active_cells(), numbers::invalid_unsigned_int);
3579 active_fe_index_transfer->cell_data_transfer->deserialize(
3580 active_fe_index_transfer->active_fe_indices);
3581
3582 // Update all locally owned active FE indices.
3583 set_active_fe_indices(active_fe_index_transfer->active_fe_indices);
3584
3585 // Update active FE indices on ghost cells.
3586 ::internal::hp::DoFHandlerImplementation::Implementation::
3587 communicate_active_fe_indices(*this);
3588
3589 // Free memory.
3590 active_fe_index_transfer.reset();
3591#endif
3592}
3593
3594
3595
3596template <int dim, int spacedim>
3597DoFHandler<dim, spacedim>::MGVertexDoFs::MGVertexDoFs()
3598 : coarsest_level(numbers::invalid_unsigned_int)
3599 , finest_level(0)
3600{}
3601
3602
3603
3604template <int dim, int spacedim>
3605void
3606DoFHandler<dim, spacedim>::MGVertexDoFs::init(
3607 const unsigned int cl,
3608 const unsigned int fl,
3609 const unsigned int dofs_per_vertex)
3610{
3611 coarsest_level = cl;
3612 finest_level = fl;
3613
3614 if (coarsest_level <= finest_level)
3615 {
3616 const unsigned int n_levels = finest_level - coarsest_level + 1;
3617 const unsigned int n_indices = n_levels * dofs_per_vertex;
3618
3619 indices = std::make_unique<types::global_dof_index[]>(n_indices);
3620 std::fill(indices.get(),
3621 indices.get() + n_indices,
3622 numbers::invalid_dof_index);
3623 }
3624 else
3625 indices.reset();
3626}
3627
3628
3629
3630template <int dim, int spacedim>
3631unsigned int
3632DoFHandler<dim, spacedim>::MGVertexDoFs::get_coarsest_level() const
3633{
3634 return coarsest_level;
3635}
3636
3637
3638
3639template <int dim, int spacedim>
3640unsigned int
3641DoFHandler<dim, spacedim>::MGVertexDoFs::get_finest_level() const
3642{
3643 return finest_level;
3644}
3645
3646/*-------------- Explicit Instantiations -------------------------------*/
3647#include "dof_handler.inst"
3648
3649
3650
3651DEAL_II_NAMESPACE_CLOSE
DoFHandler::MGVertexDoFs::init
void init(const unsigned int coarsest_level, const unsigned int finest_level, const unsigned int dofs_per_vertex)
DoFHandler::MGVertexDoFs::get_finest_level
unsigned int get_finest_level() const
DoFHandler::MGVertexDoFs::get_coarsest_level
unsigned int get_coarsest_level() const
DoFHandler::end
cell_iterator end() const
DoFHandler::pre_distributed_transfer_action
void pre_distributed_transfer_action()
DoFHandler::post_transfer_action
void post_transfer_action()
DoFHandler::set_dof_index
void set_dof_index(const unsigned int obj_level, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index) const
DoFHandler::memory_consumption
virtual std::size_t memory_consumption() const
DoFHandler::hp_cell_active_fe_indices
std::vector< std::vector< active_fe_index_type > > hp_cell_active_fe_indices
Definition: dof_handler.h:1626
DoFHandler::max_couplings_between_dofs
unsigned int max_couplings_between_dofs() const
DoFHandler::mg_levels
std::vector< std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel< dim > > > mg_levels
Definition: dof_handler.h:1646
DoFHandler::get_fe
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
DoFHandler::active_fe_index_transfer
std::unique_ptr< ActiveFEIndexTransfer > active_fe_index_transfer
Definition: dof_handler.h:1658
DoFHandler::fe_collection
hp::FECollection< dim, spacedim > fe_collection
Definition: dof_handler.h:1550
DoFHandler::create_active_fe_table
void create_active_fe_table()
DoFHandler::tria
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
Definition: dof_handler.h:1543
DoFHandler::set_active_fe_indices
void set_active_fe_indices(const std::vector< unsigned int > &active_fe_indices)
DoFHandler::end_mg
level_cell_iterator end_mg() const
DoFHandler::renumber_dofs
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
DoFHandler::begin_mg
level_cell_iterator begin_mg(const unsigned int level=0) const
DoFHandler::distribute_dofs
void distribute_dofs(const FiniteElement< dim, spacedim > &fe)
DoFHandler::n_boundary_dofs
types::global_dof_index n_boundary_dofs() const
DoFHandler::clear_mg_space
void clear_mg_space()
DoFHandler::connect_to_triangulation_signals
void connect_to_triangulation_signals()
DoFHandler::cell_dof_cache_indices
std::vector< std::vector< types::global_dof_index > > cell_dof_cache_indices
Definition: dof_handler.h:1581
DoFHandler::mg_vertex_dofs
std::vector< MGVertexDoFs > mg_vertex_dofs
Definition: dof_handler.h:1639
DoFHandler::hp_cell_future_fe_indices
std::vector< std::vector< active_fe_index_type > > hp_cell_future_fe_indices
Definition: dof_handler.h:1633
DoFHandler::post_distributed_transfer_action
void post_distributed_transfer_action()
DoFHandler::object_dof_indices
std::vector< std::array< std::vector< types::global_dof_index >, dim+1 > > object_dof_indices
Definition: dof_handler.h:1595
DoFHandler::get_triangulation
const Triangulation< dim, spacedim > & get_triangulation() const
DoFHandler::get_dof_index
types::global_dof_index get_dof_index(const unsigned int obj_level, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index) const
DoFHandler::pre_transfer_action
void pre_transfer_action()
DoFHandler::clear_space
void clear_space()
DoFHandler::hp_object_fe_indices
std::array< std::vector< active_fe_index_type >, dim+1 > hp_object_fe_indices
Definition: dof_handler.h:1614
DoFHandler::reinit
void reinit(const Triangulation< dim, spacedim > &tria)
DoFHandler::object_dof_ptr
std::vector< std::array< std::vector< offset_type >, dim+1 > > object_dof_ptr
Definition: dof_handler.h:1606
DoFHandler::mg_faces
std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces< dim > > mg_faces
Definition: dof_handler.h:1652
DoFHandler::begin_active
active_cell_iterator begin_active(const unsigned int level=0) const
DoFHandler::distribute_mg_dofs
void distribute_mg_dofs()
DoFHandler::end_active
active_cell_iterator end_active(const unsigned int level) const
DoFHandler::set_fe
void set_fe(const FiniteElement< dim, spacedim > &fe)
DoFHandler::hp_capability_enabled
bool hp_capability_enabled
Definition: dof_handler.h:1537
DoFHandler::initialize_local_block_info
void initialize_local_block_info()
DoFHandler::n_dofs
types::global_dof_index n_dofs() const
DoFHandler::update_active_fe_table
void update_active_fe_table()
DoFHandler::cell_dof_cache_ptr
std::vector< std::vector< offset_type > > cell_dof_cache_ptr
Definition: dof_handler.h:1587
DoFHandler::prepare_for_serialization_of_active_fe_indices
void prepare_for_serialization_of_active_fe_indices()
DoFHandler::begin
cell_iterator begin(const unsigned int level=0) const
DoFHandler::hp_object_fe_ptr
std::array< std::vector< offset_type >, dim+1 > hp_object_fe_ptr
Definition: dof_handler.h:1619
DoFHandler::clear
void clear()
DoFHandler::get_active_fe_indices
void get_active_fe_indices(std::vector< unsigned int > &active_fe_indices) const
DoFHandler::max_couplings_between_boundary_dofs
unsigned int max_couplings_between_boundary_dofs() const
DoFHandler::setup_policy
void setup_policy()
DoFHandler::initialize
void initialize(const Triangulation< dim, spacedim > &tria, const FiniteElement< dim, spacedim > &fe)
DoFHandler::deserialize_active_fe_indices
void deserialize_active_fe_indices()
DoFHandler::~DoFHandler
virtual ~DoFHandler() override
FiniteElementData::n_dofs_per_vertex
unsigned int n_dofs_per_vertex() const
FiniteElementData::n_dofs_per_cell
unsigned int n_dofs_per_cell() const
FiniteElementData::n_dofs_per_line
unsigned int n_dofs_per_line() const
IndexSet::is_element
bool is_element(const size_type index) const
Definition: index_set.h:1765
Triangulation::get_mesh_smoothing
virtual const MeshSmoothing & get_mesh_smoothing() const
Triangulation::begin
cell_iterator begin(const unsigned int level=0) const
Triangulation::n_raw_lines
unsigned int n_raw_lines() const
Triangulation::n_raw_faces
unsigned int n_raw_faces() const
Triangulation::n_levels
unsigned int n_levels() const
Triangulation::end
cell_iterator end() const
Triangulation::n_raw_cells
unsigned int n_raw_cells(const unsigned int level) const
Triangulation::max_adjacent_cells
unsigned int max_adjacent_cells() const
Triangulation::vertex_used
bool vertex_used(const unsigned int index) const
Triangulation::end_active
active_cell_iterator end_active(const unsigned int level) const
Triangulation::n_raw_quads
unsigned int n_raw_quads() const
Triangulation::n_cells
unsigned int n_cells() const
Triangulation::signals
Signals signals
Definition: tria.h:2295
Triangulation::n_vertices
unsigned int n_vertices() const
hp::Collection::size
unsigned int size() const
Definition: collection.h:109
hp::FECollection::max_dofs_per_line
unsigned int max_dofs_per_line() const
Definition: fe_collection.h:865
hp::FECollection::max_dofs_per_hex
unsigned int max_dofs_per_hex() const
Definition: fe_collection.h:895
hp::FECollection::max_dofs_per_vertex
unsigned int max_dofs_per_vertex() const
Definition: fe_collection.h:850
hp::FECollection::max_dofs_per_quad
unsigned int max_dofs_per_quad() const
Definition: fe_collection.h:880
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:402
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:403
dof_handler_policy.h
dof_handler.h
fully_distributed_tria.h
geometry_info.h
level
unsigned int level
Definition: grid_out.cc:4590
grid_tools.h
DoFHandler::active_cell_iterators
IteratorRange< active_cell_iterator > active_cell_iterators() const
DoFHandler::mg_cell_iterators
IteratorRange< level_cell_iterator > mg_cell_iterators() const
DoFHandler::cell_iterators_on_level
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
DoFHandler::cell_iterators
IteratorRange< cell_iterator > cell_iterators() const
DoFHandler::active_cell_iterators_on_level
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
Triangulation::cell_iterators_on_level
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
Triangulation::cell_iterators
IteratorRange< cell_iterator > cell_iterators() const
DoFHandler::mg_cell_iterators_on_level
IteratorRange< level_cell_iterator > mg_cell_iterators_on_level(const unsigned int level) const
Threads::TaskGroup::join_all
void join_all() const
Definition: thread_management.h:1524
DoFTools::ExcInvalidBoundaryIndicator
static ::ExceptionBase & ExcInvalidBoundaryIndicator()
StandardExceptions::ExcNotImplemented
static ::ExceptionBase & ExcNotImplemented()
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1465
AssertDimension
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1622
AssertIndexRange
#define AssertIndexRange(index, range)
Definition: exceptions.h:1690
internal::hp::DoFHandlerImplementation::ExcNoDominatedFiniteElementOnChildren
static ::ExceptionBase & ExcNoDominatedFiniteElementOnChildren()
StandardExceptions::ExcInternalError
static ::ExceptionBase & ExcInternalError()
StandardExceptions::ExcDimensionMismatch
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
DoFTools::ExcNoFESelected
static ::ExceptionBase & ExcNoFESelected()
StandardExceptions::ExcMessage
static ::ExceptionBase & ExcMessage(std::string arg1)
AssertThrow
#define AssertThrow(cond, exc)
Definition: exceptions.h:1575
GridTools::exchange_cell_data_to_ghosts
void exchange_cell_data_to_ghosts(const MeshType &mesh, const std::function< std_cxx17::optional< DataType >(const typename MeshType::active_cell_iterator &)> &pack, const std::function< void(const typename MeshType::active_cell_iterator &, const DataType &)> &unpack, const std::function< bool(const typename MeshType::active_cell_iterator &)> &cell_filter=always_return< typename MeshType::active_cell_iterator, bool >{true})
Threads::new_task
Task< RT > new_task(const std::function< RT()> &function)
Definition: thread_management.h:1328
memory_consumption.h
DoFTools::get_active_fe_indices
void get_active_fe_indices(const DoFHandler< dim, spacedim > &dof_handler, std::vector< unsigned int > &active_fe_indices)
Definition: dof_tools.cc:1399
GridRefinement::refine
void refine(Triangulation< dim, spacedim > &tria, const Vector< Number > &criteria, const double threshold, const unsigned int max_to_mark=numbers::invalid_unsigned_int)
Definition: grid_refinement.cc:168
GridRefinement::coarsen
void coarsen(Triangulation< dim, spacedim > &tria, const Vector< Number > &criteria, const double threshold)
Definition: grid_refinement.cc:218
IteratorState::valid
@ valid
Iterator points to a valid object.
Definition: tria_iterator_base.h:38
LAPACKSupport::T
static const char T
Definition: lapack_support.h:165
MemoryConsumption::memory_consumption
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
Definition: memory_consumption.h:266
MemorySpace::swap
void swap(MemorySpaceData< Number, MemorySpace > &, MemorySpaceData< Number, MemorySpace > &)
Definition: memory_space_data.h:92
Physics::Elasticity::Kinematics::l
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
Physics::Elasticity::Kinematics::d
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
TrilinosWrappers::internal::end
VectorType::value_type * end(VectorType &V)
Definition: trilinos_sparse_matrix.cc:67
TrilinosWrappers::internal::begin
VectorType::value_type * begin(VectorType &V)
Definition: trilinos_sparse_matrix.cc:53
TrilinosWrappers::global_index
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
Definition: trilinos_index_access.h:82
Utilities::MPI::sum
T sum(const T &t, const MPI_Comm &mpi_communicator)
Utilities::pack
size_t pack(const T &object, std::vector< char > &dest_buffer, const bool allow_compression=true)
Definition: utilities.h:1218
Utilities::int_to_string
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:473
Utilities::unpack
T unpack(const std::vector< char > &buffer, const bool allow_compression=true)
Definition: utilities.h:1321
hp
Definition: hp.h:118
internal::TriangulationImplementation::n_active_cells
unsigned int n_active_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:12594
internal::hp::DoFHandlerImplementation::dominated_future_fe_on_children
unsigned int dominated_future_fe_on_children(const typename DoFHandler< dim, spacedim >::cell_iterator &parent)
Definition: dof_handler.cc:2093
internal::hp::DoFHandlerImplementation::communicate_future_fe_indices
void communicate_future_fe_indices(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:2080
internal::reinit
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
internal::policy_to_string
std::string policy_to_string(const ::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > &policy)
Definition: dof_handler.cc:51
numbers::internal_face_boundary_id
const types::boundary_id internal_face_boundary_id
Definition: types.h:255
numbers::invalid_unsigned_int
static const unsigned int invalid_unsigned_int
Definition: types.h:196
numbers::invalid_dof_index
const types::global_dof_index invalid_dof_index
Definition: types.h:211
parallel::Triangulation
TriangulationBase< dim, spacedim > Triangulation
Definition: tria_base.h:396
std::min
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
Definition: vectorization.h:5501
types::boundary_id
unsigned int boundary_id
Definition: types.h:129
shared_tria.h
internal::DoFHandlerImplementation::Implementation::get_dof_index
static types::global_dof_index get_dof_index(const DoFHandler< 3, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 3 > > &mg_level, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 3 > > &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant< int, 3 >)
Definition: dof_handler.cc:779
internal::DoFHandlerImplementation::Implementation::reserve_space_mg
static void reserve_space_mg(DoFHandler< 3, spacedim > &dof_handler)
Definition: dof_handler.cc:587
internal::DoFHandlerImplementation::Implementation::set_dof_index
static void set_dof_index(const DoFHandler< 3, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 3 > > &, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 3 > > &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant< int, 2 >)
Definition: dof_handler.cc:897
internal::DoFHandlerImplementation::Implementation::get_dof_index
static types::global_dof_index get_dof_index(const DoFHandler< 3, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 3 > > &, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 3 > > &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant< int, 2 >)
Definition: dof_handler.cc:757
internal::DoFHandlerImplementation::Implementation::reserve_space
static void reserve_space(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:393
internal::DoFHandlerImplementation::Implementation::get_dof_index
static types::global_dof_index get_dof_index(const DoFHandler< 3, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 3 > > &, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 3 > > &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant< int, 1 >)
Definition: dof_handler.cc:735
internal::DoFHandlerImplementation::Implementation::set_dof_index
static void set_dof_index(const DoFHandler< 2, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 2 > > &, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 2 > > &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant< int, 1 >)
Definition: dof_handler.cc:825
internal::DoFHandlerImplementation::Implementation::set_dof_index
static void set_dof_index(const DoFHandler< 2, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 2 > > &mg_level, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 2 > > &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant< int, 2 >)
Definition: dof_handler.cc:849
internal::DoFHandlerImplementation::Implementation::reserve_subentities
static void reserve_subentities(DoFHandler< dim, spacedim > &dof_handler, const unsigned int structdim, const unsigned int n_raw_entities, const T &cell_process)
Definition: dof_handler.cc:334
internal::DoFHandlerImplementation::Implementation::get_dof_index
static types::global_dof_index get_dof_index(const DoFHandler< 2, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 2 > > &, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 2 > > &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant< int, 1 >)
Definition: dof_handler.cc:693
internal::DoFHandlerImplementation::Implementation::set_dof_index
static void set_dof_index(const DoFHandler< 1, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 1 > > &mg_level, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 1 > > &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant< int, 1 >)
Definition: dof_handler.cc:801
internal::DoFHandlerImplementation::Implementation::max_couplings_between_dofs
static unsigned int max_couplings_between_dofs(const DoFHandler< 2, spacedim > &dof_handler)
Definition: dof_handler.cc:106
internal::DoFHandlerImplementation::Implementation::set_dof_index
static void set_dof_index(const DoFHandler< 3, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 3 > > &, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 3 > > &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant< int, 1 >)
Definition: dof_handler.cc:873
internal::DoFHandlerImplementation::Implementation::set_dof_index
static void set_dof_index(const DoFHandler< 3, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 3 > > &mg_level, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 3 > > &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant< int, 3 >)
Definition: dof_handler.cc:921
internal::DoFHandlerImplementation::Implementation::reset_to_empty_objects
static void reset_to_empty_objects(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:261
internal::DoFHandlerImplementation::Implementation::max_couplings_between_dofs
static unsigned int max_couplings_between_dofs(const DoFHandler< 3, spacedim > &dof_handler)
Definition: dof_handler.cc:225
internal::DoFHandlerImplementation::Implementation::reserve_space_mg
static void reserve_space_mg(DoFHandler< 1, spacedim > &dof_handler)
Definition: dof_handler.cc:442
internal::DoFHandlerImplementation::Implementation::get_dof_index
static types::global_dof_index get_dof_index(const DoFHandler< 2, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 2 > > &mg_level, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 2 > > &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant< int, 2 >)
Definition: dof_handler.cc:713
internal::DoFHandlerImplementation::Implementation::max_couplings_between_dofs
static unsigned int max_couplings_between_dofs(const DoFHandler< 1, spacedim > &dof_handler)
Definition: dof_handler.cc:96
internal::DoFHandlerImplementation::Implementation::reserve_space_mg
static void reserve_space_mg(DoFHandler< 2, spacedim > &dof_handler)
Definition: dof_handler.cc:511
internal::DoFHandlerImplementation::Implementation::get_dof_index
static types::global_dof_index get_dof_index(const DoFHandler< 1, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 1 > > &mg_level, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 1 > > &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant< int, 1 >)
Definition: dof_handler.cc:670
internal::DoFHandlerImplementation::Implementation::reserve_cells
static void reserve_cells(DoFHandler< dim, spacedim > &dof_handler, const unsigned int n_inner_dofs_per_cell)
Definition: dof_handler.cc:285
internal::hp::DoFHandlerImplementation::Implementation::collect_fe_indices_on_cells_to_be_refined
static void collect_fe_indices_on_cells_to_be_refined(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:1864
internal::hp::DoFHandlerImplementation::Implementation::determine_fe_from_children
static unsigned int determine_fe_from_children(const typename Triangulation< dim, spacedim >::cell_iterator &, const std::vector< unsigned int > &children_fe_indices, const ::hp::FECollection< dim, spacedim > &fe_collection)
Definition: dof_handler.cc:1993
internal::hp::DoFHandlerImplementation::Implementation::reserve_space_cells
static void reserve_space_cells(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:1117
internal::hp::DoFHandlerImplementation::Implementation::dominated_future_fe_on_children
static unsigned int dominated_future_fe_on_children(const typename DoFHandler< dim, spacedim >::cell_iterator &parent)
Definition: dof_handler.cc:2024
internal::hp::DoFHandlerImplementation::Implementation::distribute_fe_indices_on_refined_cells
static void distribute_fe_indices_on_refined_cells(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:1935
internal::hp::DoFHandlerImplementation::Implementation::ensure_absence_of_future_fe_indices
static void ensure_absence_of_future_fe_indices(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:964
internal::hp::DoFHandlerImplementation::Implementation::reserve_space
static void reserve_space(::DoFHandler< 1, spacedim > &dof_handler)
Definition: dof_handler.cc:1449
internal::hp::DoFHandlerImplementation::Implementation::reserve_space_faces
static void reserve_space_faces(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:1176
internal::hp::DoFHandlerImplementation::Implementation::communicate_future_fe_indices
static void communicate_future_fe_indices(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:1769
internal::hp::DoFHandlerImplementation::Implementation::communicate_active_fe_indices
static void communicate_active_fe_indices(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:1657
internal::hp::DoFHandlerImplementation::Implementation::reserve_space_vertices
static void reserve_space_vertices(DoFHandler< dim, spacedim > &dof_handler)
Definition: dof_handler.cc:984
internal::hp::DoFHandlerImplementation::Implementation::reserve_space
static void reserve_space(::DoFHandler< 2, spacedim > &dof_handler)
Definition: dof_handler.cc:1473
internal::hp::DoFHandlerImplementation::Implementation::reserve_space
static void reserve_space(::DoFHandler< 3, spacedim > &dof_handler)
Definition: dof_handler.cc:1499
tria_accessor.h
tria_iterator.h
tria_levels.h
Generated by doxygen 1.9.4