libMesh
quadrature_gauss.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_QUADRATURE_GAUSS_H
21 #define LIBMESH_QUADRATURE_GAUSS_H
22 
23 // Local includes
24 #include "libmesh/quadrature.h"
25 
26 namespace libMesh
27 {
28 
39 class QGauss : public QBase
40 {
41 public:
42 
52  QGauss (const unsigned int _dim,
53  const Order _order=INVALID_ORDER) :
54  QBase(_dim, _order)
55  {
56  if (_dim == 1)
57  init(EDGE2);
58  }
59 
63  ~QGauss() {}
64 
68  virtual QuadratureType type() const libmesh_override { return QGAUSS; }
69 
70 
71 private:
72 
73  virtual void init_1D (const ElemType _type=INVALID_ELEM,
74  unsigned int p_level=0) libmesh_override;
75  virtual void init_2D (const ElemType _type=INVALID_ELEM,
76  unsigned int p_level=0) libmesh_override;
77  virtual void init_3D (const ElemType _type=INVALID_ELEM,
78  unsigned int p_level=0) libmesh_override;
79 
85  void dunavant_rule(const Real rule_data[][4],
86  const unsigned int n_pts);
87 
88  void dunavant_rule2(const Real * wts,
89  const Real * a,
90  const Real * b,
91  const unsigned int * permutation_ids,
92  const unsigned int n_wts);
93 
99  void keast_rule(const Real rule_data[][4],
100  const unsigned int n_pts);
101 };
102 
103 } // namespace libMesh
104 
105 #endif // LIBMESH_QUADRATURE_GAUSS_H
virtual void init(const ElemType type=INVALID_ELEM, unsigned int p_level=0)
Initializes the data structures for a quadrature rule for an element of type type.
Definition: quadrature.C:28
QuadratureType
Defines an enum for currently available quadrature rules.
const unsigned int _dim
The spatial dimension of the quadrature rule.
Definition: quadrature.h:311
virtual void init_1D(const ElemType _type=INVALID_ELEM, unsigned int p_level=0) libmesh_override
Initializes the 1D quadrature rule by filling the points and weights vectors with the appropriate val...
ElemType
Defines an enum for geometric element types.
ElemType _type
The type of element for which the current values have been computed.
Definition: quadrature.h:323
The libMesh namespace provides an interface to certain functionality in the library.
virtual QuadratureType type() const libmesh_override
virtual void init_3D(const ElemType _type=INVALID_ELEM, unsigned int p_level=0) libmesh_override
Initializes the 3D quadrature rule by filling the points and weights vectors with the appropriate val...
void dunavant_rule2(const Real *wts, const Real *a, const Real *b, const unsigned int *permutation_ids, const unsigned int n_wts)
void dunavant_rule(const Real rule_data[][4], const unsigned int n_pts)
The Dunavant rules are for triangles.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void keast_rule(const Real rule_data[][4], const unsigned int n_pts)
The Keast rules are for tets.
~QGauss()
Destructor.
Order
defines an enum for polynomial orders.
Definition: enum_order.h:32
This class implements specific orders of Gauss quadrature.
virtual void init_2D(const ElemType _type=INVALID_ELEM, unsigned int p_level=0) libmesh_override
Initializes the 2D quadrature rule by filling the points and weights vectors with the appropriate val...
QGauss(const unsigned int _dim, const Order _order=INVALID_ORDER)
Constructor.
The QBase class provides the basic functionality from which various quadrature rules can be derived...
Definition: quadrature.h:53
const Order _order
The polynomial order which the quadrature rule is capable of integrating exactly. ...
Definition: quadrature.h:317