From: maier Date: Tue, 22 May 2012 15:42:07 +0000 (+0000) Subject: A small test for DoFTools::make_periodicity_constraints X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=5841ff3a58bfcfbb6d120eeb9b82a04957ada576;p=dealii-svn.git A small test for DoFTools::make_periodicity_constraints git-svn-id: https://svn.dealii.org/trunk@25539 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/bits/dof_tools_21.cc b/tests/bits/dof_tools_21.cc new file mode 100644 index 0000000000..f840571896 --- /dev/null +++ b/tests/bits/dof_tools_21.cc @@ -0,0 +1,83 @@ +//---------------------------- grid_tools.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2008, 2011, 2012 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. +// +//---------------------------- dof_tools_21.cc --------------------------- + +#include "../tests.h" +#include "dof_tools_periodic.h" +#include +#include +#include +#include +#include +#include + +#include +#include +#include + +// A simple test for +// DoFTools:: +// make_periodicity_constraints (const FaceIterator &, +// const FaceIterator &, +// dealii::ConstraintMatrix &, +// const std::vector &) +// +// We project an already periodic function onto the FE space of +// periodic functions and store the resulting L2 difference +// between the projected and unprojected function. +// This should reveal any errors in the constraint_matrix: + +std::string output_file_name("dof_tools_21/output"); + + +template +void +check_this (const DoFHandler &dof_handler) +{ + Functions::CosineFunction test_func (dof_handler.get_fe().n_components ()); + + ConstraintMatrix cm; + + // Apply periodic boundary conditions only in the one direction where + // we can match the (locally refined) faces: + DoFTools::make_periodicity_constraints (dof_handler.begin(0)->face(0), + dof_handler.begin(0)->face(1), + cm); + cm.close (); + + deallog << cm.n_constraints () << std::endl; + deallog << cm.max_constraint_indirections () << std::endl; + + QGauss quadrature (6); + + Vector unconstrained (dof_handler.n_dofs ()); + Vector constrained (dof_handler.n_dofs ()); + + // Ensure that we can interpolate: + if (dof_handler.get_fe().get_unit_support_points().size() == 0) + return; + + VectorTools::interpolate (dof_handler, + test_func, + unconstrained); + + constrained = unconstrained; + cm.distribute (constrained); + + constrained -= unconstrained; + + const double p_l2_error = constrained.l2_norm(); + + Assert (p_l2_error < 1e-11, ExcInternalError()); + + deallog << "L2_Error : " << p_l2_error << std::endl; +} diff --git a/tests/bits/dof_tools_21/cmp/generic b/tests/bits/dof_tools_21/cmp/generic new file mode 100644 index 0000000000..7ac94c50bf --- /dev/null +++ b/tests/bits/dof_tools_21/cmp/generic @@ -0,0 +1,159 @@ + +DEAL::Checking Q1 in 2d: +DEAL::7 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking Q1 in 3d: +DEAL::41 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking Q2 in 2d: +DEAL::13 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking Q2 in 3d: +DEAL::141 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking Q3 in 2d: +DEAL::19 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking Q3 in 3d: +DEAL::297 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking DGQ0 in 2d: +DEAL::0 +DEAL::0 +DEAL::L2_Error : 0 +DEAL::Checking DGQ0 in 3d: +DEAL::0 +DEAL::0 +DEAL::L2_Error : 0 +DEAL::Checking DGQ1 in 2d: +DEAL::0 +DEAL::0 +DEAL::L2_Error : 0 +DEAL::Checking DGQ1 in 3d: +DEAL::0 +DEAL::0 +DEAL::L2_Error : 0 +DEAL::Checking DGQ3 in 2d: +DEAL::0 +DEAL::0 +DEAL::L2_Error : 0 +DEAL::Checking DGQ3 in 3d: +DEAL::0 +DEAL::0 +DEAL::L2_Error : 0 +DEAL::Checking DGP0 in 2d: +DEAL::0 +DEAL::0 +DEAL::Checking DGP0 in 3d: +DEAL::0 +DEAL::0 +DEAL::Checking DGP1 in 2d: +DEAL::0 +DEAL::0 +DEAL::Checking DGP1 in 3d: +DEAL::0 +DEAL::0 +DEAL::Checking DGP3 in 2d: +DEAL::0 +DEAL::0 +DEAL::Checking DGP3 in 3d: +DEAL::0 +DEAL::0 +DEAL::Checking Nedelec0 in 2d: +DEAL::6 +DEAL::1 +DEAL::Checking Nedelec0 in 3d: +DEAL::72 +DEAL::1 +DEAL::Checking RaviartThomas0 in 2d: +DEAL::6 +DEAL::1 +DEAL::Checking RaviartThomas1 in 2d: +DEAL::12 +DEAL::1 +DEAL::Checking RaviartThomas2 in 2d: +DEAL::18 +DEAL::1 +DEAL::Checking RaviartThomasNodal0 in 2d: +DEAL::6 +DEAL::1 +DEAL::Checking RaviartThomasNodal1 in 2d: +DEAL::12 +DEAL::1 +DEAL::Checking RaviartThomasNodal2 in 2d: +DEAL::18 +DEAL::1 +DEAL::Checking RaviartThomasNodal0 in 3d: +DEAL::28 +DEAL::1 +DEAL::Checking RaviartThomasNodal1 in 3d: +DEAL::112 +DEAL::1 +DEAL::Checking RaviartThomasNodal2 in 3d: +DEAL::252 +DEAL::1 +DEAL::Checking FE_Q<2>(1)3 in 2d: +DEAL::21 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking FE_DGQ<2>(2)2 in 2d: +DEAL::0 +DEAL::0 +DEAL::L2_Error : 0 +DEAL::Checking FE_DGP<2>(3)1 in 2d: +DEAL::0 +DEAL::0 +DEAL::Checking FE_Q<3>(1)3 in 3d: +DEAL::123 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking FE_DGQ<3>(2)2 in 3d: +DEAL::0 +DEAL::0 +DEAL::L2_Error : 0 +DEAL::Checking FE_DGP<3>(3)1 in 3d: +DEAL::0 +DEAL::0 +DEAL::Checking FE_Q<2>(1)3FE_DGQ<2>(2)2 in 2d: +DEAL::21 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking FE_DGQ<2>(2)2FE_DGP<2>(3)1 in 2d: +DEAL::0 +DEAL::0 +DEAL::Checking FE_DGP<2>(3)1FE_DGQ<2>(2)2 in 2d: +DEAL::0 +DEAL::0 +DEAL::Checking FE_Q<2>(1)3FE_DGP<2>(3)1FE_Q<2>(1)3 in 2d: +DEAL::42 +DEAL::1 +DEAL::Checking FE_DGQ<2>(2)2FE_DGQ<2>(2)2FE_Q<2>(3)3 in 2d: +DEAL::57 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking FE_DGP<2>(3)1FE_DGP<2>(3)1FE_Q<2>(2)3 in 2d: +DEAL::39 +DEAL::1 +DEAL::Checking FE_DGQ<3>(1)3FE_DGP<3>(3)1FE_Q<3>(1)3 in 3d: +DEAL::123 +DEAL::1 +DEAL::Checking (FESystem<2>(FE_Q<2>(1),3))3FE_DGQ<2>(0)1FE_Q<2>(1)3 in 2d: +DEAL::84 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking FE_DGQ<2>(3)1FESystem<2>(FE_DGQ<2>(0),3)1FESystem<2>(FE_Q<2>(2),1, FE_DGQ<2>(0),1)2 in 2d: +DEAL::26 +DEAL::1 +DEAL::L2_Error : 0 +DEAL::Checking FE_DGQ<2>(3)1FE_Nedelec<2>(0)2 in 2d: +DEAL::12 +DEAL::1 +DEAL::Checking FE_Nedelec<2>(0)1FESystem<2>(FE_DGQ<2>(1),2)1FESystem<2>(FE_Q<2>(2),1, FE_Nedelec<2>(0),2)2 in 2d: +DEAL::56 +DEAL::1 diff --git a/tests/bits/dof_tools_periodic.h b/tests/bits/dof_tools_periodic.h new file mode 100644 index 0000000000..db2eb4c153 --- /dev/null +++ b/tests/bits/dof_tools_periodic.h @@ -0,0 +1,243 @@ +//---------------------------- dof_tools_common.h --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003, 2004, 2005, 2006, 2008, 2010, 2011, 2012 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. +// +//---------------------------- dof_tools_periodic.h --------------------------- + + +// Code duplication: +// This is a copy of dof_tools_common.h with a slightly different mesh +// generation to fit the assumptions made by make_periodicity_constraints + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + + +// forward declaration of the function that must be provided in the +// .cc files +template +void +check_this (const DoFHandler &dof_handler); + +// forward declaration of a variable with the name of the output file +extern std::string output_file_name; + + +void +output_bool_vector (std::vector &v) +{ +for (unsigned int i=0; i &) +{} + + + +template +void +check (const FiniteElement &fe, + const std::string &name) +{ + deallog << "Checking " << name + << " in " << dim << "d:" + << std::endl; + + // create tria and dofhandler + // objects. set different boundary + // and sub-domain ids + Triangulation tria; + GridGenerator::hyper_cube(tria, -2., 2.); + tria.refine_global (1); + // refine the first two cells: + for (int i=0; i<2; ++i) + { + tria.begin_active()->set_refine_flag(); + tria.execute_coarsening_and_refinement (); + } + tria.refine_global (1); + for (typename Triangulation::active_cell_iterator + cell=tria.begin_active(); + cell!=tria.end(); ++cell) + cell->set_subdomain_id (cell->level()); + DoFHandler dof_handler (tria); + dof_handler.distribute_dofs (fe); + + // call main function in .cc files + check_this (dof_handler); +} + + + + + +#define CHECK(EL,deg,dim)\ + { FE_ ## EL EL(deg); \ + check(EL, #EL #deg); } + +#define CHECK_SYS1(sub1,N1,dim) \ + { FESystem q(sub1, N1); \ + check(q, #sub1 #N1); } + +#define CHECK_SYS2(sub1,N1,sub2,N2,dim) \ + { FESystem q(sub1, N1, sub2, N2); \ + check(q, #sub1 #N1 #sub2 #N2); } + +#define CHECK_SYS3(sub1,N1,sub2,N2,sub3,N3,dim) \ + { FESystem q(sub1, N1, sub2, N2, sub3, N3); \ + check(q, #sub1 #N1 #sub2 #N2 #sub3 #N3); } + + +#define CHECK_ALL(EL,deg)\ + { /*CHECK(EL,deg,1);*/ \ + CHECK(EL,deg,2); \ + CHECK(EL,deg,3); \ + } + + +int +main() +{ + try + { + std::ofstream logfile(output_file_name.c_str()); + logfile << std::setprecision (2); + deallog << std::setprecision (2); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + CHECK_ALL(Q,1); + CHECK_ALL(Q,2); + CHECK_ALL(Q,3); + + CHECK_ALL(DGQ,0); + CHECK_ALL(DGQ,1); + CHECK_ALL(DGQ,3); + + CHECK_ALL(DGP,0); + CHECK_ALL(DGP,1); + CHECK_ALL(DGP,3); + + CHECK(Nedelec, 0, 2); + CHECK(Nedelec, 0, 3); + + CHECK(RaviartThomas, 0, 2); + CHECK(RaviartThomas, 1, 2); + CHECK(RaviartThomas, 2, 2); + + CHECK(RaviartThomasNodal, 0, 2); + CHECK(RaviartThomasNodal, 1, 2); + CHECK(RaviartThomasNodal, 2, 2); + CHECK(RaviartThomasNodal, 0, 3); + CHECK(RaviartThomasNodal, 1, 3); + CHECK(RaviartThomasNodal, 2, 3); + + /*CHECK_SYS1(FE_Q<1>(1), 3,1); + CHECK_SYS1(FE_DGQ<1>(2),2,1); + CHECK_SYS1(FE_DGP<1>(3),1,1);*/ + + CHECK_SYS1(FE_Q<2>(1), 3,2); + CHECK_SYS1(FE_DGQ<2>(2),2,2); + CHECK_SYS1(FE_DGP<2>(3),1,2); + + CHECK_SYS1(FE_Q<3>(1), 3,3); + CHECK_SYS1(FE_DGQ<3>(2),2,3); + CHECK_SYS1(FE_DGP<3>(3),1,3); + + + /*CHECK_SYS2(FE_Q<1>(1), 3,FE_DGQ<1>(2),2,1); + CHECK_SYS2(FE_DGQ<1>(2),2,FE_DGP<1>(3),1,1); + CHECK_SYS2(FE_DGP<1>(3),1,FE_DGQ<1>(2),2,1);*/ + + CHECK_SYS2(FE_Q<2>(1), 3,FE_DGQ<2>(2),2,2); + CHECK_SYS2(FE_DGQ<2>(2),2,FE_DGP<2>(3),1,2); + CHECK_SYS2(FE_DGP<2>(3),1,FE_DGQ<2>(2),2,2); + + /*CHECK_SYS3(FE_Q<1>(1), 3,FE_DGP<1>(3),1,FE_Q<1>(1),3,1); + CHECK_SYS3(FE_DGQ<1>(2),2,FE_DGQ<1>(2),2,FE_Q<1>(3),3,1); + CHECK_SYS3(FE_DGP<1>(3),1,FE_DGP<1>(3),1,FE_Q<1>(2),3,1);*/ + + CHECK_SYS3(FE_Q<2>(1), 3,FE_DGP<2>(3),1,FE_Q<2>(1),3,2); + CHECK_SYS3(FE_DGQ<2>(2),2,FE_DGQ<2>(2),2,FE_Q<2>(3),3,2); + CHECK_SYS3(FE_DGP<2>(3),1,FE_DGP<2>(3),1,FE_Q<2>(2),3,2); + + CHECK_SYS3(FE_DGQ<3>(1), 3,FE_DGP<3>(3),1,FE_Q<3>(1),3,3); + + // systems of systems + CHECK_SYS3((FESystem<2>(FE_Q<2>(1),3)), 3, + FE_DGQ<2>(0), 1, + FE_Q<2>(1), 3, + 2); + CHECK_SYS3(FE_DGQ<2>(3), 1, + FESystem<2>(FE_DGQ<2>(0),3), 1, + FESystem<2>(FE_Q<2>(2),1, + FE_DGQ<2>(0),1),2, + 2); + + // systems with Nedelec elements + CHECK_SYS2 (FE_DGQ<2>(3), 1, + FE_Nedelec<2>(0), 2, + 2); + CHECK_SYS3(FE_Nedelec<2>(0), 1, + FESystem<2>(FE_DGQ<2>(1),2), 1, + FESystem<2>(FE_Q<2>(2),1, + FE_Nedelec<2>(0),2),2, + 2); + + return 0; + } + catch (std::exception &exc) + { + deallog << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + deallog << "Exception on processing: " << std::endl + << exc.what() << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + } + catch (...) + { + deallog << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + deallog << "Unknown exception!" << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + }; +} +