From 98bacfe6882089bf6f6fa4beab167b930c649175 Mon Sep 17 00:00:00 2001 From: Sean Ingimarson Date: Mon, 14 Feb 2022 11:44:54 -0500 Subject: [PATCH] Allow loading Gmsh files with some cells inverted doc fix --- include/deal.II/grid/grid_tools.h | 15 ++++ source/grid/grid_in.cc | 2 +- source/grid/grid_tools.cc | 18 +++- source/grid/grid_tools.inst.in | 5 ++ tests/grid/grid_in_gmsh_02.cc | 86 ++++++++++++++++++ tests/grid/grid_in_gmsh_02.geo | 26 ++++++ tests/grid/grid_in_gmsh_02.msh | 141 ++++++++++++++++++++++++++++++ tests/grid/grid_in_gmsh_02.output | 7 ++ 8 files changed, 296 insertions(+), 4 deletions(-) create mode 100644 tests/grid/grid_in_gmsh_02.cc create mode 100644 tests/grid/grid_in_gmsh_02.geo create mode 100644 tests/grid/grid_in_gmsh_02.msh create mode 100644 tests/grid/grid_in_gmsh_02.output diff --git a/include/deal.II/grid/grid_tools.h b/include/deal.II/grid/grid_tools.h index 5f5a5b9b91..d02bd22642 100644 --- a/include/deal.II/grid/grid_tools.h +++ b/include/deal.II/grid/grid_tools.h @@ -454,6 +454,21 @@ namespace GridTools const std::vector> &all_vertices, std::vector> & cells); + /** + * Check the given cells and inverts any cell that is considered to have + * negative measure/volume in the orientation required by deal.II. + * + * This function is identical to invert_all_negative_measure_cells() except it + * does not throw an error if only some of the cells are inverted. Instead, + * this function returns how many cells were inverted. Additionally, it will + * always throw an exception outside of codimension 0. + */ + template + std::size_t + invert_cells_with_negative_measure( + const std::vector> &all_vertices, + std::vector> & cells); + /** * Given a vector of CellData objects describing a mesh, reorder their * vertices so that all lines are consistently oriented. diff --git a/source/grid/grid_in.cc b/source/grid/grid_in.cc index afa6c45ac4..69203c35f4 100644 --- a/source/grid/grid_in.cc +++ b/source/grid/grid_in.cc @@ -2688,7 +2688,7 @@ GridIn::read_msh(std::istream &in) GridTools::delete_unused_vertices(vertices, cells, subcelldata); // ... and cells if (dim == spacedim) - GridTools::invert_all_negative_measure_cells(vertices, cells); + GridTools::invert_cells_with_negative_measure(vertices, cells); GridTools::consistently_order_cells(cells); } tria->create_triangulation(vertices, cells, subcelldata); diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index 9cd09727d8..8f6e7d832d 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -859,13 +859,13 @@ namespace GridTools template - void - invert_all_negative_measure_cells( + std::size_t + invert_cells_with_negative_measure( const std::vector> &all_vertices, std::vector> & cells) { if (dim == 1) - return; + return 0; if (dim == 2 && spacedim == 3) Assert(false, ExcNotImplemented()); @@ -902,6 +902,18 @@ namespace GridTools ExcInternalError()); } } + return n_negative_cells; + } + + + template + void + invert_all_negative_measure_cells( + const std::vector> &all_vertices, + std::vector> & cells) + { + const std::size_t n_negative_cells = + invert_cells_with_negative_measure(all_vertices, cells); // We assume that all cells of a grid have // either positive or negative volumes but diff --git a/source/grid/grid_tools.inst.in b/source/grid/grid_tools.inst.in index 7b01f27fe7..1ec72b814a 100644 --- a/source/grid/grid_tools.inst.in +++ b/source/grid/grid_tools.inst.in @@ -256,6 +256,11 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS) const std::vector> &, std::vector> &); + template std::size_t + invert_cells_with_negative_measure( + const std::vector> &, + std::vector> &); + # if deal_II_dimension == deal_II_space_dimension template void consistently_order_cells(std::vector> &); diff --git a/tests/grid/grid_in_gmsh_02.cc b/tests/grid/grid_in_gmsh_02.cc new file mode 100644 index 0000000000..333e71d154 --- /dev/null +++ b/tests/grid/grid_in_gmsh_02.cc @@ -0,0 +1,86 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2016 - 2022 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.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +// Test loading a .msh file from gmsh with some cells having negative volume + +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include + +using namespace dealii; + + +template +void +print_mesh_info(const Triangulation &triangulation, + const std::string & filename) +{ + deallog << "Mesh info:" << std::endl + << " dimension: " << dim << std::endl + << " no. of cells: " << triangulation.n_active_cells() << std::endl; + + { + std::map boundary_count; + for (const auto &face : triangulation.active_face_iterators()) + if (face->at_boundary()) + boundary_count[face->boundary_id()]++; + + deallog << " boundary indicators: "; + for (const std::pair &pair : + boundary_count) + { + deallog << pair.first << "(" << pair.second << " times) "; + } + deallog << std::endl; + } + + // Finally, produce a graphical representation of the mesh to an output + // file: + std::ofstream out(filename); + GridOut grid_out; + grid_out.write_vtu(triangulation, out); + deallog << " written to " << filename << std::endl << std::endl; +} + + +void +grid() +{ + Triangulation<2> triangulation; + + GridIn<2> gridin; + gridin.attach_triangulation(triangulation); + std::ifstream f(SOURCE_DIR "/grid_in_gmsh_02.msh"); + gridin.read_msh(f); + + print_mesh_info(triangulation, "grid.vtu"); +} + + +int +main() +{ + initlog(); + grid(); +} diff --git a/tests/grid/grid_in_gmsh_02.geo b/tests/grid/grid_in_gmsh_02.geo new file mode 100644 index 0000000000..f642887d3a --- /dev/null +++ b/tests/grid/grid_in_gmsh_02.geo @@ -0,0 +1,26 @@ +meshsize = 2; + +Point(1) = {-1, 0, 0, meshsize}; +Point(2) = {1, 0, 0, meshsize}; +Point(3) = {1, 2, 0, meshsize}; +Point(4) = {-1, 2, 0, meshsize}; +Point(5) = {-1, 4, 0, meshsize}; +Point(6) = {1, 4, 0, meshsize}; +Line(1) = {5, 6}; +Line(2) = {4, 3}; +Line(3) = {1, 2}; +Line(4) = {5, 4}; +Line(5) = {4, 1}; +Line(6) = {6, 3}; +Line(7) = {3, 2}; + +Line Loop(9) = {1, 6, -2, -4}; +Plane Surface(9) = {9}; +Line Loop(11) = {3, -7, -2, 5}; +Plane Surface(11) = {11}; + +Physical Curve("1") = {1}; +Physical Curve("2") = {3}; +Physical Surface(1) = {9, 11}; + +Mesh.SubdivisionAlgorithm = 1; diff --git a/tests/grid/grid_in_gmsh_02.msh b/tests/grid/grid_in_gmsh_02.msh new file mode 100644 index 0000000000..da1df887fc --- /dev/null +++ b/tests/grid/grid_in_gmsh_02.msh @@ -0,0 +1,141 @@ +$MeshFormat +4.1 0 8 +$EndMeshFormat +$PhysicalNames +2 +1 1 "1" +1 2 "2" +$EndPhysicalNames +$Entities +6 7 2 0 +1 -1 0 0 0 +2 1 0 0 0 +3 1 2 0 0 +4 -1 2 0 0 +5 -1 4 0 0 +6 1 4 0 0 +1 -1 4 0 1 4 0 1 1 2 5 -6 +2 -1 2 0 1 2 0 0 2 4 -3 +3 -1 0 0 1 0 0 1 2 2 1 -2 +4 -1 2 0 -1 4 0 0 2 5 -4 +5 -1 0 0 -1 2 0 0 2 4 -1 +6 1 2 0 1 4 0 0 2 6 -3 +7 1 0 0 1 2 0 0 2 3 -2 +9 -1 2 0 1 4 0 1 1 4 1 6 -2 -4 +11 -1 0 0 1 2 0 1 1 4 3 -7 -2 5 +$EndEntities +$Nodes +15 31 1 31 +0 1 0 1 +1 +-1 0 0 +0 2 0 1 +2 +1 0 0 +0 3 0 1 +3 +1 2 0 +0 4 0 1 +4 +-1 2 0 +0 5 0 1 +5 +-1 4 0 +0 6 0 1 +6 +1 4 0 +1 1 0 1 +7 +-2.625677453238495e-12 4 0 +1 2 0 1 +8 +-2.625677453238495e-12 2 0 +1 3 0 1 +9 +-2.625677453238495e-12 0 0 +1 4 0 1 +10 +-1 3.000000000002758 0 +1 5 0 1 +11 +-1 1.000000000002661 0 +1 6 0 1 +12 +1 3.000000000002758 0 +1 7 0 1 +13 +1 1.000000000002661 0 +2 9 0 9 +14 +15 +16 +17 +18 +19 +20 +21 +22 +0 3 0 +-0.5 2.5 0 +0.5 2.5 0 +0.5 3.5 0 +-0.5 3.5 0 +-4.376082829438133e-13 2.333333333333333 0 +0.6666666666666667 3.00000000000046 0 +-0.6666666666666667 3.00000000000046 0 +-4.376082829438133e-13 3.666666666666667 0 +2 11 0 9 +23 +24 +25 +26 +27 +28 +29 +30 +31 +0 1 0 +0.5 0.5 0 +-0.5 0.5 0 +-0.5 1.5 0 +0.5 1.5 0 +-4.376082829438133e-13 0.3333333333333333 0 +-0.6666666666666667 1.000000000000444 0 +0.6666666666666667 1.000000000000443 0 +-4.376082829438133e-13 1.666666666666667 0 +$EndNodes +$Elements +4 28 1 28 +1 1 1 2 +1 5 7 +2 7 6 +1 3 1 2 +3 1 9 +4 9 2 +2 9 3 12 +5 3 8 19 16 +6 8 4 15 19 +7 16 19 15 14 +8 6 12 20 17 +9 12 3 16 20 +10 17 20 16 14 +11 4 10 21 15 +12 10 5 18 21 +13 15 21 18 14 +14 5 7 22 18 +15 7 6 17 22 +16 18 22 17 14 +2 11 3 12 +17 1 9 28 25 +18 9 2 24 28 +19 25 28 24 23 +20 4 11 29 26 +21 11 1 25 29 +22 26 29 25 23 +23 2 13 30 24 +24 13 3 27 30 +25 24 30 27 23 +26 3 8 31 27 +27 8 4 26 31 +28 27 31 26 23 +$EndElements diff --git a/tests/grid/grid_in_gmsh_02.output b/tests/grid/grid_in_gmsh_02.output new file mode 100644 index 0000000000..724c6b60cf --- /dev/null +++ b/tests/grid/grid_in_gmsh_02.output @@ -0,0 +1,7 @@ + +DEAL::Mesh info: +DEAL:: dimension: 2 +DEAL:: no. of cells: 24 +DEAL:: boundary indicators: 0(8 times) 1(2 times) 2(2 times) +DEAL:: written to grid.vtu +DEAL:: -- 2.39.5