libMesh
Public Member Functions | Static Public Member Functions | Public Attributes | Protected Types | Protected Member Functions | Protected Attributes | Static Protected Attributes | Private Member Functions | List of all members
libMesh::QMonomial Class Referencefinal

This class defines alternate quadrature rules on "tensor-product" elements (quadrilaterals and hexahedra) which can be useful when integrating monomial finite element bases. More...

#include <quadrature_monomial.h>

Inheritance diagram for libMesh::QMonomial:
[legend]

Public Member Functions

 QMonomial (unsigned int dim, Order order=INVALID_ORDER)
 Constructor. More...
 
 QMonomial (const QMonomial &)=default
 Copy/move ctor, copy/move assignment operator, and destructor are all explicitly defaulted for this simple class. More...
 
 QMonomial (QMonomial &&)=default
 
QMonomial & operator= (const QMonomial &)=default
 
QMonomial & operator= (QMonomial &&)=default
 
virtual ~QMonomial ()=default
 
virtual QuadratureType type () const override
 
ElemType get_elem_type () const
 
unsigned int get_p_level () const
 
unsigned int n_points () const
 
unsigned int get_dim () const
 
const std::vector< Point > & get_points () const
 
std::vector< Point > & get_points ()
 
const std::vector< Real > & get_weights () const
 
std::vector< Real > & get_weights ()
 
Point qp (const unsigned int i) const
 
Real w (const unsigned int i) const
 
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. More...
 
virtual void init (const Elem &elem, const std::vector< Real > &vertex_distance_func, unsigned int p_level=0)
 Initializes the data structures for an element potentially "cut" by a signed distance function. More...
 
Order get_order () const
 
void print_info (std::ostream &os=libMesh::out) const
 Prints information relevant to the quadrature rule, by default to libMesh::out. More...
 
void scale (std::pair< Real, Real > old_range, std::pair< Real, Real > new_range)
 Maps the points of a 1D quadrature rule defined by "old_range" to another 1D interval defined by "new_range" and scales the weights accordingly. More...
 
virtual bool shapes_need_reinit ()
 

Static Public Member Functions

static std::unique_ptr< QBase > build (const std::string &name, const unsigned int dim, const Order order=INVALID_ORDER)
 Builds a specific quadrature rule based on the name string. More...
 
static std::unique_ptr< QBase > build (const QuadratureType qt, const unsigned int dim, const Order order=INVALID_ORDER)
 Builds a specific quadrature rule based on the QuadratureType. More...
 
static void print_info (std::ostream &out=libMesh::out)
 Prints the reference information, by default to libMesh::out. More...
 
static std::string get_info ()
 Gets a string containing the reference information. More...
 
static unsigned int n_objects ()
 Prints the number of outstanding (created, but not yet destroyed) objects. More...
 
static void enable_print_counter_info ()
 Methods to enable/disable the reference counter output from print_info() More...
 
static void disable_print_counter_info ()
 

Public Attributes

bool allow_rules_with_negative_weights
 Flag (default true) controlling the use of quadrature rules with negative weights. More...
 

Protected Types

typedef std::map< std::string, std::pair< unsigned int, unsigned int > > Counts
 Data structure to log the information. More...
 

Protected Member Functions

virtual void init_0D (const ElemType type=INVALID_ELEM, unsigned int p_level=0)
 Initializes the 0D quadrature rule by filling the points and weights vectors with the appropriate values. More...
 
void tensor_product_quad (const QBase &q1D)
 Constructs a 2D rule from the tensor product of q1D with itself. More...
 
void tensor_product_hex (const QBase &q1D)
 Computes the tensor product quadrature rule [q1D x q1D x q1D] from the 1D rule q1D. More...
 
void tensor_product_prism (const QBase &q1D, const QBase &q2D)
 Computes the tensor product of a 1D quadrature rule and a 2D quadrature rule. More...
 
void increment_constructor_count (const std::string &name)
 Increments the construction counter. More...
 
void increment_destructor_count (const std::string &name)
 Increments the destruction counter. More...
 

Protected Attributes

unsigned int _dim
 The spatial dimension of the quadrature rule. More...
 
Order _order
 The polynomial order which the quadrature rule is capable of integrating exactly. More...
 
ElemType _type
 The type of element for which the current values have been computed. More...
 
unsigned int _p_level
 The p-level of the element for which the current values have been computed. More...
 
std::vector< Point > _points
 The locations of the quadrature points in reference element space. More...
 
std::vector< Real > _weights
 The quadrature weights. More...
 

Static Protected Attributes

static Counts _counts
 Actually holds the data. More...
 
static Threads::atomic< unsigned int > _n_objects
 The number of objects. More...
 
static Threads::spin_mutex _mutex
 Mutual exclusion object to enable thread-safe reference counting. More...
 
static bool _enable_print_counter = true
 Flag to control whether reference count information is printed when print_info is called. More...
 

Private Member Functions

virtual void init_1D (const ElemType, unsigned int) override
 Uses a Gauss rule in 1D. More...
 
virtual void init_2D (const ElemType, unsigned int) override
 Initializes the 2D quadrature rule by filling the points and weights vectors with the appropriate values. More...
 
virtual void init_3D (const ElemType, unsigned int) override
 Initializes the 3D quadrature rule by filling the points and weights vectors with the appropriate values. More...
 
void wissmann_rule (const Real rule_data[][3], const unsigned int n_pts)
 Wissmann published three interesting "partially symmetric" rules for integrating degree 4, 6, and 8 polynomials exactly on QUADs. More...
 
void stroud_rule (const Real rule_data[][3], const unsigned int *rule_symmetry, const unsigned int n_pts)
 Stroud's rules for quads and hexes can have one of several different types of symmetry. More...
 
void kim_rule (const Real rule_data[][4], const unsigned int *rule_id, const unsigned int n_pts)
 Rules from Kim and Song, Comm. More...
 

Detailed Description

This class defines alternate quadrature rules on "tensor-product" elements (quadrilaterals and hexahedra) which can be useful when integrating monomial finite element bases.

While tensor product rules are optimal for integrating bi/tri-linear, bi/tri-quadratic, etc. (i.e. tensor product) bases (which consist of incomplete polynomials up to degree=dim*p) they are not optimal for the MONOMIAL or FEXYZ bases, which consist of complete polynomials of degree=p.

This class provides quadrature rules which are more efficient than tensor product rules when they are available, and falls back on Gaussian quadrature rules otherwise.

A number of these rules have been helpfully collected in electronic form by: Prof. Ronald Cools Katholieke Universiteit Leuven, Dept. Computerwetenschappen http://www.cs.kuleuven.ac.be/~nines/research/ecf/ecf.html A username and password to access the tables is available by request.

We also provide the original reference for each rule when it is available.

Author
John W. Peterson
Date
2008

Implements quadrature rules for non-tensor polynomials.

Definition at line 58 of file quadrature_monomial.h.

Member Typedef Documentation

◆ Counts

typedef std::map<std::string, std::pair<unsigned int, unsigned int> > libMesh::ReferenceCounter::Counts
protectedinherited

Data structure to log the information.

The log is identified by the class name.

Definition at line 117 of file reference_counter.h.

Constructor & Destructor Documentation

◆ QMonomial() [1/3]

libMesh::QMonomial::QMonomial ( unsigned int  dim,
Order  order = INVALID_ORDER 
)
inline

Constructor.

Declares the order of the quadrature rule.

Definition at line 65 of file quadrature_monomial.h.

66  :
67  QBase(dim,order)
68  {}

◆ QMonomial() [2/3]

libMesh::QMonomial::QMonomial ( const QMonomial &  )
default

Copy/move ctor, copy/move assignment operator, and destructor are all explicitly defaulted for this simple class.

◆ QMonomial() [3/3]

libMesh::QMonomial::QMonomial ( QMonomial &&  )
default

◆ ~QMonomial()

virtual libMesh::QMonomial::~QMonomial ( )
virtualdefault

Member Function Documentation

◆ build() [1/2]

std::unique_ptr< QBase > libMesh::QBase::build ( const QuadratureType  qt,
const unsigned int  dim,
const Order  order = INVALID_ORDER 
)
staticinherited

Builds a specific quadrature rule based on the QuadratureType.

This enables selection of the quadrature rule at run-time.

This function allocates memory, therefore a std::unique_ptr<QBase> is returned so that the user does not accidentally leak it.

Definition at line 52 of file quadrature_build.C.

55 {
56  switch (_qt)
57  {
58 
59  case QCLOUGH:
60  {
61 #ifdef DEBUG
62  if (_order > TWENTYTHIRD)
63  {
64  libMesh::out << "WARNING: Clough quadrature implemented" << std::endl
65  << " up to TWENTYTHIRD order." << std::endl;
66  }
67 #endif
68 
69  return libmesh_make_unique<QClough>(_dim, _order);
70  }
71 
72  case QGAUSS:
73  {
74 
75 #ifdef DEBUG
76  if (_order > FORTYTHIRD)
77  {
78  libMesh::out << "WARNING: Gauss quadrature implemented" << std::endl
79  << " up to FORTYTHIRD order." << std::endl;
80  }
81 #endif
82 
83  return libmesh_make_unique<QGauss>(_dim, _order);
84  }
85 
86  case QJACOBI_1_0:
87  {
88 
89 #ifdef DEBUG
90  if (_order > FORTYTHIRD)
91  {
92  libMesh::out << "WARNING: Jacobi(1,0) quadrature implemented" << std::endl
93  << " up to FORTYTHIRD order." << std::endl;
94  }
95 
96  if (_dim > 1)
97  {
98  libMesh::out << "WARNING: Jacobi(1,0) quadrature implemented" << std::endl
99  << " in 1D only." << std::endl;
100  }
101 #endif
102 
103  return libmesh_make_unique<QJacobi>(_dim, _order, 1, 0);
104  }
105 
106  case QJACOBI_2_0:
107  {
108 
109 #ifdef DEBUG
110  if (_order > FORTYTHIRD)
111  {
112  libMesh::out << "WARNING: Jacobi(2,0) quadrature implemented" << std::endl
113  << " up to FORTYTHIRD order." << std::endl;
114  }
115 
116  if (_dim > 1)
117  {
118  libMesh::out << "WARNING: Jacobi(2,0) quadrature implemented" << std::endl
119  << " in 1D only." << std::endl;
120  }
121 #endif
122 
123  return libmesh_make_unique<QJacobi>(_dim, _order, 2, 0);
124  }
125 
126  case QSIMPSON:
127  {
128 
129 #ifdef DEBUG
130  if (_order > THIRD)
131  {
132  libMesh::out << "WARNING: Simpson rule provides only" << std::endl
133  << " THIRD order!" << std::endl;
134  }
135 #endif
136 
137  return libmesh_make_unique<QSimpson>(_dim);
138  }
139 
140  case QTRAP:
141  {
142 
143 #ifdef DEBUG
144  if (_order > FIRST)
145  {
146  libMesh::out << "WARNING: Trapezoidal rule provides only" << std::endl
147  << " FIRST order!" << std::endl;
148  }
149 #endif
150 
151  return libmesh_make_unique<QTrap>(_dim);
152  }
153 
154  case QGRID:
155  return libmesh_make_unique<QGrid>(_dim, _order);
156 
157  case QGRUNDMANN_MOLLER:
158  return libmesh_make_unique<QGrundmann_Moller>(_dim, _order);
159 
160  case QMONOMIAL:
161  return libmesh_make_unique<QMonomial>(_dim, _order);
162 
163  case QGAUSS_LOBATTO:
164  return libmesh_make_unique<QGaussLobatto>(_dim, _order);
165 
166  case QCONICAL:
167  return libmesh_make_unique<QConical>(_dim, _order);
168 
169  case QNODAL:
170  return libmesh_make_unique<QNodal>(_dim, _order);
171 
172  default:
173  libmesh_error_msg("ERROR: Bad qt=" << _qt);
174  }
175 }

References libMesh::QBase::_dim, libMesh::QBase::_order, libMesh::FIRST, libMesh::FORTYTHIRD, libMesh::out, libMesh::QCLOUGH, libMesh::QCONICAL, libMesh::QGAUSS, libMesh::QGAUSS_LOBATTO, libMesh::QGRID, libMesh::QGRUNDMANN_MOLLER, libMesh::QJACOBI_1_0, libMesh::QJACOBI_2_0, libMesh::QMONOMIAL, libMesh::QNODAL, libMesh::QSIMPSON, libMesh::QTRAP, libMesh::THIRD, and libMesh::TWENTYTHIRD.

◆ build() [2/2]

std::unique_ptr< QBase > libMesh::QBase::build ( const std::string &  name,
const unsigned int  dim,
const Order  order = INVALID_ORDER 
)
staticinherited

Builds a specific quadrature rule based on the name string.

This enables selection of the quadrature rule at run-time. The input parameter name must be mappable through the Utility::string_to_enum<>() function.

This function allocates memory, therefore a std::unique_ptr<QBase> is returned so that the user does not accidentally leak it.

Definition at line 41 of file quadrature_build.C.

44 {
45  return QBase::build (Utility::string_to_enum<QuadratureType> (type),
46  _dim,
47  _order);
48 }

References libMesh::QBase::_dim, libMesh::QBase::_order, and libMesh::QBase::type().

Referenced by assemble_poisson(), libMesh::InfFE< Dim, T_radial, T_map >::attach_quadrature_rule(), libMesh::InfFE< Dim, T_radial, T_map >::reinit(), QuadratureTest::test1DWeights(), QuadratureTest::test2DWeights(), QuadratureTest::test3DWeights(), QuadratureTest::testBuild(), QuadratureTest::testJacobi(), QuadratureTest::testMonomialQuadrature(), QuadratureTest::testTetQuadrature(), and QuadratureTest::testTriQuadrature().

◆ disable_print_counter_info()

void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited

Definition at line 106 of file reference_counter.C.

107 {
108  _enable_print_counter = false;
109  return;
110 }

References libMesh::ReferenceCounter::_enable_print_counter.

Referenced by libMesh::LibMeshInit::LibMeshInit().

