Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
dof_faces.h
1 // ---------------------------------------------------------------------
2 // @f$Id: dof_faces.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 2006 - 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_faces_h
18 #define __deal2__hp_dof_faces_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/hp/fe_collection.h>
23 
24 #include <vector>
25 
27 
28 namespace hp
29 {
30  template <int dim, int spacedim>
31  class FECollection;
32 }
33 
34 
35 namespace internal
36 {
37  namespace hp
38  {
103  template <int structdim>
105  {
106  public:
115  std::vector<unsigned int> dof_offsets;
116 
123  std::vector<types::global_dof_index> dofs;
124 
149  template <int dim, int spacedim>
150  void
151  set_dof_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
152  const unsigned int obj_index,
153  const unsigned int fe_index,
154  const unsigned int local_index,
155  const types::global_dof_index global_index,
156  const unsigned int obj_level);
157 
181  template <int dim, int spacedim>
183  get_dof_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
184  const unsigned int obj_index,
185  const unsigned int fe_index,
186  const unsigned int local_index,
187  const unsigned int obj_level) const;
188 
212  template <int dim, int spacedim>
213  unsigned int
214  n_active_fe_indices (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
215  const unsigned int obj_index) const;
216 
222  template <int dim, int spacedim>
224  nth_active_fe_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
225  const unsigned int obj_level,
226  const unsigned int obj_index,
227  const unsigned int n) const;
228 
235  template <int dim, int spacedim>
236  bool
237  fe_index_is_active (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
238  const unsigned int obj_index,
239  const unsigned int fe_index,
240  const unsigned int obj_level) const;
241 
247  std::size_t memory_consumption () const;
248  };
249 
250 
251 
276  template<int dim>
278 
279 
287  template<>
289  {
290  public:
296  std::size_t memory_consumption () const;
297  };
298 
305  template<>
307  {
308  public:
313 
319  std::size_t memory_consumption () const;
320  };
321 
329  template<>
331  {
332  public:
337 
342 
348  std::size_t memory_consumption () const;
349  };
350 
351 
352  // --------------------- inline and template functions ------------------
353 
354  template <int structdim>
355  template <int dim, int spacedim>
356  inline
359  get_dof_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
360  const unsigned int obj_index,
361  const unsigned int fe_index,
362  const unsigned int local_index,
363  const unsigned int obj_level) const
364  {
366  ExcMessage ("You need to specify a FE index when working "
367  "with hp DoFHandlers"));
368  Assert (&dof_handler != 0,
369  ExcMessage ("No DoFHandler is specified for this iterator"));
370  Assert (&dof_handler.get_fe() != 0,
371  ExcMessage ("No finite element collection is associated with "
372  "this DoFHandler"));
373  Assert (fe_index < dof_handler.get_fe().size(),
374  ExcIndexRange (fe_index, 0, dof_handler.get_fe().size()));
375  Assert (local_index <
376  dof_handler.get_fe()[fe_index].template n_dofs_per_object<structdim>(),
377  ExcIndexRange(local_index, 0,
378  dof_handler.get_fe()[fe_index]
379  .template n_dofs_per_object<structdim>()));
380  Assert (obj_index < dof_offsets.size(),
381  ExcIndexRange (obj_index, 0, dof_offsets.size()));
382 
383  // make sure we are on an
384  // object for which DoFs have
385  // been allocated at all
386  Assert (dof_offsets[obj_index] != numbers::invalid_unsigned_int,
387  ExcMessage ("You are trying to access degree of freedom "
388  "information for an object on which no such "
389  "information is available"));
390 
391  Assert (structdim<dim, ExcMessage ("This object can not be used for cells."));
392 
393  // there may be multiple finite elements associated with
394  // this object. hop along the list of index sets until we
395  // find the one with the correct fe_index, and then poke
396  // into that part. trigger an exception if we can't find a
397  // set for this particular fe_index
398  const types::global_dof_index starting_offset = dof_offsets[obj_index];
399  const types::global_dof_index *pointer = &dofs[starting_offset];
400  while (true)
401  {
402  Assert (*pointer != numbers::invalid_dof_index,
403  ExcInternalError());
404  if (*pointer == fe_index)
405  return *(pointer + 1 + local_index);
406  else
407  pointer += static_cast<types::global_dof_index>(
408  dof_handler.get_fe()[*pointer]
409  .template n_dofs_per_object<structdim>() + 1);
410  }
411  }
412 
413 
414 
415  template <int structdim>
416  template <int dim, int spacedim>
417  inline
418  void
420  set_dof_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
421  const unsigned int obj_index,
422  const unsigned int fe_index,
423  const unsigned int local_index,
424  const types::global_dof_index global_index,
425  const unsigned int obj_level)
426  {
428  ExcMessage ("You need to specify a FE index when working "
429  "with hp DoFHandlers"));
430  Assert (&dof_handler != 0,
431  ExcMessage ("No DoFHandler is specified for this iterator"));
432  Assert (&dof_handler.get_fe() != 0,
433  ExcMessage ("No finite element collection is associated with "
434  "this DoFHandler"));
435  Assert (fe_index < dof_handler.get_fe().size(),
436  ExcIndexRange (fe_index, 0, dof_handler.get_fe().size()));
437  Assert (local_index <
438  dof_handler.get_fe()[fe_index].template n_dofs_per_object<structdim>(),
439  ExcIndexRange(local_index, 0,
440  dof_handler.get_fe()[fe_index]
441  .template n_dofs_per_object<structdim>()));
442  Assert (obj_index < dof_offsets.size(),
443  ExcIndexRange (obj_index, 0, dof_offsets.size()));
444 
445  // make sure we are on an
446  // object for which DoFs have
447  // been allocated at all
448  Assert (dof_offsets[obj_index] != numbers::invalid_unsigned_int,
449  ExcMessage ("You are trying to access degree of freedom "
450  "information for an object on which no such "
451  "information is available"));
452 
453  Assert (structdim<dim, ExcMessage ("This object can not be used for cells."));
454 
455  // there may be multiple finite elements associated with
456  // this object. hop along the list of index sets until we
457  // find the one with the correct fe_index, and then poke
458  // into that part. trigger an exception if we can't find a
459  // set for this particular fe_index
460  const types::global_dof_index starting_offset = dof_offsets[obj_index];
461  types::global_dof_index *pointer = &dofs[starting_offset];
462  while (true)
463  {
464  Assert (*pointer != numbers::invalid_dof_index,
465  ExcInternalError());
466  if (*pointer == fe_index)
467  {
468  *(pointer + 1 + local_index) = global_index;
469  return;
470  }
471  else
472  pointer += dof_handler.get_fe()[*pointer]
473  .template n_dofs_per_object<structdim>() + 1;
474  }
475  }
476 
477 
478 
479  template <int structdim>
480  template <int dim, int spacedim>
481  inline
482  unsigned int
484  n_active_fe_indices (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
485  const unsigned int obj_index) const
486  {
487  Assert (&dof_handler != 0,
488  ExcMessage ("No DoFHandler is specified for this iterator"));
489  Assert (&dof_handler.get_fe() != 0,
490  ExcMessage ("No finite element collection is associated with "
491  "this DoFHandler"));
492  Assert (obj_index < dof_offsets.size(),
493  ExcIndexRange (obj_index, 0, dof_offsets.size()));
494 
495  // make sure we are on an
496  // object for which DoFs have
497  // been allocated at all
498  if (dof_offsets[obj_index] == numbers::invalid_unsigned_int)
499  return 0;
500 
501  Assert (structdim<dim, ExcMessage ("This object can not be used for cells."));
502 
503  // there may be multiple finite elements associated with this
504  // object. hop along the list of index sets until we find the
505  // one with the correct fe_index, and then poke into that
506  // part. trigger an exception if we can't find a set for this
507  // particular fe_index
508  const unsigned int starting_offset = dof_offsets[obj_index];
509  const types::global_dof_index *pointer = &dofs[starting_offset];
510  unsigned int counter = 0;
511  while (true)
512  {
513  if (*pointer == numbers::invalid_dof_index)
514  // end of list reached
515  return counter;
516  else
517  {
518  ++counter;
519  pointer += dof_handler.get_fe()[*pointer]
520  .template n_dofs_per_object<structdim>() + 1;
521  }
522  }
523  }
524 
525 
526 
527  template <int structdim>
528  template <int dim, int spacedim>
529  inline
532  nth_active_fe_index (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
533  const unsigned int obj_level,
534  const unsigned int obj_index,
535  const unsigned int n) const
536  {
537  Assert (&dof_handler != 0,
538  ExcMessage ("No DoFHandler is specified for this iterator"));
539  Assert (&dof_handler.get_fe() != 0,
540  ExcMessage ("No finite element collection is associated with "
541  "this DoFHandler"));
542  Assert (obj_index < dof_offsets.size(),
543  ExcIndexRange (obj_index, 0, dof_offsets.size()));
544 
545  // make sure we are on an
546  // object for which DoFs have
547  // been allocated at all
548  Assert (dof_offsets[obj_index] != numbers::invalid_unsigned_int,
549  ExcMessage ("You are trying to access degree of freedom "
550  "information for an object on which no such "
551  "information is available"));
552 
553  Assert (structdim<dim, ExcMessage ("This object can not be used for cells."));
554 
555  Assert (n < n_active_fe_indices(dof_handler, obj_index),
556  ExcIndexRange (n, 0,
557  n_active_fe_indices(dof_handler, obj_index)));
558 
559  // there may be multiple finite elements associated with
560  // this object. hop along the list of index sets until we
561  // find the one with the correct fe_index, and then poke
562  // into that part. trigger an exception if we can't find a
563  // set for this particular fe_index
564  const unsigned int starting_offset = dof_offsets[obj_index];
565  const types::global_dof_index *pointer = &dofs[starting_offset];
566  unsigned int counter = 0;
567  while (true)
568  {
569  Assert (*pointer != numbers::invalid_dof_index,
570  ExcInternalError());
571 
572  const unsigned int fe_index = *pointer;
573 
574  Assert (fe_index < dof_handler.get_fe().size(),
575  ExcInternalError());
576 
577  if (counter == n)
578  return fe_index;
579 
580  ++counter;
581  pointer += dof_handler.get_fe()[fe_index]
582  .template n_dofs_per_object<structdim>() + 1;
583  }
584  }
585 
586 
587 
588  template <int structdim>
589  template <int dim, int spacedim>
590  inline
591  bool
593  fe_index_is_active (const ::hp::DoFHandler<dim,spacedim> &dof_handler,
594  const unsigned int obj_index,
595  const unsigned int fe_index,
596  const unsigned int obj_level) const
597  {
598  Assert (&dof_handler != 0,
599  ExcMessage ("No DoFHandler is specified for this iterator"));
600  Assert (&dof_handler.get_fe() != 0,
601  ExcMessage ("No finite element collection is associated with "
602  "this DoFHandler"));
603  Assert (obj_index < dof_offsets.size(),
604  ExcIndexRange (obj_index, 0, static_cast<unsigned int>(dof_offsets.size())));
606  ExcMessage ("You need to specify a FE index when working "
607  "with hp DoFHandlers"));
608  Assert (fe_index < dof_handler.get_fe().size(),
609  ExcIndexRange (fe_index, 0, dof_handler.get_fe().size()));
610 
611  // make sure we are on an
612  // object for which DoFs have
613  // been allocated at all
614  Assert (dof_offsets[obj_index] != numbers::invalid_unsigned_int,
615  ExcMessage ("You are trying to access degree of freedom "
616  "information for an object on which no such "
617  "information is available"));
618 
619  Assert (structdim<dim, ExcMessage ("This object can not be used for cells."));
620 
621  // there may be multiple finite elements associated with
622  // this object. hop along the list of index sets until we
623  // find the one with the correct fe_index, and then poke
624  // into that part. trigger an exception if we can't find a
625  // set for this particular fe_index
626  const types::global_dof_index starting_offset = dof_offsets[obj_index];
627  const types::global_dof_index *pointer = &dofs[starting_offset];
628  while (true)
629  {
630  if (*pointer == numbers::invalid_dof_index)
631  // end of list reached
632  return false;
633  else if (*pointer == fe_index)
634  return true;
635  else
636  pointer += static_cast<types::global_dof_index>(
637  dof_handler.get_fe()[*pointer]
638  .template n_dofs_per_object<structdim>()+1);
639  }
640  }
641 
642 
643  } // namespace hp
644 
645 } // namespace internal
646 
647 DEAL_II_NAMESPACE_CLOSE
648 
649 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:191
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:420
::ExceptionBase & ExcMessage(std::string arg1)
internal::hp::DoFIndicesOnFacesOrEdges< 1 > lines
Definition: dof_faces.h:312
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:359
unsigned int global_dof_index
Definition: types.h:100
unsigned int n_active_fe_indices(const ::hp::DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_index) const
Definition: dof_faces.h:484
#define Assert(cond, exc)
Definition: exceptions.h:299
std::vector< types::global_dof_index > dofs
Definition: dof_faces.h:123
std::size_t memory_consumption() const
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
std::vector< unsigned int > dof_offsets
Definition: dof_faces.h:115
Definition: hp.h:103
internal::hp::DoFIndicesOnFacesOrEdges< 2 > quads
Definition: dof_faces.h:341
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:532
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:593
const types::global_dof_index invalid_dof_index
Definition: types.h:217
::ExceptionBase & ExcInternalError()
internal::hp::DoFIndicesOnFacesOrEdges< 1 > lines
Definition: dof_faces.h:336