]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove the QGauss2...QGauss7 classes after 6.5 years of deprecation.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 19 Sep 2010 00:36:17 +0000 (00:36 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 19 Sep 2010 00:36:17 +0000 (00:36 +0000)
git-svn-id: https://svn.dealii.org/trunk@22041 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/quadrature_lib.h
deal.II/base/source/quadrature_lib.cc
deal.II/base/source/quadrature_selector.cc
deal.II/doc/news/changes.h

index c51533b26c80f240f445ccbf18a1a838c159b7f4..21f2e3241e0ad90660fa9cd76eb3501f44615fde 100644 (file)
@@ -139,118 +139,6 @@ class QGaussLobatto : public Quadrature<dim>
 
 
 
-/**
- * @deprecated Use QGauss for arbitrary order Gauss formulae instead!
- *
- *  2-Point-Gauss quadrature formula, exact for polynomials of degree 3.
- *
- *  Reference: Ward Cheney, David Kincaid: "Numerical Mathematics and Computing".
- *  For a comprehensive list of Gaussian quadrature formulae, see also:
- *  A. H. Strout, D. Secrest: "Gaussian Quadrature Formulas"
- */
-template <int dim>
-class QGauss2 : public Quadrature<dim>
-{
-  public:
-    QGauss2 ();
-};
-
-
-/**
- * @deprecated Use QGauss for arbitrary order Gauss formulae instead!
- *
- *  3-Point-Gauss quadrature formula, exact for polynomials of degree 5.
- *
- *  Reference: Ward Cheney, David Kincaid: "Numerical Mathematics and Computing".
- *  For a comprehensive list of Gaussian quadrature formulae, see also:
- *  A. H. Strout, D. Secrest: "Gaussian Quadrature Formulas"
- */
-template <int dim>
-class QGauss3 : public Quadrature<dim>
-{
-  public:
-    QGauss3 ();
-};
-
-
-/**
- * @deprecated Use QGauss for arbitrary order Gauss formulae instead!
- *
- * 4-Point-Gauss quadrature formula, exact for polynomials of degree 7.
- *
- *  Reference: Ward Cheney, David Kincaid: "Numerical Mathematics and Computing".
- *  For a comprehensive list of Gaussian quadrature formulae, see also:
- *  A. H. Strout, D. Secrest: "Gaussian Quadrature Formulas"
- */
-template <int dim>
-class QGauss4 : public Quadrature<dim>
-{
-  public:
-    QGauss4 ();
-};
-
-
-/**
- * @deprecated Use QGauss for arbitrary order Gauss formulae instead!
- *
- *  5-Point-Gauss quadrature formula, exact for polynomials of degree 9.
- *
- *  Reference: Ward Cheney, David Kincaid: "Numerical Mathematics and Computing".
- *  For a comprehensive list of Gaussian quadrature formulae, see also:
- *  A. H. Strout, D. Secrest: "Gaussian Quadrature Formulas"
- */
-template <int dim>
-class QGauss5 : public Quadrature<dim>
-{
-  public:
-    QGauss5 ();
-};
-
-
-/**
- * @deprecated Use QGauss for arbitrary order Gauss formulae instead!
- *
- *  6-Point-Gauss quadrature formula, exact for polynomials of degree 11.
- *  We have not found explicit
- *  representations of the zeros of the Legendre functions of sixth
- *  and higher degree. If anyone finds them, please replace the existing
- *  numbers by these expressions.
- *
- *  Reference: J. E. Akin: "Application and Implementation of Finite
- *  Element Methods"
- *  For a comprehensive list of Gaussian quadrature formulae, see also:
- *  A. H. Strout, D. Secrest: "Gaussian Quadrature Formulas"
- */
-template <int dim>
-class QGauss6 : public Quadrature<dim>
-{
-  public:
-    QGauss6 ();
-};
-
-
-/**
- * @deprecated Use QGauss for arbitrary order Gauss formulae instead!
- *
- *  7-Point-Gauss quadrature formula, exact for polynomials of degree 13.
- *  We have not found explicit
- *  representations of the zeros of the Legendre functions of sixth
- *  and higher degree. If anyone finds them, please replace the existing
- *  numbers by these expressions.
- *
- *  Reference: J. E. Akin: "Application and Implementation of Finite
- *  Element Methods"
- *  For a comprehensive list of Gaussian quadrature formulae, see also:
- *  A. H. Strout, D. Secrest: "Gaussian Quadrature Formulas"
- */
-template <int dim>
-class QGauss7 : public Quadrature<dim>
-{
-  public:
-    QGauss7 ();
-};
-
-
 /**
  * Midpoint quadrature rule, exact for linear polynomials.
  */
@@ -565,12 +453,6 @@ QGaussLobatto<1>::gamma(const unsigned int n) const;
 template <> std::vector<double> QGaussLog<1>::set_quadrature_points(const unsigned int) const;
 template <> std::vector<double> QGaussLog<1>::set_quadrature_weights(const unsigned int) const;
 
-template <> QGauss2<1>::QGauss2 ();
-template <> QGauss3<1>::QGauss3 ();
-template <> QGauss4<1>::QGauss4 ();
-template <> QGauss5<1>::QGauss5 ();
-template <> QGauss6<1>::QGauss6 ();
-template <> QGauss7<1>::QGauss7 ();
 template <> QMidpoint<1>::QMidpoint ();
 template <> QTrapez<1>::QTrapez ();
 template <> QSimpson<1>::QSimpson ();
index caf20b2350034e70bd008553ae231d32ee49b6b1..231abf308a1e7302f586e7a16ada07387b9ae368 100644 (file)
@@ -347,225 +347,6 @@ long double QGaussLobatto<1>::gamma(const unsigned int n) const
 
 
 
-template <>
-QGauss2<1>::QGauss2 ()
-                :
-               Quadrature<1> (2)
-{
-                                  // points on [-1,1]
-  static const double xpts_normal[] = { -std::sqrt(1./3.), std::sqrt(1./3.) };
-                                  // weights on [-1,1]
-  static const double wts_normal[]  = { 1., 1. };
-
-                                  // points and weights on [0,1]
-  static const double xpts[] = { (xpts_normal[0]+1)/2.,
-                                (xpts_normal[1]+1)/2.  };
-  static const double wts[]  = { wts_normal[0]/2.,
-                                wts_normal[1]/2.  };
-
-  for (unsigned int i=0; i<this->size(); ++i)
-    {
-      this->quadrature_points[i] = Point<1>(xpts[i]);
-      this->weights[i] = wts[i];
-    };
-}
-
-
-
-template <>
-QGauss3<1>::QGauss3 ()
-                :
-               Quadrature<1> (3)
-{
-                                  // points on [-1,1]
-  static const double xpts_normal[] = { -std::sqrt(3./5.),
-                                       0.,
-                                       std::sqrt(3./5.) };
-                                  // weights on [-1,1]
-  static const double wts_normal[]  = { 5./9.,
-                                       8./9.,
-                                       5./9. };
-
-                                  // points and weights on [0,1]
-  static const double xpts[] = { (xpts_normal[0]+1)/2.,
-                                (xpts_normal[1]+1)/2.,
-                                (xpts_normal[2]+1)/2. };
-  static const double wts[]  = { wts_normal[0]/2.,
-                                wts_normal[1]/2.,
-                                wts_normal[2]/2. };
-
-  for (unsigned int i=0; i<this->size(); ++i)
-    {
-      this->quadrature_points[i] = Point<1>(xpts[i]);
-      this->weights[i] = wts[i];
-    };
-}
-
-
-
-template <>
-QGauss4<1>::QGauss4 ()
-                :
-               Quadrature<1> (4)
-{
-                                  // points on [-1,1]
-  static const double xpts_normal[] = { -std::sqrt(1./7.*(3+4*std::sqrt(0.3))),
-                                       -std::sqrt(1./7.*(3-4*std::sqrt(0.3))),
-                                       +std::sqrt(1./7.*(3-4*std::sqrt(0.3))),
-                                       +std::sqrt(1./7.*(3+4*std::sqrt(0.3)))  };
-                                  // weights on [-1,1]
-  static const double wts_normal[]  = { 1./2. - 1./12.*std::sqrt(10./3.),
-                                       1./2. + 1./12.*std::sqrt(10./3.),
-                                       1./2. + 1./12.*std::sqrt(10./3.),
-                                       1./2. - 1./12.*std::sqrt(10./3.)  };
-
-                                  // points and weights on [0,1]
-  static const double xpts[] = { (xpts_normal[0]+1)/2.,
-                                (xpts_normal[1]+1)/2.,
-                                (xpts_normal[2]+1)/2.,
-                                (xpts_normal[3]+1)/2. };
-  static const double wts[]  = { wts_normal[0]/2.,
-                                wts_normal[1]/2.,
-                                wts_normal[2]/2.,
-                                wts_normal[3]/2. };
-
-  for (unsigned int i=0; i<this->size(); ++i)
-    {
-      this->quadrature_points[i] = Point<1>(xpts[i]);
-      this->weights[i] = wts[i];
-    };
-}
-
-
-template <>
-QGauss5<1>::QGauss5 ()
-                :
-               Quadrature<1> (5)
-{
-                                  // points on [-1,1]
-  static const double xpts_normal[] = { -std::sqrt(1./9.*(5.+2*std::sqrt(10./7.))),
-                                       -std::sqrt(1./9.*(5.-2*std::sqrt(10./7.))),
-                                       0,
-                                       +std::sqrt(1./9.*(5.-2*std::sqrt(10./7.))),
-                                       +std::sqrt(1./9.*(5.+2*std::sqrt(10./7.)))  };
-                                  // weights on [-1,1]
-  static const double wts_normal[]  = { 0.3*(+0.7+5.*std::sqrt(0.7))/(+2.+5.*std::sqrt(0.7)),
-                                       0.3*(-0.7+5.*std::sqrt(0.7))/(-2.+5.*std::sqrt(0.7)),
-                                       128./225.,
-                                       0.3*(-0.7+5.*std::sqrt(0.7))/(-2.+5.*std::sqrt(0.7)),
-                                       0.3*(+0.7+5.*std::sqrt(0.7))/(+2.+5.*std::sqrt(0.7)) };
-
-                                  // points and weights on [0,1]
-  static const double xpts[] = { (xpts_normal[0]+1)/2.,
-                                (xpts_normal[1]+1)/2.,
-                                (xpts_normal[2]+1)/2.,
-                                (xpts_normal[3]+1)/2.,
-                                (xpts_normal[4]+1)/2. };
-  static const double wts[]  = { wts_normal[0]/2.,
-                                wts_normal[1]/2.,
-                                wts_normal[2]/2.,
-                                wts_normal[3]/2.,
-                                wts_normal[4]/2. };
-
-  for (unsigned int i=0; i<this->size(); ++i)
-    {
-      this->quadrature_points[i] = Point<1>(xpts[i]);
-      this->weights[i] = wts[i];
-    };
-}
-
-
-
-template <>
-QGauss6<1>::QGauss6 ()
-                :
-               Quadrature<1> (6)
-{
-                                  // points on [-1,1]
-  static const double xpts_normal[] = { -0.932469514203152,
-                                       -0.661209386466265,
-                                       -0.238619186083197,
-                                       +0.238619186083197,
-                                       +0.661209386466265,
-                                       +0.932469514203152  };
-                                  // weights on [-1,1]
-  static const double wts_normal[]  = { 0.171324492379170,
-                                       0.360761573048139,
-                                       0.467913934572691,
-                                       0.467913934572691,
-                                       0.360761573048139,
-                                       0.171324492379170  };
-
-                                  // points and weights on [0,1]
-  static const double xpts[] = { (xpts_normal[0]+1)/2.,
-                                (xpts_normal[1]+1)/2.,
-                                (xpts_normal[2]+1)/2.,
-                                (xpts_normal[3]+1)/2.,
-                                (xpts_normal[4]+1)/2.,
-                                (xpts_normal[5]+1)/2. };
-  static const double wts[]  = { wts_normal[0]/2.,
-                                wts_normal[1]/2.,
-                                wts_normal[2]/2.,
-                                wts_normal[3]/2.,
-                                wts_normal[4]/2.,
-                                wts_normal[5]/2. };
-
-  for (unsigned int i=0; i<this->size(); ++i)
-    {
-      this->quadrature_points[i] = Point<1>(xpts[i]);
-      this->weights[i] = wts[i];
-    };
-}
-
-
-
-template <>
-QGauss7<1>::QGauss7 ()
-                :
-               Quadrature<1> (7)
-{
-                                  // points on [-1,1]
-  static const double xpts_normal[] = { -0.949107912342759,
-                                       -0.741531185599394,
-                                       -0.405845151377397,
-                                       0,
-                                       +0.405845151377397,
-                                       +0.741531185599394,
-                                       +0.949107912342759 };
-                                  // weights on [-1,1]
-  static const double wts_normal[]  = { 0.129484966168870,
-                                       0.279705391489277,
-                                       0.381830050505119,
-                                       0.417959183673469,
-                                       0.381830050505119,
-                                       0.279705391489277,
-                                       0.129484966168870  };
-
-                                  // points and weights on [0,1]
-  static const double xpts[] = { (xpts_normal[0]+1)/2.,
-                                (xpts_normal[1]+1)/2.,
-                                (xpts_normal[2]+1)/2.,
-                                (xpts_normal[3]+1)/2.,
-                                (xpts_normal[4]+1)/2.,
-                                (xpts_normal[5]+1)/2.,
-                                (xpts_normal[6]+1)/2. };
-  static const double wts[]  = { wts_normal[0]/2.,
-                                wts_normal[1]/2.,
-                                wts_normal[2]/2.,
-                                wts_normal[3]/2.,
-                                wts_normal[4]/2.,
-                                wts_normal[5]/2.,
-                                wts_normal[6]/2. };
-
-  for (unsigned int i=0; i<this->size(); ++i)
-    {
-      this->quadrature_points[i] = Point<1>(xpts[i]);
-      this->weights[i] = wts[i];
-    };
-}
-
-
-
 template <>
 QMidpoint<1>::QMidpoint ()
                 :
@@ -1173,54 +954,6 @@ QGaussLobatto<dim>::QGaussLobatto (const unsigned int n)
 
 
 
-template <int dim>
-QGauss2<dim>::QGauss2 ()
-                :
-                Quadrature<dim> (QGauss2<dim-1>(), QGauss2<1>())
-{}
-
-
-
-template <int dim>
-QGauss3<dim>::QGauss3 ()
-                :
-               Quadrature<dim> (QGauss3<dim-1>(), QGauss3<1>())
-{}
-
-
-
-template <int dim>
-QGauss4<dim>::QGauss4 ()
-                :
-               Quadrature<dim> (QGauss4<dim-1>(), QGauss4<1>())
-{}
-
-
-
-template <int dim>
-QGauss5<dim>::QGauss5 ()
-                :
-               Quadrature<dim> (QGauss5<dim-1>(), QGauss5<1>())
-{}
-
-
-
-template <int dim>
-QGauss6<dim>::QGauss6 ()
-                :
-               Quadrature<dim> (QGauss6<dim-1>(), QGauss6<1>())
-{}
-
-
-
-template <int dim>
-QGauss7<dim>::QGauss7 ()
-                :
-               Quadrature<dim> (QGauss7<dim-1>(), QGauss7<1>())
-{}
-
-
-
 template <int dim>
 QMidpoint<dim>::QMidpoint ()
                 :
@@ -1263,12 +996,6 @@ QWeddle<dim>::QWeddle ()
 // note that 1d formulae are specialized by implementation above
 template class QGauss<2>;
 template class QGaussLobatto<2>;
-template class QGauss2<2>;
-template class QGauss3<2>;
-template class QGauss4<2>;
-template class QGauss5<2>;
-template class QGauss6<2>;
-template class QGauss7<2>;
 template class QMidpoint<2>;
 template class QTrapez<2>;
 template class QSimpson<2>;
@@ -1277,12 +1004,6 @@ template class QWeddle<2>;
 
 template class QGauss<3>;
 template class QGaussLobatto<3>;
-template class QGauss2<3>;
-template class QGauss3<3>;
-template class QGauss4<3>;
-template class QGauss5<3>;
-template class QGauss6<3>;
-template class QGauss7<3>;
 template class QMidpoint<3>;
 template class QTrapez<3>;
 template class QSimpson<3>;
index 0be5da66d4ccfbe85b510e45a8ad1f5e4dba243d..7d8bab0095e08a02077bd74f8c5896816427840d 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2003, 2005, 2006 by the deal.II authors
+//    Copyright (C) 2003, 2005, 2006, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -27,16 +27,7 @@ create_quadrature (const std::string &s,
   if(s == "gauss")
     {
       AssertThrow(order >= 2, ExcInvalidQGaussOrder(order));
-      switch(order)
-       {
-         case  2: return QGauss2<dim>();
-         case  3: return QGauss3<dim>();
-         case  4: return QGauss4<dim>();
-         case  5: return QGauss5<dim>();
-         case  6: return QGauss6<dim>();
-         case  7: return QGauss7<dim>();
-         default: return QGauss<dim>(order);
-       }
+      return QGauss<dim>(order);
     }
   else
     {
@@ -54,7 +45,7 @@ create_quadrature (const std::string &s,
                                   // return something to suppress
                                   // stupid warnings by some
                                   // compilers
-  return QGauss2<dim>();
+  return Quadrature<dim>();
 }
 
 
index ceadb44286638836d8c877ec6288bb00d9d00f87..0efd08050bfba3e7309d503a0ef82d6e927b3be2 100644 (file)
@@ -30,6 +30,12 @@ inconvenience this causes.
 </p>
 
 <ol>
+<li> The <code>QGauss2, QGauss3, ..., QGauss7</code> classes &mdash; deprecated
+for more than 6 years already &mdash; have finally been removed.
+You should use code like <code>QGauss@<dim@>(4)</code> instead.
+<br>
+(WB 2010/09/19)</li>
+
 
 <li> The FE_Nedelec class had previously implemented the lowest order
 when the value 1 was passed to the constructor. This has now been

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.