libMesh::QConical Class Referencefinal

Conical product quadrature rules for Tri and Tet elements. More...

#include <quadrature_conical.h>

Inheritance diagram for libMesh::QConical:

Public Member Functions

 QConical (unsigned int d, Order o=INVALID_ORDER)
 
 QConical (const QConical &)=default
 
 QConical (QConical &&)=default
 
QConicaloperator= (const QConical &)=default
 
QConicaloperator= (QConical &&)=default
 
virtual ~QConical ()=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)
 
virtual void init (const Elem &elem, const std::vector< Real > &vertex_distance_func, unsigned int p_level=0)
 
Order get_order () const
 
void print_info (std::ostream &os=libMesh::out) const
 
void scale (std::pair< Real, Real > old_range, std::pair< Real, Real > new_range)
 
virtual bool shapes_need_reinit ()
 

Static Public Member Functions

static std::unique_ptr< QBasebuild (const std::string &name, const unsigned int dim, const Order order=INVALID_ORDER)
 
static std::unique_ptr< QBasebuild (const QuadratureType qt, const unsigned int dim, const Order order=INVALID_ORDER)
 
static void print_info (std::ostream &out=libMesh::out)
 
static std::string get_info ()
 
static unsigned int n_objects ()
 
static void enable_print_counter_info ()
 
static void disable_print_counter_info ()
 

Public Attributes

bool allow_rules_with_negative_weights
 

Protected Types

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

Protected Member Functions

virtual void init_0D (const ElemType type=INVALID_ELEM, unsigned int p_level=0)
 
void tensor_product_quad (const QBase &q1D)
 
void tensor_product_hex (const QBase &q1D)
 
void tensor_product_prism (const QBase &q1D, const QBase &q2D)
 
void increment_constructor_count (const std::string &name)
 
void increment_destructor_count (const std::string &name)
 

Protected Attributes

unsigned int _dim
 
Order _order
 
ElemType _type
 
unsigned int _p_level
 
std::vector< Point_points
 
std::vector< Real_weights
 

Static Protected Attributes

static Counts _counts
 
static Threads::atomic< unsigned int > _n_objects
 
static Threads::spin_mutex _mutex
 
static bool _enable_print_counter = true
 

Private Member Functions

virtual void init_1D (const ElemType, unsigned int=0) override
 
virtual void init_2D (const ElemType _type=INVALID_ELEM, unsigned int p_level=0) override
 
virtual void init_3D (const ElemType _type=INVALID_ELEM, unsigned int p_level=0) override
 
void conical_product_tri (unsigned int p)
 
void conical_product_tet (unsigned int p)
 
void conical_product_pyramid (unsigned int p)
 

Detailed Description

Conical product quadrature rules for Tri and Tet elements.

This class implements the so-called conical product quadrature rules for Tri and Tet elements. These rules are generally non-optimal in the number of evaluation points, but have the following nice properties: .) All positive weights. .) Easily constructed for any order given the underlying 1D Gauss and Jacobi quadrature rules. The construction of these rules is given by e.g. A. H. Stroud, "Approximate Calculation of Multiple Integrals.", 1972.

Author
John W. Peterson
Date
2008

Definition at line 43 of file quadrature_conical.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

◆ QConical() [1/3]

libMesh::QConical::QConical ( unsigned int  d,
Order  o = INVALID_ORDER 
)
inline

Constructor. Declares the order of the quadrature rule.

Definition at line 50 of file quadrature_conical.h.

51  :
52  QBase(d,o)
53  {}
QBase(unsigned int dim, Order order=INVALID_ORDER)
Definition: quadrature.h:364

◆ QConical() [2/3]

libMesh::QConical::QConical ( const QConical )
default

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

◆ QConical() [3/3]

libMesh::QConical::QConical ( QConical &&  )
default

◆ ~QConical()

virtual libMesh::QConical::~QConical ( )
virtualdefault

Member Function Documentation

