Dune Core Modules (2.7.1)

quadraturerules.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 
4 #ifndef DUNE_GEOMETRY_QUADRATURERULES_HH
5 #define DUNE_GEOMETRY_QUADRATURERULES_HH
6 
7 #include <algorithm>
8 #include <iostream>
9 #include <limits>
10 #include <mutex>
11 #include <utility>
12 #include <vector>
13 
14 #include <dune/common/fvector.hh>
17 #include <dune/common/stdthread.hh>
19 
20 #include <dune/geometry/type.hh>
22 
28 namespace Dune {
29 
35 
41  template<typename ct, int dim>
43  public:
45  enum { dimension = dim };
46 
48  typedef ct Field;
49 
52 
54  QuadraturePoint (const Vector& x, ct w) : local(x)
55  {
56  weight_ = w;
57  }
58 
60  const Vector& position () const
61  {
62  return local;
63  }
64 
66  const ct &weight () const
67  {
68  return weight_;
69  }
70 
71  protected:
73  ct weight_;
74  };
75 
79  namespace QuadratureType {
80  enum Enum {
89 
96 
103 
116  GaussLobatto = 4,
117  size
118  };
119  }
120 
124  template<typename ct, int dim>
125  class QuadratureRule : public std::vector<QuadraturePoint<ct,dim> >
126  {
127  public:
133  QuadratureRule() : delivered_order(-1) {}
134 
135  protected:
137  QuadratureRule(GeometryType t) : geometry_type(t), delivered_order(-1) {}
138 
140  QuadratureRule(GeometryType t, int order) : geometry_type(t), delivered_order(order) {}
141  public:
143  enum { d=dim };
144 
146  typedef ct CoordType;
147 
149  virtual int order () const { return delivered_order; }
150 
152  virtual GeometryType type () const { return geometry_type; }
153  virtual ~QuadratureRule(){}
154 
157  typedef typename std::vector<QuadraturePoint<ct,dim> >::const_iterator iterator;
158 
159  protected:
160  GeometryType geometry_type;
161  int delivered_order;
162  };
163 
164  // Forward declaration of the factory class,
165  // needed internally by the QuadratureRules container class.
166  template<typename ctype, int dim> class QuadratureRuleFactory;
167 
171  template<typename ctype, int dim>
173 
178  static void initQuadratureRule(QuadratureRule *qr, QuadratureType::Enum qt,
179  const GeometryType &t, int p)
180  {
182  }
183 
184  typedef std::vector<std::pair<std::once_flag, QuadratureRule> >
185  QuadratureOrderVector; // indexed by quadrature order
188  static void initQuadratureOrderVector(QuadratureOrderVector *qov,
190  const GeometryType &t)
191  {
192  if(dim == 0)
193  // we only need one quadrature rule for points, not maxint
194  *qov = QuadratureOrderVector(1);
195  else
196  *qov = QuadratureOrderVector(QuadratureRuleFactory<ctype,dim>::maxOrder(t,qt)+1);
197  }
198 
199  typedef std::vector<std::pair<std::once_flag, QuadratureOrderVector> >
200  GeometryTypeVector; // indexed by geometry type
203  static void initGeometryTypeVector(GeometryTypeVector *gtv)
204  {
205  *gtv = GeometryTypeVector(LocalGeometryTypeIndex::size(dim));
206  }
207 
210  {
211  assert(t.dim()==dim);
212 
213  DUNE_ASSERT_CALL_ONCE();
214 
215  static std::vector<std::pair< // indexed by quadrature type
216  std::once_flag,
217  GeometryTypeVector
218  > > quadratureCache(QuadratureType::size);
219 
220  auto & quadratureTypeLevel = quadratureCache[qt];
221  std::call_once(quadratureTypeLevel.first, initGeometryTypeVector,
222  &quadratureTypeLevel.second);
223 
224  auto & geometryTypeLevel =
225  quadratureTypeLevel.second[LocalGeometryTypeIndex::index(t)];
226  std::call_once(geometryTypeLevel.first, initQuadratureOrderVector,
227  &geometryTypeLevel.second, qt, t);
228 
229  // we only have one quadrature rule for points
230  auto & quadratureOrderLevel = geometryTypeLevel.second[dim == 0 ? 0 : p];
231  std::call_once(quadratureOrderLevel.first, initQuadratureRule,
232  &quadratureOrderLevel.second, qt, t, p);
233 
234  return quadratureOrderLevel.second;
235  }
237  DUNE_EXPORT static QuadratureRules& instance()
238  {
239  static QuadratureRules instance;
240  return instance;
241  }
243  QuadratureRules () {}
244  public:
246  static unsigned
249  {
251  }
252 
255  {
256  return instance()._rule(t,p,qt);
257  }
258 
262  {
263  GeometryType gt(t,dim);
264  return instance()._rule(gt,p,qt);
265  }
267  };
268 
269 } // end namespace Dune
270 
271 #include "quadraturerules/pointquadrature.hh"
272 
273 namespace Dune {
274 
276  template<typename ct, bool fundamental = std::numeric_limits<ct>::is_specialized>
278  template<typename ct>
279  struct GaussQuadratureInitHelper<ct, true> {
280  static void init(int p,
281  std::vector< FieldVector<ct, 1> > & _points,
282  std::vector< ct > & _weight,
283  int & delivered_order);
284  };
285  template<typename ct>
286  struct GaussQuadratureInitHelper<ct, false> {
287  static void init(int p,
288  std::vector< FieldVector<ct, 1> > & _points,
289  std::vector< ct > & _weight,
290  int & delivered_order);
291  };
292 
294  template<typename ct>
296  public QuadratureRule<ct,1>
297  {
298  public:
299  // compile time parameters
300  enum { dim=1 };
301  enum { highest_order=61 };
302 
304  private:
305  friend class QuadratureRuleFactory<ct,dim>;
306  GaussQuadratureRule1D (int p)
308  {
310  std::vector< FieldVector<ct, dim> > _points;
311  std::vector< ct > _weight;
312 
314  (p, _points, _weight, this->delivered_order);
315 
316  assert(_points.size() == _weight.size());
317  for (size_t i = 0; i < _points.size(); i++)
318  this->push_back(QuadraturePoint<ct,dim>(_points[i], _weight[i]));
319  }
320  };
321 
324 
325 } // namespace Dune
326 
327 #define DUNE_INCLUDING_IMPLEMENTATION
328 #include "quadraturerules/gauss_imp.hh"
329 
330 namespace Dune {
331 
333  template<typename ct,
334  bool fundamental = std::numeric_limits<ct>::is_specialized>
336  template<typename ct>
337  struct Jacobi1QuadratureInitHelper<ct, true> {
338  static void init(int p,
339  std::vector< FieldVector<ct, 1> > & _points,
340  std::vector< ct > & _weight,
341  int & delivered_order);
342  };
343  template<typename ct>
344  struct Jacobi1QuadratureInitHelper<ct, false> {
345  static void init(int p,
346  std::vector< FieldVector<ct, 1> > & _points,
347  std::vector< ct > & _weight,
348  int & delivered_order);
349  };
350 
354  template<typename ct>
356  public QuadratureRule<ct,1>
357  {
358  public:
360  enum { dim=1 };
361 
363  enum { highest_order=61 };
364 
366  private:
367  friend class QuadratureRuleFactory<ct,dim>;
370  {
372  std::vector< FieldVector<ct, dim> > _points;
373  std::vector< ct > _weight;
374 
375  int deliveredOrder_;
376 
378  (p, _points, _weight, deliveredOrder_);
379  this->delivered_order = deliveredOrder_;
380  assert(_points.size() == _weight.size());
381  for (size_t i = 0; i < _points.size(); i++)
382  this->push_back(QuadraturePoint<ct,dim>(_points[i], _weight[i]));
383  }
384  };
385 
386 #ifndef DOXYGEN
389 #endif // !DOXYGEN
390 
391 } // namespace Dune
392 
393 #define DUNE_INCLUDING_IMPLEMENTATION
394 #include "quadraturerules/jacobi_1_0_imp.hh"
395 
396 namespace Dune {
397 
399  template<typename ct,
400  bool fundamental = std::numeric_limits<ct>::is_specialized>
402  template<typename ct>
403  struct Jacobi2QuadratureInitHelper<ct, true> {
404  static void init(int p,
405  std::vector< FieldVector<ct, 1> > & _points,
406  std::vector< ct > & _weight,
407  int & delivered_order);
408  };
409  template<typename ct>
410  struct Jacobi2QuadratureInitHelper<ct, false> {
411  static void init(int p,
412  std::vector< FieldVector<ct, 1> > & _points,
413  std::vector< ct > & _weight,
414  int & delivered_order);
415  };
416 
420  template<typename ct>
422  public QuadratureRule<ct,1>
423  {
424  public:
426  enum { dim=1 };
427 
429  enum { highest_order=61 };
430 
432  private:
433  friend class QuadratureRuleFactory<ct,dim>;
436  {
438  std::vector< FieldVector<ct, dim> > _points;
439  std::vector< ct > _weight;
440 
441  int deliveredOrder_;
442 
444  (p, _points, _weight, deliveredOrder_);
445 
446  this->delivered_order = deliveredOrder_;
447  assert(_points.size() == _weight.size());
448  for (size_t i = 0; i < _points.size(); i++)
449  this->push_back(QuadraturePoint<ct,dim>(_points[i], _weight[i]));
450  }
451  };
452 
453 #ifndef DOXYGEN
456 #endif // !DOXYGEN
457 
458 } // namespace Dune
459 
460 #define DUNE_INCLUDING_IMPLEMENTATION
461 #include "quadraturerules/jacobi_2_0_imp.hh"
462 
463 namespace Dune {
464 
466  template<typename ct,
467  bool fundamental = std::numeric_limits<ct>::is_specialized>
469  template<typename ct>
470  struct GaussLobattoQuadratureInitHelper<ct, true> {
471  static void init(int p,
472  std::vector< FieldVector<ct, 1> > & _points,
473  std::vector< ct > & _weight,
474  int & delivered_order);
475  };
476  template<typename ct>
477  struct GaussLobattoQuadratureInitHelper<ct, false> {
478  static void init(int p,
479  std::vector< FieldVector<ct, 1> > & _points,
480  std::vector< ct > & _weight,
481  int & delivered_order);
482  };
483 
487  template<typename ct>
489  public QuadratureRule<ct,1>
490  {
491  public:
493  enum { dim=1 };
494 
496  enum { highest_order=31 };
497 
499  private:
500  friend class QuadratureRuleFactory<ct,dim>;
503  {
505  std::vector< FieldVector<ct, dim> > _points;
506  std::vector< ct > _weight;
507 
508  int deliveredOrder_;
509 
511  (p, _points, _weight, deliveredOrder_);
512 
513  this->delivered_order = deliveredOrder_;
514  assert(_points.size() == _weight.size());
515  for (size_t i = 0; i < _points.size(); i++)
516  this->push_back(QuadraturePoint<ct,dim>(_points[i], _weight[i]));
517  }
518  };
519 
520 #ifndef DOXYGEN
523 #endif // !DOXYGEN
524 
525 } // namespace Dune
526 
527 #define DUNE_INCLUDING_IMPLEMENTATION
528 #include "quadraturerules/gausslobatto_imp.hh"
529 
530 #include "quadraturerules/tensorproductquadrature.hh"
531 
532 #include "quadraturerules/simplexquadrature.hh"
533 
534 namespace Dune {
535 
536  /***********************************
537  * quadrature for Prism
538  **********************************/
539 
541  template<int dim>
543 
545  template<>
547  {
548  public:
549  enum { MAXP=6};
550  enum { highest_order=2 };
551 
554  {
555  int m = 0;
556  O[m] = 0;
557 
558  // polynom degree 0 ???
559  m = 6;
560  G[m][0][0] = 0.0;
561  G[m][0][1] = 0.0;
562  G[m][0][2] = 0.0;
563 
564  G[m][1][0] = 1.0;
565  G[m][1][1] = 0.0;
566  G[m][1][2] = 0.0;
567 
568  G[m][2][0] = 0.0;
569  G[m][2][1] = 1.0;
570  G[m][2][2] = 0.0;
571 
572  G[m][3][0] = 0.0;
573  G[m][3][1] = 0.0;
574  G[m][3][2] = 1.0;
575 
576  G[m][4][0] = 1.0;
577  G[m][4][1] = 0.0;
578  G[m][4][2] = 1.0;
579 
580  G[m][5][0] = 0.0;
581  G[m][5][1] = 0.1;
582  G[m][5][2] = 1.0;
583 
584  W[m][0] = 0.16666666666666666 / 2.0;
585  W[m][1] = 0.16666666666666666 / 2.0;
586  W[m][2] = 0.16666666666666666 / 2.0;
587  W[m][3] = 0.16666666666666666 / 2.0;
588  W[m][4] = 0.16666666666666666 / 2.0;
589  W[m][5] = 0.16666666666666666 / 2.0;
590 
591  O[m] = 0; // verify ????????
592 
593 
594  // polynom degree 2 ???
595  m = 6;
596  G[m][0][0] =0.66666666666666666 ;
597  G[m][0][1] =0.16666666666666666 ;
598  G[m][0][2] =0.211324865405187 ;
599 
600  G[m][1][0] = 0.16666666666666666;
601  G[m][1][1] =0.66666666666666666 ;
602  G[m][1][2] = 0.211324865405187;
603 
604  G[m][2][0] = 0.16666666666666666;
605  G[m][2][1] = 0.16666666666666666;
606  G[m][2][2] = 0.211324865405187;
607 
608  G[m][3][0] = 0.66666666666666666;
609  G[m][3][1] = 0.16666666666666666;
610  G[m][3][2] = 0.788675134594813;
611 
612  G[m][4][0] = 0.16666666666666666;
613  G[m][4][1] = 0.66666666666666666;
614  G[m][4][2] = 0.788675134594813;
615 
616  G[m][5][0] = 0.16666666666666666;
617  G[m][5][1] = 0.16666666666666666;
618  G[m][5][2] = 0.788675134594813;
619 
620  W[m][0] = 0.16666666666666666 / 2.0;
621  W[m][1] = 0.16666666666666666 / 2.0;
622  W[m][2] = 0.16666666666666666 / 2.0;
623  W[m][3] = 0.16666666666666666 / 2.0;
624  W[m][4] = 0.16666666666666666 / 2.0;
625  W[m][5] = 0.16666666666666666 / 2.0;
626 
627  O[m] = 2; // verify ????????
628 
629  }
630 
633  {
634  return G[m][i];
635  }
636 
638  double weight (int m, int i)
639  {
640  return W[m][i];
641  }
642 
644  int order (int m)
645  {
646  return O[m];
647  }
648 
649  private:
650  FieldVector<double, 3> G[MAXP+1][MAXP]; //positions
651 
652  double W[MAXP+1][MAXP]; // weights associated with points
653  int O[MAXP+1]; // order of the rule
654  };
655 
656 
660  template<int dim>
662  static PrismQuadraturePoints<3> prqp;
663  };
664 
668  template<>
670  static PrismQuadraturePoints<3> prqp;
671  };
672 
676  template<typename ct, int dim>
678 
682  template<typename ct>
684  {
685  public:
686 
688  enum { d = 3 };
689 
691  enum { highest_order = 2 };
692 
694  private:
695  friend class QuadratureRuleFactory<ct,d>;
697  {
698  int m=6;
699  this->delivered_order = PrismQuadraturePointsSingleton<3>::prqp.order(m);
700  for(int i=0; i<m; ++i)
701  {
702  FieldVector<ct,3> local;
703  for (int k=0; k<d; k++)
704  local[k] = PrismQuadraturePointsSingleton<3>::prqp.point(m,i)[k];
705  double weight =
707  // put in container
708  this->push_back(QuadraturePoint<ct,d>(local,weight));
709  }
710  }
711  };
712 
719  template<typename ctype, int dim>
721  private:
722  friend class QuadratureRules<ctype, dim>;
723  static unsigned maxOrder(const GeometryType &t, QuadratureType::Enum qt)
724  {
725  return TensorProductQuadratureRule<ctype,dim>::maxOrder(t.id(), qt);
726  }
727  static QuadratureRule<ctype, dim> rule(const GeometryType& t, int p, QuadratureType::Enum qt)
728  {
729  return TensorProductQuadratureRule<ctype,dim>(t.id(), p, qt);
730  }
731  };
732 
733  template<typename ctype>
735  private:
736  enum { dim = 0 };
737  friend class QuadratureRules<ctype, dim>;
738  static unsigned maxOrder(const GeometryType &t, QuadratureType::Enum qt)
739  {
740  if (t.isVertex())
741  {
743  }
744  DUNE_THROW(Exception, "Unknown GeometryType");
745  }
746  static QuadratureRule<ctype, dim> rule(const GeometryType& t, int p, QuadratureType::Enum qt)
747  {
748  if (t.isVertex())
749  {
750  return PointQuadratureRule<ctype>();
751  }
752  DUNE_THROW(Exception, "Unknown GeometryType");
753  }
754  };
755 
756  template<typename ctype>
757  class QuadratureRuleFactory<ctype, 1> {
758  private:
759  enum { dim = 1 };
760  friend class QuadratureRules<ctype, dim>;
761  static unsigned maxOrder(const GeometryType &t, QuadratureType::Enum qt)
762  {
763  if (t.isLine())
764  {
765  switch (qt) {
767  return GaussQuadratureRule1D<ctype>::highest_order;
769  return Jacobi1QuadratureRule1D<ctype>::highest_order;
771  return Jacobi2QuadratureRule1D<ctype>::highest_order;
772  case QuadratureType::GaussLobatto :
773  return GaussLobattoQuadratureRule1D<ctype>::highest_order;
775  return JacobiNQuadratureRule1D<ctype>::maxOrder();
776  default :
777  DUNE_THROW(Exception, "Unknown QuadratureType");
778  }
779  }
780  DUNE_THROW(Exception, "Unknown GeometryType");
781  }
782  static QuadratureRule<ctype, dim> rule(const GeometryType& t, int p, QuadratureType::Enum qt)
783  {
784  if (t.isLine())
785  {
786  switch (qt) {
788  return GaussQuadratureRule1D<ctype>(p);
790  return Jacobi1QuadratureRule1D<ctype>(p);
792  return Jacobi2QuadratureRule1D<ctype>(p);
793  case QuadratureType::GaussLobatto :
794  return GaussLobattoQuadratureRule1D<ctype>(p);
796  return JacobiNQuadratureRule1D<ctype>(p);
797  default :
798  DUNE_THROW(Exception, "Unknown QuadratureType");
799  }
800  }
801  DUNE_THROW(Exception, "Unknown GeometryType");
802  }
803  };
804 
805  template<typename ctype>
806  class QuadratureRuleFactory<ctype, 2> {
807  private:
808  enum { dim = 2 };
809  friend class QuadratureRules<ctype, dim>;
810  static unsigned maxOrder(const GeometryType &t, QuadratureType::Enum qt)
811  {
812  unsigned order =
813  TensorProductQuadratureRule<ctype,dim>::maxOrder(t.id(), qt);
814  if (t.isSimplex())
815  order = std::max
816  (order, unsigned(SimplexQuadratureRule<ctype,dim>::highest_order));
817  return order;
818  }
819  static QuadratureRule<ctype, dim> rule(const GeometryType& t, int p, QuadratureType::Enum qt)
820  {
821  if (t.isSimplex()
823  && p <= SimplexQuadratureRule<ctype,dim>::highest_order)
824  {
825  return SimplexQuadratureRule<ctype,dim>(p);
826  }
827  return TensorProductQuadratureRule<ctype,dim>(t.id(), p, qt);
828  }
829  };
830 
831  template<typename ctype>
832  class QuadratureRuleFactory<ctype, 3> {
833  private:
834  enum { dim = 3 };
835  friend class QuadratureRules<ctype, dim>;
836  static unsigned maxOrder(const GeometryType &t, QuadratureType::Enum qt)
837  {
838  unsigned order =
839  TensorProductQuadratureRule<ctype,dim>::maxOrder(t.id(), qt);
840  if (t.isSimplex())
841  order = std::max
842  (order, unsigned(SimplexQuadratureRule<ctype,dim>::highest_order));
843  if (t.isPrism())
844  order = std::max
845  (order, unsigned(PrismQuadratureRule<ctype,dim>::highest_order));
846  return order;
847  }
848  static QuadratureRule<ctype, dim> rule(const GeometryType& t, int p, QuadratureType::Enum qt)
849  {
850 
851  if (t.isSimplex()
853  && p <= SimplexQuadratureRule<ctype,dim>::highest_order)
854  {
855  return SimplexQuadratureRule<ctype,dim>(p);
856  }
857  if (t.isPrism()
859  && p <= PrismQuadratureRule<ctype,dim>::highest_order)
860  {
861  return PrismQuadratureRule<ctype,dim>(p);
862  }
863  return TensorProductQuadratureRule<ctype,dim>(t.id(), p, qt);
864  }
865  };
866 
867 } // end namespace
868 
869 #endif // DUNE_GEOMETRY_QUADRATURERULES_HH
Jacobi-Gauss quadrature for alpha=2, beta=0.
Definition: quadraturerules.hh:490
Gauss quadrature rule in 1D.
Definition: quadraturerules.hh:297
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:279
constexpr unsigned int dim() const
Return dimension of the type.
Definition: type.hh:644
BasicType
Each entity can be tagged by one of these basic types plus its space dimension.
Definition: type.hh:285
constexpr unsigned int id() const
Return the topology id of the type.
Definition: type.hh:649
Jacobi-Gauss quadrature for alpha=1, beta=0.
Definition: quadraturerules.hh:357
Jacobi-Gauss quadrature for alpha=2, beta=0.
Definition: quadraturerules.hh:423
static constexpr std::size_t size(std::size_t dim)
Compute total number of geometry types for the given dimension.
Definition: typeindex.hh:56
static constexpr std::size_t index(const GeometryType &gt)
Compute the index for the given geometry type within its dimension.
Definition: typeindex.hh:68
Default exception for dummy implementations.
Definition: exceptions.hh:261
Definition: quadraturerules.hh:547
double weight(int m, int i)
Definition: quadraturerules.hh:638
int order(int m)
Definition: quadraturerules.hh:644
PrismQuadraturePoints()
initialize quadrature points on the interval for all orders
Definition: quadraturerules.hh:553
FieldVector< double, 3 > point(int m, int i)
Definition: quadraturerules.hh:632
Definition: quadraturerules.hh:542
Quadrature rules for prisms.
Definition: quadraturerules.hh:684
Quadrature rules for prisms.
Definition: quadraturerules.hh:677
Exception thrown if a desired QuadratureRule is not available, because the requested order is to high...
Definition: quadraturerules.hh:34
Single evaluation point in a quadrature rule.
Definition: quadraturerules.hh:42
Dune::FieldVector< ct, dim > Vector
Type used for the position of a quadrature point.
Definition: quadraturerules.hh:51
ct Field
Number type used for coordinates and quadrature weights.
Definition: quadraturerules.hh:48
const Vector & position() const
return local coordinates of integration point i
Definition: quadraturerules.hh:60
const ct & weight() const
return weight associated with integration point i
Definition: quadraturerules.hh:66
QuadraturePoint(const Vector &x, ct w)
set up quadrature of given order in d dimensions
Definition: quadraturerules.hh:54
Factory class for creation of quadrature rules, depending on GeometryType, order and QuadratureType.
Definition: quadraturerules.hh:720
Abstract base class for quadrature rules.
Definition: quadraturerules.hh:126
virtual GeometryType type() const
return type of element
Definition: quadraturerules.hh:152
QuadratureRule(GeometryType t, int order)
Constructor for a given geometry type and a given quadrature order.
Definition: quadraturerules.hh:140
ct CoordType
The type used for coordinates.
Definition: quadraturerules.hh:146
QuadratureRule()
Default constructor.
Definition: quadraturerules.hh:133
virtual int order() const
return order
Definition: quadraturerules.hh:149
QuadratureRule(GeometryType t)
Constructor for a given geometry type. Leaves the quadrature order invalid
Definition: quadraturerules.hh:137
std::vector< QuadraturePoint< ct, dim > >::const_iterator iterator
Definition: quadraturerules.hh:157
A container for all quadrature rules of dimension dim
Definition: quadraturerules.hh:172
static unsigned maxOrder(const GeometryType &t, QuadratureType::Enum qt=QuadratureType::GaussLegendre)
maximum quadrature order for given geometry type and quadrature type
Definition: quadraturerules.hh:247
static const QuadratureRule & rule(const GeometryType &t, int p, QuadratureType::Enum qt=QuadratureType::GaussLegendre)
select the appropriate QuadratureRule for GeometryType t and order p
Definition: quadraturerules.hh:254
static DUNE_NO_DEPRECATED_BEGIN const QuadratureRule & rule(const GeometryType::BasicType t, int p, QuadratureType::Enum qt=QuadratureType::GaussLegendre)
select the appropriate QuadratureRule for GeometryType t and order p
Definition: quadraturerules.hh:261
A few common exception classes.
Implements a vector constructed from a given type representing a field and a compile-time given size.
#define DUNE_NO_DEPRECATED_END
Ignore deprecation warnings (end)
Definition: deprecated.hh:198
#define DUNE_NO_DEPRECATED_BEGIN
Ignore deprecation warnings (start)
Definition: deprecated.hh:192
#define DUNE_THROW(E, m)
Definition: exceptions.hh:216
bool gt(const T &first, const T &second, typename EpsilonType< T >::Type epsilon)
test if first greater than second
Definition: float_cmp.cc:156
constexpr GeometryType line
GeometryType representing a line.
Definition: type.hh:803
constexpr GeometryType prism
GeometryType representing a 3D prism.
Definition: type.hh:833
auto max(ADLTag< 0 >, const V &v1, const V &v2)
implements binary Simd::max()
Definition: defaults.hh:79
Enum
Definition: quadraturerules.hh:80
@ GaussJacobi_n_0
Gauss-Legendre rules with .
Definition: quadraturerules.hh:115
@ GaussJacobi_2_0
Gauss-Legendre rules with .
Definition: quadraturerules.hh:102
@ GaussJacobi_1_0
Gauss-Jacobi rules with .
Definition: quadraturerules.hh:95
@ GaussLegendre
Gauss-Legendre rules (default)
Definition: quadraturerules.hh:88
Dune namespace.
Definition: alignedallocator.hh:14
Standard Dune debug streams.
Definition: quadraturerules.hh:468
Definition: quadraturerules.hh:277
Definition: quadraturerules.hh:335
Definition: quadraturerules.hh:401
Singleton holding the Prism Quadrature points.
Definition: quadraturerules.hh:669
Singleton holding the Prism Quadrature points.
Definition: quadraturerules.hh:661
A unique label for each type of element that can occur in a grid.
Helper classes to provide indices for geometrytypes for use in a vector.
Definition of macros controlling symbol visibility at the ABI level.
#define DUNE_EXPORT
Export a symbol as part of the public ABI.
Definition: visibility.hh:18
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.80.0 (May 15, 22:30, 2024)