From 97596dea7552c448ab59396d8adda3b0abdcdb8f Mon Sep 17 00:00:00 2001 From: Tobias Leicht Date: Thu, 16 Aug 2007 07:43:30 +0000 Subject: [PATCH] New test, that currently fails: The neighbor_of_neighbor function does not work properly if two cells are neighbors at more than one face. git-svn-id: https://svn.dealii.org/trunk@14963 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/bits/Makefile | 3 +- tests/bits/neighboring_cells_at_two_faces.cc | 101 ++++++++++++++++++ .../cmp/generic | 3 + 3 files changed, 106 insertions(+), 1 deletion(-) create mode 100644 tests/bits/neighboring_cells_at_two_faces.cc create mode 100644 tests/bits/neighboring_cells_at_two_faces/cmp/generic diff --git a/tests/bits/Makefile b/tests/bits/Makefile index a4d698ccab..164ba8b138 100644 --- a/tests/bits/Makefile +++ b/tests/bits/Makefile @@ -89,7 +89,8 @@ tests_x = grid_generator_?? \ tria_crash_* \ face_orientation_* \ max_n_cells_* \ - joa_* + joa_* \ + neighboring_cells_at_two_faces # tests for the hp branch: # fe_collection_* diff --git a/tests/bits/neighboring_cells_at_two_faces.cc b/tests/bits/neighboring_cells_at_two_faces.cc new file mode 100644 index 0000000000..fa34fffa35 --- /dev/null +++ b/tests/bits/neighboring_cells_at_two_faces.cc @@ -0,0 +1,101 @@ +//--------------- neighboring_cells_at_two_faces.cc ---------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 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. +// +//--------------- neighboring_cells_at_two_faces.cc ---------------- + +#include "../tests.h" +#include +#include +#include +#include +#include +#include + +#include +#include + + +// generate two degenerate quads which reduce to the form of triangles, where +// two sides of the quad form one side of the triangle. the neighboring quads +// are neighboring along both of these sides. +// +// *-------* +// |1 /| +// | a/ | +// | / | +// | * | +// | / | +// | /b | +// |/ 2| +// *-------* +// +// looking from cell 1 at face a and asking for the corresponding face (number) +// of cell 2 neighbor_of_neighbor returns something corrwesponding to face b, +// not a, as it should. Simply looking at the identity of the neighboring cell +// is not enough, we have to look at the (index of the) face instead. + +void create_grid (Triangulation<2> &tria) +{ + const unsigned int n_points=5; + + const Point<2> points[n_points] = { Point<2>(0.0,0.0), + Point<2>(1.0,1.0), + Point<2>(2.0,2.0), + Point<2>(0.0,2.0), + Point<2>(2.0,0.0)}; + std::vector > vertices(n_points); + vertices.assign(points,points+n_points); + + std::vector > cells(2); + cells[0].vertices[0] = 0; + cells[0].vertices[1] = 1; + cells[0].vertices[2] = 2; + cells[0].vertices[3] = 3; + cells[1].vertices[0] = 0; + cells[1].vertices[1] = 4; + cells[1].vertices[2] = 2; + cells[1].vertices[3] = 1; + + // generate a triangulation + // out of this + GridReordering<2>::reorder_cells (cells); + tria.create_triangulation_compatibility (vertices, cells, SubCellData()); +} + + +void check_neighbors (const Triangulation<2> &tria) +{ + Triangulation<2>::cell_iterator cell = tria.begin(); + for (unsigned int f=0; f::faces_per_cell; ++f) + if (cell->neighbor(f).state()==IteratorState::valid) + { + const unsigned int neighbor_neighbor=cell->neighbor_of_neighbor(f); + deallog << "At face " << f << ": neighbor_of_neighbor=" + << neighbor_neighbor << std::endl; + Assert(cell->face(f)==cell->neighbor(f)->face(neighbor_neighbor), + ExcMessage("Error in neighbor_of_neighbor() function!")); + } +} + + + +int main() +{ + std::ofstream logfile("neighboring_cells_at_two_faces/output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + Triangulation<2> tria; + create_grid(tria); + check_neighbors(tria); +} + diff --git a/tests/bits/neighboring_cells_at_two_faces/cmp/generic b/tests/bits/neighboring_cells_at_two_faces/cmp/generic new file mode 100644 index 0000000000..576ed308b8 --- /dev/null +++ b/tests/bits/neighboring_cells_at_two_faces/cmp/generic @@ -0,0 +1,3 @@ + +DEAL::At face 1: neighbor_of_neighbor=3 +DEAL::At face 2: neighbor_of_neighbor=0 -- 2.39.5