dune-functions  2.6-dev
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, 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  quadrilateralOffset_ = dofsPerTriangle * gridView_.size(Dune::GeometryTypes::triangle);
91  break;
92  }
93  case 3:
94  {
95  prismOffset_ = dofsPerTetrahedron * gridView_.size(Dune::GeometryTypes::tetrahedron);
96 
97  hexahedronOffset_ = prismOffset_ + dofsPerPrism * gridView_.size(Dune::GeometryTypes::prism);
98 
99  pyramidOffset_ = hexahedronOffset_ + dofsPerHexahedron * gridView_.size(Dune::GeometryTypes::hexahedron);
100  break;
101  }
102  }
103  }
104 
107  const GridView& gridView() const
108  {
109  return gridView_;
110  }
111 
112  void update(const GridView& gv)
113  {
114  gridView_ = gv;
115  }
116 
117  template<class TP>
118  Node<TP> node(const TP& tp) const
119  {
120  return Node<TP>{tp};
121  }
122 
123  template<class TP>
125  {
126  return IndexSet<TP>{*this};
127  }
128 
129  size_type size() const
130  {
131  switch (dim)
132  {
133  case 1:
134  return dofsPerEdge*gridView_.size(0);
135  case 2:
136  {
137  return dofsPerTriangle*gridView_.size(Dune::GeometryTypes::triangle) + dofsPerQuad*gridView_.size(Dune::GeometryTypes::quadrilateral);
138  }
139  case 3:
140  {
141  return dofsPerTetrahedron*gridView_.size(Dune::GeometryTypes::tetrahedron)
142  + dofsPerPyramid*gridView_.size(Dune::GeometryTypes::pyramid)
143  + dofsPerPrism*gridView_.size(Dune::GeometryTypes::prism)
144  + dofsPerHexahedron*gridView_.size(Dune::GeometryTypes::hexahedron);
145  }
146  }
147  DUNE_THROW(Dune::NotImplemented, "No size method for " << dim << "d grids available yet!");
148  }
149 
151  size_type size(const SizePrefix prefix) const
152  {
153  assert(prefix.size() == 0 || prefix.size() == 1);
154  return (prefix.size() == 0) ? size() : 0;
155  }
156 
159  {
160  return size();
161  }
162 
164  {
165  return StaticPower<(k+1),GV::dimension>::power;
166  }
167 
168 //protected:
170 
173  size_t prismOffset_;
175 };
176 
177 
178 
179 template<typename GV, int k, class MI, class TP>
181 {
182  // Cannot be an enum -- otherwise the switch statement below produces compiler warnings
183  static const int dim = GV::dimension;
184 
185 public:
186 
187  using size_type = std::size_t;
188 
190  using MultiIndex = MI;
191 
193 
194  using Node = typename PreBasis::template Node<TP>;
195 
196  LagrangeDGNodeIndexSet(const PreBasis& preBasis) :
197  preBasis_(&preBasis)
198  {}
199 
205  void bind(const Node& node)
206  {
207  node_ = &node;
208  }
209 
212  void unbind()
213  {
214  node_ = nullptr;
215  }
216 
219  size_type size() const
220  {
221  return node_->finiteElement().size();
222  }
223 
225  template<typename It>
226  It indices(It it) const
227  {
228  const auto& gridIndexSet = preBasis_->gridView().indexSet();
229  const auto& element = node_->element();
230 
231  for (size_type i = 0, end = size() ; i < end ; ++i, ++it)
232  {
233  switch (dim)
234  {
235  case 1:
236  {
237  *it = {preBasis_->dofsPerEdge*gridIndexSet.subIndex(element,0,0) + i};
238  continue;
239  }
240  case 2:
241  {
242  if (element.type().isTriangle())
243  {
244  *it = {preBasis_->dofsPerTriangle*gridIndexSet.subIndex(element,0,0) + i};
245  continue;
246  }
247  else if (element.type().isQuadrilateral())
248  {
249  *it = { preBasis_->quadrilateralOffset_ + preBasis_->dofsPerQuad*gridIndexSet.subIndex(element,0,0) + i};
250  continue;
251  }
252  else
253  DUNE_THROW(Dune::NotImplemented, "2d elements have to be triangles or quadrilaterals");
254  }
255  case 3:
256  {
257  if (element.type().isTetrahedron())
258  {
259  *it = {preBasis_->dofsPerTetrahedron*gridIndexSet.subIndex(element,0,0) + i};
260  continue;
261  }
262  else if (element.type().isPrism())
263  {
264  *it = { preBasis_->prismOffset_ + preBasis_->dofsPerPrism*gridIndexSet.subIndex(element,0,0) + i};
265  continue;
266  }
267  else if (element.type().isHexahedron())
268  {
269  *it = { preBasis_->hexahedronOffset_ + preBasis_->dofsPerHexahedron*gridIndexSet.subIndex(element,0,0) + i};
270  continue;
271  }
272  else if (element.type().isPyramid())
273  {
274  *it = { preBasis_->pyramidOffset_ + preBasis_->dofsPerPyramid*gridIndexSet.subIndex(element,0,0) + i};
275  continue;
276  }
277  else
278  DUNE_THROW(Dune::NotImplemented, "3d elements have to be tetrahedrons, prisms, hexahedrons or pyramids");
279  }
280  }
281  DUNE_THROW(Dune::NotImplemented, "No index method for " << dim << "d grids available yet!");
282  }
283  return it;
284  }
285 
286 protected:
288 
289  const Node* node_;
290 };
291 
292 
293 
294 // *****************************************************************************
295 // This is the actual global basis implementation based on the reusable parts.
296 // *****************************************************************************
297 
305 template<typename GV, int k>
307 
308 
309 
310 } // end namespace Functions
311 } // end namespace Dune
312 
313 
314 #endif // DUNE_FUNCTIONS_FUNCTIONSPACEBASES_LAGRANGEDGBASIS_HH
Dune::Functions::LagrangeDGNodeIndexSet
Definition: lagrangedgbasis.hh:38
Dune::Functions::LagrangeDGPreBasis::prismOffset_
size_t prismOffset_
Definition: lagrangedgbasis.hh:173
Dune::Functions::LagrangeDGPreBasis::size
size_type size() const
Definition: lagrangedgbasis.hh:129
Dune
Definition: polynomial.hh:7
Dune::Functions::LagrangeDGNodeIndexSet::size
size_type size() const
Size of subtree rooted in this node (element-local)
Definition: lagrangedgbasis.hh:219
Dune::Functions::LagrangeDGPreBasis::hexahedronOffset_
size_t hexahedronOffset_
Definition: lagrangedgbasis.hh:174
pqknodalbasis.hh
Dune::Functions::LagrangeDGNodeIndexSet::LagrangeDGNodeIndexSet
LagrangeDGNodeIndexSet(const PreBasis &preBasis)
Definition: lagrangedgbasis.hh:196
flatmultiindex.hh
Dune::Functions::LagrangeDGNodeIndexSet::unbind
void unbind()
Unbind the view.
Definition: lagrangedgbasis.hh:212
Dune::Functions::LagrangeDGPreBasis::size_type
std::size_t size_type
Definition: lagrangedgbasis.hh:50
Dune::Functions::LagrangeDGNodeIndexSet::size_type
std::size_t size_type
Definition: lagrangedgbasis.hh:187
Dune::Functions::LagrangeDGPreBasis::gridView
const GridView & gridView() const
Obtain the grid view that the basis is defined on.
Definition: lagrangedgbasis.hh:107
Dune::Functions::LagrangeDGNodeIndexSet::node_
const Node * node_
Definition: lagrangedgbasis.hh:289
Dune::Functions::LagrangeDGPreBasis::dofsPerQuad
const static int dofsPerQuad
Definition: lagrangedgbasis.hh:56
Dune::Functions::LagrangeDGNodeIndexSet::MultiIndex
MI MultiIndex
Type used for global numbering of the basis vectors.
Definition: lagrangedgbasis.hh:190
Dune::Functions::LagrangeDGPreBasis::dofsPerEdge
const static int dofsPerEdge
Definition: lagrangedgbasis.hh:54
Dune::Functions::LagrangeDGPreBasis::size
size_type size(const SizePrefix prefix) const
Return number possible values for next position in multi index.
Definition: lagrangedgbasis.hh:151
Dune::Functions::LagrangeDGPreBasis::dofsPerPrism
const static int dofsPerPrism
Definition: lagrangedgbasis.hh:58
nodes.hh
Dune::Functions::LagrangeDGPreBasis::MultiIndex
MI MultiIndex
Type used for global numbering of the basis vectors.
Definition: lagrangedgbasis.hh:70
defaultglobalbasis.hh
Dune::Functions::LagrangeDGPreBasis::quadrilateralOffset_
size_t quadrilateralOffset_
Definition: lagrangedgbasis.hh:171
Dune::Functions::LagrangeDGPreBasis::SizePrefix
Dune::ReservedVector< size_type, 2 > SizePrefix
Definition: lagrangedgbasis.hh:72
Dune::Functions::LagrangeDGPreBasis::maxNodeSize
size_type maxNodeSize() const
Definition: lagrangedgbasis.hh:163
Dune::Functions::LagrangeDGNodeIndexSet::preBasis_
const PreBasis * preBasis_
Definition: lagrangedgbasis.hh:287
Dune::Functions::LagrangeDGPreBasis::initializeIndices
void initializeIndices()
Definition: lagrangedgbasis.hh:80
Dune::Functions::PQkNode
Definition: pqknodalbasis.hh:36
Dune::Functions::LagrangeDGPreBasis
Definition: lagrangedgbasis.hh:42
Dune::Functions::LagrangeDGPreBasis::pyramidOffset_
size_t pyramidOffset_
Definition: lagrangedgbasis.hh:172
Dune::Functions::LagrangeDGPreBasis::GridView
GV GridView
The grid view that the FE space is defined on.
Definition: lagrangedgbasis.hh:49
Dune::Functions::LagrangeDGNodeIndexSet::Node
typename PreBasis::template Node< TP > Node
Definition: lagrangedgbasis.hh:194
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:226
Dune::Functions::LagrangeDGPreBasis::update
void update(const GridView &gv)
Definition: lagrangedgbasis.hh:112
Dune::Functions::LagrangeDGPreBasis::dimension
size_type dimension() const
Definition: lagrangedgbasis.hh:158
Dune::Functions::LagrangeDGPreBasis::LagrangeDGPreBasis
LagrangeDGPreBasis(const GridView &gv)
Constructor for a given grid view object.
Definition: lagrangedgbasis.hh:75
Dune::Functions::LagrangeDGPreBasis::dofsPerTriangle
const static int dofsPerTriangle
Definition: lagrangedgbasis.hh:55
Dune::Functions::LagrangeDGPreBasis::indexSet
IndexSet< TP > indexSet() const
Definition: lagrangedgbasis.hh:124
Dune::Functions::LagrangeDGPreBasis::dofsPerPyramid
const static int dofsPerPyramid
Definition: lagrangedgbasis.hh:60
Dune::Functions::LagrangeDGNodeIndexSet::bind
void bind(const Node &node)
Bind the view to a grid element.
Definition: lagrangedgbasis.hh:205
Dune::Functions::LagrangeDGPreBasis::dofsPerTetrahedron
const static int dofsPerTetrahedron
Definition: lagrangedgbasis.hh:57
Dune::Functions::LagrangeDGPreBasis::node
Node< TP > node(const TP &tp) const
Definition: lagrangedgbasis.hh:118
Dune::Functions::LagrangeDGPreBasis::gridView_
GridView gridView_
Definition: lagrangedgbasis.hh:169
Dune::Functions::BasisBuilder::power
auto power(ChildPreBasisFactory &&childPreBasisFactory, const IndexMergingStrategy &ims)
Create a pre-basis factory that can build a PowerPreBasis.
Definition: powerbasis.hh:493
Dune::Functions::DefaultGlobalBasis
Global basis for given pre-basis.
Definition: defaultglobalbasis.hh:42
Dune::Functions::LagrangeDGPreBasis::dofsPerHexahedron
const static int dofsPerHexahedron
Definition: lagrangedgbasis.hh:59