Reference documentation for deal.II version Git 689de043d4 2020-08-10 16:46:15 +0200
\(\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 
119 template <int dim>
120 class QTrapez : public Quadrature<dim>
121 {
122 public:
123  QTrapez();
124 };
125 
126 
127 
134 template <int dim>
135 class QMilne : public Quadrature<dim>
136 {
137 public:
138  QMilne();
139 };
140 
141 
148 template <int dim>
149 class QWeddle : public Quadrature<dim>
150 {
151 public:
152  QWeddle();
153 };
154 
155 
156 
170 template <int dim>
171 class QGaussLog : public Quadrature<dim>
172 {
173 public:
177  QGaussLog(const unsigned int n, const bool revert = false);
178 
179 private:
183  static std::vector<double>
184  get_quadrature_points(const unsigned int n);
185 
189  static std::vector<double>
190  get_quadrature_weights(const unsigned int n);
191 };
192 
193 
194 
231 template <int dim>
232 class QGaussLogR : public Quadrature<dim>
233 {
234 public:
242  QGaussLogR(const unsigned int n,
243  const Point<dim> x0 = Point<dim>(),
244  const double alpha = 1,
245  const bool factor_out_singular_weight = false);
246 
252  QGaussLogR(QGaussLogR<dim> &&) noexcept = default;
253 
254 protected:
259  const double fraction;
260 };
261 
262 
286 template <int dim>
287 class QGaussOneOverR : public Quadrature<dim>
288 {
289 public:
322  QGaussOneOverR(const unsigned int n,
323  const Point<dim> singularity,
324  const bool factor_out_singular_weight = false);
359  QGaussOneOverR(const unsigned int n,
360  const unsigned int vertex_index,
361  const bool factor_out_singular_weight = false);
362 
363 private:
369  static unsigned int
370  quad_size(const Point<dim> singularity, const unsigned int n);
371 };
372 
373 
374 
384 template <int dim>
385 class QSorted : public Quadrature<dim>
386 {
387 public:
392  QSorted(const Quadrature<dim> &quad);
393 
394 private:
400  bool
401  compare_weights(const unsigned int a, const unsigned int b) const;
402 };
403 
458 template <int dim>
459 class QTelles : public Quadrature<dim>
460 {
461 public:
468  QTelles(const Quadrature<1> &base_quad, const Point<dim> &singularity);
474  QTelles(const unsigned int n, const Point<dim> &singularity);
475 };
476 
488 template <int dim>
489 class QGaussChebyshev : public Quadrature<dim>
490 {
491 public:
493  QGaussChebyshev(const unsigned int n);
494 };
495 
496 
511 template <int dim>
512 class QGaussRadauChebyshev : public Quadrature<dim>
513 {
514 public:
515  /* EndPoint is used to specify which of the two endpoints of the unit interval
516  * is used also as quadrature point
517  */
518  enum EndPoint
519  {
527  right
528  };
530  QGaussRadauChebyshev(const unsigned int n,
532 
537  QGaussRadauChebyshev(QGaussRadauChebyshev<dim> &&) noexcept = default;
538 
539 private:
540  const EndPoint ep;
541 };
542 
556 template <int dim>
558 {
559 public:
561  QGaussLobattoChebyshev(const unsigned int n);
562 };
563 
594 template <int dim>
595 class QSimplex : public Quadrature<dim>
596 {
597 public:
604  QSimplex(const Quadrature<dim> &quad);
605 
630  compute_affine_transformation(
631  const std::array<Point<dim>, dim + 1> &vertices) const;
632 };
633 
652 class QTrianglePolar : public QSimplex<2>
653 {
654 public:
662  QTrianglePolar(const Quadrature<1> &radial_quadrature,
663  const Quadrature<1> &angular_quadrature);
664 
671  QTrianglePolar(const unsigned int n);
672 };
673 
706 class QDuffy : public QSimplex<2>
707 {
708 public:
723  QDuffy(const Quadrature<1> &radial_quadrature,
724  const Quadrature<1> &angular_quadrature,
725  const double beta = 1.0);
726 
734  QDuffy(const unsigned int n, const double beta);
735 };
736 
741 template <int dim>
742 class QSplit : public Quadrature<dim>
743 {
744 public:
779  QSplit(const QSimplex<dim> &base, const Point<dim> &split_point);
780 };
781 
784 /* -------------- declaration of explicit specializations ------------- */
785 
786 #ifndef DOXYGEN
787 template <>
788 QGauss<1>::QGauss(const unsigned int n);
789 template <>
790 QGaussLobatto<1>::QGaussLobatto(const unsigned int n);
791 
792 template <>
793 std::vector<double>
794 QGaussLog<1>::get_quadrature_points(const unsigned int);
795 template <>
796 std::vector<double>
797 QGaussLog<1>::get_quadrature_weights(const unsigned int);
798 
799 template <>
801 template <>
803 template <>
805 template <>
807 template <>
809 template <>
810 QGaussLog<1>::QGaussLog(const unsigned int n, const bool revert);
811 template <>
812 QGaussLogR<1>::QGaussLogR(const unsigned int n,
813  const Point<1> x0,
814  const double alpha,
815  const bool flag);
816 template <>
817 QGaussOneOverR<2>::QGaussOneOverR(const unsigned int n,
818  const unsigned int index,
819  const bool flag);
820 template <>
821 QTelles<1>::QTelles(const Quadrature<1> &base_quad,
822  const Point<1> & singularity);
823 #endif // DOXYGEN
824 
825 
826 
828 #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:363
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:362
QGaussLogR(const unsigned int n, const Point< dim > x0=Point< dim >(), const double alpha=1, const bool factor_out_singular_weight=false)