From: Wolfgang Bangerth Date: Mon, 30 Jun 2014 10:23:46 +0000 (+0000) Subject: Add another testcase that currently fails, see issue 197 (https://code.google.com... X-Git-Tag: v8.2.0-rc1~348 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b50a3b8d605c7809dc072af7bf08f258c9453b98;p=dealii.git Add another testcase that currently fails, see issue 197 (https://code.google.com/p/dealii/issues/detail?id=197). git-svn-id: https://svn.dealii.org/trunk@33100 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/fail/merge_triangulations_03.cc b/tests/fail/merge_triangulations_03.cc new file mode 100644 index 0000000000..f263da0f40 --- /dev/null +++ b/tests/fail/merge_triangulations_03.cc @@ -0,0 +1,172 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2014 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +// This script reproduces an issue with merging a hyper ball within a +// hyper shell using dealII-8.0.0. +// +// Test by Bamdad Hosseini (see mail to mailing +// list on 6/20/2014) + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +std::ofstream logfile("output"); + + +template +void flatten_triangulation( Triangulation &tria_in, Triangulation &tria_out ) +// takes a possibly refined triangulation and returns a new coarse triangulation +// created from the vertices of the input +{ + + // we start by extracting the number of vertices and cell is + // the original triangulation + + unsigned int n_vertices = tria_in.n_vertices(); + unsigned int n_cells = tria_in.n_cells(); + unsigned int n_active_cells = tria_in.n_active_cells(); + + // also get all vertices + const std::vector< Point > vertices = tria_in.get_vertices(); + + // now we need to loop over all cells and extract the cell data. + typename Triangulation::active_cell_iterator + cell=tria_in.begin_active(), + endc=tria_in.end(); + + std::vector< CellData > cells(n_active_cells, CellData()); + + unsigned int cell_counter = 0; + for(; cell!=endc; ++cell) + { + for (unsigned int j=0; + j < GeometryInfo::vertices_per_cell; + ++j) { + cells[cell_counter].vertices[j] = cell->vertex_index(j); + } + cells[cell_counter].material_id=0; + ++cell_counter; + } + + // pass the info to create_triangulation to create the flat mesh + tria_out.create_triangulation(vertices, cells, SubCellData() ); +} + +// this is a function to provide some basic information about the +// triangulations. The code inside can be ignored. It produces an +// eps image of the mesh so we can see the nodes are at the right +// place. + +template +void mesh_info(const Triangulation &tria, + const std::string &filename) +{ + deallog << "Mesh info for " << filename << ":" << std::endl + << " dimension: " << dim << std::endl + << " no. of cells: " << tria.n_active_cells() << std::endl; + + { + std::map boundary_count; + typename Triangulation::active_cell_iterator + cell = tria.begin_active(), + endc = tria.end(); + for (; cell!=endc; ++cell) + { + for (unsigned int face=0; + face::faces_per_cell; + ++face) + { + if (cell->face(face)->at_boundary()) + boundary_count[cell->face(face)->boundary_indicator()]++; + } + } + deallog << " boundary indicators: "; + for (std::map::iterator it=boundary_count.begin(); + it!=boundary_count.end(); + ++it) + { + deallog << it->first << "(" << it->second << " times) "; + } + deallog << std::endl; + } +} + + +void test () +{ + Point<2> center; + + Triangulation<2> ball; + Triangulation<2> flat_ball; + Triangulation<2> shell; + Triangulation<2> tria_out; + + // create the two meshes + GridGenerator::hyper_ball( ball, center, 0.5 ); + GridGenerator::hyper_shell(shell, center, 0.5, 1, 8, true); //colorize flag set to true so outer bnd is 1 and inner is 0 + + // set boundaries + static const HyperBallBoundary<2> inner_bnd(center, 0.5); + static const HyperBallBoundary<2> outer_bnd(center, 1); + + ball.set_boundary(0, inner_bnd); + shell.set_boundary(0, inner_bnd); + shell.set_boundary(1, outer_bnd); + + // first we need to refine the ball once to get the right + // nodes to merge the meshes + ball.refine_global(1); + + // flatten the refined ball so we have a coarse mesh to merge + flatten_triangulation( ball, flat_ball ); + + mesh_info( shell, "shell.eps" ); + mesh_info( ball, "ball.eps" ); + mesh_info( flat_ball, "flat_ball.eps" ); + + // now we merge (this throws an exception with no error messages + GridGenerator::merge_triangulations( flat_ball, shell, tria_out ); + mesh_info(tria_out, "tria_out.eps"); +} + + +int main () +{ + deallog << std::setprecision(2); + logfile << std::setprecision(2); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + test (); + + return 0; +}