◆ build() [1/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 40 of file quadrature_build.C.

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

Referenced by libMesh::InfFE< Dim, T_radial, T_map >::attach_quadrature_rule().

43 {
44  return QBase::build (Utility::string_to_enum<QuadratureType> (type),
45  _dim,
46  _order);
47 }
unsigned int _dim
Definition: quadrature.h:325
virtual QuadratureType type() const =0
static std::unique_ptr< QBase > build(const std::string &name, const unsigned int dim, const Order order=INVALID_ORDER)

◆ build() [2/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 51 of file quadrature_build.C.

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::QSIMPSON, libMesh::QTRAP, libMesh::THIRD, and libMesh::TWENTYTHIRD.

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

◆ conical_product_pyramid()

void libMesh::QConical::conical_product_pyramid ( unsigned int  p)
private

Implementation of conical product rule for a Pyramid in 3D of order = _order+2*p.

Definition at line 181 of file quadrature_conical.C.

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

Referenced by init_3D().

182 {
183  // Be sure the underlying rule object was built with the same dimension as the
184  // rule we are about to construct.
185  libmesh_assert_equal_to (this->get_dim(), 3);
186 
187  QGauss gauss1D(1,static_cast<Order>(_order+2*p));
188  QJacobi jac1D(1,static_cast<Order>(_order+2*p),2,0);
189 
190  // These rules should have the same number of points
191  libmesh_assert_equal_to (gauss1D.n_points(), jac1D.n_points());
192 
193  // Save the number of points as a convenient variable
194  const unsigned int np = gauss1D.n_points();
195 
196  // Resize the points and weights vectors
197  _points.resize(np * np * np);
198  _weights.resize(np * np * np);
199 
200  // Compute the conical product
201  unsigned int q = 0;
202  for (unsigned int i=0; i<np; ++i)
203  for (unsigned int j=0; j<np; ++j)
204  for (unsigned int k=0; k<np; ++k, ++q)
205  {
206  const Real xi=gauss1D.qp(i)(0);
207  const Real yj=gauss1D.qp(j)(0);
208  const Real zk=jac1D.qp(k)(0);
209 
210  _points[q](0) = (1.-zk) * xi;
211  _points[q](1) = (1.-zk) * yj;
212  _points[q](2) = zk;
213  _weights[q] = gauss1D.w(i) * gauss1D.w(j) * jac1D.w(k);
214  }
215 }
std::vector< Point > _points
Definition: quadrature.h:349
std::vector< Real > _weights
Definition: quadrature.h:355
unsigned int get_dim() const
Definition: quadrature.h:136
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ conical_product_tet()

void libMesh::QConical::conical_product_tet ( unsigned int  p)
private

Implementation of conical product rule for a Tet in 3D of order = _order+2*p.

Definition at line 104 of file quadrature_conical.C.

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

Referenced by init_3D().

105 {
106  // Be sure the underlying rule object was built with the same dimension as the
107  // rule we are about to construct.
108  libmesh_assert_equal_to (this->get_dim(), 3);
109 
110  QGauss gauss1D(1,static_cast<Order>(_order+2*p));
111  QJacobi jacA1D(1,static_cast<Order>(_order+2*p),1,0);
112  QJacobi jacB1D(1,static_cast<Order>(_order+2*p),2,0);
113 
114  // The Gauss rule needs to be scaled to [0,1]
115  std::pair<Real, Real> old_range(-1.0L, 1.0L);
116  std::pair<Real, Real> new_range( 0.0L, 1.0L);
117  gauss1D.scale(old_range,
118  new_range);
119 
120  // Now construct the points and weights for the conical product rule.
121 
122  // All rules should have the same number of points
123  libmesh_assert_equal_to (gauss1D.n_points(), jacA1D.n_points());
124  libmesh_assert_equal_to (jacA1D.n_points(), jacB1D.n_points());
125 
126  // Save the number of points as a convenient variable
127  const unsigned int np = gauss1D.n_points();
128 
129  // All rules should be between x=0 and x=1
130  libmesh_assert_greater_equal (gauss1D.qp(0)(0), 0.0);
131  libmesh_assert_less_equal (gauss1D.qp(np-1)(0), 1.0);
132  libmesh_assert_greater_equal (jacA1D.qp(0)(0), 0.0);
133  libmesh_assert_less_equal (jacA1D.qp(np-1)(0), 1.0);
134  libmesh_assert_greater_equal (jacB1D.qp(0)(0), 0.0);
135  libmesh_assert_less_equal (jacB1D.qp(np-1)(0), 1.0);
136 
137  // Resize the points and weights vectors
138  _points.resize(np * np * np);
139  _weights.resize(np * np * np);
140 
141  // Compute the conical product
142  unsigned int gp = 0;
143  for (unsigned int i=0; i<np; i++)
144  for (unsigned int j=0; j<np; j++)
145  for (unsigned int k=0; k<np; k++)
146  {
147  _points[gp](0) = jacB1D.qp(k)(0); //t[k];
148  _points[gp](1) = jacA1D.qp(j)(0) * (1.-jacB1D.qp(k)(0)); //s[j]*(1.-t[k]);
149  _points[gp](2) = gauss1D.qp(i)(0) * (1.-jacA1D.qp(j)(0)) * (1.-jacB1D.qp(k)(0)); //r[i]*(1.-s[j])*(1.-t[k]);
150  _weights[gp] = gauss1D.w(i) * jacA1D.w(j) * jacB1D.w(k); //A[i]*B[j]*C[k];
151  gp++;
152  }
153 }
std::vector< Point > _points
Definition: quadrature.h:349
std::vector< Real > _weights
Definition: quadrature.h:355
unsigned int get_dim() const
Definition: quadrature.h:136

◆ conical_product_tri()

void libMesh::QConical::conical_product_tri ( unsigned int  p)
private

Implementation of conical product rule for a Tri in 2D of order = _order+2*p.

Definition at line 53 of file quadrature_conical.C.

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

Referenced by init_2D().

54 {
55  // Be sure the underlying rule object was built with the same dimension as the
56  // rule we are about to construct.
57  libmesh_assert_equal_to (this->get_dim(), 2);
58 
59  QGauss gauss1D(1,static_cast<Order>(_order+2*p));
60  QJacobi jac1D(1,static_cast<Order>(_order+2*p),1,0);
61 
62  // The Gauss rule needs to be scaled to [0,1]
63  std::pair<Real, Real> old_range(-1.0L, 1.0L);
64  std::pair<Real, Real> new_range( 0.0L, 1.0L);
65  gauss1D.scale(old_range,
66  new_range);
67 
68  // Now construct the points and weights for the conical product rule.
69 
70  // Both rules should have the same number of points.
71  libmesh_assert_equal_to (gauss1D.n_points(), jac1D.n_points());
72 
73  // Save the number of points as a convenient variable
74  const unsigned int np = gauss1D.n_points();
75 
76  // Both rules should be between x=0 and x=1
77  libmesh_assert_greater_equal (gauss1D.qp(0)(0), 0.0);
78  libmesh_assert_less_equal (gauss1D.qp(np-1)(0), 1.0);
79  libmesh_assert_greater_equal (jac1D.qp(0)(0), 0.0);
80  libmesh_assert_less_equal (jac1D.qp(np-1)(0), 1.0);
81 
82  // Resize the points and weights vectors
83  _points.resize(np * np);
84  _weights.resize(np * np);
85 
86  // Compute the conical product
87  unsigned int gp = 0;
88  for (unsigned int i=0; i<np; i++)
89  for (unsigned int j=0; j<np; j++)
90  {
91  _points[gp](0) = jac1D.qp(j)(0); //s[j];
92  _points[gp](1) = gauss1D.qp(i)(0) * (1.-jac1D.qp(j)(0)); //r[i]*(1.-s[j]);
93  _weights[gp] = gauss1D.w(i) * jac1D.w(j); //A[i]*B[j];
94  gp++;
95  }
96 }
std::vector< Point > _points
Definition: quadrature.h:349
std::vector< Real > _weights
Definition: quadrature.h:355
unsigned int get_dim() const
Definition: quadrature.h:136

◆ disable_print_counter_info()

void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited

Definition at line 106 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

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

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

◆ 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.

References libMesh::ReferenceCounter::_enable_print_counter.

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

◆ get_dim()

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

Definition at line 136 of file quadrature.h.

References libMesh::QBase::_dim.

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

136 { return _dim; }
unsigned int _dim
Definition: quadrature.h:325

◆ get_elem_type()

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

Definition at line 117 of file quadrature.h.

References libMesh::QBase::_type.

117 { return _type; }
ElemType _type
Definition: quadrature.h:337

◆ 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.

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

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

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 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:42

◆ get_order()

Order libMesh::QBase::get_order ( ) const
inlineinherited
Returns
The order of the quadrature rule.

Definition at line 203 of file quadrature.h.

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

Referenced by libMesh::InfFE< Dim, T_radial, T_map >::attach_quadrature_rule().

203 { return static_cast<Order>(_order + _p_level); }
unsigned int _p_level
Definition: quadrature.h:343

◆ get_p_level()

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

Definition at line 122 of file quadrature.h.

References libMesh::QBase::_p_level.

122 { return _p_level; }
unsigned int _p_level
Definition: quadrature.h:343

◆ get_points() [1/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 142 of file quadrature.h.

References libMesh::QBase::_points.

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

142 { return _points; }
std::vector< Point > _points
Definition: quadrature.h:349

◆ get_points() [2/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 148 of file quadrature.h.

References libMesh::QBase::_points.

148 { return _points; }
std::vector< Point > _points
Definition: quadrature.h:349

◆ get_weights() [1/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 154 of file quadrature.h.

References libMesh::QBase::_weights.

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

154 { return _weights; }
std::vector< Real > _weights
Definition: quadrature.h:355

◆ get_weights() [2/2]

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

Definition at line 160 of file quadrature.h.

References libMesh::QBase::_weights.

160 { return _weights; }
std::vector< Real > _weights
Definition: quadrature.h:355

◆ 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.

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

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

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 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:42
spin_mutex spin_mtx
Definition: threads.C:29

◆ 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.

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

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

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 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:42
spin_mutex spin_mtx
Definition: threads.C:29

◆ init() [1/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 28 of file quadrature.C.

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::QMonomial::init_1D(), libMesh::QClough::init_2D(), libMesh::QGaussLobatto::init_2D(), libMesh::QGrid::init_2D(), libMesh::QTrap::init_2D(), libMesh::QSimpson::init_2D(), libMesh::QGauss::init_2D(), libMesh::QMonomial::init_2D(), libMesh::QGaussLobatto::init_3D(), libMesh::QTrap::init_3D(), libMesh::QGrid::init_3D(), libMesh::QSimpson::init_3D(), libMesh::QGauss::init_3D(), libMesh::QMonomial::init_3D(), libMesh::QGauss::QGauss(), libMesh::QGaussLobatto::QGaussLobatto(), libMesh::QJacobi::QJacobi(), libMesh::QSimpson::QSimpson(), and libMesh::QTrap::QTrap().

30 {
31  // check to see if we have already
32  // done the work for this quadrature rule
33  if (t == _type && p == _p_level)
34  return;
35  else
36  {
37  _type = t;
38  _p_level = p;
39  }
40 
41 
42 
43  switch(_dim)
44  {
45  case 0:
46  this->init_0D(_type,_p_level);
47 
48  return;
49 
50  case 1:
51  this->init_1D(_type,_p_level);
52 
53  return;
54 
55  case 2:
56  this->init_2D(_type,_p_level);
57 
58  return;
59 
60  case 3:
61  this->init_3D(_type,_p_level);
62 
63  return;
64 
65  default:
66  libmesh_error_msg("Invalid dimension _dim = " << _dim);
67  }
68 }
ElemType _type
Definition: quadrature.h:337
unsigned int _dim
Definition: quadrature.h:325
unsigned int _p_level
Definition: quadrature.h:343
virtual void init_2D(const ElemType, unsigned int=0)
Definition: quadrature.h:277
virtual void init_1D(const ElemType type=INVALID_ELEM, unsigned int p_level=0)=0
virtual void init_0D(const ElemType type=INVALID_ELEM, unsigned int p_level=0)
Definition: quadrature.C:82
virtual void init_3D(const ElemType, unsigned int=0)
Definition: quadrature.h:293

◆ init() [2/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 72 of file quadrature.C.

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

75 {
76  // dispatch generic implementation
77  this->init(elem.type(), p_level);
78 }
virtual void init(const ElemType type=INVALID_ELEM, unsigned int p_level=0)
Definition: quadrature.C:28

◆ 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.

Definition at line 82 of file quadrature.C.

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

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

84 {
85  _points.resize(1);
86  _weights.resize(1);
87  _points[0] = Point(0.);
88  _weights[0] = 1.0;
89 }
std::vector< Point > _points
Definition: quadrature.h:349
std::vector< Real > _weights
Definition: quadrature.h:355

◆ init_1D()

void libMesh::QConical::init_1D ( const ElemType  ,
unsigned int  p = 0 
)
overrideprivatevirtual

The optimal "conical product" rule in 1D is simply Gauss.

Implements libMesh::QBase.

Definition at line 38 of file quadrature_conical.C.

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

40 {
41  QGauss gauss1D(1, static_cast<Order>(_order+2*p));
42 
43  // Swap points and weights with the about-to-be destroyed rule.
44  _points.swap(gauss1D.get_points());
45  _weights.swap(gauss1D.get_weights());
46 }
std::vector< Point > _points
Definition: quadrature.h:349
std::vector< Real > _weights
Definition: quadrature.h:355

◆ init_2D()

void libMesh::QConical::init_2D ( const ElemType  _type = INVALID_ELEM,
unsigned int  p_level = 0 
)
overrideprivatevirtual

The conical product rules are defined in 2D only for Tris.

Reimplemented from libMesh::QBase.

Definition at line 27 of file quadrature_conical_2D.C.

References conical_product_tri(), libMesh::TRI3, and libMesh::TRI6.

29 {
30  switch (type_in)
31  {
32  case TRI3:
33  case TRI6:
34  {
35  this->conical_product_tri(p);
36  return;
37  } // end case TRI3, TRI6
38 
39  //---------------------------------------------
40  // Unsupported element type
41  default:
42  libmesh_error_msg("ERROR: Unsupported element type: " << type_in);
43  } // end switch (type_in)
44 }
void conical_product_tri(unsigned int p)

◆ init_3D()

void libMesh::QConical::init_3D ( const ElemType  _type = INVALID_ELEM,
unsigned int  p_level = 0 
)
overrideprivatevirtual

The conical product rules are defined in 3D only for Tets.

Reimplemented from libMesh::QBase.

Definition at line 27 of file quadrature_conical_3D.C.

References conical_product_pyramid(), conical_product_tet(), libMesh::PYRAMID13, libMesh::PYRAMID14, libMesh::PYRAMID5, libMesh::TET10, and libMesh::TET4.

29 {
30  switch (type_in)
31  {
32  case TET4:
33  case TET10:
34  {
35  this->conical_product_tet(p);
36  return;
37  } // end case TET4, TET10
38 
39  case PYRAMID5:
40  case PYRAMID13:
41  case PYRAMID14:
42  {
43  this->conical_product_pyramid(p);
44  return;
45  } // end case PYRAMID5
46 
47 
48  //---------------------------------------------
49  // Unsupported element type
50  default:
51  libmesh_error_msg("ERROR: Unsupported element type: " << type_in);
52  } // end switch (type_in)
53 }
void conical_product_tet(unsigned int p)
void conical_product_pyramid(unsigned int p)

◆ 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.

References libMesh::ReferenceCounter::_n_objects.

84  { return _n_objects; }
static Threads::atomic< unsigned int > _n_objects

◆ n_points()

unsigned int libMesh::QBase::n_points ( ) const
inlineinherited

◆ operator=() [1/2]

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

◆ operator=() [2/2]

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

◆ print_info() [1/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.

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

88 {
90  out_stream << ReferenceCounter::get_info();
91 }
static std::string get_info()

◆ print_info() [2/2]

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

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

Definition at line 378 of file quadrature.h.

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

Referenced by libMesh::operator<<().

379 {
380  libmesh_assert(!_points.empty());
381  libmesh_assert(!_weights.empty());
382 
383  Real summed_weights=0;
384  os << "N_Q_Points=" << this->n_points() << std::endl << std::endl;
385  for (unsigned int qpoint=0; qpoint<this->n_points(); qpoint++)
386  {
387  os << " Point " << qpoint << ":\n"
388  << " "
389  << _points[qpoint]
390  << "\n Weight:\n "
391  << " w=" << _weights[qpoint] << "\n" << std::endl;
392 
393  summed_weights += _weights[qpoint];
394  }
395  os << "Summed Weights: " << summed_weights << std::endl;
396 }
std::vector< Point > _points
Definition: quadrature.h:349
std::vector< Real > _weights
Definition: quadrature.h:355
unsigned int n_points() const
Definition: quadrature.h:127
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ qp()

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

Definition at line 165 of file quadrature.h.

References libMesh::QBase::_points.

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

166  {
167  libmesh_assert_less (i, _points.size());
168  return _points[i];
169  }
std::vector< Point > _points
Definition: quadrature.h:349

◆ 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 93 of file quadrature.C.

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

Referenced by conical_product_tet(), and conical_product_tri().

95 {
96  // Make sure we are in 1D
97  libmesh_assert_equal_to (_dim, 1);
98 
99  Real
100  h_new = new_range.second - new_range.first,
101  h_old = old_range.second - old_range.first;
102 
103  // Make sure that we have sane ranges
104  libmesh_assert_greater (h_new, 0.);
105  libmesh_assert_greater (h_old, 0.);
106 
107  // Make sure there are some points
108  libmesh_assert_greater (_points.size(), 0);
109 
110  // Compute the scale factor
111  Real scfact = h_new/h_old;
112 
113  // We're mapping from old_range -> new_range
114  for (std::size_t i=0; i<_points.size(); i++)
115  {
116  _points[i](0) = new_range.first +
117  (_points[i](0) - old_range.first) * scfact;
118 
119  // Scale the weights
120  _weights[i] *= scfact;
121  }
122 }
unsigned int _dim
Definition: quadrature.h:325
std::vector< Point > _points
Definition: quadrature.h:349
std::vector< Real > _weights
Definition: quadrature.h:355
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ 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 231 of file quadrature.h.

231 { return false; }

◆ 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 154 of file quadrature.C.

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().

155 {
156  const unsigned int np = q1D.n_points();
157 
158  _points.resize(np * np * np);
159 
160  _weights.resize(np * np * np);
161 
162  unsigned int q=0;
163 
164  for (unsigned int k=0; k<np; k++)
165  for (unsigned int j=0; j<np; j++)
166  for (unsigned int i=0; i<np; i++)
167  {
168  _points[q](0) = q1D.qp(i)(0);
169  _points[q](1) = q1D.qp(j)(0);
170  _points[q](2) = q1D.qp(k)(0);
171 
172  _weights[q] = q1D.w(i) * q1D.w(j) * q1D.w(k);
173 
174  q++;
175  }
176 }
std::vector< Point > _points
Definition: quadrature.h:349
std::vector< Real > _weights
Definition: quadrature.h:355

◆ 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 181 of file quadrature.C.

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().

182 {
183  const unsigned int n_points1D = q1D.n_points();
184  const unsigned int n_points2D = q2D.n_points();
185 
186  _points.resize (n_points1D * n_points2D);
187  _weights.resize (n_points1D * n_points2D);
188 
189  unsigned int q=0;
190 
191  for (unsigned int j=0; j<n_points1D; j++)
192  for (unsigned int i=0; i<n_points2D; i++)
193  {
194  _points[q](0) = q2D.qp(i)(0);
195  _points[q](1) = q2D.qp(i)(1);
196  _points[q](2) = q1D.qp(j)(0);
197 
198  _weights[q] = q2D.w(i) * q1D.w(j);
199 
200  q++;
201  }
202 
203 }
std::vector< Point > _points
Definition: quadrature.h:349
std::vector< Real > _weights
Definition: quadrature.h:355

◆ 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 127 of file quadrature.C.

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().

128 {
129 
130  const unsigned int np = q1D.n_points();
131 
132  _points.resize(np * np);
133 
134  _weights.resize(np * np);
135 
136  unsigned int q=0;
137 
138  for (unsigned int j=0; j<np; j++)
139  for (unsigned int i=0; i<np; i++)
140  {
141  _points[q](0) = q1D.qp(i)(0);
142  _points[q](1) = q1D.qp(j)(0);
143 
144  _weights[q] = q1D.w(i)*q1D.w(j);
145 
146  q++;
147  }
148 }
std::vector< Point > _points
Definition: quadrature.h:349
std::vector< Real > _weights
Definition: quadrature.h:355

◆ type()

QuadratureType libMesh::QConical::type ( ) const
overridevirtual
Returns
The QuadratureType for this class.

Implements libMesh::QBase.

Definition at line 33 of file quadrature_conical.C.

References libMesh::QCONICAL.

34 {
35  return QCONICAL;
36 }

◆ w()

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

Definition at line 174 of file quadrature.h.

References libMesh::QBase::_weights.

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

175  {
176  libmesh_assert_less (i, _weights.size());
177  return _weights[i];
178  }
std::vector< Real > _weights
Definition: quadrature.h:355

Member Data Documentation

◆ _counts

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

◆ _dim

unsigned int libMesh::QBase::_dim
protectedinherited

◆ _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

◆ _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 343 of file quadrature.h.

Referenced by libMesh::QBase::get_order(), libMesh::QBase::get_p_level(), and libMesh::QBase::init().

◆ _points

◆ _type

ElemType libMesh::QBase::_type
protectedinherited

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

Definition at line 337 of file quadrature.h.

Referenced by libMesh::QBase::get_elem_type(), and libMesh::QBase::init().

◆ _weights

◆ 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 246 of file quadrature.h.

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


The documentation for this class was generated from the following files: