]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Allow loading Gmsh files with some cells inverted
authorSean Ingimarson <singima@clemson.edu>
Mon, 14 Feb 2022 16:44:54 +0000 (11:44 -0500)
committerSean Ingimarson <singima@clemson.edu>
Fri, 25 Feb 2022 18:19:34 +0000 (13:19 -0500)
doc fix

include/deal.II/grid/grid_tools.h
source/grid/grid_in.cc
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/grid/grid_in_gmsh_02.cc [new file with mode: 0644]
tests/grid/grid_in_gmsh_02.geo [new file with mode: 0644]
tests/grid/grid_in_gmsh_02.msh [new file with mode: 0644]
tests/grid/grid_in_gmsh_02.output [new file with mode: 0644]

index 5f5a5b9b9174daade7256dcabf2ef7bef11acfd8..d02bd22642d1a1474f4243188c4ffb53108ac077 100644 (file)
@@ -454,6 +454,21 @@ namespace GridTools
     const std::vector<Point<spacedim>> &all_vertices,
     std::vector<CellData<dim>> &        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 <int dim, int spacedim>
+  std::size_t
+  invert_cells_with_negative_measure(
+    const std::vector<Point<spacedim>> &all_vertices,
+    std::vector<CellData<dim>> &        cells);
+
   /**
    * Given a vector of CellData objects describing a mesh, reorder their
    * vertices so that all lines are consistently oriented.
index afa6c45ac4fcd1f9c6539f19a34e1f53edb74261..69203c35f4a68b0ac75f052df7a2064ba90b62ef 100644 (file)
@@ -2688,7 +2688,7 @@ GridIn<dim, spacedim>::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);
index 9cd09727d8d16cee168ca8b868789461a0226806..8f6e7d832dfb778b89c1a686a71a467ee1061c7f 100644 (file)
@@ -859,13 +859,13 @@ namespace GridTools
 
 
   template <int dim, int spacedim>
-  void
-  invert_all_negative_measure_cells(
+  std::size_t
+  invert_cells_with_negative_measure(
     const std::vector<Point<spacedim>> &all_vertices,
     std::vector<CellData<dim>> &        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 <int dim, int spacedim>
+  void
+  invert_all_negative_measure_cells(
+    const std::vector<Point<spacedim>> &all_vertices,
+    std::vector<CellData<dim>> &        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
index 7b01f27fe78a7571febb2f6d39312b9c112cd326..1ec72b814a07b400e6532979e5de3dfac80f57a4 100644 (file)
@@ -256,6 +256,11 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
         const std::vector<Point<deal_II_space_dimension>> &,
         std::vector<CellData<deal_II_dimension>> &);
 
+      template std::size_t
+      invert_cells_with_negative_measure(
+        const std::vector<Point<deal_II_space_dimension>> &,
+        std::vector<CellData<deal_II_dimension>> &);
+
 #  if deal_II_dimension == deal_II_space_dimension
       template void
       consistently_order_cells(std::vector<CellData<deal_II_dimension>> &);
diff --git a/tests/grid/grid_in_gmsh_02.cc b/tests/grid/grid_in_gmsh_02.cc
new file mode 100644 (file)
index 0000000..333e71d
--- /dev/null
@@ -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 <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_in.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/manifold_lib.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include <fstream>
+#include <iostream>
+#include <map>
+
+using namespace dealii;
+
+
+template <int dim>
+void
+print_mesh_info(const Triangulation<dim> &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<types::boundary_id, unsigned int> 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<const types::boundary_id, unsigned int> &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 (file)
index 0000000..f642887
--- /dev/null
@@ -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 (file)
index 0000000..da1df88
--- /dev/null
@@ -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 (file)
index 0000000..724c6b6
--- /dev/null
@@ -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::

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.