From: Wolfgang Bangerth Date: Mon, 29 Sep 2003 19:35:39 +0000 (+0000) Subject: New test. X-Git-Tag: v8.0.0~16147 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=3b7b36eb3c6340ce1bef29bec6b69b95dff4c292;p=dealii.git New test. git-svn-id: https://svn.dealii.org/trunk@8063 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/bits/q_points.cc b/tests/bits/q_points.cc new file mode 100644 index 0000000000..3c826dd382 --- /dev/null +++ b/tests/bits/q_points.cc @@ -0,0 +1,135 @@ +//---------------------------- q_points.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003 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. +// +//---------------------------- q_points.cc --------------------------- + + +// make sure that if we loop over all quadrature points on a face and +// over the same quadrature points on the same face of the neighboring +// cell, that we then visit them in the same order. this test is +// basically the same as the mesh_3d_7 test, except that we leave out +// the complication of mis-oriented faces + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include + + +void create_two_cubes (Triangulation<3> &coarse_grid) +{ + const Point<3> points[6] = { Point<3>(0,0,0), + Point<3>(1,0,0), + Point<3>(1,1,0), + Point<3>(0,1,0), + Point<3>(2,0,0), + Point<3>(2,1,0)}; + std::vector > vertices; + for (unsigned int i=0; i<6; ++i) + vertices.push_back (points[i]); + for (unsigned int i=0; i<6; ++i) + vertices.push_back (points[i] + Point<3>(0,0,-1)); + + std::vector > cells; + + const unsigned int connectivity[2][4] + = { { 0,1,2,3 }, + { 1,4,5,2 } }; + for (unsigned int i=0; i<2; ++i) + { + CellData<3> cell; + for (unsigned int j=0; j<4; ++j) + { + cell.vertices[j] = connectivity[i][j]; + cell.vertices[j+4] = connectivity[i][j]+6; + } + cells.push_back (cell); + } + + // finally generate a triangulation + // out of this + GridReordering<3>::reorder_cells (cells); + coarse_grid.create_triangulation (vertices, cells, SubCellData()); +} + + +void check (Triangulation<3> &tria) +{ + QGauss3<2> quadrature; + FE_Q<3> fe(1); + FEFaceValues<3> fe_face_values1 (fe, quadrature, + update_q_points | update_JxW_values); + FEFaceValues<3> fe_face_values2 (fe, quadrature, + update_q_points | update_JxW_values); + + DoFHandler<3> dof_handler (tria); + dof_handler.distribute_dofs (fe); + + for (DoFHandler<3>::cell_iterator cell=dof_handler.begin(); + cell != dof_handler.end(); ++cell) + for (unsigned int f=0; f::faces_per_cell; ++f) + if (!cell->at_boundary(f)) + { + const unsigned int nn + = cell->neighbor_of_neighbor (f); + + // this test is about + // properly oriented + // faces. mesh_3d_7 does it + // for mis-oriented faces + Assert (cell->get_face_orientation(f) == true, + ExcInternalError()); + Assert (cell->neighbor(f)->get_face_orientation(nn) == true, + ExcInternalError()); + + fe_face_values1.reinit (cell, f); + fe_face_values2.reinit (cell->neighbor(f), nn); + + deallog << "Cell " << cell << ", face " << f + << std::endl; + + for (unsigned int q=0; q coarse_grid; + create_two_cubes (coarse_grid); + check (coarse_grid); + + coarse_grid.refine_global (1); + check (coarse_grid); +} + + +