Reference documentation for deal.II version 8.1.0
dof_level.h
1 // ---------------------------------------------------------------------
2 // @f$Id: dof_level.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 2003 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__hp_dof_level_h
18 #define __deal2__hp_dof_level_h
19 
20 
21 #include <deal.II/base/config.h>
23 
24 #include <vector>
25 
26 
28 
29 namespace hp
30 {
31  template <int, int> class DoFHandler;
32  template <int, int> class FECollection;
33 }
34 
35 
36 namespace internal
37 {
38  namespace hp
39  {
40  namespace DoFHandler
41  {
42  struct Implementation;
43  }
44  }
45  namespace DoFCellAccessor
46  {
47  struct Implementation;
48  }
49 }
50 
51 
52 namespace internal
53 {
54  namespace hp
55  {
103  class DoFLevel
104  {
105  private:
109  typedef unsigned int offset_type;
110 
114  typedef unsigned short int active_fe_index_type;
115 
120  typedef signed short int signed_active_fe_index_type;
121 
133  std::vector<active_fe_index_type> active_fe_indices;
134 
147  std::vector<offset_type> dof_offsets;
148 
154  std::vector<types::global_dof_index> dof_indices;
155 
159  std::vector<offset_type> cell_cache_offsets;
160 
166  std::vector<types::global_dof_index> cell_dof_indices_cache;
167 
168  public:
169 
194  void
195  set_dof_index (const unsigned int obj_index,
196  const unsigned int fe_index,
197  const unsigned int local_index,
198  const types::global_dof_index global_index);
199 
224  get_dof_index (const unsigned int obj_index,
225  const unsigned int fe_index,
226  const unsigned int local_index) const;
227 
233  unsigned int
234  active_fe_index (const unsigned int obj_index) const;
235 
242  bool
243  fe_index_is_active (const unsigned int obj_index,
244  const unsigned int fe_index) const;
245 
251  void
252  set_active_fe_index (const unsigned int obj_index,
253  const unsigned int fe_index);
254 
267  get_cell_cache_start (const unsigned int obj_index,
268  const unsigned int dofs_per_cell) const;
269 
275  std::size_t memory_consumption () const;
276 
277  private:
287  template <int dim, int spacedim>
288  void compress_data (const ::hp::FECollection<dim,spacedim> &fe_collection);
289 
299  template <int dim, int spacedim>
300  void uncompress_data (const ::hp::FECollection<dim,spacedim> &fe_collection);
301 
306  template <int, int> friend class ::hp::DoFHandler;
307  friend struct ::internal::hp::DoFHandler::Implementation;
308  friend struct ::internal::DoFCellAccessor::Implementation;
309  };
310 
311 
312  // -------------------- template functions --------------------------------
313 
314  inline
316  DoFLevel::
317  get_dof_index (const unsigned int obj_index,
318  const unsigned int fe_index,
319  const unsigned int local_index) const
320  {
321  Assert (obj_index < dof_offsets.size(),
322  ExcIndexRange (obj_index, 0, dof_offsets.size()));
323 
324  // make sure we are on an
325  // object for which DoFs have
326  // been allocated at all
327  Assert (dof_offsets[obj_index] != (offset_type)(-1),
328  ExcMessage ("You are trying to access degree of freedom "
329  "information for an object on which no such "
330  "information is available"));
331 
332  Assert (fe_index == active_fe_indices[obj_index],
333  ExcMessage ("FE index does not match that of the present cell"));
334 
335  // see if the dof_indices array has been compressed for this
336  // particular cell
337  if ((signed_active_fe_index_type)active_fe_indices[obj_index]>=0)
338  return dof_indices[dof_offsets[obj_index]+local_index];
339  else
340  return dof_indices[dof_offsets[obj_index]]+local_index;
341  }
342 
343 
344 
345  inline
346  void
347  DoFLevel::
348  set_dof_index (const unsigned int obj_index,
349  const unsigned int fe_index,
350  const unsigned int local_index,
351  const types::global_dof_index global_index)
352  {
353  Assert (obj_index < dof_offsets.size(),
354  ExcIndexRange (obj_index, 0, dof_offsets.size()));
355 
356  // make sure we are on an
357  // object for which DoFs have
358  // been allocated at all
359  Assert (dof_offsets[obj_index] != (offset_type)(-1),
360  ExcMessage ("You are trying to access degree of freedom "
361  "information for an object on which no such "
362  "information is available"));
363  Assert ((signed_active_fe_index_type)active_fe_indices[obj_index]>=0,
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;
368  }
369 
370 
371 
372  inline
373  unsigned int
374  DoFLevel::
375  active_fe_index (const unsigned int obj_index) const
376  {
377  Assert (obj_index < active_fe_indices.size(),
378  ExcIndexRange (obj_index, 0, active_fe_indices.size()));
379 
380  if (((signed_active_fe_index_type)active_fe_indices[obj_index]) >= 0)
381  return active_fe_indices[obj_index];
382  else
383  return (active_fe_index_type)~(signed_active_fe_index_type)active_fe_indices[obj_index];
384  }
385 
386 
387 
388  inline
389  bool
390  DoFLevel::
391  fe_index_is_active (const unsigned int obj_index,
392  const unsigned int fe_index) const
393  {
394  return (fe_index == active_fe_index(obj_index));
395  }
396 
397 
398  inline
399  void
400  DoFLevel::
401  set_active_fe_index (const unsigned int obj_index,
402  const unsigned int fe_index)
403  {
404  Assert (obj_index < active_fe_indices.size(),
405  ExcIndexRange (obj_index, 0, active_fe_indices.size()));
406 
407  active_fe_indices[obj_index] = fe_index;
408  }
409 
410 
411 
412  inline
414  DoFLevel::get_cell_cache_start (const unsigned int obj_index,
415  const unsigned int dofs_per_cell) const
416  {
417  Assert (obj_index < cell_cache_offsets.size(),
418  ExcInternalError());
419  Assert (cell_cache_offsets[obj_index]+dofs_per_cell
420  <=
421  cell_dof_indices_cache.size(),
422  ExcInternalError());
423 
424  return &cell_dof_indices_cache[cell_cache_offsets[obj_index]];
425  }
426  } // namespace hp
427 
428 } // namespace internal
429 
430 DEAL_II_NAMESPACE_CLOSE
431 
432 #endif
unsigned int offset_type
Definition: dof_level.h:109
::ExceptionBase & ExcMessage(std::string arg1)
std::vector< offset_type > cell_cache_offsets
Definition: dof_level.h:159
unsigned short int active_fe_index_type
Definition: dof_level.h:114
unsigned int global_dof_index
Definition: types.h:100
#define Assert(cond, exc)
Definition: exceptions.h:299
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
Definition: hp.h:103
std::vector< offset_type > dof_offsets
Definition: dof_level.h:147
std::vector< types::global_dof_index > dof_indices
Definition: dof_level.h:154
std::vector< types::global_dof_index > cell_dof_indices_cache
Definition: dof_level.h:166
signed short int signed_active_fe_index_type
Definition: dof_level.h:120
::ExceptionBase & ExcInternalError()
std::vector< active_fe_index_type > active_fe_indices
Definition: dof_level.h:133