◆ enable_print_counter_info()

void libMesh::ReferenceCounter::enable_print_counter_info ( )
staticinherited

Methods to enable/disable the reference counter output from print_info()

Definition at line 100 of file reference_counter.C.

101 {
102  _enable_print_counter = true;
103  return;
104 }

References libMesh::ReferenceCounter::_enable_print_counter.

◆ get_dim()

unsigned int libMesh::QBase::get_dim ( ) const
inlineinherited
Returns
The spatial dimension of the quadrature rule.

Definition at line 135 of file quadrature.h.

135 { return _dim; }

References libMesh::QBase::_dim.

Referenced by libMesh::InfFE< Dim, T_radial, T_map >::attach_quadrature_rule(), libMesh::QConical::conical_product_pyramid(), libMesh::QConical::conical_product_tet(), and libMesh::QConical::conical_product_tri().

◆ get_elem_type()

ElemType libMesh::QBase::get_elem_type ( ) const
inlineinherited
Returns
The element type we're currently using.

Definition at line 116 of file quadrature.h.

116 { return _type; }

References libMesh::QBase::_type.

◆ get_info()

std::string libMesh::ReferenceCounter::get_info ( )
staticinherited

Gets a string containing the reference information.

Definition at line 47 of file reference_counter.C.

48 {
49 #if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
50 
51  std::ostringstream oss;
52 
53  oss << '\n'
54  << " ---------------------------------------------------------------------------- \n"
55  << "| Reference count information |\n"
56  << " ---------------------------------------------------------------------------- \n";
57 
58  for (const auto & pr : _counts)
59  {
60  const std::string name(pr.first);
61  const unsigned int creations = pr.second.first;
62  const unsigned int destructions = pr.second.second;
63 
64  oss << "| " << name << " reference count information:\n"
65  << "| Creations: " << creations << '\n'
66  << "| Destructions: " << destructions << '\n';
67  }
68 
69  oss << " ---------------------------------------------------------------------------- \n";
70 
71  return oss.str();
72 
73 #else
74 
75  return "";
76 
77 #endif
78 }

References libMesh::ReferenceCounter::_counts, and libMesh::Quality::name().

Referenced by libMesh::ReferenceCounter::print_info().

◆ get_order()

Order libMesh::QBase::get_order ( ) const
inlineinherited
Returns
The current "total" order of the quadrature rule which can vary element by element, depending on the Elem::p_level(), which gets passed to us during init().

Each additional power of p increases the quadrature order required to integrate the mass matrix by 2, hence the formula below.

Definition at line 211 of file quadrature.h.

211 { return static_cast<Order>(_order + 2 * _p_level); }

References libMesh::QBase::_order, and libMesh::QBase::_p_level.

Referenced by libMesh::InfFE< Dim, T_radial, T_map >::attach_quadrature_rule(), libMesh::QConical::conical_product_pyramid(), libMesh::QConical::conical_product_tet(), libMesh::QConical::conical_product_tri(), libMesh::QGaussLobatto::init_1D(), libMesh::QConical::init_1D(), libMesh::QGauss::init_1D(), libMesh::QJacobi::init_1D(), libMesh::QGrundmann_Moller::init_1D(), libMesh::QGauss::init_2D(), init_2D(), libMesh::QGrundmann_Moller::init_2D(), libMesh::QGauss::init_3D(), init_3D(), and libMesh::QGrundmann_Moller::init_3D().

◆ get_p_level()

unsigned int libMesh::QBase::get_p_level ( ) const
inlineinherited
Returns
The p-refinement level we're currently using.

Definition at line 121 of file quadrature.h.

121 { return _p_level; }

References libMesh::QBase::_p_level.

◆ get_points() [1/2]

std::vector<Point>& libMesh::QBase::get_points ( )
inlineinherited
Returns
A std::vector containing the quadrature point locations in reference element space as a writable reference.

Definition at line 147 of file quadrature.h.

147 { return _points; }

References libMesh::QBase::_points.

◆ get_points() [2/2]

const std::vector<Point>& libMesh::QBase::get_points ( ) const
inlineinherited
Returns
A std::vector containing the quadrature point locations in reference element space.

Definition at line 141 of file quadrature.h.

141 { return _points; }

References libMesh::QBase::_points.

Referenced by assemble_ellipticdg(), libMesh::QClough::init_1D(), libMesh::QConical::init_1D(), libMesh::QNodal::init_1D(), init_1D(), libMesh::QGrundmann_Moller::init_1D(), libMesh::QClough::init_2D(), libMesh::QGauss::init_2D(), libMesh::QNodal::init_2D(), init_2D(), libMesh::QGauss::init_3D(), libMesh::QNodal::init_3D(), init_3D(), EIM_F::interior_assembly(), and AssemblyEIM::interior_assembly().

◆ get_weights() [1/2]

std::vector<Real>& libMesh::QBase::get_weights ( )
inlineinherited
Returns
A writable references to a std::vector containing the quadrature weights.

Definition at line 159 of file quadrature.h.

159 { return _weights; }

References libMesh::QBase::_weights.

◆ get_weights() [2/2]

const std::vector<Real>& libMesh::QBase::get_weights ( ) const
inlineinherited
Returns
A constant reference to a std::vector containing the quadrature weights.

Definition at line 153 of file quadrature.h.

153 { return _weights; }

References libMesh::QBase::_weights.

Referenced by libMesh::QClough::init_1D(), libMesh::QConical::init_1D(), libMesh::QNodal::init_1D(), init_1D(), libMesh::QGrundmann_Moller::init_1D(), libMesh::QClough::init_2D(), libMesh::QGauss::init_2D(), libMesh::QNodal::init_2D(), init_2D(), libMesh::QGauss::init_3D(), libMesh::QNodal::init_3D(), and init_3D().

◆ increment_constructor_count()

void libMesh::ReferenceCounter::increment_constructor_count ( const std::string &  name)
inlineprotectedinherited

Increments the construction counter.

Should be called in the constructor of any derived class that will be reference counted.

Definition at line 181 of file reference_counter.h.

182 {
183  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
184  std::pair<unsigned int, unsigned int> & p = _counts[name];
185 
186  p.first++;
187 }

References libMesh::ReferenceCounter::_counts, libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCountedObject< RBParametrized >::ReferenceCountedObject().

◆ increment_destructor_count()

void libMesh::ReferenceCounter::increment_destructor_count ( const std::string &  name)
inlineprotectedinherited

Increments the destruction counter.

Should be called in the destructor of any derived class that will be reference counted.

Definition at line 194 of file reference_counter.h.

195 {
196  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
197  std::pair<unsigned int, unsigned int> & p = _counts[name];
198 
199  p.second++;
200 }

References libMesh::ReferenceCounter::_counts, libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCountedObject< RBParametrized >::~ReferenceCountedObject().

◆ init() [1/2]

void libMesh::QBase::init ( const Elem &  elem,
const std::vector< Real > &  vertex_distance_func,
unsigned int  p_level = 0 
)
virtualinherited

Initializes the data structures for an element potentially "cut" by a signed distance function.

The array vertex_distance_func contains vertex values of the signed distance function. If the signed distance function changes sign on the vertices, then the element is considered to be cut.) This interface can be extended by derived classes in order to subdivide the element and construct a composite quadrature rule.

Definition at line 103 of file quadrature.C.

106 {
107  // dispatch generic implementation
108  this->init(elem.type(), p_level);
109 }

References libMesh::QBase::init(), and libMesh::Elem::type().

◆ init() [2/2]

void libMesh::QBase::init ( const ElemType  type = INVALID_ELEM,
unsigned int  p_level = 0 
)
virtualinherited

Initializes the data structures for a quadrature rule for an element of type type.

Definition at line 59 of file quadrature.C.

61 {
62  // check to see if we have already
63  // done the work for this quadrature rule
64  if (t == _type && p == _p_level)
65  return;
66  else
67  {
68  _type = t;
69  _p_level = p;
70  }
71 
72 
73 
74  switch(_dim)
75  {
76  case 0:
77  this->init_0D();
78 
79  return;
80 
81  case 1:
82  this->init_1D();
83 
84  return;
85 
86  case 2:
87  this->init_2D();
88 
89  return;
90 
91  case 3:
92  this->init_3D();
93 
94  return;
95 
96  default:
97  libmesh_error_msg("Invalid dimension _dim = " << _dim);
98  }
99 }

References libMesh::QBase::_dim, libMesh::QBase::_p_level, libMesh::QBase::_type, libMesh::QBase::init_0D(), libMesh::QBase::init_1D(), libMesh::QBase::init_2D(), and libMesh::QBase::init_3D().

Referenced by libMesh::QBase::init(), libMesh::QClough::init_1D(), libMesh::QNodal::init_1D(), init_1D(), libMesh::QClough::init_2D(), libMesh::QGaussLobatto::init_2D(), libMesh::QTrap::init_2D(), libMesh::QGrid::init_2D(), libMesh::QGauss::init_2D(), libMesh::QSimpson::init_2D(), libMesh::QNodal::init_2D(), init_2D(), libMesh::QGaussLobatto::init_3D(), libMesh::QTrap::init_3D(), libMesh::QGrid::init_3D(), libMesh::QSimpson::init_3D(), libMesh::QGauss::init_3D(), libMesh::QNodal::init_3D(), init_3D(), libMesh::QGauss::QGauss(), libMesh::QGaussLobatto::QGaussLobatto(), libMesh::QJacobi::QJacobi(), libMesh::QNodal::QNodal(), libMesh::QSimpson::QSimpson(), and libMesh::QTrap::QTrap().

◆ init_0D()

void libMesh::QBase::init_0D ( const ElemType  type = INVALID_ELEM,
unsigned int  p_level = 0 
)
protectedvirtualinherited

Initializes the 0D quadrature rule by filling the points and weights vectors with the appropriate values.

Generally this is just one point with weight 1.

Note
The arguments should no longer be used for anything in derived classes, they are only maintained for backwards compatibility and will eventually be removed.

Definition at line 113 of file quadrature.C.

114 {
115  _points.resize(1);
116  _weights.resize(1);
117  _points[0] = Point(0.);
118  _weights[0] = 1.0;
119 }

References libMesh::QBase::_points, and libMesh::QBase::_weights.

Referenced by libMesh::QBase::init().

◆ init_1D()

void libMesh::QMonomial::init_1D ( const  ElemType,
unsigned int   
)
overrideprivatevirtual

Uses a Gauss rule in 1D.

More efficient rules for non tensor product bases on quadrilaterals and hexahedra.

Implements libMesh::QBase.

Definition at line 29 of file quadrature_monomial_1D.C.

30 {
31  QGauss gauss_rule(1, _order);
32  gauss_rule.init(_type, _p_level);
33 
34  _points.swap(gauss_rule.get_points());
35  _weights.swap(gauss_rule.get_weights());
36 }

References libMesh::QBase::_order, libMesh::QBase::_p_level, libMesh::QBase::_points, libMesh::QBase::_type, libMesh::QBase::_weights, libMesh::QBase::get_points(), libMesh::QBase::get_weights(), and libMesh::QBase::init().

◆ init_2D()

void libMesh::QMonomial::init_2D ( const  type,
unsigned int  p_level 
)
overrideprivatevirtual

Initializes the 2D quadrature rule by filling the points and weights vectors with the appropriate values.

The order of the rule will be defined by the implementing class. Should not be pure virtual since a derived quadrature rule may only be defined in 1D. If not overridden, throws an error.

Note
The arguments should no longer be used for anything in derived classes, they are only maintained for backwards compatibility and will eventually be removed.

Reimplemented from libMesh::QBase.

Definition at line 28 of file quadrature_monomial_2D.C.

