Reference documentation for deal.II version 8.1.0
dof_handler.h
1 // ---------------------------------------------------------------------
2 // @f$Id: dof_handler.h 31932 2013-12-08 02:15:54Z heister @f$
3 //
4 // Copyright (C) 2005 - 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_handler_h
18 #define __deal2__hp_dof_handler_h
19 
20 
21 
22 #include <deal.II/base/config.h>
24 #include <deal.II/base/template_constraints.h>
25 #include <deal.II/base/smartpointer.h>
26 #include <deal.II/dofs/function_map.h>
27 #include <deal.II/dofs/dof_iterator_selector.h>
28 #include <deal.II/dofs/number_cache.h>
29 #include <deal.II/hp/fe_collection.h>
30 
31 #include <vector>
32 #include <map>
33 #include <set>
34 
36 
37 namespace internal
38 {
39  namespace hp
40  {
41  class DoFLevel;
42  template <int> class DoFIndicesOnFaces;
43  template <int> class DoFIndicesOnFacesOrEdges;
44 
45  namespace DoFHandler
46  {
47  struct Implementation;
48  }
49  }
50 }
51 
52 namespace internal
53 {
54  namespace DoFAccessor
55  {
56  struct Implementation;
57  }
58 
59  namespace DoFCellAccessor
60  {
61  struct Implementation;
62  }
63 }
64 
65 
66 
67 namespace hp
68 {
69 
80  template <int dim, int spacedim=dim>
81  class DoFHandler : public Subscriptor
82  {
83  typedef ::internal::DoFHandler::Iterators<DoFHandler<dim,spacedim>, false> ActiveSelector;
84  typedef ::internal::DoFHandler::Iterators<DoFHandler<dim,spacedim>, true> LevelSelector;
85  public:
86  typedef typename ActiveSelector::CellAccessor cell_accessor;
87  typedef typename ActiveSelector::FaceAccessor face_accessor;
88 
89  typedef typename ActiveSelector::line_iterator line_iterator;
90  typedef typename ActiveSelector::active_line_iterator active_line_iterator;
91 
92  typedef typename ActiveSelector::quad_iterator quad_iterator;
93  typedef typename ActiveSelector::active_quad_iterator active_quad_iterator;
94 
95  typedef typename ActiveSelector::hex_iterator hex_iterator;
96  typedef typename ActiveSelector::active_hex_iterator active_hex_iterator;
97 
98  typedef typename ActiveSelector::active_cell_iterator active_cell_iterator;
99 
100  typedef typename ActiveSelector::face_iterator face_iterator;
101  typedef typename ActiveSelector::active_face_iterator active_face_iterator;
102 
103  typedef typename LevelSelector::CellAccessor level_cell_accessor;
104  typedef typename LevelSelector::FaceAccessor level_face_accessor;
105 
106  typedef typename LevelSelector::cell_iterator level_cell_iterator;
107  typedef typename LevelSelector::face_iterator level_face_iterator;
108 
109  typedef level_cell_iterator cell_iterator;
110 
116 
121  static const unsigned int dimension = dim;
122 
127  static const unsigned int space_dimension = spacedim;
128 
144  static const types::global_dof_index invalid_dof_index = numbers::invalid_dof_index;
145 
168  static const unsigned int default_fe_index = numbers::invalid_unsigned_int;
169 
170 
176 
180  virtual ~DoFHandler ();
181 
204  virtual void distribute_dofs (const hp::FECollection<dim,spacedim> &fe);
205 
212  void set_active_fe_indices (const std::vector<unsigned int> &active_fe_indices);
213 
221  void get_active_fe_indices (std::vector<unsigned int> &active_fe_indices) const;
222 
229  virtual void clear ();
230 
262  void renumber_dofs (const std::vector<types::global_dof_index> &new_numbers);
263 
290  unsigned int max_couplings_between_dofs () const;
291 
305  unsigned int max_couplings_between_boundary_dofs () const;
306 
315  cell_iterator begin (const unsigned int level = 0) const;
316 
329  active_cell_iterator begin_active(const unsigned int level = 0) const;
330 
338  cell_iterator end () const;
339 
347  cell_iterator end (const unsigned int level) const;
348 
354  active_cell_iterator end_active (const unsigned int level) const;
355 
358  /*---------------------------------------*/
359 
360 
392  types::global_dof_index n_dofs () const;
393 
403  types::global_dof_index n_dofs(const unsigned int level) const;
404 
409  types::global_dof_index n_boundary_dofs () const;
410 
425  n_boundary_dofs (const FunctionMap &boundary_indicators) const;
426 
436  n_boundary_dofs (const std::set<types::boundary_id> &boundary_indicators) const;
437 
468  types::global_dof_index n_locally_owned_dofs() const;
469 
478  const IndexSet &locally_owned_dofs() const;
479 
480 
497  const std::vector<IndexSet> &
498  locally_owned_dofs_per_processor () const;
499 
526  const std::vector<types::global_dof_index> &
527  n_locally_owned_dofs_per_processor () const;
528 
535  const hp::FECollection<dim,spacedim> &get_fe () const;
536 
541  const Triangulation<dim,spacedim> &get_tria () const;
542 
555  virtual std::size_t memory_consumption () const;
556 
560  DeclException0 (ExcInvalidTriangulation);
564  DeclException0 (ExcNoFESelected);
568  DeclException0 (ExcRenumberingIncomplete);
572  DeclException0 (ExcGridsDoNotMatch);
576  DeclException0 (ExcInvalidBoundaryIndicator);
580  DeclException1 (ExcMatrixHasWrongSize,
581  int,
582  << "The matrix has the wrong dimension " << arg1);
586  DeclException0 (ExcFunctionNotUseful);
590  DeclException1 (ExcNewNumbersNotConsecutive,
592  << "The given list of new dof indices is not consecutive: "
593  << "the index " << arg1 << " does not exist.");
597  DeclException2 (ExcInvalidFEIndex,
598  int, int,
599  << "The mesh contains a cell with an active_fe_index of "
600  << arg1 << ", but the finite element collection only has "
601  << arg2 << " elements");
605  DeclException1 (ExcInvalidLevel,
606  int,
607  << "The given level " << arg1
608  << " is not in the valid range!");
612  DeclException0 (ExcFacesHaveNoLevel);
617  DeclException1 (ExcEmptyLevel,
618  int,
619  << "You tried to do something on level " << arg1
620  << ", but this level is empty.");
621 
622  protected:
623 
629 
646 
647  private:
648 
658  DoFHandler (const DoFHandler &);
659 
669  DoFHandler &operator = (const DoFHandler &);
670 
672  {
673  public:
674  MGVertexDoFs ();
675  ~MGVertexDoFs ();
676  types::global_dof_index get_index (const unsigned int level, const unsigned int dof_number) const;
677  void set_index (const unsigned int level, const unsigned int dof_number, const types::global_dof_index index);
678  };
679 
683  void clear_space ();
684 
685  template<int structdim>
686  types::global_dof_index get_dof_index (const unsigned int obj_level, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index) const;
687 
688  template<int structdim>
689  void set_dof_index (const unsigned int obj_level, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index) const;
690 
707  void create_active_fe_table ();
708 
719  void pre_refinement_action ();
720  void post_refinement_action ();
721 
728  void
729  compute_vertex_dof_identities (std::vector<types::global_dof_index> &new_dof_indices) const;
730 
737  void
738  compute_line_dof_identities (std::vector<types::global_dof_index> &new_dof_indices) const;
739 
746  void
747  compute_quad_dof_identities (std::vector<types::global_dof_index> &new_dof_indices) const;
748 
774  void renumber_dofs_internal (const std::vector<types::global_dof_index> &new_numbers,
776 
777  void renumber_dofs_internal (const std::vector<types::global_dof_index> &new_numbers,
779 
780  void renumber_dofs_internal (const std::vector<types::global_dof_index> &new_numbers,
782 
783  void renumber_dofs_internal (const std::vector<types::global_dof_index> &new_numbers,
785 
793  std::vector<::internal::hp::DoFLevel *> levels;
794 
803 
816 
840  std::vector<types::global_dof_index> vertex_dofs;
841 
863  std::vector<types::global_dof_index> vertex_dofs_offsets;
864 
865  std::vector<MGVertexDoFs> mg_vertex_dofs;
866 
879  std::vector<std::vector<bool> *> has_children;
880 
886  std::vector<boost::signals2::connection> tria_listeners;
887 
891  template <int, class, bool> friend class ::DoFAccessor;
892  template <class, bool> friend class ::DoFCellAccessor;
893  friend struct ::internal::DoFAccessor::Implementation;
894  friend struct ::internal::DoFCellAccessor::Implementation;
895 
903  template <int> friend class ::internal::hp::DoFIndicesOnFacesOrEdges;
904  friend struct ::internal::hp::DoFHandler::Implementation;
905  };
906 
907 
908 
909 #ifndef DOXYGEN
910 
911 
912  /* ----------------------- Inline functions ---------------------------------- */
913 
914  template <int dim, int spacedim>
915  inline
918  {
919  return number_cache.n_global_dofs;
920  }
921 
922 
923  template <int dim, int spacedim>
924  inline
926  DoFHandler<dim,spacedim>::n_dofs (const unsigned int) const
927  {
929  }
930 
931 
932  template <int dim, int spacedim>
935  {
936  return number_cache.n_locally_owned_dofs;
937  }
938 
939 
940  template <int dim, int spacedim>
941  const IndexSet &
943  {
944  return number_cache.locally_owned_dofs;
945  }
946 
947 
948  template <int dim, int spacedim>
949  const std::vector<types::global_dof_index> &
951  {
952  return number_cache.n_locally_owned_dofs_per_processor;
953  }
954 
955 
956  template <int dim, int spacedim>
957  const std::vector<IndexSet> &
959  {
960  return number_cache.locally_owned_dofs_per_processor;
961  }
962 
963 
964 
965  template<int dim, int spacedim>
966  inline
969  {
970  Assert (finite_elements != 0,
971  ExcMessage ("No finite element collection is associated with "
972  "this DoFHandler"));
973  return *finite_elements;
974  }
975 
976 
977  template<int dim, int spacedim>
978  inline
981  {
982  return *tria;
983  }
984 
985  template<int dim, int spacedim>
986  inline
988  {
989  Assert (false, ExcNotImplemented ());
990  }
991 
992  template<int dim, int spacedim>
993  inline
995  {
996  Assert (false, ExcNotImplemented ());
997  }
998 
999  template<int dim, int spacedim>
1000  inline
1002  const unsigned int) const
1003  {
1004  Assert (false, ExcNotImplemented ());
1005  return invalid_dof_index;
1006  }
1007 
1008  template<int dim, int spacedim>
1009  inline
1010  void DoFHandler<dim, spacedim>::MGVertexDoFs::set_index (const unsigned int,
1011  const unsigned int,
1013  {
1014  Assert (false, ExcNotImplemented ());
1015  }
1016 
1017 
1018 #endif
1019 
1020 }
1021 
1022 DEAL_II_NAMESPACE_CLOSE
1023 
1024 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:191
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:525
SmartPointer< const hp::FECollection< dim, spacedim >, hp::DoFHandler< dim, spacedim > > finite_elements
Definition: dof_handler.h:645
void set_index(const unsigned int level, const unsigned int dof_number, const types::global_dof_index index)
std::map< types::boundary_id, const Function< dim > * > type
Definition: function_map.h:56
::ExceptionBase & ExcMessage(std::string arg1)
::internal::hp::DoFIndicesOnFaces< dim > * faces
Definition: dof_handler.h:802
::internal::DoFHandler::NumberCache number_cache
Definition: dof_handler.h:815
std::vector< boost::signals2::connection > tria_listeners
Definition: dof_handler.h:886
const std::vector< types::global_dof_index > & n_locally_owned_dofs_per_processor() const
std::vector<::internal::hp::DoFLevel * > levels
Definition: dof_handler.h:793
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
unsigned int global_dof_index
Definition: types.h:100
#define Assert(cond, exc)
Definition: exceptions.h:299
types::global_dof_index n_dofs() const
unsigned int n_locally_owned_dofs() const
#define DeclException0(Exception0)
Definition: exceptions.h:505
types::global_dof_index get_index(const unsigned int level, const unsigned int dof_number) const
Definition: hp.h:103
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
Definition: dof_handler.h:628
FunctionMap< spacedim >::type FunctionMap
Definition: dof_handler.h:115
std::vector< std::vector< bool > * > has_children
Definition: dof_handler.h:879
::ExceptionBase & ExcNotImplemented()
std::vector< types::global_dof_index > vertex_dofs_offsets
Definition: dof_handler.h:863
const types::global_dof_index invalid_dof_index
Definition: types.h:217
const Triangulation< dim, spacedim > & get_tria() const
const IndexSet & locally_owned_dofs() const
const std::vector< IndexSet > & locally_owned_dofs_per_processor() const
std::vector< types::global_dof_index > vertex_dofs
Definition: dof_handler.h:840
const FiniteElement< dim, spacedim > & get_fe() const