16 #ifndef dealii_hp_dof_level_h 17 #define dealii_hp_dof_level_h 20 #include <deal.II/base/config.h> 22 #include <deal.II/base/exceptions.h> 27 DEAL_II_NAMESPACE_OPEN
42 namespace DoFHandlerImplementation
44 struct Implementation;
47 namespace DoFCellAccessorImplementation
49 struct Implementation;
201 set_dof_index(
const unsigned int obj_index,
202 const unsigned int fe_index,
203 const unsigned int local_index,
218 get_dof_index(
const unsigned int obj_index,
219 const unsigned int fe_index,
220 const unsigned int local_index)
const;
226 active_fe_index(
const unsigned int obj_index)
const;
233 fe_index_is_active(
const unsigned int obj_index,
234 const unsigned int fe_index)
const;
240 set_active_fe_index(
const unsigned int obj_index,
241 const unsigned int fe_index);
255 get_cell_cache_start(
const unsigned int obj_index,
256 const unsigned int dofs_per_cell)
const;
263 memory_consumption()
const;
269 template <
class Archive>
271 serialize(Archive &ar,
const unsigned int version);
282 template <
int dim,
int spacedim>
285 const ::hp::FECollection<dim, spacedim> &fe_collection);
295 template <
int dim,
int spacedim>
298 const ::hp::FECollection<dim, spacedim> &fe_collection);
319 normalize_active_fe_indices();
327 friend class ::hp::DoFHandler;
328 friend struct ::internal::hp::DoFHandlerImplementation::
330 friend struct ::internal::DoFCellAccessorImplementation::
347 DoFLevel::get_toggled_compression_state(
359 DoFLevel::get_dof_index(
const unsigned int obj_index,
360 const unsigned int fe_index,
361 const unsigned int local_index)
const 364 Assert(obj_index < dof_offsets.size(),
370 ExcMessage(
"You are trying to access degree of freedom " 371 "information for an object on which no such " 372 "information is available"));
375 (is_compressed_entry(active_fe_indices[obj_index]) ==
false ?
376 active_fe_indices[obj_index] :
377 get_toggled_compression_state(active_fe_indices[obj_index])),
378 ExcMessage(
"FE index does not match that of the present cell"));
382 if (is_compressed_entry(active_fe_indices[obj_index]) ==
false)
383 return dof_indices[dof_offsets[obj_index] + local_index];
385 return dof_indices[dof_offsets[obj_index]] + local_index;
391 DoFLevel::set_dof_index(
const unsigned int obj_index,
392 const unsigned int fe_index,
393 const unsigned int local_index,
397 Assert(obj_index < dof_offsets.size(),
404 ExcMessage(
"You are trying to access degree of freedom " 405 "information for an object on which no such " 406 "information is available"));
408 is_compressed_entry(active_fe_indices[obj_index]) ==
false,
410 "This function can no longer be called after compressing the dof_indices array"));
411 Assert(fe_index == active_fe_indices[obj_index],
412 ExcMessage(
"FE index does not match that of the present cell"));
413 dof_indices[dof_offsets[obj_index] + local_index] = global_index;
419 DoFLevel::active_fe_index(
const unsigned int obj_index)
const 421 Assert(obj_index < active_fe_indices.size(),
424 if (is_compressed_entry(active_fe_indices[obj_index]) ==
false)
425 return active_fe_indices[obj_index];
427 return get_toggled_compression_state(active_fe_indices[obj_index]);
433 DoFLevel::fe_index_is_active(
const unsigned int obj_index,
434 const unsigned int fe_index)
const 436 return (fe_index == active_fe_index(obj_index));
442 DoFLevel::set_active_fe_index(
const unsigned int obj_index,
443 const unsigned int fe_index)
445 Assert(obj_index < active_fe_indices.size(),
453 Assert(is_compressed_entry(fe_index) ==
false,
455 "You are using an active_fe_index that is larger than an " 456 "internal limitation for these objects. Try to work with " 457 "hp::FECollection objects that have a more modest size."));
459 active_fe_indices[obj_index] = fe_index;
465 DoFLevel::get_cell_cache_start(
const unsigned int obj_index,
466 const unsigned int dofs_per_cell)
const 470 (obj_index < cell_cache_offsets.size()) &&
471 (cell_cache_offsets[obj_index] + dofs_per_cell <=
472 cell_dof_indices_cache.size()),
474 "You are trying to access an element of the cache that stores " 475 "the indices of all degrees of freedom that live on one cell. " 476 "However, this element does not exist. Did you forget to call " 477 "DoFHandler::distribute_dofs(), or did you forget to call it " 478 "again after changing the active_fe_index of one of the cells?"));
480 return &cell_dof_indices_cache[cell_cache_offsets[obj_index]];
485 template <
class Archive>
487 DoFLevel::serialize(Archive &ar,
const unsigned int)
489 ar & this->active_fe_indices;
490 ar & this->cell_cache_offsets;
491 ar & this->cell_dof_indices_cache;
492 ar & this->dof_indices;
493 ar & this->dof_offsets;
499 DEAL_II_NAMESPACE_CLOSE
std::vector< offset_type > cell_cache_offsets
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
unsigned long long int global_dof_index
unsigned short int active_fe_index_type
static::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
std::vector< offset_type > dof_offsets
std::vector< types::global_dof_index > dof_indices
signed short int signed_active_fe_index_type
std::vector< types::global_dof_index > cell_dof_indices_cache
std::vector< active_fe_index_type > active_fe_indices