29 {
30 
31  switch (_type)
32  {
33  //---------------------------------------------
34  // Quadrilateral quadrature rules
35  case QUAD4:
36  case QUADSHELL4:
37  case QUAD8:
38  case QUADSHELL8:
39  case QUAD9:
40  {
41  switch(get_order())
42  {
43  case SECOND:
44  {
45  // A degree=2 rule for the QUAD with 3 points.
46  // A tensor product degree-2 Gauss would have 4 points.
47  // This rule (or a variation on it) is probably available in
48  //
49  // A.H. Stroud, Approximate calculation of multiple integrals,
50  // Prentice-Hall, Englewood Cliffs, N.J., 1971.
51  //
52  // though I have never actually seen a reference for it.
53  // Luckily it's fairly easy to derive, which is what I've done
54  // here [JWP].
55  const Real
56  s=std::sqrt(Real(1)/3),
57  t=std::sqrt(Real(2)/3);
58 
59  const Real data[2][3] =
60  {
61  {0.0, s, 2.0},
62  { t, -s, 1.0}
63  };
64 
65  _points.resize(3);
66  _weights.resize(3);
67 
68  wissmann_rule(data, 2);
69 
70  return;
71  } // end case SECOND
72 
73 
74 
75  // For third-order, fall through to default case, use 2x2 Gauss product rule.
76  // case THIRD:
77  // {
78  // } // end case THIRD
79 
80  // Tabulated-in-double-precision rules aren't accurate enough for
81  // higher precision, so fall back on Gauss
82 #if !defined(LIBMESH_DEFAULT_TRIPLE_PRECISION) && !defined(LIBMESH_DEFAULT_QUADRUPLE_PRECISION)
83  case FOURTH:
84  {
85  // A pair of degree=4 rules for the QUAD "C2" due to
86  // Wissmann and Becker. These rules both have six points.
87  // A tensor product degree-4 Gauss would have 9 points.
88  //
89  // J. W. Wissmann and T. Becker, Partially symmetric cubature
90  // formulas for even degrees of exactness, SIAM J. Numer. Anal. 23
91  // (1986), 676--685.
92  const Real data[4][3] =
93  {
94  // First of 2 degree-4 rules given by Wissmann
95  {Real(0.0000000000000000e+00), Real(0.0000000000000000e+00), Real(1.1428571428571428e+00)},
96  {Real(0.0000000000000000e+00), Real(9.6609178307929590e-01), Real(4.3956043956043956e-01)},
97  {Real(8.5191465330460049e-01), Real(4.5560372783619284e-01), Real(5.6607220700753210e-01)},
98  {Real(6.3091278897675402e-01), Real(-7.3162995157313452e-01), Real(6.4271900178367668e-01)}
99  //
100  // Second of 2 degree-4 rules given by Wissmann. These both
101  // yield 4th-order accurate rules, I just chose the one that
102  // happened to contain the origin.
103  // {0.000000000000000, -0.356822089773090, 1.286412084888852},
104  // {0.000000000000000, 0.934172358962716, 0.491365692888926},
105  // {0.774596669241483, 0.390885162530071, 0.761883709085613},
106  // {0.774596669241483, -0.852765377881771, 0.349227402025498}
107  };
108 
109  _points.resize(6);
110  _weights.resize(6);
111 
112  wissmann_rule(data, 4);
113 
114  return;
115  } // end case FOURTH
116 #endif
117 
118 
119 
120 
121  case FIFTH:
122  {
123  // A degree 5, 7-point rule due to Stroud.
124  //
125  // A.H. Stroud, Approximate calculation of multiple integrals,
126  // Prentice-Hall, Englewood Cliffs, N.J., 1971.
127  //
128  // This rule is provably minimal in the number of points.
129  // A tensor-product rule accurate for "bi-quintic" polynomials would have 9 points.
130  const Real data[3][3] =
131  {
132  { 0.L, 0.L, Real(8)/7 }, // 1
133  { 0.L, std::sqrt(Real(14)/15), Real(20)/63}, // 2
134  {std::sqrt(Real(3)/5), std::sqrt(Real(1)/3), Real(20)/36} // 4
135  };
136 
137  const unsigned int symmetry[3] = {
138  0, // Origin
139  7, // Central Symmetry
140  6 // Rectangular
141  };
142 
143  _points.resize (7);
144  _weights.resize(7);
145 
146  stroud_rule(data, symmetry, 3);
147 
148  return;
149  } // end case FIFTH
150 
151 
152 
153 
154  // Tabulated-in-double-precision rules aren't accurate enough for
155  // higher precision, so fall back on Gauss
156 #if !defined(LIBMESH_DEFAULT_TRIPLE_PRECISION) && !defined(LIBMESH_DEFAULT_QUADRUPLE_PRECISION)
157  case SIXTH:
158  {
159  // A pair of degree=6 rules for the QUAD "C2" due to
160  // Wissmann and Becker. These rules both have 10 points.
161  // A tensor product degree-6 Gauss would have 16 points.
162  //
163  // J. W. Wissmann and T. Becker, Partially symmetric cubature
164  // formulas for even degrees of exactness, SIAM J. Numer. Anal. 23
165  // (1986), 676--685.
166  const Real data[6][3] =
167  {
168  // First of 2 degree-6, 10 point rules given by Wissmann
169  // {0.000000000000000, 0.836405633697626, 0.455343245714174},
170  // {0.000000000000000, -0.357460165391307, 0.827395973202966},
171  // {0.888764014654765, 0.872101531193131, 0.144000884599645},
172  // {0.604857639464685, 0.305985162155427, 0.668259104262665},
173  // {0.955447506641064, -0.410270899466658, 0.225474004890679},
174  // {0.565459993438754, -0.872869311156879, 0.320896396788441}
175  //
176  // Second of 2 degree-6, 10 point rules given by Wissmann.
177  // Either of these will work, I just chose the one with points
178  // slightly further into the element interior.
179  {Real(0.0000000000000000e+00), Real(8.6983337525005900e-01), Real(3.9275059096434794e-01)},
180  {Real(0.0000000000000000e+00), Real(-4.7940635161211124e-01), Real(7.5476288124261053e-01)},
181  {Real(8.6374282634615388e-01), Real(8.0283751620765670e-01), Real(2.0616605058827902e-01)},
182  {Real(5.1869052139258234e-01), Real(2.6214366550805818e-01), Real(6.8999213848986375e-01)},
183  {Real(9.3397254497284950e-01), Real(-3.6309658314806653e-01), Real(2.6051748873231697e-01)},
184  {Real(6.0897753601635630e-01), Real(-8.9660863276245265e-01), Real(2.6956758608606100e-01)}
185  };
186 
187  _points.resize(10);
188  _weights.resize(10);
189 
190  wissmann_rule(data, 6);
191 
192  return;
193  } // end case SIXTH
194 #endif
195 
196 
197 
198 
199  case SEVENTH:
200  {
201  // A degree 7, 12-point rule due to Tyler, can be found in Stroud's book
202  //
203  // A.H. Stroud, Approximate calculation of multiple integrals,
204  // Prentice-Hall, Englewood Cliffs, N.J., 1971.
205  //
206  // This rule is fully-symmetric and provably minimal in the number of points.
207  // A tensor-product rule accurate for "bi-septic" polynomials would have 16 points.
208  const Real
209  r = std::sqrt(Real(6)/7),
210  s = std::sqrt( (114 - 3*std::sqrt(Real(583))) / 287 ),
211  t = std::sqrt( (114 + 3*std::sqrt(Real(583))) / 287 ),
212  B1 = Real(196)/810,
213  B2 = 4 * (178981 + 2769*std::sqrt(Real(583))) / 1888920,
214  B3 = 4 * (178981 - 2769*std::sqrt(Real(583))) / 1888920;
215 
216  const Real data[3][3] =
217  {
218  {r, 0.0, B1}, // 4
219  {s, 0.0, B2}, // 4
220  {t, 0.0, B3} // 4
221  };
222 
223  const unsigned int symmetry[3] = {
224  3, // Full Symmetry, (x,0)
225  2, // Full Symmetry, (x,x)
226  2 // Full Symmetry, (x,x)
227  };
228 
229  _points.resize (12);
230  _weights.resize(12);
231 
232  stroud_rule(data, symmetry, 3);
233 
234  return;
235  } // end case SEVENTH
236 
237 
238 
239 
240  // Tabulated-in-double-precision rules aren't accurate enough for
241  // higher precision, so fall back on Gauss
242 #if !defined(LIBMESH_DEFAULT_TRIPLE_PRECISION) && !defined(LIBMESH_DEFAULT_QUADRUPLE_PRECISION)
243  case EIGHTH:
244  {
245  // A pair of degree=8 rules for the QUAD "C2" due to
246  // Wissmann and Becker. These rules both have 16 points.
247  // A tensor product degree-6 Gauss would have 25 points.
248  //
249  // J. W. Wissmann and T. Becker, Partially symmetric cubature
250  // formulas for even degrees of exactness, SIAM J. Numer. Anal. 23
251  // (1986), 676--685.
252  const Real data[10][3] =
253  {
254  // First of 2 degree-8, 16 point rules given by Wissmann
255  // {0.000000000000000, 0.000000000000000, 0.055364705621440},
256  // {0.000000000000000, 0.757629177660505, 0.404389368726076},
257  // {0.000000000000000, -0.236871842255702, 0.533546604952635},
258  // {0.000000000000000, -0.989717929044527, 0.117054188786739},
259  // {0.639091304900370, 0.950520955645667, 0.125614417613747},
260  // {0.937069076924990, 0.663882736885633, 0.136544584733588},
261  // {0.537083530541494, 0.304210681724104, 0.483408479211257},
262  // {0.887188506449625, -0.236496718536120, 0.252528506429544},
263  // {0.494698820670197, -0.698953476086564, 0.361262323882172},
264  // {0.897495818279768, -0.900390774211580, 0.085464254086247}
265  //
266  // Second of 2 degree-8, 16 point rules given by Wissmann.
267  // Either of these will work, I just chose the one with points
268  // further into the element interior.
269  {Real(0.0000000000000000e+00), Real(6.5956013196034176e-01), Real(4.5027677630559029e-01)},
270  {Real(0.0000000000000000e+00), Real(-9.4914292304312538e-01), Real(1.6657042677781274e-01)},
271  {Real(9.5250946607156228e-01), Real(7.6505181955768362e-01), Real(9.8869459933431422e-02)},
272  {Real(5.3232745407420624e-01), Real(9.3697598108841598e-01), Real(1.5369674714081197e-01)},
273  {Real(6.8473629795173504e-01), Real(3.3365671773574759e-01), Real(3.9668697607290278e-01)},
274  {Real(2.3314324080140552e-01), Real(-7.9583272377396852e-02), Real(3.5201436794569501e-01)},
275  {Real(9.2768331930611748e-01), Real(-2.7224008061253425e-01), Real(1.8958905457779799e-01)},
276  {Real(4.5312068740374942e-01), Real(-6.1373535339802760e-01), Real(3.7510100114758727e-01)},
277  {Real(8.3750364042281223e-01), Real(-8.8847765053597136e-01), Real(1.2561879164007201e-01)}
278  };
279 
280  _points.resize(16);
281  _weights.resize(16);
282 
283  wissmann_rule(data, /*10*/ 9);
284 
285  return;
286  } // end case EIGHTH
287 
288 
289 
290 
291  case NINTH:
292  {
293  // A degree 9, 17-point rule due to Moller.
294  //
295  // H.M. Moller, Kubaturformeln mit minimaler Knotenzahl,
296  // Numer. Math. 25 (1976), 185--200.
297  //
298  // This rule is provably minimal in the number of points.
299  // A tensor-product rule accurate for "bi-ninth" degree polynomials would have 25 points.
300  const Real data[5][3] =
301  {
302  {Real(0.0000000000000000e+00), Real(0.0000000000000000e+00), Real(5.2674897119341563e-01)}, // 1
303  {Real(6.3068011973166885e-01), Real(9.6884996636197772e-01), Real(8.8879378170198706e-02)}, // 4
304  {Real(9.2796164595956966e-01), Real(7.5027709997890053e-01), Real(1.1209960212959648e-01)}, // 4
305  {Real(4.5333982113564719e-01), Real(5.2373582021442933e-01), Real(3.9828243926207009e-01)}, // 4
306  {Real(8.5261572933366230e-01), Real(7.6208328192617173e-02), Real(2.6905133763978080e-01)} // 4
307  };
308 
309  const unsigned int symmetry[5] = {
310  0, // Single point
311  4, // Rotational Invariant
312  4, // Rotational Invariant
313  4, // Rotational Invariant
314  4 // Rotational Invariant
315  };
316 
317  _points.resize (17);
318  _weights.resize(17);
319 
320  stroud_rule(data, symmetry, 5);
321 
322  return;
323  } // end case NINTH
324 
325 
326 
327 
328  case TENTH:
329  case ELEVENTH:
330  {
331  // A degree 11, 24-point rule due to Cools and Haegemans.
332  //
333  // R. Cools and A. Haegemans, Another step forward in searching for
334  // cubature formulae with a minimal number of knots for the square,
335  // Computing 40 (1988), 139--146.
336  //
337  // P. Verlinden and R. Cools, The algebraic construction of a minimal
338  // cubature formula of degree 11 for the square, Cubature Formulas
339  // and their Applications (Russian) (Krasnoyarsk) (M.V. Noskov, ed.),
340  // 1994, pp. 13--23.
341  //
342  // This rule is provably minimal in the number of points.
343  // A tensor-product rule accurate for "bi-tenth" or "bi-eleventh" degree polynomials would have 36 points.
344  const Real data[6][3] =
345  {
346  {Real(6.9807610454956756e-01), Real(9.8263922354085547e-01), Real(4.8020763350723814e-02)}, // 4
347  {Real(9.3948638281673690e-01), Real(8.2577583590296393e-01), Real(6.6071329164550595e-02)}, // 4
348  {Real(9.5353952820153201e-01), Real(1.8858613871864195e-01), Real(9.7386777358668164e-02)}, // 4
349  {Real(3.1562343291525419e-01), Real(8.1252054830481310e-01), Real(2.1173634999894860e-01)}, // 4
350  {Real(7.1200191307533630e-01), Real(5.2532025036454776e-01), Real(2.2562606172886338e-01)}, // 4
351  {Real(4.2484724884866925e-01), Real(4.1658071912022368e-02), Real(3.5115871839824543e-01)} // 4
352  };
353 
354  const unsigned int symmetry[6] = {
355  4, // Rotational Invariant
356  4, // Rotational Invariant
357  4, // Rotational Invariant
358  4, // Rotational Invariant
359  4, // Rotational Invariant
360  4 // Rotational Invariant
361  };
362 
363  _points.resize (24);
364  _weights.resize(24);
365 
366  stroud_rule(data, symmetry, 6);
367 
368  return;
369  } // end case TENTH,ELEVENTH
370 
371 
372 
373 
374  case TWELFTH:
375  case THIRTEENTH:
376  {
377  // A degree 13, 33-point rule due to Cools and Haegemans.
378  //
379  // R. Cools and A. Haegemans, Another step forward in searching for
380  // cubature formulae with a minimal number of knots for the square,
381  // Computing 40 (1988), 139--146.
382  //
383  // A tensor-product rule accurate for "bi-12" or "bi-13" degree polynomials would have 49 points.
384  const Real data[9][3] =
385  {
386  {Real(0.0000000000000000e+00), Real(0.0000000000000000e+00), Real(3.0038211543122536e-01)}, // 1
387  {Real(9.8348668243987226e-01), Real(7.7880971155441942e-01), Real(2.9991838864499131e-02)}, // 4
388  {Real(8.5955600564163892e-01), Real(9.5729769978630736e-01), Real(3.8174421317083669e-02)}, // 4
389  {Real(9.5892517028753485e-01), Real(1.3818345986246535e-01), Real(6.0424923817749980e-02)}, // 4
390  {Real(3.9073621612946100e-01), Real(9.4132722587292523e-01), Real(7.7492738533105339e-02)}, // 4
391  {Real(8.5007667369974857e-01), Real(4.7580862521827590e-01), Real(1.1884466730059560e-01)}, // 4
392  {Real(6.4782163718701073e-01), Real(7.5580535657208143e-01), Real(1.2976355037000271e-01)}, // 4
393  {Real(7.0741508996444936e-02), Real(6.9625007849174941e-01), Real(2.1334158145718938e-01)}, // 4
394  {Real(4.0930456169403884e-01), Real(3.4271655604040678e-01), Real(2.5687074948196783e-01)} // 4
395  };
396 
397  const unsigned int symmetry[9] = {
398  0, // Single point
399  4, // Rotational Invariant
400  4, // Rotational Invariant
401  4, // Rotational Invariant
402  4, // Rotational Invariant
403  4, // Rotational Invariant
404  4, // Rotational Invariant
405  4, // Rotational Invariant
406  4 // Rotational Invariant
407  };
408 
409  _points.resize (33);
410  _weights.resize(33);
411 
412  stroud_rule(data, symmetry, 9);
413 
414  return;
415  } // end case TWELFTH,THIRTEENTH
416 
417 
418 
419 
420  case FOURTEENTH:
421  case FIFTEENTH:
422  {
423  // A degree-15, 48 point rule originally due to Rabinowitz and Richter,
424  // can be found in Cools' 1971 book.
425  //
426  // A.H. Stroud, Approximate calculation of multiple integrals,
427  // Prentice-Hall, Englewood Cliffs, N.J., 1971.
428  //
429  // The product Gauss rule for this order has 8^2=64 points.
430  const Real data[9][3] =
431  {
432  {Real(9.915377816777667e-01L), Real(0.0000000000000000e+00), Real(3.01245207981210e-02L)}, // 4
433  {Real(8.020163879230440e-01L), Real(0.0000000000000000e+00), Real(8.71146840209092e-02L)}, // 4
434  {Real(5.648674875232742e-01L), Real(0.0000000000000000e+00), Real(1.250080294351494e-01L)}, // 4
435  {Real(9.354392392539896e-01L), Real(0.0000000000000000e+00), Real(2.67651407861666e-02L)}, // 4
436  {Real(7.624563338825799e-01L), Real(0.0000000000000000e+00), Real(9.59651863624437e-02L)}, // 4
437  {Real(2.156164241427213e-01L), Real(0.0000000000000000e+00), Real(1.750832998343375e-01L)}, // 4
438  {Real(9.769662659711761e-01L), Real(6.684480048977932e-01L), Real(2.83136372033274e-02L)}, // 4
439  {Real(8.937128379503403e-01L), Real(3.735205277617582e-01L), Real(8.66414716025093e-02L)}, // 4
440  {Real(6.122485619312083e-01L), Real(4.078983303613935e-01L), Real(1.150144605755996e-01L)} // 4
441  };
442 
443  const unsigned int symmetry[9] = {
444  3, // Full Symmetry, (x,0)
445  3, // Full Symmetry, (x,0)
446  3, // Full Symmetry, (x,0)
447  2, // Full Symmetry, (x,x)
448  2, // Full Symmetry, (x,x)
449  2, // Full Symmetry, (x,x)
450  1, // Full Symmetry, (x,y)
451  1, // Full Symmetry, (x,y)
452  1, // Full Symmetry, (x,y)
453  };
454 
455  _points.resize (48);
456  _weights.resize(48);
457 
458  stroud_rule(data, symmetry, 9);
459 
460  return;
461  } // case FOURTEENTH, FIFTEENTH:
462 
463 
464 
465 
466  case SIXTEENTH:
467  case SEVENTEENTH:
468  {
469  // A degree 17, 60-point rule due to Cools and Haegemans.
470  //
471  // R. Cools and A. Haegemans, Another step forward in searching for
472  // cubature formulae with a minimal number of knots for the square,
473  // Computing 40 (1988), 139--146.
474  //
475  // A tensor-product rule accurate for "bi-14" or "bi-15" degree polynomials would have 64 points.
476  // A tensor-product rule accurate for "bi-16" or "bi-17" degree polynomials would have 81 points.
477  const Real data[10][3] =
478  {
479  {Real(9.8935307451260049e-01), Real(0.0000000000000000e+00), Real(2.0614915919990959e-02)}, // 4
480  {Real(3.7628520715797329e-01), Real(0.0000000000000000e+00), Real(1.2802571617990983e-01)}, // 4
481  {Real(9.7884827926223311e-01), Real(0.0000000000000000e+00), Real(5.5117395340318905e-03)}, // 4
482  {Real(8.8579472916411612e-01), Real(0.0000000000000000e+00), Real(3.9207712457141880e-02)}, // 4
483  {Real(1.7175612383834817e-01), Real(0.0000000000000000e+00), Real(7.6396945079863302e-02)}, // 4
484  {Real(5.9049927380600241e-01), Real(3.1950503663457394e-01), Real(1.4151372994997245e-01)}, // 8
485  {Real(7.9907913191686325e-01), Real(5.9797245192945738e-01), Real(8.3903279363797602e-02)}, // 8
486  {Real(8.0374396295874471e-01), Real(5.8344481776550529e-02), Real(6.0394163649684546e-02)}, // 8
487  {Real(9.3650627612749478e-01), Real(3.4738631616620267e-01), Real(5.7387752969212695e-02)}, // 8
488  {Real(9.8132117980545229e-01), Real(7.0600028779864611e-01), Real(2.1922559481863763e-02)}, // 8
489  };
490 
491  const unsigned int symmetry[10] = {
492  3, // Fully symmetric (x,0)
493  3, // Fully symmetric (x,0)
494  2, // Fully symmetric (x,x)
495  2, // Fully symmetric (x,x)
496  2, // Fully symmetric (x,x)
497  1, // Fully symmetric (x,y)
498  1, // Fully symmetric (x,y)
499  1, // Fully symmetric (x,y)
500  1, // Fully symmetric (x,y)
501  1 // Fully symmetric (x,y)
502  };
503 
504  _points.resize (60);
505  _weights.resize(60);
506 
507  stroud_rule(data, symmetry, 10);
508 
509  return;
510  } // end case FOURTEENTH through SEVENTEENTH
511 #endif
512 
513 
514 
515  // By default: construct and use a Gauss quadrature rule
516  default:
517  {
518  // Break out and fall down into the default: case for the
519  // outer switch statement.
520  break;
521  }
522 
523  } // end switch(_order + 2*p)
524  } // end case QUAD4/8/9
525 
526  libmesh_fallthrough();
527 
528  // By default: construct and use a Gauss quadrature rule
529  default:
530  {
531  QGauss gauss_rule(2, _order);
532  gauss_rule.init(_type, _p_level);
533 
534  // Swap points and weights with the about-to-be destroyed rule.
535  _points.swap (gauss_rule.get_points() );
536  _weights.swap(gauss_rule.get_weights());
537 
538  return;
539  }
540  } // end switch (_type)
541 }

