quadrature.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_QUADRATURE_H
21 #define LIBMESH_QUADRATURE_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
26 #include "libmesh/point.h"
27 #include "libmesh/enum_elem_type.h"
28 #include "libmesh/enum_order.h"
30 #include "libmesh/auto_ptr.h"
31 
32 // C++ includes
33 #include <vector>
34 #include <string>
35 #include <utility>
36 
37 namespace libMesh
38 {
39 
40 
41 
51 // ------------------------------------------------------------
52 // QBase class definition
53 
54 class QBase : public ReferenceCountedObject<QBase>
55 {
56 protected:
57 
62  QBase (const unsigned int _dim,
63  const Order _order=INVALID_ORDER);
64 
65 public:
66 
70  virtual ~QBase() {}
71 
75  virtual QuadratureType type() const = 0;
76 
86  static AutoPtr<QBase> build (const std::string &name,
87  const unsigned int _dim,
88  const Order _order=INVALID_ORDER);
89 
97  static AutoPtr<QBase> build (const QuadratureType _qt,
98  const unsigned int _dim,
99  const Order _order=INVALID_ORDER);
100 
105  { return _type; }
106 
110  unsigned int get_p_level() const
111  { return _p_level; }
112 
116  unsigned int n_points() const
117  { libmesh_assert (!_points.empty());
118  return libmesh_cast_int<unsigned int>(_points.size()); }
119 
123  unsigned int get_dim() const { return _dim; }
124 
129  const std::vector<Point>& get_points() const { return _points; }
130 
135  std::vector<Point>& get_points() { return _points; }
136 
140  const std::vector<Real>& get_weights() const { return _weights; }
141 
145  std::vector<Real>& get_weights() { return _weights; }
146 
150  Point qp(const unsigned int i) const
151  { libmesh_assert_less (i, _points.size()); return _points[i]; }
152 
156  Real w(const unsigned int i) const
157  { libmesh_assert_less (i, _weights.size()); return _weights[i]; }
158 
163  void init (const ElemType type=INVALID_ELEM,
164  unsigned int p_level=0);
165 
169  Order get_order() const { return static_cast<Order>(_order + _p_level); }
170 
175  void print_info(std::ostream& os=libMesh::out) const;
176 
183  void scale(std::pair<Real, Real> old_range,
184  std::pair<Real, Real> new_range);
185 
189  friend std::ostream& operator << (std::ostream& os, const QBase& q);
190 
198  virtual bool shapes_need_reinit() { return false; }
199 
216 
217 protected:
218 
219 
225  virtual void init_0D (const ElemType type=INVALID_ELEM,
226  unsigned int p_level=0);
227 
235  virtual void init_1D (const ElemType type=INVALID_ELEM,
236  unsigned int p_level=0) = 0;
237 
246  virtual void init_2D (const ElemType,
247  unsigned int =0)
248 #ifndef DEBUG
249  {}
250 #else
251  {
252  libMesh::err << "ERROR: Seems as if this quadrature rule" << std::endl
253  << " is not implemented for 2D." << std::endl;
254  libmesh_error();
255  }
256 #endif
257 
266  virtual void init_3D (const ElemType,
267  unsigned int =0)
268 #ifndef DEBUG
269  {}
270 #else
271  {
272  libMesh::err << "ERROR: Seems as if this quadrature rule" << std::endl
273  << " is not implemented for 3D." << std::endl;
274  libmesh_error();
275  }
276 #endif
277 
278 
285  void tensor_product_quad (const QBase& q1D);
286 
293  void tensor_product_hex (const QBase& q1D);
294 
302  void tensor_product_prism (const QBase& q1D, const QBase& q2D);
303 
304 
305 
309  const unsigned int _dim;
310 
314  const Order _order;
315 
320  ElemType _type;
321 
326  unsigned int _p_level;
327 
332  std::vector<Point> _points;
333 
337  std::vector<Real> _weights;
338 };
339 
340 
341 
342 
343 
344 // ------------------------------------------------------------
345 // QBase class members
346 
347 inline
348 QBase::QBase(const unsigned int d,
349  const Order o) :
350  allow_rules_with_negative_weights(true),
351  _dim(d),
352  _order(o),
353  _type(INVALID_ELEM),
354  _p_level(0)
355 {
356 }
357 
358 
359 
360 
361 inline
362 void QBase::print_info(std::ostream& os) const
363 {
364  libmesh_assert(!_points.empty());
365  libmesh_assert(!_weights.empty());
366 
367  os << "N_Q_Points=" << this->n_points() << std::endl << std::endl;
368  for (unsigned int qpoint=0; qpoint<this->n_points(); qpoint++)
369  {
370  os << " Point " << qpoint << ":\n"
371  << " "
372  << _points[qpoint]
373  << " Weight:\n "
374  << " w=" << _weights[qpoint] << "\n" << std::endl;
375  }
376 }
377 
378 } // namespace libMesh
379 
380 #endif // LIBMESH_QUADRATURE_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:06 UTC

Hosted By:
SourceForge.net Logo