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