References libMesh::QBase::_order, libMesh::QBase::_p_level, libMesh::QBase::_points, libMesh::QBase::_type, libMesh::QBase::_weights, data, libMesh::EIGHTH, libMesh::ELEVENTH, libMesh::FIFTEENTH, libMesh::FIFTH, libMesh::FOURTEENTH, libMesh::FOURTH, libMesh::QBase::get_order(), libMesh::QBase::get_points(), libMesh::QBase::get_weights(), libMesh::QBase::init(), libMesh::NINTH, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::QUADSHELL4, libMesh::QUADSHELL8, libMesh::Real, libMesh::SECOND, libMesh::SEVENTEENTH, libMesh::SEVENTH, libMesh::SIXTEENTH, libMesh::SIXTH, std::sqrt(), stroud_rule(), libMesh::TENTH, libMesh::THIRTEENTH, libMesh::TWELFTH, and wissmann_rule().

◆ init_3D()

void libMesh::QMonomial::init_3D ( const  type,
unsigned int  p_level 
)
overrideprivatevirtual

Initializes the 3D quadrature rule by filling the points and weights vectors with the appropriate values.

The order of the rule will be defined by the implementing class. Should not be pure virtual since a derived quadrature rule may only be defined in 1D. If not overridden, throws an error.

Note
The arguments should no longer be used for anything in derived classes, they are only maintained for backwards compatibility and will eventually be removed.

Reimplemented from libMesh::QBase.

Definition at line 28 of file quadrature_monomial_3D.C.

