Reference documentation for deal.II version Git e159ac89de 2020-06-06 19:38:41 +0200
\(\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 - 2020 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 
18 
21 
23 #include <deal.II/dofs/dof_faces.h>
27 
28 #include <deal.II/fe/fe.h>
29 
30 #include <deal.II/grid/tria.h>
34 
35 #include <algorithm>
36 #include <memory>
37 #include <set>
38 #include <unordered_set>
39 
41 
42 
43 template <int dim, int spacedim>
44 const unsigned int DoFHandler<dim, spacedim>::dimension;
45 
46 template <int dim, int spacedim>
48 
49 template <int dim, int spacedim>
51 
52 
53 // reference the invalid_dof_index variable explicitly to work around
54 // a bug in the icc8 compiler
55 namespace internal
56 {
57  template <int dim, int spacedim>
60  {
62  }
63 } // namespace internal
64 
65 
66 
67 namespace internal
68 {
69  template <int dim, int spacedim>
70  std::string
71  policy_to_string(const ::internal::DoFHandlerImplementation::Policy::
72  PolicyBase<dim, spacedim> &policy)
73  {
74  std::string policy_name;
75  if (dynamic_cast<const typename ::internal::DoFHandlerImplementation::
76  Policy::Sequential<::DoFHandler<dim, spacedim>> *>(
77  &policy) ||
78  dynamic_cast<
79  const typename ::internal::DoFHandlerImplementation::Policy::
80  Sequential<::hp::DoFHandler<dim, spacedim>> *>(&policy))
81  policy_name = "Policy::Sequential<";
82  else if (dynamic_cast<
83  const typename ::internal::DoFHandlerImplementation::
84  Policy::ParallelDistributed<::DoFHandler<dim, spacedim>>
85  *>(&policy) ||
86  dynamic_cast<
87  const typename ::internal::DoFHandlerImplementation::
88  Policy::ParallelDistributed<
90  policy_name = "Policy::ParallelDistributed<";
91  else if (dynamic_cast<
92  const typename ::internal::DoFHandlerImplementation::
93  Policy::ParallelShared<::DoFHandler<dim, spacedim>> *>(
94  &policy) ||
95  dynamic_cast<const typename ::internal::
98  &policy))
99  policy_name = "Policy::ParallelShared<";
100  else
101  AssertThrow(false, ExcNotImplemented());
102  policy_name += Utilities::int_to_string(dim) + "," +
103  Utilities::int_to_string(spacedim) + ">";
104  return policy_name;
105  }
106 
107 
108  namespace DoFHandlerImplementation
109  {
110  // access class
111  // ::DoFHandler instead of
112  // namespace internal::DoFHandler
113  using ::DoFHandler;
114 
115 
121  {
126  template <int spacedim>
127  static unsigned int
129  {
130  return std::min(static_cast<types::global_dof_index>(
131  3 * dof_handler.get_fe().dofs_per_vertex +
132  2 * dof_handler.get_fe().dofs_per_line),
133  dof_handler.n_dofs());
134  }
135 
136 
137 
138  template <int spacedim>
139  static unsigned int
141  {
142  // get these numbers by drawing pictures
143  // and counting...
144  // example:
145  // | | |
146  // --x-----x--x--X--
147  // | | | |
148  // | x--x--x
149  // | | | |
150  // --x--x--*--x--x--
151  // | | | |
152  // x--x--x |
153  // | | | |
154  // --X--x--x-----x--
155  // | | |
156  // x = vertices connected with center vertex *;
157  // = total of 19
158  // (the X vertices are connected with * if
159  // the vertices adjacent to X are hanging
160  // nodes)
161  // count lines -> 28 (don't forget to count
162  // mother and children separately!)
163  types::global_dof_index max_couplings;
164  switch (dof_handler.tria->max_adjacent_cells())
165  {
166  case 4:
167  max_couplings = 19 * dof_handler.get_fe().dofs_per_vertex +
168  28 * dof_handler.get_fe().dofs_per_line +
169  8 * dof_handler.get_fe().dofs_per_quad;
170  break;
171  case 5:
172  max_couplings = 21 * dof_handler.get_fe().dofs_per_vertex +
173  31 * dof_handler.get_fe().dofs_per_line +
174  9 * dof_handler.get_fe().dofs_per_quad;
175  break;
176  case 6:
177  max_couplings = 28 * dof_handler.get_fe().dofs_per_vertex +
178  42 * dof_handler.get_fe().dofs_per_line +
179  12 * dof_handler.get_fe().dofs_per_quad;
180  break;
181  case 7:
182  max_couplings = 30 * dof_handler.get_fe().dofs_per_vertex +
183  45 * dof_handler.get_fe().dofs_per_line +
184  13 * dof_handler.get_fe().dofs_per_quad;
185  break;
186  case 8:
187  max_couplings = 37 * dof_handler.get_fe().dofs_per_vertex +
188  56 * dof_handler.get_fe().dofs_per_line +
189  16 * dof_handler.get_fe().dofs_per_quad;
190  break;
191 
192  // the following numbers are not based on actual counting but by
193  // extrapolating the number sequences from the previous ones (for
194  // example, for dofs_per_vertex, the sequence above is 19, 21, 28,
195  // 30, 37, and is continued as follows):
196  case 9:
197  max_couplings = 39 * dof_handler.get_fe().dofs_per_vertex +
198  59 * dof_handler.get_fe().dofs_per_line +
199  17 * dof_handler.get_fe().dofs_per_quad;
200  break;
201  case 10:
202  max_couplings = 46 * dof_handler.get_fe().dofs_per_vertex +
203  70 * dof_handler.get_fe().dofs_per_line +
204  20 * dof_handler.get_fe().dofs_per_quad;
205  break;
206  case 11:
207  max_couplings = 48 * dof_handler.get_fe().dofs_per_vertex +
208  73 * dof_handler.get_fe().dofs_per_line +
209  21 * dof_handler.get_fe().dofs_per_quad;
210  break;
211  case 12:
212  max_couplings = 55 * dof_handler.get_fe().dofs_per_vertex +
213  84 * dof_handler.get_fe().dofs_per_line +
214  24 * dof_handler.get_fe().dofs_per_quad;
215  break;
216  case 13:
217  max_couplings = 57 * dof_handler.get_fe().dofs_per_vertex +
218  87 * dof_handler.get_fe().dofs_per_line +
219  25 * dof_handler.get_fe().dofs_per_quad;
220  break;
221  case 14:
222  max_couplings = 63 * dof_handler.get_fe().dofs_per_vertex +
223  98 * dof_handler.get_fe().dofs_per_line +
224  28 * dof_handler.get_fe().dofs_per_quad;
225  break;
226  case 15:
227  max_couplings = 65 * dof_handler.get_fe().dofs_per_vertex +
228  103 * dof_handler.get_fe().dofs_per_line +
229  29 * dof_handler.get_fe().dofs_per_quad;
230  break;
231  case 16:
232  max_couplings = 72 * dof_handler.get_fe().dofs_per_vertex +
233  114 * dof_handler.get_fe().dofs_per_line +
234  32 * dof_handler.get_fe().dofs_per_quad;
235  break;
236 
237  default:
238  Assert(false, ExcNotImplemented());
239  max_couplings = 0;
240  }
241  return std::min(max_couplings, dof_handler.n_dofs());
242  }
243 
244 
245  template <int spacedim>
246  static unsigned int
248  {
249  // TODO:[?] Invent significantly better estimates than the ones in this
250  // function
251 
252  // doing the same thing here is a
253  // rather complicated thing, compared
254  // to the 2d case, since it is hard
255  // to draw pictures with several
256  // refined hexahedra :-) so I
257  // presently only give a coarse
258  // estimate for the case that at most
259  // 8 hexes meet at each vertex
260  //
261  // can anyone give better estimate
262  // here?
263  const unsigned int max_adjacent_cells =
264  dof_handler.tria->max_adjacent_cells();
265 
266  types::global_dof_index max_couplings;
267  if (max_adjacent_cells <= 8)
268  max_couplings = 7 * 7 * 7 * dof_handler.get_fe().dofs_per_vertex +
269  7 * 6 * 7 * 3 * dof_handler.get_fe().dofs_per_line +
270  9 * 4 * 7 * 3 * dof_handler.get_fe().dofs_per_quad +
271  27 * dof_handler.get_fe().dofs_per_hex;
272  else
273  {
274  Assert(false, ExcNotImplemented());
275  max_couplings = 0;
276  }
277 
278  return std::min(max_couplings, dof_handler.n_dofs());
279  }
280 
281 
291  template <int spacedim>
292  static void reserve_space(DoFHandler<1, spacedim> &dof_handler)
293  {
294  dof_handler.vertex_dofs.resize(dof_handler.tria->n_vertices() *
295  dof_handler.get_fe().dofs_per_vertex,
297 
298  for (unsigned int i = 0; i < dof_handler.tria->n_levels(); ++i)
299  {
300  dof_handler.levels.emplace_back(
302 
303  dof_handler.levels.back()->dof_object.dofs.resize(
304  dof_handler.tria->n_raw_cells(i) *
305  dof_handler.get_fe().dofs_per_line,
307 
308  dof_handler.levels.back()->cell_dof_indices_cache.resize(
309  dof_handler.tria->n_raw_cells(i) *
310  dof_handler.get_fe().dofs_per_cell,
312  }
313  }
314 
315 
316  template <int spacedim>
317  static void reserve_space(DoFHandler<2, spacedim> &dof_handler)
318  {
319  dof_handler.vertex_dofs.resize(dof_handler.tria->n_vertices() *
320  dof_handler.get_fe().dofs_per_vertex,
322 
323  for (unsigned int i = 0; i < dof_handler.tria->n_levels(); ++i)
324  {
325  dof_handler.levels.emplace_back(
327 
328  dof_handler.levels.back()->dof_object.dofs.resize(
329  dof_handler.tria->n_raw_cells(i) *
330  dof_handler.get_fe().dofs_per_quad,
332 
333  dof_handler.levels.back()->cell_dof_indices_cache.resize(
334  dof_handler.tria->n_raw_cells(i) *
335  dof_handler.get_fe().dofs_per_cell,
337  }
338 
339  dof_handler.faces =
340  std::make_unique<internal::DoFHandlerImplementation::DoFFaces<2>>();
341  // avoid access to n_raw_lines when there are no cells
342  if (dof_handler.tria->n_cells() > 0)
343  {
344  dof_handler.faces->lines.dofs.resize(
345  dof_handler.tria->n_raw_lines() *
346  dof_handler.get_fe().dofs_per_line,
348  }
349  }
350 
351 
352  template <int spacedim>
353  static void reserve_space(DoFHandler<3, spacedim> &dof_handler)
354  {
355  dof_handler.vertex_dofs.resize(dof_handler.tria->n_vertices() *
356  dof_handler.get_fe().dofs_per_vertex,
358 
359  for (unsigned int i = 0; i < dof_handler.tria->n_levels(); ++i)
360  {
361  dof_handler.levels.emplace_back(
363 
364  dof_handler.levels.back()->dof_object.dofs.resize(
365  dof_handler.tria->n_raw_cells(i) *
366  dof_handler.get_fe().dofs_per_hex,
368 
369  dof_handler.levels.back()->cell_dof_indices_cache.resize(
370  dof_handler.tria->n_raw_cells(i) *
371  dof_handler.get_fe().dofs_per_cell,
373  }
374  dof_handler.faces =
375  std::make_unique<internal::DoFHandlerImplementation::DoFFaces<3>>();
376 
377  // avoid access to n_raw_lines when there are no cells
378  if (dof_handler.tria->n_cells() > 0)
379  {
380  dof_handler.faces->lines.dofs.resize(
381  dof_handler.tria->n_raw_lines() *
382  dof_handler.get_fe().dofs_per_line,
384  dof_handler.faces->quads.dofs.resize(
385  dof_handler.tria->n_raw_quads() *
386  dof_handler.get_fe().dofs_per_quad,
388  }
389  }
390 
391  template <int spacedim>
392  static void reserve_space_mg(DoFHandler<1, spacedim> &dof_handler)
393  {
394  Assert(dof_handler.get_triangulation().n_levels() > 0,
395  ExcMessage("Invalid triangulation"));
396  dof_handler.clear_mg_space();
397 
398  const ::Triangulation<1, spacedim> &tria =
399  dof_handler.get_triangulation();
400  const unsigned int dofs_per_line = dof_handler.get_fe().dofs_per_line;
401  const unsigned int n_levels = tria.n_levels();
402 
403  for (unsigned int i = 0; i < n_levels; ++i)
404  {
405  dof_handler.mg_levels.emplace_back(
407  dof_handler.mg_levels.back()->dof_object.dofs =
408  std::vector<types::global_dof_index>(tria.n_raw_lines(i) *
409  dofs_per_line,
411  }
412 
413  const unsigned int n_vertices = tria.n_vertices();
414 
415  dof_handler.mg_vertex_dofs.resize(n_vertices);
416 
417  std::vector<unsigned int> max_level(n_vertices, 0);
418  std::vector<unsigned int> min_level(n_vertices, n_levels);
419 
420  for (typename ::Triangulation<1, spacedim>::cell_iterator cell =
421  tria.begin();
422  cell != tria.end();
423  ++cell)
424  {
425  const unsigned int level = cell->level();
426 
427  for (const unsigned int vertex : GeometryInfo<1>::vertex_indices())
428  {
429  const unsigned int vertex_index = cell->vertex_index(vertex);
430 
431  if (min_level[vertex_index] > level)
432  min_level[vertex_index] = level;
433 
434  if (max_level[vertex_index] < level)
435  max_level[vertex_index] = level;
436  }
437  }
438 
439  for (unsigned int vertex = 0; vertex < n_vertices; ++vertex)
440  if (tria.vertex_used(vertex))
441  {
442  Assert(min_level[vertex] < n_levels, ExcInternalError());
443  Assert(max_level[vertex] >= min_level[vertex],
444  ExcInternalError());
445  dof_handler.mg_vertex_dofs[vertex].init(
446  min_level[vertex],
447  max_level[vertex],
448  dof_handler.get_fe().dofs_per_vertex);
449  }
450 
451  else
452  {
453  Assert(min_level[vertex] == n_levels, ExcInternalError());
454  Assert(max_level[vertex] == 0, ExcInternalError());
455  dof_handler.mg_vertex_dofs[vertex].init(1, 0, 0);
456  }
457  }
458 
459  template <int spacedim>
460  static void reserve_space_mg(DoFHandler<2, spacedim> &dof_handler)
461  {
462  Assert(dof_handler.get_triangulation().n_levels() > 0,
463  ExcMessage("Invalid triangulation"));
464  dof_handler.clear_mg_space();
465 
466  const ::FiniteElement<2, spacedim> &fe = dof_handler.get_fe();
467  const ::Triangulation<2, spacedim> &tria =
468  dof_handler.get_triangulation();
469  const unsigned int n_levels = tria.n_levels();
470 
471  for (unsigned int i = 0; i < n_levels; ++i)
472  {
473  dof_handler.mg_levels.emplace_back(
474  std::make_unique<
476  dof_handler.mg_levels.back()->dof_object.dofs =
477  std::vector<types::global_dof_index>(tria.n_raw_quads(i) *
478  fe.dofs_per_quad,
480  }
481 
482  dof_handler.mg_faces =
483  std::make_unique<internal::DoFHandlerImplementation::DoFFaces<2>>();
484  dof_handler.mg_faces->lines.dofs = std::vector<types::global_dof_index>(
485  tria.n_raw_lines() * fe.dofs_per_line, numbers::invalid_dof_index);
486 
487  const unsigned int n_vertices = tria.n_vertices();
488 
489  dof_handler.mg_vertex_dofs.resize(n_vertices);
490 
491  std::vector<unsigned int> max_level(n_vertices, 0);
492  std::vector<unsigned int> min_level(n_vertices, n_levels);
493 
494  for (typename ::Triangulation<2, spacedim>::cell_iterator cell =
495  tria.begin();
496  cell != tria.end();
497  ++cell)
498  {
499  const unsigned int level = cell->level();
500 
501  for (const unsigned int vertex : GeometryInfo<2>::vertex_indices())
502  {
503  const unsigned int vertex_index = cell->vertex_index(vertex);
504 
505  if (min_level[vertex_index] > level)
506  min_level[vertex_index] = level;
507 
508  if (max_level[vertex_index] < level)
509  max_level[vertex_index] = level;
510  }
511  }
512 
513  for (unsigned int vertex = 0; vertex < n_vertices; ++vertex)
514  if (tria.vertex_used(vertex))
515  {
516  Assert(min_level[vertex] < n_levels, ExcInternalError());
517  Assert(max_level[vertex] >= min_level[vertex],
518  ExcInternalError());
519  dof_handler.mg_vertex_dofs[vertex].init(min_level[vertex],
520  max_level[vertex],
521  fe.dofs_per_vertex);
522  }
523 
524  else
525  {
526  Assert(min_level[vertex] == n_levels, ExcInternalError());
527  Assert(max_level[vertex] == 0, ExcInternalError());
528  dof_handler.mg_vertex_dofs[vertex].init(1, 0, 0);
529  }
530  }
531 
532  template <int spacedim>
533  static void reserve_space_mg(DoFHandler<3, spacedim> &dof_handler)
534  {
535  Assert(dof_handler.get_triangulation().n_levels() > 0,
536  ExcMessage("Invalid triangulation"));
537  dof_handler.clear_mg_space();
538 
539  const ::FiniteElement<3, spacedim> &fe = dof_handler.get_fe();
540  const ::Triangulation<3, spacedim> &tria =
541  dof_handler.get_triangulation();
542  const unsigned int n_levels = tria.n_levels();
543 
544  for (unsigned int i = 0; i < n_levels; ++i)
545  {
546  dof_handler.mg_levels.emplace_back(
547  std::make_unique<
549  dof_handler.mg_levels.back()->dof_object.dofs =
550  std::vector<types::global_dof_index>(tria.n_raw_hexs(i) *
551  fe.dofs_per_hex,
553  }
554 
555  dof_handler.mg_faces =
556  std::make_unique<internal::DoFHandlerImplementation::DoFFaces<3>>();
557  dof_handler.mg_faces->lines.dofs = std::vector<types::global_dof_index>(
558  tria.n_raw_lines() * fe.dofs_per_line, numbers::invalid_dof_index);
559  dof_handler.mg_faces->quads.dofs = std::vector<types::global_dof_index>(
560  tria.n_raw_quads() * fe.dofs_per_quad, numbers::invalid_dof_index);
561 
562  const unsigned int n_vertices = tria.n_vertices();
563 
564  dof_handler.mg_vertex_dofs.resize(n_vertices);
565 
566  std::vector<unsigned int> max_level(n_vertices, 0);
567  std::vector<unsigned int> min_level(n_vertices, n_levels);
568 
569  for (typename ::Triangulation<3, spacedim>::cell_iterator cell =
570  tria.begin();
571  cell != tria.end();
572  ++cell)
573  {
574  const unsigned int level = cell->level();
575 
576  for (const unsigned int vertex : GeometryInfo<3>::vertex_indices())
577  {
578  const unsigned int vertex_index = cell->vertex_index(vertex);
579 
580  if (min_level[vertex_index] > level)
581  min_level[vertex_index] = level;
582 
583  if (max_level[vertex_index] < level)
584  max_level[vertex_index] = level;
585  }
586  }
587 
588  for (unsigned int vertex = 0; vertex < n_vertices; ++vertex)
589  if (tria.vertex_used(vertex))
590  {
591  Assert(min_level[vertex] < n_levels, ExcInternalError());
592  Assert(max_level[vertex] >= min_level[vertex],
593  ExcInternalError());
594  dof_handler.mg_vertex_dofs[vertex].init(min_level[vertex],
595  max_level[vertex],
596  fe.dofs_per_vertex);
597  }
598 
599  else
600  {
601  Assert(min_level[vertex] == n_levels, ExcInternalError());
602  Assert(max_level[vertex] == 0, ExcInternalError());
603  dof_handler.mg_vertex_dofs[vertex].init(1, 0, 0);
604  }
605  }
606 
607 
608 
609  template <int spacedim>
612  const DoFHandler<1, spacedim> &dof_handler,
614  &mg_level,
616  &,
617  const unsigned int obj_index,
618  const unsigned int fe_index,
619  const unsigned int local_index,
620  const std::integral_constant<int, 1>)
621  {
622  return mg_level->dof_object.get_dof_index(dof_handler,
623  obj_index,
624  fe_index,
625  local_index);
626  }
627 
628  template <int spacedim>
631  const DoFHandler<2, spacedim> &dof_handler,
633  &,
635  & mg_faces,
636  const unsigned int obj_index,
637  const unsigned int fe_index,
638  const unsigned int local_index,
639  const std::integral_constant<int, 1>)
640  {
641  return mg_faces->lines.get_dof_index(dof_handler,
642  obj_index,
643  fe_index,
644  local_index);
645  }
646 
647  template <int spacedim>
650  const DoFHandler<2, spacedim> &dof_handler,
652  &mg_level,
654  &,
655  const unsigned int obj_index,
656  const unsigned int fe_index,
657  const unsigned int local_index,
658  const std::integral_constant<int, 2>)
659  {
660  return mg_level->dof_object.get_dof_index(dof_handler,
661  obj_index,
662  fe_index,
663  local_index);
664  }
665 
666  template <int spacedim>
669  const DoFHandler<3, spacedim> &dof_handler,
671  &,
673  & mg_faces,
674  const unsigned int obj_index,
675  const unsigned int fe_index,
676  const unsigned int local_index,
677  const std::integral_constant<int, 1>)
678  {
679  return mg_faces->lines.get_dof_index(dof_handler,
680  obj_index,
681  fe_index,
682  local_index);
683  }
684 
685  template <int spacedim>
688  const DoFHandler<3, spacedim> &dof_handler,
690  &,
692  & mg_faces,
693  const unsigned int obj_index,
694  const unsigned int fe_index,
695  const unsigned int local_index,
696  const std::integral_constant<int, 2>)
697  {
698  return mg_faces->quads.get_dof_index(dof_handler,
699  obj_index,
700  fe_index,
701  local_index);
702  }
703 
704  template <int spacedim>
707  const DoFHandler<3, spacedim> &dof_handler,
709  &mg_level,
711  &,
712  const unsigned int obj_index,
713  const unsigned int fe_index,
714  const unsigned int local_index,
715  const std::integral_constant<int, 3>)
716  {
717  return mg_level->dof_object.get_dof_index(dof_handler,
718  obj_index,
719  fe_index,
720  local_index);
721  }
722 
723  template <int spacedim>
724  static void
726  const DoFHandler<1, spacedim> &dof_handler,
728  &mg_level,
730  &,
731  const unsigned int obj_index,
732  const unsigned int fe_index,
733  const unsigned int local_index,
735  const std::integral_constant<int, 1>)
736  {
737  mg_level->dof_object.set_dof_index(
738  dof_handler, obj_index, fe_index, local_index, global_index);
739  }
740 
741  template <int spacedim>
742  static void
744  const DoFHandler<2, spacedim> &dof_handler,
746  &,
748  & mg_faces,
749  const unsigned int obj_index,
750  const unsigned int fe_index,
751  const unsigned int local_index,
753  const std::integral_constant<int, 1>)
754  {
755  mg_faces->lines.set_dof_index(
756  dof_handler, obj_index, fe_index, local_index, global_index);
757  }
758 
759  template <int spacedim>
760  static void
762  const DoFHandler<2, spacedim> &dof_handler,
764  &mg_level,
766  &,
767  const unsigned int obj_index,
768  const unsigned int fe_index,
769  const unsigned int local_index,
771  const std::integral_constant<int, 2>)
772  {
773  mg_level->dof_object.set_dof_index(
774  dof_handler, obj_index, fe_index, local_index, global_index);
775  }
776 
777  template <int spacedim>
778  static void
780  const DoFHandler<3, spacedim> &dof_handler,
782  &,
784  & mg_faces,
785  const unsigned int obj_index,
786  const unsigned int fe_index,
787  const unsigned int local_index,
789  const std::integral_constant<int, 1>)
790  {
791  mg_faces->lines.set_dof_index(
792  dof_handler, obj_index, fe_index, local_index, global_index);
793  }
794 
795  template <int spacedim>
796  static void
798  const DoFHandler<3, spacedim> &dof_handler,
800  &,
802  & mg_faces,
803  const unsigned int obj_index,
804  const unsigned int fe_index,
805  const unsigned int local_index,
807  const std::integral_constant<int, 2>)
808  {
809  mg_faces->quads.set_dof_index(
810  dof_handler, obj_index, fe_index, local_index, global_index);
811  }
812 
813  template <int spacedim>
814  static void
816  const DoFHandler<3, spacedim> &dof_handler,
818  &mg_level,
820  &,
821  const unsigned int obj_index,
822  const unsigned int fe_index,
823  const unsigned int local_index,
825  const std::integral_constant<int, 3>)
826  {
827  mg_level->dof_object.set_dof_index(
828  dof_handler, obj_index, fe_index, local_index, global_index);
829  }
830  };
831  } // namespace DoFHandlerImplementation
832 } // namespace internal
833 
834 
835 
836 template <int dim, int spacedim>
838  : tria(&tria, typeid(*this).name())
839  , faces(nullptr)
840  , mg_faces(nullptr)
841 {
842  // decide whether we need a sequential or a parallel distributed policy
843  if (dynamic_cast<const parallel::shared::Triangulation<dim, spacedim> *>(
844  &tria) != nullptr)
845  policy =
846  std::make_unique<internal::DoFHandlerImplementation::Policy::
847  ParallelShared<DoFHandler<dim, spacedim>>>(*this);
848  else if (dynamic_cast<
850  &tria) == nullptr)
851  policy =
853  DoFHandler<dim, spacedim>>>(*this);
854  else
855  policy =
856  std::make_unique<internal::DoFHandlerImplementation::Policy::
857  ParallelDistributed<DoFHandler<dim, spacedim>>>(*this);
858 }
859 
860 
861 template <int dim, int spacedim>
863  : tria(nullptr, typeid(*this).name())
864 {}
865 
866 
867 template <int dim, int spacedim>
869 {
870  // release allocated memory
871  // virtual functions called in constructors and destructors never use the
872  // override in a derived class
873  // for clarity be explicit on which function is called
875 
876  // also release the policy. this needs to happen before the
877  // current object disappears because the policy objects
878  // store references to the DoFhandler object they work on
879  policy.reset();
880 }
881 
882 
883 template <int dim, int spacedim>
884 void
887 {
888  tria = &t;
889  faces = nullptr;
891 
892  // decide whether we need a sequential or a parallel distributed policy
893  if (dynamic_cast<const parallel::shared::Triangulation<dim, spacedim> *>(
894  &t) != nullptr)
895  policy =
896  std::make_unique<internal::DoFHandlerImplementation::Policy::
897  ParallelShared<DoFHandler<dim, spacedim>>>(*this);
898  else if (dynamic_cast<
900  &t) != nullptr)
901  policy =
902  std::make_unique<internal::DoFHandlerImplementation::Policy::
903  ParallelDistributed<DoFHandler<dim, spacedim>>>(*this);
904  else
905  policy =
907  DoFHandler<dim, spacedim>>>(*this);
908 
909  distribute_dofs(fe);
910 }
911 
912 
913 
914 /*------------------------ Cell iterator functions ------------------------*/
915 
916 template <int dim, int spacedim>
918 DoFHandler<dim, spacedim>::begin(const unsigned int level) const
919 {
921  this->get_triangulation().begin(level);
922  if (cell == this->get_triangulation().end(level))
923  return end(level);
924  return cell_iterator(*cell, this);
925 }
926 
927 
928 
929 template <int dim, int spacedim>
932 {
933  // level is checked in begin
934  cell_iterator i = begin(level);
935  if (i.state() != IteratorState::valid)
936  return i;
937  while (i->has_children())
938  if ((++i).state() != IteratorState::valid)
939  return i;
940  return i;
941 }
942 
943 
944 
945 template <int dim, int spacedim>
948 {
949  return cell_iterator(&this->get_triangulation(), -1, -1, this);
950 }
951 
952 
953 template <int dim, int spacedim>
955 DoFHandler<dim, spacedim>::end(const unsigned int level) const
956 {
958  this->get_triangulation().end(level);
959  if (cell.state() != IteratorState::valid)
960  return end();
961  return cell_iterator(*cell, this);
962 }
963 
964 
965 template <int dim, int spacedim>
968 {
970  this->get_triangulation().end_active(level);
971  if (cell.state() != IteratorState::valid)
972  return active_cell_iterator(end());
973  return active_cell_iterator(*cell, this);
974 }
975 
976 
977 
978 template <int dim, int spacedim>
981 {
982  // Assert(this->has_level_dofs(), ExcMessage("You can only iterate over mg "
983  // "levels if mg dofs got distributed."));
985  this->get_triangulation().begin(level);
986  if (cell == this->get_triangulation().end(level))
987  return end_mg(level);
988  return level_cell_iterator(*cell, this);
989 }
990 
991 
992 template <int dim, int spacedim>
994 DoFHandler<dim, spacedim>::end_mg(const unsigned int level) const
995 {
996  // Assert(this->has_level_dofs(), ExcMessage("You can only iterate over mg "
997  // "levels if mg dofs got distributed."));
999  this->get_triangulation().end(level);
1000  if (cell.state() != IteratorState::valid)
1001  return end();
1002  return level_cell_iterator(*cell, this);
1003 }
1004 
1005 
1006 template <int dim, int spacedim>
1009 {
1010  return level_cell_iterator(&this->get_triangulation(), -1, -1, this);
1011 }
1012 
1013 
1014 
1015 template <int dim, int spacedim>
1018 {
1020  begin(), end());
1021 }
1022 
1023 
1024 template <int dim, int spacedim>
1027 {
1028  return IteratorRange<
1030  end());
1031 }
1032 
1033 
1034 
1035 template <int dim, int spacedim>
1038 {
1040  begin_mg(), end_mg());
1041 }
1042 
1043 
1044 
1045 template <int dim, int spacedim>
1048  const unsigned int level) const
1049 {
1051  begin(level), end(level));
1052 }
1053 
1054 
1055 
1056 template <int dim, int spacedim>
1059  const unsigned int level) const
1060 {
1061  return IteratorRange<
1063  begin_active(level), end_active(level));
1064 }
1065 
1066 
1067 
1068 template <int dim, int spacedim>
1071  const unsigned int level) const
1072 {
1074  begin_mg(level), end_mg(level));
1075 }
1076 
1077 
1078 
1079 //---------------------------------------------------------------------------
1080 
1081 
1082 
1083 template <int dim, int spacedim>
1086 {
1087  const FiniteElement<dim, spacedim> &fe = get_fe();
1088  const unsigned int dofs_per_face = fe.n_dofs_per_face();
1089  std::vector<::types::global_dof_index> dofs_on_face(dofs_per_face);
1090 
1091  const IndexSet &owned_dofs = locally_owned_dofs();
1092 
1093  std::unordered_set<unsigned int> boundary_dof_indices;
1094  // Loop over all faces of all cells and see whether they are at a boundary.
1095  // note (i) that we visit interior faces twice (which we don't care about) but
1096  // exterior faces only once as is appropriate, and (ii) that we need not take
1097  // special care of single lines (using @p{cell->has_boundary_lines}), since we
1098  // do not support boundaries of dimension dim-2, and so every boundary line is
1099  // also part of a boundary face.
1100  for (const auto &cell : this->active_cell_iterators())
1101  {
1102  if (cell->is_locally_owned() && cell->at_boundary())
1103  {
1104  for (const unsigned int iface :
1106  {
1107  const auto face = cell->face(iface);
1108  if (face->at_boundary())
1109  {
1110  cell->face(iface)->get_dof_indices(dofs_on_face);
1111  for (unsigned int idof_face = 0; idof_face < dofs_per_face;
1112  ++idof_face)
1113  {
1114  const unsigned int global_idof_index =
1115  dofs_on_face[idof_face];
1116 
1117  // If dof is not already in our hashmap, then insert it
1118  // into our set of surface indices. However, do not add if
1119  // it is not locally owned, it will get picked up by
1120  // another processor
1121  if (owned_dofs.is_element(global_idof_index))
1122  {
1123  boundary_dof_indices.insert(global_idof_index);
1124  }
1125  }
1126  }
1127  }
1128  }
1129  }
1130  return boundary_dof_indices.size();
1131 }
1132 
1133 
1134 
1135 template <int dim, int spacedim>
1138  const std::set<types::boundary_id> &boundary_ids) const
1139 {
1140  Assert(boundary_ids.find(numbers::internal_face_boundary_id) ==
1141  boundary_ids.end(),
1143 
1144  const FiniteElement<dim, spacedim> &fe = get_fe();
1145  const unsigned int dofs_per_face = fe.n_dofs_per_face();
1146  std::vector<::types::global_dof_index> dofs_on_face(dofs_per_face);
1147 
1148  const IndexSet &owned_dofs = locally_owned_dofs();
1149 
1150  std::unordered_set<unsigned int> boundary_dof_indices;
1151  // Loop over all faces of all cells and see whether they are at a boundary.
1152  // note (i) that we visit interior faces twice (which we don't care about) but
1153  // exterior faces only once as is appropriate, and (ii) that we need not take
1154  // special care of single lines (using @p{cell->has_boundary_lines}), since we
1155  // do not support boundaries of dimension dim-2, and so every boundary line is
1156  // also part of a boundary face.
1157  for (const auto &cell : this->active_cell_iterators())
1158  {
1159  if (cell->is_locally_owned() && cell->at_boundary())
1160  {
1161  for (const unsigned int iface :
1163  {
1164  const auto face = cell->face(iface);
1165  const unsigned int boundary_id = face->boundary_id();
1166 
1167  if (face->at_boundary() &&
1168  boundary_ids.find(boundary_id) != boundary_ids.end())
1169  {
1170  cell->face(iface)->get_dof_indices(dofs_on_face);
1171  for (unsigned int idof_face = 0; idof_face < dofs_per_face;
1172  ++idof_face)
1173  {
1174  const unsigned int global_idof_index =
1175  dofs_on_face[idof_face];
1176 
1177  // If dof is not already in our hashmap, then insert it
1178  // into our set of surface indices. However, do not add if
1179  // it is not locally owned, it will get picked up by
1180  // another processor
1181  if (owned_dofs.is_element(global_idof_index))
1182  {
1183  boundary_dof_indices.insert(global_idof_index);
1184  }
1185  }
1186  }
1187  }
1188  }
1189  }
1190  return boundary_dof_indices.size();
1191 }
1192 
1193 
1194 
1195 template <int dim, int spacedim>
1196 std::size_t
1198 {
1199  std::size_t mem =
1208  for (unsigned int i = 0; i < levels.size(); ++i)
1210 
1211  for (unsigned int level = 0; level < mg_levels.size(); ++level)
1212  mem += mg_levels[level]->memory_consumption();
1213 
1214  if (mg_faces != nullptr)
1216 
1217  for (unsigned int i = 0; i < mg_vertex_dofs.size(); ++i)
1218  mem += sizeof(MGVertexDoFs) + (1 + mg_vertex_dofs[i].get_finest_level() -
1219  mg_vertex_dofs[i].get_coarsest_level()) *
1220  sizeof(types::global_dof_index);
1221 
1222  return mem;
1223 }
1224 
1225 
1226 
1227 template <int dim, int spacedim>
1228 void
1230 {
1231  // Only recreate the FECollection if we don't already store
1232  // the exact same FiniteElement object.
1233  if (fe_collection.size() == 0 || fe_collection[0] != ff)
1235 }
1236 
1237 
1238 
1239 template <int dim, int spacedim>
1240 void
1242  const FiniteElement<dim, spacedim> &ff)
1243 {
1244  Assert(
1245  tria != nullptr,
1246  ExcMessage(
1247  "You need to set the Triangulation in the DoFHandler using initialize() or "
1248  "in the constructor before you can distribute DoFs."));
1249  Assert(tria->n_levels() > 0,
1250  ExcMessage("The Triangulation you are using is empty!"));
1251 
1252  // first, assign the finite_element
1253  set_fe(ff);
1254 
1255  // delete all levels and set them
1256  // up newly. note that we still
1257  // have to allocate space for all
1258  // degrees of freedom on this mesh
1259  // (including ghost and cells that
1260  // are entirely stored on different
1261  // processors), though we may not
1262  // assign numbers to some of them
1263  // (i.e. they will remain at
1264  // invalid_dof_index). We need to
1265  // allocate the space because we
1266  // will want to be able to query
1267  // the dof_indices on each cell,
1268  // and simply be told that we don't
1269  // know them on some cell (i.e. get
1270  // back invalid_dof_index)
1271  clear_space();
1273 
1274  // hand things off to the policy
1275  number_cache = policy->distribute_dofs();
1276 
1277  // initialize the block info object
1278  // only if this is a sequential
1279  // triangulation. it doesn't work
1280  // correctly yet if it is parallel
1282  *>(&*tria) == nullptr)
1283  block_info_object.initialize(*this, false, true);
1284 }
1285 
1286 
1287 
1288 template <int dim, int spacedim>
1289 void
1291 {
1292  Assert(
1293  levels.size() > 0,
1294  ExcMessage(
1295  "Distribute active DoFs using distribute_dofs() before calling distribute_mg_dofs()."));
1296 
1297  Assert(
1298  ((tria->get_mesh_smoothing() &
1301  ExcMessage(
1302  "The mesh smoothing requirement 'limit_level_difference_at_vertices' has to be set for using multigrid!"));
1303 
1304  clear_mg_space();
1305 
1307  mg_number_cache = policy->distribute_mg_dofs();
1308 
1309  // initialize the block info object
1310  // only if this is a sequential
1311  // triangulation. it doesn't work
1312  // correctly yet if it is parallel
1313  if (dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
1314  &*tria) == nullptr)
1315  block_info_object.initialize(*this, true, false);
1316 }
1317 
1318 
1319 
1320 template <int dim, int spacedim>
1321 void
1323 {
1324  mg_levels.clear();
1325  mg_faces.reset();
1326 
1327  std::vector<MGVertexDoFs> tmp;
1328 
1329  std::swap(mg_vertex_dofs, tmp);
1330 
1331  mg_number_cache.clear();
1332 }
1333 
1334 
1335 template <int dim, int spacedim>
1336 void
1338 {
1340 }
1341 
1342 
1343 
1344 template <int dim, int spacedim>
1345 void
1347 {
1348  // release memory
1349  clear_space();
1350  clear_mg_space();
1351 }
1352 
1353 
1354 
1355 template <int dim, int spacedim>
1356 void
1358  const std::vector<types::global_dof_index> &new_numbers)
1359 {
1360  Assert(levels.size() > 0,
1361  ExcMessage(
1362  "You need to distribute DoFs before you can renumber them."));
1363 
1364 #ifdef DEBUG
1365  if (dynamic_cast<const parallel::shared::Triangulation<dim, spacedim> *>(
1366  &*tria) != nullptr)
1367  {
1368  Assert(new_numbers.size() == n_dofs() ||
1369  new_numbers.size() == n_locally_owned_dofs(),
1370  ExcMessage("Incorrect size of the input array."));
1371  }
1372  else if (dynamic_cast<
1374  &*tria) != nullptr)
1375  {
1376  AssertDimension(new_numbers.size(), n_locally_owned_dofs());
1377  }
1378  else
1379  {
1380  AssertDimension(new_numbers.size(), n_dofs());
1381  }
1382 
1383  // assert that the new indices are
1384  // consecutively numbered if we are
1385  // working on a single
1386  // processor. this doesn't need to
1387  // hold in the case of a parallel
1388  // mesh since we map the interval
1389  // [0...n_dofs()) into itself but
1390  // only globally, not on each
1391  // processor
1392  if (n_locally_owned_dofs() == n_dofs())
1393  {
1394  std::vector<types::global_dof_index> tmp(new_numbers);
1395  std::sort(tmp.begin(), tmp.end());
1396  std::vector<types::global_dof_index>::const_iterator p = tmp.begin();
1398  for (; p != tmp.end(); ++p, ++i)
1399  Assert(*p == i, ExcNewNumbersNotConsecutive(i));
1400  }
1401  else
1402  for (const auto new_number : new_numbers)
1403  Assert(new_number < n_dofs(),
1404  ExcMessage(
1405  "New DoF index is not less than the total number of dofs."));
1406 #endif
1407 
1408  number_cache = policy->renumber_dofs(new_numbers);
1409 }
1410 
1411 
1412 template <int dim, int spacedim>
1413 void
1415  const unsigned int level,
1416  const std::vector<types::global_dof_index> &new_numbers)
1417 {
1418  Assert(
1419  mg_levels.size() > 0 && levels.size() > 0,
1420  ExcMessage(
1421  "You need to distribute active and level DoFs before you can renumber level DoFs."));
1422  AssertIndexRange(level, get_triangulation().n_global_levels());
1423  AssertDimension(new_numbers.size(),
1425 
1426 #ifdef DEBUG
1427  // assert that the new indices are consecutively numbered if we are working
1428  // on a single processor. this doesn't need to hold in the case of a
1429  // parallel mesh since we map the interval [0...n_dofs(level)) into itself
1430  // but only globally, not on each processor
1431  if (n_locally_owned_dofs() == n_dofs())
1432  {
1433  std::vector<types::global_dof_index> tmp(new_numbers);
1434  std::sort(tmp.begin(), tmp.end());
1435  std::vector<types::global_dof_index>::const_iterator p = tmp.begin();
1437  for (; p != tmp.end(); ++p, ++i)
1438  Assert(*p == i, ExcNewNumbersNotConsecutive(i));
1439  }
1440  else
1441  for (const auto new_number : new_numbers)
1442  Assert(new_number < n_dofs(level),
1443  ExcMessage(
1444  "New DoF index is not less than the total number of dofs."));
1445 #endif
1446 
1447  mg_number_cache[level] = policy->renumber_mg_dofs(level, new_numbers);
1448 }
1449 
1450 
1451 
1452 template <int dim, int spacedim>
1453 unsigned int
1455 {
1458 }
1459 
1460 
1461 
1462 template <int dim, int spacedim>
1463 unsigned int
1465 {
1466  switch (dim)
1467  {
1468  case 1:
1469  return get_fe().dofs_per_vertex;
1470  case 2:
1471  return (3 * get_fe().dofs_per_vertex + 2 * get_fe().dofs_per_line);
1472  case 3:
1473  // we need to take refinement of
1474  // one boundary face into
1475  // consideration here; in fact,
1476  // this function returns what
1477  // #max_coupling_between_dofs<2>
1478  // returns
1479  //
1480  // we assume here, that only four
1481  // faces meet at the boundary;
1482  // this assumption is not
1483  // justified and needs to be
1484  // fixed some time. fortunately,
1485  // omitting it for now does no
1486  // harm since the matrix will cry
1487  // foul if its requirements are
1488  // not satisfied
1489  return (19 * get_fe().dofs_per_vertex + 28 * get_fe().dofs_per_line +
1490  8 * get_fe().dofs_per_quad);
1491  default:
1492  Assert(false, ExcNotImplemented());
1494  }
1495 }
1496 
1497 
1498 
1499 template <int dim, int spacedim>
1500 void
1502 {
1503  levels.clear();
1504  faces.reset();
1505 
1506  std::vector<types::global_dof_index> tmp;
1507  std::swap(vertex_dofs, tmp);
1508 
1509  number_cache.clear();
1510 }
1511 
1512 
1513 
1514 template <int dim, int spacedim>
1515 template <int structdim>
1517 DoFHandler<dim, spacedim>::get_dof_index(const unsigned int obj_level,
1518  const unsigned int obj_index,
1519  const unsigned int fe_index,
1520  const unsigned int local_index) const
1521 {
1523  *this,
1524  this->mg_levels[obj_level],
1525  this->mg_faces,
1526  obj_index,
1527  fe_index,
1528  local_index,
1529  std::integral_constant<int, structdim>());
1530 }
1531 
1532 
1533 
1534 template <int dim, int spacedim>
1535 template <int structdim>
1536 void
1538  const unsigned int obj_level,
1539  const unsigned int obj_index,
1540  const unsigned int fe_index,
1541  const unsigned int local_index,
1543 {
1545  *this,
1546  this->mg_levels[obj_level],
1547  this->mg_faces,
1548  obj_index,
1549  fe_index,
1550  local_index,
1551  global_index,
1552  std::integral_constant<int, structdim>());
1553 }
1554 
1555 
1556 
1557 template <int dim, int spacedim>
1559  : coarsest_level(numbers::invalid_unsigned_int)
1560  , finest_level(0)
1561 {}
1562 
1563 
1564 
1565 template <int dim, int spacedim>
1566 void
1568  const unsigned int cl,
1569  const unsigned int fl,
1570  const unsigned int dofs_per_vertex)
1571 {
1572  coarsest_level = cl;
1573  finest_level = fl;
1574 
1576  {
1577  const unsigned int n_levels = finest_level - coarsest_level + 1;
1578  const unsigned int n_indices = n_levels * dofs_per_vertex;
1579 
1580  indices = std::make_unique<types::global_dof_index[]>(n_indices);
1581  std::fill(indices.get(),
1582  indices.get() + n_indices,
1584  }
1585  else
1586  indices.reset();
1587 }
1588 
1589 
1590 
1591 template <int dim, int spacedim>
1592 unsigned int
1594 {
1595  return coarsest_level;
1596 }
1597 
1598 
1599 
1600 template <int dim, int spacedim>
1601 unsigned int
1603 {
1604  return finest_level;
1605 }
1606 
1607 
1608 /*-------------- Explicit Instantiations -------------------------------*/
1609 #include "dof_handler.inst"
1610 
1611 
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
unsigned int max_couplings_between_dofs() const
unsigned int get_coarsest_level() const
static void reserve_space_mg(DoFHandler< 2, spacedim > &dof_handler)
Definition: dof_handler.cc:460
std::vector< MGVertexDoFs > mg_vertex_dofs
Definition: dof_handler.h:1326
const Triangulation< dim, spacedim > & get_triangulation() const
static void reserve_space(DoFHandler< 3, spacedim > &dof_handler)
Definition: dof_handler.cc:353
level_cell_iterator end_mg() const
static const unsigned int invalid_unsigned_int
Definition: types.h:191
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1560
static void reserve_space_mg(DoFHandler< 1, spacedim > &dof_handler)
Definition: dof_handler.cc:392
cell_iterator begin(const unsigned int level=0) const
Definition: dof_handler.cc:918
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:779
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:668
cell_iterator end() const
Definition: dof_handler.cc:947
::internal::DoFHandlerImplementation::NumberCache number_cache
Definition: dof_handler.h:1204
#define AssertIndexRange(index, range)
Definition: exceptions.h:1628
virtual void clear()
static unsigned int max_couplings_between_dofs(const DoFHandler< 1, spacedim > &dof_handler)
Definition: dof_handler.cc:128
static unsigned int max_couplings_between_dofs(const DoFHandler< 2, spacedim > &dof_handler)
Definition: dof_handler.cc:140
IteratorRange< level_cell_iterator > mg_cell_iterators_on_level(const unsigned int level) const
void clear_mg_space()
#define AssertThrow(cond, exc)
Definition: exceptions.h:1513
virtual std::size_t memory_consumption() const
void initialize(const Triangulation< dim, spacedim > &tria, const FiniteElement< dim, spacedim > &fe)
Definition: dof_handler.cc:885
std::vector< std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel< dim > > > mg_levels
Definition: dof_handler.h:1338
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
unsigned int n_dofs_per_face() const
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:761
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:815
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
static ::ExceptionBase & ExcInvalidBoundaryIndicator()
active_cell_iterator begin_active(const unsigned int level=0) const
Definition: dof_handler.cc:931
types::global_dof_index n_locally_owned_dofs() const
BlockInfo block_info_object
Definition: dof_handler.h:1174
static ::ExceptionBase & ExcMessage(std::string arg1)
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:743
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:630
std::vector<::internal::DoFHandlerImplementation::NumberCache > mg_number_cache
Definition: dof_handler.h:1210
virtual void distribute_mg_dofs()
#define Assert(cond, exc)
Definition: exceptions.h:1403
hp::FECollection< dim, spacedim > fe_collection
Definition: dof_handler.h:1187
IteratorRange< active_cell_iterator > active_cell_iterators() const
active_cell_iterator end_active(const unsigned int level) const
Definition: dof_handler.cc:967
std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces< dim > > mg_faces
Definition: dof_handler.h:1349
types::global_dof_index n_boundary_dofs() const
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:361
unsigned int level
Definition: grid_out.cc:4355
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
Definition: dof_handler.h:1180
types::global_dof_index n_dofs() const
typename LevelSelector::cell_iterator level_cell_iterator
Definition: dof_handler.h:391
level_cell_iterator begin_mg(const unsigned int level=0) const
Definition: dof_handler.cc:980
std::unique_ptr<::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > > policy
Definition: dof_handler.h:1195
static unsigned int max_couplings_between_dofs(const DoFHandler< 3, spacedim > &dof_handler)
Definition: dof_handler.cc:247
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:706
virtual void set_fe(const FiniteElement< dim, spacedim > &fe)
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
static void reserve_space_mg(DoFHandler< 3, spacedim > &dof_handler)
Definition: dof_handler.cc:533
std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces< dim > > faces
Definition: dof_handler.h:1346
unsigned int get_finest_level() const
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:474
const types::global_dof_index * dummy()
Definition: dof_handler.cc:59
std::vector< std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel< dim > > > levels
Definition: dof_handler.h:1334
void initialize(const DoFHandler< dim, spacedim > &, bool levels_only=false, bool active_only=false)
Fill the object with values describing block structure of the DoFHandler.
Definition: block_info.cc:30
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
std::unique_ptr< types::global_dof_index[]> indices
Definition: dof_handler.h:1287
void initialize_local_block_info()
static void reserve_space(DoFHandler< 1, spacedim > &dof_handler)
Definition: dof_handler.cc:292
std::string policy_to_string(const ::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > &policy)
Definition: dof_handler.cc:71
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:687
void swap(MemorySpaceData< Number, MemorySpace > &, MemorySpaceData< Number, MemorySpace > &)
Definition: memory_space.h:102
void initialize_local(const DoFHandler< dim, spacedim > &)
Initialize block structure on cells and compute renumbering between cell dofs and block cell dofs...
Definition: block_info.cc:73
IteratorRange< level_cell_iterator > mg_cell_iterators() const
static void reserve_space(DoFHandler< 2, spacedim > &dof_handler)
Definition: dof_handler.cc:317
unsigned int coarsest_level
Definition: dof_handler.h:1271
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:360
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:797
T min(const T &t, const MPI_Comm &mpi_communicator)
void clear_space()
void init(const unsigned int coarsest_level, const unsigned int finest_level, const unsigned int dofs_per_vertex)
virtual void distribute_dofs(const FiniteElement< dim, spacedim > &fe)
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
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
static ::ExceptionBase & ExcNotImplemented()
Iterator points to a valid object.
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
bool is_element(const size_type index) const
Definition: index_set.h:1763
static ::ExceptionBase & ExcNewNumbersNotConsecutive(types::global_dof_index arg1)
typename ActiveSelector::cell_iterator cell_iterator
Definition: dof_handler.h:355
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:725
const types::boundary_id internal_face_boundary_id
Definition: types.h:250
unsigned int max_couplings_between_boundary_dofs() const
typename ActiveSelector::active_cell_iterator active_cell_iterator
Definition: dof_handler.h:327
const types::global_dof_index invalid_dof_index
Definition: types.h:206
IteratorState::IteratorStates state() const
virtual ~DoFHandler() override
Definition: dof_handler.cc:868
const IndexSet & locally_owned_dofs() const
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:611
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:649
unsigned int boundary_id
Definition: types.h:129
size_type n_elements() const
Definition: index_set.h:1830
IteratorRange< cell_iterator > cell_iterators() const
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
std::vector< types::global_dof_index > vertex_dofs
Definition: dof_handler.h:1320
static ::ExceptionBase & ExcInternalError()