Reference documentation for deal.II version Git 6d63218887 2020-10-30 17:17:53 -0400
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
quadrature_lib.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2020 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_lib_h
17 #define dealii_quadrature_lib_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 
25 
28 
37 template <int dim>
38 class QGauss : public Quadrature<dim>
39 {
40 public:
45  QGauss(const unsigned int n);
46 };
47 
48 
71 template <int dim>
72 class QGaussLobatto : public Quadrature<dim>
73 {
74 public:
79  QGaussLobatto(const unsigned int n);
80 };
81 
82 
83 
88 template <int dim>
89 class QMidpoint : public Quadrature<dim>
90 {
91 public:
92  QMidpoint();
93 };
94 
95 
100 template <int dim>
101 class QSimpson : public Quadrature<dim>
102 {
103 public:
104  QSimpson();
105 };
106 
107 
108 
117 template <int dim>
118 class QTrapezoid : public Quadrature<dim>
119 {
120 public:
121  QTrapezoid();
122 };
123 
124 
136 template <int dim>
138 
139 
140 
147 template <int dim>
148 class QMilne : public Quadrature<dim>
149 {
150 public:
151  QMilne();
152 };
153 
154 
161 template <int dim>
162 class QWeddle : public Quadrature<dim>
163 {
164 public:
165  QWeddle();
166 };
167 
168 
169 
183 template <int dim>
184 class QGaussLog : public Quadrature<dim>
185 {
186 public:
190  QGaussLog(const unsigned int n, const bool revert = false);
191 
192 private:
196  static std::vector<double>
197  get_quadrature_points(const unsigned int n);
198 
202  static std::vector<double>
203  get_quadrature_weights(const unsigned int n);
204 };
205 
206 
207 
244 template <int dim>
245 class QGaussLogR : public Quadrature<dim>
246 {
247 public:
255  QGaussLogR(const unsigned int n,
256  const Point<dim> x0 = Point<dim>(),
257  const double alpha = 1,
258  const bool factor_out_singular_weight = false);
259 
265  QGaussLogR(QGaussLogR<dim> &&) noexcept = default;
266 
267 protected:
272  const double fraction;
273 };
274 
275 
299 template <int dim>
300 class QGaussOneOverR : public Quadrature<dim>
301 {
302 public:
335  QGaussOneOverR(const unsigned int n,
336  const Point<dim> singularity,
337  const bool factor_out_singular_weight = false);
372  QGaussOneOverR(const unsigned int n,
373  const unsigned int vertex_index,
374  const bool factor_out_singular_weight = false);
375 
376 private:
382  static unsigned int
383  quad_size(const Point<dim> singularity, const unsigned int n);
384 };
385 
386 
387 
397 template <int dim>
398 class QSorted : public Quadrature<dim>
399 {
400 public:
405  QSorted(const Quadrature<dim> &quad);
406 
407 private:
413  bool
414  compare_weights(const unsigned int a, const unsigned int b) const;
415 };
416 
471 template <int dim>
472 class QTelles : public Quadrature<dim>
473 {
474 public:
481  QTelles(const Quadrature<1> &base_quad, const Point<dim> &singularity);
487  QTelles(const unsigned int n, const Point<dim> &singularity);
488 };
489 
501 template <int dim>
502 class QGaussChebyshev : public Quadrature<dim>
503 {
504 public:
506  QGaussChebyshev(const unsigned int n);
507 };
508 
509 
524 template <int dim>
525 class QGaussRadauChebyshev : public Quadrature<dim>
526 {
527 public:
528  /* EndPoint is used to specify which of the two endpoints of the unit interval
529  * is used also as quadrature point
530  */
531  enum EndPoint
532  {
540  right
541  };
543  QGaussRadauChebyshev(const unsigned int n,
545 
550  QGaussRadauChebyshev(QGaussRadauChebyshev<dim> &&) noexcept = default;
551 
552 private:
553  const EndPoint ep;
554 };
555 
569 template <int dim>
571 {
572 public:
574  QGaussLobattoChebyshev(const unsigned int n);
575 };
576 
607 template <int dim>
608 class QSimplex : public Quadrature<dim>
609 {
610 public:
617  QSimplex(const Quadrature<dim> &quad);
618 
643  compute_affine_transformation(
644  const std::array<Point<dim>, dim + 1> &vertices) const;
645 };
646 
665 class QTrianglePolar : public QSimplex<2>
666 {
667 public:
675  QTrianglePolar(const Quadrature<1> &radial_quadrature,
676  const Quadrature<1> &angular_quadrature);
677 
684  QTrianglePolar(const unsigned int n);
685 };
686 
719 class QDuffy : public QSimplex<2>
720 {
721 public:
736  QDuffy(const Quadrature<1> &radial_quadrature,
737  const Quadrature<1> &angular_quadrature,
738  const double beta = 1.0);
739 
747  QDuffy(const unsigned int n, const double beta);
748 };
749 
754 template <int dim>
755 class QSplit : public Quadrature<dim>
756 {
757 public:
792  QSplit(const QSimplex<dim> &base, const Point<dim> &split_point);
793 };
794 
797 /* -------------- declaration of explicit specializations ------------- */
798 
799 #ifndef DOXYGEN
800 template <>
801 QGauss<1>::QGauss(const unsigned int n);
802 template <>
803 QGaussLobatto<1>::QGaussLobatto(const unsigned int n);
804 
805 template <>
806 std::vector<double>
807 QGaussLog<1>::get_quadrature_points(const unsigned int);
808 template <>
809 std::vector<double>
810 QGaussLog<1>::get_quadrature_weights(const unsigned int);
811 
812 template <>
814 template <>
816 template <>
818 template <>
820 template <>
822 template <>
823 QGaussLog<1>::QGaussLog(const unsigned int n, const bool revert);
824 template <>
825 QGaussLogR<1>::QGaussLogR(const unsigned int n,
826  const Point<1> x0,
827  const double alpha,
828  const bool flag);
829 template <>
830 QGaussOneOverR<2>::QGaussOneOverR(const unsigned int n,
831  const unsigned int index,
832  const bool flag);
833 template <>
834 QTelles<1>::QTelles(const Quadrature<1> &base_quad,
835  const Point<1> & singularity);
836 #endif // DOXYGEN
837 
838 
839 
841 #endif
QGaussLog(const unsigned int n, const bool revert=false)
QGaussOneOverR(const unsigned int n, const Point< dim > singularity, const bool factor_out_singular_weight=false)
const double fraction
std::vector< double > get_quadrature_weights(const unsigned int n)
static std::vector< double > get_quadrature_points(const unsigned int n)
std::vector< double > get_quadrature_points(const unsigned int n)
QGauss(const unsigned int n)
void split_point(const Point< dim1+dim2 > &source, Point< dim1 > &p1, Point< dim2 > &p2)
static std::vector< double > get_quadrature_weights(const unsigned int n)
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:369
Point< 3 > vertices[4]
QGaussLobatto(const unsigned int n)
QTelles(const Quadrature< 1 > &base_quad, const Point< dim > &singularity)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:368
QGaussLogR(const unsigned int n, const Point< dim > x0=Point< dim >(), const double alpha=1, const bool factor_out_singular_weight=false)
#define DEAL_II_DEPRECATED
Definition: config.h:152