Reference documentation for deal.II version 8.5.1
quadrature_lib.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2017 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 at
12 // the top level of the deal.II distribution.
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 #include <deal.II/base/quadrature.h>
22 
23 DEAL_II_NAMESPACE_OPEN
24 
27 
38 template <int dim>
39 class QGauss : public Quadrature<dim>
40 {
41 public:
46  QGauss (const unsigned int n);
47 };
48 
49 
74 template <int dim>
75 class QGaussLobatto : public Quadrature<dim>
76 {
77 public:
82  QGaussLobatto(const unsigned int n);
83 
84 protected:
92  std::vector<long double>
93  compute_quadrature_points (const unsigned int q,
94  const int alpha,
95  const int beta) const;
96 
104  std::vector<long double>
105  compute_quadrature_weights (const std::vector<long double> &x,
106  const int alpha,
107  const int beta) const;
108 
115  long double JacobiP(const long double x,
116  const int alpha,
117  const int beta,
118  const unsigned int n) const;
119 
124  static long double gamma(const unsigned int n);
125 };
126 
127 
128 
133 template <int dim>
134 class QMidpoint : public Quadrature<dim>
135 {
136 public:
137  QMidpoint ();
138 };
139 
140 
145 template <int dim>
146 class QSimpson : public Quadrature<dim>
147 {
148 public:
149  QSimpson ();
150 };
151 
152 
153 
166 template <int dim>
167 class QTrapez : public Quadrature<dim>
168 {
169 public:
170  QTrapez ();
171 };
172 
173 
174 
181 template <int dim>
182 class QMilne : public Quadrature<dim>
183 {
184 public:
185  QMilne ();
186 };
187 
188 
195 template <int dim>
196 class QWeddle : public Quadrature<dim>
197 {
198 public:
199  QWeddle ();
200 };
201 
202 
203 
217 template <int dim>
218 class QGaussLog : public Quadrature<dim>
219 {
220 public:
224  QGaussLog(const unsigned int n,
225  const bool revert=false);
226 
227 private:
231  static
232  std::vector<double>
233  get_quadrature_points(const unsigned int n);
234 
238  static
239  std::vector<double>
240  get_quadrature_weights(const unsigned int n);
241 };
242 
243 
244 
245 
282 template <int dim>
283 class QGaussLogR : public Quadrature<dim>
284 {
285 public:
293  QGaussLogR(const unsigned int n,
294  const Point<dim> x0 = Point<dim>(),
295  const double alpha = 1,
296  const bool factor_out_singular_weight=false);
297 
298 #ifdef DEAL_II_WITH_CXX11
299 
303  QGaussLogR(QGaussLogR<dim> &&) = default;
304 #endif
305 
306 protected:
311  const double fraction;
312 };
313 
314 
338 template <int dim>
339 class QGaussOneOverR : public Quadrature<dim>
340 {
341 public:
374  QGaussOneOverR(const unsigned int n,
375  const Point<dim> singularity,
376  const bool factor_out_singular_weight=false);
411  QGaussOneOverR(const unsigned int n,
412  const unsigned int vertex_index,
413  const bool factor_out_singular_weight=false);
414 private:
420  static unsigned int quad_size(const Point<dim> singularity,
421  const unsigned int n);
422 };
423 
424 
425 
435 template <int dim>
436 class QSorted : public Quadrature<dim>
437 {
438 public:
443  QSorted (const Quadrature<dim> &quad);
444 
445 private:
451  bool compare_weights(const unsigned int a,
452  const unsigned int b) const;
453 };
454 
511 template <int dim>
512 class QTelles: public Quadrature<dim>
513 {
514 public:
521  QTelles (const Quadrature<1> &base_quad, const Point<dim> &singularity);
527  QTelles (const unsigned int n, const Point<dim> &singularity);
528 
529 };
530 
546 template <int dim>
547 class QGaussChebyshev : public Quadrature<dim>
548 {
549 public:
551  QGaussChebyshev(const unsigned int n);
552 
553 private:
555  static std::vector<double>
556  get_quadrature_points(const unsigned int n);
557 
559  static std::vector<double>
560  get_quadrature_weights(const unsigned int n);
561 
562 };
563 
564 
581 template <int dim>
582 class QGaussRadauChebyshev : public Quadrature<dim>
583 {
584 public:
585  /* EndPoint is used to specify which of the two endpoints of the unit interval
586  * is used also as quadrature point
587  */
588  enum EndPoint
589  {
597  right
598  };
600  QGaussRadauChebyshev(const unsigned int n,
602 
603 #ifdef DEAL_II_WITH_CXX11
604 
609 #endif
610 
611 private:
612  const EndPoint ep;
614  static std::vector<double>
615  get_quadrature_points(const unsigned int n, EndPoint ep);
616 
618  static std::vector<double>
619  get_quadrature_weights(const unsigned int n, EndPoint ep);
620 
621 };
622 
638 template <int dim>
640 {
641 public:
643  QGaussLobattoChebyshev(const unsigned int n);
644 
645 private:
647  static std::vector<double>
648  get_quadrature_points(const unsigned int n);
649 
651  static std::vector<double>
652  get_quadrature_weights(const unsigned int n);
653 
654 };
655 
656 /* -------------- declaration of explicit specializations ------------- */
657 
658 template <> QGauss<1>::QGauss (const unsigned int n);
659 template <> QGaussLobatto<1>::QGaussLobatto (const unsigned int n);
660 template <>
661 std::vector<long double> QGaussLobatto<1>::
662 compute_quadrature_points(const unsigned int, const int, const int) const;
663 template <>
664 std::vector<long double> QGaussLobatto<1>::
665 compute_quadrature_weights(const std::vector<long double> &, const int, const int) const;
666 template <>
667 long double QGaussLobatto<1>::
668 JacobiP(const long double, const int, const int, const unsigned int) const;
669 template <>
670 long double
671 QGaussLobatto<1>::gamma(const unsigned int n);
672 
673 template <> std::vector<double> QGaussLog<1>::get_quadrature_points(const unsigned int);
674 template <> std::vector<double> QGaussLog<1>::get_quadrature_weights(const unsigned int);
675 
676 template <> QMidpoint<1>::QMidpoint ();
677 template <> QTrapez<1>::QTrapez ();
678 template <> QSimpson<1>::QSimpson ();
679 template <> QMilne<1>::QMilne ();
680 template <> QWeddle<1>::QWeddle ();
681 template <> QGaussLog<1>::QGaussLog (const unsigned int n, const bool revert);
682 template <> QGaussLogR<1>::QGaussLogR (const unsigned int n, const Point<1> x0, const double alpha, const bool flag);
683 template <> QGaussOneOverR<2>::QGaussOneOverR (const unsigned int n, const unsigned int index, const bool flag);
684 template <> QTelles<1>::QTelles(const Quadrature<1> &base_quad, const Point<1> &singularity);
685 
686 
687 
688 DEAL_II_NAMESPACE_CLOSE
689 
690 #endif
static long double gamma(const unsigned int n)
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
static std::vector< double > get_quadrature_points(const unsigned int n)
QGauss(const unsigned int n)
long double JacobiP(const long double x, const int alpha, const int beta, const unsigned int n) const
std::vector< long double > compute_quadrature_points(const unsigned int q, const int alpha, const int beta) const
static std::vector< double > get_quadrature_weights(const unsigned int n)
QGaussLobatto(const unsigned int n)
std::vector< long double > compute_quadrature_weights(const std::vector< long double > &x, const int alpha, const int beta) const
QTelles(const Quadrature< 1 > &base_quad, const Point< dim > &singularity)
QGaussLogR(const unsigned int n, const Point< dim > x0=Point< dim >(), const double alpha=1, const bool factor_out_singular_weight=false)