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