Reference documentation for deal.II version Git 840703ecf3 2019-10-14 17:15:20 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
dof_faces.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2006 - 2018 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_faces_h
17 #define dealii_hp_dof_faces_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/hp/fe_collection.h>
23 
24 #include <vector>
25 
26 DEAL_II_NAMESPACE_OPEN
27 
28 namespace internal
29 {
30  namespace hp
31  {
94  template <int structdim>
96  {
97  public:
105  std::vector<unsigned int> dof_offsets;
106 
111  std::vector<types::global_dof_index> dofs;
112 
125  template <int dim, int spacedim>
126  void
127  set_dof_index(const ::hp::DoFHandler<dim, spacedim> &dof_handler,
128  const unsigned int obj_index,
129  const unsigned int fe_index,
130  const unsigned int local_index,
131  const types::global_dof_index global_index,
132  const unsigned int obj_level);
133 
145  template <int dim, int spacedim>
147  get_dof_index(const ::hp::DoFHandler<dim, spacedim> &dof_handler,
148  const unsigned int obj_index,
149  const unsigned int fe_index,
150  const unsigned int local_index,
151  const unsigned int obj_level) const;
152 
164  template <int dim, int spacedim>
165  unsigned int
167  const ::hp::DoFHandler<dim, spacedim> &dof_handler,
168  const unsigned int obj_index) const;
169 
173  template <int dim, int spacedim>
176  const ::hp::DoFHandler<dim, spacedim> &dof_handler,
177  const unsigned int obj_level,
178  const unsigned int obj_index,
179  const unsigned int n) const;
180 
185  template <int dim, int spacedim>
186  bool
188  const ::hp::DoFHandler<dim, spacedim> &dof_handler,
189  const unsigned int obj_index,
190  const unsigned int fe_index,
191  const unsigned int obj_level) const;
192 
197  std::size_t
198  memory_consumption() const;
199 
204  template <class Archive>
205  void
206  serialize(Archive &ar, const unsigned int version);
207  };
208 
209 
210 
238  template <int dim>
240 
241 
249  template <>
251  {
252  public:
257  std::size_t
258  memory_consumption() const;
259 
264  template <class Archive>
265  void
266  serialize(Archive &ar, const unsigned int version);
267  };
268 
276  template <>
278  {
279  public:
284 
289  std::size_t
290  memory_consumption() const;
291 
296  template <class Archive>
297  void
298  serialize(Archive &ar, const unsigned int version);
299  };
300 
308  template <>
310  {
311  public:
316 
321 
326  std::size_t
327  memory_consumption() const;
328 
333  template <class Archive>
334  void
335  serialize(Archive &ar, const unsigned int version);
336  };
337 
338 
339  // --------------------- inline and template functions ------------------
340  template <class Archive>
341  void
342  DoFIndicesOnFaces<1>::serialize(Archive &, const unsigned int)
343  {}
344 
345 
346 
347  inline std::size_t
349  {
350  return 0;
351  }
352 
353 
354 
355  template <class Archive>
356  void
357  DoFIndicesOnFaces<2>::serialize(Archive &ar, const unsigned int)
358  {
359  ar &lines;
360  }
361 
362 
363 
364  inline std::size_t
366  {
368  }
369 
370 
371 
372  template <class Archive>
373  void
374  DoFIndicesOnFaces<3>::serialize(Archive &ar, const unsigned int)
375  {
376  ar &lines &quads;
377  }
378 
379 
380 
381  inline std::size_t
383  {
386  }
387 
388  template <int structdim>
389  template <int dim, int spacedim>
392  const ::hp::DoFHandler<dim, spacedim> &dof_handler,
393  const unsigned int obj_index,
394  const unsigned int fe_index,
395  const unsigned int local_index,
396  const unsigned int /*obj_level*/) const
397  {
398  Assert((fe_index !=
400  ExcMessage("You need to specify a FE index when working "
401  "with hp DoFHandlers"));
402  Assert(fe_index < dof_handler.get_fe_collection().size(),
403  ExcIndexRange(fe_index,
404  0,
405  dof_handler.get_fe_collection().size()));
406  Assert(local_index < dof_handler.get_fe(fe_index)
407  .template n_dofs_per_object<structdim>(),
408  ExcIndexRange(local_index,
409  0,
410  dof_handler.get_fe(fe_index)
411  .template n_dofs_per_object<structdim>()));
412  Assert(obj_index < dof_offsets.size(),
413  ExcIndexRange(obj_index, 0, dof_offsets.size()));
414 
415  // make sure we are on an
416  // object for which DoFs have
417  // been allocated at all
419  ExcMessage("You are trying to access degree of freedom "
420  "information for an object on which no such "
421  "information is available"));
422 
423  Assert(structdim < dim,
424  ExcMessage("This object can not be used for cells."));
425 
426  // there may be multiple finite elements associated with
427  // this object. hop along the list of index sets until we
428  // find the one with the correct fe_index, and then poke
429  // into that part. trigger an exception if we can't find a
430  // set for this particular fe_index
431  const types::global_dof_index starting_offset = dof_offsets[obj_index];
432  const types::global_dof_index *pointer = &dofs[starting_offset];
433  while (true)
434  {
436  if (*pointer == fe_index)
437  return *(pointer + 1 + local_index);
438  else
439  pointer += static_cast<types::global_dof_index>(
440  dof_handler.get_fe(*pointer)
441  .template n_dofs_per_object<structdim>() +
442  1);
443  }
444  }
445 
446 
447 
448  template <int structdim>
449  template <int dim, int spacedim>
450  inline void
452  const ::hp::DoFHandler<dim, spacedim> &dof_handler,
453  const unsigned int obj_index,
454  const unsigned int fe_index,
455  const unsigned int local_index,
456  const types::global_dof_index global_index,
457  const unsigned int /*obj_level*/)
458  {
459  Assert((fe_index !=
461  ExcMessage("You need to specify a FE index when working "
462  "with hp DoFHandlers"));
463  Assert(fe_index < dof_handler.get_fe_collection().size(),
464  ExcIndexRange(fe_index,
465  0,
466  dof_handler.get_fe_collection().size()));
467  Assert(local_index < dof_handler.get_fe(fe_index)
468  .template n_dofs_per_object<structdim>(),
469  ExcIndexRange(local_index,
470  0,
471  dof_handler.get_fe(fe_index)
472  .template n_dofs_per_object<structdim>()));
473  Assert(obj_index < dof_offsets.size(),
474  ExcIndexRange(obj_index, 0, dof_offsets.size()));
475 
476  // make sure we are on an
477  // object for which DoFs have
478  // been allocated at all
480  ExcMessage("You are trying to access degree of freedom "
481  "information for an object on which no such "
482  "information is available"));
483 
484  Assert(structdim < dim,
485  ExcMessage("This object can not be used for cells."));
486 
487  // there may be multiple finite elements associated with
488  // this object. hop along the list of index sets until we
489  // find the one with the correct fe_index, and then poke
490  // into that part. trigger an exception if we can't find a
491  // set for this particular fe_index
492  const types::global_dof_index starting_offset = dof_offsets[obj_index];
493  types::global_dof_index * pointer = &dofs[starting_offset];
494  while (true)
495  {
497  if (*pointer == fe_index)
498  {
499  *(pointer + 1 + local_index) = global_index;
500  return;
501  }
502  else
503  pointer += dof_handler.get_fe(*pointer)
504  .template n_dofs_per_object<structdim>() +
505  1;
506  }
507  }
508 
509 
510 
511  template <int structdim>
512  template <int dim, int spacedim>
513  inline unsigned int
515  const ::hp::DoFHandler<dim, spacedim> &dof_handler,
516  const unsigned int obj_index) const
517  {
518  Assert(obj_index < dof_offsets.size(),
519  ExcIndexRange(obj_index, 0, dof_offsets.size()));
520 
521  // make sure we are on an
522  // object for which DoFs have
523  // been allocated at all
524  if (dof_offsets[obj_index] == numbers::invalid_unsigned_int)
525  return 0;
526 
527  Assert(structdim < dim,
528  ExcMessage("This object can not be used for cells."));
529 
530  // there may be multiple finite elements associated with this
531  // object. hop along the list of index sets until we find the
532  // one with the correct fe_index, and then poke into that
533  // part. trigger an exception if we can't find a set for this
534  // particular fe_index
535  const unsigned int starting_offset = dof_offsets[obj_index];
536  const types::global_dof_index *pointer = &dofs[starting_offset];
537  unsigned int counter = 0;
538  while (true)
539  {
540  if (*pointer == numbers::invalid_dof_index)
541  // end of list reached
542  return counter;
543  else
544  {
545  ++counter;
546  pointer += dof_handler.get_fe(*pointer)
547  .template n_dofs_per_object<structdim>() +
548  1;
549  }
550  }
551  }
552 
553 
554 
555  template <int structdim>
556  template <int dim, int spacedim>
559  const ::hp::DoFHandler<dim, spacedim> &dof_handler,
560  const unsigned int /*obj_level*/,
561  const unsigned int obj_index,
562  const unsigned int n) const
563  {
564  Assert(obj_index < dof_offsets.size(),
565  ExcIndexRange(obj_index, 0, dof_offsets.size()));
566 
567  // make sure we are on an
568  // object for which DoFs have
569  // been allocated at all
571  ExcMessage("You are trying to access degree of freedom "
572  "information for an object on which no such "
573  "information is available"));
574 
575  Assert(structdim < dim,
576  ExcMessage("This object can not be used for cells."));
577 
578  Assert(n < n_active_fe_indices(dof_handler, obj_index),
579  ExcIndexRange(n, 0, n_active_fe_indices(dof_handler, obj_index)));
580 
581  // there may be multiple finite elements associated with
582  // this object. hop along the list of index sets until we
583  // find the one with the correct fe_index, and then poke
584  // into that part. trigger an exception if we can't find a
585  // set for this particular fe_index
586  const unsigned int starting_offset = dof_offsets[obj_index];
587  const types::global_dof_index *pointer = &dofs[starting_offset];
588  unsigned int counter = 0;
589  while (true)
590  {
592 
593  const unsigned int fe_index = *pointer;
594 
595  Assert(fe_index < dof_handler.get_fe_collection().size(),
596  ExcInternalError());
597 
598  if (counter == n)
599  return fe_index;
600 
601  ++counter;
602  pointer += dof_handler.get_fe(fe_index)
603  .template n_dofs_per_object<structdim>() +
604  1;
605  }
606  }
607 
608 
609 
610  template <int structdim>
611  template <int dim, int spacedim>
612  inline bool
614  const ::hp::DoFHandler<dim, spacedim> &dof_handler,
615  const unsigned int obj_index,
616  const unsigned int fe_index,
617  const unsigned int /*obj_level*/) const
618  {
619  Assert(obj_index < dof_offsets.size(),
620  ExcIndexRange(obj_index,
621  0,
622  static_cast<unsigned int>(dof_offsets.size())));
623  Assert((fe_index !=
625  ExcMessage("You need to specify a FE index when working "
626  "with hp DoFHandlers"));
627  Assert(fe_index < dof_handler.get_fe_collection().size(),
628  ExcIndexRange(fe_index,
629  0,
630  dof_handler.get_fe_collection().size()));
631 
632  // make sure we are on an
633  // object for which DoFs have
634  // been allocated at all
636  ExcMessage("You are trying to access degree of freedom "
637  "information for an object on which no such "
638  "information is available"));
639 
640  Assert(structdim < dim,
641  ExcMessage("This object can not be used for cells."));
642 
643  // there may be multiple finite elements associated with
644  // this object. hop along the list of index sets until we
645  // find the one with the correct fe_index, and then poke
646  // into that part. trigger an exception if we can't find a
647  // set for this particular fe_index
648  const types::global_dof_index starting_offset = dof_offsets[obj_index];
649  const types::global_dof_index *pointer = &dofs[starting_offset];
650  while (true)
651  {
652  if (*pointer == numbers::invalid_dof_index)
653  // end of list reached
654  return false;
655  else if (*pointer == fe_index)
656  return true;
657  else
658  pointer += static_cast<types::global_dof_index>(
659  dof_handler.get_fe(*pointer)
660  .template n_dofs_per_object<structdim>() +
661  1);
662  }
663  }
664 
665  template <int structdim>
666  template <class Archive>
667  void
669  const unsigned int)
670  {
671  ar &dofs;
672  ar &dof_offsets;
673  }
674 
675 
676  } // namespace hp
677 
678 } // namespace internal
679 
680 DEAL_II_NAMESPACE_CLOSE
681 
682 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:187
void set_dof_index(const ::hp::DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index, const unsigned int obj_level)
Definition: dof_faces.h:451
types::global_dof_index get_dof_index(const ::hp::DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const unsigned int obj_level) const
Definition: dof_faces.h:391
unsigned int n_active_fe_indices(const ::hp::DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_index) const
Definition: dof_faces.h:514
internal::hp::DoFIndicesOnFacesOrEdges< 1 > lines
Definition: dof_faces.h:283
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1407
void serialize(Archive &ar, const unsigned int version)
Definition: dof_faces.h:668
std::size_t memory_consumption() const
Definition: dof_faces.cc:30
std::vector< types::global_dof_index > dofs
Definition: dof_faces.h:111
std::vector< unsigned int > dof_offsets
Definition: dof_faces.h:105
Definition: hp.h:117
internal::hp::DoFIndicesOnFacesOrEdges< 2 > quads
Definition: dof_faces.h:320
bool fe_index_is_active(const ::hp::DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_index, const unsigned int fe_index, const unsigned int obj_level) const
Definition: dof_faces.h:613
unsigned int global_dof_index
Definition: types.h:89
types::global_dof_index nth_active_fe_index(const ::hp::DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const unsigned int n) const
Definition: dof_faces.h:558
const types::global_dof_index invalid_dof_index
Definition: types.h:202
internal::hp::DoFIndicesOnFacesOrEdges< 1 > lines
Definition: dof_faces.h:315
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
static ::ExceptionBase & ExcInternalError()