17 #ifndef __deal2__hp_dof_faces_h 18 #define __deal2__hp_dof_faces_h 21 #include <deal.II/base/config.h> 22 #include <deal.II/hp/fe_collection.h> 30 template <
int dim,
int spacedim>
103 template <
int structdim>
104 class DoFIndicesOnFacesOrEdges
123 std::vector<types::global_dof_index>
dofs;
149 template <
int dim,
int spacedim>
151 set_dof_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
152 const unsigned int obj_index,
153 const unsigned int fe_index,
154 const unsigned int local_index,
156 const unsigned int obj_level);
181 template <
int dim,
int spacedim>
183 get_dof_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
184 const unsigned int obj_index,
185 const unsigned int fe_index,
186 const unsigned int local_index,
187 const unsigned int obj_level)
const;
212 template <
int dim,
int spacedim>
214 n_active_fe_indices (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
215 const unsigned int obj_index)
const;
222 template <
int dim,
int spacedim>
224 nth_active_fe_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
225 const unsigned int obj_level,
226 const unsigned int obj_index,
227 const unsigned int n)
const;
235 template <
int dim,
int spacedim>
237 fe_index_is_active (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
238 const unsigned int obj_index,
239 const unsigned int fe_index,
240 const unsigned int obj_level)
const;
247 std::size_t memory_consumption ()
const;
296 std::size_t memory_consumption ()
const;
319 std::size_t memory_consumption ()
const;
348 std::size_t memory_consumption ()
const;
354 template <
int structdim>
355 template <
int dim,
int spacedim>
360 const unsigned int obj_index,
361 const unsigned int fe_index,
362 const unsigned int local_index,
363 const unsigned int obj_level)
const 366 ExcMessage (
"You need to specify a FE index when working " 367 "with hp DoFHandlers"));
368 Assert (&dof_handler != 0,
369 ExcMessage (
"No DoFHandler is specified for this iterator"));
370 Assert (&dof_handler.get_fe() != 0,
371 ExcMessage (
"No finite element collection is associated with " 373 Assert (fe_index < dof_handler.get_fe().size(),
376 dof_handler.get_fe()[fe_index].template n_dofs_per_object<structdim>(),
378 dof_handler.get_fe()[fe_index]
379 .template n_dofs_per_object<structdim>()));
380 Assert (obj_index < dof_offsets.size(),
387 ExcMessage (
"You are trying to access degree of freedom " 388 "information for an object on which no such " 389 "information is available"));
391 Assert (structdim<dim,
ExcMessage (
"This object can not be used for cells."));
404 if (*pointer == fe_index)
405 return *(pointer + 1 + local_index);
408 dof_handler.get_fe()[*pointer]
409 .template n_dofs_per_object<structdim>() + 1);
415 template <
int structdim>
416 template <
int dim,
int spacedim>
421 const unsigned int obj_index,
422 const unsigned int fe_index,
423 const unsigned int local_index,
425 const unsigned int obj_level)
428 ExcMessage (
"You need to specify a FE index when working " 429 "with hp DoFHandlers"));
430 Assert (&dof_handler != 0,
431 ExcMessage (
"No DoFHandler is specified for this iterator"));
432 Assert (&dof_handler.get_fe() != 0,
433 ExcMessage (
"No finite element collection is associated with " 435 Assert (fe_index < dof_handler.get_fe().size(),
438 dof_handler.get_fe()[fe_index].template n_dofs_per_object<structdim>(),
440 dof_handler.get_fe()[fe_index]
441 .template n_dofs_per_object<structdim>()));
442 Assert (obj_index < dof_offsets.size(),
449 ExcMessage (
"You are trying to access degree of freedom " 450 "information for an object on which no such " 451 "information is available"));
453 Assert (structdim<dim,
ExcMessage (
"This object can not be used for cells."));
466 if (*pointer == fe_index)
468 *(pointer + 1 + local_index) = global_index;
472 pointer += dof_handler.get_fe()[*pointer]
473 .template n_dofs_per_object<structdim>() + 1;
479 template <
int structdim>
480 template <
int dim,
int spacedim>
485 const unsigned int obj_index)
const 487 Assert (&dof_handler != 0,
488 ExcMessage (
"No DoFHandler is specified for this iterator"));
489 Assert (&dof_handler.get_fe() != 0,
490 ExcMessage (
"No finite element collection is associated with " 492 Assert (obj_index < dof_offsets.size(),
501 Assert (structdim<dim,
ExcMessage (
"This object can not be used for cells."));
508 const unsigned int starting_offset = dof_offsets[obj_index];
510 unsigned int counter = 0;
519 pointer += dof_handler.get_fe()[*pointer]
520 .template n_dofs_per_object<structdim>() + 1;
527 template <
int structdim>
528 template <
int dim,
int spacedim>
533 const unsigned int obj_level,
534 const unsigned int obj_index,
535 const unsigned int n)
const 537 Assert (&dof_handler != 0,
538 ExcMessage (
"No DoFHandler is specified for this iterator"));
539 Assert (&dof_handler.get_fe() != 0,
540 ExcMessage (
"No finite element collection is associated with " 542 Assert (obj_index < dof_offsets.size(),
549 ExcMessage (
"You are trying to access degree of freedom " 550 "information for an object on which no such " 551 "information is available"));
553 Assert (structdim<dim,
ExcMessage (
"This object can not be used for cells."));
555 Assert (n < n_active_fe_indices(dof_handler, obj_index),
557 n_active_fe_indices(dof_handler, obj_index)));
564 const unsigned int starting_offset = dof_offsets[obj_index];
566 unsigned int counter = 0;
572 const unsigned int fe_index = *pointer;
574 Assert (fe_index < dof_handler.get_fe().size(),
581 pointer += dof_handler.get_fe()[fe_index]
582 .template n_dofs_per_object<structdim>() + 1;
588 template <
int structdim>
589 template <
int dim,
int spacedim>
594 const unsigned int obj_index,
595 const unsigned int fe_index,
596 const unsigned int obj_level)
const 598 Assert (&dof_handler != 0,
599 ExcMessage (
"No DoFHandler is specified for this iterator"));
600 Assert (&dof_handler.get_fe() != 0,
601 ExcMessage (
"No finite element collection is associated with " 603 Assert (obj_index < dof_offsets.size(),
604 ExcIndexRange (obj_index, 0, static_cast<unsigned int>(dof_offsets.size())));
606 ExcMessage (
"You need to specify a FE index when working " 607 "with hp DoFHandlers"));
608 Assert (fe_index < dof_handler.get_fe().size(),
615 ExcMessage (
"You are trying to access degree of freedom " 616 "information for an object on which no such " 617 "information is available"));
619 Assert (structdim<dim,
ExcMessage (
"This object can not be used for cells."));
633 else if (*pointer == fe_index)
637 dof_handler.get_fe()[*pointer]
638 .template n_dofs_per_object<structdim>()+1);
647 DEAL_II_NAMESPACE_CLOSE
static const unsigned int invalid_unsigned_int
::ExceptionBase & ExcMessage(std::string arg1)
internal::hp::DoFIndicesOnFacesOrEdges< 1 > lines
unsigned int global_dof_index
#define Assert(cond, exc)
std::vector< types::global_dof_index > dofs
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
std::vector< unsigned int > dof_offsets
internal::hp::DoFIndicesOnFacesOrEdges< 2 > quads
const types::global_dof_index invalid_dof_index
::ExceptionBase & ExcInternalError()
internal::hp::DoFIndicesOnFacesOrEdges< 1 > lines