From: bangerth Date: Tue, 16 Nov 2010 20:34:14 +0000 (+0000) Subject: Minor cleanups. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0308f4f5428e7139fccf2efc233b19964e07d034;p=dealii-svn.git Minor cleanups. git-svn-id: https://svn.dealii.org/trunk@22771 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/source/grid/grid_tools.cc b/deal.II/source/grid/grid_tools.cc index 862fe7a252..57e526cefb 100644 --- a/deal.II/source/grid/grid_tools.cc +++ b/deal.II/source/grid/grid_tools.cc @@ -1813,42 +1813,26 @@ fix_up_distorted_child_cells (const typename Triangulation::Distor template -void +void GridTools::extract_boundary_mesh (const Triangulation &volume_mesh, Triangulation &surface_mesh, std::map::cell_iterator, - typename Triangulation::face_iterator> + typename Triangulation::face_iterator> &surface_to_volume_mapping, const std::set &boundary_ids) { +// Assumption: +// We are relying below on the fact that Triangulation::create_triangulation(...) will keep the order +// pass by CellData and that it will not reorder the vertices. - //

Assumptions

- // We are relying heavily on that - // Triangulation::create_triangulation(...) will keep the order - // pass by CellData and that it will not reorder the vertices. - - //

TODO

- //
    - //
  • std::set &boundary_ids needs to be set - // to a default value of {0}, no is set to empty - //
  • Add mapping from volume to surface - //
  • Add some assertions - //
  • Create a function (reconcile?) to reconcile the triangulations - // after refinement/coarsening is performed in one of them - //
  • Create a class that manages the triangulations, as direct - // access to refining the triangulations will most likely lead - // to inconsitencies - //
- - surface_to_volume_mapping.clear(); // Is this reasonable? + surface_to_volume_mapping.clear(); const unsigned int boundary_dim = dim-1; //dimension of the boundary mesh - // First create surface mesh and mapping from only level(0) cells of volume_mesh - - std::vector::face_iterator> - mapping; // temporary map for level==0 + // First create surface mesh and mapping from only level(0) cells of volume_mesh + std::vector::face_iterator> + mapping; // temporary map for level==0 typename Triangulation::cell_iterator @@ -1861,87 +1845,81 @@ GridTools::extract_boundary_mesh (const Triangulation &volume_mesh std::vector< Point > vertices; std::map map_vert_index; //volume vertex indices to surf ones - + unsigned int v_index; CellData< boundary_dim > c_data; for (; cell!=endc; ++cell) - - for (unsigned int i=0; - i < GeometryInfo::faces_per_cell; ++i){ - - face = cell->face(i); - - if ( face->at_boundary() - && - (boundary_ids.empty() || - ( boundary_ids.find(face->boundary_indicator()) != boundary_ids.end())) - ) { - for (unsigned int j=0; - j::vertices_per_cell; ++j) { - - v_index = face->vertex_index(j); - - if ( !touched[v_index] ){ - vertices.push_back(face->vertex(j)); - map_vert_index[v_index] = vertices.size() - 1; - touched[v_index] = true; - } - - c_data.vertices[j] = map_vert_index[v_index]; - c_data.material_id = 0; - - } - - cells.push_back(c_data); - mapping.push_back(face); - - } - } + for (unsigned int i=0; i < GeometryInfo::faces_per_cell; ++i) + { + face = cell->face(i); - // create level 0 surface triangulation - surface_mesh.create_triangulation (vertices, cells, SubCellData()); + if ( face->at_boundary() + && + (boundary_ids.empty() || + ( boundary_ids.find(face->boundary_indicator()) != boundary_ids.end())) ) + { + for (unsigned int j=0; + j::vertices_per_cell; ++j) + { + v_index = face->vertex_index(j); - { // Make the actual mapping - typename Triangulation::active_cell_iterator - cell = surface_mesh.begin(0), - endc = surface_mesh.end(0); - for (; cell!=endc; ++cell) - surface_to_volume_mapping[cell] = mapping.at(cell->index()); - } + if ( !touched[v_index] ) + { + vertices.push_back(face->vertex(j)); + map_vert_index[v_index] = vertices.size() - 1; + touched[v_index] = true; + } + + c_data.vertices[j] = map_vert_index[v_index]; + c_data.material_id = 0; + } - do { - bool changed = false; - typename Triangulation::active_cell_iterator - cell = surface_mesh.begin_active(), - endc = surface_mesh.end(); - - for (; cell!=endc; ++cell) - if (surface_to_volume_mapping[cell]->has_children() == true ) { - cell->set_refine_flag (); - changed = true; + cells.push_back(c_data); + mapping.push_back(face); + } } - + // create level 0 surface triangulation + surface_mesh.create_triangulation (vertices, cells, SubCellData()); + + // Make the actual mapping + for (typename Triangulation::active_cell_iterator + cell = surface_mesh.begin(0); + cell!=surface_mesh.end(0); ++cell) + surface_to_volume_mapping[cell] = mapping.at(cell->index()); - if (changed){ - - surface_mesh.execute_coarsening_and_refinement (); + do + { + bool changed = false; + typename Triangulation::active_cell_iterator + cell = surface_mesh.begin_active(), + endc = surface_mesh.end(); - typename Triangulation::cell_iterator - cell = surface_mesh.begin(), - endc = surface_mesh.end(); for (; cell!=endc; ++cell) - for (unsigned int c=0; cn_children(); c++) - if (surface_to_volume_mapping.find(cell->child(c)) == surface_to_volume_mapping.end()) - surface_to_volume_mapping[cell->child(c)] - = surface_to_volume_mapping[cell]->child(c); + if (surface_to_volume_mapping[cell]->has_children() == true ) + { + cell->set_refine_flag (); + changed = true; + } + + if (changed) + { + surface_mesh.execute_coarsening_and_refinement (); + + typename Triangulation::cell_iterator + cell = surface_mesh.begin(), + endc = surface_mesh.end(); + for (; cell!=endc; ++cell) + for (unsigned int c=0; cn_children(); c++) + if (surface_to_volume_mapping.find(cell->child(c)) == surface_to_volume_mapping.end()) + surface_to_volume_mapping[cell->child(c)] + = surface_to_volume_mapping[cell]->child(c); + } + else + break; } - else - break; - } while (true); - } // explicit instantiations