Reference documentation for deal.II version Git 04ded8f 2017-09-19 10:11:38 +0200
dof_handler.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2016 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 at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #include <deal.II/base/memory_consumption.h>
17 #include <deal.II/dofs/dof_handler.h>
18 #include <deal.II/dofs/dof_handler_policy.h>
19 #include <deal.II/dofs/dof_levels.h>
20 #include <deal.II/dofs/dof_faces.h>
21 #include <deal.II/dofs/dof_accessor.h>
22 #include <deal.II/grid/tria_accessor.h>
23 #include <deal.II/grid/tria_iterator.h>
24 #include <deal.II/grid/tria_levels.h>
25 #include <deal.II/grid/tria.h>
26 #include <deal.II/base/geometry_info.h>
27 #include <deal.II/fe/fe.h>
28 #include <deal.II/distributed/shared_tria.h>
29 #include <deal.II/distributed/tria.h>
30 #include <deal.II/base/std_cxx14/memory.h>
31 
32 #include <set>
33 #include <algorithm>
34 
35 DEAL_II_NAMESPACE_OPEN
36 
37 
38 template <int dim, int spacedim>
39 const unsigned int DoFHandler<dim,spacedim>::dimension;
40 
41 template <int dim, int spacedim>
43 
44 template <int dim, int spacedim>
46 
47 template <int dim, int spacedim>
49 
50 
51 // reference the invalid_dof_index variable explicitly to work around
52 // a bug in the icc8 compiler
53 namespace internal
54 {
55  template <int dim, int spacedim>
56  const types::global_dof_index *dummy ()
57  {
58  return &::numbers::invalid_dof_index;
59  }
60 }
61 
62 
63 
64 namespace internal
65 {
66  template <int dim, int spacedim>
67  std::string policy_to_string(const ::internal::DoFHandler::Policy::PolicyBase<dim,spacedim> &policy)
68  {
69  std::string policy_name;
70  if (dynamic_cast<const typename ::internal::DoFHandler::Policy::Sequential<::DoFHandler<dim,spacedim> >*>(&policy)
71  || dynamic_cast<const typename ::internal::DoFHandler::Policy::Sequential<::hp::DoFHandler<dim,spacedim> >*>(&policy))
72  policy_name = "Policy::Sequential<";
73  else if (dynamic_cast<const typename ::internal::DoFHandler::Policy::ParallelDistributed<::DoFHandler<dim,spacedim> >*>(&policy)
74  || dynamic_cast<const typename ::internal::DoFHandler::Policy::ParallelDistributed<::hp::DoFHandler<dim,spacedim> >*>(&policy))
75  policy_name = "Policy::ParallelDistributed<";
76  else if (dynamic_cast<const typename ::internal::DoFHandler::Policy::ParallelShared<::DoFHandler<dim,spacedim> >*>(&policy)
77  || dynamic_cast<const typename ::internal::DoFHandler::Policy::ParallelShared<::hp::DoFHandler<dim,spacedim> >*>(&policy))
78  policy_name = "Policy::ParallelShared<";
79  else
81  policy_name += Utilities::int_to_string(dim)+
82  ","+Utilities::int_to_string(spacedim)+">";
83  return policy_name;
84  }
85 
86 
87  namespace DoFHandler
88  {
89  // access class
90  // ::DoFHandler instead of
91  // namespace internal::DoFHandler
92  using ::DoFHandler;
93 
94 
100  {
105  template <int spacedim>
106  static
107  unsigned int
109  {
110  return std::min(static_cast<types::global_dof_index>(3*dof_handler.get_fe().dofs_per_vertex +
111  2*dof_handler.get_fe().dofs_per_line),
112  dof_handler.n_dofs());
113  }
114 
115 
116 
117  template <int spacedim>
118  static
119  unsigned int
121  {
122 
123  // get these numbers by drawing pictures
124  // and counting...
125  // example:
126  // | | |
127  // --x-----x--x--X--
128  // | | | |
129  // | x--x--x
130  // | | | |
131  // --x--x--*--x--x--
132  // | | | |
133  // x--x--x |
134  // | | | |
135  // --X--x--x-----x--
136  // | | |
137  // x = vertices connected with center vertex *;
138  // = total of 19
139  // (the X vertices are connected with * if
140  // the vertices adjacent to X are hanging
141  // nodes)
142  // count lines -> 28 (don't forget to count
143  // mother and children separately!)
144  types::global_dof_index max_couplings;
145  switch (dof_handler.tria->max_adjacent_cells())
146  {
147  case 4:
148  max_couplings=19*dof_handler.get_fe().dofs_per_vertex +
149  28*dof_handler.get_fe().dofs_per_line +
150  8*dof_handler.get_fe().dofs_per_quad;
151  break;
152  case 5:
153  max_couplings=21*dof_handler.get_fe().dofs_per_vertex +
154  31*dof_handler.get_fe().dofs_per_line +
155  9*dof_handler.get_fe().dofs_per_quad;
156  break;
157  case 6:
158  max_couplings=28*dof_handler.get_fe().dofs_per_vertex +
159  42*dof_handler.get_fe().dofs_per_line +
160  12*dof_handler.get_fe().dofs_per_quad;
161  break;
162  case 7:
163  max_couplings=30*dof_handler.get_fe().dofs_per_vertex +
164  45*dof_handler.get_fe().dofs_per_line +
165  13*dof_handler.get_fe().dofs_per_quad;
166  break;
167  case 8:
168  max_couplings=37*dof_handler.get_fe().dofs_per_vertex +
169  56*dof_handler.get_fe().dofs_per_line +
170  16*dof_handler.get_fe().dofs_per_quad;
171  break;
172 
173  // the following numbers are not based on actual counting but by
174  // extrapolating the number sequences from the previous ones (for
175  // example, for dofs_per_vertex, the sequence above is 19, 21, 28,
176  // 30, 37, and is continued as follows):
177  case 9:
178  max_couplings=39*dof_handler.get_fe().dofs_per_vertex +
179  59*dof_handler.get_fe().dofs_per_line +
180  17*dof_handler.get_fe().dofs_per_quad;
181  break;
182  case 10:
183  max_couplings=46*dof_handler.get_fe().dofs_per_vertex +
184  70*dof_handler.get_fe().dofs_per_line +
185  20*dof_handler.get_fe().dofs_per_quad;
186  break;
187  case 11:
188  max_couplings=48*dof_handler.get_fe().dofs_per_vertex +
189  73*dof_handler.get_fe().dofs_per_line +
190  21*dof_handler.get_fe().dofs_per_quad;
191  break;
192  case 12:
193  max_couplings=55*dof_handler.get_fe().dofs_per_vertex +
194  84*dof_handler.get_fe().dofs_per_line +
195  24*dof_handler.get_fe().dofs_per_quad;
196  break;
197  case 13:
198  max_couplings=57*dof_handler.get_fe().dofs_per_vertex +
199  87*dof_handler.get_fe().dofs_per_line +
200  25*dof_handler.get_fe().dofs_per_quad;
201  break;
202  case 14:
203  max_couplings=63*dof_handler.get_fe().dofs_per_vertex +
204  98*dof_handler.get_fe().dofs_per_line +
205  28*dof_handler.get_fe().dofs_per_quad;
206  break;
207  case 15:
208  max_couplings=65*dof_handler.get_fe().dofs_per_vertex +
209  103*dof_handler.get_fe().dofs_per_line +
210  29*dof_handler.get_fe().dofs_per_quad;
211  break;
212  case 16:
213  max_couplings=72*dof_handler.get_fe().dofs_per_vertex +
214  114*dof_handler.get_fe().dofs_per_line +
215  32*dof_handler.get_fe().dofs_per_quad;
216  break;
217 
218  default:
219  Assert (false, ExcNotImplemented());
220  max_couplings=0;
221  }
222  return std::min(max_couplings,dof_handler.n_dofs());
223  }
224 
225 
226  template <int spacedim>
227  static
228  unsigned int
230  {
231 //TODO:[?] Invent significantly better estimates than the ones in this function
232 
233  // doing the same thing here is a
234  // rather complicated thing, compared
235  // to the 2d case, since it is hard
236  // to draw pictures with several
237  // refined hexahedra :-) so I
238  // presently only give a coarse
239  // estimate for the case that at most
240  // 8 hexes meet at each vertex
241  //
242  // can anyone give better estimate
243  // here?
244  const unsigned int max_adjacent_cells
245  = dof_handler.tria->max_adjacent_cells();
246 
247  types::global_dof_index max_couplings;
248  if (max_adjacent_cells <= 8)
249  max_couplings=7*7*7*dof_handler.get_fe().dofs_per_vertex +
250  7*6*7*3*dof_handler.get_fe().dofs_per_line +
251  9*4*7*3*dof_handler.get_fe().dofs_per_quad +
252  27*dof_handler.get_fe().dofs_per_hex;
253  else
254  {
255  Assert (false, ExcNotImplemented());
256  max_couplings=0;
257  }
258 
259  return std::min(max_couplings,dof_handler.n_dofs());
260  }
261 
262 
272  template <int spacedim>
273  static
275  {
276  dof_handler.vertex_dofs
277  .resize(dof_handler.tria->n_vertices() *
278  dof_handler.get_fe().dofs_per_vertex,
280 
281  for (unsigned int i=0; i<dof_handler.tria->n_levels(); ++i)
282  {
283  dof_handler.levels.emplace_back (new internal::DoFHandler::DoFLevel<1>);
284 
285  dof_handler.levels.back()->dof_object.dofs
286  .resize (dof_handler.tria->n_raw_cells(i) *
287  dof_handler.get_fe().dofs_per_line,
289 
290  dof_handler.levels.back()->cell_dof_indices_cache
291  .resize (dof_handler.tria->n_raw_cells(i) *
292  dof_handler.get_fe().dofs_per_cell,
294  }
295  }
296 
297 
298  template <int spacedim>
299  static
300  void reserve_space (DoFHandler<2,spacedim> &dof_handler)
301  {
302  dof_handler.vertex_dofs
303  .resize(dof_handler.tria->n_vertices() *
304  dof_handler.get_fe().dofs_per_vertex,
306 
307  for (unsigned int i=0; i<dof_handler.tria->n_levels(); ++i)
308  {
309  dof_handler.levels.emplace_back (new internal::DoFHandler::DoFLevel<2>);
310 
311  dof_handler.levels.back()->dof_object.dofs
312  .resize (dof_handler.tria->n_raw_cells(i) *
313  dof_handler.get_fe().dofs_per_quad,
315 
316  dof_handler.levels.back()->cell_dof_indices_cache
317  .resize (dof_handler.tria->n_raw_cells(i) *
318  dof_handler.get_fe().dofs_per_cell,
320  }
321 
322  dof_handler.faces.reset (new internal::DoFHandler::DoFFaces<2>);
323  // avoid access to n_raw_lines when there are no cells
324  if (dof_handler.tria->n_cells() > 0)
325  {
326  dof_handler.faces->lines.dofs
327  .resize (dof_handler.tria->n_raw_lines() *
328  dof_handler.get_fe().dofs_per_line,
330  }
331  }
332 
333 
334  template <int spacedim>
335  static
336  void reserve_space (DoFHandler<3,spacedim> &dof_handler)
337  {
338  dof_handler.vertex_dofs
339  .resize(dof_handler.tria->n_vertices() *
340  dof_handler.get_fe().dofs_per_vertex,
342 
343  for (unsigned int i=0; i<dof_handler.tria->n_levels(); ++i)
344  {
345  dof_handler.levels.emplace_back (new internal::DoFHandler::DoFLevel<3>);
346 
347  dof_handler.levels.back()->dof_object.dofs
348  .resize (dof_handler.tria->n_raw_cells(i) *
349  dof_handler.get_fe().dofs_per_hex,
351 
352  dof_handler.levels.back()->cell_dof_indices_cache
353  .resize (dof_handler.tria->n_raw_cells(i) *
354  dof_handler.get_fe().dofs_per_cell,
356  }
357  dof_handler.faces.reset (new internal::DoFHandler::DoFFaces<3>);
358 
359  // avoid access to n_raw_lines when there are no cells
360  if (dof_handler.tria->n_cells() > 0)
361  {
362  dof_handler.faces->lines.dofs
363  .resize (dof_handler.tria->n_raw_lines() *
364  dof_handler.get_fe().dofs_per_line,
366  dof_handler.faces->quads.dofs
367  .resize (dof_handler.tria->n_raw_quads() *
368  dof_handler.get_fe().dofs_per_quad,
370  }
371  }
372 
373  template <int spacedim>
374  static
375  void reserve_space_mg (DoFHandler<1, spacedim> &dof_handler)
376  {
377  Assert (dof_handler.get_triangulation().n_levels () > 0, ExcMessage ("Invalid triangulation"));
378  dof_handler.clear_mg_space ();
379 
380  const ::Triangulation<1, spacedim> &tria = dof_handler.get_triangulation();
381  const unsigned int &dofs_per_line = dof_handler.get_fe ().dofs_per_line;
382  const unsigned int &n_levels = tria.n_levels ();
383 
384  for (unsigned int i = 0; i < n_levels; ++i)
385  {
386  dof_handler.mg_levels.emplace_back (new internal::DoFHandler::DoFLevel<1>);
387  dof_handler.mg_levels.back ()->dof_object.dofs = std::vector<types::global_dof_index> (tria.n_raw_lines (i) * dofs_per_line, numbers::invalid_dof_index);
388  }
389 
390  const unsigned int &n_vertices = tria.n_vertices ();
391 
392  dof_handler.mg_vertex_dofs.resize (n_vertices);
393 
394  std::vector<unsigned int> max_level (n_vertices, 0);
395  std::vector<unsigned int> min_level (n_vertices, n_levels);
396 
397  for (typename ::Triangulation<1, spacedim>::cell_iterator cell = tria.begin (); cell != tria.end (); ++cell)
398  {
399  const unsigned int level = cell->level ();
400 
401  for (unsigned int vertex = 0; vertex < GeometryInfo<1>::vertices_per_cell; ++vertex)
402  {
403  const unsigned int vertex_index = cell->vertex_index (vertex);
404 
405  if (min_level[vertex_index] > level)
406  min_level[vertex_index] = level;
407 
408  if (max_level[vertex_index] < level)
409  max_level[vertex_index] = level;
410  }
411  }
412 
413  for (unsigned int vertex = 0; vertex < n_vertices; ++vertex)
414  if (tria.vertex_used (vertex))
415  {
416  Assert (min_level[vertex] < n_levels, ExcInternalError ());
417  Assert (max_level[vertex] >= min_level[vertex], ExcInternalError ());
418  dof_handler.mg_vertex_dofs[vertex].init (min_level[vertex], max_level[vertex], dof_handler.get_fe ().dofs_per_vertex);
419  }
420 
421  else
422  {
423  Assert (min_level[vertex] == n_levels, ExcInternalError ());
424  Assert (max_level[vertex] == 0, ExcInternalError ());
425  dof_handler.mg_vertex_dofs[vertex].init (1, 0, 0);
426  }
427  }
428 
429  template <int spacedim>
430  static
431  void reserve_space_mg (DoFHandler<2, spacedim> &dof_handler)
432  {
433  Assert (dof_handler.get_triangulation().n_levels () > 0, ExcMessage ("Invalid triangulation"));
434  dof_handler.clear_mg_space ();
435 
436  const ::FiniteElement<2, spacedim> &fe = dof_handler.get_fe ();
437  const ::Triangulation<2, spacedim> &tria = dof_handler.get_triangulation();
438  const unsigned int &n_levels = tria.n_levels ();
439 
440  for (unsigned int i = 0; i < n_levels; ++i)
441  {
442  dof_handler.mg_levels.emplace_back (new internal::DoFHandler::DoFLevel<2>);
443  dof_handler.mg_levels.back ()->dof_object.dofs = std::vector<types::global_dof_index> (tria.n_raw_quads (i) * fe.dofs_per_quad, numbers::invalid_dof_index);
444  }
445 
446  dof_handler.mg_faces.reset (new internal::DoFHandler::DoFFaces<2>);
447  dof_handler.mg_faces->lines.dofs = std::vector<types::global_dof_index> (tria.n_raw_lines () * fe.dofs_per_line, numbers::invalid_dof_index);
448 
449  const unsigned int &n_vertices = tria.n_vertices ();
450 
451  dof_handler.mg_vertex_dofs.resize (n_vertices);
452 
453  std::vector<unsigned int> max_level (n_vertices, 0);
454  std::vector<unsigned int> min_level (n_vertices, n_levels);
455 
456  for (typename ::Triangulation<2, spacedim>::cell_iterator cell = tria.begin (); cell != tria.end (); ++cell)
457  {
458  const unsigned int level = cell->level ();
459 
460  for (unsigned int vertex = 0; vertex < GeometryInfo<2>::vertices_per_cell; ++vertex)
461  {
462  const unsigned int vertex_index = cell->vertex_index (vertex);
463 
464  if (min_level[vertex_index] > level)
465  min_level[vertex_index] = level;
466 
467  if (max_level[vertex_index] < level)
468  max_level[vertex_index] = level;
469  }
470  }
471 
472  for (unsigned int vertex = 0; vertex < n_vertices; ++vertex)
473  if (tria.vertex_used (vertex))
474  {
475  Assert (min_level[vertex] < n_levels, ExcInternalError ());
476  Assert (max_level[vertex] >= min_level[vertex], ExcInternalError ());
477  dof_handler.mg_vertex_dofs[vertex].init (min_level[vertex], max_level[vertex], fe.dofs_per_vertex);
478  }
479 
480  else
481  {
482  Assert (min_level[vertex] == n_levels, ExcInternalError ());
483  Assert (max_level[vertex] == 0, ExcInternalError ());
484  dof_handler.mg_vertex_dofs[vertex].init (1, 0, 0);
485  }
486  }
487 
488  template <int spacedim>
489  static
490  void reserve_space_mg (DoFHandler<3, spacedim> &dof_handler)
491  {
492  Assert (dof_handler.get_triangulation().n_levels () > 0, ExcMessage ("Invalid triangulation"));
493  dof_handler.clear_mg_space ();
494 
495  const ::FiniteElement<3, spacedim> &fe = dof_handler.get_fe ();
496  const ::Triangulation<3, spacedim> &tria = dof_handler.get_triangulation();
497  const unsigned int &n_levels = tria.n_levels ();
498 
499  for (unsigned int i = 0; i < n_levels; ++i)
500  {
501  dof_handler.mg_levels.emplace_back (new internal::DoFHandler::DoFLevel<3>);
502  dof_handler.mg_levels.back ()->dof_object.dofs = std::vector<types::global_dof_index> (tria.n_raw_hexs (i) * fe.dofs_per_hex, numbers::invalid_dof_index);
503  }
504 
505  dof_handler.mg_faces.reset (new internal::DoFHandler::DoFFaces<3>);
506  dof_handler.mg_faces->lines.dofs = std::vector<types::global_dof_index> (tria.n_raw_lines () * fe.dofs_per_line, numbers::invalid_dof_index);
507  dof_handler.mg_faces->quads.dofs = std::vector<types::global_dof_index> (tria.n_raw_quads () * fe.dofs_per_quad, numbers::invalid_dof_index);
508 
509  const unsigned int &n_vertices = tria.n_vertices ();
510 
511  dof_handler.mg_vertex_dofs.resize (n_vertices);
512 
513  std::vector<unsigned int> max_level (n_vertices, 0);
514  std::vector<unsigned int> min_level (n_vertices, n_levels);
515 
516  for (typename ::Triangulation<3, spacedim>::cell_iterator cell = tria.begin (); cell != tria.end (); ++cell)
517  {
518  const unsigned int level = cell->level ();
519 
520  for (unsigned int vertex = 0; vertex < GeometryInfo<3>::vertices_per_cell; ++vertex)
521  {
522  const unsigned int vertex_index = cell->vertex_index (vertex);
523 
524  if (min_level[vertex_index] > level)
525  min_level[vertex_index] = level;
526 
527  if (max_level[vertex_index] < level)
528  max_level[vertex_index] = level;
529  }
530  }
531 
532  for (unsigned int vertex = 0; vertex < n_vertices; ++vertex)
533  if (tria.vertex_used (vertex))
534  {
535  Assert (min_level[vertex] < n_levels, ExcInternalError ());
536  Assert (max_level[vertex] >= min_level[vertex], ExcInternalError ());
537  dof_handler.mg_vertex_dofs[vertex].init (min_level[vertex], max_level[vertex], fe.dofs_per_vertex);
538  }
539 
540  else
541  {
542  Assert (min_level[vertex] == n_levels, ExcInternalError ());
543  Assert (max_level[vertex] == 0, ExcInternalError ());
544  dof_handler.mg_vertex_dofs[vertex].init (1, 0, 0);
545  }
546  }
547 
548 
549 
550  template <int spacedim>
551  static
553  get_dof_index (
554  const DoFHandler<1, spacedim> &dof_handler,
557  const unsigned int obj_index,
558  const unsigned int fe_index,
559  const unsigned int local_index,
560  const std::integral_constant<int, 1>)
561  {
562  return mg_level.dof_object.get_dof_index (dof_handler, obj_index, fe_index, local_index);
563  }
564 
565  template <int spacedim>
566  static
568  get_dof_index (const DoFHandler<2, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<2> &, internal::DoFHandler::DoFFaces<2> &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant<int, 1>)
569  {
570  return mg_faces.lines.get_dof_index (dof_handler, obj_index, fe_index, local_index);
571  }
572 
573  template <int spacedim>
574  static
576  get_dof_index (const DoFHandler<2, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<2> &mg_level, internal::DoFHandler::DoFFaces<2> &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant<int, 2>)
577  {
578  return mg_level.dof_object.get_dof_index (dof_handler, obj_index, fe_index, local_index);
579  }
580 
581  template <int spacedim>
582  static
584  get_dof_index (const DoFHandler<3, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<3> &, internal::DoFHandler::DoFFaces<3> &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant<int, 1>)
585  {
586  return mg_faces.lines.get_dof_index (dof_handler, obj_index, fe_index, local_index);
587  }
588 
589  template <int spacedim>
590  static
592  get_dof_index (const DoFHandler<3, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<3> &, internal::DoFHandler::DoFFaces<3> &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant<int, 2>)
593  {
594  return mg_faces.quads.get_dof_index (dof_handler, obj_index, fe_index, local_index);
595  }
596 
597  template <int spacedim>
598  static
600  get_dof_index (const DoFHandler<3, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<3> &mg_level, internal::DoFHandler::DoFFaces<3> &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const std::integral_constant<int, 3>)
601  {
602  return mg_level.dof_object.get_dof_index (dof_handler, obj_index, fe_index, local_index);
603  }
604 
605  template <int spacedim>
606  static
607  void set_dof_index (const DoFHandler<1, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<1> &mg_level, internal::DoFHandler::DoFFaces<1> &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant<int, 1>)
608  {
609  mg_level.dof_object.set_dof_index (dof_handler, obj_index, fe_index, local_index, global_index);
610  }
611 
612  template <int spacedim>
613  static
614  void set_dof_index (const DoFHandler<2, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<2> &, internal::DoFHandler::DoFFaces<2> &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant<int, 1>)
615  {
616  mg_faces.lines.set_dof_index (dof_handler, obj_index, fe_index, local_index, global_index);
617  }
618 
619  template <int spacedim>
620  static
621  void set_dof_index (const DoFHandler<2, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<2> &mg_level, internal::DoFHandler::DoFFaces<2> &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant<int, 2>)
622  {
623  mg_level.dof_object.set_dof_index (dof_handler, obj_index, fe_index, local_index, global_index);
624  }
625 
626  template <int spacedim>
627  static
628  void set_dof_index (const DoFHandler<3, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<3> &, internal::DoFHandler::DoFFaces<3> &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant<int, 1>)
629  {
630  mg_faces.lines.set_dof_index (dof_handler, obj_index, fe_index, local_index, global_index);
631  }
632 
633  template <int spacedim>
634  static
635  void set_dof_index (const DoFHandler<3, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<3> &, internal::DoFHandler::DoFFaces<3> &mg_faces, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant<int, 2>)
636  {
637  mg_faces.quads.set_dof_index (dof_handler, obj_index, fe_index, local_index, global_index);
638  }
639 
640  template <int spacedim>
641  static
642  void set_dof_index (const DoFHandler<3, spacedim> &dof_handler, internal::DoFHandler::DoFLevel<3> &mg_level, internal::DoFHandler::DoFFaces<3> &, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const std::integral_constant<int, 3>)
643  {
644  mg_level.dof_object.set_dof_index (dof_handler, obj_index, fe_index, local_index, global_index);
645  }
646  };
647  }
648 }
649 
650 
651 
652 template <int dim, int spacedim>
654  :
655  tria(&tria, typeid(*this).name()),
656  fe_collection(nullptr),
657  faces(nullptr),
658  mg_faces (nullptr)
659 {
660  // decide whether we need a sequential or a parallel distributed policy
661  if (dynamic_cast<const parallel::shared::Triangulation< dim, spacedim>*>
662  (&tria)
663  != nullptr)
665  else if (dynamic_cast<const parallel::distributed::Triangulation< dim, spacedim >*>
666  (&tria)
667  == nullptr)
669  else
671 }
672 
673 
674 template <int dim, int spacedim>
676  :
677  tria(nullptr, typeid(*this).name()),
678  fe_collection(nullptr)
679 {}
680 
681 
682 template <int dim, int spacedim>
684 {
685  // release allocated memory
686  clear ();
687 
688  // also release the policy. this needs to happen before the
689  // current object disappears because the policy objects
690  // store references to the DoFhandler object they work on
691  policy.reset ();
692 }
693 
694 
695 template <int dim, int spacedim>
696 void
699  const FiniteElement<dim,spacedim> &fe)
700 {
701  tria = &t;
702  faces = nullptr;
703  number_cache.n_global_dofs = 0;
704 
705  // decide whether we need a sequential or a parallel distributed policy
706  if (dynamic_cast<const parallel::shared::Triangulation< dim, spacedim>*> (&t) != nullptr)
708  else if (dynamic_cast<const parallel::distributed::Triangulation< dim, spacedim >*> (&t) != nullptr)
710  else
712 
713  distribute_dofs(fe);
714 }
715 
716 
717 
718 /*------------------------ Cell iterator functions ------------------------*/
719 
720 template <int dim, int spacedim>
722 DoFHandler<dim,spacedim>::begin (const unsigned int level) const
723 {
724  typename Triangulation<dim,spacedim>::cell_iterator cell = this->get_triangulation().begin(level);
725  if (cell == this->get_triangulation().end(level))
726  return end(level);
727  return cell_iterator (*cell, this);
728 }
729 
730 
731 
732 template <int dim, int spacedim>
734 DoFHandler<dim,spacedim>::begin_active (const unsigned int level) const
735 {
736  // level is checked in begin
737  cell_iterator i = begin (level);
738  if (i.state() != IteratorState::valid)
739  return i;
740  while (i->has_children())
741  if ((++i).state() != IteratorState::valid)
742  return i;
743  return i;
744 }
745 
746 
747 
748 template <int dim, int spacedim>
751 {
752  return cell_iterator (&this->get_triangulation(),
753  -1,
754  -1,
755  this);
756 }
757 
758 
759 template <int dim, int spacedim>
761 DoFHandler<dim,spacedim>::end (const unsigned int level) const
762 {
763  typename Triangulation<dim,spacedim>::cell_iterator cell = this->get_triangulation().end(level);
764  if (cell.state() != IteratorState::valid)
765  return end();
766  return cell_iterator (*cell, this);
767 }
768 
769 
770 template <int dim, int spacedim>
772 DoFHandler<dim, spacedim>::end_active (const unsigned int level) const
773 {
774  typename Triangulation<dim,spacedim>::cell_iterator cell = this->get_triangulation().end_active(level);
775  if (cell.state() != IteratorState::valid)
776  return active_cell_iterator(end());
777  return active_cell_iterator (*cell, this);
778 }
779 
780 
781 
782 template <int dim, int spacedim>
783 typename DoFHandler<dim, spacedim>::level_cell_iterator
784 DoFHandler<dim, spacedim>::begin_mg (const unsigned int level) const
785 {
786  // Assert(this->has_level_dofs(), ExcMessage("You can only iterate over mg "
787  // "levels if mg dofs got distributed."));
788  typename Triangulation<dim,spacedim>::cell_iterator cell = this->get_triangulation().begin(level);
789  if (cell == this->get_triangulation().end(level))
790  return end_mg(level);
791  return level_cell_iterator (*cell, this);
792 }
793 
794 
795 template <int dim, int spacedim>
796 typename DoFHandler<dim, spacedim>::level_cell_iterator
797 DoFHandler<dim, spacedim>::end_mg (const unsigned int level) const
798 {
799  // Assert(this->has_level_dofs(), ExcMessage("You can only iterate over mg "
800  // "levels if mg dofs got distributed."));
801  typename Triangulation<dim,spacedim>::cell_iterator cell = this->get_triangulation().end(level);
802  if (cell.state() != IteratorState::valid)
803  return end();
804  return level_cell_iterator (*cell, this);
805 }
806 
807 
808 template <int dim, int spacedim>
809 typename DoFHandler<dim, spacedim>::level_cell_iterator
811 {
812  return level_cell_iterator (&this->get_triangulation(), -1, -1, this);
813 }
814 
815 
816 
817 template <int dim, int spacedim>
820 {
821  return
823  (begin(), end());
824 }
825 
826 
827 template <int dim, int spacedim>
830 {
831  return
833  (begin_active(), end());
834 }
835 
836 
837 
838 template <int dim, int spacedim>
841 {
842  return
844  (begin_mg(), end_mg());
845 }
846 
847 
848 
849 
850 template <int dim, int spacedim>
853 {
854  return
856  (begin(level), end(level));
857 }
858 
859 
860 
861 template <int dim, int spacedim>
864 {
865  return
867  (begin_active(level), end_active(level));
868 }
869 
870 
871 
872 template <int dim, int spacedim>
875 {
876  return
878  (begin_mg(level), end_mg(level));
879 }
880 
881 
882 
883 //---------------------------------------------------------------------------
884 
885 
886 
887 template <int dim, int spacedim>
889 {
890  std::set<int> boundary_dofs;
891 
892  const unsigned int dofs_per_face = get_fe().dofs_per_face;
893  std::vector<types::global_dof_index> dofs_on_face(dofs_per_face);
894 
895  // loop over all faces of all cells
896  // and see whether they are at a
897  // boundary. note (i) that we visit
898  // interior faces twice (which we
899  // don't care about) but exterior
900  // faces only once as is
901  // appropriate, and (ii) that we
902  // need not take special care of
903  // single lines (using
904  // @p{cell->has_boundary_lines}),
905  // since we do not support
906  // boundaries of dimension dim-2,
907  // and so every boundary line is
908  // also part of a boundary face.
909  active_cell_iterator cell = begin_active (),
910  endc = end();
911  for (; cell!=endc; ++cell)
912  for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
913  if (cell->at_boundary(f))
914  {
915  cell->face(f)->get_dof_indices (dofs_on_face);
916  for (unsigned int i=0; i<dofs_per_face; ++i)
917  boundary_dofs.insert(dofs_on_face[i]);
918  }
919 
920  return boundary_dofs.size();
921 }
922 
923 
924 
925 template <int dim, int spacedim>
927 DoFHandler<dim,spacedim>::n_boundary_dofs (const std::set<types::boundary_id> &boundary_ids) const
928 {
929  Assert (boundary_ids.find (numbers::internal_face_boundary_id) == boundary_ids.end(),
931 
932  std::set<types::global_dof_index> boundary_dofs;
933 
934  const unsigned int dofs_per_face = get_fe().dofs_per_face;
935  std::vector<types::global_dof_index> dofs_on_face(dofs_per_face);
936 
937  // same as in the previous
938  // function, but with a different
939  // check for the boundary indicator
940  active_cell_iterator cell = begin_active (),
941  endc = end();
942  for (; cell!=endc; ++cell)
943  for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
944  if (cell->at_boundary(f)
945  &&
946  (std::find (boundary_ids.begin(),
947  boundary_ids.end(),
948  cell->face(f)->boundary_id()) !=
949  boundary_ids.end()))
950  {
951  cell->face(f)->get_dof_indices (dofs_on_face);
952  for (unsigned int i=0; i<dofs_per_face; ++i)
953  boundary_dofs.insert(dofs_on_face[i]);
954  }
955 
956  return boundary_dofs.size();
957 }
958 
959 
960 
961 template <int dim, int spacedim>
962 std::size_t
964 {
965  std::size_t mem = (MemoryConsumption::memory_consumption (tria) +
966  MemoryConsumption::memory_consumption (fe_collection) +
967  MemoryConsumption::memory_consumption (block_info_object) +
971  sizeof (number_cache) +
972  MemoryConsumption::memory_consumption (mg_number_cache) +
974  for (unsigned int i=0; i<levels.size(); ++i)
975  mem += MemoryConsumption::memory_consumption (*levels[i]);
976 
977  for (unsigned int level = 0; level < mg_levels.size (); ++level)
978  mem += mg_levels[level]->memory_consumption ();
979 
980  if (mg_faces != nullptr)
981  mem += MemoryConsumption::memory_consumption (*mg_faces);
982 
983  for (unsigned int i = 0; i < mg_vertex_dofs.size (); ++i)
984  mem += sizeof (MGVertexDoFs) + (1 + mg_vertex_dofs[i].get_finest_level () - mg_vertex_dofs[i].get_coarsest_level ()) * sizeof (types::global_dof_index);
985 
986  return mem;
987 }
988 
989 
990 
991 template <int dim, int spacedim>
993 {
994  fe_collection = std_cxx14::make_unique<hp::FECollection<dim, spacedim>>(ff);
995 
996  // delete all levels and set them
997  // up newly. note that we still
998  // have to allocate space for all
999  // degrees of freedom on this mesh
1000  // (including ghost and cells that
1001  // are entirely stored on different
1002  // processors), though we may not
1003  // assign numbers to some of them
1004  // (i.e. they will remain at
1005  // invalid_dof_index). We need to
1006  // allocate the space because we
1007  // will want to be able to query
1008  // the dof_indices on each cell,
1009  // and simply be told that we don't
1010  // know them on some cell (i.e. get
1011  // back invalid_dof_index)
1012  clear_space ();
1014 
1015  // hand things off to the policy
1016  number_cache = policy->distribute_dofs ();
1017 
1018  // initialize the block info object
1019  // only if this is a sequential
1020  // triangulation. it doesn't work
1021  // correctly yet if it is parallel
1022  if (dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim>*>(&*tria) == nullptr)
1023  block_info_object.initialize(*this, false, true);
1024 }
1025 
1026 
1027 
1028 template <int dim, int spacedim>
1030 {
1031  this->distribute_mg_dofs();
1032 }
1033 
1034 
1035 
1036 template <int dim, int spacedim>
1038 {
1039  Assert(levels.size()>0, ExcMessage("Distribute active DoFs using distribute_dofs() before calling distribute_mg_dofs()."));
1040 
1043  ExcMessage("The mesh smoothing requirement 'limit_level_difference_at_vertices' has to be set for using multigrid!"));
1044 
1045  clear_mg_space();
1046 
1047  internal::DoFHandler::Implementation::reserve_space_mg (*this);
1048  mg_number_cache = policy->distribute_mg_dofs ();
1049 
1050  // initialize the block info object
1051  // only if this is a sequential
1052  // triangulation. it doesn't work
1053  // correctly yet if it is parallel
1054  if (dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim>*>(&*tria) == nullptr)
1055  block_info_object.initialize (*this, true, false);
1056 }
1057 
1058 
1059 
1060 template <int dim, int spacedim>
1062 {
1063  //TODO: move this to distribute_mg_dofs and remove function
1064 }
1065 
1066 
1067 
1068 template <int dim, int spacedim>
1070 {
1071  mg_levels.clear ();
1072  mg_faces.reset ();
1073 
1074  std::vector<MGVertexDoFs> tmp;
1075 
1076  std::swap (mg_vertex_dofs, tmp);
1077 
1078  mg_number_cache.clear();
1079 }
1080 
1081 
1082 template <int dim, int spacedim>
1084 {
1085  block_info_object.initialize_local(*this);
1086 }
1087 
1088 
1089 
1090 template <int dim, int spacedim>
1092 {
1093  // release lock to old fe
1094  fe_collection.release();
1095 
1096  // release memory
1097  clear_space ();
1098  clear_mg_space ();
1099 }
1100 
1101 
1102 
1103 template <int dim, int spacedim>
1104 void
1105 DoFHandler<dim,spacedim>::renumber_dofs (const std::vector<types::global_dof_index> &new_numbers)
1106 {
1107  Assert(levels.size()>0, ExcMessage("You need to distribute DoFs before you can renumber them."));
1108 
1109  Assert (new_numbers.size() == n_locally_owned_dofs(),
1111 
1112 #ifdef DEBUG
1113  // assert that the new indices are
1114  // consecutively numbered if we are
1115  // working on a single
1116  // processor. this doesn't need to
1117  // hold in the case of a parallel
1118  // mesh since we map the interval
1119  // [0...n_dofs()) into itself but
1120  // only globally, not on each
1121  // processor
1122  if (n_locally_owned_dofs() == n_dofs())
1123  {
1124  std::vector<types::global_dof_index> tmp(new_numbers);
1125  std::sort (tmp.begin(), tmp.end());
1126  std::vector<types::global_dof_index>::const_iterator p = tmp.begin();
1128  for (; p!=tmp.end(); ++p, ++i)
1129  Assert (*p == i, ExcNewNumbersNotConsecutive(i));
1130  }
1131  else
1132  for (types::global_dof_index i=0; i<new_numbers.size(); ++i)
1133  Assert (new_numbers[i] < n_dofs(),
1134  ExcMessage ("New DoF index is not less than the total number of dofs."));
1135 #endif
1136 
1137  number_cache = policy->renumber_dofs (new_numbers);
1138 }
1139 
1140 
1141 template <int dim, int spacedim>
1142 void
1144  const std::vector<types::global_dof_index> &new_numbers)
1145 {
1146  Assert(mg_levels.size()>0 && levels.size()>0,
1147  ExcMessage("You need to distribute active and level DoFs before you can renumber level DoFs."));
1148  AssertIndexRange(level, levels.size());
1149  Assert (new_numbers.size() == n_dofs(level), ExcRenumberingIncomplete());
1150 
1151  mg_number_cache[level] = policy->renumber_mg_dofs (level, new_numbers);
1152 }
1153 
1154 
1155 
1156 template <int dim, int spacedim>
1157 unsigned int
1159 {
1161 }
1162 
1163 
1164 
1165 template <int dim, int spacedim>
1166 unsigned int
1168 {
1169  switch (dim)
1170  {
1171  case 1:
1172  return get_fe().dofs_per_vertex;
1173  case 2:
1174  return (3*get_fe().dofs_per_vertex +
1175  2*get_fe().dofs_per_line);
1176  case 3:
1177  // we need to take refinement of
1178  // one boundary face into
1179  // consideration here; in fact,
1180  // this function returns what
1181  // #max_coupling_between_dofs<2>
1182  // returns
1183  //
1184  // we assume here, that only four
1185  // faces meet at the boundary;
1186  // this assumption is not
1187  // justified and needs to be
1188  // fixed some time. fortunately,
1189  // omitting it for now does no
1190  // harm since the matrix will cry
1191  // foul if its requirements are
1192  // not satisfied
1193  return (19*get_fe().dofs_per_vertex +
1194  28*get_fe().dofs_per_line +
1195  8*get_fe().dofs_per_quad);
1196  default:
1197  Assert (false, ExcNotImplemented());
1199  }
1200 }
1201 
1202 
1203 
1204 template <int dim, int spacedim>
1206 {
1207  levels.clear ();
1208  faces.reset ();
1209 
1210  std::vector<types::global_dof_index> tmp;
1211  std::swap (vertex_dofs, tmp);
1212 
1213  number_cache.clear ();
1214 }
1215 
1216 
1217 
1218 template <int dim, int spacedim>
1219 template <int structdim>
1221 DoFHandler<dim, spacedim>::get_dof_index (const unsigned int obj_level,
1222  const unsigned int obj_index,
1223  const unsigned int fe_index,
1224  const unsigned int local_index) const
1225 {
1226  return internal::DoFHandler::Implementation::get_dof_index (*this, *this->mg_levels[obj_level],
1227  *this->mg_faces, obj_index,
1228  fe_index, local_index,
1229  std::integral_constant<int, structdim> ());
1230 }
1231 
1232 
1233 
1234 template <int dim, int spacedim>
1235 template <int structdim>
1236 void DoFHandler<dim, spacedim>::set_dof_index (const unsigned int obj_level,
1237  const unsigned int obj_index,
1238  const unsigned int fe_index,
1239  const unsigned int local_index,
1240  const types::global_dof_index global_index) const
1241 {
1242  internal::DoFHandler::Implementation::set_dof_index (*this,
1243  *this->mg_levels[obj_level],
1244  *this->mg_faces,
1245  obj_index,
1246  fe_index,
1247  local_index,
1248  global_index,
1249  std::integral_constant<int, structdim> ());
1250 }
1251 
1252 
1253 
1254 template <int dim, int spacedim>
1256  :
1257  coarsest_level (numbers::invalid_unsigned_int),
1258  finest_level (0)
1259 {}
1260 
1261 
1262 
1263 template <int dim, int spacedim>
1265  const unsigned int fl,
1266  const unsigned int dofs_per_vertex)
1267 {
1268  coarsest_level = cl;
1269  finest_level = fl;
1270 
1271  if (coarsest_level <= finest_level)
1272  {
1273  const unsigned int n_levels = finest_level - coarsest_level + 1;
1274  const unsigned int n_indices = n_levels * dofs_per_vertex;
1275 
1276  indices.reset (new types::global_dof_index[n_indices]);
1277  std::fill (indices.get(), indices.get()+n_indices,
1279  }
1280  else
1281  indices.reset ();
1282 }
1283 
1284 
1285 
1286 template <int dim, int spacedim>
1288 {
1289  return coarsest_level;
1290 }
1291 
1292 
1293 
1294 template <int dim, int spacedim>
1296 {
1297  return finest_level;
1298 }
1299 
1300 
1301 /*-------------- Explicit Instantiations -------------------------------*/
1302 #include "dof_handler.inst"
1303 
1304 
1305 DEAL_II_NAMESPACE_CLOSE
unsigned int max_couplings_between_boundary_dofs() const
std::vector< MGVertexDoFs > mg_vertex_dofs
Definition: dof_handler.h:1079
active_cell_iterator end_active(const unsigned int level) const
Definition: dof_handler.cc:772
static const unsigned int invalid_unsigned_int
Definition: types.h:173
active_cell_iterator begin_active(const unsigned int level=0) const
Definition: dof_handler.cc:734
Iterator points to a valid object.
level_cell_iterator begin_mg(const unsigned int level=0) const
Definition: dof_handler.cc:784
#define AssertIndexRange(index, range)
Definition: exceptions.h:1199
virtual void clear()
const Triangulation< dim, spacedim > & get_triangulation() const
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
void clear_mg_space()
ActiveSelector::active_cell_iterator active_cell_iterator
Definition: dof_handler.h:229
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
Definition: dof_handler.h:939
#define AssertThrow(cond, exc)
Definition: exceptions.h:398
unsigned int max_couplings_between_dofs() const
void initialize(const Triangulation< dim, spacedim > &tria, const FiniteElement< dim, spacedim > &fe)
Definition: dof_handler.cc:698
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
Definition: dof_handler.cc:863
static::ExceptionBase & ExcMessage(std::string arg1)
level_cell_iterator end_mg() const
Definition: dof_handler.cc:810
ActiveSelector::cell_iterator cell_iterator
Definition: dof_handler.h:257
void reserve_space()
unsigned int global_dof_index
Definition: types.h:88
std::vector< std::unique_ptr<::internal::DoFHandler::DoFLevel< dim > > > levels
Definition: dof_handler.h:1085
virtual void distribute_mg_dofs()
#define Assert(cond, exc)
Definition: exceptions.h:337
types::global_dof_index n_dofs() const
void set_dof_index(const ::DoFHandler< dh_dim, spacedim > &dof_handler, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index)
Definition: dof_objects.cc:41
unsigned int get_finest_level() const
types::global_dof_index get_dof_index(const ::DoFHandler< dh_dim, spacedim > &dof_handler, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index) const
Definition: dof_objects.h:189
unsigned int get_coarsest_level() const
static::ExceptionBase & ExcRenumberingIncomplete()
IteratorState::IteratorStates state() const
internal::DoFHandler::DoFObjects< 1 > lines
Definition: dof_faces.h:116
DoFObjects< dim > dof_object
Definition: dof_levels.h:82
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:93
cell_iterator end() const
Definition: dof_handler.cc:750
virtual ~DoFHandler()
Definition: dof_handler.cc:683
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
Definition: dof_handler.cc:852
IteratorRange< cell_iterator > cell_iterators() const
Definition: dof_handler.cc:819
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
void initialize_local_block_info()
static unsigned int max_couplings_between_dofs(const DoFHandler< 1, spacedim > &dof_handler)
Definition: dof_handler.cc:108
internal::DoFHandler::DoFObjects< 1 > lines
Definition: dof_faces.h:146
cell_iterator begin(const unsigned int level=0) const
Definition: dof_handler.cc:722
IteratorRange< level_cell_iterator > mg_cell_iterators() const
Definition: dof_handler.cc:840
types::global_dof_index n_boundary_dofs() const
Definition: dof_handler.cc:888
void clear_space()
internal::DoFHandler::DoFObjects< 2 > quads
Definition: dof_faces.h:151
IteratorRange< level_cell_iterator > mg_cell_iterators_on_level(const unsigned int level) const
Definition: dof_handler.cc:874
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)
Definition: dof_handler.cc:992
static::ExceptionBase & ExcInvalidBoundaryIndicator()
static::ExceptionBase & ExcNotImplemented()
const types::boundary_id internal_face_boundary_id
Definition: types.h:219
virtual std::size_t memory_consumption() const
Definition: dof_handler.cc:963
const types::global_dof_index invalid_dof_index
Definition: types.h:187
static void reserve_space(DoFHandler< 1, spacedim > &dof_handler)
Definition: dof_handler.cc:274
IteratorRange< active_cell_iterator > active_cell_iterators() const
Definition: dof_handler.cc:829
std::unique_ptr<::internal::DoFHandler::DoFFaces< dim > > faces
Definition: dof_handler.h:1094
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:1073
static::ExceptionBase & ExcInternalError()