17 #ifndef __deal2__hp_dof_level_h 18 #define __deal2__hp_dof_level_h 21 #include <deal.II/base/config.h> 32 template <
int,
int>
class FECollection;
42 struct Implementation;
47 struct Implementation;
195 set_dof_index (
const unsigned int obj_index,
196 const unsigned int fe_index,
197 const unsigned int local_index,
224 get_dof_index (
const unsigned int obj_index,
225 const unsigned int fe_index,
226 const unsigned int local_index)
const;
234 active_fe_index (
const unsigned int obj_index)
const;
243 fe_index_is_active (
const unsigned int obj_index,
244 const unsigned int fe_index)
const;
252 set_active_fe_index (
const unsigned int obj_index,
253 const unsigned int fe_index);
267 get_cell_cache_start (
const unsigned int obj_index,
268 const unsigned int dofs_per_cell)
const;
275 std::size_t memory_consumption ()
const;
287 template <
int dim,
int spacedim>
288 void compress_data (const ::hp::FECollection<dim,spacedim> &fe_collection);
299 template <
int dim,
int spacedim>
300 void uncompress_data (const ::hp::FECollection<dim,spacedim> &fe_collection);
306 template <
int,
int>
friend class ::hp::DoFHandler;
307 friend struct ::internal::hp::DoFHandler::Implementation;
308 friend struct ::internal::DoFCellAccessor::Implementation;
317 get_dof_index (
const unsigned int obj_index,
318 const unsigned int fe_index,
319 const unsigned int local_index)
const 321 Assert (obj_index < dof_offsets.size(),
328 ExcMessage (
"You are trying to access degree of freedom " 329 "information for an object on which no such " 330 "information is available"));
332 Assert (fe_index == active_fe_indices[obj_index],
333 ExcMessage (
"FE index does not match that of the present cell"));
338 return dof_indices[dof_offsets[obj_index]+local_index];
340 return dof_indices[dof_offsets[obj_index]]+local_index;
348 set_dof_index (
const unsigned int obj_index,
349 const unsigned int fe_index,
350 const unsigned int local_index,
353 Assert (obj_index < dof_offsets.size(),
360 ExcMessage (
"You are trying to access degree of freedom " 361 "information for an object on which no such " 362 "information is available"));
364 ExcMessage (
"This function can no longer be called after compressing the dof_indices array"));
365 Assert (fe_index == active_fe_indices[obj_index],
366 ExcMessage (
"FE index does not match that of the present cell"));
367 dof_indices[dof_offsets[obj_index]+local_index] = global_index;
375 active_fe_index (
const unsigned int obj_index)
const 377 Assert (obj_index < active_fe_indices.size(),
381 return active_fe_indices[obj_index];
391 fe_index_is_active (
const unsigned int obj_index,
392 const unsigned int fe_index)
const 394 return (fe_index == active_fe_index(obj_index));
401 set_active_fe_index (
const unsigned int obj_index,
402 const unsigned int fe_index)
404 Assert (obj_index < active_fe_indices.size(),
407 active_fe_indices[obj_index] = fe_index;
414 DoFLevel::get_cell_cache_start (
const unsigned int obj_index,
415 const unsigned int dofs_per_cell)
const 417 Assert (obj_index < cell_cache_offsets.size(),
419 Assert (cell_cache_offsets[obj_index]+dofs_per_cell
421 cell_dof_indices_cache.size(),
424 return &cell_dof_indices_cache[cell_cache_offsets[obj_index]];
430 DEAL_II_NAMESPACE_CLOSE
::ExceptionBase & ExcMessage(std::string arg1)
std::vector< offset_type > cell_cache_offsets
unsigned short int active_fe_index_type
unsigned int global_dof_index
#define Assert(cond, exc)
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
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
::ExceptionBase & ExcInternalError()
std::vector< active_fe_index_type > active_fe_indices