From 2265620119cf04fcab63bbc190f55678187c6f05 Mon Sep 17 00:00:00 2001
From: Emilie Marchandise <emilie.marchandise@uclouvain.be>
Date: Wed, 25 Apr 2012 09:44:14 +0000
Subject: [PATCH] centerline field detects now inlet and outlets separately

---
 Geo/GFaceCompound.cpp    |  4 ++--
 Mesh/CenterlineField.cpp | 29 +++++++++++++++++++++++------
 Mesh/CenterlineField.h   |  3 +++
 Mesh/meshGFace.cpp       |  7 +++----
 Mesh/meshGFace.h         |  2 +-
 5 files changed, 32 insertions(+), 13 deletions(-)

diff --git a/Geo/GFaceCompound.cpp b/Geo/GFaceCompound.cpp
index dbba693362..6ff9949aae 100644
--- a/Geo/GFaceCompound.cpp
+++ b/Geo/GFaceCompound.cpp
@@ -511,7 +511,7 @@ void GFaceCompound::printFillTris() const{
 void GFaceCompound::fillNeumannBCS_Plane() const
 {
 
-  Msg::Info("Meshing %d interior holes with planes ", _interior_loops.size()-1);
+  Msg::Debug("Meshing %d interior holes with planes ", _interior_loops.size()-1);
 
   fillTris.clear();
   fillNodes.clear();
@@ -542,7 +542,7 @@ void GFaceCompound::fillNeumannBCS_Plane() const
       opt_mesh_algo2d(0, GMSH_SET, 1.0); //mesh adapt
       opt_mesh_recombine_all(0, GMSH_SET, 0.0); //no recombination
       meshGFace mgf;
-      mgf(newFace);
+      mgf(newFace,false);
       opt_mesh_algo2d(0, GMSH_SET, meshingAlgo);
       opt_mesh_recombine_all(0, GMSH_SET, recombine);
       for(unsigned int i = 0; i < newFace->triangles.size(); ++i){
diff --git a/Mesh/CenterlineField.cpp b/Mesh/CenterlineField.cpp
index 5897f27cd0..8eb07acf3e 100644
--- a/Mesh/CenterlineField.cpp
+++ b/Mesh/CenterlineField.cpp
@@ -405,7 +405,9 @@ void Centerline::importFile(std::string fileName){
 
   int maxN = 0.0;
   std::vector<GEdge*> modEdges = mod->bindingsGetEdges();
- for (unsigned int i = 0; i < modEdges.size(); i++){
+  MVertex *vin = modEdges[0]->lines[0]->getVertex(0);
+  ptin = SPoint3(vin->x(), vin->y(), vin->z());
+  for (unsigned int i = 0; i < modEdges.size(); i++){
     GEdge *ge = modEdges[i];
     for(unsigned int j = 0; j < ge->lines.size(); j++){
       MLine *l = ge->lines[j];
@@ -660,7 +662,7 @@ void Centerline::createSplitCompounds(){
 
     gfc->meshAttributes.recombine = recombine;
     gfc->addPhysicalEntity(100);
-    current->setPhysicalName("newsurf", 2, 100);
+    current->setPhysicalName("wall", 2, 100);
 
   }
 
@@ -798,6 +800,8 @@ void Centerline::createClosedVolume(){
 
   //identify the boundary edges by looping over all discreteFaces
   std::vector<GEdge*> boundEdges;
+  double dist_inlet = 1.e6;
+  GEdge *gin = NULL;
   for (int i= 0; i< NF; i++){
     GFace *gf = current->getFaceByTag(i+1);
     std::list<GEdge*> l_edges = gf->edges();
@@ -805,9 +809,16 @@ void Centerline::createClosedVolume(){
       std::vector<GEdge*>::iterator ite = std::find(boundEdges.begin(), boundEdges.end(), *it);
       if (ite != boundEdges.end()) boundEdges.erase(ite);
       else boundEdges.push_back(*it);
+      GVertex *gv = (*it)->getBeginVertex();
+      SPoint3 pt(gv->x(), gv->y(), gv->z());
+      double dist = pt.distance(ptin);
+      if(dist < dist_inlet){
+	dist_inlet = dist;
+	gin = *it;
+      }			      
     }
   }
-
+ 
   current->setFactory("Gmsh");
   std::vector<std::vector<GFace *> > myFaceLoops;
   std::vector<GFace *> myFaces;
@@ -820,8 +831,14 @@ void Centerline::createClosedVolume(){
     myEdges.push_back(gec);
     myEdgeLoops.push_back(myEdges);
     GFace *newFace = current->addPlanarFace(myEdgeLoops);
-    newFace->addPhysicalEntity(200);
-    current->setPhysicalName("in/out", 2, 200);
+    if (gin==boundEdges[i]) {
+      newFace->addPhysicalEntity(200);
+      current->setPhysicalName("inlet", 2, 200);
+    }
+    else{
+      newFace->addPhysicalEntity(300);
+      current->setPhysicalName("outlets", 2, 300);
+    }
     myFaces.push_back(newFace);
   }
 
@@ -834,7 +851,7 @@ void Centerline::createClosedVolume(){
   myFaceLoops.push_back(myFaces);
   GRegion *reg = current->addVolume(myFaceLoops);
   reg->addPhysicalEntity(reg->tag());
-  current->setPhysicalName("newvol", 3, reg->tag());
+  current->setPhysicalName("lumenVolume", 3, reg->tag());
 
   Msg::Info("Centerline action (closeVolume) has created volume %d ", reg->tag());
 
diff --git a/Mesh/CenterlineField.h b/Mesh/CenterlineField.h
index 7fd5ff5ee2..532d4749a1 100644
--- a/Mesh/CenterlineField.h
+++ b/Mesh/CenterlineField.h
@@ -25,6 +25,7 @@ class MTriangle;
 class discreteEdge;
 class discreteFace;
 class MElement;
+class SPoint3;
 
 // A branch of a 1D tree
 struct Branch{
@@ -67,6 +68,8 @@ class Centerline : public Field{
   double hLayer;
   int nbElemLayer;
 
+  //inlet point
+  SPoint3 ptin;
   //all (unique) lines of centerlines
   std::vector<MLine*> lines;
   //the stuctured tree of the centerlines
diff --git a/Mesh/meshGFace.cpp b/Mesh/meshGFace.cpp
index 99bf99e2ab..d5a5c6b885 100644
--- a/Mesh/meshGFace.cpp
+++ b/Mesh/meshGFace.cpp
@@ -1790,7 +1790,7 @@ void deMeshGFace::operator() (GFace *gf)
 // for debugging, change value from -1 to -100;
 int debugSurface = -1; //-1;
 
-void meshGFace::operator() (GFace *gf)
+void meshGFace::operator() (GFace *gf, bool print)
 {
 
   gf->model()->setCurrentMeshEntity(gf);
@@ -1846,7 +1846,8 @@ void meshGFace::operator() (GFace *gf)
     algo = "MeshAdapt";
   }
 
-  Msg::Info("Meshing surface %d (%s, %s)", gf->tag(), gf->getTypeString().c_str(), algo);
+  if (print)
+    Msg::Info("Meshing surface %d (%s, %s)", gf->tag(), gf->getTypeString().c_str(), algo);
 
   // compute loops on the fly (indices indicate start and end points
   // of a loop; loops are not yet oriented)
@@ -2127,8 +2128,6 @@ void orientMeshGFace::operator()(GFace *gf)
   //do sthg for compound face
   if(gf->geomType() == GEntity::CompoundSurface ) {
     GFaceCompound *gfc = (GFaceCompound*) gf;
-    //if (gfc->getCompounds().size() != 1) return;
-    //else{printf("compound face %d orient \n", gfc->tag());}
   
     std::list<GFace*> comp = gfc->getCompounds(); 
     MTriangle *lt = (*comp.begin())->triangles[0];
diff --git a/Mesh/meshGFace.h b/Mesh/meshGFace.h
index dd31a952c6..e939d2d6d5 100644
--- a/Mesh/meshGFace.h
+++ b/Mesh/meshGFace.h
@@ -25,7 +25,7 @@ class meshGFace {
     : repairSelfIntersecting1dMesh(r), twoPassesMesh(t), onlyInitialMesh(false)
   {
   }
-  void operator()(GFace *);
+  void operator()(GFace *, bool print=true);
   void setOnlyInitial(){ onlyInitialMesh = true; }
 };
 
-- 
GitLab