From: pauletti Date: Thu, 18 Nov 2010 17:37:27 +0000 (+0000) Subject: Extracting a piece of the boundary of the sphere does not seem X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=bc6b832cebac0ef2f26b1233e7ba3781e9f5b36f;p=dealii-svn.git Extracting a piece of the boundary of the sphere does not seem to project some nodes properly git-svn-id: https://svn.dealii.org/trunk@22809 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/codim_one/extract_boundary_mesh_04.cc b/tests/codim_one/extract_boundary_mesh_04.cc new file mode 100644 index 0000000000..6a03470313 --- /dev/null +++ b/tests/codim_one/extract_boundary_mesh_04.cc @@ -0,0 +1,138 @@ +//---------------------------- extract_boundary_mesh_03.cc --------------------------- +// $Id: bem.cc 22693 2010-11-11 20:11:47Z kanschat $ +// Version: $Name$ +// +// Copyright (C) 2010 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. +// +//---------------------------- extract_boundary_mesh_03.cc --------------------------- + +/* + Code for testing the function + GridTools::extract_boundary_mesh (...). + We test that the order of cells and the orientation + of the vertices is consistent between the two meshes. + + Like _03, but set volume boundary indicator to 1 only on the upper + half of the sphere . +*/ + +#include "../tests.h" + +#include +#include +#include +#include +#include + +using namespace std; + + +template +void test_vertices_orientation(const Triangulation &boundary_mesh, + map< typename Triangulation::cell_iterator, + typename Triangulation::face_iterator > + &surface_to_volume_mapping) +{ + typename Triangulation::active_cell_iterator + cell = boundary_mesh.begin_active(), + endc = boundary_mesh.end(); + typename Triangulation::face_iterator face; + + for (; cell!=endc; ++cell){ + + face = surface_to_volume_mapping [cell]; + Assert (face->at_boundary(), ExcInternalError()); + + deallog << "Surface cell: " << cell << " with vertices:" << std::endl; + for (unsigned int k=0; k::vertices_per_cell; ++k) + { + deallog << " " << cell->vertex(k) << std::endl; + Assert (std::fabs(cell->vertex(k).distance (Point()) - 1) < 1e-12, + ExcInternalError()); + } + + deallog << "Volume face: " << face << " with vertices:" << std::endl; + for (unsigned int k=0; k::vertices_per_cell; ++k) + { + deallog << " " << face->vertex(k) << std::endl; + Assert (std::fabs(face->vertex(k).distance (Point()) - 1) < 1e-12, + ExcInternalError()); + } + + + for (unsigned int k=0; k::vertices_per_cell; ++k) + { + Point diff(face->vertex(k)); + diff -= cell->vertex(k); + Assert (diff.square() == 0, ExcInternalError()); + } + } +} + +template +void save_mesh(const Triangulation& tria) +{ + GridOut grid_out; + grid_out.write_gnuplot (tria, deallog.get_file_stream()); +} + + +int main () +{ + + ofstream logfile("extract_boundary_mesh_04/output"); + deallog.attach(logfile); + deallog.depth_console(0); + + + { // Extract the boundary of a hyper-sphere + + const int dim = 3; + deallog << "Testing hyper_cube in dim: " << dim << "..."<< endl; + + map< Triangulation::cell_iterator, + Triangulation::face_iterator> + surface_to_volume_mapping; + const HyperBallBoundary boundary_description; + Triangulation volume_mesh; + GridGenerator::hyper_ball(volume_mesh); + for (Triangulation::active_cell_iterator + cell = volume_mesh.begin_active(); + cell != volume_mesh.end(); ++cell) + for (unsigned int f=0; f::faces_per_cell; ++f) + if (cell->at_boundary(f) && (cell->center()[dim-1]>0)) + cell->face(f)->set_all_boundary_indicators (1); + + volume_mesh.set_boundary (1, boundary_description); + volume_mesh.set_boundary (0, boundary_description); + volume_mesh.refine_global (3); + + const HyperBallBoundary surface_description; + Triangulation boundary_mesh; + boundary_mesh.set_boundary (1, surface_description); + boundary_mesh.set_boundary (0, surface_description); + + set boundary_ids; + boundary_ids.insert(1); + + GridTools::extract_boundary_mesh (volume_mesh, boundary_mesh, + surface_to_volume_mapping, + boundary_ids); + + deallog << volume_mesh.n_active_cells () << std::endl; + deallog << boundary_mesh.n_active_cells () << std::endl; + //save_mesh(boundary_mesh); + + + test_vertices_orientation(boundary_mesh, surface_to_volume_mapping); + //save_mesh(boundary_mesh); + } + + + return 0; + }