dune-functions  2.7.0
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 // LagrangeDGPreBasis
26 // LagrangeDGNodeIndexSet
27 // LagrangeDGNode
28 //
29 // The pre-basis 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>
36 
37 template<typename GV, int k, class MI>
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 
64 
66 
68  using MultiIndex = MI;
69 
70  using SizePrefix = Dune::ReservedVector<size_type, 1>;
71 
74  gridView_(gv)
75  {}
76 
77 
79  {
80  switch (dim)
81  {
82  case 1:
83  {
84  break;
85  }
86  case 2:
87  {
88  quadrilateralOffset_ = dofsPerTriangle * gridView_.size(Dune::GeometryTypes::triangle);
89  break;
90  }
91  case 3:
92  {
93  prismOffset_ = dofsPerTetrahedron * gridView_.size(Dune::GeometryTypes::tetrahedron);
94 
95  hexahedronOffset_ = prismOffset_ + dofsPerPrism * gridView_.size(Dune::GeometryTypes::prism);
96 
97  pyramidOffset_ = hexahedronOffset_ + dofsPerHexahedron * gridView_.size(Dune::GeometryTypes::hexahedron);
98  break;
99  }
100  }
101  }
102 
105  const GridView& gridView() const
106  {
107  return gridView_;
108  }
109 
110  void update(const GridView& gv)
111  {
112  gridView_ = gv;
113  }
114 
118  Node makeNode() const
119  {
120  return Node{};
121  }
122 
130  {
131  return IndexSet{*this};
132  }
133 
134  size_type size() const
135  {
136  switch (dim)
137  {
138  case 1:
139  return dofsPerEdge*gridView_.size(0);
140  case 2:
141  {
142  return dofsPerTriangle*gridView_.size(Dune::GeometryTypes::triangle) + dofsPerQuad*gridView_.size(Dune::GeometryTypes::quadrilateral);
143  }
144  case 3:
145  {
146  return dofsPerTetrahedron*gridView_.size(Dune::GeometryTypes::tetrahedron)
147  + dofsPerPyramid*gridView_.size(Dune::GeometryTypes::pyramid)
148  + dofsPerPrism*gridView_.size(Dune::GeometryTypes::prism)
149  + dofsPerHexahedron*gridView_.size(Dune::GeometryTypes::hexahedron);
150  }
151  }
152  DUNE_THROW(Dune::NotImplemented, "No size method for " << dim << "d grids available yet!");
153  }
154 
156  size_type size(const SizePrefix prefix) const
157  {
158  assert(prefix.size() == 0 || prefix.size() == 1);
159  return (prefix.size() == 0) ? size() : 0;
160  }
161 
164  {
165  return size();
166  }
167 
169  {
170  return StaticPower<(k+1),GV::dimension>::power;
171  }
172 
173 //protected:
175 
178  size_t prismOffset_;
180 };
181 
182 
183 
184 template<typename GV, int k, class MI>
186 {
187  // Cannot be an enum -- otherwise the switch statement below produces compiler warnings
188  static const int dim = GV::dimension;
189 
190 public:
191 
192  using size_type = std::size_t;
193 
195  using MultiIndex = MI;
196 
198 
200 
201  LagrangeDGNodeIndexSet(const PreBasis& preBasis) :
202  preBasis_(&preBasis)
203  {}
204 
210  void bind(const Node& node)
211  {
212  node_ = &node;
213  }
214 
217  void unbind()
218  {
219  node_ = nullptr;
220  }
221 
224  size_type size() const
225  {
226  return node_->finiteElement().size();
227  }
228 
230  template<typename It>
231  It indices(It it) const
232  {
233  const auto& gridIndexSet = preBasis_->gridView().indexSet();
234  const auto& element = node_->element();
235 
236  for (size_type i = 0, end = size() ; i < end ; ++i, ++it)
237  {
238  switch (dim)
239  {
240  case 1:
241  {
242  *it = {preBasis_->dofsPerEdge*gridIndexSet.subIndex(element,0,0) + i};
243  continue;
244  }
245  case 2:
246  {
247  if (element.type().isTriangle())
248  {
249  *it = {preBasis_->dofsPerTriangle*gridIndexSet.subIndex(element,0,0) + i};
250  continue;
251  }
252  else if (element.type().isQuadrilateral())
253  {
254  *it = { preBasis_->quadrilateralOffset_ + preBasis_->dofsPerQuad*gridIndexSet.subIndex(element,0,0) + i};
255  continue;
256  }
257  else
258  DUNE_THROW(Dune::NotImplemented, "2d elements have to be triangles or quadrilaterals");
259  }
260  case 3:
261  {
262  if (element.type().isTetrahedron())
263  {
264  *it = {preBasis_->dofsPerTetrahedron*gridIndexSet.subIndex(element,0,0) + i};
265  continue;
266  }
267  else if (element.type().isPrism())
268  {
269  *it = { preBasis_->prismOffset_ + preBasis_->dofsPerPrism*gridIndexSet.subIndex(element,0,0) + i};
270  continue;
271  }
272  else if (element.type().isHexahedron())
273  {
274  *it = { preBasis_->hexahedronOffset_ + preBasis_->dofsPerHexahedron*gridIndexSet.subIndex(element,0,0) + i};
275  continue;
276  }
277  else if (element.type().isPyramid())
278  {
279  *it = { preBasis_->pyramidOffset_ + preBasis_->dofsPerPyramid*gridIndexSet.subIndex(element,0,0) + i};
280  continue;
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  return it;
289  }
290 
291 protected:
293 
294  const Node* node_;
295 };
296 
297 
298 
299 
300 namespace BasisFactory {
301 
302 namespace Imp {
303 
304 template<std::size_t k>
305 class LagrangeDGPreBasisFactory
306 {
307 public:
308  static const std::size_t requiredMultiIndexSize = 1;
309 
310  template<class MultiIndex, class GridView>
311  auto makePreBasis(const GridView& gridView) const
312  {
314  }
315 
316 };
317 
318 } // end namespace BasisFactory::Imp
319 
320 
321 
329 template<std::size_t k>
331 {
332  return Imp::LagrangeDGPreBasisFactory<k>();
333 }
334 
335 } // end namespace BasisFactory
336 
337 
338 
339 // *****************************************************************************
340 // This is the actual global basis implementation based on the reusable parts.
341 // *****************************************************************************
342 
350 template<typename GV, int k>
352 
353 
354 
355 } // end namespace Functions
356 } // end namespace Dune
357 
358 
359 #endif // DUNE_FUNCTIONS_FUNCTIONSPACEBASES_LAGRANGEDGBASIS_HH
Dune::Functions::LagrangeDGPreBasis::pyramidOffset_
size_t pyramidOffset_
Definition: lagrangedgbasis.hh:177
Dune::Functions::LagrangeNode::element
const Element & element() const
Return current element, throw if unbound.
Definition: lagrangebasis.hh:351
Dune::Functions::LagrangeDGPreBasis::gridView
const GridView & gridView() const
Obtain the grid view that the basis is defined on.
Definition: lagrangedgbasis.hh:105
Dune::Functions::LagrangeDGPreBasis::maxNodeSize
size_type maxNodeSize() const
Definition: lagrangedgbasis.hh:168
Dune::Functions::LagrangeDGNodeIndexSet::size_type
std::size_t size_type
Definition: lagrangedgbasis.hh:192
Dune::Functions::LagrangeDGPreBasis::dofsPerQuad
const static int dofsPerQuad
Definition: lagrangedgbasis.hh:56
Dune::Functions::LagrangeDGPreBasis::dofsPerTetrahedron
const static int dofsPerTetrahedron
Definition: lagrangedgbasis.hh:57
Dune::Functions::LagrangeDGPreBasis::SizePrefix
Dune::ReservedVector< size_type, 1 > SizePrefix
Definition: lagrangedgbasis.hh:70
Dune::Functions::LagrangeDGNodeIndexSet::indices
It indices(It it) const
Maps from subtree index set [0..size-1] to a globally unique multi index in global basis.
Definition: lagrangedgbasis.hh:231
Dune::Functions::LagrangeDGNodeIndexSet::unbind
void unbind()
Unbind the view.
Definition: lagrangedgbasis.hh:217
Dune::Functions::LagrangeDGPreBasis::makeNode
Node makeNode() const
Create tree node.
Definition: lagrangedgbasis.hh:118
Dune::Functions::LagrangeDGPreBasis::dofsPerTriangle
const static int dofsPerTriangle
Definition: lagrangedgbasis.hh:55
Dune::Functions::LagrangeDGPreBasis::makeIndexSet
IndexSet makeIndexSet() const
Create tree node index set.
Definition: lagrangedgbasis.hh:129
Dune::Functions::LagrangeDGNodeIndexSet::MultiIndex
MI MultiIndex
Type used for global numbering of the basis vectors.
Definition: lagrangedgbasis.hh:195
Dune::Functions::LagrangeDGNodeIndexSet::bind
void bind(const Node &node)
Bind the view to a grid element.
Definition: lagrangedgbasis.hh:210
Dune::Functions::LagrangeDGPreBasis::dimension
size_type dimension() const
Definition: lagrangedgbasis.hh:163
Dune::Functions::LagrangeDGNodeIndexSet::size
size_type size() const
Size of subtree rooted in this node (element-local)
Definition: lagrangedgbasis.hh:224
lagrangebasis.hh
Dune::Functions::LagrangeDGPreBasis::hexahedronOffset_
size_t hexahedronOffset_
Definition: lagrangedgbasis.hh:179
Dune::Functions::LagrangeDGPreBasis::dofsPerPrism
const static int dofsPerPrism
Definition: lagrangedgbasis.hh:58
Dune::Functions::LagrangeDGPreBasis::dofsPerEdge
const static int dofsPerEdge
Definition: lagrangedgbasis.hh:54
Dune::Functions::LagrangeDGPreBasis
Definition: lagrangedgbasis.hh:42
Dune::Functions::DefaultGlobalBasis
Global basis for given pre-basis.
Definition: defaultglobalbasis.hh:46
Dune::Functions::LagrangeDGPreBasis::prismOffset_
size_t prismOffset_
Definition: lagrangedgbasis.hh:178
Dune::Functions::LagrangeDGPreBasis::dofsPerPyramid
const static int dofsPerPyramid
Definition: lagrangedgbasis.hh:60
Dune::Functions::LagrangeDGPreBasis::initializeIndices
void initializeIndices()
Definition: lagrangedgbasis.hh:78
Dune::Functions::LagrangeDGPreBasis::quadrilateralOffset_
size_t quadrilateralOffset_
Definition: lagrangedgbasis.hh:176
Dune::Functions::LagrangeDGNodeIndexSet
Definition: lagrangedgbasis.hh:38
Dune::Functions::LagrangeDGPreBasis::size
size_type size() const
Definition: lagrangedgbasis.hh:134
nodes.hh
Dune::Functions::BasisFactory::lagrangeDG
auto lagrangeDG()
Create a pre-basis factory that can create a LagrangeDG pre-basis.
Definition: lagrangedgbasis.hh:330
defaultglobalbasis.hh
Dune
Definition: polynomial.hh:10
Dune::Functions::LagrangeNode
Definition: lagrangebasis.hh:33
Dune::Functions::LagrangeDGPreBasis::LagrangeDGPreBasis
LagrangeDGPreBasis(const GridView &gv)
Constructor for a given grid view object.
Definition: lagrangedgbasis.hh:73
Dune::Functions::LagrangeDGNodeIndexSet::preBasis_
const PreBasis * preBasis_
Definition: lagrangedgbasis.hh:292
Dune::Functions::LagrangeDGPreBasis::update
void update(const GridView &gv)
Definition: lagrangedgbasis.hh:110
Dune::Functions::LagrangeDGPreBasis::GridView
GV GridView
The grid view that the FE space is defined on.
Definition: lagrangedgbasis.hh:49
Dune::Functions::LagrangeDGPreBasis::MultiIndex
MI MultiIndex
Type used for global numbering of the basis vectors.
Definition: lagrangedgbasis.hh:68
Dune::Functions::LagrangeDGPreBasis::size_type
std::size_t size_type
Definition: lagrangedgbasis.hh:50
Dune::Functions::LagrangeDGPreBasis::gridView_
GridView gridView_
Definition: lagrangedgbasis.hh:174
Dune::Functions::BasisFactory::power
auto power(ChildPreBasisFactory &&childPreBasisFactory, const IndexMergingStrategy &ims)
Create a pre-basis factory that can build a PowerPreBasis.
Definition: powerbasis.hh:483
Dune::Functions::LagrangeDGNodeIndexSet::LagrangeDGNodeIndexSet
LagrangeDGNodeIndexSet(const PreBasis &preBasis)
Definition: lagrangedgbasis.hh:201
flatmultiindex.hh
Dune::Functions::LagrangeDGPreBasis::size
size_type size(const SizePrefix prefix) const
Return number possible values for next position in multi index.
Definition: lagrangedgbasis.hh:156
Dune::Functions::LagrangeDGPreBasis::dofsPerHexahedron
const static int dofsPerHexahedron
Definition: lagrangedgbasis.hh:59
Dune::Functions::LagrangeNode::finiteElement
const FiniteElement & finiteElement() const
Return the LocalFiniteElement for the element we are bound to.
Definition: lagrangebasis.hh:360
Dune::Functions::LagrangeDGNodeIndexSet::node_
const Node * node_
Definition: lagrangedgbasis.hh:294