dune-fem  2.6-git
hpdg/anisotropic.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_HPDG_SPACE_DISCONTINUOUSGALERKIN_ANISOTROPIC_HH
2 #define DUNE_FEM_HPDG_SPACE_DISCONTINUOUSGALERKIN_ANISOTROPIC_HH
3 
4 #include <algorithm>
5 #include <utility>
6 
7 #include <dune/grid/common/gridenums.hh>
8 
14 
16 
17 #include "blockmapper.hh"
18 #include "space.hh"
19 
20 namespace Dune
21 {
22 
23  namespace Fem
24  {
25 
26  namespace hpDG
27  {
28 
29  // Internal forward declaration
30  // ----------------------------
31 
32  template< class FunctionSpace, class GridPart, int order, bool caching = true >
33  class AnisotropicDiscontinuousGalerkinSpace;
34 
35 
36 
37 #ifndef DOXYGEN
38 
39  // AnisotropicDiscontinuousGalerkinSpaceTraits
40  // -------------------------------------------
41 
42  template< class FunctionSpace, class GridPart, int order, bool caching >
43  struct AnisotropicDiscontinuousGalerkinSpaceTraits
44  {
45  using DiscreteFunctionSpaceType = hpDG::AnisotropicDiscontinuousGalerkinSpace< FunctionSpace, GridPart, order >;
46 
47  using FunctionSpaceType = FunctionSpace;
48 
49  using GridPartType = GridPart;
50 
51  using BasisFunctionSetsType = hpDG::AnisotropicBasisFunctionSets< FunctionSpaceType, GridPartType, order, caching >;
52  using BasisFunctionSetType = typename BasisFunctionSetsType::BasisFunctionSetType;
53 
54  static const int codimension = BasisFunctionSetType::EntityType::codimension;
55 
56  using BlockMapperType = hpDG::DiscontinuousGalerkinBlockMapper< GridPartType, BasisFunctionSetsType >;
57  static const int localBlockSize = BasisFunctionSetsType::localBlockSize;
58 
59  typedef Hybrid::IndexRange< int, localBlockSize > LocalBlockIndices;
60 
61  template< class DiscreteFunction, class Operation = Dune::Fem::DFCommunicationOperation::Copy >
62  struct CommDataHandle
63  {
64  using OperationType = Operation;
66  };
67  };
68 
69 #endif // #ifndef DOXYGEN
70 
71 
72 
73  // AnisotropicDiscontinuousGalerkinSpace
74  // --------------------------------
75 
85  template< class FunctionSpace, class GridPart, int order, bool caching >
87  : public hpDG::DiscontinuousGalerkinSpace< AnisotropicDiscontinuousGalerkinSpaceTraits< FunctionSpace, GridPart, order, caching > >
88  {
90 
91  public:
94 
96  const Dune::InterfaceType interface = Dune::InteriorBorder_All_Interface,
97  const Dune::CommunicationDirection direction = Dune::ForwardCommunication )
98  : BaseType( gridPart, BasisFunctionSetsType(), defaultKey(), interface, direction )
99  {}
100 
101  template< class Function >
103  const Dune::InterfaceType interface = Dune::InteriorBorder_All_Interface,
104  const Dune::CommunicationDirection direction = Dune::ForwardCommunication )
105  : BaseType( gridPart, BasisFunctionSetsType(), defaultKey(), function, interface, direction )
106  {}
107 
108  private:
109  static typename BaseType::KeyType defaultKey ()
110  {
111  typename BaseType::KeyType key;
112  std::fill( key.begin(), key.end(), order );
113  return std::move( key );
114  }
115  };
116 
117  } // namespace hpDG
118 
119 
120 
121 #ifndef DOXYGEN
122 
123  // DefaultLocalRestrictProlong
124  // ---------------------------
125 
126  template< class FunctionSpace, class GridPart, int order, bool caching >
127  class DefaultLocalRestrictProlong< hpDG::AnisotropicDiscontinuousGalerkinSpace< FunctionSpace, GridPart, order, caching > >
128  : public DiscontinuousGalerkinLocalRestrictProlong< hpDG::AnisotropicDiscontinuousGalerkinSpace< FunctionSpace, GridPart, order, caching >, false >
129  {
130  using BaseType = DiscontinuousGalerkinLocalRestrictProlong< hpDG::AnisotropicDiscontinuousGalerkinSpace< FunctionSpace, GridPart, order, caching >, false >;
131 
132  public:
133  explicit DefaultLocalRestrictProlong ( const typename BaseType::DiscreteFunctionSpaceType &space )
134  : BaseType( space )
135  {}
136  };
137 
138 
139 
140  // ISTLParallelMatrixAdapter
141  // -------------------------
142 
143  template< class Matrix, class Traits >
144  struct ISTLParallelMatrixAdapter< Matrix, hpDG::DiscontinuousGalerkinSpace< Traits > >
145  {
146  using Type = DGParallelMatrixAdapter< Matrix >;
147  };
148 
149 #endif //#ifndef DOXYGEN
150 
151  } // namespace Fem
152 
153 } // namespace Dune
154 
155 #endif // #ifndef DUNE_FEM_HPDG_SPACE_DISCONTINUOUSGALERKIN_ANISOTROPIC_HH
Definition: bindguard.hh:11
Abstract class representing a function.
Definition: common/function.hh:50
Default communication handler for discrete functions.
Definition: defaultcommhandler.hh:29
Implementation of an -adaptive discrete function space using anisotropic product Legendre polynomials...
Definition: hpdg/anisotropic.hh:88
AnisotropicDiscontinuousGalerkinSpace(GridPartType &gridPart, Function function, const Dune::InterfaceType interface=Dune::InteriorBorder_All_Interface, const Dune::CommunicationDirection direction=Dune::ForwardCommunication)
Definition: hpdg/anisotropic.hh:102
typename BaseType::BasisFunctionSetsType BasisFunctionSetsType
Definition: hpdg/anisotropic.hh:93
AnisotropicDiscontinuousGalerkinSpace(GridPartType &gridPart, const Dune::InterfaceType interface=Dune::InteriorBorder_All_Interface, const Dune::CommunicationDirection direction=Dune::ForwardCommunication)
Definition: hpdg/anisotropic.hh:95
typename BaseType::GridPartType GridPartType
Definition: hpdg/anisotropic.hh:92
Generic implementation of a -adaptive discontinuous finite element space.
Definition: hpdg/space.hh:45
typename Traits::BasisFunctionSetsType BasisFunctionSetsType
basis function sets type
Definition: hpdg/space.hh:55
typename BasisFunctionSetsType::KeyType KeyType
key type identifying a basis function set
Definition: hpdg/space.hh:57
const KeyType & key(const EntityType &entity) const
get identifiying basis function set key assigned to given entity
Definition: hpdg/space.hh:204