29 {
30 
31  switch (_type)
32  {
33  //---------------------------------------------
34  // Hex quadrature rules
35  case HEX8:
36  case HEX20:
37  case HEX27:
38  {
39  switch(get_order())
40  {
41 
42  // The CONSTANT/FIRST rule is the 1-point Gauss "product" rule...we fall
43  // through to the default case for this rule.
44 
45  case SECOND:
46  case THIRD:
47  {
48  // A degree 3, 6-point, "rotationally-symmetric" rule by
49  // Kim and Song, Comm. Korean Math. Soc vol. 13, no. 4, 1998, pp. 913-931.
50  //
51  // Warning: this rule contains points on the boundary of the reference
52  // element, and therefore may be unsuitable for some problems. The alternative
53  // would be a 2x2x2 Gauss product rule.
54  const Real data[1][4] =
55  {
56  {1.0L, 0.0L, 0.0L, Real(4)/3}
57  };
58 
59  const unsigned int rule_id[1] = {
60  1 // (x,0,0) -> 6 permutations
61  };
62 
63  _points.resize(6);
64  _weights.resize(6);
65 
66  kim_rule(data, rule_id, 1);
67  return;
68  } // end case SECOND,THIRD
69 
70  case FOURTH:
71  case FIFTH:
72  {
73  // A degree 5, 13-point rule by Stroud,
74  // AH Stroud, "Some Fifth Degree Integration Formulas for Symmetric Regions II.",
75  // Numerische Mathematik 9, pp. 460-468 (1967).
76  //
77  // This rule is provably minimal in the number of points. The equations given for
78  // the n-cube on pg. 466 of the paper for mu/gamma and gamma are wrong, at least for
79  // the n=3 case. The analytical values given here were computed by me [JWP] in Maple.
80 
81  // Convenient intermediate values.
82  const Real sqrt19 = Real(std::sqrt(19.L));
83  const Real tp = Real(std::sqrt(71440 + 6802*sqrt19));
84 
85  // Point data for permutations.
86  const Real eta = 0;
87 
88  const Real lambda = Real(std::sqrt(1919.L/3285.L - 148.L*sqrt19/3285.L + 4.L*tp/3285.L));
89  // 8.8030440669930978047737818209860e-01L;
90 
91  const Real xi = Real(-std::sqrt(1121.L/3285.L + 74.L*sqrt19/3285.L - 2.L*tp/3285.L));
92  // -4.9584817142571115281421242364290e-01L;
93 
94  const Real mu = Real(std::sqrt(1121.L/3285.L + 74.L*sqrt19/3285.L + 2.L*tp/3285.L));
95  // 7.9562142216409541542982482567580e-01L;
96 
97  const Real gamma = Real(std::sqrt(1919.L/3285.L - 148.L*sqrt19/3285.L - 4.L*tp/3285.L));
98  // 2.5293711744842581347389255929324e-02L;
99 
100  // Weights: the centroid weight is given analytically. Weight B (resp C) goes
101  // with the {lambda,xi} (resp {gamma,mu}) permutation. The single-precision
102  // results reported by Stroud are given for reference.
103 
104  const Real A = Real(32)/19;
105  // Stroud: 0.21052632 * 8.0 = 1.684210560;
106 
107  const Real B = Real(1) / (Real(260072)/133225 - 1520*sqrt19/133225 + (133 - 37*sqrt19)*tp/133225);
108  // 5.4498735127757671684690782180890e-01L; // Stroud: 0.068123420 * 8.0 = 0.544987360;
109 
110  const Real C = Real(1) / (Real(260072)/133225 - 1520*sqrt19/133225 - (133 - 37*sqrt19)*tp/133225);
111  // 5.0764422766979170420572375713840e-01L; // Stroud: 0.063455527 * 8.0 = 0.507644216;
112 
113  _points.resize(13);
114  _weights.resize(13);
115 
116  unsigned int c=0;
117 
118  // Point with weight A (origin)
119  _points[c] = Point(eta, eta, eta);
120  _weights[c++] = A;
121 
122  // Points with weight B
123  _points[c] = Point(lambda, xi, xi);
124  _weights[c++] = B;
125  _points[c] = -_points[c-1];
126  _weights[c++] = B;
127 
128  _points[c] = Point(xi, lambda, xi);
129  _weights[c++] = B;
130  _points[c] = -_points[c-1];
131  _weights[c++] = B;
132 
133  _points[c] = Point(xi, xi, lambda);
134  _weights[c++] = B;
135  _points[c] = -_points[c-1];
136  _weights[c++] = B;
137 
138  // Points with weight C
139  _points[c] = Point(mu, mu, gamma);
140  _weights[c++] = C;
141  _points[c] = -_points[c-1];
142  _weights[c++] = C;
143 
144  _points[c] = Point(mu, gamma, mu);
145  _weights[c++] = C;
146  _points[c] = -_points[c-1];
147  _weights[c++] = C;
148 
149  _points[c] = Point(gamma, mu, mu);
150  _weights[c++] = C;
151  _points[c] = -_points[c-1];
152  _weights[c++] = C;
153 
154  return;
155 
156 
157  // // A degree 5, 14-point, "rotationally-symmetric" rule by
158  // // Kim and Song, Comm. Korean Math. Soc vol. 13, no. 4, 1998, pp. 913-931.
159  // // Was also reported in Stroud's 1971 book.
160  // const Real data[2][4] =
161  // {
162  // {7.95822425754221463264548820476135e-01L, 0.00000000000000000000000000000000e+00L, 0.00000000000000000000000000000000e+00L, 8.86426592797783933518005540166204e-01L},
163  // {7.58786910639328146269034278112267e-01L, 7.58786910639328146269034278112267e-01L, 7.58786910639328146269034278112267e-01L, 3.35180055401662049861495844875346e-01L}
164  // };
165 
166  // const unsigned int rule_id[2] = {
167  // 1, // (x,0,0) -> 6 permutations
168  // 4 // (x,x,x) -> 8 permutations
169  // };
170 
171  // _points.resize(14);
172  // _weights.resize(14);
173 
174  // kim_rule(data, rule_id, 2);
175  // return;
176  } // end case FOURTH,FIFTH
177 
178 
179  case SIXTH:
180  case SEVENTH:
181  {
182  if (allow_rules_with_negative_weights)
183  {
184  // A degree 7, 31-point, "rotationally-symmetric" rule by
185  // Kim and Song, Comm. Korean Math. Soc vol. 13, no. 4, 1998, pp. 913-931.
186  // This rule contains a negative weight, so only use it if such type of
187  // rules are allowed.
188  const Real data[3][4] =
189  {
190  {Real(0.00000000000000000000000000000000e+00L), Real(0.00000000000000000000000000000000e+00L), Real(0.00000000000000000000000000000000e+00L), Real(-1.27536231884057971014492753623188e+00L)},
191  {Real(5.85540043769119907612630781744060e-01L), Real(0.00000000000000000000000000000000e+00L), Real(0.00000000000000000000000000000000e+00L), Real(8.71111111111111111111111111111111e-01L)},
192  {Real(6.94470135991704766602025803883310e-01L), Real(9.37161638568208038511047377665396e-01L), Real(4.15659267604065126239606672567031e-01L), Real(1.68695652173913043478260869565217e-01L)}
193  };
194 
195  const unsigned int rule_id[3] = {
196  0, // (0,0,0) -> 1 permutation
197  1, // (x,0,0) -> 6 permutations
198  6 // (x,y,z) -> 24 permutations
199  };
200 
201  _points.resize(31);
202  _weights.resize(31);
203 
204  kim_rule(data, rule_id, 3);
205  return;
206  } // end if (allow_rules_with_negative_weights)
207 
208 
209  // A degree 7, 34-point, "fully-symmetric" rule, first published in
210  // P.C. Hammer and A.H. Stroud, "Numerical Evaluation of Multiple Integrals II",
211  // Mathematical Tables and Other Aids to Computation, vol 12., no 64, 1958, pp. 272-280
212  //
213  // This rule happens to fall under the same general
214  // construction as the Kim rules, so we've re-used
215  // that code here. Stroud gives 16 digits for his rule,
216  // and this is the most accurate version I've found.
217  //
218  // For comparison, a SEVENTH-order Gauss product rule
219  // (which integrates tri-7th order polynomials) would
220  // have 4^3=64 points.
221  const Real
222  r = Real(std::sqrt(6.L/7.L)),
223  s = Real(std::sqrt((960.L - 3.L*std::sqrt(28798.L)) / 2726.L)),
224  t = Real(std::sqrt((960.L + 3.L*std::sqrt(28798.L)) / 2726.L)),
225  B1 = Real(8624)/29160,
226  B2 = Real(2744)/29160,
227  B3 = 8*(774*t*t - 230)/(9720*(t*t-s*s)),
228  B4 = 8*(230 - 774*s*s)/(9720*(t*t-s*s));
229 
230  const Real data[4][4] =
231  {
232  {r, 0, 0, B1},
233  {r, r, 0, B2},
234  {s, s, s, B3},
235  {t, t, t, B4}
236  };
237 
238  const unsigned int rule_id[4] = {
239  1, // (x,0,0) -> 6 permutations
240  2, // (x,x,0) -> 12 permutations
241  4, // (x,x,x) -> 8 permutations
242  4 // (x,x,x) -> 8 permutations
243  };
244 
245  _points.resize(34);
246  _weights.resize(34);
247 
248  kim_rule(data, rule_id, 4);
249  return;
250 
251 
252  // // A degree 7, 38-point, "rotationally-symmetric" rule by
253  // // Kim and Song, Comm. Korean Math. Soc vol. 13, no. 4, 1998, pp. 913-931.
254  // //
255  // // This rule is obviously inferior to the 34-point rule above...
256  // const Real data[3][4] =
257  //{
258  // {9.01687807821291289082811566285950e-01L, 0.00000000000000000000000000000000e+00L, 0.00000000000000000000000000000000e+00L, 2.95189738262622903181631100062774e-01L},
259  // {4.08372221499474674069588900002128e-01L, 4.08372221499474674069588900002128e-01L, 4.08372221499474674069588900002128e-01L, 4.04055417266200582425904380777126e-01L},
260  // {8.59523090201054193116477875786220e-01L, 8.59523090201054193116477875786220e-01L, 4.14735913727987720499709244748633e-01L, 1.24850759678944080062624098058597e-01L}
261  //};
262  //
263  // const unsigned int rule_id[3] = {
264  //1, // (x,0,0) -> 6 permutations
265  //4, // (x,x,x) -> 8 permutations
266  //5 // (x,x,z) -> 24 permutations
267  // };
268  //
269  // _points.resize(38);
270  // _weights.resize(38);
271  //
272  // kim_rule(data, rule_id, 3);
273  // return;
274  } // end case SIXTH,SEVENTH
275 
276  case EIGHTH:
277  {
278  // A degree 8, 47-point, "rotationally-symmetric" rule by
279  // Kim and Song, Comm. Korean Math. Soc vol. 13, no. 4, 1998, pp. 913-931.
280  //
281  // A EIGHTH-order Gauss product rule (which integrates tri-8th order polynomials)
282  // would have 5^3=125 points.
283  const Real data[5][4] =
284  {
285  {Real(0.00000000000000000000000000000000e+00L), Real(0.00000000000000000000000000000000e+00L), Real(0.00000000000000000000000000000000e+00L), Real(4.51903714875199690490763818699555e-01L)},
286  {Real(7.82460796435951590652813975429717e-01L), Real(0.00000000000000000000000000000000e+00L), Real(0.00000000000000000000000000000000e+00L), Real(2.99379177352338919703385618576171e-01L)},
287  {Real(4.88094669706366480526729301468686e-01L), Real(4.88094669706366480526729301468686e-01L), Real(4.88094669706366480526729301468686e-01L), Real(3.00876159371240019939698689791164e-01L)},
288  {Real(8.62218927661481188856422891110042e-01L), Real(8.62218927661481188856422891110042e-01L), Real(8.62218927661481188856422891110042e-01L), Real(4.94843255877038125738173175714853e-02L)},
289  {Real(2.81113909408341856058098281846420e-01L), Real(9.44196578292008195318687494773744e-01L), Real(6.97574833707236996779391729948984e-01L), Real(1.22872389222467338799199767122592e-01L)}
290  };
291 
292  const unsigned int rule_id[5] = {
293  0, // (0,0,0) -> 1 permutation
294  1, // (x,0,0) -> 6 permutations
295  4, // (x,x,x) -> 8 permutations
296  4, // (x,x,x) -> 8 permutations
297  6 // (x,y,z) -> 24 permutations
298  };
299 
300  _points.resize(47);
301  _weights.resize(47);
302 
303  kim_rule(data, rule_id, 5);
304  return;
305  } // end case EIGHTH
306 
307 
308  // By default: construct and use a Gauss quadrature rule
309  default:
310  {
311  // Break out and fall down into the default: case for the
312  // outer switch statement.
313  break;
314  }
315 
316  } // end switch(_order + 2*p)
317  } // end case HEX8/20/27
318 
319  libmesh_fallthrough();
320 
321  // By default: construct and use a Gauss quadrature rule
322  default:
323  {
324  QGauss gauss_rule(3, _order);
325  gauss_rule.init(_type, _p_level);
326 
327  // Swap points and weights with the about-to-be destroyed rule.
328  _points.swap (gauss_rule.get_points() );
329  _weights.swap(gauss_rule.get_weights());
330 
331  return;
332  }
333  } // end switch (_type)
334 }

References libMesh::QBase::_order, libMesh::QBase::_p_level, libMesh::QBase::_points, libMesh::QBase::_type, libMesh::QBase::_weights, A, libMesh::QBase::allow_rules_with_negative_weights, data, libMesh::EIGHTH, libMesh::FIFTH, libMesh::FOURTH, libMesh::QBase::get_order(), libMesh::QBase::get_points(), libMesh::QBase::get_weights(), libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::QBase::init(), kim_rule(), libMesh::Real, libMesh::SECOND, libMesh::SEVENTH, libMesh::SIXTH, std::sqrt(), and libMesh::THIRD.

◆ kim_rule()

void libMesh::QMonomial::kim_rule ( const Real  rule_data[][4],
const unsigned int *  rule_id,
const unsigned int  n_pts 
)
private

Rules from Kim and Song, Comm.

Korean Math. Soc vol. 13, no. 4, 1998, pp. 913-931. The rules are obtained by considering the group G^{rot} of rotations of the reference hex, and the invariant polynomials of this group.

In Kim and Song's rules, quadrature points are described by the following points and their unique permutations under the G^{rot} group:

0.) (0,0,0) ( 1 perm ) -> [0, 0, 0] 1.) (x,0,0) ( 6 perms) -> [x, 0, 0], [0, -x, 0], [-x, 0, 0], [0, x, 0], [0, 0, -x], [0, 0, x] 2.) (x,x,0) (12 perms) -> [x, x, 0], [x, -x, 0], [-x, -x, 0], [-x, x, 0], [x, 0, -x], [x, 0, x], [0, x, -x], [0, x, x], [0, -x, -x], [-x, 0, -x], [0, -x, x], [-x, 0, x] 3.) (x,y,0) (24 perms) -> [x, y, 0], [y, -x, 0], [-x, -y, 0], [-y, x, 0], [x, 0, -y], [x, -y, 0], [x, 0, y], [0, y, -x], [-x, y, 0], [0, y, x], [y, 0, -x], [0, -y, -x], [-y, 0, -x], [y, x, 0], [-y, -x, 0], [y, 0, x], [0, -y, x], [-y, 0, x], [-x, 0, y], [0, -x, -y], [0, -x, y], [-x, 0, -y], [0, x, y], [0, x, -y] 4.) (x,x,x) ( 8 perms) -> [x, x, x], [x, -x, x], [-x, -x, x], [-x, x, x], [x, x, -x], [x, -x, -x], [-x, x, -x], [-x, -x, -x] 5.) (x,x,z) (24 perms) -> [x, x, z], [x, -x, z], [-x, -x, z], [-x, x, z], [x, z, -x], [x, -x, -z], [x, -z, x], [z, x, -x], [-x, x, -z], [-z, x, x], [x, -z, -x], [-z, -x, -x], [-x, z, -x], [x, x, -z], [-x, -x, -z], [x, z, x], [z, -x, x], [-x, -z, x], [-x, z, x], [z, -x, -x], [-z, -x, x], [-x, -z, -x], [z, x, x], [-z, x, -x] 6.) (x,y,z) (24 perms) -> [x, y, z], [y, -x, z], [-x, -y, z], [-y, x, z], [x, z, -y], [x, -y, -z], [x, -z, y], [z, y, -x], [-x, y, -z], [-z, y, x], [y, -z, -x], [-z, -y, -x], [-y, z, -x], [y, x, -z], [-y, -x, -z], [y, z, x], [z, -y, x], [-y, -z, x], [-x, z, y], [z, -x, -y], [-z, -x, y], [-x, -z, -y], [z, x, y], [-z, x, -y]

