From 362ad234c1fe52b6d911045613b0b108d0218b6e Mon Sep 17 00:00:00 2001 From: wolf Date: Mon, 21 Apr 2003 16:21:57 +0000 Subject: [PATCH] Add some new tests. git-svn-id: https://svn.dealii.org/trunk@7423 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/base/Makefile | 7 +- tests/base/anisotropic_1.cc | 174 ++++++++++++++++++++++++++++++++++++ tests/base/anisotropic_2.cc | 166 ++++++++++++++++++++++++++++++++++ tests/base/hierarchical.cc | 48 ++++++++++ 4 files changed, 393 insertions(+), 2 deletions(-) create mode 100644 tests/base/anisotropic_1.cc create mode 100644 tests/base/anisotropic_2.cc create mode 100644 tests/base/hierarchical.cc diff --git a/tests/base/Makefile b/tests/base/Makefile index 2747de2d9a..dc044d7627 100644 --- a/tests/base/Makefile +++ b/tests/base/Makefile @@ -34,12 +34,15 @@ tensor.exe : tensor.g.$(OBJEXT) $(libraries) timer.exe : timer.g.$(OBJEXT) $(libraries) threads.exe : threads.g.$(OBJEXT) $(libraries) auto_derivative_function.exe : auto_derivative_function.g.$(OBJEXT) $(libraries) - +anisotropic_1.exe : anisotropic_1.g.$(OBJEXT) $(lib-base.g) +anisotropic_2.exe : anisotropic_2.g.$(OBJEXT) $(lib-base.g) +hierarchical.exe : hierarchical.g.$(OBJEXT) $(lib-base.g) tests = $(sort \ logtest reference quadrature_test table tensor \ timer threads polynomial1d polynomial_test \ - auto_derivative_function) + auto_derivative_function anisotropic_1 anisotropic_2 \ + hierarchical) ############################################################ diff --git a/tests/base/anisotropic_1.cc b/tests/base/anisotropic_1.cc new file mode 100644 index 0000000000..9c3f727f2f --- /dev/null +++ b/tests/base/anisotropic_1.cc @@ -0,0 +1,174 @@ +//---------------------------- anisotropic_1.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003 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. +// +//---------------------------- anisotropic_1.cc --------------------------- + + +// check that TensorProducPolynomials and AnisotropicPolynomials +// compute the same thing if the basis polynomials for the anisotropic +// class happen to be the same for all coordinate directions + + +#include +#include + +#include +#include + + +using namespace Polynomials; + + +template +void check_poly(const Point &x, + const POLY1 &p, + const POLY2 &q) +{ + const unsigned int n = p.n(); + + std::vector values1 (n), values2 (n); + std::vector > gradients1(n), gradients2(n); + std::vector > second1(n), second2(n); + + p.compute (x, values1, gradients1, second1); + q.compute (x, values2, gradients2, second2); + + for (unsigned int k=0; k 5.0E-16) + deallog << 'P' << k << ": values differ " << val1 << " != " + << values1[k] << std::endl; + double val2 = q.compute_value(k,x); + if (fabs(val2 - values2[k]) > 5.0E-16) + deallog << 'Q' << k << ": values differ " << val2 << " != " + << values2[k] << std::endl; + if (fabs(val2 - val1) > 5.0E-16) + deallog << "PQ" << k << ": values differ " << val1 << " != " + << val2 << std::endl; + + // Check if compute_grad is ok + Tensor<1,dim> grad1 = p.compute_grad(k,x); + if (grad1 != gradients1[k]) + deallog << 'P' << k << ": gradients differ " << grad1 << " != " + << gradients1[k] << std::endl; + Tensor<1,dim> grad2 = q.compute_grad(k,x); + if (grad2 != gradients2[k]) + deallog << 'Q' << k << ": gradients differ " << grad1 << " != " + << gradients2[k] << std::endl; + if (grad2 != grad1) + deallog << "PQ" << k << ": gradients differ " << grad1 << " != " + << grad2 << std::endl; + + // Check if compute_grad_grad is ok + Tensor<2,dim> grad_grad1 = p.compute_grad_grad(k,x); + if (grad_grad1 != second1[k]) + deallog << 'P' << k << ": second derivatives differ " << grad_grad1 << " != " + << second1[k] << std::endl; + Tensor<2,dim> grad_grad2 = q.compute_grad_grad(k,x); + if (grad_grad2 != second2[k]) + deallog << 'Q' << k << ": second derivatives differ " << grad_grad2 << " != " + << second2[k] << std::endl; + if (grad_grad2 != grad_grad1) + deallog << "PQ" << k << ": second derivatives differ " << grad_grad1 << " != " + << grad_grad2 << std::endl; + + + // finally output values, + // gradients, etc, to make sure + // that they are not only + // consistent, but also + // correct. Multiply them + // somewhat to make them + // significant despite our + // two-post-dot-digits limit + values1[k] *= pow(10, dim); + gradients1[k] *= pow(10, dim); + + deallog << 'P' << k << "\t= " << values1[k] + << "\tgradient\t"; + for (unsigned int d=0;d +void +check_tensor (const std::vector >& v, + const Point& x) +{ + TensorProductPolynomials p(v); + + std::vector > > pols (dim, v); + AnisotropicPolynomials q(pols); + + check_poly (x, p, q); +} + + +void +check_dimensions (const std::vector >& p) +{ + deallog.push("1d"); + check_tensor(p, Point<1>(.5)); + deallog.pop(); + + deallog.push("2d"); + check_tensor(p, Point<2>(.5, .2)); + deallog.pop(); + + deallog.push("3d"); + check_tensor(p, Point<3>(.5, .2, .3)); + deallog.pop(); +} + +int main() +{ + std::ofstream logfile("anisotropic_1.output"); + logfile.precision(2); + deallog.attach(logfile); + deallog.depth_console(0); + + deallog.push("Lagrange"); + std::vector > p; + for (unsigned int i=0;i<3;++i) + p.push_back (LagrangeEquidistant(3, i)); + + check_dimensions(p); + + deallog.pop(); + deallog.push("Legendre"); + + p.clear (); + for (unsigned int i=0;i<3;++i) + p.push_back (Legendre(i)); + + check_dimensions(p); + + deallog.pop(); + deallog.push("Hierarchical"); + + p.clear (); + for (unsigned int i=0;i<3;++i) + p.push_back (Hierarchical(i)); + + check_dimensions(p); + + deallog.pop(); +} diff --git a/tests/base/anisotropic_2.cc b/tests/base/anisotropic_2.cc new file mode 100644 index 0000000000..6fd4d51846 --- /dev/null +++ b/tests/base/anisotropic_2.cc @@ -0,0 +1,166 @@ +//---------------------------- anisotropic_2.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003 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. +// +//---------------------------- anisotropic_2.cc --------------------------- + + +// check AnisotropicPolynomials + + +#include +#include + +#include +#include + + +using namespace Polynomials; + +typedef std::vector > PolVector; + + +void print_2d (const AnisotropicPolynomials<2> &aniso) +{ + const unsigned int N=10, M=13; + for (unsigned int i=0; i<=N; ++i) + { + deallog << std::endl; + for (unsigned int j=0; j<=M; ++j) + { + deallog << 1.*i/N << " " + << 1.*j/M << " "; + for (unsigned int k=0; k(1.*i/N, 1.*j/M)) + << " "; + deallog << std::endl; + for (unsigned int k=0; k(1.*i/N, 1.*j/M)) + << " "; + deallog << std::endl; + } + } +} + + + +template +void check_2d () +{ + // two checks with higher degree in + // x or y direction + { + PolVector pols[2] = { Pol::generate_complete_basis (3), + Pol::generate_complete_basis (1) }; + AnisotropicPolynomials<2> aniso (std::vector(&pols[0], + &pols[2])); + print_2d (aniso); + } + { + PolVector pols[2] = { Pol::generate_complete_basis (2), + Pol::generate_complete_basis (3) }; + AnisotropicPolynomials<2> aniso (std::vector(&pols[0], + &pols[2])); + + print_2d (aniso); + } +} + + +void print_3d (const AnisotropicPolynomials<3> &aniso) +{ + const unsigned int N=4, M=3, P=5; + for (unsigned int i=0; i<=N; ++i) + { + deallog << std::endl; + for (unsigned int j=0; j<=M; ++j) + { + deallog << std::endl; + for (unsigned int k=0; k<=P; ++k) + { + deallog << 1.*i/N << " " + << 1.*j/M << " " + << 1.*k/P; + for (unsigned int k=0; k(1.*i/N, 1.*j/M, 1.*k/P)) + << " "; + deallog << std::endl; + for (unsigned int k=0; k(1.*i/N, 1.*j/M, 1.*k/P)) + << " "; + deallog << std::endl; + } + } + } +} + + + +template +void check_3d () +{ + // three checks with higher degree + // in x, y or z direction + { + PolVector pols[3] = { Pol::generate_complete_basis (3), + Pol::generate_complete_basis (1), + Pol::generate_complete_basis (1) }; + AnisotropicPolynomials<3> aniso (std::vector(&pols[0], + &pols[3])); + print_3d (aniso); + } + { + PolVector pols[3] = { Pol::generate_complete_basis (1), + Pol::generate_complete_basis (3), + Pol::generate_complete_basis (1) }; + AnisotropicPolynomials<3> aniso (std::vector(&pols[0], + &pols[3])); + print_3d (aniso); + } + { + PolVector pols[3] = { Pol::generate_complete_basis (1), + Pol::generate_complete_basis (2), + Pol::generate_complete_basis (3) }; + AnisotropicPolynomials<3> aniso (std::vector(&pols[0], + &pols[3])); + print_3d (aniso); + } +} + + + +template +void check () +{ + check_2d (); + check_3d (); +} + + + +int main() +{ + std::ofstream logfile("anisotropic_2.output"); + logfile.precision(2); + deallog.attach(logfile); + deallog.depth_console(0); + +// deallog.push("Lagrange"); +// check (); +// deallog.pop(); + +// deallog.push("Legendre"); +// check > (); +// deallog.pop(); + + deallog.push("Hierarchical"); + check > (); + deallog.pop(); +} diff --git a/tests/base/hierarchical.cc b/tests/base/hierarchical.cc new file mode 100644 index 0000000000..fdab9246eb --- /dev/null +++ b/tests/base/hierarchical.cc @@ -0,0 +1,48 @@ +//----------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003 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. +// +//----------------------------------------------------------------------------- + +// generate a hierarchical basis and display the values of the shape +// functions at equidistant points. (I needed this output at one point +// in time, so why not make it a testcase -- WB) + +#include +#include + +#include +#include + + +using namespace Polynomials; + + +int main () +{ + std::ofstream logfile("hierarchical.output"); + logfile.precision(3); + deallog.attach(logfile); + deallog.depth_console(0); + + const std::vector > + p = Hierarchical::generate_complete_basis (10); + + const unsigned int div=30; + for (unsigned int i=0; i<=div; ++i) + { + const double x = 1.*i/div; + deallog << x << " "; + for (unsigned int j=0; j