Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
quadrature.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2018 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.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_quadrature_h
17 #define dealii_quadrature_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/point.h>
23 #include <deal.II/base/subscriptor.h>
24 
25 #include <array>
26 #include <memory>
27 #include <vector>
28 
29 DEAL_II_NAMESPACE_OPEN
30 
33 
84 template <int dim>
85 class Quadrature : public Subscriptor
86 {
87 public:
92  using SubQuadrature = Quadrature<dim - 1>;
93 
102  explicit Quadrature(const unsigned int n_quadrature_points = 0);
103 
112  Quadrature(const SubQuadrature &, const Quadrature<1> &);
113 
130  explicit Quadrature(const Quadrature<dim != 1 ? 1 : 0> &quadrature_1d);
131 
135  Quadrature(const Quadrature<dim> &q);
136 
141  Quadrature(Quadrature<dim> &&) noexcept = default;
142 
148  Quadrature(const std::vector<Point<dim>> &points,
149  const std::vector<double> & weights);
150 
158  Quadrature(const std::vector<Point<dim>> &points);
159 
164  Quadrature(const Point<dim> &point);
165 
169  virtual ~Quadrature() override = default;
170 
175  Quadrature &
176  operator=(const Quadrature<dim> &);
177 
182  Quadrature &
183  operator=(Quadrature<dim> &&) = default; // NOLINT
184 
188  bool
189  operator==(const Quadrature<dim> &p) const;
190 
195  void
196  initialize(const std::vector<Point<dim>> &points,
197  const std::vector<double> & weights);
198 
202  unsigned int
203  size() const;
204 
208  const Point<dim> &
209  point(const unsigned int i) const;
210 
214  const std::vector<Point<dim>> &
215  get_points() const;
216 
220  double
221  weight(const unsigned int i) const;
222 
226  const std::vector<double> &
227  get_weights() const;
228 
233  std::size_t
234  memory_consumption() const;
235 
240  template <class Archive>
241  void
242  serialize(Archive &ar, const unsigned int version);
243 
249  bool
250  is_tensor_product() const;
251 
270 #ifndef DOXYGEN
271  typename std::conditional<dim == 1,
272  std::array<Quadrature<1>, dim>,
273  const std::array<Quadrature<1>, dim> &>::type
274 #else
275  const std::array<Quadrature<1>, dim> &
276 #endif
277  get_tensor_basis() const;
278 
279 protected:
284  std::vector<Point<dim>> quadrature_points;
285 
290  std::vector<double> weights;
291 
300 
305  std::unique_ptr<std::array<Quadrature<1>, dim>> tensor_basis;
306 };
307 
308 
319 template <int dim>
320 class QAnisotropic : public Quadrature<dim>
321 {
322 public:
327  QAnisotropic(const Quadrature<1> &qx);
328 
332  QAnisotropic(const Quadrature<1> &qx, const Quadrature<1> &qy);
333 
337  QAnisotropic(const Quadrature<1> &qx,
338  const Quadrature<1> &qy,
339  const Quadrature<1> &qz);
340 };
341 
342 
368 template <int dim>
369 class QIterated : public Quadrature<dim>
370 {
371 public:
376  QIterated(const Quadrature<1> &base_quadrature, const unsigned int n_copies);
377 
382  "The quadrature formula you provided cannot be used "
383  "as the basis for iteration.");
384 };
385 
386 
387 
390 #ifndef DOXYGEN
391 
392 // ------------------- inline and template functions ----------------
393 
394 
395 template <int dim>
396 inline unsigned int
397 Quadrature<dim>::size() const
398 {
399  return weights.size();
400 }
401 
402 
403 template <int dim>
404 inline const Point<dim> &
405 Quadrature<dim>::point(const unsigned int i) const
406 {
407  AssertIndexRange(i, size());
408  return quadrature_points[i];
409 }
410 
411 
412 
413 template <int dim>
414 double
415 Quadrature<dim>::weight(const unsigned int i) const
416 {
417  AssertIndexRange(i, size());
418  return weights[i];
419 }
420 
421 
422 
423 template <int dim>
424 inline const std::vector<Point<dim>> &
426 {
427  return quadrature_points;
428 }
429 
430 
431 
432 template <int dim>
433 inline const std::vector<double> &
435 {
436  return weights;
437 }
438 
439 
440 
441 template <int dim>
442 inline bool
444 {
445  return is_tensor_product_flag;
446 }
447 
448 
449 
450 template <int dim>
451 template <class Archive>
452 inline void
453 Quadrature<dim>::serialize(Archive &ar, const unsigned int)
454 {
455  // forward to serialization
456  // function in the base class.
457  ar &static_cast<Subscriptor &>(*this);
458 
459  ar &quadrature_points &weights;
460 }
461 
462 
463 
464 /* -------------- declaration of explicit specializations ------------- */
465 
466 template <>
467 Quadrature<0>::Quadrature(const unsigned int);
468 template <>
470 template <>
472 
473 template <>
475 
476 template <>
478 
479 #endif // DOXYGEN
480 DEAL_II_NAMESPACE_CLOSE
481 
482 #endif
std::vector< double > weights
Definition: quadrature.h:290
const std::vector< Point< dim > > & get_points() const
const std::vector< double > & get_weights() const
Quadrature(const unsigned int n_quadrature_points=0)
Definition: quadrature.cc:40
#define AssertIndexRange(index, range)
Definition: exceptions.h:1637
const Point< dim > & point(const unsigned int i) const
STL namespace.
const std::array< Quadrature< 1 >, dim > & get_tensor_basis() const
Definition: quadrature.cc:316
Definition: point.h:110
QAnisotropic(const Quadrature< 1 > &qx)
Definition: quadrature.cc:343
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
std::size_t memory_consumption() const
Definition: quadrature.cc:304
static ::ExceptionBase & ExcInvalidQuadratureFormula()
bool is_tensor_product_flag
Definition: quadrature.h:299
std::vector< Point< dim > > quadrature_points
Definition: quadrature.h:284
unsigned int size() const
QIterated(const Quadrature< 1 > &base_quadrature, const unsigned int n_copies)
Definition: quadrature.cc:1777
std::unique_ptr< std::array< Quadrature< 1 >, dim > > tensor_basis
Definition: quadrature.h:305
void initialize(const std::vector< Point< dim >> &points, const std::vector< double > &weights)
Definition: quadrature.cc:50
void serialize(Archive &ar, const unsigned int version)
bool is_tensor_product() const
double weight(const unsigned int i) const