libMesh
quadrature.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 // C++ includes
20 
21 // Local includes
22 #include "libmesh/elem.h"
23 #include "libmesh/quadrature.h"
24 
25 namespace libMesh
26 {
27 
28 void QBase::init(const ElemType t,
29  unsigned int p)
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 }
69 
70 
71 
72 void QBase::init (const Elem & elem,
73  const std::vector<Real> & /* vertex_distance_func */,
74  unsigned int p_level)
75 {
76  // dispatch generic implementation
77  this->init(elem.type(), p_level);
78 }
79 
80 
81 
83  unsigned int)
84 {
85  _points.resize(1);
86  _weights.resize(1);
87  _points[0] = Point(0.);
88  _weights[0] = 1.0;
89 }
90 
91 
92 
93 void QBase::scale(std::pair<Real, Real> old_range,
94  std::pair<Real, Real> new_range)
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 }
123 
124 
125 
126 
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 }
149 
150 
151 
152 
153 
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 }
177 
178 
179 
180 
181 void QBase::tensor_product_prism(const QBase & q1D, const QBase & q2D)
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 }
204 
205 
206 
207 
208 std::ostream & operator << (std::ostream & os, const QBase & q)
209 {
210  q.print_info(os);
211  return os;
212 }
213 
214 } // namespace libMesh
virtual void init(const ElemType type=INVALID_ELEM, unsigned int p_level=0)
Initializes the data structures for a quadrature rule for an element of type type.
Definition: quadrature.C:28
virtual ElemType type() const =0
const unsigned int _dim
The spatial dimension of the quadrature rule.
Definition: quadrature.h:311
ElemType
Defines an enum for geometric element types.
Real w(const unsigned int i) const
Definition: quadrature.h:160
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
ElemType _type
The type of element for which the current values have been computed.
Definition: quadrature.h:323
friend std::ostream & operator<<(std::ostream &os, const QBase &q)
Same as above, but allows you to use the stream syntax.
Definition: quadrature.C:208
The libMesh namespace provides an interface to certain functionality in the library.
Point qp(const unsigned int i) const
Definition: quadrature.h:151
std::vector< Point > _points
The locations of the quadrature points in reference element space.
Definition: quadrature.h:335
std::vector< Real > _weights
The quadrature weights.
Definition: quadrature.h:341
void tensor_product_prism(const QBase &q1D, const QBase &q2D)
Computes the tensor product of a 1D quadrature rule and a 2D quadrature rule.
Definition: quadrature.C:181
unsigned int _p_level
The p-level of the element for which the current values have been computed.
Definition: quadrature.h:329
virtual void init_2D(const ElemType, unsigned int=0)
Initializes the 2D quadrature rule by filling the points and weights vectors with the appropriate val...
Definition: quadrature.h:263
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...
void tensor_product_hex(const QBase &q1D)
Computes the tensor product quadrature rule [q1D x q1D x q1D] from the 1D rule q1D.
Definition: quadrature.C:154
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:82
void print_info(std::ostream &os=libMesh::out) const
Prints information relevant to the quadrature rule, by default to libMesh::out.
Definition: quadrature.h:364
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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...
Definition: quadrature.C:93
void tensor_product_quad(const QBase &q1D)
Constructs a 2D rule from the tensor product of q1D with itself.
Definition: quadrature.C:127
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
virtual void init_3D(const ElemType, unsigned int=0)
Initializes the 3D quadrature rule by filling the points and weights vectors with the appropriate val...
Definition: quadrature.h:279
unsigned int n_points() const
Definition: quadrature.h:113
The QBase class provides the basic functionality from which various quadrature rules can be derived...
Definition: quadrature.h:53