Reference documentation for deal.II version GIT 6da2e5d553 2022-07-01 18: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 
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 
93 template <int dim>
94 class QMidpoint : public Quadrature<dim>
95 {
96 public:
97  QMidpoint();
98 };
99 
100 
105 template <int dim>
106 class QSimpson : public Quadrature<dim>
107 {
108 public:
109  QSimpson();
110 };
111 
112 
113 
122 template <int dim>
123 class QTrapezoid : public Quadrature<dim>
124 {
125 public:
126  QTrapezoid();
127 };
128 
129 
136 template <int dim>
137 class QMilne : public Quadrature<dim>
138 {
139 public:
140  QMilne();
141 };
142 
143 
150 template <int dim>
151 class QWeddle : public Quadrature<dim>
152 {
153 public:
154  QWeddle();
155 };
156 
157 
158 
172 template <int dim>
173 class QGaussLog : public Quadrature<dim>
174 {
175 public:
179  QGaussLog(const unsigned int n, const bool revert = false);
180 
181 private:
185  static std::vector<double>
186  get_quadrature_points(const unsigned int n);
187 
191  static std::vector<double>
192  get_quadrature_weights(const unsigned int n);
193 };
194 
195 
196 
233 template <int dim>
234 class QGaussLogR : public Quadrature<dim>
235 {
236 public:
244  QGaussLogR(const unsigned int n,
245  const Point<dim> & x0 = Point<dim>(),
246  const double alpha = 1,
247  const bool factor_out_singular_weight = false);
248 
254  QGaussLogR(QGaussLogR<dim> &&) noexcept = default;
255 
256 protected:
261  const double fraction;
262 };
263 
264 
288 template <int dim>
289 class QGaussOneOverR : public Quadrature<dim>
290 {
291 public:
324  QGaussOneOverR(const unsigned int n,
325  const Point<dim> & singularity,
326  const bool factor_out_singular_weight = false);
361  QGaussOneOverR(const unsigned int n,
362  const unsigned int vertex_index,
363  const bool factor_out_singular_weight = false);
364 
365 private:
371  static unsigned int
372  quad_size(const Point<dim> &singularity, const unsigned int n);
373 };
374 
375 
376 
386 template <int dim>
387 class QSorted : public Quadrature<dim>
388 {
389 public:
394  QSorted(const Quadrature<dim> &quad);
395 
396 private:
402  bool
403  compare_weights(const unsigned int a, const unsigned int b) const;
404 };
405 
460 template <int dim>
461 class QTelles : public Quadrature<dim>
462 {
463 public:
470  QTelles(const Quadrature<1> &base_quad, const Point<dim> &singularity);
476  QTelles(const unsigned int n, const Point<dim> &singularity);
477 };
478 
490 template <int dim>
491 class QGaussChebyshev : public Quadrature<dim>
492 {
493 public:
495  QGaussChebyshev(const unsigned int n);
496 };
497 
498 
513 template <int dim>
514 class QGaussRadauChebyshev : public Quadrature<dim>
515 {
516 public:
517  /* EndPoint is used to specify which of the two endpoints of the unit interval
518  * is used also as quadrature point
519  */
520  enum EndPoint
521  {
529  right
530  };
532  QGaussRadauChebyshev(const unsigned int n,
533  EndPoint ep = QGaussRadauChebyshev::left);
534 
540 
541 private:
542  const EndPoint ep;
543 };
544 
558 template <int dim>
560 {
561 public:
563  QGaussLobattoChebyshev(const unsigned int n);
564 };
565 
599 template <int dim>
600 class QSimplex : public Quadrature<dim>
601 {
602 public:
609  QSimplex(const Quadrature<dim> &quad);
610 
635  compute_affine_transformation(
636  const std::array<Point<dim>, dim + 1> &vertices) const;
637 };
638 
657 class QTrianglePolar : public QSimplex<2>
658 {
659 public:
667  QTrianglePolar(const Quadrature<1> &radial_quadrature,
668  const Quadrature<1> &angular_quadrature);
669 
676  QTrianglePolar(const unsigned int n);
677 };
678 
711 class QDuffy : public QSimplex<2>
712 {
713 public:
728  QDuffy(const Quadrature<1> &radial_quadrature,
729  const Quadrature<1> &angular_quadrature,
730  const double beta = 1.0);
731 
739  QDuffy(const unsigned int n, const double beta);
740 };
741 
746 template <int dim>
747 class QSplit : public Quadrature<dim>
748 {
749 public:
784  QSplit(const QSimplex<dim> &base, const Point<dim> &split_point);
785 };
786 
811 template <int dim>
812 class QGaussSimplex : public QSimplex<dim>
813 {
814 public:
819  explicit QGaussSimplex(const unsigned int n_points_1D);
820 };
821 
852 template <int dim>
854 {
855 public:
862  explicit QWitherdenVincentSimplex(const unsigned int n_points_1D,
863  const bool use_odd_order = true);
864 };
865 
873 template <int dim>
874 class QIteratedSimplex : public Quadrature<dim>
875 {
876 public:
877  QIteratedSimplex(const Quadrature<dim> &base_quadrature,
878  const unsigned int n_copies);
879 };
880 
884 template <int dim>
885 class QGaussWedge : public Quadrature<dim>
886 {
887 public:
893  explicit QGaussWedge(const unsigned int n_points_1D);
894 };
895 
899 template <int dim>
900 class QGaussPyramid : public Quadrature<dim>
901 {
902 public:
908  explicit QGaussPyramid(const unsigned int n_points_1D);
909 };
910 
913 /* -------------- declaration of explicit specializations ------------- */
914 
915 #ifndef DOXYGEN
916 template <>
917 QGauss<1>::QGauss(const unsigned int n);
918 template <>
919 QGaussLobatto<1>::QGaussLobatto(const unsigned int n);
920 
921 template <>
922 std::vector<double>
923 QGaussLog<1>::get_quadrature_points(const unsigned int);
924 template <>
925 std::vector<double>
926 QGaussLog<1>::get_quadrature_weights(const unsigned int);
927 
928 template <>
930 template <>
932 template <>
934 template <>
936 template <>
938 template <>
939 QGaussLog<1>::QGaussLog(const unsigned int n, const bool revert);
940 template <>
941 QGaussLogR<1>::QGaussLogR(const unsigned int n,
942  const Point<1> & x0,
943  const double alpha,
944  const bool flag);
945 template <>
946 QGaussOneOverR<2>::QGaussOneOverR(const unsigned int n,
947  const unsigned int index,
948  const bool flag);
949 template <>
950 QTelles<1>::QTelles(const Quadrature<1> &base_quad,
951  const Point<1> & singularity);
952 #endif // DOXYGEN
953 
954 
955 
957 #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:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
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)