dune-fem  2.8-git
hpdg/orthogonal.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_HPDG_SPACE_DISCONTINUOUSGALERKIN_ORTHOGONAL_HH
2 #define DUNE_FEM_HPDG_SPACE_DISCONTINUOUSGALERKIN_ORTHOGONAL_HH
3 
4 #include <dune/grid/common/gridenums.hh>
5 
7 
13 
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, class Storage = Fem::CachingStorage >
33  class OrthogonalDiscontinuousGalerkinSpace;
34 
35 
36 
37 #ifndef DOXYGEN
38 
39  // OrthogonalDiscontinuousGalerkinSpaceTraits
40  // ------------------------------------------
41 
42  template< class FunctionSpace, class GridPart, int order, class Storage >
43  struct OrthogonalDiscontinuousGalerkinSpaceTraits
44  {
45  using DiscreteFunctionSpaceType = hpDG::OrthogonalDiscontinuousGalerkinSpace< FunctionSpace, GridPart, order, Storage >;
46 
47  using FunctionSpaceType = FunctionSpace;
48 
49  using GridPartType = GridPart;
50 
51  using BasisFunctionSetsType = hpDG::OrthogonalBasisFunctionSets< FunctionSpaceType, GridPartType, order, Storage >;
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  static_assert( localBlockSize == FunctionSpace::dimRange, " dimRange prob ");
59 
60  typedef Hybrid::IndexRange< int, localBlockSize > LocalBlockIndices;
61 
62  template< class DiscreteFunction, class Operation = Dune::Fem::DFCommunicationOperation::Copy >
63  struct CommDataHandle
64  {
65  using OperationType = Operation;
67  };
68  };
69 
70 #endif // #ifndef DOXYGEN
71 
72 
73 
74  // OrthogonalDiscontinuousGalerkinSpace
75  // ------------------------------------
76 
86  template< class FunctionSpace, class GridPart, int order, class Storage >
88  : public hpDG::DiscontinuousGalerkinSpace< OrthogonalDiscontinuousGalerkinSpaceTraits< FunctionSpace, GridPart, order, Storage > >
89  {
91 
92  public:
93 
94  static const int polynomialOrder = order ;
95 
97  using EntityType = typename BaseType::EntityType;
99 
101  const Dune::InterfaceType interface = Dune::InteriorBorder_All_Interface,
102  const Dune::CommunicationDirection direction = Dune::ForwardCommunication )
103  : BaseType( gridPart, BasisFunctionSetsType{}, order, interface, direction )
104  {}
105 
106  template <class Function,
107  std::enable_if_t<
108  std::is_arithmetic<
109  decltype(Function(std::declval<const EntityType>()))>::value,int> i=0>
111  const Dune::InterfaceType interface = Dune::InteriorBorder_All_Interface,
112  const Dune::CommunicationDirection direction = Dune::ForwardCommunication )
113  : BaseType( gridPart, BasisFunctionSetsType{}, order, function, interface, direction )
114  {}
115  };
116 
117  } // namespace hpDG
118 
120  template <class FunctionSpaceImp,
121  class GridPartImp,
122  int polOrd,
123  class Storage,
124  class VolumeQuadratureImp>
126  hpDG::OrthogonalDiscontinuousGalerkinSpace< FunctionSpaceImp, GridPartImp, polOrd, Storage >, VolumeQuadratureImp >
128  hpDG::OrthogonalDiscontinuousGalerkinSpace< FunctionSpaceImp, GridPartImp, polOrd, Storage >, VolumeQuadratureImp >
129  {
132  public:
133  using BaseType::BaseType;
134  };
135 
136 
137 #ifndef DOXYGEN
138 
139  // DefaultLocalRestrictProlong
140  // ---------------------------
141 
142  template< class FunctionSpace, class GridPart, int order, class Storage >
143  class DefaultLocalRestrictProlong< hpDG::OrthogonalDiscontinuousGalerkinSpace< FunctionSpace, GridPart, order, Storage > >
144  : public DiscontinuousGalerkinLocalRestrictProlong< hpDG::OrthogonalDiscontinuousGalerkinSpace< FunctionSpace, GridPart, order, Storage >, false >
145  {
147 
148  public:
149  explicit DefaultLocalRestrictProlong ( const typename BaseType::DiscreteFunctionSpaceType &space )
150  : BaseType( space )
151  {}
152  };
153 
154 
155 #endif // #ifndef DOXYGEN
156 
157 
158  namespace Capabilities
159  {
161  // hpDG::OrthogonalDiscontinuousGalerkinSpace
163 
164  template< class FunctionSpace, class GridPart, int polOrder, class Storage >
165  struct hasStaticPolynomialOrder< hpDG::OrthogonalDiscontinuousGalerkinSpace< FunctionSpace, GridPart, polOrder, Storage > >
166  {
167  static const bool v = true;
168  static const int order = polOrder;
169  };
170 
171  template< class FunctionSpace, class GridPart, int polOrder, class Storage >
172  struct isLocalized< hpDG::OrthogonalDiscontinuousGalerkinSpace< FunctionSpace, GridPart, polOrder, Storage > >
173  {
174  static const bool v = true;
175  };
176 
177  template< class FunctionSpace, class GridPart, int polOrder, class Storage >
178  struct isAdaptive< hpDG::OrthogonalDiscontinuousGalerkinSpace< FunctionSpace, GridPart, polOrder, Storage > >
179  {
180  static const bool v = true;
181  };
182 
183  template< class FunctionSpace, class GridPart, int polOrder, class Storage >
184  struct viewThreadSafe< hpDG::OrthogonalDiscontinuousGalerkinSpace< FunctionSpace, GridPart, polOrder, Storage > >
185  {
186  static const bool v = true;
187  };
188 
189  template< class FunctionSpace, class GridPart, int polOrder, class Storage >
190  struct isHierarchic< hpDG::OrthogonalDiscontinuousGalerkinSpace< FunctionSpace, GridPart, polOrder, Storage > >
191  {
192  static const bool v = true;
193  };
194 
195  } // namespace Capabilities
196 
197 
198  } // namespace Fem
199 
200 } // namespace Dune
201 
202 #endif // #ifndef DUNE_FEM_HPDG_SPACE_DISCONTINUOUSGALERKIN_ORTHOGONAL_HH
Definition: bindguard.hh:11
Abstract class representing a function.
Definition: common/function.hh:50
Local Mass Matrix inversion implementation, select the correct method in your implementation.
Definition: localmassmatrix.hh:34
Local Mass Matrix for arbitrary spaces.
Definition: localmassmatrix.hh:909
DG Local Mass Matrix for arbitrary spaces.
Definition: localmassmatrix.hh:929
specialize with true if polynomial order fixed and compile time static
Definition: space/common/capabilities.hh:37
static const bool v
Definition: space/common/capabilities.hh:38
static const int order
Definition: space/common/capabilities.hh:39
specialize with true if the space is localized, * i.e., the basis function set is based on a shape fu...
Definition: space/common/capabilities.hh:68
static const bool v
Definition: space/common/capabilities.hh:69
specialize with true if space can be used with AdaptiveDiscreteFunction
Definition: space/common/capabilities.hh:81
static const bool v
Definition: space/common/capabilities.hh:82
specialize with true if the space implementation is thread safe, while it is not modified
Definition: space/common/capabilities.hh:108
static const bool v
Definition: space/common/capabilities.hh:109
specialize with true if for a space the basis functions are sorted by the polynomial order,...
Definition: space/common/capabilities.hh:120
static const bool v
Definition: space/common/capabilities.hh:121
Default communication handler for discrete functions.
Definition: defaultcommhandler.hh:29
@ dimRange
dimension of range vector space
Definition: functionspaceinterface.hh:48
Definition: common/localrestrictprolong.hh:16
Definition: discontinuousgalerkin/localrestrictprolong.hh:31
Implementation of an -adaptive discrete function space using orthogonal polynomials.
Definition: hpdg/orthogonal.hh:89
typename BaseType::BasisFunctionSetsType BasisFunctionSetsType
Definition: hpdg/orthogonal.hh:98
typename BaseType::GridPartType GridPartType
Definition: hpdg/orthogonal.hh:96
typename BaseType::EntityType EntityType
Definition: hpdg/orthogonal.hh:97
OrthogonalDiscontinuousGalerkinSpace(GridPartType &gridPart, const Dune::InterfaceType interface=Dune::InteriorBorder_All_Interface, const Dune::CommunicationDirection direction=Dune::ForwardCommunication)
Definition: hpdg/orthogonal.hh:100
OrthogonalDiscontinuousGalerkinSpace(GridPartType &gridPart, Function function, const Dune::InterfaceType interface=Dune::InteriorBorder_All_Interface, const Dune::CommunicationDirection direction=Dune::ForwardCommunication)
Definition: hpdg/orthogonal.hh:110
static const int polynomialOrder
Definition: hpdg/orthogonal.hh:94
Generic implementation of a -adaptive discontinuous finite element space.
Definition: hpdg/space.hh:46
typename Traits::BasisFunctionSetsType BasisFunctionSetsType
basis function sets type
Definition: hpdg/space.hh:56