Only two of Kim and Song's rules are particularly useful for FEM calculations: the degree 7, 38-point rule and their degree 8, 47-point rule. The others either contain negative weights or points outside the reference interval. The points and weights, to 32 digits, were obtained from: Ronald Cools' website (http://www.cs.kuleuven.ac.be/~nines/research/ecf/ecf.html) and the unique permutations of G^{rot} were computed by me [JWP] using Maple.

Definition at line 205 of file quadrature_monomial.C.

208 {
209  for (unsigned int i=0, c=0; i<n_pts; ++i)
210  {
211  const Real
212  x=rule_data[i][0],
213  y=rule_data[i][1],
214  z=rule_data[i][2],
215  wt=rule_data[i][3];
216 
217  switch(rule_id[i])
218  {
219  case 0: // (0,0,0) 1 permutation
220  {
221  _points[c] = Point( x, y, z); _weights[c++] = wt;
222 
223  break;
224  }
225  case 1: // (x,0,0) 6 permutations
226  {
227  _points[c] = Point( x, 0., 0.); _weights[c++] = wt;
228  _points[c] = Point(0., -x, 0.); _weights[c++] = wt;
229  _points[c] = Point(-x, 0., 0.); _weights[c++] = wt;
230  _points[c] = Point(0., x, 0.); _weights[c++] = wt;
231  _points[c] = Point(0., 0., -x); _weights[c++] = wt;
232  _points[c] = Point(0., 0., x); _weights[c++] = wt;
233 
234  break;
235  }
236  case 2: // (x,x,0) 12 permutations
237  {
238  _points[c] = Point( x, x, 0.); _weights[c++] = wt;
239  _points[c] = Point( x, -x, 0.); _weights[c++] = wt;
240  _points[c] = Point(-x, -x, 0.); _weights[c++] = wt;
241  _points[c] = Point(-x, x, 0.); _weights[c++] = wt;
242  _points[c] = Point( x, 0., -x); _weights[c++] = wt;
243  _points[c] = Point( x, 0., x); _weights[c++] = wt;
244  _points[c] = Point(0., x, -x); _weights[c++] = wt;
245  _points[c] = Point(0., x, x); _weights[c++] = wt;
246  _points[c] = Point(0., -x, -x); _weights[c++] = wt;
247  _points[c] = Point(-x, 0., -x); _weights[c++] = wt;
248  _points[c] = Point(0., -x, x); _weights[c++] = wt;
249  _points[c] = Point(-x, 0., x); _weights[c++] = wt;
250 
251  break;
252  }
253  case 3: // (x,y,0) 24 permutations
254  {
255  _points[c] = Point( x, y, 0.); _weights[c++] = wt;
256  _points[c] = Point( y, -x, 0.); _weights[c++] = wt;
257  _points[c] = Point(-x, -y, 0.); _weights[c++] = wt;
258  _points[c] = Point(-y, x, 0.); _weights[c++] = wt;
259  _points[c] = Point( x, 0., -y); _weights[c++] = wt;
260  _points[c] = Point( x, -y, 0.); _weights[c++] = wt;
261  _points[c] = Point( x, 0., y); _weights[c++] = wt;
262  _points[c] = Point(0., y, -x); _weights[c++] = wt;
263  _points[c] = Point(-x, y, 0.); _weights[c++] = wt;
264  _points[c] = Point(0., y, x); _weights[c++] = wt;
265  _points[c] = Point( y, 0., -x); _weights[c++] = wt;
266  _points[c] = Point(0., -y, -x); _weights[c++] = wt;
267  _points[c] = Point(-y, 0., -x); _weights[c++] = wt;
268  _points[c] = Point( y, x, 0.); _weights[c++] = wt;
269  _points[c] = Point(-y, -x, 0.); _weights[c++] = wt;
270  _points[c] = Point( y, 0., x); _weights[c++] = wt;
271  _points[c] = Point(0., -y, x); _weights[c++] = wt;
272  _points[c] = Point(-y, 0., x); _weights[c++] = wt;
273  _points[c] = Point(-x, 0., y); _weights[c++] = wt;
274  _points[c] = Point(0., -x, -y); _weights[c++] = wt;
275  _points[c] = Point(0., -x, y); _weights[c++] = wt;
276  _points[c] = Point(-x, 0., -y); _weights[c++] = wt;
277  _points[c] = Point(0., x, y); _weights[c++] = wt;
278  _points[c] = Point(0., x, -y); _weights[c++] = wt;
279 
280  break;
281  }
282  case 4: // (x,x,x) 8 permutations
283  {
284  _points[c] = Point( x, x, x); _weights[c++] = wt;
285  _points[c] = Point( x, -x, x); _weights[c++] = wt;
286  _points[c] = Point(-x, -x, x); _weights[c++] = wt;
287  _points[c] = Point(-x, x, x); _weights[c++] = wt;
288  _points[c] = Point( x, x, -x); _weights[c++] = wt;
289  _points[c] = Point( x, -x, -x); _weights[c++] = wt;
290  _points[c] = Point(-x, x, -x); _weights[c++] = wt;
291  _points[c] = Point(-x, -x, -x); _weights[c++] = wt;
292 
293  break;
294  }
295  case 5: // (x,x,z) 24 permutations
296  {
297  _points[c] = Point( x, x, z); _weights[c++] = wt;
298  _points[c] = Point( x, -x, z); _weights[c++] = wt;
299  _points[c] = Point(-x, -x, z); _weights[c++] = wt;
300  _points[c] = Point(-x, x, z); _weights[c++] = wt;
301  _points[c] = Point( x, z, -x); _weights[c++] = wt;
302  _points[c] = Point( x, -x, -z); _weights[c++] = wt;
303  _points[c] = Point( x, -z, x); _weights[c++] = wt;
304  _points[c] = Point( z, x, -x); _weights[c++] = wt;
305  _points[c] = Point(-x, x, -z); _weights[c++] = wt;
306  _points[c] = Point(-z, x, x); _weights[c++] = wt;
307  _points[c] = Point( x, -z, -x); _weights[c++] = wt;
308  _points[c] = Point(-z, -x, -x); _weights[c++] = wt;
309  _points[c] = Point(-x, z, -x); _weights[c++] = wt;
310  _points[c] = Point( x, x, -z); _weights[c++] = wt;
311  _points[c] = Point(-x, -x, -z); _weights[c++] = wt;
312  _points[c] = Point( x, z, x); _weights[c++] = wt;
313  _points[c] = Point( z, -x, x); _weights[c++] = wt;
314  _points[c] = Point(-x, -z, x); _weights[c++] = wt;
315  _points[c] = Point(-x, z, x); _weights[c++] = wt;
316  _points[c] = Point( z, -x, -x); _weights[c++] = wt;
317  _points[c] = Point(-z, -x, x); _weights[c++] = wt;
318  _points[c] = Point(-x, -z, -x); _weights[c++] = wt;
319  _points[c] = Point( z, x, x); _weights[c++] = wt;
320  _points[c] = Point(-z, x, -x); _weights[c++] = wt;
321 
322  break;
323  }
324  case 6: // (x,y,z) 24 permutations
325  {
326  _points[c] = Point( x, y, z); _weights[c++] = wt;
327  _points[c] = Point( y, -x, z); _weights[c++] = wt;
328  _points[c] = Point(-x, -y, z); _weights[c++] = wt;
329  _points[c] = Point(-y, x, z); _weights[c++] = wt;
330  _points[c] = Point( x, z, -y); _weights[c++] = wt;
331  _points[c] = Point( x, -y, -z); _weights[c++] = wt;
332  _points[c] = Point( x, -z, y); _weights[c++] = wt;
333  _points[c] = Point( z, y, -x); _weights[c++] = wt;
334  _points[c] = Point(-x, y, -z); _weights[c++] = wt;
335  _points[c] = Point(-z, y, x); _weights[c++] = wt;
336  _points[c] = Point( y, -z, -x); _weights[c++] = wt;
337  _points[c] = Point(-z, -y, -x); _weights[c++] = wt;
338  _points[c] = Point(-y, z, -x); _weights[c++] = wt;
339  _points[c] = Point( y, x, -z); _weights[c++] = wt;
340  _points[c] = Point(-y, -x, -z); _weights[c++] = wt;
341  _points[c] = Point( y, z, x); _weights[c++] = wt;
342  _points[c] = Point( z, -y, x); _weights[c++] = wt;
343  _points[c] = Point(-y, -z, x); _weights[c++] = wt;
344  _points[c] = Point(-x, z, y); _weights[c++] = wt;
345  _points[c] = Point( z, -x, -y); _weights[c++] = wt;
346  _points[c] = Point(-z, -x, y); _weights[c++] = wt;
347  _points[c] = Point(-x, -z, -y); _weights[c++] = wt;
348  _points[c] = Point( z, x, y); _weights[c++] = wt;
349  _points[c] = Point(-z, x, -y); _weights[c++] = wt;
350 
351  break;
352  }
353  default:
354  libmesh_error_msg("Unknown rule ID: " << rule_id[i] << "!");
355  } // end switch(rule_id[i])
356  }
357 }

References libMesh::QBase::_points, libMesh::QBase::_weights, and libMesh::Real.

Referenced by init_3D().

◆ n_objects()

static unsigned int libMesh::ReferenceCounter::n_objects ( )
inlinestaticinherited

Prints the number of outstanding (created, but not yet destroyed) objects.

Definition at line 83 of file reference_counter.h.

84  { return _n_objects; }

References libMesh::ReferenceCounter::_n_objects.

◆ n_points()

unsigned int libMesh::QBase::n_points ( ) const
inlineinherited
Returns
The number of points associated with the quadrature rule.

Definition at line 126 of file quadrature.h.

127  {
128  libmesh_assert (!_points.empty());
129  return cast_int<unsigned int>(_points.size());
130  }

References libMesh::QBase::_points, and libMesh::libmesh_assert().

Referenced by libMesh::ExactSolution::_compute_error(), assemble(), assemble_1D(), assemble_ellipticdg(), assemble_helmholtz(), assemble_poisson(), assemble_shell(), assemble_wave(), assembly_with_dg_fem_context(), AssemblyA0::boundary_assembly(), AssemblyF0::boundary_assembly(), AssemblyA1::boundary_assembly(), AssemblyF1::boundary_assembly(), AssemblyF2::boundary_assembly(), A2::boundary_assembly(), AssemblyA2::boundary_assembly(), A3::boundary_assembly(), F0::boundary_assembly(), Output0::boundary_assembly(), libMesh::System::calculate_norm(), libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), libMesh::QConical::conical_product_pyramid(), libMesh::QConical::conical_product_tet(), libMesh::QConical::conical_product_tri(), SecondOrderScalarSystemSecondOrderTimeSolverBase::damping_residual(), SecondOrderScalarSystemFirstOrderTimeSolverBase::damping_residual(), NavierSystem::element_constraint(), CoupledSystem::element_constraint(), PoissonSystem::element_postprocess(), LaplaceSystem::element_postprocess(), LaplaceQoI::element_qoi(), LaplaceQoI::element_qoi_derivative(), LaplaceSystem::element_qoi_derivative(), HeatSystem::element_qoi_derivative(), NavierSystem::element_time_derivative(), SolidSystem::element_time_derivative(), PoissonSystem::element_time_derivative(), LaplaceSystem::element_time_derivative(), CurlCurlSystem::element_time_derivative(), ElasticitySystem::element_time_derivative(), L2System::element_time_derivative(), CoupledSystem::element_time_derivative(), FirstOrderScalarSystemBase::element_time_derivative(), SecondOrderScalarSystemFirstOrderTimeSolverBase::element_time_derivative(), libMesh::RBEIMConstruction::enrich_RB_space(), A0::interior_assembly(), B::interior_assembly(), M0::interior_assembly(), A1::interior_assembly(), AssemblyA0::interior_assembly(), EIM_IP_assembly::interior_assembly(), AcousticsInnerProduct::interior_assembly(), AssemblyA1::interior_assembly(), A2::interior_assembly(), AssemblyA2::interior_assembly(), EIM_F::interior_assembly(), F0::interior_assembly(), OutputAssembly::interior_assembly(), InnerProductAssembly::interior_assembly(), AssemblyEIM::interior_assembly(), AssemblyF0::interior_assembly(), AssemblyF1::interior_assembly(), Ex6InnerProduct::interior_assembly(), Ex6EIMInnerProduct::interior_assembly(), LaplaceYoung::jacobian(), NavierSystem::mass_residual(), ElasticitySystem::mass_residual(), libMesh::FEMPhysics::mass_residual(), FirstOrderScalarSystemBase::mass_residual(), SecondOrderScalarSystemSecondOrderTimeSolverBase::mass_residual(), SecondOrderScalarSystemFirstOrderTimeSolverBase::mass_residual(), libMesh::QBase::print_info(), LaplaceYoung::residual(), LaplaceSystem::side_constraint(), LaplaceSystem::side_postprocess(), CoupledSystemQoI::side_qoi(), CoupledSystemQoI::side_qoi_derivative(), LaplaceSystem::side_qoi_derivative(), SolidSystem::side_time_derivative(), CurlCurlSystem::side_time_derivative(), ElasticitySystem::side_time_derivative(), libMesh::QBase::tensor_product_hex(), libMesh::QBase::tensor_product_prism(), libMesh::QBase::tensor_product_quad(), and libMesh::RBEIMConstruction::truth_solve().

◆ operator=() [1/2]

QMonomial& libMesh::QMonomial::operator= ( const QMonomial &  )
default

◆ operator=() [2/2]

QMonomial& libMesh::QMonomial::operator= ( QMonomial &&  )
default

◆ print_info() [1/2]

void libMesh::QBase::print_info ( std::ostream &  os = libMesh::out) const
inherited

Prints information relevant to the quadrature rule, by default to libMesh::out.

Definition at line 37 of file quadrature.C.

38 {
39  libmesh_assert(!_points.empty());
40  libmesh_assert(!_weights.empty());
41 
42  Real summed_weights=0;
43  os << "N_Q_Points=" << this->n_points() << std::endl << std::endl;
44  for (auto qpoint: index_range(_points))
45  {
46  os << " Point " << qpoint << ":\n"
47  << " "
48  << _points[qpoint]
49  << "\n Weight:\n "
50  << " w=" << _weights[qpoint] << "\n" << std::endl;
51 
52  summed_weights += _weights[qpoint];
53  }
54  os << "Summed Weights: " << summed_weights << std::endl;
55 }

References libMesh::QBase::_points, libMesh::QBase::_weights, libMesh::index_range(), libMesh::libmesh_assert(), libMesh::QBase::n_points(), and libMesh::Real.

Referenced by libMesh::operator<<().

◆ print_info() [2/2]

void libMesh::ReferenceCounter::print_info ( std::ostream &  out = libMesh::out)
staticinherited

Prints the reference information, by default to libMesh::out.

Definition at line 87 of file reference_counter.C.

88 {
89  if (_enable_print_counter)
90  out_stream << ReferenceCounter::get_info();
91 }

References libMesh::ReferenceCounter::_enable_print_counter, and libMesh::ReferenceCounter::get_info().

◆ qp()

Point libMesh::QBase::qp ( const unsigned int  i) const
inlineinherited
Returns
The \( i^{th} \) quadrature point in reference element space.

Definition at line 164 of file quadrature.h.

165  {
166  libmesh_assert_less (i, _points.size());
167  return _points[i];
168  }

References libMesh::QBase::_points.

Referenced by libMesh::QConical::conical_product_pyramid(), libMesh::QConical::conical_product_tet(), libMesh::QConical::conical_product_tri(), libMesh::QBase::tensor_product_hex(), libMesh::QBase::tensor_product_prism(), and libMesh::QBase::tensor_product_quad().

◆ scale()

void libMesh::QBase::scale ( std::pair< Real, Real >  old_range,
std::pair< Real, Real >  new_range 
)
inherited

Maps the points of a 1D quadrature rule defined by "old_range" to another 1D interval defined by "new_range" and scales the weights accordingly.

Definition at line 137 of file quadrature.C.

139 {
140  // Make sure we are in 1D
141  libmesh_assert_equal_to (_dim, 1);
142 
143  Real
144  h_new = new_range.second - new_range.first,
145  h_old = old_range.second - old_range.first;
146 
147  // Make sure that we have sane ranges
148  libmesh_assert_greater (h_new, 0.);
149  libmesh_assert_greater (h_old, 0.);
150 
151  // Make sure there are some points
152  libmesh_assert_greater (_points.size(), 0);
153 
154  // Compute the scale factor
155  Real scfact = h_new/h_old;
156 
157  // We're mapping from old_range -> new_range
158  for (auto i : index_range(_points))
159  {
160  _points[i](0) = new_range.first +
161  (_points[i](0) - old_range.first) * scfact;
162 
163  // Scale the weights
164  _weights[i] *= scfact;
165  }
166 }

References libMesh::QBase::_dim, libMesh::QBase::_points, libMesh::QBase::_weights, libMesh::index_range(), and libMesh::Real.

Referenced by libMesh::QConical::conical_product_tet(), and libMesh::QConical::conical_product_tri().

◆ shapes_need_reinit()

