Reference documentation for deal.II version Git 7026f387cc 2019-10-15 14:19:01 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
dof_level.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2003 - 2019 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_hp_dof_level_h
17 #define dealii_hp_dof_level_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/exceptions.h>
23 
24 #include <vector>
25 
26 
27 DEAL_II_NAMESPACE_OPEN
28 
29 // Forward declarations
30 #ifndef DOXYGEN
31 namespace hp
32 {
33  template <int, int>
34  class DoFHandler;
35  template <int, int>
36  class FECollection;
37 } // namespace hp
38 
39 
40 namespace internal
41 {
42  namespace hp
43  {
44  namespace DoFHandlerImplementation
45  {
46  struct Implementation;
47  }
48  } // namespace hp
49  namespace DoFCellAccessorImplementation
50  {
51  struct Implementation;
52  }
53 } // namespace internal
54 #endif
55 
56 
57 namespace internal
58 {
59  namespace hp
60  {
109  class DoFLevel
110  {
111  private:
115  using offset_type = unsigned int;
116 
120  using active_fe_index_type = unsigned short int;
121 
126  using signed_active_fe_index_type = signed short int;
127 
133 
139  static bool
140  is_compressed_entry(const active_fe_index_type active_fe_index);
141 
148  static active_fe_index_type
149  get_toggled_compression_state(const active_fe_index_type active_fe_index);
150 
161  std::vector<active_fe_index_type> active_fe_indices;
162 
175  std::vector<active_fe_index_type> future_fe_indices;
176 
189  std::vector<offset_type> dof_offsets;
190 
196  std::vector<types::global_dof_index> dof_indices;
197 
201  std::vector<offset_type> cell_cache_offsets;
202 
208  std::vector<types::global_dof_index> cell_dof_indices_cache;
209 
210  public:
223  void
224  set_dof_index(const unsigned int obj_index,
225  const unsigned int fe_index,
226  const unsigned int local_index,
227  const types::global_dof_index global_index);
228 
241  get_dof_index(const unsigned int obj_index,
242  const unsigned int fe_index,
243  const unsigned int local_index) const;
244 
248  unsigned int
249  active_fe_index(const unsigned int obj_index) const;
250 
255  bool
256  fe_index_is_active(const unsigned int obj_index,
257  const unsigned int fe_index) const;
258 
262  void
263  set_active_fe_index(const unsigned int obj_index,
264  const unsigned int fe_index);
265 
270  unsigned int
271  future_fe_index(const unsigned int obj_index) const;
272 
276  void
277  set_future_fe_index(const unsigned int obj_index,
278  const unsigned int fe_index);
279 
283  bool
284  future_fe_index_set(const unsigned int obj_index) const;
285 
289  void
290  clear_future_fe_index(const unsigned int obj_index);
291 
304  get_cell_cache_start(const unsigned int obj_index,
305  const unsigned int dofs_per_cell) const;
306 
311  std::size_t
312  memory_consumption() const;
313 
318  template <class Archive>
319  void
320  serialize(Archive &ar, const unsigned int version);
321 
322  private:
331  template <int dim, int spacedim>
332  void
333  compress_data(
334  const ::hp::FECollection<dim, spacedim> &fe_collection);
335 
344  template <int dim, int spacedim>
345  void
346  uncompress_data(
347  const ::hp::FECollection<dim, spacedim> &fe_collection);
348 
349 
367  void
368  normalize_active_fe_indices();
369 
370 
371  // Make hp::DoFHandler and its auxiliary class a friend since it is the
372  // class that needs to create these data structures.
373  template <int, int>
374  friend class ::hp::DoFHandler;
375  friend struct ::internal::hp::DoFHandlerImplementation::
376  Implementation;
377  friend struct ::internal::DoFCellAccessorImplementation::
378  Implementation;
379  };
380 
381 
382  // -------------------- template functions --------------------------------
383 
384 
385  inline bool
386  DoFLevel::is_compressed_entry(const active_fe_index_type active_fe_index)
387  {
388  return (static_cast<signed_active_fe_index_type>(active_fe_index) < 0);
389  }
390 
391 
392 
394  DoFLevel::get_toggled_compression_state(
395  const active_fe_index_type active_fe_index)
396  {
397  // convert the active_fe_index into a signed type, flip all
398  // bits, and get the unsigned representation back
399  return static_cast<active_fe_index_type>(
400  ~(static_cast<signed_active_fe_index_type>(active_fe_index)));
401  }
402 
403 
404 
406  DoFLevel::get_dof_index(const unsigned int obj_index,
407  const unsigned int fe_index,
408  const unsigned int local_index) const
409  {
410  (void)fe_index;
411  Assert(obj_index < dof_offsets.size(),
412  ExcIndexRange(obj_index, 0, dof_offsets.size()));
413 
414  // make sure we are on an object for which DoFs have been
415  // allocated at all
416  Assert(dof_offsets[obj_index] != static_cast<offset_type>(-1),
417  ExcMessage("You are trying to access degree of freedom "
418  "information for an object on which no such "
419  "information is available"));
420 
421  Assert(fe_index ==
422  (is_compressed_entry(active_fe_indices[obj_index]) == false ?
423  active_fe_indices[obj_index] :
424  get_toggled_compression_state(active_fe_indices[obj_index])),
425  ExcMessage("FE index does not match that of the present cell"));
426 
427  // see if the dof_indices array has been compressed for this
428  // particular cell
429  if (is_compressed_entry(active_fe_indices[obj_index]) == false)
430  return dof_indices[dof_offsets[obj_index] + local_index];
431  else
432  return dof_indices[dof_offsets[obj_index]] + local_index;
433  }
434 
435 
436 
437  inline void
438  DoFLevel::set_dof_index(const unsigned int obj_index,
439  const unsigned int fe_index,
440  const unsigned int local_index,
441  const types::global_dof_index global_index)
442  {
443  (void)fe_index;
444  Assert(obj_index < dof_offsets.size(),
445  ExcIndexRange(obj_index, 0, dof_offsets.size()));
446 
447  // make sure we are on an
448  // object for which DoFs have
449  // been allocated at all
450  Assert(dof_offsets[obj_index] != static_cast<offset_type>(-1),
451  ExcMessage("You are trying to access degree of freedom "
452  "information for an object on which no such "
453  "information is available"));
454  Assert(
455  is_compressed_entry(active_fe_indices[obj_index]) == false,
456  ExcMessage(
457  "This function can no longer be called after compressing the dof_indices array"));
458  Assert(fe_index == active_fe_indices[obj_index],
459  ExcMessage("FE index does not match that of the present cell"));
460  dof_indices[dof_offsets[obj_index] + local_index] = global_index;
461  }
462 
463 
464 
465  inline unsigned int
466  DoFLevel::active_fe_index(const unsigned int obj_index) const
467  {
468  Assert(obj_index < active_fe_indices.size(),
469  ExcIndexRange(obj_index, 0, active_fe_indices.size()));
470 
471  if (is_compressed_entry(active_fe_indices[obj_index]) == false)
472  return active_fe_indices[obj_index];
473  else
474  return get_toggled_compression_state(active_fe_indices[obj_index]);
475  }
476 
477 
478 
479  inline bool
480  DoFLevel::fe_index_is_active(const unsigned int obj_index,
481  const unsigned int fe_index) const
482  {
483  return (fe_index == active_fe_index(obj_index));
484  }
485 
486 
487 
488  inline void
489  DoFLevel::set_active_fe_index(const unsigned int obj_index,
490  const unsigned int fe_index)
491  {
492  Assert(obj_index < active_fe_indices.size(),
493  ExcIndexRange(obj_index, 0, active_fe_indices.size()));
494 
495  // check whether the given fe_index is within the range of
496  // values that we interpret as "not compressed". if not, then
497  // the index is so large that we cannot accept it. (but this
498  // will not likely happen because it requires someone using an
499  // FECollection that has more than 32k entries.)
500  Assert(is_compressed_entry(fe_index) == false,
501  ExcMessage(
502  "You are using an active_fe_index that is larger than an "
503  "internal limitation for these objects. Try to work with "
504  "hp::FECollection objects that have a more modest size."));
505  Assert(fe_index != invalid_active_fe_index,
506  ExcMessage(
507  "You are using an active_fe_index that is reserved for "
508  "internal purposes for these objects. Try to work with "
509  "hp::FECollection objects that have a more modest size."));
510 
511  active_fe_indices[obj_index] = fe_index;
512  }
513 
514 
515 
516  inline unsigned int
517  DoFLevel::future_fe_index(const unsigned int obj_index) const
518  {
519  Assert(obj_index < future_fe_indices.size(),
520  ExcIndexRange(obj_index, 0, future_fe_indices.size()));
521 
522  if (future_fe_index_set(obj_index))
523  return future_fe_indices[obj_index];
524 
525  return active_fe_index(obj_index);
526  }
527 
528 
529 
530  inline void
531  DoFLevel::set_future_fe_index(const unsigned int obj_index,
532  const unsigned int fe_index)
533  {
534  Assert(obj_index < future_fe_indices.size(),
535  ExcIndexRange(obj_index, 0, future_fe_indices.size()));
536 
537  // check whether the given fe_index is within the range of
538  // values that we interpret as "not compressed". if not, then
539  // the index is so large that we cannot accept it. (but this
540  // will not likely happen because it requires someone using an
541  // FECollection that has more than 32k entries.)
542  Assert(is_compressed_entry(fe_index) == false,
543  ExcMessage(
544  "You are using a future_fe_index that is larger than an "
545  "internal limitation for these objects. Try to work with "
546  "hp::FECollection objects that have a more modest size."));
547  Assert(fe_index != invalid_active_fe_index,
548  ExcMessage(
549  "You are using a future_fe_index that is reserved for "
550  "internal purposes for these objects. Try to work with "
551  "hp::FECollection objects that have a more modest size."));
552 
553  future_fe_indices[obj_index] = fe_index;
554  }
555 
556 
557 
558  inline bool
559  DoFLevel::future_fe_index_set(const unsigned int obj_index) const
560  {
561  Assert(obj_index < future_fe_indices.size(),
562  ExcIndexRange(obj_index, 0, future_fe_indices.size()));
563 
564  return (future_fe_indices[obj_index] != invalid_active_fe_index);
565  }
566 
567 
568 
569  inline void
570  DoFLevel::clear_future_fe_index(const unsigned int obj_index)
571  {
572  Assert(obj_index < future_fe_indices.size(),
573  ExcIndexRange(obj_index, 0, future_fe_indices.size()));
574 
575  future_fe_indices[obj_index] = invalid_active_fe_index;
576  }
577 
578 
579 
580  inline const types::global_dof_index *
581  DoFLevel::get_cell_cache_start(const unsigned int obj_index,
582  const unsigned int dofs_per_cell) const
583  {
584  (void)dofs_per_cell;
585  Assert(
586  (obj_index < cell_cache_offsets.size()) &&
587  (cell_cache_offsets[obj_index] + dofs_per_cell <=
588  cell_dof_indices_cache.size()),
589  ExcMessage(
590  "You are trying to access an element of the cache that stores "
591  "the indices of all degrees of freedom that live on one cell. "
592  "However, this element does not exist. Did you forget to call "
593  "DoFHandler::distribute_dofs(), or did you forget to call it "
594  "again after changing the active_fe_index of one of the cells?"));
595 
596  return &cell_dof_indices_cache[cell_cache_offsets[obj_index]];
597  }
598 
599 
600 
601  template <class Archive>
602  inline void
603  DoFLevel::serialize(Archive &ar, const unsigned int)
604  {
605  ar & this->active_fe_indices;
606  ar & this->cell_cache_offsets;
607  ar & this->cell_dof_indices_cache;
608  ar & this->dof_indices;
609  ar & this->dof_offsets;
610  ar & this->future_fe_indices;
611  }
612  } // namespace hp
613 
614 } // namespace internal
615 
616 DEAL_II_NAMESPACE_CLOSE
617 
618 #endif
unsigned int offset_type
Definition: dof_level.h:115
std::vector< offset_type > cell_cache_offsets
Definition: dof_level.h:201
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
unsigned short int active_fe_index_type
Definition: dof_level.h:120
static ::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1407
Definition: hp.h:117
std::vector< offset_type > dof_offsets
Definition: dof_level.h:189
unsigned int global_dof_index
Definition: types.h:89
std::vector< types::global_dof_index > dof_indices
Definition: dof_level.h:196
signed short int signed_active_fe_index_type
Definition: dof_level.h:126
std::vector< types::global_dof_index > cell_dof_indices_cache
Definition: dof_level.h:208
static const active_fe_index_type invalid_active_fe_index
Definition: dof_level.h:132
std::vector< active_fe_index_type > future_fe_indices
Definition: dof_level.h:175
std::vector< active_fe_index_type > active_fe_indices
Definition: dof_level.h:161