]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Testcase for new classes Polynomial and LagrangeEquidistant.
authorRalf Hartmann <Ralf.Hartmann@dlr.de>
Fri, 15 Dec 2000 07:02:41 +0000 (07:02 +0000)
committerRalf Hartmann <Ralf.Hartmann@dlr.de>
Fri, 15 Dec 2000 07:02:41 +0000 (07:02 +0000)
git-svn-id: https://svn.dealii.org/trunk@3535 0785d39b-7218-0410-832d-ea1e28bc413d

tests/base/Makefile.in
tests/base/polynomial_test.cc [new file with mode: 0644]
tests/base/polynomial_test.checked [new file with mode: 0644]

index 4de80563b9f87881596f937accac74284e037dfd..d0c0a946a85f9de0ab7292cc2b3ac53c8c5910e5 100644 (file)
@@ -137,6 +137,8 @@ quadrature_test.testcase: $(quadrature_test-o-files) $(libraries)
        @$(CXX) $(flags) -o $@ $^
 
 
+############################################################
+
 
 timer-cc-files = timer.cc
 
@@ -150,6 +152,23 @@ timer.testcase: $(timer-o-files) $(libraries)
        @echo =====linking======= $<
        @$(CXX) $(flags) -o $@ $^
 
+
+############################################################
+
+
+polynomial_test-cc-files = polynomial_test.cc
+
+ifeq ($(debug-mode),on)
+polynomial_test-o-files = $(polynomial_test-cc-files:.cc=.go)
+else
+polynomial_test-o-files = $(polynomial_test-cc-files:.cc=.o)
+endif
+
+polynomial_test.testcase: $(polynomial_test-o-files) $(libraries)
+       @echo =====linking======= $<
+       @$(CXX) $(flags) -o $@ $^
+
+
 ############################################################
 # Continue with other targets if needed
 ############################################################
