From: Wolfgang Bangerth Date: Thu, 9 Nov 2017 14:52:31 +0000 (-0700) Subject: Add tests for VectorTools::integrate_difference() using complex data types. X-Git-Tag: v9.0.0-rc1~795^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ca6da4aa963f3f47f131a20a2680ff850bf096d0;p=dealii.git Add tests for VectorTools::integrate_difference() using complex data types. --- diff --git a/tests/vector_tools/integrate_difference_01_complex_01.cc b/tests/vector_tools/integrate_difference_01_complex_01.cc new file mode 100644 index 0000000000..ecaf6feb26 --- /dev/null +++ b/tests/vector_tools/integrate_difference_01_complex_01.cc @@ -0,0 +1,116 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2003 - 2015 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// Test integrate_difference for complex-valued vectors. This test is +// like integrate_difference_01, using a real-valued function but +// using complex data types + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +using namespace dealii; + + +// x+y(+z), x^2+y^2 (, z+xy) +// div = 1+2y (+1) +template +class Ref : public Function > +{ +public: + Ref() + :Function >(dim) + {} + + std::complex value (const Point &p, const unsigned int c) const + { + if (c==0) + return p[0]+p[1]+((dim==3)?p[2]:0.0); + if (c==1) + return p[0]*p[0]+p[1]*p[1]; + if (c==2) + return p[2]+p[0]*p[1]; + else + return numbers::signaling_nan(); + } +}; + + + +template +void test(VectorTools::NormType norm, double value) +{ + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(1); + + FESystem fe(FE_Q(4),dim); + DoFHandler dofh(tria); + dofh.distribute_dofs(fe); + + Vector > solution (dofh.n_dofs ()); + VectorTools::interpolate(dofh, Ref(), solution); + + Vector cellwise_errors (tria.n_active_cells()); + VectorTools::integrate_difference (dofh, + solution, + Functions::ZeroFunction(dim), + cellwise_errors, + QGauss(5), + norm); + + const double error = cellwise_errors.l2_norm(); + + const double difference = std::abs(error-value); + deallog << "computed: " << error + << " expected: " << value + << " difference: " << difference + << std::endl; + Assert(difference<1e-10, ExcMessage("Error in integrate_difference")); +} + +template +void test() +{ + deallog << "Hdiv_seminorm:" << std::endl; + // sqrt(\int (div f)^2 = sqrt(\int (1+2y)^2) + test(VectorTools::Hdiv_seminorm, std::sqrt(13.0/3.0)); + deallog << "L2_norm:" << std::endl; + // sqrt(\int_\Omega f^2) = sqrt(\int (x+y)^2+(x^2+y^2)^2) + test(VectorTools::L2_norm, std::sqrt(161.0/90.0)); + deallog << "H1_seminorm:" << std::endl; + // sqrt( sum | d/dxi f |_0^2 ) = sqrt( sum \int ) + test(VectorTools::H1_seminorm, std::sqrt(14.0/3.0)); + + deallog << "OK" << std::endl; +} + + +int main (int argc, char **argv) +{ + initlog(); + test<2>(); +} diff --git a/tests/vector_tools/integrate_difference_01_complex_01.output b/tests/vector_tools/integrate_difference_01_complex_01.output new file mode 100644 index 0000000000..9706d0feec --- /dev/null +++ b/tests/vector_tools/integrate_difference_01_complex_01.output @@ -0,0 +1,8 @@ + +DEAL::Hdiv_seminorm: +DEAL::computed: 2.08167 expected: 2.08167 difference: 4.44089e-16 +DEAL::L2_norm: +DEAL::computed: 1.33749 expected: 1.33749 difference: 2.22045e-16 +DEAL::H1_seminorm: +DEAL::computed: 2.16025 expected: 2.16025 difference: 4.44089e-16 +DEAL::OK diff --git a/tests/vector_tools/integrate_difference_01_complex_02.cc b/tests/vector_tools/integrate_difference_01_complex_02.cc new file mode 100644 index 0000000000..7b20a50926 --- /dev/null +++ b/tests/vector_tools/integrate_difference_01_complex_02.cc @@ -0,0 +1,116 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2003 - 2015 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// Test integrate_difference for complex-valued vectors. This test is +// like integrate_difference_01, but multiplies the real-valued +// function by sqrt(-1) and stores everything using complex data types + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +using namespace dealii; + + +// x+y(+z), x^2+y^2 (, z+xy) times i +// div = (1+2y (+1))i +template +class Ref : public Function > +{ +public: + Ref() + :Function >(dim) + {} + + std::complex value (const Point &p, const unsigned int c) const + { + if (c==0) + return std::complex(0, p[0]+p[1]+((dim==3)?p[2]:0.0)); + if (c==1) + return std::complex(0, p[0]*p[0]+p[1]*p[1]); + if (c==2) + return std::complex(0, p[2]+p[0]*p[1]); + else + return numbers::signaling_nan(); + } +}; + + + +template +void test(VectorTools::NormType norm, double value) +{ + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(1); + + FESystem fe(FE_Q(4),dim); + DoFHandler dofh(tria); + dofh.distribute_dofs(fe); + + Vector > solution (dofh.n_dofs ()); + VectorTools::interpolate(dofh, Ref(), solution); + + Vector cellwise_errors (tria.n_active_cells()); + VectorTools::integrate_difference (dofh, + solution, + Functions::ZeroFunction(dim), + cellwise_errors, + QGauss(5), + norm); + + const double error = cellwise_errors.l2_norm(); + + const double difference = std::abs(error-value); + deallog << "computed: " << error + << " expected: " << value + << " difference: " << difference + << std::endl; + Assert(difference<1e-10, ExcMessage("Error in integrate_difference")); +} + +template +void test() +{ + deallog << "Hdiv_seminorm:" << std::endl; + // sqrt(\int (div f)^2 = sqrt(\int (1+2y)^2) + test(VectorTools::Hdiv_seminorm, std::sqrt(13.0/3.0)); + deallog << "L2_norm:" << std::endl; + // sqrt(\int_\Omega f^2) = sqrt(\int (x+y)^2+(x^2+y^2)^2) + test(VectorTools::L2_norm, std::sqrt(161.0/90.0)); + deallog << "H1_seminorm:" << std::endl; + // sqrt( sum | d/dxi f |_0^2 ) = sqrt( sum \int ) + test(VectorTools::H1_seminorm, std::sqrt(14.0/3.0)); + + deallog << "OK" << std::endl; +} + + +int main (int argc, char **argv) +{ + initlog(); + test<2>(); +} diff --git a/tests/vector_tools/integrate_difference_01_complex_02.output b/tests/vector_tools/integrate_difference_01_complex_02.output new file mode 100644 index 0000000000..9706d0feec --- /dev/null +++ b/tests/vector_tools/integrate_difference_01_complex_02.output @@ -0,0 +1,8 @@ + +DEAL::Hdiv_seminorm: +DEAL::computed: 2.08167 expected: 2.08167 difference: 4.44089e-16 +DEAL::L2_norm: +DEAL::computed: 1.33749 expected: 1.33749 difference: 2.22045e-16 +DEAL::H1_seminorm: +DEAL::computed: 2.16025 expected: 2.16025 difference: 4.44089e-16 +DEAL::OK diff --git a/tests/vector_tools/integrate_difference_01_complex_03.cc b/tests/vector_tools/integrate_difference_01_complex_03.cc new file mode 100644 index 0000000000..7662c67d80 --- /dev/null +++ b/tests/vector_tools/integrate_difference_01_complex_03.cc @@ -0,0 +1,116 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2003 - 2015 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// Test integrate_difference for complex-valued vectors. This test is +// like integrate_difference_01, but multiplies the real-valued +// function by (1/sqrt(2)+1/sqrt(2)i) and stores everything using complex data types + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +using namespace dealii; + + +// x+y(+z), x^2+y^2 (, z+xy) times (1+i)/sqrt(2) +// div = (1+2y (+1)) times (1+i)/sqrt(2) +template +class Ref : public Function > +{ +public: + Ref() + :Function >(dim) + {} + + std::complex value (const Point &p, const unsigned int c) const + { + if (c==0) + return (p[0]+p[1]+((dim==3)?p[2]:0.0)) * std::complex(1./std::sqrt(2.), 1./std::sqrt(2.)); + if (c==1) + return (p[0]*p[0]+p[1]*p[1]) * std::complex(1./std::sqrt(2.), 1./std::sqrt(2.)); + if (c==2) + return (p[2]+p[0]*p[1]) * std::complex(1./std::sqrt(2.), 1./std::sqrt(2.)); + else + return numbers::signaling_nan(); + } +}; + + + +template +void test(VectorTools::NormType norm, double value) +{ + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global(1); + + FESystem fe(FE_Q(4),dim); + DoFHandler dofh(tria); + dofh.distribute_dofs(fe); + + Vector > solution (dofh.n_dofs ()); + VectorTools::interpolate(dofh, Ref(), solution); + + Vector cellwise_errors (tria.n_active_cells()); + VectorTools::integrate_difference (dofh, + solution, + Functions::ZeroFunction(dim), + cellwise_errors, + QGauss(5), + norm); + + const double error = cellwise_errors.l2_norm(); + + const double difference = std::abs(error-value); + deallog << "computed: " << error + << " expected: " << value + << " difference: " << difference + << std::endl; + Assert(difference<1e-10, ExcMessage("Error in integrate_difference")); +} + +template +void test() +{ + deallog << "Hdiv_seminorm:" << std::endl; + // sqrt(\int (div f)^2 = sqrt(\int (1+2y)^2) + test(VectorTools::Hdiv_seminorm, std::sqrt(13.0/3.0)); + deallog << "L2_norm:" << std::endl; + // sqrt(\int_\Omega f^2) = sqrt(\int (x+y)^2+(x^2+y^2)^2) + test(VectorTools::L2_norm, std::sqrt(161.0/90.0)); + deallog << "H1_seminorm:" << std::endl; + // sqrt( sum | d/dxi f |_0^2 ) = sqrt( sum \int ) + test(VectorTools::H1_seminorm, std::sqrt(14.0/3.0)); + + deallog << "OK" << std::endl; +} + + +int main (int argc, char **argv) +{ + initlog(); + test<2>(); +} diff --git a/tests/vector_tools/integrate_difference_01_complex_03.output b/tests/vector_tools/integrate_difference_01_complex_03.output new file mode 100644 index 0000000000..9706d0feec --- /dev/null +++ b/tests/vector_tools/integrate_difference_01_complex_03.output @@ -0,0 +1,8 @@ + +DEAL::Hdiv_seminorm: +DEAL::computed: 2.08167 expected: 2.08167 difference: 4.44089e-16 +DEAL::L2_norm: +DEAL::computed: 1.33749 expected: 1.33749 difference: 2.22045e-16 +DEAL::H1_seminorm: +DEAL::computed: 2.16025 expected: 2.16025 difference: 4.44089e-16 +DEAL::OK