Reference documentation for deal.II version 8.1.0
quadrature_lib.h
1 // ---------------------------------------------------------------------
2 // @f$Id: quadrature_lib.h 31792 2013-11-25 10:36:42Z felix.gruber @f$
3 //
4 // Copyright (C) 1998 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__quadrature_lib_h
18 #define __deal2__quadrature_lib_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/quadrature.h>
23 
25 
28 
37 template <int dim>
38 class QGauss : public Quadrature<dim>
39 {
40 public:
48  QGauss (const unsigned int n);
49 };
50 
51 
78 template<int dim>
79 class QGaussLobatto : public Quadrature<dim>
80 {
81 public:
87  QGaussLobatto(const unsigned int n);
88 
89 protected:
102  std::vector<long double>
103  compute_quadrature_points (const unsigned int q,
104  const int alpha,
105  const int beta) const;
106 
116  std::vector<long double>
117  compute_quadrature_weights (const std::vector<long double> &x,
118  const int alpha,
119  const int beta) const;
120 
131  long double JacobiP(const long double x,
132  const int alpha,
133  const int beta,
134  const unsigned int n) const;
135 
141  long double gamma(const unsigned int n) const;
142 };
143 
144 
145 
149 template <int dim>
150 class QMidpoint : public Quadrature<dim>
151 {
152 public:
153  QMidpoint ();
154 };
155 
156 
160 template <int dim>
161 class QSimpson : public Quadrature<dim>
162 {
163 public:
164  QSimpson ();
165 };
166 
167 
171 template <int dim>
172 class QTrapez : public Quadrature<dim>
173 {
174 public:
175  QTrapez ();
176 };
177 
182 template <int dim>
183 class QMilne : public Quadrature<dim>
184 {
185 public:
186  QMilne ();
187 };
188 
189 
194 template <int dim>
195 class QWeddle : public Quadrature<dim>
196 {
197 public:
198  QWeddle ();
199 };
200 
201 
202 
218 template <int dim>
219 class QGaussLog : public Quadrature<dim>
220 {
221 public:
226  QGaussLog(const unsigned int n,
227  const bool revert=false);
228 
229 protected:
234  std::vector<double>
235  set_quadrature_points(const unsigned int n) const;
236 
241  std::vector<double>
242  set_quadrature_weights(const unsigned int n) const;
243 
244 };
245 
246 
247 
248 
292 template<int dim>
293 class QGaussLogR : public Quadrature<dim>
294 {
295 public:
308  QGaussLogR(const unsigned int n,
309  const Point<dim> x0 = Point<dim>(),
310  const double alpha = 1,
311  const bool factor_out_singular_weight=false);
312 
313 protected:
319  const double fraction;
320 };
321 
322 
347 template<int dim>
348 class QGaussOneOverR : public Quadrature<dim>
349 {
350 public:
385  QGaussOneOverR(const unsigned int n,
386  const Point<dim> singularity,
387  const bool factor_out_singular_weight=false);
424  QGaussOneOverR(const unsigned int n,
425  const unsigned int vertex_index,
426  const bool factor_out_singular_weight=false);
427 private:
432  static unsigned int quad_size(const Point<dim> singularity,
433  const unsigned int n);
434 };
435 
436 
437 
440 /* -------------- declaration of explicit specializations ------------- */
441 
442 template <> QGauss<1>::QGauss (const unsigned int n);
443 template <> QGaussLobatto<1>::QGaussLobatto (const unsigned int n);
444 template <>
445 std::vector<long double> QGaussLobatto<1>::
446 compute_quadrature_points(const unsigned int, const int, const int) const;
447 template <>
448 std::vector<long double> QGaussLobatto<1>::
449 compute_quadrature_weights(const std::vector<long double> &, const int, const int) const;
450 template <>
451 long double QGaussLobatto<1>::
452 JacobiP(const long double, const int, const int, const unsigned int) const;
453 template <>
454 long double
455 QGaussLobatto<1>::gamma(const unsigned int n) const;
456 
457 template <> std::vector<double> QGaussLog<1>::set_quadrature_points(const unsigned int) const;
458 template <> std::vector<double> QGaussLog<1>::set_quadrature_weights(const unsigned int) const;
459 
460 template <> QMidpoint<1>::QMidpoint ();
461 template <> QTrapez<1>::QTrapez ();
462 template <> QSimpson<1>::QSimpson ();
463 template <> QMilne<1>::QMilne ();
464 template <> QWeddle<1>::QWeddle ();
465 template <> QGaussLog<1>::QGaussLog (const unsigned int n, const bool revert);
466 template <> QGaussLogR<1>::QGaussLogR (const unsigned int n, const Point<1> x0, const double alpha, const bool flag);
467 template <> QGaussOneOverR<2>::QGaussOneOverR (const unsigned int n, const unsigned int index, const bool flag);
468 
469 
470 
471 
472 DEAL_II_NAMESPACE_CLOSE
473 
474 #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)
long double gamma(const unsigned int n) const
const double fraction
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
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
std::vector< double > set_quadrature_weights(const unsigned int n) const
QGaussLogR(const unsigned int n, const Point< dim > x0=Point< dim >(), const double alpha=1, const bool factor_out_singular_weight=false)
std::vector< double > set_quadrature_points(const unsigned int n) const