From 9f58e5c38ed1004e6135f253670d9d0794d35974 Mon Sep 17 00:00:00 2001
From: Christophe Geuzaine <cgeuzaine@ulg.ac.be>
Date: Thu, 3 Oct 2019 12:52:05 +0200
Subject: [PATCH] When OpenBoundaryPhysicalGroup is given and is an open curve,
 we should not duplicate its boundary nodes. Thanks @bsalves!

---
 CREDITS.txt      |  2 +-
 Plugin/Crack.cpp | 53 ++++++++++++++++++++++++++++++++++++++++--------
 2 files changed, 46 insertions(+), 9 deletions(-)

diff --git a/CREDITS.txt b/CREDITS.txt
index 235eb5e586..eb5d989bb2 100644
--- a/CREDITS.txt
+++ b/CREDITS.txt
@@ -50,7 +50,7 @@ Almeida, Guillaume Demesy, Wendy Merks-Swolfs, Cosmin Stefan Deaconu, Nigel
 Nunn, Serban Georgescu, Julien Troufflard, Michele Mocciola, Matthijs Sypkens
 Smit, Sauli Ruuska, Romain Boman, Fredrik Ekre, Mark Burton, Max Orok, Paul
 Cristini, Isuru Fernando, Jose Paulo Moitinho de Almeida, Sophie Le Bras,
-Alberto Escrig, Samy Mukadi, Peter Johnston.
+Alberto Escrig, Samy Mukadi, Peter Johnston, Bruno de Sousa Alves.
 
 Special thanks to Bill Spitzak, Michael Sweet, Matthias Melcher, Greg Ercolano
 and others for the Fast Light Tool Kit on which Gmsh's GUI is based. See
diff --git a/Plugin/Crack.cpp b/Plugin/Crack.cpp
index e4a20d7326..d169a10e3b 100644
--- a/Plugin/Crack.cpp
+++ b/Plugin/Crack.cpp
@@ -31,15 +31,15 @@ std::string GMSH_CrackPlugin::getHelp() const
   return "Plugin(Crack) creates a crack around the physical "
          "group `PhysicalGroup' of dimension `Dimension' (1 or 2), "
          "embedded in a mesh of dimension `Dimension' + 1. "
-         "The plugin duplicates the vertices and the elements on "
+         "The plugin duplicates the nodes and the elements on "
          "the crack and stores them in a new discrete curve "
          "(`Dimension' = 1) or surface (`Dimension' = 2). The "
          "elements touching the crack on the ``negative'' side "
-         "are modified to use the newly generated vertices."
+         "are modified to use the newly generated nodes."
          "If `OpenBoundaryPhysicalGroup' is given (> 0), its "
-         "vertices are duplicated and the crack will be left "
+         "nodes are duplicated and the crack will be left "
          "open on that (part of the) boundary. Otherwise, the "
-         "lips of the crack are sealed, i.e., its vertices are "
+         "lips of the crack are sealed, i.e., its nodes are "
          "not duplicated. For 1D cracks, `NormalX', `NormalY' and "
          "`NormalZ' provide the reference normal of the surface "
          "in which the crack is supposed to be embedded.";
@@ -115,7 +115,7 @@ PView *GMSH_CrackPlugin::execute(PView *view)
     for(std::size_t j = 0; j < entities[i]->getNumMeshElements(); j++)
       crackElements.push_back(entities[i]->getMeshElement(j));
 
-  // get internal crack vertices and boundary vertices
+  // get internal crack nodes and boundary nodes
   std::set<MVertex *> crackVertices, bndVertices;
   if(dim == 1) {
     for(std::size_t i = 0; i < crackElements.size(); i++) {
@@ -154,13 +154,34 @@ PView *GMSH_CrackPlugin::execute(PView *view)
       bndVertices.insert(it->data.begin(), it->data.end());
   }
 
-  // get (forced) open boundary vertices and remove them from boundary vertices
+  // compute the boundary nodes (if any) of the "OpenBoundary" physical group
+  std::set<MVertex *> bndVerticesFromOpenBoundary;
   for(std::size_t i = 0; i < openEntities.size(); i++) {
     for(std::size_t j = 0; j < openEntities[i]->getNumMeshElements(); j++) {
       MElement *e = openEntities[i]->getMeshElement(j);
       for(std::size_t k = 0; k < e->getNumVertices(); k++) {
         MVertex *v = e->getVertex(k);
-        bndVertices.erase(v);
+        if(bndVerticesFromOpenBoundary.find(v) == bndVerticesFromOpenBoundary.end())
+          bndVerticesFromOpenBoundary.insert(v);
+        else
+          bndVerticesFromOpenBoundary.erase(v);
+      }
+    }
+  }
+
+  if(bndVerticesFromOpenBoundary.size())
+    Msg::Info("%u nodes on boundary of OpenBoundaryPhysicalGroup",
+              bndVerticesFromOpenBoundary.size());
+
+  // get open boundary nodes and remove them from boundary nodes (if they are
+  // not on the "boundary of the open boundary" ;-)
+  for(std::size_t i = 0; i < openEntities.size(); i++) {
+    for(std::size_t j = 0; j < openEntities[i]->getNumMeshElements(); j++) {
+      MElement *e = openEntities[i]->getMeshElement(j);
+      for(std::size_t k = 0; k < e->getNumVertices(); k++) {
+        MVertex *v = e->getVertex(k);
+        if(bndVerticesFromOpenBoundary.find(v) == bndVerticesFromOpenBoundary.end())
+          bndVertices.erase(v);
       }
     }
   }
@@ -216,6 +237,22 @@ PView *GMSH_CrackPlugin::execute(PView *view)
   */
 
   // create new crack entity
+
+  // TODO: the new discrete entities do not have a consistent topology: we don't
+  // specify their bounding points/curves
+  //   a) This is easy to fix if there's no OpenBoundaryPhysicalGroup and
+  //      we crack a *single* elementary entity
+  //   b) If there is an open boundary, we need to create a new elementary
+  //      entity on the boundary, and correctly classify the nodes on it...
+  //      and we also need to create boundary elements
+  //   c) If we crack a group made of multiple elementary entities we might
+  //      want to create multiple crackes entities, and do the same as (b)
+  //      for all internal seams
+  //
+  // In practice, c) is not crucial - the current approach simply creates a
+  // single new surface/curve, which is probably fine as in solvers we won't use
+  // the internal seams.
+
   GEdge *crackEdge = 0;
   GFace *crackFace = 0;
   if(dim == 1) {
@@ -230,7 +267,7 @@ PView *GMSH_CrackPlugin::execute(PView *view)
     crackEdge ? (GEntity *)crackEdge : (GEntity *)crackFace;
   crackEntity->physicals.push_back(physical);
 
-  // duplicate internal crack vertices
+  // duplicate internal crack nodes
   std::map<MVertex *, MVertex *> vxv;
   for(std::set<MVertex *>::iterator it = crackVertices.begin();
       it != crackVertices.end(); it++) {
-- 
GitLab