virtual bool libMesh::QBase::shapes_need_reinit ( )
inlinevirtualinherited
Returns
true if the shape functions need to be recalculated, false otherwise.

This may be required if the number of quadrature points or their position changes.

Definition at line 239 of file quadrature.h.

239 { return false; }

◆ stroud_rule()

void libMesh::QMonomial::stroud_rule ( const Real  rule_data[][3],
const unsigned int *  rule_symmetry,
const unsigned int  n_pts 
)
private

Stroud's rules for quads and hexes can have one of several different types of symmetry.

The rule_symmetry array describes how the different lines of the rule_data array are to be applied. The different rule_symmetry possibilities are: 0) Origin or single-point: (x,y) Fully-symmetric, 3 cases: 1) (x,y) -> (x,y), (-x,y), (x,-y), (-x,-y) (y,x), (-y,x), (y,-x), (-y,-x) 2) (x,x) -> (x,x), (-x,x), (x,-x), (-x,-x) 3) (x,0) -> (x,0), (-x,0), (0, x), ( 0,-x) 4) Rotational Invariant, (x,y) -> (x,y), (-x,-y), (-y, x), (y,-x) 5) Partial Symmetry, (x,y) -> (x,y), (-x, y) [x!=0] 6) Rectangular Symmetry, (x,y) -> (x,y), (-x, y), (-x,-y), (x,-y) 7) Central Symmetry, (0,y) -> (0,y), ( 0,-y)

Not all rules with these symmetries are due to Stroud, however, his book is probably the most frequently-cited compendium of quadrature rules and later authors certainly built upon his work.

Definition at line 57 of file quadrature_monomial.C.

60 {
61  for (unsigned int i=0, c=0; i<n_pts; ++i)
62  {
63  const Real
64  x=rule_data[i][0],
65  y=rule_data[i][1],
66  wt=rule_data[i][2];
67 
68  switch(rule_symmetry[i])
69  {
70  case 0: // Single point (no symmetry)
71  {
72  _points[c] = Point( x, y);
73  _weights[c++] = wt;
74 
75  break;
76  }
77  case 1: // Fully-symmetric (x,y)
78  {
79  _points[c] = Point( x, y);
80  _weights[c++] = wt;
81 
82  _points[c] = Point(-x, y);
83  _weights[c++] = wt;
84 
85  _points[c] = Point( x,-y);
86  _weights[c++] = wt;
87 
88  _points[c] = Point(-x,-y);
89  _weights[c++] = wt;
90 
91  _points[c] = Point( y, x);
92  _weights[c++] = wt;
93 
94  _points[c] = Point(-y, x);
95  _weights[c++] = wt;
96 
97  _points[c] = Point( y,-x);
98  _weights[c++] = wt;
99 
100  _points[c] = Point(-y,-x);
101  _weights[c++] = wt;
102 
103  break;
104  }
105  case 2: // Fully-symmetric (x,x)
106  {
107  _points[c] = Point( x, x);
108  _weights[c++] = wt;
109 
110  _points[c] = Point(-x, x);
111  _weights[c++] = wt;
112 
113  _points[c] = Point( x,-x);
114  _weights[c++] = wt;
115 
116  _points[c] = Point(-x,-x);
117  _weights[c++] = wt;
118 
119  break;
120  }
121  case 3: // Fully-symmetric (x,0)
122  {
123  libmesh_assert_equal_to (y, 0.0);
124 
125  _points[c] = Point( x,0.);
126  _weights[c++] = wt;
127 
128  _points[c] = Point(-x,0.);
129  _weights[c++] = wt;
130 
131  _points[c] = Point(0., x);
132  _weights[c++] = wt;
133 
134  _points[c] = Point(0.,-x);
135  _weights[c++] = wt;
136 
137  break;
138  }
139  case 4: // Rotational invariant
140  {
141  _points[c] = Point( x, y);
142  _weights[c++] = wt;
143 
144  _points[c] = Point(-x,-y);
145  _weights[c++] = wt;
146 
147  _points[c] = Point(-y, x);
148  _weights[c++] = wt;
149 
150  _points[c] = Point( y,-x);
151  _weights[c++] = wt;
152 
153  break;
154  }
155  case 5: // Partial symmetry (Wissman's rules)
156  {
157  libmesh_assert_not_equal_to (x, 0.0);
158 
159  _points[c] = Point( x, y);
160  _weights[c++] = wt;
161 
162  _points[c] = Point(-x, y);
163  _weights[c++] = wt;
164 
165  break;
166  }
167  case 6: // Rectangular symmetry
168  {
169  _points[c] = Point( x, y);
170  _weights[c++] = wt;
171 
172  _points[c] = Point(-x, y);
173  _weights[c++] = wt;
174 
175  _points[c] = Point(-x,-y);
176  _weights[c++] = wt;
177 
178  _points[c] = Point( x,-y);
179  _weights[c++] = wt;
180 
181  break;
182  }
183  case 7: // Central symmetry
184  {
185  libmesh_assert_equal_to (x, 0.0);
186  libmesh_assert_not_equal_to (y, 0.0);
187 
188  _points[c] = Point(0., y);
189  _weights[c++] = wt;
190 
191  _points[c] = Point(0.,-y);
192  _weights[c++] = wt;
193 
194  break;
195  }
196  default:
197  libmesh_error_msg("Unknown symmetry!");
198  } // end switch(rule_symmetry[i])
199  }
200 }

References libMesh::QBase::_points, libMesh::QBase::_weights, and libMesh::Real.

Referenced by init_2D().

◆ tensor_product_hex()

void libMesh::QBase::tensor_product_hex ( const QBase &  q1D)
protectedinherited

Computes the tensor product quadrature rule [q1D x q1D x q1D] from the 1D rule q1D.

Used in the init_3D routines for hexahedral element types.

Definition at line 198 of file quadrature.C.

199 {
200  const unsigned int np = q1D.n_points();
201 
202  _points.resize(np * np * np);
203 
204  _weights.resize(np * np * np);
205 
206  unsigned int q=0;
207 
208  for (unsigned int k=0; k<np; k++)
209  for (unsigned int j=0; j<np; j++)
210  for (unsigned int i=0; i<np; i++)
211  {
212  _points[q](0) = q1D.qp(i)(0);
213  _points[q](1) = q1D.qp(j)(0);
214  _points[q](2) = q1D.qp(k)(0);
215 
216  _weights[q] = q1D.w(i) * q1D.w(j) * q1D.w(k);
217 
218  q++;
219  }
220 }

References libMesh::QBase::_points, libMesh::QBase::_weights, libMesh::QBase::n_points(), libMesh::QBase::qp(), and libMesh::QBase::w().

Referenced by libMesh::QGaussLobatto::init_3D(), libMesh::QTrap::init_3D(), libMesh::QGrid::init_3D(), libMesh::QSimpson::init_3D(), and libMesh::QGauss::init_3D().

◆ tensor_product_prism()

void libMesh::QBase::tensor_product_prism ( const QBase &  q1D,
const QBase &  q2D 
)
protectedinherited

Computes the tensor product of a 1D quadrature rule and a 2D quadrature rule.

Used in the init_3D routines for prismatic element types.

Definition at line 225 of file quadrature.C.

226 {
227  const unsigned int n_points1D = q1D.n_points();
228  const unsigned int n_points2D = q2D.n_points();
229 
230  _points.resize (n_points1D * n_points2D);
231  _weights.resize (n_points1D * n_points2D);
232 
233  unsigned int q=0;
234 
235  for (unsigned int j=0; j<n_points1D; j++)
236  for (unsigned int i=0; i<n_points2D; i++)
237  {
238  _points[q](0) = q2D.qp(i)(0);
239  _points[q](1) = q2D.qp(i)(1);
240  _points[q](2) = q1D.qp(j)(0);
241 
242  _weights[q] = q2D.w(i) * q1D.w(j);
243 
244  q++;
245  }
246 
247 }

References libMesh::QBase::_points, libMesh::QBase::_weights, libMesh::QBase::n_points(), libMesh::QBase::qp(), and libMesh::QBase::w().

Referenced by libMesh::QGrid::init_3D(), libMesh::QTrap::init_3D(), libMesh::QSimpson::init_3D(), and libMesh::QGauss::init_3D().

◆ tensor_product_quad()

void libMesh::QBase::tensor_product_quad ( const QBase &  q1D)
protectedinherited

Constructs a 2D rule from the tensor product of q1D with itself.

Used in the init_2D() routines for quadrilateral element types.

Definition at line 171 of file quadrature.C.

172 {
173 
174  const unsigned int np = q1D.n_points();
175 
176  _points.resize(np * np);
177 
178  _weights.resize(np * np);
179 
180  unsigned int q=0;
181 
182  for (unsigned int j=0; j<np; j++)
183  for (unsigned int i=0; i<np; i++)
184  {
185  _points[q](0) = q1D.qp(i)(0);
186  _points[q](1) = q1D.qp(j)(0);
187 
188  _weights[q] = q1D.w(i)*q1D.w(j);
189 
190  q++;
191  }
192 }

References libMesh::QBase::_points, libMesh::QBase::_weights, libMesh::QBase::n_points(), libMesh::QBase::qp(), and libMesh::QBase::w().

Referenced by libMesh::QGaussLobatto::init_2D(), libMesh::QTrap::init_2D(), libMesh::QGrid::init_2D(), libMesh::QSimpson::init_2D(), and libMesh::QGauss::init_2D().

◆ type()

QuadratureType libMesh::QMonomial::type ( ) const
overridevirtual
Returns
QMONOMIAL.

Implements libMesh::QBase.

Definition at line 32 of file quadrature_monomial.C.

33 {
34  return QMONOMIAL;
35 }

References libMesh::QMONOMIAL.

◆ w()

Real libMesh::QBase::w ( const unsigned int  i) const
inlineinherited
Returns
The \( i^{th} \) quadrature weight.

Definition at line 173 of file quadrature.h.

174  {
175  libmesh_assert_less (i, _weights.size());
176  return _weights[i];
177  }

References libMesh::QBase::_weights.

Referenced by libMesh::QConical::conical_product_pyramid(), libMesh::QConical::conical_product_tet(), libMesh::QConical::conical_product_tri(), libMesh::QBase::tensor_product_hex(), libMesh::QBase::tensor_product_prism(), and libMesh::QBase::tensor_product_quad().

◆ wissmann_rule()

void libMesh::QMonomial::wissmann_rule ( const Real  rule_data[][3],
const unsigned int  n_pts 
)
private

Wissmann published three interesting "partially symmetric" rules for integrating degree 4, 6, and 8 polynomials exactly on QUADs.

These rules have all positive weights, all points inside the reference element, and have fewer points than tensor-product rules of equivalent order, making them superior to those rules for monomial bases.

J. W. Wissman and T. Becker, Partially symmetric cubature formulas for even degrees of exactness, SIAM J. Numer. Anal. 23 (1986), 676–685.

Definition at line 37 of file quadrature_monomial.C.

39 {
40  for (unsigned int i=0, c=0; i<n_pts; ++i)
41  {
42  _points[c] = Point( rule_data[i][0], rule_data[i][1] );
43  _weights[c++] = rule_data[i][2];
44 
45  // This may be an (x1,x2) -> (-x1,x2) point, in which case
46  // we will also generate the mirror point using the same weight.
47  if (rule_data[i][0] != static_cast<Real>(0.0))
48  {
49  _points[c] = Point( -rule_data[i][0], rule_data[i][1] );
50  _weights[c++] = rule_data[i][2];
51  }
52  }
53 }

References libMesh::QBase::_points, and libMesh::QBase::_weights.

Referenced by init_2D().

Member Data Documentation

◆ _counts

ReferenceCounter::Counts libMesh::ReferenceCounter::_counts
staticprotectedinherited

Actually holds the data.

Definition at line 122 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::get_info(), libMesh::ReferenceCounter::increment_constructor_count(), and libMesh::ReferenceCounter::increment_destructor_count().

◆ _dim

unsigned int libMesh::QBase::_dim
protectedinherited

The spatial dimension of the quadrature rule.

Definition at line 339 of file quadrature.h.

Referenced by libMesh::QBase::build(), libMesh::QBase::get_dim(), libMesh::QBase::init(), libMesh::QGaussLobatto::QGaussLobatto(), libMesh::QJacobi::QJacobi(), libMesh::QNodal::QNodal(), libMesh::QSimpson::QSimpson(), libMesh::QTrap::QTrap(), and libMesh::QBase::scale().

◆ _enable_print_counter

bool libMesh::ReferenceCounter::_enable_print_counter = true
staticprotectedinherited

Flag to control whether reference count information is printed when print_info is called.

Definition at line 141 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::disable_print_counter_info(), libMesh::ReferenceCounter::enable_print_counter_info(), and libMesh::ReferenceCounter::print_info().

◆ _mutex

Threads::spin_mutex libMesh::ReferenceCounter::_mutex
staticprotectedinherited

Mutual exclusion object to enable thread-safe reference counting.

Definition at line 135 of file reference_counter.h.

◆ _n_objects

Threads::atomic< unsigned int > libMesh::ReferenceCounter::_n_objects
staticprotectedinherited

The number of objects.

Print the reference count information when the number returns to 0.

Definition at line 130 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::n_objects(), libMesh::ReferenceCounter::ReferenceCounter(), and libMesh::ReferenceCounter::~ReferenceCounter().

◆ _order

Order libMesh::QBase::_order
protectedinherited

The polynomial order which the quadrature rule is capable of integrating exactly.

Definition at line 345 of file quadrature.h.

Referenced by libMesh::QBase::build(), libMesh::QBase::get_order(), libMesh::QClough::init_1D(), libMesh::QGaussLobatto::init_1D(), libMesh::QGrid::init_1D(), libMesh::QGauss::init_1D(), libMesh::QNodal::init_1D(), libMesh::QJacobi::init_1D(), init_1D(), libMesh::QClough::init_2D(), libMesh::QGaussLobatto::init_2D(), libMesh::QGrid::init_2D(), libMesh::QGauss::init_2D(), libMesh::QNodal::init_2D(), init_2D(), libMesh::QGaussLobatto::init_3D(), libMesh::QGrid::init_3D(), libMesh::QGauss::init_3D(), libMesh::QNodal::init_3D(), and init_3D().

