From: Timo Heister Date: Wed, 20 Jan 2016 15:24:19 +0000 (-0500) Subject: add failing test X-Git-Tag: v8.4.0-rc2~66^2~9 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=fd542267b86ac31c10b2e10963bd37b6483ba4ab;p=dealii.git add failing test --- diff --git a/tests/codim_one/error_estimator_02.cc b/tests/codim_one/error_estimator_02.cc new file mode 100644 index 0000000000..665a550139 --- /dev/null +++ b/tests/codim_one/error_estimator_02.cc @@ -0,0 +1,160 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2004 - 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. +// +// --------------------------------------------------------------------- + + +// check that the kelly error estimator uses correct normals on surfaces +// with kinks by interpolating a function that is inside the FE space +// and should produce zero errors. + + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + + +template +class MyFunction : public Function +{ +public: + MyFunction () : Function(1) {} + + virtual double value (const Point &p, + const unsigned int component) const + { + return p(0)*p(0)+p(0)*p(1); + } + + virtual void vector_value (const Point &p, + Vector &values) const + { + values(0) = value(p,0); + values(1) = value(p,1); + } +}; + + + +template +Quadrature & +get_q_face () +{ + static QGauss q(4); + return q; +} + +template<> +Quadrature<0> & +get_q_face <1>() +{ + Quadrature<0> *q = 0; + return *q; +} + +template +void make_mesh (Triangulation &tria) +{ + // two faces of a hyper_cube + Triangulation volume_mesh; + GridGenerator::hyper_cube(volume_mesh, 0, 1, true); + std::set boundary_ids; + boundary_ids.insert (1); + boundary_ids.insert (2); + GridGenerator::extract_boundary_mesh (volume_mesh, tria, + boundary_ids); + tria.refine_global (2); +} + + + + +template +void +check () +{ + MyFunction function; + + Triangulation tria; + make_mesh (tria); + + FE_Q element(2); + DoFHandler dof(tria); + dof.distribute_dofs(element); + + MappingQ mapping(3); + Quadrature &q_face = get_q_face(); + + std::map*> neumann_bc; + neumann_bc[0] = &function; + + Vector v (dof.n_dofs()); + VectorTools::interpolate (mapping, dof, function, v); + + Vector error (tria.n_active_cells()); + + KellyErrorEstimator::estimate (mapping, dof, q_face, + typename FunctionMap::type(), + v, error); + deallog << "Estimated error indicators:" << std::endl; + for (unsigned int i=0; i > data_out; + data_out.attach_dof_handler (dof); + data_out.add_data_vector (v, + "solution", + DataOut >::type_dof_data); + data_out.add_data_vector (error, + "error"); + data_out.build_patches (); + std::string filename = spacedim == 2 ? + "solution-2d-" : + "solution-3d-"; + filename += Utilities::int_to_string(0,2) + ".vtk"; + std::ofstream output (filename.c_str()); + data_out.write_vtk (output); + + } + + deallog << "OK" << std::endl; +} + + +int main () +{ + initlog(); + + deallog.push ("1d"); + check<1,2> (); + deallog.pop(); + deallog.push ("2d"); + check<2,3> (); + deallog.pop(); +} diff --git a/tests/codim_one/error_estimator_02.output b/tests/codim_one/error_estimator_02.output new file mode 100644 index 0000000000..e69de29bb2