From cb884d453d2c30f3eeebeb044773a13106ddf4d0 Mon Sep 17 00:00:00 2001 From: bangerth Date: Mon, 8 Jul 2013 03:51:51 +0000 Subject: [PATCH] Yet another version of this test. It should be producing the same result every time, but it doesn't with my current version of work_stream.h. git-svn-id: https://svn.dealii.org/trunk@29948 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/base/work_stream_03.cc | 68 ++++++++---------------------------- 1 file changed, 15 insertions(+), 53 deletions(-) diff --git a/tests/base/work_stream_03.cc b/tests/base/work_stream_03.cc index 618700df8d..f912c0fc46 100644 --- a/tests/base/work_stream_03.cc +++ b/tests/base/work_stream_03.cc @@ -57,11 +57,9 @@ template class F : public Function { public: - F (const unsigned int q, - const unsigned int n_components) + F () : - Function(n_components), - q(q) + q(1) {} virtual double value (const Point &p, @@ -76,19 +74,6 @@ class F : public Function return val; } - - virtual void vector_value (const Point &p, - Vector &v) const - { - for (unsigned int c=0; c template void do_project (const Triangulation &triangulation, - const FiniteElement &fe, - const unsigned int p, - const unsigned int order_difference) + const FiniteElement &fe) { DoFHandler dof_handler(triangulation); dof_handler.distribute_dofs (fe); - std::cout << "n_dofs=" << dof_handler.n_dofs() << std::endl; - ConstraintMatrix constraints; DoFTools::make_hanging_node_constraints (dof_handler, constraints); constraints.close (); Vector projection (dof_handler.n_dofs()); - Vector error (triangulation.n_active_cells()); - for (unsigned int q=0; q<=p+2-order_difference; ++q) - { - // project the function + + // project the function VectorTools::project (dof_handler, constraints, - QGauss(p+2), - F (q, fe.n_components()), + QGauss(3), + F (), projection); // just to make sure it doesn't get // forgotten: handle hanging node // constraints constraints.distribute (projection); - // then compute the interpolation error - VectorTools::integrate_difference (dof_handler, - projection, - F (q, fe.n_components()), - error, - QGauss(std::max(p,q)+1), - VectorTools::L2_norm); - std::cout << fe.get_name() << ", P_" << q - << ", rel. error=" << error.l2_norm() / projection.l2_norm() - << std::endl; - - if (q<=p-order_difference) - if (error.l2_norm() > 1e-10*projection.l2_norm()) - { - std::cout << "Projection failed with relative error " - << error.l2_norm() / projection.l2_norm() - << std::endl; - Assert (false, ExcInternalError()); - } - } + double sum=0; + for (unsigned int i=0; i &triangulation, // of polynomial degree p has normal components of degree p-1 and therefore // can only represent polynomials of degree p-1 exactly. the gap is then 1. template -void test_no_hanging_nodes (const FiniteElement &fe, - const unsigned int p, - const unsigned int order_difference = 0) +void test_no_hanging_nodes (const FiniteElement &fe) { Triangulation triangulation; GridGenerator::hyper_cube (triangulation); triangulation.refine_global (3); - do_project (triangulation, fe, p, order_difference); + for (unsigned int i=0; i<12; ++i) + do_project (triangulation, fe); } @@ -170,7 +132,7 @@ void test_no_hanging_nodes (const FiniteElement &fe, template void test () { - test_no_hanging_nodes (FE_Q(1), 1); + test_no_hanging_nodes (FE_Q(1)); } -- 2.39.5