dune-functions  2.5.1
lagrangedgbasis.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_FUNCTIONS_FUNCTIONSPACEBASES_LAGRANGEDGBASIS_HH
4 #define DUNE_FUNCTIONS_FUNCTIONSPACEBASES_LAGRANGEDGBASIS_HH
5 
6 #include <array>
7 #include <dune/common/exceptions.hh>
8 
13 
14 
15 
16 
17 namespace Dune {
18 namespace Functions {
19 
20 
21 
22 // *****************************************************************************
23 // This is the reusable part of the basis. It contains
24 //
25 // LagrangeDGNodeFactory
26 // LagrangeDGNodeIndexSet
27 // LagrangeDGNode
28 //
29 // The factory allows to create the others and is the owner of possible shared
30 // state. These three components do _not_ depend on the global basis or index
31 // set and can be used without a global basis.
32 // *****************************************************************************
33 
34 template<typename GV, int k, typename TP>
36 
37 template<typename GV, int k, class MI, class TP>
39 
40 
41 template<typename GV, int k, class MI>
43 {
44  static const int dim = GV::dimension;
45 
46 public:
47 
49  using GridView = GV;
50  using size_type = std::size_t;
51 
52 
53  // Precompute the number of dofs per entity type
54  const static int dofsPerEdge = k+1;
55  const static int dofsPerTriangle = (k+1)*(k+2)/2;
56  const static int dofsPerQuad = (k+1)*(k+1);
57  const static int dofsPerTetrahedron = (k+1)*(k+2)*(k+3)/6;
58  const static int dofsPerPrism = (k+1)*(k+1)*(k+2)/2;
59  const static int dofsPerHexahedron = (k+1)*(k+1)*(k+1);
60  const static int dofsPerPyramid = (k+1)*(k+2)*(2*k+3)/6;
61 
62 
63  template<class TP>
65 
66  template<class TP>
68 
70  using MultiIndex = MI;
71 
72  using SizePrefix = Dune::ReservedVector<size_type, 2>;
73 
76  gridView_(gv)
77  {}
78 
79 
81  {
82  switch (dim)
83  {
84  case 1:
85  {
86  break;
87  }
88  case 2:
89  {
90  GeometryType triangle;
91  triangle.makeTriangle();
92  quadrilateralOffset_ = dofsPerTriangle * gridView_.size(triangle);
93  break;
94  }
95  case 3:
96  {
97  GeometryType tetrahedron;
98  tetrahedron.makeSimplex(3);
99  prismOffset_ = dofsPerTetrahedron * gridView_.size(tetrahedron);
100 
101  GeometryType prism;
102  prism.makePrism();
103  hexahedronOffset_ = prismOffset_ + dofsPerPrism * gridView_.size(prism);
104 
105  GeometryType hexahedron;
106  hexahedron.makeCube(3);
107  pyramidOffset_ = hexahedronOffset_ + dofsPerHexahedron * gridView_.size(hexahedron);
108  break;
109  }
110  }
111  }
112 
115  const GridView& gridView() const
116  {
117  return gridView_;
118  }
119 
120  void update(const GridView& gv)
121  {
122  gridView_ = gv;
123  }
124 
125  template<class TP>
126  Node<TP> node(const TP& tp) const
127  {
128  return Node<TP>{tp};
129  }
130 
131  template<class TP>
133  {
134  return IndexSet<TP>{*this};
135  }
136 
137  size_type size() const
138  {
139  switch (dim)
140  {
141  case 1:
142  return dofsPerEdge*gridView_.size(0);
143  case 2:
144  {
145  GeometryType triangle, quad;
146  triangle.makeTriangle();
147  quad.makeQuadrilateral();
148  return dofsPerTriangle*gridView_.size(triangle) + dofsPerQuad*gridView_.size(quad);
149  }
150  case 3:
151  {
152  GeometryType tetrahedron, pyramid, prism, hexahedron;
153  tetrahedron.makeTetrahedron();
154  pyramid.makePyramid();
155  prism.makePrism();
156  hexahedron.makeCube(3);
157  return dofsPerTetrahedron*gridView_.size(tetrahedron) + dofsPerPyramid*gridView_.size(pyramid)
158  + dofsPerPrism*gridView_.size(prism) + dofsPerHexahedron*gridView_.size(hexahedron);
159  }
160  }
161  DUNE_THROW(Dune::NotImplemented, "No size method for " << dim << "d grids available yet!");
162  }
163 
165  size_type size(const SizePrefix prefix) const
166  {
167  if (prefix.size() == 0)
168  return size();
169  assert(false);
170  }
171 
174  {
175  return size();
176  }
177 
179  {
181  }
182 
183 //protected:
185 
188  size_t prismOffset_;
190 };
191 
192 
193 
194 template<typename GV, int k, class MI, class TP>
196 {
197  // Cannot be an enum -- otherwise the switch statement below produces compiler warnings
198  static const int dim = GV::dimension;
199 
200 public:
201 
202  using size_type = std::size_t;
203 
205  using MultiIndex = MI;
206 
208 
209  using Node = typename NodeFactory::template Node<TP>;
210 
211  LagrangeDGNodeIndexSet(const NodeFactory& nodeFactory) :
212  nodeFactory_(&nodeFactory)
213  {}
214 
220  void bind(const Node& node)
221  {
222  node_ = &node;
223  }
224 
227  void unbind()
228  {
229  node_ = nullptr;
230  }
231 
234  size_type size() const
235  {
236  return node_->finiteElement().size();
237  }
238 
241  {
242  const auto& gridIndexSet = nodeFactory_->gridView().indexSet();
243  const auto& element = node_->element();
244 
245  switch (dim)
246  {
247  case 1:
248  {
249  return {nodeFactory_->dofsPerEdge*gridIndexSet.subIndex(element,0,0) + i};
250  }
251  case 2:
252  {
253  if (element.type().isTriangle())
254  {
255  return {nodeFactory_->dofsPerTriangle*gridIndexSet.subIndex(element,0,0) + i};
256  }
257  else if (element.type().isQuadrilateral())
258  {
259  return { nodeFactory_->quadrilateralOffset_ + nodeFactory_->dofsPerQuad*gridIndexSet.subIndex(element,0,0) + i};
260  }
261  else
262  DUNE_THROW(Dune::NotImplemented, "2d elements have to be triangles or quadrilaterals");
263  }
264  case 3:
265  {
266  if (element.type().isTetrahedron())
267  {
268  return {nodeFactory_->dofsPerTetrahedron*gridIndexSet.subIndex(element,0,0) + i};
269  }
270  else if (element.type().isPrism())
271  {
272  return { nodeFactory_->prismOffset_ + nodeFactory_->dofsPerPrism*gridIndexSet.subIndex(element,0,0) + i};
273  }
274  else if (element.type().isHexahedron())
275  {
276  return { nodeFactory_->hexahedronOffset_ + nodeFactory_->dofsPerHexahedron*gridIndexSet.subIndex(element,0,0) + i};
277  }
278  else if (element.type().isPyramid())
279  {
280  return { nodeFactory_->pyramidOffset_ + nodeFactory_->dofsPerPyramid*gridIndexSet.subIndex(element,0,0) + i};
281  }
282  else
283  DUNE_THROW(Dune::NotImplemented, "3d elements have to be tetrahedrons, prisms, hexahedrons or pyramids");
284  }
285  }
286  DUNE_THROW(Dune::NotImplemented, "No index method for " << dim << "d grids available yet!");
287  }
288 
289 protected:
291 
292  const Node* node_;
293 };
294 
295 
296 
297 // *****************************************************************************
298 // This is the actual global basis implementation based on the reusable parts.
299 // *****************************************************************************
300 
308 template<typename GV, int k>
310 
311 
312 
313 } // end namespace Functions
314 } // end namespace Dune
315 
316 
317 #endif // DUNE_FUNCTIONS_FUNCTIONSPACEBASES_LAGRANGEDGBASIS_HH
size_type maxNodeSize() const
Definition: lagrangedgbasis.hh:178
void initializeIndices()
Definition: lagrangedgbasis.hh:80
GridView gridView_
Definition: lagrangedgbasis.hh:184
Definition: lagrangedgbasis.hh:42
static const int dofsPerPrism
Definition: lagrangedgbasis.hh:58
typename NodeFactory::template Node< TP > Node
Definition: lagrangedgbasis.hh:209
static const int dofsPerTriangle
Definition: lagrangedgbasis.hh:55
std::size_t size_type
Definition: lagrangedgbasis.hh:202
size_t prismOffset_
Definition: lagrangedgbasis.hh:188
Global basis for given node factory.
Definition: defaultglobalbasis.hh:42
Definition: polynomial.hh:7
const NodeFactory * nodeFactory_
Definition: lagrangedgbasis.hh:290
size_t pyramidOffset_
Definition: lagrangedgbasis.hh:187
static const int dofsPerEdge
Definition: lagrangedgbasis.hh:54
void bind(const Node &node)
Bind the view to a grid element.
Definition: lagrangedgbasis.hh:220
IndexSet< TP > indexSet() const
Definition: lagrangedgbasis.hh:132
const Node * node_
Definition: lagrangedgbasis.hh:292
size_type size() const
Size of subtree rooted in this node (element-local)
Definition: lagrangedgbasis.hh:234
size_type size(const SizePrefix prefix) const
Return number possible values for next position in multi index.
Definition: lagrangedgbasis.hh:165
Node< TP > node(const TP &tp) const
Definition: lagrangedgbasis.hh:126
void update(const GridView &gv)
Definition: lagrangedgbasis.hh:120
MI MultiIndex
Type used for global numbering of the basis vectors.
Definition: lagrangedgbasis.hh:70
std::size_t size_type
Definition: lagrangedgbasis.hh:50
static const int dofsPerTetrahedron
Definition: lagrangedgbasis.hh:57
Definition: pqknodalbasis.hh:36
const GridView & gridView() const
Obtain the grid view that the basis is defined on.
Definition: lagrangedgbasis.hh:115
size_type dimension() const
Definition: lagrangedgbasis.hh:173
static const int dofsPerHexahedron
Definition: lagrangedgbasis.hh:59
void unbind()
Unbind the view.
Definition: lagrangedgbasis.hh:227
static const int dofsPerQuad
Definition: lagrangedgbasis.hh:56
LagrangeDGNodeIndexSet(const NodeFactory &nodeFactory)
Definition: lagrangedgbasis.hh:211
Imp::PowerNodeFactoryBuilder< k, IndexMergingStrategy, SubFactoryTag > power(SubFactoryTag &&tag, const IndexMergingStrategy &ims)
Create a factory builder that can build a PowerNodeFactory.
Definition: powerbasis.hh:424
Definition: lagrangedgbasis.hh:38
size_type size() const
Definition: lagrangedgbasis.hh:137
MI MultiIndex
Type used for global numbering of the basis vectors.
Definition: lagrangedgbasis.hh:205
MultiIndex index(size_type i) const
Maps from subtree index set [0..size-1] to a globally unique multi index in global basis...
Definition: lagrangedgbasis.hh:240
static const int dofsPerPyramid
Definition: lagrangedgbasis.hh:60
Dune::ReservedVector< size_type, 2 > SizePrefix
Definition: lagrangedgbasis.hh:72
size_t hexahedronOffset_
Definition: lagrangedgbasis.hh:189
GV GridView
The grid view that the FE space is defined on.
Definition: lagrangedgbasis.hh:49
size_t quadrilateralOffset_
Definition: lagrangedgbasis.hh:186
LagrangeDGNodeFactory(const GridView &gv)
Constructor for a given grid view object.
Definition: lagrangedgbasis.hh:75