From 3b576589cff7325a341392e9aca8acd53b901994 Mon Sep 17 00:00:00 2001
From: Martin Kronbichler <kronbichler@lnm.mw.tum.de>
Date: Mon, 22 Jul 2013 06:03:33 +0000
Subject: [PATCH] Simplify test case

git-svn-id: https://svn.dealii.org/trunk@30087 0785d39b-7218-0410-832d-ea1e28bc413d
---
 tests/base/gauss_lobatto.cc | 175 +-----------------------------------
 1 file changed, 3 insertions(+), 172 deletions(-)

diff --git a/tests/base/gauss_lobatto.cc b/tests/base/gauss_lobatto.cc
index 69c56353da..466d3dfa36 100644
--- a/tests/base/gauss_lobatto.cc
+++ b/tests/base/gauss_lobatto.cc
@@ -15,11 +15,7 @@
 // ---------------------------------------------------------------------
 
 
-
-// check accuracy of various quadrature formulas by using them to
-// integrate polynomials of increasing degree, and finding the degree
-// until which they integrate exactly
-
+// check points and weights for Gauss-Lobatto quadrature formula
 
 #include "../tests.h"
 #include <iomanip>
@@ -27,170 +23,8 @@
 
 #include <deal.II/base/logstream.h>
 #include <deal.II/base/quadrature_lib.h>
-#include <deal.II/base/qprojector.h>
 #include <cmath>
 
-template <int dim>
-void
-fill_vector (std::vector<Quadrature<dim> *> &quadratures)
-{
-  quadratures.push_back (new QMidpoint<dim>());
-  quadratures.push_back (new QTrapez<dim>());
-  quadratures.push_back (new QSimpson<dim>());
-  quadratures.push_back (new QMilne<dim>());
-  quadratures.push_back (new QWeddle<dim>());
-  for (unsigned int i=0; i<9; ++i)
-    {
-      quadratures.push_back (new QGauss<dim>(i));
-    }
-  QMilne<1> q1d;
-  quadratures.push_back (new Quadrature<dim>(q1d));
-  for (unsigned int i=2; i<8; ++i)
-    {
-      quadratures.push_back (new QGaussLobatto<dim>(i));
-    }
-}
-
-template <int dim>
-void
-check_cells (std::vector<Quadrature<dim>*> &quadratures)
-{
-  Quadrature<dim> quadrature;
-  for (unsigned int n=0; n<quadratures.size(); ++n)
-    {
-      quadrature = *quadratures[n];
-      const std::vector<Point<dim> > &points=quadrature.get_points();
-      const std::vector<double> &weights=quadrature.get_weights();
-
-      deallog << "Quadrature no." << n;
-
-      unsigned int i=0;
-      double quadrature_int=0;
-      double exact_int=0;
-      double err = 0;
-
-      do
-        {
-          ++i;
-
-          quadrature_int=0;
-          // Check the polynomial x^i*y^i
-
-          for (unsigned int x=0; x<quadrature.size(); ++x)
-            {
-              double f=1.;
-              switch (dim)
-                {
-                case 3:
-                  f *= std::pow(static_cast<double>(points[x](2)), i*1.0);
-                case 2:
-                  f *= std::pow(static_cast<double>(points[x](1)), i*1.0);
-                case 1:
-                  f *= std::pow(static_cast<double>(points[x](0)), i*1.0);
-                }
-              quadrature_int+=f*weights[x];
-            }
-
-          // the exact integral is 1/(i+1)
-          exact_int=1./std::pow(static_cast<double>(i+1),dim);
-          err = std::fabs(quadrature_int-exact_int);
-        }
-      while (err<1e-14);
-      // Uncomment here for testing
-//      deallog << " (Int " << quadrature_int << ',' << exact_int << ")";
-      deallog << " is exact for polynomials of degree " << i-1 << std::endl;
-
-      if (dim==1)
-        {
-          // check the ordering of
-          // the quadrature points
-          bool in_order=true;
-          for (unsigned int x=1; x<quadrature.size(); ++x)
-            {
-              if (points[x](0)<=points[x-1](0))
-                in_order=false;
-            }
-          if (!in_order)
-            for (unsigned int x=0; x<quadrature.size(); ++x)
-              deallog << points[x] << std::endl;
-        }
-    }
-}
-
-
-template <int dim>
-void
-check_faces (const std::vector<Quadrature<dim-1>*>& quadratures, const bool sub)
-{
-  if (sub)
-    deallog.push("subfaces");
-  else
-    deallog.push("faces");
-
-  for (unsigned int n=0; n<quadratures.size(); ++n)
-    {
-      Quadrature<dim> quadrature (sub == false?
-                                  QProjector<dim>::project_to_all_faces(*quadratures[n]) :
-                                  QProjector<dim>::project_to_all_subfaces(*quadratures[n]));
-      const std::vector<Point<dim> > &points=quadrature.get_points();
-      const std::vector<double> &weights=quadrature.get_weights();
-
-      deallog << "Quadrature no." << n;
-
-      unsigned int i=0;
-      long double quadrature_int=0;
-      double exact_int=0;
-      double err = 0;
-
-      do
-        {
-          ++i;
-
-          quadrature_int=0;
-          // Check the polynomial
-          // x^i*y^i*z^i
-
-          for (unsigned int x=0; x<quadrature.size(); ++x)
-            {
-              long double f=1.;
-              switch (dim)
-                {
-                case 3:
-                  f *= std::pow((long double) points[x](2), i*1.0L);
-                case 2:
-                  f *= std::pow((long double) points[x](1), i*1.0L);
-                case 1:
-                  f *= std::pow((long double) points[x](0), i*1.0L);
-                }
-              quadrature_int+=f*weights[x];
-            }
-
-          // the exact integral is
-          // 1/(i+1)^(dim-1)
-          switch (dim)
-            {
-            case 2:
-              exact_int = 2 * (sub ? 2:1) / (double) (i+1);
-              break;
-            case 3:
-              exact_int = 3 * (sub ? (4+2+2):1)*8 / (double) (i+1)/(i+1);
-              break;
-            }
-
-          err = std::fabs(quadrature_int-exact_int);
-        }
-      // for comparison: use factor 8 in case
-      // of dim==3, as we integrate 8 times
-      // over the whole surface (all
-      // combinations of face_orientation,
-      // face_flip and face_rotation)
-      while (err < (dim==3 ? 8 : 1) * 2e-14);
-      // Uncomment here for testing
-      //      deallog << " (Int " << quadrature_int << '-' << exact_int << '=' << err << ")";
-      deallog << " is exact for polynomials of degree " << i-1 << std::endl;
-    }
-  deallog.pop();
-}
 
 int main()
 {
@@ -207,16 +41,13 @@ int main()
       for (unsigned int i=0; i<q.size(); ++i)
         deallog << q.point(i) << ' ' << q.weight(i) << std::endl;
 
-      // the points must be
-      // symmetrically located around
-      // 0.5
+      // the points must be symmetrically located around 0.5
       double p = 0;
       for (unsigned int i=0; i<q.size(); ++i)
         p += (q.point(i)[0] - 0.5);
       Assert (std::fabs(p) < 1e-12, ExcInternalError());
 
-      // the sum of weights must be
-      // one
+      // the sum of weights must be one
       double w = 0;
       for (unsigned int i=0; i<q.size(); ++i)
         w += q.weight(i);
-- 
2.39.5