Reference documentation for deal.II version 8.4.1
quadrature.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2015 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii__quadrature_h
17 #define dealii__quadrature_h
18 
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/base/point.h>
22 #include <deal.II/base/subscriptor.h>
23 #include <vector>
24 
25 DEAL_II_NAMESPACE_OPEN
26 
29 
80 template <int dim>
81 class Quadrature : public Subscriptor
82 {
83 public:
88  typedef Quadrature<dim-1> SubQuadrature;
89 
98  explicit Quadrature (const unsigned int n_quadrature_points = 0);
99 
106  Quadrature (const SubQuadrature &,
107  const Quadrature<1> &);
108 
121  explicit Quadrature (const Quadrature<dim != 1 ? 1 : 0> &quadrature_1d);
122 
126  Quadrature (const Quadrature<dim> &q);
127 
133  Quadrature (const std::vector<Point<dim> > &points,
134  const std::vector<double> &weights);
135 
143  Quadrature (const std::vector<Point<dim> > &points);
144 
149  Quadrature (const Point<dim> &point);
150 
154  virtual ~Quadrature ();
155 
161 
165  bool operator == (const Quadrature<dim> &p) const;
166 
171  void initialize(const std::vector<Point<dim> > &points,
172  const std::vector<double> &weights);
173 
177  unsigned int size () const;
178 
182  const Point<dim> &point (const unsigned int i) const;
183 
187  const std::vector<Point<dim> > &get_points () const;
188 
192  double weight (const unsigned int i) const;
193 
197  const std::vector<double> &get_weights () const;
198 
203  std::size_t memory_consumption () const;
204 
209  template <class Archive>
210  void serialize (Archive &ar, const unsigned int version);
211 
212 protected:
217  std::vector<Point<dim> > quadrature_points;
218 
223  std::vector<double> weights;
224 };
225 
226 
237 template <int dim>
238 class QAnisotropic : public Quadrature<dim>
239 {
240 public:
245  QAnisotropic(const Quadrature<1> &qx);
246 
250  QAnisotropic(const Quadrature<1> &qx,
251  const Quadrature<1> &qy);
252 
256  QAnisotropic(const Quadrature<1> &qx,
257  const Quadrature<1> &qy,
258  const Quadrature<1> &qz);
259 };
260 
261 
287 template <int dim>
288 class QIterated : public Quadrature<dim>
289 {
290 public:
295  QIterated (const Quadrature<1> &base_quadrature,
296  const unsigned int n_copies);
297 
301  DeclExceptionMsg (ExcInvalidQuadratureFormula,
302  "The quadrature formula you provided cannot be used "
303  "as the basis for iteration.");
304 private:
309  static bool
310  uses_both_endpoints (const Quadrature<1> &base_quadrature);
311 };
312 
313 
314 
317 #ifndef DOXYGEN
318 
319 // ------------------- inline and template functions ----------------
320 
321 
322 template<int dim>
323 inline
324 unsigned int
325 Quadrature<dim>::size () const
326 {
327  return weights.size();
328 }
329 
330 
331 template <int dim>
332 inline
333 const Point<dim> &
334 Quadrature<dim>::point (const unsigned int i) const
335 {
336  AssertIndexRange(i, size());
337  return quadrature_points[i];
338 }
339 
340 
341 
342 template <int dim>
343 double
344 Quadrature<dim>::weight (const unsigned int i) const
345 {
346  AssertIndexRange(i, size());
347  return weights[i];
348 }
349 
350 
351 
352 template <int dim>
353 inline
354 const std::vector<Point<dim> > &
356 {
357  return quadrature_points;
358 }
359 
360 
361 
362 template <int dim>
363 inline
364 const std::vector<double> &
366 {
367  return weights;
368 }
369 
370 
371 
372 template <int dim>
373 template <class Archive>
374 inline
375 void
376 Quadrature<dim>::serialize (Archive &ar, const unsigned int)
377 {
378  // forward to serialization
379  // function in the base class.
380  ar &static_cast<Subscriptor &>(*this);
381 
382  ar &quadrature_points &weights;
383 }
384 
385 
386 
387 /* -------------- declaration of explicit specializations ------------- */
388 
389 template <>
390 Quadrature<0>::Quadrature (const unsigned int);
391 template <>
393  const Quadrature<1> &);
394 template <>
396 template <>
398 
399 template <>
401  const Quadrature<1> &);
402 
403 template <>
405 
406 #endif // DOXYGEN
407 DEAL_II_NAMESPACE_CLOSE
408 
409 #endif
std::vector< double > weights
Definition: quadrature.h:223
Quadrature(const unsigned int n_quadrature_points=0)
Definition: quadrature.cc:45
DeclExceptionMsg(ExcInvalidQuadratureFormula,"The quadrature formula you provided cannot be used ""as the basis for iteration.")
#define AssertIndexRange(index, range)
Definition: exceptions.h:1081
static bool uses_both_endpoints(const Quadrature< 1 > &base_quadrature)
const Point< dim > & point(const unsigned int i) const
Quadrature & operator=(const Quadrature< dim > &)
Definition: quadrature.cc:252
QAnisotropic(const Quadrature< 1 > &qx)
Definition: quadrature.cc:289
const std::vector< Point< dim > > & get_points() const
unsigned int size() const
void initialize(const std::vector< Point< dim > > &points, const std::vector< double > &weights)
Definition: quadrature.cc:55
virtual ~Quadrature()
Definition: quadrature.cc:273
double weight(const unsigned int i) const
Quadrature< dim-1 > SubQuadrature
Definition: quadrature.h:88
std::vector< Point< dim > > quadrature_points
Definition: quadrature.h:217
const std::vector< double > & get_weights() const
std::size_t memory_consumption() const
Definition: quadrature.cc:280
QIterated(const Quadrature< 1 > &base_quadrature, const unsigned int n_copies)
Definition: quadrature.cc:1758
bool operator==(const Quadrature< dim > &p) const
Definition: quadrature.cc:263
void serialize(Archive &ar, const unsigned int version)