16 #ifndef dealii__hp_dof_level_h 17 #define dealii__hp_dof_level_h 20 #include <deal.II/base/config.h> 21 #include <deal.II/base/exceptions.h> 26 DEAL_II_NAMESPACE_OPEN
31 template <
int,
int>
class FECollection;
41 struct Implementation;
46 struct Implementation;
129 is_compressed_entry (
const active_fe_index_type active_fe_index);
139 get_toggled_compression_state (
const active_fe_index_type active_fe_index);
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;
262 std::size_t memory_consumption ()
const;
268 template <
class Archive>
269 void serialize(Archive &ar,
270 const unsigned int version);
281 template <
int dim,
int spacedim>
282 void compress_data (const ::hp::FECollection<dim,spacedim> &fe_collection);
292 template <
int dim,
int spacedim>
293 void uncompress_data (const ::hp::FECollection<dim,spacedim> &fe_collection);
313 void normalize_active_fe_indices ();
320 template <
int,
int>
friend class ::hp::DoFHandler;
321 friend struct ::internal::hp::DoFHandler::Implementation;
322 friend struct ::internal::DoFCellAccessor::Implementation;
351 get_dof_index (
const unsigned int obj_index,
352 const unsigned int fe_index,
353 const unsigned int local_index)
const 356 Assert (obj_index < dof_offsets.size(),
362 ExcMessage (
"You are trying to access degree of freedom " 363 "information for an object on which no such " 364 "information is available"));
366 Assert (fe_index == (is_compressed_entry(active_fe_indices[obj_index]) ==
false ?
367 active_fe_indices[obj_index] :
368 get_toggled_compression_state(active_fe_indices[obj_index])),
369 ExcMessage (
"FE index does not match that of the present cell"));
373 if (is_compressed_entry(active_fe_indices[obj_index]) ==
false)
374 return dof_indices[dof_offsets[obj_index]+local_index];
376 return dof_indices[dof_offsets[obj_index]]+local_index;
384 set_dof_index (
const unsigned int obj_index,
385 const unsigned int fe_index,
386 const unsigned int local_index,
390 Assert (obj_index < dof_offsets.size(),
397 ExcMessage (
"You are trying to access degree of freedom " 398 "information for an object on which no such " 399 "information is available"));
400 Assert (is_compressed_entry(active_fe_indices[obj_index]) ==
false,
401 ExcMessage (
"This function can no longer be called after compressing the dof_indices array"));
402 Assert (fe_index == active_fe_indices[obj_index],
403 ExcMessage (
"FE index does not match that of the present cell"));
404 dof_indices[dof_offsets[obj_index]+local_index] = global_index;
412 active_fe_index (
const unsigned int obj_index)
const 414 Assert (obj_index < active_fe_indices.size(),
417 if (is_compressed_entry(active_fe_indices[obj_index]) ==
false)
418 return active_fe_indices[obj_index];
420 return get_toggled_compression_state(active_fe_indices[obj_index]);
428 fe_index_is_active (
const unsigned int obj_index,
429 const unsigned int fe_index)
const 431 return (fe_index == active_fe_index(obj_index));
439 set_active_fe_index (
const unsigned int obj_index,
440 const unsigned int fe_index)
442 Assert (obj_index < active_fe_indices.size(),
450 Assert (is_compressed_entry (fe_index) ==
false,
451 ExcMessage (
"You are using an active_fe_index that is larger than an " 452 "internal limitation for these objects. Try to work with " 453 "hp::FECollection objects that have a more modest size."));
455 active_fe_indices[obj_index] = fe_index;
462 DoFLevel::get_cell_cache_start (
const unsigned int obj_index,
463 const unsigned int dofs_per_cell)
const 466 Assert ((obj_index < cell_cache_offsets.size())
468 (cell_cache_offsets[obj_index]+dofs_per_cell
470 cell_dof_indices_cache.size()),
471 ExcMessage(
"You are trying to access an element of the cache that stores " 472 "the indices of all degrees of freedom that live on one cell. " 473 "However, this element does not exist. Did you forget to call " 474 "DoFHandler::distribute_dofs(), or did you forget to call it " 475 "again after changing the active_fe_index of one of the cells?"));
477 return &cell_dof_indices_cache[cell_cache_offsets[obj_index]];
482 template <
class Archive>
485 DoFLevel::serialize(Archive &ar,
488 ar &this->active_fe_indices;
489 ar &this->cell_cache_offsets;
490 ar &this->cell_dof_indices_cache;
491 ar &this->dof_indices;
492 ar &this->dof_offsets;
498 DEAL_II_NAMESPACE_CLOSE
std::vector< offset_type > cell_cache_offsets
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
static::ExceptionBase & ExcMessage(std::string arg1)
unsigned short int active_fe_index_type
unsigned int global_dof_index
#define Assert(cond, exc)
std::vector< offset_type > dof_offsets
std::vector< types::global_dof_index > dof_indices
std::vector< types::global_dof_index > cell_dof_indices_cache
signed short int signed_active_fe_index_type
std::vector< active_fe_index_type > active_fe_indices