Reference documentation for deal.II version GIT c00406db16 2023-09-28 18:00: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 - 2023 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 
79 template <int dim>
80 class QGaussRadau : public Quadrature<dim>
81 {
82 public:
87  enum EndPoint
88  {
96  right
97  };
104  QGaussRadau(const unsigned int n,
105  const EndPoint end_point = QGaussRadau::EndPoint::left);
106 
110  QGaussRadau(QGaussRadau<dim> &&) noexcept = default;
111 
112 private:
114 };
115 
116 
139 template <int dim>
140 class QGaussLobatto : public Quadrature<dim>
141 {
142 public:
147  QGaussLobatto(const unsigned int n);
148 };
149 
150 
151 
161 template <int dim>
162 class QMidpoint : public Quadrature<dim>
163 {
164 public:
165  QMidpoint();
166 };
167 
168 
173 template <int dim>
174 class QSimpson : public Quadrature<dim>
175 {
176 public:
177  QSimpson();
178 };
179 
180 
181 
190 template <int dim>
191 class QTrapezoid : public Quadrature<dim>
192 {
193 public:
194  QTrapezoid();
195 };
196 
197 
204 template <int dim>
205 class QMilne : public Quadrature<dim>
206 {
207 public:
208  QMilne();
209 };
210 
211 
218 template <int dim>
219 class QWeddle : public Quadrature<dim>
220 {
221 public:
222  QWeddle();
223 };
224 
225 
226 
240 template <int dim>
241 class QGaussLog : public Quadrature<dim>
242 {
243 public:
247  QGaussLog(const unsigned int n, const bool revert = false);
248 
249 private:
253  static std::vector<double>
254  get_quadrature_points(const unsigned int n);
255 
259  static std::vector<double>
260  get_quadrature_weights(const unsigned int n);
261 };
262 
263 
264 
301 template <int dim>
302 class QGaussLogR : public Quadrature<dim>
303 {
304 public:
312  QGaussLogR(const unsigned int n,
313  const Point<dim> &x0 = Point<dim>(),
314  const double alpha = 1,
315  const bool factor_out_singular_weight = false);
316 
322  QGaussLogR(QGaussLogR<dim> &&) noexcept = default;
323 
324 protected:
329  const double fraction;
330 };
331 
332 
356 template <int dim>
357 class QGaussOneOverR : public Quadrature<dim>
358 {
359 public:
392  QGaussOneOverR(const unsigned int n,
393  const Point<dim> &singularity,
394  const bool factor_out_singular_weight = false);
429  QGaussOneOverR(const unsigned int n,
430  const unsigned int vertex_index,
431  const bool factor_out_singular_weight = false);
432 
433 private:
439  static unsigned int
440  quad_size(const Point<dim> &singularity, const unsigned int n);
441 };
442 
443 
444 
454 template <int dim>
455 class QSorted : public Quadrature<dim>
456 {
457 public:
462  QSorted(const Quadrature<dim> &quad);
463 
464 private:
470  bool
471  compare_weights(const unsigned int a, const unsigned int b) const;
472 };
473 
528 template <int dim>
529 class QTelles : public Quadrature<dim>
530 {
531 public:
538  QTelles(const Quadrature<1> &base_quad, const Point<dim> &singularity);
544  QTelles(const unsigned int n, const Point<dim> &singularity);
545 };
546 
558 template <int dim>
559 class QGaussChebyshev : public Quadrature<dim>
560 {
561 public:
563  QGaussChebyshev(const unsigned int n);
564 };
565 
566 
581 template <int dim>
582 class QGaussRadauChebyshev : public Quadrature<dim>
583 {
584 public:
589  enum EndPoint
590  {
598  right
599  };
602  const unsigned int n,
603  const EndPoint end_point = QGaussRadauChebyshev::EndPoint::left);
604 
609 
610 private:
612 };
613 
627 template <int dim>
629 {
630 public:
632  QGaussLobattoChebyshev(const unsigned int n);
633 };
634 
668 template <int dim>
669 class QSimplex : public Quadrature<dim>
670 {
671 public:
678  QSimplex(const Quadrature<dim> &quad);
679 
708  template <int spacedim = dim>
710  compute_affine_transformation(
711  const std::array<Point<spacedim>, dim + 1> &vertices) const;
712 
725  template <int spacedim = dim>
727  mapped_quadrature(
728  const std::vector<std::array<Point<spacedim>, dim + 1>> &simplices) const;
729 };
730 
749 class QTrianglePolar : public QSimplex<2>
750 {
751 public:
759  QTrianglePolar(const Quadrature<1> &radial_quadrature,
760  const Quadrature<1> &angular_quadrature);
761 
768  QTrianglePolar(const unsigned int n);
769 };
770 
803 class QDuffy : public QSimplex<2>
804 {
805 public:
820  QDuffy(const Quadrature<1> &radial_quadrature,
821  const Quadrature<1> &angular_quadrature,
822  const double beta = 1.0);
823 
831  QDuffy(const unsigned int n, const double beta);
832 };
833 
838 template <int dim>
839 class QSplit : public Quadrature<dim>
840 {
841 public:
876  QSplit(const QSimplex<dim> &base, const Point<dim> &split_point);
877 };
878 
904 template <int dim>
905 class QGaussSimplex : public QSimplex<dim>
906 {
907 public:
912  explicit QGaussSimplex(const unsigned int n_points_1D);
913 };
914 
946 template <int dim>
948 {
949 public:
956  explicit QWitherdenVincentSimplex(const unsigned int n_points_1D,
957  const bool use_odd_order = true);
958 };
959 
967 template <int dim>
968 class QIteratedSimplex : public Quadrature<dim>
969 {
970 public:
971  QIteratedSimplex(const Quadrature<dim> &base_quadrature,
972  const unsigned int n_copies);
973 };
974 
978 template <int dim>
979 class QGaussWedge : public Quadrature<dim>
980 {
981 public:
987  explicit QGaussWedge(const unsigned int n_points_1D);
988 };
989 
993 template <int dim>
994 class QGaussPyramid : public Quadrature<dim>
995 {
996 public:
1002  explicit QGaussPyramid(const unsigned int n_points_1D);
1003 };
1004 
1007 /* -------------- declaration of explicit specializations ------------- */
1008 
1009 #ifndef DOXYGEN
1010 template <>
1011 QGauss<1>::QGauss(const unsigned int n);
1012 template <>
1013 QGaussRadau<1>::QGaussRadau(const unsigned int n,
1015 template <>
1016 QGaussLobatto<1>::QGaussLobatto(const unsigned int n);
1017 
1018 template <>
1019 std::vector<double>
1020 QGaussLog<1>::get_quadrature_points(const unsigned int);
1021 template <>
1022 std::vector<double>
1023 QGaussLog<1>::get_quadrature_weights(const unsigned int);
1024 
1025 template <>
1027 template <>
1029 template <>
1031 template <>
1033 template <>
1035 template <>
1036 QGaussLog<1>::QGaussLog(const unsigned int n, const bool revert);
1037 template <>
1038 QGaussLogR<1>::QGaussLogR(const unsigned int n,
1039  const Point<1> &x0,
1040  const double alpha,
1041  const bool flag);
1042 template <>
1043 QGaussOneOverR<2>::QGaussOneOverR(const unsigned int n,
1044  const unsigned int index,
1045  const bool flag);
1046 template <>
1047 QTelles<1>::QTelles(const Quadrature<1> &base_quad,
1048  const Point<1> &singularity);
1049 #endif // DOXYGEN
1050 
1051 
1052 
1054 #endif
Definition: point.h:112
QGaussLobatto(const unsigned int n)
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
const EndPoint end_point
QGaussRadau(QGaussRadau< dim > &&) noexcept=default
QGaussRadau(const unsigned int n, const EndPoint end_point=QGaussRadau::EndPoint::left)
QGauss(const unsigned int n)
QTelles(const Quadrature< 1 > &base_quad, const Point< dim > &singularity)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
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)