diff --git a/tests/base/polynomial_test.cc b/tests/base/polynomial_test.cc
new file mode 100644 (file)
index 0000000..49f6559
--- /dev/null
@@ -0,0 +1,118 @@
+//----------------------------  polynomial_test.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  polynomial_test.cc  ---------------------------
+
+#include <iostream>
+#include <fstream>
+
+#include <base/logstream.h>
+#include <base/polynomial.h>
+
+
+
+bool equals_delta_ij(double value, unsigned int i, unsigned int j)
+{
+  double eps=1e-14;
+  if ((i==j && fabs(value-1)<eps) || (i!=j && fabs(value)<eps))
+    return true;
+  else
+    return false;
+}
+
+
+
+
+int main(int, char)
+{
+  ofstream logfile("polynomial_test.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  
+  vector<double> values(1);
+  deallog << "LagrangeEquidistant polynoms:" << endl;
+  for (unsigned int order=1; order<=4; ++order)
+    {
+      deallog << "Polynomial p of order " << order << endl;
+      for (unsigned int s_point=0; s_point<=order; ++s_point)
+       {
+         LagrangeEquidistant polynom(order, s_point);
+
+                                          // support points in vertices
+         if (order>0)
+           for (unsigned int i=0; i<=1; ++i)
+             {
+               double x=i;
+               polynom.value(x, values);
+               deallog << " p_" << s_point << "(" << x << ")";
+//             deallog << "=" << values[0];
+               if (equals_delta_ij(values[0], s_point, i))
+                 deallog << "   ok";
+               else
+                 deallog << "   false";
+               deallog << endl;
+             }
+                                          // support points on line
+         if (order>1)
+           for (unsigned int i=1; i<order; ++i)
+             {
+               double x=static_cast<double>(i)/order;
+               polynom.value(x, values);
+               deallog << " p_" << s_point << "(" << x << ")";
+//             deallog << "=" << values[0];
+               if (equals_delta_ij(values[0], s_point, i+1))
+                 deallog << "   ok";
+               else
+                 deallog << "   false";
+               deallog << endl;
+             }
+         deallog << endl;
+       }
+    }
+
+  deallog << endl << "Test derivatives computed by the Horner scheme:" << endl;
+  LagrangeEquidistant pol(4, 3);
+  vector<double> v_horner(6);
+  for (unsigned int i=0; i<=10; ++i)
+    {
+      double xi=i*0.1;
+      deallog << "x=" << xi << ",    all derivatives: ";
+      vector<double> v_exact(6);
+      
+      v_exact[0]=64.0*xi*xi*xi*xi-128.0*xi*xi*xi+76.0*xi*xi-12.0*xi;
+      v_exact[1]=256.0*xi*xi*xi-384.0*xi*xi+152.0*xi-12.0;
+      v_exact[2]=768.0*xi*xi-768.0*xi+152.0;
+      v_exact[3]=1536*xi-768;
+      v_exact[4]=1536;
+      v_exact[5]=0;
+
+      pol.value(xi, v_horner);
+
+      bool ok=true;
+      for (unsigned int i=0; i<v_exact.size(); ++i)
+       {
+//       deallog << "v_horner[i]=" << v_horner[i]
+//            << "   v_exact[i]=" << v_exact[i] << endl;
+         if (fabs(v_horner[i]-v_exact[i])>1e-12)
+           ok=false;
+       }
+
+      if (ok)
+       deallog << "ok";
+      else
+       deallog << "false";
+
+      deallog << endl;
+    }
+}
+
+
+
diff --git a/tests/base/polynomial_test.checked b/tests/base/polynomial_test.checked
new file mode 100644 (file)
index 0000000..06dbd84
--- /dev/null
@@ -0,0 +1,87 @@
+
+DEAL::LagrangeEquidistant polynoms:
+DEAL::Polynomial p of order 1
+DEAL:: p_0(0.00000)   ok
+DEAL:: p_0(1.00000)   ok
+
+DEAL:: p_1(0.00000)   ok
+DEAL:: p_1(1.00000)   ok
+
+DEAL::Polynomial p of order 2
+DEAL:: p_0(0.00000)   ok
+DEAL:: p_0(1.00000)   ok
+DEAL:: p_0(0.500000)   ok
+
+DEAL:: p_1(0.00000)   ok
+DEAL:: p_1(1.00000)   ok
+DEAL:: p_1(0.500000)   ok
+
+DEAL:: p_2(0.00000)   ok
+DEAL:: p_2(1.00000)   ok
+DEAL:: p_2(0.500000)   ok
+
+DEAL::Polynomial p of order 3
+DEAL:: p_0(0.00000)   ok
+DEAL:: p_0(1.00000)   ok
+DEAL:: p_0(0.333333)   ok
+DEAL:: p_0(0.666667)   ok
+
+DEAL:: p_1(0.00000)   ok
+DEAL:: p_1(1.00000)   ok
+DEAL:: p_1(0.333333)   ok
+DEAL:: p_1(0.666667)   ok
+
+DEAL:: p_2(0.00000)   ok
+DEAL:: p_2(1.00000)   ok
+DEAL:: p_2(0.333333)   ok
+DEAL:: p_2(0.666667)   ok
+
+DEAL:: p_3(0.00000)   ok
+DEAL:: p_3(1.00000)   ok
+DEAL:: p_3(0.333333)   ok
+DEAL:: p_3(0.666667)   ok
+
+DEAL::Polynomial p of order 4
+DEAL:: p_0(0.00000)   ok
+DEAL:: p_0(1.00000)   ok
+DEAL:: p_0(0.250000)   ok
+DEAL:: p_0(0.500000)   ok
+DEAL:: p_0(0.750000)   ok
+
+DEAL:: p_1(0.00000)   ok
+DEAL:: p_1(1.00000)   ok
+DEAL:: p_1(0.250000)   ok
+DEAL:: p_1(0.500000)   ok
+DEAL:: p_1(0.750000)   ok
+
+DEAL:: p_2(0.00000)   ok
+DEAL:: p_2(1.00000)   ok
+DEAL:: p_2(0.250000)   ok
+DEAL:: p_2(0.500000)   ok
+DEAL:: p_2(0.750000)   ok
+
+DEAL:: p_3(0.00000)   ok
+DEAL:: p_3(1.00000)   ok
+DEAL:: p_3(0.250000)   ok
+DEAL:: p_3(0.500000)   ok
+DEAL:: p_3(0.750000)   ok
+
+DEAL:: p_4(0.00000)   ok
+DEAL:: p_4(1.00000)   ok
+DEAL:: p_4(0.250000)   ok
+DEAL:: p_4(0.500000)   ok
+DEAL:: p_4(0.750000)   ok
+
+
+DEAL::Test derivatives computed by the Horner scheme:
+DEAL::x=0.00000,    all derivatives: ok
+DEAL::x=0.100000,    all derivatives: ok
+DEAL::x=0.200000,    all derivatives: ok
+DEAL::x=0.300000,    all derivatives: ok
+DEAL::x=0.400000,    all derivatives: ok
+DEAL::x=0.500000,    all derivatives: ok
+DEAL::x=0.600000,    all derivatives: ok
+DEAL::x=0.700000,    all derivatives: ok
+DEAL::x=0.800000,    all derivatives: ok
+DEAL::x=0.900000,    all derivatives: ok
+DEAL::x=1.00000,    all derivatives: ok

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.