From: leicht Date: Fri, 26 Jan 2007 09:09:55 +0000 (+0000) Subject: New test: test dof reordering if the face is rotated against the standard orientation... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=385f5b01be6cfce3dde2400eb2f34539f196be36;p=dealii-svn.git New test: test dof reordering if the face is rotated against the standard orientation. Test this via projection of a function to a mesh using FE_Q<3>(3) elements. Currently only works on branch_general_meshes. Should be moved to tests/bits eventually. git-svn-id: https://svn.dealii.org/trunk@14383 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/fail/face_orientation_and_fe_q_02.cc b/tests/fail/face_orientation_and_fe_q_02.cc new file mode 100644 index 0000000000..bc1ba72b4c --- /dev/null +++ b/tests/fail/face_orientation_and_fe_q_02.cc @@ -0,0 +1,213 @@ +//---------------------------- face_orientation_and_fe_q_02.cc --------------------------- +// $Id: face_orientation_and_fe_q_01.cc 12732 2006-03-28 23:15:45Z wolf $ +// Version: $Name$ +// +// Copyright (C) 2006, 2007 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- face_orientation_and_fe_q_02.cc --------------------------- + + +// make sure that we treat FE_Q elements correctly if face_flip==true || +// face_rotation==true. for p>=3, we need to revert the order of dofs somehow +// between the two sides of the face. The same applies for lines in non-standard +// orientation. this test is derived from deal.II/project_q_03 + +char logname[] = "face_orientation_and_fe_q_02/output"; + + +#include "../tests.h" +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + + +template +void test (); + + +template +class F : public Function +{ + public: + F (const unsigned int q, + const unsigned int n_components) + : + Function(n_components), + q(q) + {} + + virtual double value (const Point &p, + const unsigned int component) const + { + Assert ((component == 0) && (this->n_components == 1), + ExcInternalError()); + double val = 0; + for (unsigned int d=0; d &p, + Vector &v) const + { + for (unsigned int c=0; c +void do_project (const Triangulation &triangulation, + const FiniteElement &fe, + const unsigned int p, + const unsigned int order_difference) +{ + DoFHandler dof_handler(triangulation); + dof_handler.distribute_dofs (fe); + + deallog << "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 + VectorTools::project (dof_handler, + constraints, + QGauss(p+2), + F (q, fe.n_components()), + 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); + deallog << fe.get_name() << ", P_" << q + << ", rel. error=" << error.l2_norm() / projection.l2_norm() + << std::endl; + + if (q<=p-order_difference) + Assert (error.l2_norm() <= 1e-10*projection.l2_norm(), + ExcFailedProjection(error.l2_norm() / projection.l2_norm())); + } +} + + + +// test with a 3d grid that has cells with face_rotation==false || +// face_flip==false and hanging nodes. maybe it triggers all sorts of +// assumptions that may be hidden in places +// +// the mesh we use is a 7 cell moebius-type in 3d, with the first and last +// cell refined in turn. that then makes 2 meshes with 14 active cells +// each. this also cycles through all possibilities of coarser or finer cell +// having face_rotation==false || face_flip==false +// +// the whole procedure is repeated for three combinations of face_rotation and +// face_flip (the standard case in which evereything is as usual is left out) +template +void test_with_wrong_face_orientation (const FiniteElement &fe, + const unsigned int p, + const unsigned int order_difference = 0) +{ + if (dim != 3) + return; + for (unsigned int j=1; j<4; ++j) + // j=1: face_rotation=true, face_flip=true + // j=2: face_rotation=false, face_flip=true + // j=3: face_rotation=true, face_flip=false + for (unsigned int i=0; i<2; ++i) + { + Triangulation triangulation; + GridGenerator::moebius (triangulation,7,j,1.0,0.2); + typename Triangulation::active_cell_iterator + cell = triangulation.begin_active(); + if (i==0) + { + std::advance (cell,6); + deallog<<"face_rotation="<face_rotation(2) + <<", face_flip="<face_flip(2)<set_refine_flag (); + triangulation.execute_coarsening_and_refinement (); + + do_project (triangulation, fe, p, order_difference); + } +} + + + + +int main () +{ + std::ofstream logfile(logname); + logfile.precision (3); + + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + test<3>(); +} + + + +template +void test () +{ + test_with_wrong_face_orientation (FE_Q(3), 3); +} diff --git a/tests/fail/face_orientation_and_fe_q_02/cmp/generic b/tests/fail/face_orientation_and_fe_q_02/cmp/generic new file mode 100644 index 0000000000..83ed06a012 --- /dev/null +++ b/tests/fail/face_orientation_and_fe_q_02/cmp/generic @@ -0,0 +1,46 @@ + +DEAL::face_rotation=1, face_flip=1 +DEAL::n_dofs=639 +DEAL::FE_Q<3>(3), P_0, rel. error=0 +DEAL::FE_Q<3>(3), P_1, rel. error=0 +DEAL::FE_Q<3>(3), P_2, rel. error=0 +DEAL::FE_Q<3>(3), P_3, rel. error=0 +DEAL::FE_Q<3>(3), P_4, rel. error=2.53e-05 +DEAL::FE_Q<3>(3), P_5, rel. error=4.19e-05 +DEAL::n_dofs=639 +DEAL::FE_Q<3>(3), P_0, rel. error=0 +DEAL::FE_Q<3>(3), P_1, rel. error=0 +DEAL::FE_Q<3>(3), P_2, rel. error=0 +DEAL::FE_Q<3>(3), P_3, rel. error=0 +DEAL::FE_Q<3>(3), P_4, rel. error=2.93e-05 +DEAL::FE_Q<3>(3), P_5, rel. error=5.91e-05 +DEAL::face_rotation=0, face_flip=1 +DEAL::n_dofs=639 +DEAL::FE_Q<3>(3), P_0, rel. error=0 +DEAL::FE_Q<3>(3), P_1, rel. error=0 +DEAL::FE_Q<3>(3), P_2, rel. error=0 +DEAL::FE_Q<3>(3), P_3, rel. error=0 +DEAL::FE_Q<3>(3), P_4, rel. error=2.50e-05 +DEAL::FE_Q<3>(3), P_5, rel. error=4.14e-05 +DEAL::n_dofs=639 +DEAL::FE_Q<3>(3), P_0, rel. error=0 +DEAL::FE_Q<3>(3), P_1, rel. error=0 +DEAL::FE_Q<3>(3), P_2, rel. error=0 +DEAL::FE_Q<3>(3), P_3, rel. error=0 +DEAL::FE_Q<3>(3), P_4, rel. error=2.89e-05 +DEAL::FE_Q<3>(3), P_5, rel. error=5.85e-05 +DEAL::face_rotation=1, face_flip=0 +DEAL::n_dofs=639 +DEAL::FE_Q<3>(3), P_0, rel. error=0 +DEAL::FE_Q<3>(3), P_1, rel. error=0 +DEAL::FE_Q<3>(3), P_2, rel. error=0 +DEAL::FE_Q<3>(3), P_3, rel. error=0 +DEAL::FE_Q<3>(3), P_4, rel. error=2.43e-05 +DEAL::FE_Q<3>(3), P_5, rel. error=4.07e-05 +DEAL::n_dofs=639 +DEAL::FE_Q<3>(3), P_0, rel. error=0 +DEAL::FE_Q<3>(3), P_1, rel. error=0 +DEAL::FE_Q<3>(3), P_2, rel. error=0 +DEAL::FE_Q<3>(3), P_3, rel. error=0 +DEAL::FE_Q<3>(3), P_4, rel. error=2.81e-05 +DEAL::FE_Q<3>(3), P_5, rel. error=5.76e-05