Reference documentation for deal.II version GIT 20f059c89a 2022-12-04 14:55:02+00:00
\(\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 - 2022 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 
39 template <int dim>
40 class QGauss : public Quadrature<dim>
41 {
42 public:
47  QGauss(const unsigned int n);
48 };
49 
50 
73 template <int dim>
74 class QGaussLobatto : public Quadrature<dim>
75 {
76 public:
81  QGaussLobatto(const unsigned int n);
82 };
83 
84 
85 
95 template <int dim>
96 class QMidpoint : public Quadrature<dim>
97 {
98 public:
99  QMidpoint();
100 };
101 
102 
107 template <int dim>
108 class QSimpson : public Quadrature<dim>
109 {
110 public:
111  QSimpson();
112 };
113 
114 
115 
124 template <int dim>
125 class QTrapezoid : public Quadrature<dim>
126 {
127 public:
128  QTrapezoid();
129 };
130 
131 
138 template <int dim>
139 class QMilne : public Quadrature<dim>
140 {
141 public:
142  QMilne();
143 };
144 
145 
152 template <int dim>
153 class QWeddle : public Quadrature<dim>
154 {
155 public:
156  QWeddle();
157 };
158 
159 
160 
174 template <int dim>
175 class QGaussLog : public Quadrature<dim>
176 {
177 public:
181  QGaussLog(const unsigned int n, const bool revert = false);
182 
183 private:
187  static std::vector<double>
188  get_quadrature_points(const unsigned int n);
189 
193  static std::vector<double>
194  get_quadrature_weights(const unsigned int n);
195 };
196 
197 
198 
235 template <int dim>
236 class QGaussLogR : public Quadrature<dim>
237 {
238 public:
246  QGaussLogR(const unsigned int n,
247  const Point<dim> & x0 = Point<dim>(),
248  const double alpha = 1,
249  const bool factor_out_singular_weight = false);
250 
256  QGaussLogR(QGaussLogR<dim> &&) noexcept = default;
257 
258 protected:
263  const double fraction;
264 };
265 
266 
290 template <int dim>
291 class QGaussOneOverR : public Quadrature<dim>
292 {
293 public:
326  QGaussOneOverR(const unsigned int n,
327  const Point<dim> & singularity,
328  const bool factor_out_singular_weight = false);
363  QGaussOneOverR(const unsigned int n,
364  const unsigned int vertex_index,
365  const bool factor_out_singular_weight = false);
366 
367 private:
373  static unsigned int
374  quad_size(const Point<dim> &singularity, const unsigned int n);
375 };
376 
377 
378 
388 template <int dim>
389 class QSorted : public Quadrature<dim>
390 {
391 public:
396  QSorted(const Quadrature<dim> &quad);
397 
398 private:
404  bool
405  compare_weights(const unsigned int a, const unsigned int b) const;
406 };
407 
462 template <int dim>
463 class QTelles : public Quadrature<dim>
464 {
465 public:
472  QTelles(const Quadrature<1> &base_quad, const Point<dim> &singularity);
478  QTelles(const unsigned int n, const Point<dim> &singularity);
479 };
480 
492 template <int dim>
493 class QGaussChebyshev : public Quadrature<dim>
494 {
495 public:
497  QGaussChebyshev(const unsigned int n);
498 };
499 
500 
515 template <int dim>
516 class QGaussRadauChebyshev : public Quadrature<dim>
517 {
518 public:
519  /* EndPoint is used to specify which of the two endpoints of the unit interval
520  * is used also as quadrature point
521  */
522  enum EndPoint
523  {
531  right
532  };
534  QGaussRadauChebyshev(const unsigned int n,
535  EndPoint ep = QGaussRadauChebyshev::left);
536 
542 
543 private:
544  const EndPoint ep;
545 };
546 
560 template <int dim>
562 {
563 public:
565  QGaussLobattoChebyshev(const unsigned int n);
566 };
567 
601 template <int dim>
602 class QSimplex : public Quadrature<dim>
603 {
604 public:
611  QSimplex(const Quadrature<dim> &quad);
612 
641  template <int spacedim = dim>
643  compute_affine_transformation(
644  const std::array<Point<spacedim>, dim + 1> &vertices) const;
645 
658  template <int spacedim = dim>
660  mapped_quadrature(
661  const std::vector<std::array<Point<spacedim>, dim + 1>> &simplices) const;
662 };
663 
682 class QTrianglePolar : public QSimplex<2>
683 {
684 public:
692  QTrianglePolar(const Quadrature<1> &radial_quadrature,
693  const Quadrature<1> &angular_quadrature);
694 
701  QTrianglePolar(const unsigned int n);
702 };
703 
736 class QDuffy : public QSimplex<2>
737 {
738 public:
753  QDuffy(const Quadrature<1> &radial_quadrature,
754  const Quadrature<1> &angular_quadrature,
755  const double beta = 1.0);
756 
764  QDuffy(const unsigned int n, const double beta);
765 };
766 
771 template <int dim>
772 class QSplit : public Quadrature<dim>
773 {
774 public:
809  QSplit(const QSimplex<dim> &base, const Point<dim> &split_point);
810 };
811 
836 template <int dim>
837 class QGaussSimplex : public QSimplex<dim>
838 {
839 public:
844  explicit QGaussSimplex(const unsigned int n_points_1D);
845 };
846 
877 template <int dim>
879 {
880 public:
887  explicit QWitherdenVincentSimplex(const unsigned int n_points_1D,
888  const bool use_odd_order = true);
889 };
890 
898 template <int dim>
899 class QIteratedSimplex : public Quadrature<dim>
900 {
901 public:
902  QIteratedSimplex(const Quadrature<dim> &base_quadrature,
903  const unsigned int n_copies);
904 };
905 
909 template <int dim>
910 class QGaussWedge : public Quadrature<dim>
911 {
912 public:
918  explicit QGaussWedge(const unsigned int n_points_1D);
919 };
920 
924 template <int dim>
925 class QGaussPyramid : public Quadrature<dim>
926 {
927 public:
933  explicit QGaussPyramid(const unsigned int n_points_1D);
934 };
935 
938 /* -------------- declaration of explicit specializations ------------- */
939 
940 #ifndef DOXYGEN
941 template <>
942 QGauss<1>::QGauss(const unsigned int n);
943 template <>
944 QGaussLobatto<1>::QGaussLobatto(const unsigned int n);
945 
946 template <>
947 std::vector<double>
948 QGaussLog<1>::get_quadrature_points(const unsigned int);
949 template <>
950 std::vector<double>
951 QGaussLog<1>::get_quadrature_weights(const unsigned int);
952 
953 template <>
955 template <>
957 template <>
959 template <>
961 template <>
963 template <>
964 QGaussLog<1>::QGaussLog(const unsigned int n, const bool revert);
965 template <>
966 QGaussLogR<1>::QGaussLogR(const unsigned int n,
967  const Point<1> & x0,
968  const double alpha,
969  const bool flag);
970 template <>
971 QGaussOneOverR<2>::QGaussOneOverR(const unsigned int n,
972  const unsigned int index,
973  const bool flag);
974 template <>
975 QTelles<1>::QTelles(const Quadrature<1> &base_quad,
976  const Point<1> & singularity);
977 #endif // DOXYGEN
978 
979 
980 
982 #endif
Definition: point.h:111
QGaussLobatto(const unsigned int n)
const double fraction
QGaussLogR(const unsigned int n, const Point< dim > &x0=Point< dim >(), const double alpha=1, const bool factor_out_singular_weight=false)
QGaussLogR(QGaussLogR< dim > &&) noexcept=default
QGaussLog(const unsigned int n, const bool revert=false)
static std::vector< double > get_quadrature_points(const unsigned int n)
static std::vector< double > get_quadrature_weights(const unsigned int n)
QGaussOneOverR(const unsigned int n, const unsigned int vertex_index, const bool factor_out_singular_weight=false)
QGaussOneOverR(const unsigned int n, const Point< dim > &singularity, const bool factor_out_singular_weight=false)
static unsigned int quad_size(const Point< dim > &singularity, const unsigned int n)
QGaussRadauChebyshev(QGaussRadauChebyshev< dim > &&) noexcept=default
QGauss(const unsigned int n)
QTelles(const Quadrature< 1 > &base_quad, const Point< dim > &singularity)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:458
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:459
Point< 3 > vertices[4]
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
void split_point(const Point< dim1+dim2 > &source, Point< dim1 > &p1, Point< dim2 > &p2)