From: Wolfgang Bangerth Date: Sat, 18 Apr 2015 02:45:12 +0000 (-0500) Subject: Add a test. X-Git-Tag: v8.3.0-rc1~233^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=9b8fa53aaa8e2d3a6dc7684b08461281a6f81e1a;p=dealii.git Add a test. --- diff --git a/tests/hp/integrate_difference_03.cc b/tests/hp/integrate_difference_03.cc new file mode 100644 index 0000000000..1fe3e337e3 --- /dev/null +++ b/tests/hp/integrate_difference_03.cc @@ -0,0 +1,136 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2005 - 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. +// +// --------------------------------------------------------------------- + + + +// like _02, but test the Hdiv seminorm. this requires a vector-valued +// element and a function we interpolate that is vector-valued. we +// simply extend the function from _02 by zeros + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + + +template +double f (const Point &p) +{ + return p[0]; +} + + + +template +void test () +{ + deallog << "dim=" << dim << std::endl; + + Triangulation tria; + GridGenerator::hyper_cube(tria); + tria.refine_global (2); + tria.begin_active()->set_refine_flag (); + tria.execute_coarsening_and_refinement (); + tria.refine_global (4-dim); + + hp::FECollection fe_collection; + hp::QCollection q_collection; + for (unsigned int i=1; i<=4; ++i) + { + fe_collection.push_back(FESystem(FE_Q (i), dim)); + q_collection.push_back (QGauss (i+2)); + } + + + hp::DoFHandler dof_handler(tria); + + for (typename hp::DoFHandler::active_cell_iterator + cell = dof_handler.begin_active(); + cell != dof_handler.end(); ++cell) + cell->set_active_fe_index (Testing::rand() % fe_collection.size()); + + dof_handler.distribute_dofs(fe_collection); + + // interpolate a linear function + Vector vec (dof_handler.n_dofs()); + VectorTools::interpolate (dof_handler, + VectorFunctionFromScalarFunctionObject (&f, 0, dim), + vec); + + Vector diff (tria.n_active_cells()); + + // H1 seminorm. the function is u(x)=x, so + // its H1 seminorm should be equal to 1/2 + { + VectorTools::integrate_difference (dof_handler, + vec, + ZeroFunction(dim), + diff, + q_collection, + VectorTools::H1_seminorm); + deallog << "H1 seminorm, diff=" << diff.l2_norm() << std::endl; + } + + // Hdiv seminorm. the function is + // u(x)=x, so the norm must be equal to 1 + // on every cell + { + VectorTools::integrate_difference (dof_handler, + vec, + ZeroFunction(dim), + diff, + q_collection, + VectorTools::W1infty_seminorm); + deallog << "Hdiv semi, diff=" << diff.linfty_norm() << std::endl; + // also ensure that we indeed get the + // same value on every cell + diff. add(-1); + AssertThrow (diff.l2_norm() == 0, ExcInternalError()); + } +} + + +int main () +{ + std::ofstream logfile("output"); + logfile.precision(2); + deallog << std::setprecision(2); + + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + test<1> (); + test<2> (); + test<3> (); + + deallog << "OK" << std::endl; +} diff --git a/tests/hp/integrate_difference_03.output b/tests/hp/integrate_difference_03.output new file mode 100644 index 0000000000..1d59a57064 --- /dev/null +++ b/tests/hp/integrate_difference_03.output @@ -0,0 +1,11 @@ + +DEAL::dim=1 +DEAL::H1 seminorm, diff=1.0 +DEAL::Hdiv semi, diff=1.0 +DEAL::dim=2 +DEAL::H1 seminorm, diff=1.0 +DEAL::Hdiv semi, diff=1.0 +DEAL::dim=3 +DEAL::H1 seminorm, diff=1.0 +DEAL::Hdiv semi, diff=1.0 +DEAL::OK