From: Daniel Garcia-Sanchez Date: Wed, 9 Oct 2019 10:42:12 +0000 (+0200) Subject: Add test integrate_difference_01_complex_04 X-Git-Tag: v9.2.0-rc1~984^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=52e8c78dc4ba049b757aff2347a296d2bd330bf9;p=dealii.git Add test integrate_difference_01_complex_04 --- diff --git a/tests/vector_tools/integrate_difference_01_complex_04.cc b/tests/vector_tools/integrate_difference_01_complex_04.cc new file mode 100644 index 0000000000..bb750b3f10 --- /dev/null +++ b/tests/vector_tools/integrate_difference_01_complex_04.cc @@ -0,0 +1,116 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2003 - 2018 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + +// Test integrate_difference for complex-valued vectors. This test is +// like integrate_difference_01_complex_04, but compares the interpolated +// solution to the exact solution. + +#include +#include +#include + +#include + +#include +#include + +#include +#include +#include +#include + +#include + +#include + +#include "../tests.h" + +using namespace dealii; + + +// x+y(+z), x^2+y^2 (, z+xy) 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, Ref(), 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 << "L2_norm:" << std::endl; + // sqrt(\int_\Omega f^2) = sqrt(\int (x+y)^2+(x^2+y^2)^2) + test(VectorTools::L2_norm, 0.); + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char **argv) +{ + initlog(); + test<2>(); +} diff --git a/tests/vector_tools/integrate_difference_01_complex_04.with_complex_values=on.output b/tests/vector_tools/integrate_difference_01_complex_04.with_complex_values=on.output new file mode 100644 index 0000000000..496092bb8e --- /dev/null +++ b/tests/vector_tools/integrate_difference_01_complex_04.with_complex_values=on.output @@ -0,0 +1,4 @@ + +DEAL::L2_norm: +DEAL::computed: 2.68836e-16 expected: 0.00000 difference: 2.68836e-16 +DEAL::OK