16 #include <deal.II/base/geometry_info.h> 17 #include <deal.II/base/memory_consumption.h> 18 #include <deal.II/base/std_cxx14/memory.h> 20 #include <deal.II/distributed/shared_tria.h> 21 #include <deal.II/distributed/tria.h> 23 #include <deal.II/dofs/dof_accessor.h> 24 #include <deal.II/dofs/dof_faces.h> 25 #include <deal.II/dofs/dof_handler.h> 26 #include <deal.II/dofs/dof_handler_policy.h> 27 #include <deal.II/dofs/dof_levels.h> 29 #include <deal.II/fe/fe.h> 31 #include <deal.II/grid/tria.h> 32 #include <deal.II/grid/tria_accessor.h> 33 #include <deal.II/grid/tria_iterator.h> 34 #include <deal.II/grid/tria_levels.h> 39 DEAL_II_NAMESPACE_OPEN
42 template <
int dim,
int spacedim>
45 template <
int dim,
int spacedim>
48 template <
int dim,
int spacedim>
51 template <
int dim,
int spacedim>
59 template <
int dim,
int spacedim>
63 return &::numbers::invalid_dof_index;
71 template <
int dim,
int spacedim>
73 policy_to_string(const ::internal::DoFHandlerImplementation::Policy::
74 PolicyBase<dim, spacedim> &policy)
76 std::string policy_name;
77 if (
dynamic_cast<const typename ::internal::DoFHandlerImplementation::
78 Policy::Sequential<::DoFHandler<dim, spacedim>
> *>(
81 const typename ::internal::DoFHandlerImplementation::Policy::
82 Sequential<::hp::DoFHandler<dim, spacedim>
> *>(&policy))
83 policy_name =
"Policy::Sequential<";
84 else if (
dynamic_cast< 85 const typename ::internal::DoFHandlerImplementation::
86 Policy::ParallelDistributed<::DoFHandler<dim, spacedim>
> 90 Policy::ParallelDistributed<
92 policy_name =
"Policy::ParallelDistributed<";
93 else if (
dynamic_cast< 94 const typename ::internal::DoFHandlerImplementation::
95 Policy::ParallelShared<::DoFHandler<dim, spacedim>
> *>(
97 dynamic_cast<const typename ::
internal::
98 DoFHandlerImplementation::Policy::ParallelShared<
101 policy_name =
"Policy::ParallelShared<";
110 namespace DoFHandlerImplementation
128 template <
int spacedim>
132 return std::min(static_cast<types::global_dof_index>(
133 3 * dof_handler.
get_fe().dofs_per_vertex +
134 2 * dof_handler.
get_fe().dofs_per_line),
140 template <
int spacedim>
166 switch (dof_handler.
tria->max_adjacent_cells())
169 max_couplings = 19 * dof_handler.
get_fe().dofs_per_vertex +
170 28 * dof_handler.
get_fe().dofs_per_line +
171 8 * dof_handler.
get_fe().dofs_per_quad;
174 max_couplings = 21 * dof_handler.
get_fe().dofs_per_vertex +
175 31 * dof_handler.
get_fe().dofs_per_line +
176 9 * dof_handler.
get_fe().dofs_per_quad;
179 max_couplings = 28 * dof_handler.
get_fe().dofs_per_vertex +
180 42 * dof_handler.
get_fe().dofs_per_line +
181 12 * dof_handler.
get_fe().dofs_per_quad;
184 max_couplings = 30 * dof_handler.
get_fe().dofs_per_vertex +
185 45 * dof_handler.
get_fe().dofs_per_line +
186 13 * dof_handler.
get_fe().dofs_per_quad;
189 max_couplings = 37 * dof_handler.
get_fe().dofs_per_vertex +
190 56 * dof_handler.
get_fe().dofs_per_line +
191 16 * dof_handler.
get_fe().dofs_per_quad;
199 max_couplings = 39 * dof_handler.
get_fe().dofs_per_vertex +
200 59 * dof_handler.
get_fe().dofs_per_line +
201 17 * dof_handler.
get_fe().dofs_per_quad;
204 max_couplings = 46 * dof_handler.
get_fe().dofs_per_vertex +
205 70 * dof_handler.
get_fe().dofs_per_line +
206 20 * dof_handler.
get_fe().dofs_per_quad;
209 max_couplings = 48 * dof_handler.
get_fe().dofs_per_vertex +
210 73 * dof_handler.
get_fe().dofs_per_line +
211 21 * dof_handler.
get_fe().dofs_per_quad;
214 max_couplings = 55 * dof_handler.
get_fe().dofs_per_vertex +
215 84 * dof_handler.
get_fe().dofs_per_line +
216 24 * dof_handler.
get_fe().dofs_per_quad;
219 max_couplings = 57 * dof_handler.
get_fe().dofs_per_vertex +
220 87 * dof_handler.
get_fe().dofs_per_line +
221 25 * dof_handler.
get_fe().dofs_per_quad;
224 max_couplings = 63 * dof_handler.
get_fe().dofs_per_vertex +
225 98 * dof_handler.
get_fe().dofs_per_line +
226 28 * dof_handler.
get_fe().dofs_per_quad;
229 max_couplings = 65 * dof_handler.
get_fe().dofs_per_vertex +
230 103 * dof_handler.
get_fe().dofs_per_line +
231 29 * dof_handler.
get_fe().dofs_per_quad;
234 max_couplings = 72 * dof_handler.
get_fe().dofs_per_vertex +
235 114 * dof_handler.
get_fe().dofs_per_line +
236 32 * dof_handler.
get_fe().dofs_per_quad;
243 return std::min(max_couplings, dof_handler.
n_dofs());
247 template <
int spacedim>
265 const unsigned int max_adjacent_cells =
266 dof_handler.
tria->max_adjacent_cells();
269 if (max_adjacent_cells <= 8)
270 max_couplings = 7 * 7 * 7 * dof_handler.
get_fe().dofs_per_vertex +
271 7 * 6 * 7 * 3 * dof_handler.
get_fe().dofs_per_line +
272 9 * 4 * 7 * 3 * dof_handler.
get_fe().dofs_per_quad +
273 27 * dof_handler.
get_fe().dofs_per_hex;
280 return std::min(max_couplings, dof_handler.
n_dofs());
293 template <
int spacedim>
297 dof_handler.
get_fe().dofs_per_vertex,
300 for (
unsigned int i = 0; i < dof_handler.
tria->n_levels(); ++i)
302 dof_handler.
levels.emplace_back(
305 dof_handler.
levels.back()->dof_object.dofs.resize(
306 dof_handler.
tria->n_raw_cells(i) *
307 dof_handler.
get_fe().dofs_per_line,
310 dof_handler.
levels.back()->cell_dof_indices_cache.resize(
311 dof_handler.
tria->n_raw_cells(i) *
312 dof_handler.
get_fe().dofs_per_cell,
318 template <
int spacedim>
322 dof_handler.
get_fe().dofs_per_vertex,
325 for (
unsigned int i = 0; i < dof_handler.
tria->n_levels(); ++i)
327 dof_handler.
levels.emplace_back(
330 dof_handler.
levels.back()->dof_object.dofs.resize(
331 dof_handler.
tria->n_raw_cells(i) *
332 dof_handler.
get_fe().dofs_per_quad,
335 dof_handler.
levels.back()->cell_dof_indices_cache.resize(
336 dof_handler.
tria->n_raw_cells(i) *
337 dof_handler.
get_fe().dofs_per_cell,
341 dof_handler.
faces = std_cxx14::make_unique<
344 if (dof_handler.
tria->n_cells() > 0)
346 dof_handler.
faces->lines.dofs.resize(
347 dof_handler.
tria->n_raw_lines() *
348 dof_handler.
get_fe().dofs_per_line,
354 template <
int spacedim>
358 dof_handler.
get_fe().dofs_per_vertex,
361 for (
unsigned int i = 0; i < dof_handler.
tria->n_levels(); ++i)
363 dof_handler.
levels.emplace_back(
366 dof_handler.
levels.back()->dof_object.dofs.resize(
367 dof_handler.
tria->n_raw_cells(i) *
368 dof_handler.
get_fe().dofs_per_hex,
371 dof_handler.
levels.back()->cell_dof_indices_cache.resize(
372 dof_handler.
tria->n_raw_cells(i) *
373 dof_handler.
get_fe().dofs_per_cell,
376 dof_handler.
faces = std_cxx14::make_unique<
380 if (dof_handler.
tria->n_cells() > 0)
382 dof_handler.
faces->lines.dofs.resize(
383 dof_handler.
tria->n_raw_lines() *
384 dof_handler.
get_fe().dofs_per_line,
386 dof_handler.
faces->quads.dofs.resize(
387 dof_handler.
tria->n_raw_quads() *
388 dof_handler.
get_fe().dofs_per_quad,
393 template <
int spacedim>
400 const ::Triangulation<1, spacedim> &tria =
402 const unsigned int &dofs_per_line = dof_handler.
get_fe().dofs_per_line;
403 const unsigned int &n_levels = tria.n_levels();
405 for (
unsigned int i = 0; i < n_levels; ++i)
407 dof_handler.mg_levels.emplace_back(
409 dof_handler.mg_levels.back()->dof_object.dofs =
410 std::vector<types::global_dof_index>(tria.n_raw_lines(i) *
415 const unsigned int &n_vertices = tria.n_vertices();
419 std::vector<unsigned int> max_level(n_vertices, 0);
420 std::vector<unsigned int> min_level(n_vertices, n_levels);
422 for (typename ::Triangulation<1, spacedim>::cell_iterator cell =
427 const unsigned int level = cell->level();
429 for (
unsigned int vertex = 0;
430 vertex < GeometryInfo<1>::vertices_per_cell;
433 const unsigned int vertex_index = cell->vertex_index(vertex);
435 if (min_level[vertex_index] > level)
436 min_level[vertex_index] = level;
438 if (max_level[vertex_index] < level)
439 max_level[vertex_index] = level;
443 for (
unsigned int vertex = 0; vertex < n_vertices; ++vertex)
444 if (tria.vertex_used(vertex))
447 Assert(max_level[vertex] >= min_level[vertex],
452 dof_handler.
get_fe().dofs_per_vertex);
463 template <
int spacedim>
470 const ::FiniteElement<2, spacedim> &fe = dof_handler.
get_fe();
471 const ::Triangulation<2, spacedim> &tria =
473 const unsigned int &n_levels = tria.n_levels();
475 for (
unsigned int i = 0; i < n_levels; ++i)
477 dof_handler.mg_levels.emplace_back(
478 std_cxx14::make_unique<
480 dof_handler.mg_levels.back()->dof_object.dofs =
481 std::vector<types::global_dof_index>(tria.n_raw_quads(i) *
486 dof_handler.mg_faces = std_cxx14::make_unique<
488 dof_handler.mg_faces->lines.dofs = std::vector<types::global_dof_index>(
491 const unsigned int &n_vertices = tria.n_vertices();
495 std::vector<unsigned int> max_level(n_vertices, 0);
496 std::vector<unsigned int> min_level(n_vertices, n_levels);
498 for (typename ::Triangulation<2, spacedim>::cell_iterator cell =
503 const unsigned int level = cell->level();
505 for (
unsigned int vertex = 0;
506 vertex < GeometryInfo<2>::vertices_per_cell;
509 const unsigned int vertex_index = cell->vertex_index(vertex);
511 if (min_level[vertex_index] > level)
512 min_level[vertex_index] = level;
514 if (max_level[vertex_index] < level)
515 max_level[vertex_index] = level;
519 for (
unsigned int vertex = 0; vertex < n_vertices; ++vertex)
520 if (tria.vertex_used(vertex))
523 Assert(max_level[vertex] >= min_level[vertex],
538 template <
int spacedim>
545 const ::FiniteElement<3, spacedim> &fe = dof_handler.
get_fe();
546 const ::Triangulation<3, spacedim> &tria =
548 const unsigned int &n_levels = tria.n_levels();
550 for (
unsigned int i = 0; i < n_levels; ++i)
552 dof_handler.mg_levels.emplace_back(
553 std_cxx14::make_unique<
555 dof_handler.mg_levels.back()->dof_object.dofs =
556 std::vector<types::global_dof_index>(tria.n_raw_hexs(i) *
561 dof_handler.mg_faces = std_cxx14::make_unique<
563 dof_handler.mg_faces->lines.dofs = std::vector<types::global_dof_index>(
565 dof_handler.mg_faces->quads.dofs = std::vector<types::global_dof_index>(
568 const unsigned int &n_vertices = tria.n_vertices();
572 std::vector<unsigned int> max_level(n_vertices, 0);
573 std::vector<unsigned int> min_level(n_vertices, n_levels);
575 for (typename ::Triangulation<3, spacedim>::cell_iterator cell =
580 const unsigned int level = cell->level();
582 for (
unsigned int vertex = 0;
583 vertex < GeometryInfo<3>::vertices_per_cell;
586 const unsigned int vertex_index = cell->vertex_index(vertex);
588 if (min_level[vertex_index] > level)
589 min_level[vertex_index] = level;
591 if (max_level[vertex_index] < level)
592 max_level[vertex_index] = level;
596 for (
unsigned int vertex = 0; vertex < n_vertices; ++vertex)
597 if (tria.vertex_used(vertex))
600 Assert(max_level[vertex] >= min_level[vertex],
617 template <
int spacedim>
625 const unsigned int obj_index,
626 const unsigned int fe_index,
627 const unsigned int local_index,
628 const std::integral_constant<int, 1>)
630 return mg_level->dof_object.get_dof_index(dof_handler,
636 template <
int spacedim>
644 const unsigned int obj_index,
645 const unsigned int fe_index,
646 const unsigned int local_index,
647 const std::integral_constant<int, 1>)
649 return mg_faces->lines.get_dof_index(dof_handler,
655 template <
int spacedim>
663 const unsigned int obj_index,
664 const unsigned int fe_index,
665 const unsigned int local_index,
666 const std::integral_constant<int, 2>)
668 return mg_level->dof_object.get_dof_index(dof_handler,
674 template <
int spacedim>
682 const unsigned int obj_index,
683 const unsigned int fe_index,
684 const unsigned int local_index,
685 const std::integral_constant<int, 1>)
687 return mg_faces->lines.get_dof_index(dof_handler,
693 template <
int spacedim>
701 const unsigned int obj_index,
702 const unsigned int fe_index,
703 const unsigned int local_index,
704 const std::integral_constant<int, 2>)
706 return mg_faces->quads.get_dof_index(dof_handler,
712 template <
int spacedim>
720 const unsigned int obj_index,
721 const unsigned int fe_index,
722 const unsigned int local_index,
723 const std::integral_constant<int, 3>)
725 return mg_level->dof_object.get_dof_index(dof_handler,
731 template <
int spacedim>
739 const unsigned int obj_index,
740 const unsigned int fe_index,
741 const unsigned int local_index,
743 const std::integral_constant<int, 1>)
745 mg_level->dof_object.set_dof_index(
746 dof_handler, obj_index, fe_index, local_index, global_index);
749 template <
int spacedim>
757 const unsigned int obj_index,
758 const unsigned int fe_index,
759 const unsigned int local_index,
761 const std::integral_constant<int, 1>)
763 mg_faces->lines.set_dof_index(
764 dof_handler, obj_index, fe_index, local_index, global_index);
767 template <
int spacedim>
775 const unsigned int obj_index,
776 const unsigned int fe_index,
777 const unsigned int local_index,
779 const std::integral_constant<int, 2>)
781 mg_level->dof_object.set_dof_index(
782 dof_handler, obj_index, fe_index, local_index, global_index);
785 template <
int spacedim>
793 const unsigned int obj_index,
794 const unsigned int fe_index,
795 const unsigned int local_index,
797 const std::integral_constant<int, 1>)
799 mg_faces->lines.set_dof_index(
800 dof_handler, obj_index, fe_index, local_index, global_index);
803 template <
int spacedim>
811 const unsigned int obj_index,
812 const unsigned int fe_index,
813 const unsigned int local_index,
815 const std::integral_constant<int, 2>)
817 mg_faces->quads.set_dof_index(
818 dof_handler, obj_index, fe_index, local_index, global_index);
821 template <
int spacedim>
829 const unsigned int obj_index,
830 const unsigned int fe_index,
831 const unsigned int local_index,
833 const std::integral_constant<int, 3>)
835 mg_level->dof_object.set_dof_index(
836 dof_handler, obj_index, fe_index, local_index, global_index);
844 template <
int dim,
int spacedim>
846 : tria(&tria, typeid(*this).name())
857 else if (
dynamic_cast< 871 template <
int dim,
int spacedim>
873 :
tria(nullptr, typeid(*this).name())
877 template <
int dim,
int spacedim>
893 template <
int dim,
int spacedim>
909 else if (
dynamic_cast< 928 template <
int dim,
int spacedim>
941 template <
int dim,
int spacedim>
949 while (i->has_children())
957 template <
int dim,
int spacedim>
965 template <
int dim,
int spacedim>
977 template <
int dim,
int spacedim>
990 template <
int dim,
int spacedim>
991 typename DoFHandler<dim, spacedim>::level_cell_iterator
1000 return level_cell_iterator(*cell,
this);
1004 template <
int dim,
int spacedim>
1005 typename DoFHandler<dim, spacedim>::level_cell_iterator
1014 return level_cell_iterator(*cell,
this);
1018 template <
int dim,
int spacedim>
1019 typename DoFHandler<dim, spacedim>::level_cell_iterator
1027 template <
int dim,
int spacedim>
1036 template <
int dim,
int spacedim>
1047 template <
int dim,
int spacedim>
1057 template <
int dim,
int spacedim>
1060 const unsigned int level)
const 1068 template <
int dim,
int spacedim>
1071 const unsigned int level)
const 1080 template <
int dim,
int spacedim>
1083 const unsigned int level)
const 1095 template <
int dim,
int spacedim>
1099 std::set<int> boundary_dofs;
1101 const unsigned int dofs_per_face =
get_fe().dofs_per_face;
1102 std::vector<types::global_dof_index> dofs_on_face(dofs_per_face);
1119 for (; cell != endc; ++cell)
1120 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
1121 if (cell->at_boundary(f))
1123 cell->face(f)->get_dof_indices(dofs_on_face);
1124 for (
unsigned int i = 0; i < dofs_per_face; ++i)
1125 boundary_dofs.insert(dofs_on_face[i]);
1128 return boundary_dofs.size();
1133 template <
int dim,
int spacedim>
1136 const std::set<types::boundary_id> &boundary_ids)
const 1142 std::set<types::global_dof_index> boundary_dofs;
1144 const unsigned int dofs_per_face =
get_fe().dofs_per_face;
1145 std::vector<types::global_dof_index> dofs_on_face(dofs_per_face);
1151 for (; cell != endc; ++cell)
1152 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
1153 if (cell->at_boundary(f) &&
1154 (std::find(boundary_ids.begin(),
1156 cell->face(f)->boundary_id()) != boundary_ids.end()))
1158 cell->face(f)->get_dof_indices(dofs_on_face);
1159 for (
unsigned int i = 0; i < dofs_per_face; ++i)
1160 boundary_dofs.insert(dofs_on_face[i]);
1163 return boundary_dofs.size();
1168 template <
int dim,
int spacedim>
1181 for (
unsigned int i = 0; i <
levels.size(); ++i)
1184 for (
unsigned int level = 0; level < mg_levels.size(); ++level)
1187 if (mg_faces !=
nullptr)
1200 template <
int dim,
int spacedim>
1208 "You need to set the Triangulation in the DoFHandler using initialize() or " 1209 "in the constructor before you can distribute DoFs."));
1211 ExcMessage(
"The Triangulation you are using is empty!"));
1251 template <
int dim,
int spacedim>
1261 template <
int dim,
int spacedim>
1268 "Distribute active DoFs using distribute_dofs() before calling distribute_mg_dofs()."));
1271 ((
tria->get_mesh_smoothing() &
1275 "The mesh smoothing requirement 'limit_level_difference_at_vertices' has to be set for using multigrid!"));
1279 internal::DoFHandlerImplementation::Implementation::reserve_space_mg(*
this);
1293 template <
int dim,
int spacedim>
1300 std::vector<MGVertexDoFs> tmp;
1308 template <
int dim,
int spacedim>
1317 template <
int dim,
int spacedim>
1328 template <
int dim,
int spacedim>
1331 const std::vector<types::global_dof_index> &new_numbers)
1335 "You need to distribute DoFs before you can renumber them."));
1343 ExcMessage(
"Incorrect size of the input array."));
1345 else if (
dynamic_cast< 1367 std::vector<types::global_dof_index> tmp(new_numbers);
1368 std::sort(tmp.begin(), tmp.end());
1369 std::vector<types::global_dof_index>::const_iterator p = tmp.begin();
1371 for (; p != tmp.end(); ++p, ++i)
1378 "New DoF index is not less than the total number of dofs."));
1385 template <
int dim,
int spacedim>
1388 const unsigned int level,
1389 const std::vector<types::global_dof_index> &new_numbers)
1392 mg_levels.size() > 0 &&
levels.size() > 0,
1394 "You need to distribute active and level DoFs before you can renumber level DoFs."));
1406 std::vector<types::global_dof_index> tmp(new_numbers);
1407 std::sort(tmp.begin(), tmp.end());
1408 std::vector<types::global_dof_index>::const_iterator p = tmp.begin();
1410 for (; p != tmp.end(); ++p, ++i)
1417 "New DoF index is not less than the total number of dofs."));
1420 mg_number_cache[level] = policy->renumber_mg_dofs(level, new_numbers);
1425 template <
int dim,
int spacedim>
1435 template <
int dim,
int spacedim>
1442 return get_fe().dofs_per_vertex;
1444 return (3 *
get_fe().dofs_per_vertex + 2 *
get_fe().dofs_per_line);
1462 return (19 *
get_fe().dofs_per_vertex + 28 *
get_fe().dofs_per_line +
1463 8 *
get_fe().dofs_per_quad);
1472 template <
int dim,
int spacedim>
1479 std::vector<types::global_dof_index> tmp;
1487 template <
int dim,
int spacedim>
1488 template <
int structdim>
1491 const unsigned int obj_index,
1492 const unsigned int fe_index,
1493 const unsigned int local_index)
const 1495 return internal::DoFHandlerImplementation::Implementation::get_dof_index(
1497 this->mg_levels[obj_level],
1502 std::integral_constant<int, structdim>());
1507 template <
int dim,
int spacedim>
1508 template <
int structdim>
1511 const unsigned int obj_level,
1512 const unsigned int obj_index,
1513 const unsigned int fe_index,
1514 const unsigned int local_index,
1517 internal::DoFHandlerImplementation::Implementation::set_dof_index(
1519 this->mg_levels[obj_level],
1525 std::integral_constant<int, structdim>());
1530 template <
int dim,
int spacedim>
1532 : coarsest_level(
numbers::invalid_unsigned_int)
1538 template <
int dim,
int spacedim>
1541 const unsigned int cl,
1542 const unsigned int fl,
1543 const unsigned int dofs_per_vertex)
1551 const unsigned int n_indices = n_levels * dofs_per_vertex;
1553 indices = std_cxx14::make_unique<types::global_dof_index[]>(n_indices);
1564 template <
int dim,
int spacedim>
1573 template <
int dim,
int spacedim>
1582 #include "dof_handler.inst" 1585 DEAL_II_NAMESPACE_CLOSE
unsigned int max_couplings_between_boundary_dofs() const
std::vector< MGVertexDoFs > mg_vertex_dofs
active_cell_iterator end_active(const unsigned int level) const
static const unsigned int invalid_unsigned_int
active_cell_iterator begin_active(const unsigned int level=0) const
#define AssertDimension(dim1, dim2)
level_cell_iterator begin_mg(const unsigned int level=0) const
::internal::DoFHandlerImplementation::NumberCache number_cache
#define AssertIndexRange(index, range)
size_type n_elements() const
static unsigned int max_couplings_between_dofs(const DoFHandler< 1, spacedim > &dof_handler)
#define AssertThrow(cond, exc)
const Triangulation< dim, spacedim > & get_triangulation() const
unsigned int max_couplings_between_dofs() const
void initialize(const Triangulation< dim, spacedim > &tria, const FiniteElement< dim, spacedim > &fe)
unsigned long long int global_dof_index
static::ExceptionBase & ExcInvalidBoundaryIndicator()
BlockInfo block_info_object
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
static::ExceptionBase & ExcMessage(std::string arg1)
level_cell_iterator end_mg() const
unsigned int finest_level
std::vector<::internal::DoFHandlerImplementation::NumberCache > mg_number_cache
virtual void distribute_mg_dofs()
#define Assert(cond, exc)
hp::FECollection< dim, spacedim > fe_collection
types::global_dof_index n_dofs() const
unsigned int get_finest_level() const
unsigned int n_locally_owned_dofs() const
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
unsigned int get_coarsest_level() const
IteratorState::IteratorStates state() const
std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces< dim > > faces
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
cell_iterator end() const
std::vector< std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel< dim > > > levels
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
void initialize(const DoFHandler< dim, spacedim > &, bool levels_only=false, bool active_only=false)
Fill the object with values describing block structure of the DoFHandler.
IteratorRange< cell_iterator > cell_iterators() const
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
void swap(Vector< Number > &u, Vector< Number > &v)
std::unique_ptr< types::global_dof_index[]> indices
void initialize_local_block_info()
static void reserve_space(DoFHandler< 1, spacedim > &dof_handler)
cell_iterator begin(const unsigned int level=0) const
void initialize_local(const DoFHandler< dim, spacedim > &)
Initialize block structure on cells and compute renumbering between cell dofs and block cell dofs...
IteratorRange< level_cell_iterator > mg_cell_iterators() const
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
types::global_dof_index n_boundary_dofs() const
unsigned int coarsest_level
IteratorRange< level_cell_iterator > mg_cell_iterators_on_level(const unsigned int level) const
void init(const unsigned int coarsest_level, const unsigned int finest_level, const unsigned int dofs_per_vertex)
virtual void distribute_dofs(const FiniteElement< dim, spacedim > &fe)
static::ExceptionBase & ExcNotImplemented()
Iterator points to a valid object.
static::ExceptionBase & ExcNewNumbersNotConsecutive(types::global_dof_index arg1)
typename ActiveSelector::cell_iterator cell_iterator
types::global_dof_index n_global_dofs
const types::boundary_id internal_face_boundary_id
typename ActiveSelector::active_cell_iterator active_cell_iterator
virtual std::size_t memory_consumption() const
const types::global_dof_index invalid_dof_index
virtual ~DoFHandler() override
IteratorRange< active_cell_iterator > active_cell_iterators() const
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
std::vector< types::global_dof_index > vertex_dofs
static::ExceptionBase & ExcInternalError()