◆ _p_level

unsigned int libMesh::QBase::_p_level
protectedinherited

The p-level of the element for which the current values have been computed.

Definition at line 357 of file quadrature.h.

Referenced by libMesh::QBase::get_order(), libMesh::QBase::get_p_level(), libMesh::QBase::init(), libMesh::QClough::init_1D(), libMesh::QNodal::init_1D(), init_1D(), libMesh::QClough::init_2D(), libMesh::QGaussLobatto::init_2D(), libMesh::QGauss::init_2D(), libMesh::QNodal::init_2D(), init_2D(), libMesh::QGaussLobatto::init_3D(), libMesh::QGauss::init_3D(), libMesh::QNodal::init_3D(), and init_3D().

◆ _points

std::vector<Point> libMesh::QBase::_points
protectedinherited

The locations of the quadrature points in reference element space.

Definition at line 363 of file quadrature.h.

Referenced by libMesh::QConical::conical_product_pyramid(), libMesh::QConical::conical_product_tet(), libMesh::QConical::conical_product_tri(), libMesh::QGauss::dunavant_rule(), libMesh::QGauss::dunavant_rule2(), libMesh::QBase::get_points(), libMesh::QGrundmann_Moller::gm_rule(), libMesh::QBase::init_0D(), libMesh::QClough::init_1D(), libMesh::QGaussLobatto::init_1D(), libMesh::QTrap::init_1D(), libMesh::QConical::init_1D(), libMesh::QGrid::init_1D(), libMesh::QGauss::init_1D(), libMesh::QSimpson::init_1D(), libMesh::QNodal::init_1D(), libMesh::QJacobi::init_1D(), init_1D(), libMesh::QGrundmann_Moller::init_1D(), libMesh::QClough::init_2D(), libMesh::QTrap::init_2D(), libMesh::QGrid::init_2D(), libMesh::QSimpson::init_2D(), libMesh::QGauss::init_2D(), libMesh::QNodal::init_2D(), init_2D(), libMesh::QGrid::init_3D(), libMesh::QTrap::init_3D(), libMesh::QSimpson::init_3D(), libMesh::QGauss::init_3D(), libMesh::QNodal::init_3D(), init_3D(), libMesh::QGrundmann_Moller::init_3D(), libMesh::QGauss::keast_rule(), kim_rule(), libMesh::QBase::n_points(), libMesh::QBase::print_info(), libMesh::QBase::qp(), libMesh::QBase::scale(), stroud_rule(), libMesh::QBase::tensor_product_hex(), libMesh::QBase::tensor_product_prism(), libMesh::QBase::tensor_product_quad(), and wissmann_rule().

◆ _type

ElemType libMesh::QBase::_type
protectedinherited

The type of element for which the current values have been computed.

Definition at line 351 of file quadrature.h.

Referenced by libMesh::QBase::get_elem_type(), libMesh::QBase::init(), libMesh::QClough::init_1D(), libMesh::QNodal::init_1D(), init_1D(), libMesh::QClough::init_2D(), libMesh::QGaussLobatto::init_2D(), libMesh::QConical::init_2D(), libMesh::QTrap::init_2D(), libMesh::QGrid::init_2D(), libMesh::QGauss::init_2D(), libMesh::QSimpson::init_2D(), libMesh::QNodal::init_2D(), init_2D(), libMesh::QGrundmann_Moller::init_2D(), libMesh::QClough::init_3D(), libMesh::QGaussLobatto::init_3D(), libMesh::QGrid::init_3D(), libMesh::QConical::init_3D(), libMesh::QTrap::init_3D(), libMesh::QGauss::init_3D(), libMesh::QSimpson::init_3D(), libMesh::QNodal::init_3D(), init_3D(), and libMesh::QGrundmann_Moller::init_3D().

◆ _weights

std::vector<Real> libMesh::QBase::_weights
protectedinherited

The quadrature weights.

The order of the weights matches the ordering of the _points vector.

Definition at line 369 of file quadrature.h.

Referenced by libMesh::QConical::conical_product_pyramid(), libMesh::QConical::conical_product_tet(), libMesh::QConical::conical_product_tri(), libMesh::QGauss::dunavant_rule(), libMesh::QGauss::dunavant_rule2(), libMesh::QBase::get_weights(), libMesh::QGrundmann_Moller::gm_rule(), libMesh::QBase::init_0D(), libMesh::QClough::init_1D(), libMesh::QGaussLobatto::init_1D(), libMesh::QConical::init_1D(), libMesh::QTrap::init_1D(), libMesh::QGrid::init_1D(), libMesh::QGauss::init_1D(), libMesh::QSimpson::init_1D(), libMesh::QNodal::init_1D(), libMesh::QJacobi::init_1D(), init_1D(), libMesh::QGrundmann_Moller::init_1D(), libMesh::QClough::init_2D(), libMesh::QTrap::init_2D(), libMesh::QGrid::init_2D(), libMesh::QGauss::init_2D(), libMesh::QSimpson::init_2D(), libMesh::QNodal::init_2D(), init_2D(), libMesh::QGrid::init_3D(), libMesh::QTrap::init_3D(), libMesh::QSimpson::init_3D(), libMesh::QGauss::init_3D(), libMesh::QNodal::init_3D(), init_3D(), libMesh::QGrundmann_Moller::init_3D(), libMesh::QGauss::keast_rule(), kim_rule(), libMesh::QBase::print_info(), libMesh::QBase::scale(), stroud_rule(), libMesh::QBase::tensor_product_hex(), libMesh::QBase::tensor_product_prism(), libMesh::QBase::tensor_product_quad(), libMesh::QBase::w(), and wissmann_rule().

◆ allow_rules_with_negative_weights

bool libMesh::QBase::allow_rules_with_negative_weights
inherited

Flag (default true) controlling the use of quadrature rules with negative weights.

Set this to false to require rules with all positive weights.

Rules with negative weights can be unsuitable for some problems. For example, it is possible for a rule with negative weights to obtain a negative result when integrating a positive function.

A particular example: if rules with negative weights are not allowed, a request for TET,THIRD (5 points) will return the TET,FIFTH (14 points) rule instead, nearly tripling the computational effort required!

Definition at line 254 of file quadrature.h.

Referenced by libMesh::QGrundmann_Moller::init_2D(), libMesh::QGauss::init_3D(), init_3D(), and libMesh::QGrundmann_Moller::init_3D().

The documentation for this class was generated from the following files:
libMesh::HEX20
Definition: enum_elem_type.h:48
libMesh::QBase::_p_level
unsigned int _p_level
The p-level of the element for which the current values have been computed.
Definition: quadrature.h:357
libMesh::QMonomial::stroud_rule
void stroud_rule(const Real rule_data[][3], const unsigned int *rule_symmetry, const unsigned int n_pts)
Stroud's rules for quads and hexes can have one of several different types of symmetry.
Definition: quadrature_monomial.C:57
libMesh::QMONOMIAL
Definition: enum_quadrature_type.h:41
libMesh::QCLOUGH
Definition: enum_quadrature_type.h:44
libMesh::HEX8
Definition: enum_elem_type.h:47
libMesh::SIXTH
Definition: enum_order.h:47
libMesh::QBase::n_points
unsigned int n_points() const
Definition: quadrature.h:126
libMesh::FOURTEENTH
Definition: enum_order.h:55
libMesh::index_range
IntRange< std::size_t > index_range(const std::vector< T > &vec)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:106
libMesh::QMonomial::kim_rule
void kim_rule(const Real rule_data[][4], const unsigned int *rule_id, const unsigned int n_pts)
Rules from Kim and Song, Comm.
Definition: quadrature_monomial.C:205
B
Definition: assembly.h:38
libMesh::ReferenceCounter::_counts
static Counts _counts
Actually holds the data.
Definition: reference_counter.h:122
std::sqrt
MetaPhysicL::DualNumber< T, D > sqrt(const MetaPhysicL::DualNumber< T, D > &in)
libMesh::QBase::init_1D
virtual void init_1D(const ElemType type=INVALID_ELEM, unsigned int p_level=0)=0
Initializes the 1D quadrature rule by filling the points and weights vectors with the appropriate val...
libMesh::ReferenceCounter::_n_objects
static Threads::atomic< unsigned int > _n_objects
The number of objects.
Definition: reference_counter.h:130
libMesh::FIFTH
Definition: enum_order.h:46
libMesh::THIRTEENTH
Definition: enum_order.h:54
libMesh::SECOND
Definition: enum_order.h:43
libMesh::QMonomial::wissmann_rule
void wissmann_rule(const Real rule_data[][3], const unsigned int n_pts)
Wissmann published three interesting "partially symmetric" rules for integrating degree 4,...
Definition: quadrature_monomial.C:37
libMesh::QGAUSS
Definition: enum_quadrature_type.h:34
libMesh::ReferenceCounter::get_info
static std::string get_info()
Gets a string containing the reference information.
Definition: reference_counter.C:47
libMesh::QBase::init_3D
virtual void init_3D(const ElemType type=INVALID_ELEM, unsigned int p_level=0)
Initializes the 3D quadrature rule by filling the points and weights vectors with the appropriate val...
Definition: quadrature.C:130
libMesh::QBase::init
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:59
libMesh::TWELFTH
Definition: enum_order.h:53
dim
unsigned int dim
Definition: adaptivity_ex3.C:113
libMesh::QBase::build
static std::unique_ptr< QBase > build(const std::string &name, const unsigned int dim, const Order order=INVALID_ORDER)
Builds a specific quadrature rule based on the name string.
Definition: quadrature_build.C:41
libMesh::FORTYTHIRD
Definition: enum_order.h:84
libMesh::QSIMPSON
Definition: enum_quadrature_type.h:37
libMesh::QBase::_order
Order _order
The polynomial order which the quadrature rule is capable of integrating exactly.
Definition: quadrature.h:345
libMesh::libmesh_assert
libmesh_assert(ctx)
libMesh::ELEVENTH
Definition: enum_order.h:52
libMesh::HEX27
Definition: enum_elem_type.h:49
libMesh::TWENTYTHIRD
Definition: enum_order.h:64
libMesh::QTRAP
Definition: enum_quadrature_type.h:38
libMesh::QNODAL
Definition: enum_quadrature_type.h:46
libMesh::Threads::spin_mtx
spin_mutex spin_mtx
A convenient spin mutex object which can be used for obtaining locks.
Definition: threads.C:29
libMesh::QUAD4
Definition: enum_elem_type.h:41
libMesh::QBase::QBase
QBase(unsigned int dim, Order order=INVALID_ORDER)
Constructor.
Definition: quadrature.C:27
libMesh::QBase::_type
ElemType _type
The type of element for which the current values have been computed.
Definition: quadrature.h:351
libMesh::EIGHTH
Definition: enum_order.h:49
A
static PetscErrorCode Mat * A
Definition: petscdmlibmeshimpl.C:1026
libMesh::TENTH
Definition: enum_order.h:51
libMesh::FOURTH
Definition: enum_order.h:45
libMesh::QBase::_weights
std::vector< Real > _weights
The quadrature weights.
Definition: quadrature.h:369
libMesh::QGAUSS_LOBATTO
Definition: enum_quadrature_type.h:43
libMesh::QBase::type
virtual QuadratureType type() const =0
libMesh::QJACOBI_1_0
Definition: enum_quadrature_type.h:35
libMesh::SEVENTH
Definition: enum_order.h:48
libMesh::QUADSHELL8
Definition: enum_elem_type.h:73
libMesh::QUADSHELL4
Definition: enum_elem_type.h:72
libMesh::FIFTEENTH
Definition: enum_order.h:56
libMesh::QBase::allow_rules_with_negative_weights
bool allow_rules_with_negative_weights
Flag (default true) controlling the use of quadrature rules with negative weights.
Definition: quadrature.h:254
libMesh::QJACOBI_2_0
Definition: enum_quadrature_type.h:36
data
IterBase * data
Ideally this private member data should have protected access.
Definition: variant_filter_iterator.h:337
libMesh::SIXTEENTH
Definition: enum_order.h:57
libMesh::THIRD
Definition: enum_order.h:44
libMesh::QUAD9
Definition: enum_elem_type.h:43
libMesh::NINTH
Definition: enum_order.h:50
libMesh::Real
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Definition: libmesh_common.h:121
libMesh::QBase::get_order
Order get_order() const
Definition: quadrature.h:211
libMesh::QCONICAL
Definition: enum_quadrature_type.h:42
libMesh::ReferenceCounter::_enable_print_counter
static bool _enable_print_counter
Flag to control whether reference count information is printed when print_info is called.
Definition: reference_counter.h:141
libMesh::QBase::init_0D
virtual void init_0D(const ElemType type=INVALID_ELEM, unsigned int p_level=0)
Initializes the 0D quadrature rule by filling the points and weights vectors with the appropriate val...
Definition: quadrature.C:113
libMesh::QBase::_dim
unsigned int _dim
The spatial dimension of the quadrature rule.
Definition: quadrature.h:339
libMesh::out
OStreamProxy out
libMesh::QBase::_points
std::vector< Point > _points
The locations of the quadrature points in reference element space.
Definition: quadrature.h:363
libMesh::FIRST
Definition: enum_order.h:42
libMesh::Quality::name
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:42
libMesh::QUAD8
Definition: enum_elem_type.h:42
libMesh::QGRUNDMANN_MOLLER
Definition: enum_quadrature_type.h:40
libMesh::QGRID
Definition: enum_quadrature_type.h:39
libMesh::QBase::init_2D
virtual void init_2D(const ElemType type=INVALID_ELEM, unsigned int p_level=0)
Initializes the 2D quadrature rule by filling the points and weights vectors with the appropriate val...
Definition: quadrature.C:123
libMesh::SEVENTEENTH
Definition: enum_order.h:58
Generated on Sat Jan 25 2020 12:08:41 for libMesh by   doxygen 1.8.16