From 0c17bd5fde6298c4165e8facb86c6855bfe3a62b Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Mon, 31 Aug 2015 13:54:38 -0400 Subject: [PATCH] add test for H1div_seminorm --- tests/vector_tools/CMakeLists.txt | 4 + tests/vector_tools/integrate_difference_01.cc | 109 ++++++++++++++++++ .../integrate_difference_01.output | 6 + 3 files changed, 119 insertions(+) create mode 100644 tests/vector_tools/CMakeLists.txt create mode 100644 tests/vector_tools/integrate_difference_01.cc create mode 100644 tests/vector_tools/integrate_difference_01.output diff --git a/tests/vector_tools/CMakeLists.txt b/tests/vector_tools/CMakeLists.txt new file mode 100644 index 0000000000..f8ac3eb029 --- /dev/null +++ b/tests/vector_tools/CMakeLists.txt @@ -0,0 +1,4 @@ +CMAKE_MINIMUM_REQUIRED(VERSION 2.8.9) +INCLUDE(../setup_testsubproject.cmake) +PROJECT(testsuite CXX) +DEAL_II_PICKUP_TESTS() diff --git a/tests/vector_tools/integrate_difference_01.cc b/tests/vector_tools/integrate_difference_01.cc new file mode 100644 index 0000000000..a94367eb3f --- /dev/null +++ b/tests/vector_tools/integrate_difference_01.cc @@ -0,0 +1,109 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2003 - 2014 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 + +#include "../tests.h" +#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) + {} + + double 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]; + } +}; + + + +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, + 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 << "OK" << std::endl; +} + + +int main (int argc, char **argv) +{ + initlog(); + test<2>(); + +} diff --git a/tests/vector_tools/integrate_difference_01.output b/tests/vector_tools/integrate_difference_01.output new file mode 100644 index 0000000000..001f3d030a --- /dev/null +++ b/tests/vector_tools/integrate_difference_01.output @@ -0,0 +1,6 @@ + +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::OK -- 2.39.5