From 5775dc74d7771613614ce03e648fccac439f7cdb Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 16 Jul 2013 19:45:42 +0000 Subject: [PATCH] Fix a bug in VectorTools::project. git-svn-id: https://svn.dealii.org/trunk@30020 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 7 + .../deal.II/numerics/vector_tools.templates.h | 2 +- tests/deal.II/project_01_curved_boundary.cc | 171 ++++++++++++++++++ .../project_01_curved_boundary/cmp/generic | 165 +++++++++++++++++ 4 files changed, 344 insertions(+), 1 deletion(-) create mode 100644 tests/deal.II/project_01_curved_boundary.cc create mode 100644 tests/deal.II/project_01_curved_boundary/cmp/generic diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index d8442a4bc5..fae3f4bc26 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -164,6 +164,13 @@ this function.
    +
  1. Fixed: VectorTools::project has an option to first project onto the +boundary. However, the implementation of this option ignored the mapping +that is provided to the function. This is now fixed. +
    +(Wolfgang Bangerth, 2013/07/16) +
  2. +
  3. Improved: The WorkStream class used throughout deal.II is now using thread local variables and initializes temporary variables on the thread that uses them, leading to better cache locality. diff --git a/deal.II/include/deal.II/numerics/vector_tools.templates.h b/deal.II/include/deal.II/numerics/vector_tools.templates.h index 6e90a698fa..89452ceffc 100644 --- a/deal.II/include/deal.II/numerics/vector_tools.templates.h +++ b/deal.II/include/deal.II/numerics/vector_tools.templates.h @@ -575,7 +575,7 @@ namespace VectorTools typename FunctionMap::type boundary_functions; for (types::boundary_id c=0; c arguments + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + + +// define a function that is zero within the square +// [-sqrt(2)/2,sqrt(2)/2]^d +// which is the domain produced by GridGenerator::hyper_sphere +// when using a linear mapping. we will want to see a solution +// that is not zero when using a higher order mapping and a Q2 +// element +template +class F : public Function +{ + public: + virtual double value (const Point &p, + const unsigned int = 0) const + { + // compute the linfty norm of p + double m = 0; + for (unsigned int d=0; d +void test() +{ + HyperBallBoundary boundary; + + // create 2 triangulations with the + // same coarse grid, and refine + // them differently + Triangulation tria; + GridGenerator::hyper_ball (tria); + tria.set_boundary (0, boundary); + + + FE_Q fe(1); + DoFHandler dh (tria); + dh.distribute_dofs (fe); + + Vector v(dh.n_dofs()); + + ConstraintMatrix cm; + cm.close (); + + // use the implicit Q1 mapping. this will yield a zero solution + { + VectorTools::project (dh, cm, QGauss(3), F(), + v); + deallog << v.l2_norm() << std::endl; + Assert (v.l2_norm() == 0, ExcInternalError()); + } + + // use an explicit Q1 mapping. this will yield a zero solution + { + VectorTools::project (MappingQ1(), + dh, cm, QGauss(3), F(), + v); + deallog << v.l2_norm() << std::endl; + Assert (v.l2_norm() == 0, ExcInternalError()); + } + + // use an explicit Q2 mapping. this will yield a nonzero solution if it's a + // straight projection since some of the quadrature points are lying outside + // the area where the function is zero + { + VectorTools::project (MappingQ(2), + dh, cm, QGauss(3), F(), + v); + deallog << v.l2_norm() << std::endl; + Assert (v.l2_norm() != 0, ExcInternalError()); + } + + // use an explicit Q2 mapping but enforce zero boundary values. this will + // yield a nonzero solution since some of the quadrature points are lying + // outside the area where the function is zero, even though the values of + // the DoFs at the boundary are in fact zero (they are interpolated only at + // points where the function is zero) + { + VectorTools::project (MappingQ(2), + dh, cm, QGauss(3), F(), + v, true); + deallog << v.l2_norm() << std::endl; + Assert (v.l2_norm() != 0, ExcInternalError()); + } + + // use an explicit Q2 mapping and project to the boundary first. this will + // yield a nonzero solution since some of the quadrature points are lying + // outside the area where the function is zero. furthermore, the values on + // the boundary will be nonzero since we project, even though the + // *interpolation* of boundary values onto the trace of the Q1 space is zero + { + VectorTools::project (MappingQ(2), + dh, cm, QGauss(3), F(), + v, false, + QGauss(2), true); + deallog << v.l2_norm() << std::endl; + Assert (v.l2_norm() != 0, ExcInternalError()); + for (typename DoFHandler::active_cell_iterator cell=dh.begin_active(); + cell != dh.end(); ++cell) + for (unsigned int i=0; i::vertices_per_cell; ++i) + deallog << cell->vertex(i) << ' ' << v(cell->vertex_dof_index(i,0)) + << std::endl; + } + + + // same as above, but use a projection with a QTrapez formula. this happens + // to evaluate the function only at points where it is zero, and + // consequently the values at the boundary should be zero + { + VectorTools::project (MappingQ(2), + dh, cm, QGauss(3), F(), + v, false, + QTrapez(), true); + deallog << v.l2_norm() << std::endl; + Assert (v.l2_norm() != 0, ExcInternalError()); + for (typename DoFHandler::active_cell_iterator cell=dh.begin_active(); + cell != dh.end(); ++cell) + for (unsigned int i=0; i::vertices_per_cell; ++i) + deallog << cell->vertex(i) << ' ' << v(cell->vertex_dof_index(i,0)) + << std::endl; + } +} + + +int main() +{ + std::ofstream logfile ("project_01_curved_boundary/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + test<2>(); + test<3>(); +} + diff --git a/tests/deal.II/project_01_curved_boundary/cmp/generic b/tests/deal.II/project_01_curved_boundary/cmp/generic new file mode 100644 index 0000000000..218888d980 --- /dev/null +++ b/tests/deal.II/project_01_curved_boundary/cmp/generic @@ -0,0 +1,165 @@ + +DEAL::0 +DEAL::0 +DEAL::0.227756 +DEAL::0.0359680 +DEAL::0.421875 +DEAL::-0.707107 -0.707107 0.195262 +DEAL::0.707107 -0.707107 0.195262 +DEAL::-0.292893 -0.292893 -0.0797954 +DEAL::0.292893 -0.292893 -0.0797954 +DEAL::-0.707107 -0.707107 0.195262 +DEAL::-0.292893 -0.292893 -0.0797954 +DEAL::-0.707107 0.707107 0.195262 +DEAL::-0.292893 0.292893 -0.0797954 +DEAL::-0.292893 -0.292893 -0.0797954 +DEAL::0.292893 -0.292893 -0.0797954 +DEAL::-0.292893 0.292893 -0.0797954 +DEAL::0.292893 0.292893 -0.0797954 +DEAL::0.707107 -0.707107 0.195262 +DEAL::0.707107 0.707107 0.195262 +DEAL::0.292893 -0.292893 -0.0797954 +DEAL::0.292893 0.292893 -0.0797954 +DEAL::-0.707107 0.707107 0.195262 +DEAL::-0.292893 0.292893 -0.0797954 +DEAL::0.707107 0.707107 0.195262 +DEAL::0.292893 0.292893 -0.0797954 +DEAL::0.0359680 +DEAL::-0.707107 -0.707107 0 +DEAL::0.707107 -0.707107 0 +DEAL::-0.292893 -0.292893 0.0179840 +DEAL::0.292893 -0.292893 0.0179840 +DEAL::-0.707107 -0.707107 0 +DEAL::-0.292893 -0.292893 0.0179840 +DEAL::-0.707107 0.707107 0 +DEAL::-0.292893 0.292893 0.0179840 +DEAL::-0.292893 -0.292893 0.0179840 +DEAL::0.292893 -0.292893 0.0179840 +DEAL::-0.292893 0.292893 0.0179840 +DEAL::0.292893 0.292893 0.0179840 +DEAL::0.707107 -0.707107 0 +DEAL::0.707107 0.707107 0 +DEAL::0.292893 -0.292893 0.0179840 +DEAL::0.292893 0.292893 0.0179840 +DEAL::-0.707107 0.707107 0 +DEAL::-0.292893 0.292893 0.0179840 +DEAL::0.707107 0.707107 0 +DEAL::0.292893 0.292893 0.0179840 +DEAL::0 +DEAL::0 +DEAL::0.166690 +DEAL::0.0709393 +DEAL::0.440214 +DEAL::-0.211325 -0.211325 -0.211325 -0.104038 +DEAL::0.211325 -0.211325 -0.211325 -0.104038 +DEAL::-0.211325 0.211325 -0.211325 -0.104038 +DEAL::0.211325 0.211325 -0.211325 -0.104038 +DEAL::-0.211325 -0.211325 0.211325 -0.104038 +DEAL::0.211325 -0.211325 0.211325 -0.104038 +DEAL::-0.211325 0.211325 0.211325 -0.104038 +DEAL::0.211325 0.211325 0.211325 -0.104038 +DEAL::-0.577350 -0.577350 -0.577350 0.115757 +DEAL::0.577350 -0.577350 -0.577350 0.115757 +DEAL::-0.577350 0.577350 -0.577350 0.115757 +DEAL::0.577350 0.577350 -0.577350 0.115757 +DEAL::-0.211325 -0.211325 -0.211325 -0.104038 +DEAL::0.211325 -0.211325 -0.211325 -0.104038 +DEAL::-0.211325 0.211325 -0.211325 -0.104038 +DEAL::0.211325 0.211325 -0.211325 -0.104038 +DEAL::0.577350 -0.577350 -0.577350 0.115757 +DEAL::0.577350 0.577350 -0.577350 0.115757 +DEAL::0.211325 -0.211325 -0.211325 -0.104038 +DEAL::0.211325 0.211325 -0.211325 -0.104038 +DEAL::0.577350 -0.577350 0.577350 0.115757 +DEAL::0.577350 0.577350 0.577350 0.115757 +DEAL::0.211325 -0.211325 0.211325 -0.104038 +DEAL::0.211325 0.211325 0.211325 -0.104038 +DEAL::-0.577350 -0.577350 0.577350 0.115757 +DEAL::0.577350 -0.577350 0.577350 0.115757 +DEAL::-0.211325 -0.211325 0.211325 -0.104038 +DEAL::0.211325 -0.211325 0.211325 -0.104038 +DEAL::-0.577350 0.577350 0.577350 0.115757 +DEAL::0.577350 0.577350 0.577350 0.115757 +DEAL::-0.211325 0.211325 0.211325 -0.104038 +DEAL::0.211325 0.211325 0.211325 -0.104038 +DEAL::-0.577350 -0.577350 -0.577350 0.115757 +DEAL::-0.211325 -0.211325 -0.211325 -0.104038 +DEAL::-0.577350 0.577350 -0.577350 0.115757 +DEAL::-0.211325 0.211325 -0.211325 -0.104038 +DEAL::-0.577350 -0.577350 0.577350 0.115757 +DEAL::-0.211325 -0.211325 0.211325 -0.104038 +DEAL::-0.577350 0.577350 0.577350 0.115757 +DEAL::-0.211325 0.211325 0.211325 -0.104038 +DEAL::-0.577350 -0.577350 -0.577350 0.115757 +DEAL::0.577350 -0.577350 -0.577350 0.115757 +DEAL::-0.211325 -0.211325 -0.211325 -0.104038 +DEAL::0.211325 -0.211325 -0.211325 -0.104038 +DEAL::-0.577350 -0.577350 0.577350 0.115757 +DEAL::0.577350 -0.577350 0.577350 0.115757 +DEAL::-0.211325 -0.211325 0.211325 -0.104038 +DEAL::0.211325 -0.211325 0.211325 -0.104038 +DEAL::-0.577350 0.577350 -0.577350 0.115757 +DEAL::-0.211325 0.211325 -0.211325 -0.104038 +DEAL::0.577350 0.577350 -0.577350 0.115757 +DEAL::0.211325 0.211325 -0.211325 -0.104038 +DEAL::-0.577350 0.577350 0.577350 0.115757 +DEAL::-0.211325 0.211325 0.211325 -0.104038 +DEAL::0.577350 0.577350 0.577350 0.115757 +DEAL::0.211325 0.211325 0.211325 -0.104038 +DEAL::0.0709393 +DEAL::-0.211325 -0.211325 -0.211325 0.0250808 +DEAL::0.211325 -0.211325 -0.211325 0.0250808 +DEAL::-0.211325 0.211325 -0.211325 0.0250808 +DEAL::0.211325 0.211325 -0.211325 0.0250808 +DEAL::-0.211325 -0.211325 0.211325 0.0250808 +DEAL::0.211325 -0.211325 0.211325 0.0250808 +DEAL::-0.211325 0.211325 0.211325 0.0250808 +DEAL::0.211325 0.211325 0.211325 0.0250808 +DEAL::-0.577350 -0.577350 -0.577350 0 +DEAL::0.577350 -0.577350 -0.577350 0 +DEAL::-0.577350 0.577350 -0.577350 0 +DEAL::0.577350 0.577350 -0.577350 0 +DEAL::-0.211325 -0.211325 -0.211325 0.0250808 +DEAL::0.211325 -0.211325 -0.211325 0.0250808 +DEAL::-0.211325 0.211325 -0.211325 0.0250808 +DEAL::0.211325 0.211325 -0.211325 0.0250808 +DEAL::0.577350 -0.577350 -0.577350 0 +DEAL::0.577350 0.577350 -0.577350 0 +DEAL::0.211325 -0.211325 -0.211325 0.0250808 +DEAL::0.211325 0.211325 -0.211325 0.0250808 +DEAL::0.577350 -0.577350 0.577350 0 +DEAL::0.577350 0.577350 0.577350 0 +DEAL::0.211325 -0.211325 0.211325 0.0250808 +DEAL::0.211325 0.211325 0.211325 0.0250808 +DEAL::-0.577350 -0.577350 0.577350 0 +DEAL::0.577350 -0.577350 0.577350 0 +DEAL::-0.211325 -0.211325 0.211325 0.0250808 +DEAL::0.211325 -0.211325 0.211325 0.0250808 +DEAL::-0.577350 0.577350 0.577350 0 +DEAL::0.577350 0.577350 0.577350 0 +DEAL::-0.211325 0.211325 0.211325 0.0250808 +DEAL::0.211325 0.211325 0.211325 0.0250808 +DEAL::-0.577350 -0.577350 -0.577350 0 +DEAL::-0.211325 -0.211325 -0.211325 0.0250808 +DEAL::-0.577350 0.577350 -0.577350 0 +DEAL::-0.211325 0.211325 -0.211325 0.0250808 +DEAL::-0.577350 -0.577350 0.577350 0 +DEAL::-0.211325 -0.211325 0.211325 0.0250808 +DEAL::-0.577350 0.577350 0.577350 0 +DEAL::-0.211325 0.211325 0.211325 0.0250808 +DEAL::-0.577350 -0.577350 -0.577350 0 +DEAL::0.577350 -0.577350 -0.577350 0 +DEAL::-0.211325 -0.211325 -0.211325 0.0250808 +DEAL::0.211325 -0.211325 -0.211325 0.0250808 +DEAL::-0.577350 -0.577350 0.577350 0 +DEAL::0.577350 -0.577350 0.577350 0 +DEAL::-0.211325 -0.211325 0.211325 0.0250808 +DEAL::0.211325 -0.211325 0.211325 0.0250808 +DEAL::-0.577350 0.577350 -0.577350 0 +DEAL::-0.211325 0.211325 -0.211325 0.0250808 +DEAL::0.577350 0.577350 -0.577350 0 +DEAL::0.211325 0.211325 -0.211325 0.0250808 +DEAL::-0.577350 0.577350 0.577350 0 +DEAL::-0.211325 0.211325 0.211325 0.0250808 +DEAL::0.577350 0.577350 0.577350 0 +DEAL::0.211325 0.211325 0.211325 0.0250808 -- 2.39.5