Reference documentation for deal.II version Git 58efed1 2018-02-20 08:43:45 -0500
dof_handler.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2018 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 = std_cxx14::make_unique<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 = std_cxx14::make_unique<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 (std_cxx14::make_unique<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 = std_cxx14::make_unique<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 (std_cxx14::make_unique<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 = std_cxx14::make_unique<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,
555  const std::unique_ptr<internal::DoFHandler::DoFLevel<1> > &mg_level,
556  const std::unique_ptr<internal::DoFHandler::DoFFaces<1> > &,
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, const std::unique_ptr<internal::DoFHandler::DoFLevel<2> > &, const std::unique_ptr<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, const std::unique_ptr<internal::DoFHandler::DoFLevel<2> > &mg_level, const std::unique_ptr<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, const std::unique_ptr<internal::DoFHandler::DoFLevel<3> > &, const std::unique_ptr<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, const std::unique_ptr<internal::DoFHandler::DoFLevel<3> > &, const std::unique_ptr<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, const std::unique_ptr<internal::DoFHandler::DoFLevel<3> > &mg_level, const std::unique_ptr<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, const std::unique_ptr<internal::DoFHandler::DoFLevel<1> > &mg_level, const std::unique_ptr<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, const std::unique_ptr<internal::DoFHandler::DoFLevel<2> > &, const std::unique_ptr<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, const std::unique_ptr<internal::DoFHandler::DoFLevel<2> > &mg_level, const std::unique_ptr<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, const std::unique_ptr<internal::DoFHandler::DoFLevel<3> > &, const std::unique_ptr<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, const std::unique_ptr<internal::DoFHandler::DoFLevel<3> > &, const std::unique_ptr<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, const std::unique_ptr<internal::DoFHandler::DoFLevel<3> > &mg_level, const std::unique_ptr<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)
664  policy = std_cxx14::make_unique<internal::DoFHandler::Policy::ParallelShared<DoFHandler<dim,spacedim> >> (*this);
665  else if (dynamic_cast<const parallel::distributed::Triangulation< dim, spacedim >*>
666  (&tria)
667  == nullptr)
668  policy = std_cxx14::make_unique<internal::DoFHandler::Policy::Sequential<DoFHandler<dim,spacedim> >> (*this);
669  else
670  policy = std_cxx14::make_unique<internal::DoFHandler::Policy::ParallelDistributed<DoFHandler<dim,spacedim> >> (*this);
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  // virtual functions called in constructors and destructors never use the
687  // override in a derived class
688  // for clarity be explicit on which function is called
690 
691  // also release the policy. this needs to happen before the
692  // current object disappears because the policy objects
693  // store references to the DoFhandler object they work on
694  policy.reset ();
695 }
696 
697 
698 template <int dim, int spacedim>
699 void
702  const FiniteElement<dim,spacedim> &fe)
703 {
704  tria = &t;
705  faces = nullptr;
706  number_cache.n_global_dofs = 0;
707 
708  // decide whether we need a sequential or a parallel distributed policy
709  if (dynamic_cast<const parallel::shared::Triangulation< dim, spacedim>*> (&t) != nullptr)
710  policy = std_cxx14::make_unique<internal::DoFHandler::Policy::ParallelShared<DoFHandler<dim,spacedim> >> (*this);
711  else if (dynamic_cast<const parallel::distributed::Triangulation< dim, spacedim >*> (&t) != nullptr)
712  policy = std_cxx14::make_unique<internal::DoFHandler::Policy::ParallelDistributed<DoFHandler<dim,spacedim> >> (*this);
713  else
714  policy = std_cxx14::make_unique<internal::DoFHandler::Policy::Sequential<DoFHandler<dim,spacedim> >> (*this);
715 
716  distribute_dofs(fe);
717 }
718 
719 
720 
721 /*------------------------ Cell iterator functions ------------------------*/
722 
723 template <int dim, int spacedim>
725 DoFHandler<dim,spacedim>::begin (const unsigned int level) const
726 {
727  typename Triangulation<dim,spacedim>::cell_iterator cell = this->get_triangulation().begin(level);
728  if (cell == this->get_triangulation().end(level))
729  return end(level);
730  return cell_iterator (*cell, this);
731 }
732 
733 
734 
735 template <int dim, int spacedim>
737 DoFHandler<dim,spacedim>::begin_active (const unsigned int level) const
738 {
739  // level is checked in begin
740  cell_iterator i = begin (level);
741  if (i.state() != IteratorState::valid)
742  return i;
743  while (i->has_children())
744  if ((++i).state() != IteratorState::valid)
745  return i;
746  return i;
747 }
748 
749 
750 
751 template <int dim, int spacedim>
754 {
755  return cell_iterator (&this->get_triangulation(),
756  -1,
757  -1,
758  this);
759 }
760 
761 
762 template <int dim, int spacedim>
764 DoFHandler<dim,spacedim>::end (const unsigned int level) const
765 {
766  typename Triangulation<dim,spacedim>::cell_iterator cell = this->get_triangulation().end(level);
767  if (cell.state() != IteratorState::valid)
768  return end();
769  return cell_iterator (*cell, this);
770 }
771 
772 
773 template <int dim, int spacedim>
775 DoFHandler<dim, spacedim>::end_active (const unsigned int level) const
776 {
777  typename Triangulation<dim,spacedim>::cell_iterator cell = this->get_triangulation().end_active(level);
778  if (cell.state() != IteratorState::valid)
779  return active_cell_iterator(end());
780  return active_cell_iterator (*cell, this);
781 }
782 
783 
784 
785 template <int dim, int spacedim>
786 typename DoFHandler<dim, spacedim>::level_cell_iterator
787 DoFHandler<dim, spacedim>::begin_mg (const unsigned int level) const
788 {
789  // Assert(this->has_level_dofs(), ExcMessage("You can only iterate over mg "
790  // "levels if mg dofs got distributed."));
791  typename Triangulation<dim,spacedim>::cell_iterator cell = this->get_triangulation().begin(level);
792  if (cell == this->get_triangulation().end(level))
793  return end_mg(level);
794  return level_cell_iterator (*cell, this);
795 }
796 
797 
798 template <int dim, int spacedim>
799 typename DoFHandler<dim, spacedim>::level_cell_iterator
800 DoFHandler<dim, spacedim>::end_mg (const unsigned int level) const
801 {
802  // Assert(this->has_level_dofs(), ExcMessage("You can only iterate over mg "
803  // "levels if mg dofs got distributed."));
804  typename Triangulation<dim,spacedim>::cell_iterator cell = this->get_triangulation().end(level);
805  if (cell.state() != IteratorState::valid)
806  return end();
807  return level_cell_iterator (*cell, this);
808 }
809 
810 
811 template <int dim, int spacedim>
812 typename DoFHandler<dim, spacedim>::level_cell_iterator
814 {
815  return level_cell_iterator (&this->get_triangulation(), -1, -1, this);
816 }
817 
818 
819 
820 template <int dim, int spacedim>
823 {
824  return
826  (begin(), end());
827 }
828 
829 
830 template <int dim, int spacedim>
833 {
834  return
836  (begin_active(), end());
837 }
838 
839 
840 
841 template <int dim, int spacedim>
844 {
845  return
847  (begin_mg(), end_mg());
848 }
849 
850 
851 
852 
853 template <int dim, int spacedim>
856 {
857  return
859  (begin(level), end(level));
860 }
861 
862 
863 
864 template <int dim, int spacedim>
867 {
868  return
870  (begin_active(level), end_active(level));
871 }
872 
873 
874 
875 template <int dim, int spacedim>
878 {
879  return
881  (begin_mg(level), end_mg(level));
882 }
883 
884 
885 
886 //---------------------------------------------------------------------------
887 
888 
889 
890 template <int dim, int spacedim>
892 {
893  std::set<int> boundary_dofs;
894 
895  const unsigned int dofs_per_face = get_fe().dofs_per_face;
896  std::vector<types::global_dof_index> dofs_on_face(dofs_per_face);
897 
898  // loop over all faces of all cells
899  // and see whether they are at a
900  // boundary. note (i) that we visit
901  // interior faces twice (which we
902  // don't care about) but exterior
903  // faces only once as is
904  // appropriate, and (ii) that we
905  // need not take special care of
906  // single lines (using
907  // @p{cell->has_boundary_lines}),
908  // since we do not support
909  // boundaries of dimension dim-2,
910  // and so every boundary line is
911  // also part of a boundary face.
912  active_cell_iterator cell = begin_active (),
913  endc = end();
914  for (; cell!=endc; ++cell)
915  for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
916  if (cell->at_boundary(f))
917  {
918  cell->face(f)->get_dof_indices (dofs_on_face);
919  for (unsigned int i=0; i<dofs_per_face; ++i)
920  boundary_dofs.insert(dofs_on_face[i]);
921  }
922 
923  return boundary_dofs.size();
924 }
925 
926 
927 
928 template <int dim, int spacedim>
930 DoFHandler<dim,spacedim>::n_boundary_dofs (const std::set<types::boundary_id> &boundary_ids) const
931 {
932  Assert (boundary_ids.find (numbers::internal_face_boundary_id) == boundary_ids.end(),
934 
935  std::set<types::global_dof_index> boundary_dofs;
936 
937  const unsigned int dofs_per_face = get_fe().dofs_per_face;
938  std::vector<types::global_dof_index> dofs_on_face(dofs_per_face);
939 
940  // same as in the previous
941  // function, but with a different
942  // check for the boundary indicator
943  active_cell_iterator cell = begin_active (),
944  endc = end();
945  for (; cell!=endc; ++cell)
946  for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
947  if (cell->at_boundary(f)
948  &&
949  (std::find (boundary_ids.begin(),
950  boundary_ids.end(),
951  cell->face(f)->boundary_id()) !=
952  boundary_ids.end()))
953  {
954  cell->face(f)->get_dof_indices (dofs_on_face);
955  for (unsigned int i=0; i<dofs_per_face; ++i)
956  boundary_dofs.insert(dofs_on_face[i]);
957  }
958 
959  return boundary_dofs.size();
960 }
961 
962 
963 
964 template <int dim, int spacedim>
965 std::size_t
967 {
968  std::size_t mem = (MemoryConsumption::memory_consumption (tria) +
969  MemoryConsumption::memory_consumption (fe_collection) +
970  MemoryConsumption::memory_consumption (block_info_object) +
974  sizeof (number_cache) +
975  MemoryConsumption::memory_consumption (mg_number_cache) +
977  for (unsigned int i=0; i<levels.size(); ++i)
978  mem += MemoryConsumption::memory_consumption (*levels[i]);
979 
980  for (unsigned int level = 0; level < mg_levels.size (); ++level)
981  mem += mg_levels[level]->memory_consumption ();
982 
983  if (mg_faces != nullptr)
984  mem += MemoryConsumption::memory_consumption (*mg_faces);
985 
986  for (unsigned int i = 0; i < mg_vertex_dofs.size (); ++i)
987  mem += sizeof (MGVertexDoFs) + (1 + mg_vertex_dofs[i].get_finest_level () - mg_vertex_dofs[i].get_coarsest_level ()) * sizeof (types::global_dof_index);
988 
989  return mem;
990 }
991 
992 
993 
994 template <int dim, int spacedim>
996 {
997  Assert(tria!=nullptr,
998  ExcMessage("You need to set the Triangulation in the DoFHandler using initialize() or "
999  "in the constructor before you can distribute DoFs."));
1000  Assert (tria->n_levels() > 0,
1001  ExcMessage("The Triangulation you are using is empty!"));
1002 
1003  // Only recreate the FECollection if we don't already store
1004  // the exact same FiniteElement object.
1005  if (fe_collection == nullptr || &((*fe_collection)[0]) != &ff)
1006  fe_collection = std_cxx14::make_unique<hp::FECollection<dim, spacedim>>(ff);
1007 
1008  // delete all levels and set them
1009  // up newly. note that we still
1010  // have to allocate space for all
1011  // degrees of freedom on this mesh
1012  // (including ghost and cells that
1013  // are entirely stored on different
1014  // processors), though we may not
1015  // assign numbers to some of them
1016  // (i.e. they will remain at
1017  // invalid_dof_index). We need to
1018  // allocate the space because we
1019  // will want to be able to query
1020  // the dof_indices on each cell,
1021  // and simply be told that we don't
1022  // know them on some cell (i.e. get
1023  // back invalid_dof_index)
1024  clear_space ();
1026 
1027  // hand things off to the policy
1028  number_cache = policy->distribute_dofs ();
1029 
1030  // initialize the block info object
1031  // only if this is a sequential
1032  // triangulation. it doesn't work
1033  // correctly yet if it is parallel
1034  if (dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim>*>(&*tria) == nullptr)
1035  block_info_object.initialize(*this, false, true);
1036 }
1037 
1038 
1039 
1040 template <int dim, int spacedim>
1042 {
1043  this->distribute_mg_dofs();
1044 }
1045 
1046 
1047 
1048 template <int dim, int spacedim>
1050 {
1051  Assert(levels.size()>0, ExcMessage("Distribute active DoFs using distribute_dofs() before calling distribute_mg_dofs()."));
1052 
1055  ExcMessage("The mesh smoothing requirement 'limit_level_difference_at_vertices' has to be set for using multigrid!"));
1056 
1057  clear_mg_space();
1058 
1059  internal::DoFHandler::Implementation::reserve_space_mg (*this);
1060  mg_number_cache = policy->distribute_mg_dofs ();
1061 
1062  // initialize the block info object
1063  // only if this is a sequential
1064  // triangulation. it doesn't work
1065  // correctly yet if it is parallel
1066  if (dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim>*>(&*tria) == nullptr)
1067  block_info_object.initialize (*this, true, false);
1068 }
1069 
1070 
1071 
1072 template <int dim, int spacedim>
1074 {
1075  mg_levels.clear ();
1076  mg_faces.reset ();
1077 
1078  std::vector<MGVertexDoFs> tmp;
1079 
1080  std::swap (mg_vertex_dofs, tmp);
1081 
1082  mg_number_cache.clear();
1083 }
1084 
1085 
1086 template <int dim, int spacedim>
1088 {
1089  block_info_object.initialize_local(*this);
1090 }
1091 
1092 
1093 
1094 template <int dim, int spacedim>
1096 {
1097  // release lock to old fe
1098  fe_collection.reset();
1099 
1100  // release memory
1101  clear_space ();
1102  clear_mg_space ();
1103 }
1104 
1105 
1106 
1107 template <int dim, int spacedim>
1108 void
1109 DoFHandler<dim,spacedim>::renumber_dofs (const std::vector<types::global_dof_index> &new_numbers)
1110 {
1111  Assert(levels.size()>0, ExcMessage("You need to distribute DoFs before you can renumber them."));
1112 
1113 #ifdef DEBUG
1114  if (dynamic_cast<const parallel::shared::Triangulation< dim, spacedim>*> (&*tria) != nullptr)
1115  {
1116  Assert(new_numbers.size() == n_dofs() || new_numbers.size() == n_locally_owned_dofs(),
1117  ExcMessage("Incorrect size of the input array."));
1118  }
1119  else if (dynamic_cast<const parallel::distributed::Triangulation< dim, spacedim >*> (&*tria) != nullptr)
1120  {
1121  AssertDimension (new_numbers.size(), n_locally_owned_dofs());
1122  }
1123  else
1124  {
1125  AssertDimension (new_numbers.size(), n_dofs());
1126  }
1127 
1128  // assert that the new indices are
1129  // consecutively numbered if we are
1130  // working on a single
1131  // processor. this doesn't need to
1132  // hold in the case of a parallel
1133  // mesh since we map the interval
1134  // [0...n_dofs()) into itself but
1135  // only globally, not on each
1136  // processor
1137  if (n_locally_owned_dofs() == n_dofs())
1138  {
1139  std::vector<types::global_dof_index> tmp(new_numbers);
1140  std::sort (tmp.begin(), tmp.end());
1141  std::vector<types::global_dof_index>::const_iterator p = tmp.begin();
1143  for (; p!=tmp.end(); ++p, ++i)
1144  Assert (*p == i, ExcNewNumbersNotConsecutive(i));
1145  }
1146  else
1147  for (types::global_dof_index i=0; i<new_numbers.size(); ++i)
1148  Assert (new_numbers[i] < n_dofs(),
1149  ExcMessage ("New DoF index is not less than the total number of dofs."));
1150 #endif
1151 
1152  number_cache = policy->renumber_dofs (new_numbers);
1153 }
1154 
1155 
1156 template <int dim, int spacedim>
1157 void
1159  const std::vector<types::global_dof_index> &new_numbers)
1160 {
1161  Assert(mg_levels.size()>0 && levels.size()>0,
1162  ExcMessage("You need to distribute active and level DoFs before you can renumber level DoFs."));
1163  AssertIndexRange(level, get_triangulation().n_global_levels());
1164  AssertDimension (new_numbers.size(), locally_owned_mg_dofs(level).n_elements());
1165 
1166 #ifdef DEBUG
1167  // assert that the new indices are consecutively numbered if we are working
1168  // on a single processor. this doesn't need to hold in the case of a
1169  // parallel mesh since we map the interval [0...n_dofs(level)) into itself
1170  // but only globally, not on each processor
1171  if (n_locally_owned_dofs() == n_dofs())
1172  {
1173  std::vector<types::global_dof_index> tmp(new_numbers);
1174  std::sort (tmp.begin(), tmp.end());
1175  std::vector<types::global_dof_index>::const_iterator p = tmp.begin();
1177  for (; p!=tmp.end(); ++p, ++i)
1178  Assert (*p == i, ExcNewNumbersNotConsecutive(i));
1179  }
1180  else
1181  for (types::global_dof_index i=0; i<new_numbers.size(); ++i)
1182  Assert (new_numbers[i] < n_dofs(level),
1183  ExcMessage ("New DoF index is not less than the total number of dofs."));
1184 #endif
1185 
1186  mg_number_cache[level] = policy->renumber_mg_dofs (level, new_numbers);
1187 }
1188 
1189 
1190 
1191 template <int dim, int spacedim>
1192 unsigned int
1194 {
1196 }
1197 
1198 
1199 
1200 template <int dim, int spacedim>
1201 unsigned int
1203 {
1204  switch (dim)
1205  {
1206  case 1:
1207  return get_fe().dofs_per_vertex;
1208  case 2:
1209  return (3*get_fe().dofs_per_vertex +
1210  2*get_fe().dofs_per_line);
1211  case 3:
1212  // we need to take refinement of
1213  // one boundary face into
1214  // consideration here; in fact,
1215  // this function returns what
1216  // #max_coupling_between_dofs<2>
1217  // returns
1218  //
1219  // we assume here, that only four
1220  // faces meet at the boundary;
1221  // this assumption is not
1222  // justified and needs to be
1223  // fixed some time. fortunately,
1224  // omitting it for now does no
1225  // harm since the matrix will cry
1226  // foul if its requirements are
1227  // not satisfied
1228  return (19*get_fe().dofs_per_vertex +
1229  28*get_fe().dofs_per_line +
1230  8*get_fe().dofs_per_quad);
1231  default:
1232  Assert (false, ExcNotImplemented());
1234  }
1235 }
1236 
1237 
1238 
1239 template <int dim, int spacedim>
1241 {
1242  levels.clear ();
1243  faces.reset ();
1244 
1245  std::vector<types::global_dof_index> tmp;
1246  std::swap (vertex_dofs, tmp);
1247 
1248  number_cache.clear ();
1249 }
1250 
1251 
1252 
1253 template <int dim, int spacedim>
1254 template <int structdim>
1256 DoFHandler<dim, spacedim>::get_dof_index (const unsigned int obj_level,
1257  const unsigned int obj_index,
1258  const unsigned int fe_index,
1259  const unsigned int local_index) const
1260 {
1261  return internal::DoFHandler::Implementation::get_dof_index (*this, this->mg_levels[obj_level],
1262  this->mg_faces, obj_index,
1263  fe_index, local_index,
1264  std::integral_constant<int, structdim> ());
1265 }
1266 
1267 
1268 
1269 template <int dim, int spacedim>
1270 template <int structdim>
1271 void DoFHandler<dim, spacedim>::set_dof_index (const unsigned int obj_level,
1272  const unsigned int obj_index,
1273  const unsigned int fe_index,
1274  const unsigned int local_index,
1275  const types::global_dof_index global_index) const
1276 {
1277  internal::DoFHandler::Implementation::set_dof_index (*this,
1278  this->mg_levels[obj_level],
1279  this->mg_faces,
1280  obj_index,
1281  fe_index,
1282  local_index,
1283  global_index,
1284  std::integral_constant<int, structdim> ());
1285 }
1286 
1287 
1288 
1289 template <int dim, int spacedim>
1291  :
1292  coarsest_level (numbers::invalid_unsigned_int),
1293  finest_level (0)
1294 {}
1295 
1296 
1297 
1298 template <int dim, int spacedim>
1300  const unsigned int fl,
1301  const unsigned int dofs_per_vertex)
1302 {
1303  coarsest_level = cl;
1304  finest_level = fl;
1305 
1306  if (coarsest_level <= finest_level)
1307  {
1308  const unsigned int n_levels = finest_level - coarsest_level + 1;
1309  const unsigned int n_indices = n_levels * dofs_per_vertex;
1310 
1311  indices = std_cxx14::make_unique<types::global_dof_index[]> (n_indices);
1312  std::fill (indices.get(), indices.get()+n_indices,
1314  }
1315  else
1316  indices.reset ();
1317 }
1318 
1319 
1320 
1321 template <int dim, int spacedim>
1323 {
1324  return coarsest_level;
1325 }
1326 
1327 
1328 
1329 template <int dim, int spacedim>
1331 {
1332  return finest_level;
1333 }
1334 
1335 
1336 /*-------------- Explicit Instantiations -------------------------------*/
1337 #include "dof_handler.inst"
1338 
1339 
1340 DEAL_II_NAMESPACE_CLOSE
unsigned int max_couplings_between_boundary_dofs() const
std::vector< MGVertexDoFs > mg_vertex_dofs
Definition: dof_handler.h:1168
active_cell_iterator end_active(const unsigned int level) const
Definition: dof_handler.cc:775
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:737
Iterator points to a valid object.
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1202
level_cell_iterator begin_mg(const unsigned int level=0) const
Definition: dof_handler.cc:787
#define AssertIndexRange(index, range)
Definition: exceptions.h:1237
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:1033
#define AssertThrow(cond, exc)
Definition: exceptions.h:410
unsigned int max_couplings_between_dofs() const
void initialize(const Triangulation< dim, spacedim > &tria, const FiniteElement< dim, spacedim > &fe)
Definition: dof_handler.cc:701
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
Definition: dof_handler.cc:866
static::ExceptionBase & ExcMessage(std::string arg1)
level_cell_iterator end_mg() const
Definition: dof_handler.cc:813
ActiveSelector::cell_iterator cell_iterator
Definition: dof_handler.h:257
unsigned int global_dof_index
Definition: types.h:88
std::vector< std::unique_ptr<::internal::DoFHandler::DoFLevel< dim > > > levels
Definition: dof_handler.h:1174
virtual void distribute_mg_dofs()
#define Assert(cond, exc)
Definition: exceptions.h:349
types::global_dof_index n_dofs() const
unsigned int get_finest_level() const
unsigned int get_coarsest_level() const
IteratorState::IteratorStates state() const
cell_iterator end() const
Definition: dof_handler.cc:753
virtual ~DoFHandler()
Definition: dof_handler.cc:683
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
Definition: dof_handler.cc:855
IteratorRange< cell_iterator > cell_iterators() const
Definition: dof_handler.cc:822
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
cell_iterator begin(const unsigned int level=0) const
Definition: dof_handler.cc:725
IteratorRange< level_cell_iterator > mg_cell_iterators() const
Definition: dof_handler.cc:843
types::global_dof_index n_boundary_dofs() const
Definition: dof_handler.cc:891
void clear_space()
IteratorRange< level_cell_iterator > mg_cell_iterators_on_level(const unsigned int level) const
Definition: dof_handler.cc:877
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:995
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:966
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:832
std::unique_ptr<::internal::DoFHandler::DoFFaces< dim > > faces
Definition: dof_handler.h:1183
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:1162
static::ExceptionBase & ExcInternalError()