From d38b6ffa0b86e0ea252ba752aebdf19551920062 Mon Sep 17 00:00:00 2001
From: Matti Pellika <matti.pellikka@tut.fi>
Date: Fri, 1 Jan 2010 17:43:34 +0000
Subject: [PATCH] Changed homology generator chain mesh element orientations
 (determined by vertex order) so that post-processing data isn't needed to
 save the generators, instead all the elements have coefficient 1 now in the
 chain. Due to this also reverted changes made to GModel writeMSH method in my
 earlier commit (they became obsolete).

---
 Geo/ChainComplex.cpp       | 11 +++++---
 Geo/GModel.h               |  3 ++-
 Geo/GModelIO_Mesh.cpp      | 54 ++++++++++++++++++--------------------
 Geo/Homology.h             |  5 ++--
 Post/PViewDataGModelIO.cpp |  2 +-
 tutorial/t10.geo           |  8 +++---
 6 files changed, 45 insertions(+), 38 deletions(-)

diff --git a/Geo/ChainComplex.cpp b/Geo/ChainComplex.cpp
index 1e32afa678..d9a522cb7d 100644
--- a/Geo/ChainComplex.cpp
+++ b/Geo/ChainComplex.cpp
@@ -871,18 +871,22 @@ void Chain::createPView(){
   
   std::vector<MElement*> elements;
   MElementFactory factory;
-
   std::map<int, std::vector<double> > data;
   
   for(citer cit = _cells.begin(); cit != _cells.end(); cit++){
     Cell* cell = (*cit).first;
     int coeff = (*cit).second;
     std::vector<MVertex*> v = cell->getVertexVector();
+    if(cell->getDim() > 0 && coeff == -1){ // flip orientation
+      MVertex* temp = v[0];
+      v[0] = v[1];
+      v[1] = temp;
+    }
     MElement *e = factory.create(cell->getTypeForMSH(), v, cell->getNum(), cell->getPartition());
     elements.push_back(e);
     
-    std::vector<double> coeffs (1,coeff);
-    data.insert(std::make_pair(e->getNum(), coeffs));
+    std::vector<double> coeffs (1,1);
+    data[e->getNum()] = coeffs;
   }
   
   int max[4];
@@ -903,6 +907,7 @@ void Chain::createPView(){
   _model->storeChain(getDim(), entityMap, physicalMap);
   _model->setPhysicalName(getName(), getDim(), physicalNum);
   
+  // only for visualization
   if(!data.empty()) PView *chain = new PView(getName(), "ElementData", getGModel(), data, 0, 1);
   
   return;
diff --git a/Geo/GModel.h b/Geo/GModel.h
index 2346ff5c4c..a71531f51e 100644
--- a/Geo/GModel.h
+++ b/Geo/GModel.h
@@ -371,7 +371,7 @@ class GModel
   // Gmsh mesh file format
   int readMSH(const std::string &name);
   int writeMSH(const std::string &name, double version=1.0, bool binary=false,
-               bool saveAll=false, bool saveParametric=false, double scalingFactor=1.0, bool eleRenumbering=true);
+               bool saveAll=false, bool saveParametric=false, double scalingFactor=1.0);
   int writeDistanceMSH(const std::string &name, double version=1.0, bool binary=false,
 		       bool saveAll=false, bool saveParametric=false, double scalingFactor=1.0);
 
@@ -448,6 +448,7 @@ class GModel
                  std::map<int, std::map<int, std::string> > &physicalMap){
     _storeElementsInEntities(entityMap);
     _storePhysicalTagsInEntities(dim, physicalMap);
+    _associateEntityWithMeshVertices();
   }
 };
 
diff --git a/Geo/GModelIO_Mesh.cpp b/Geo/GModelIO_Mesh.cpp
index edff2d5d6d..4c14ba50ae 100644
--- a/Geo/GModelIO_Mesh.cpp
+++ b/Geo/GModelIO_Mesh.cpp
@@ -529,16 +529,14 @@ static void writeElementHeaderMSH(bool binary, FILE *fp, std::map<int,int> &elem
 template<class T>
 static void writeElementsMSH(FILE *fp, T *ele, bool saveAll, 
                              double version, bool binary, int &num, int elementary, 
-                             std::vector<int> &physicals, int parentNum, bool eleRenumbering=true)
+                             std::vector<int> &physicals, int parentNum)
 {
   if(saveAll){
-    if(eleRenumbering) ele->writeMSH(fp, version, binary, ++num, elementary, 0, parentNum);
-    else ele->writeMSH(fp, version, binary, 0, elementary, 0, parentNum);
+    ele->writeMSH(fp, version, binary, ++num, elementary, 0, parentNum);
   }
   else{
     for(unsigned int j = 0; j < physicals.size(); j++){
-      if(eleRenumbering) ele->writeMSH(fp, version, binary, ++num, elementary, physicals[j], parentNum);
-      else ele->writeMSH(fp, version, binary, 0, elementary, physicals[j], parentNum);
+      ele->writeMSH(fp, version, binary, ++num, elementary, physicals[j], parentNum);
     }
   }
 }
@@ -546,24 +544,24 @@ static void writeElementsMSH(FILE *fp, T *ele, bool saveAll,
 template<class T>
 static void writeElementsMSH(FILE *fp, std::vector<T*> &ele, bool saveAll,
                              double version, bool binary, int &num, int elementary,
-                             std::vector<int> &physicals, bool eleRenumbering=true)
+                             std::vector<int> &physicals)
 {
   for(unsigned int i = 0; i < ele.size(); i++)
     writeElementsMSH(fp, ele[i], saveAll, version, binary, num,
-                     elementary, physicals, 0, eleRenumbering);
+                     elementary, physicals, 0);
 }
 
 template<class T>
 static void writeElementsMSH(FILE *fp, std::vector<T*> &ele, bool saveAll,
                              double version, bool binary, int &num, int elementary,
-                             std::vector<int> &physicals, std::map<MElement *, int> &parentsNum, bool eleRenumbering=true)
+                             std::vector<int> &physicals, std::map<MElement *, int> &parentsNum)
 {
   for(unsigned int i = 0; i < ele.size(); i++) {
     int parentNum = (ele[i]->getParent() != NULL &&
                      parentsNum.find(ele[i]->getParent()) != parentsNum.end()) ?
       parentsNum.find(ele[i]->getParent())->second : 0;
     writeElementsMSH(fp, ele[i], saveAll, version, binary, num,
-                     elementary, physicals, parentNum, eleRenumbering);
+                     elementary, physicals, parentNum);
   }
 }
 int GModel::writeDistanceMSH(const std::string &name, double version, bool binary,
@@ -658,7 +656,7 @@ int GModel::writeDistanceMSH(const std::string &name, double version, bool binar
 
 }
 int GModel::writeMSH(const std::string &name, double version, bool binary,
-                     bool saveAll, bool saveParametric, double scalingFactor, bool eleRenumbering)
+                     bool saveAll, bool saveParametric, double scalingFactor)
 {
   FILE *fp = fopen(name.c_str(), binary ? "wb" : "w");
   if(!fp){
@@ -794,99 +792,99 @@ int GModel::writeMSH(const std::string &name, double version, bool binary,
   writeElementHeaderMSH(binary, fp, elements, MSH_PNT);
   for(viter it = firstVertex(); it != lastVertex(); ++it)
     writeElementsMSH(fp, (*it)->points, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, eleRenumbering);
+                     (*it)->tag(), (*it)->physicals);
 
   writeElementHeaderMSH(binary, fp, elements, MSH_LIN_2, MSH_LIN_3, MSH_LIN_4,
                         MSH_LIN_5);
   for(eiter it = firstEdge(); it != lastEdge(); ++it)
    writeElementsMSH(fp, (*it)->lines, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, eleRenumbering);
+                     (*it)->tag(), (*it)->physicals);
 
   writeElementHeaderMSH(binary, fp, elements, MSH_TRI_3, MSH_TRI_6, MSH_TRI_9,
                         MSH_TRI_10, MSH_TRI_12, MSH_TRI_15, MSH_TRI_15I, MSH_TRI_21);
   for(itP = parents[0].begin(); itP != parents[0].end(); itP++)
     if(itP->first->getType() == TYPE_TRI) {
       writeElementsMSH(fp, itP->first, saveAll, version, binary, num,
-                       itP->second->tag(), itP->second->physicals, 0, eleRenumbering);
+                       itP->second->tag(), itP->second->physicals, 0);
       parentsNum[itP->first] = num;
     }
   for(fiter it = firstFace(); it != lastFace(); ++it)
     writeElementsMSH(fp, (*it)->triangles, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, eleRenumbering);
+                     (*it)->tag(), (*it)->physicals);
   writeElementHeaderMSH(binary, fp, elements, MSH_QUA_4, MSH_QUA_9, MSH_QUA_8, 
                         MSH_QUA_16, MSH_QUA_25, MSH_QUA_36, MSH_QUA_12, MSH_QUA_16I,
                         MSH_QUA_20);
   for(itP = parents[0].begin(); itP != parents[0].end(); itP++)
     if(itP->first->getType() == TYPE_QUA) {
       writeElementsMSH(fp, itP->first, saveAll, version, binary, num,
-                       itP->second->tag(), itP->second->physicals, 0, eleRenumbering);
+                       itP->second->tag(), itP->second->physicals, 0);
       parentsNum[itP->first] = num;
     }
   for(fiter it = firstFace(); it != lastFace(); ++it)
     writeElementsMSH(fp, (*it)->quadrangles, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, eleRenumbering);
+                     (*it)->tag(), (*it)->physicals);
   writeElementHeaderMSH(binary, fp, elements, MSH_POLYG_);
   for(itP = parents[0].begin(); itP != parents[0].end(); itP++)
     if(itP->first->getType() == TYPE_POLYG) {
       writeElementsMSH(fp, itP->first, saveAll, version, binary, num,
-                       itP->second->tag(), itP->second->physicals, 0, eleRenumbering);
+                       itP->second->tag(), itP->second->physicals, 0);
       parentsNum[itP->first] = num;
     }
   for(fiter it = firstFace(); it != lastFace(); it++)
     writeElementsMSH(fp, (*it)->polygons, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, parentsNum, eleRenumbering);
+                     (*it)->tag(), (*it)->physicals, parentsNum);
 
   writeElementHeaderMSH(binary, fp, elements, MSH_TET_4, MSH_TET_10, MSH_TET_20, 
                         MSH_TET_35, MSH_TET_56, MSH_TET_52);
   for(itP = parents[1].begin(); itP != parents[1].end(); itP++)
     if(itP->first->getType() == TYPE_TET) {
       writeElementsMSH(fp, itP->first, saveAll, version, binary, num,
-                       itP->second->tag(), itP->second->physicals, 0, eleRenumbering);
+                       itP->second->tag(), itP->second->physicals, 0);
       parentsNum[itP->first] = num;
     }
   for(riter it = firstRegion(); it != lastRegion(); ++it)
     writeElementsMSH(fp, (*it)->tetrahedra, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, eleRenumbering);
+                     (*it)->tag(), (*it)->physicals);
   writeElementHeaderMSH(binary, fp, elements, MSH_HEX_8, MSH_HEX_27, MSH_HEX_20);
   for(itP = parents[1].begin(); itP != parents[1].end(); itP++)
     if(itP->first->getType() == TYPE_HEX) {
       writeElementsMSH(fp, itP->first, saveAll, version, binary, num,
-                       itP->second->tag(), itP->second->physicals, 0, eleRenumbering);
+                       itP->second->tag(), itP->second->physicals, 0);
       parentsNum[itP->first] = num;
     }
   for(riter it = firstRegion(); it != lastRegion(); ++it)
     writeElementsMSH(fp, (*it)->hexahedra, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, eleRenumbering);
+                     (*it)->tag(), (*it)->physicals);
   writeElementHeaderMSH(binary, fp, elements, MSH_PRI_6, MSH_PRI_18, MSH_PRI_15);
   for(itP = parents[1].begin(); itP != parents[1].end(); itP++)
     if(itP->first->getType() == TYPE_PRI) {
       writeElementsMSH(fp, itP->first, saveAll, version, binary, num,
-                       itP->second->tag(), itP->second->physicals, 0, eleRenumbering);
+                       itP->second->tag(), itP->second->physicals, 0);
       parentsNum[itP->first] = num;
     }
   for(riter it = firstRegion(); it != lastRegion(); ++it)
     writeElementsMSH(fp, (*it)->prisms, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, eleRenumbering);
+                     (*it)->tag(), (*it)->physicals);
   writeElementHeaderMSH(binary, fp, elements, MSH_PYR_5, MSH_PYR_14, MSH_PYR_13);
   for(itP = parents[1].begin(); itP != parents[1].end(); itP++)
     if(itP->first->getType() == TYPE_PYR) {
       writeElementsMSH(fp, itP->first, saveAll, version, binary, num,
-                       itP->second->tag(), itP->second->physicals, 0, eleRenumbering);
+                       itP->second->tag(), itP->second->physicals, 0);
       parentsNum[itP->first] = num;
     }
   for(riter it = firstRegion(); it != lastRegion(); ++it)
     writeElementsMSH(fp, (*it)->pyramids, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, eleRenumbering);
+                     (*it)->tag(), (*it)->physicals);
   writeElementHeaderMSH(binary, fp, elements, MSH_POLYH_);
   for(itP = parents[1].begin(); itP != parents[1].end(); itP++)
     if(itP->first->getType() == TYPE_POLYH) {
       writeElementsMSH(fp, itP->first, saveAll, version, binary, num,
-                       itP->second->tag(), itP->second->physicals, 0, eleRenumbering);
+                       itP->second->tag(), itP->second->physicals, 0);
       parentsNum[itP->first] = num;
     }
   for(riter it = firstRegion(); it != lastRegion(); ++it)
     writeElementsMSH(fp, (*it)->polyhedra, saveAll, version, binary, num,
-                     (*it)->tag(), (*it)->physicals, parentsNum, eleRenumbering);
+                     (*it)->tag(), (*it)->physicals, parentsNum);
   
   if(binary) fprintf(fp, "\n");
 
diff --git a/Geo/Homology.h b/Geo/Homology.h
index 57786ec3b2..b1a29f0db4 100644
--- a/Geo/Homology.h
+++ b/Geo/Homology.h
@@ -108,13 +108,14 @@ class Homology
    
    // write the generators to a file
    bool writeGeneratorsMSH(std::string fileName, bool binary=false){
-     if(!_model->writeMSH(fileName, 2.0, binary, false, false, 1.0, false)) return false;
+     if(!_model->writeMSH(fileName, 2.0, binary)) return false;
+     /*
      for(int i = 0; i < 4; i++){
        for(int j = 0; j < _generators[i].size(); j++){
          Chain* chain = _generators[i].at(j);
          if(!chain->writeChainMSH(fileName)) return false;
        }
-     }
+     }*/
      Msg::Info("Wrote results to %s.", fileName.c_str());
      printf("Wrote results to %s. \n", fileName.c_str());
      
diff --git a/Post/PViewDataGModelIO.cpp b/Post/PViewDataGModelIO.cpp
index 48304312c7..ac55cda4d5 100644
--- a/Post/PViewDataGModelIO.cpp
+++ b/Post/PViewDataGModelIO.cpp
@@ -117,7 +117,7 @@ bool PViewDataGModel::writeMSH(std::string fileName, bool binary)
 
   GModel *model = _steps[0]->getModel();
 
-  if(!model->writeMSH(fileName, 2.0, binary, false, false, 1.0, false)) return false;
+  if(!model->writeMSH(fileName, 2.0, binary)) return false;
 
   // append data
   FILE *fp = fopen(fileName.c_str(), binary ? "ab" : "a");
diff --git a/tutorial/t10.geo b/tutorial/t10.geo
index 0335ebc0b3..bd92364cc0 100644
--- a/tutorial/t10.geo
+++ b/tutorial/t10.geo
@@ -6,7 +6,10 @@
  *
  *********************************************************************/
  
-// Homology computation in Gmsh finds representatives of (relative) homology spaces using a mesh of a model. Those representatives generate the (relative) homology spaces of the model. Alternativelly, Gmsh can only look for the ranks of the (relative) homology spaces, the Betti numbers of the model. 
+// Homology computation in Gmsh finds representative chains of (relative) homology spaces using a mesh of a model. Those representatives generate the (relative) homology spaces of the model. Alternatively, Gmsh can only look for the ranks of the (relative) homology spaces, the Betti numbers of the model. 
+
+// The generators chains are stored in a given .msh-file as physical groups, whose mesh elements are oriented such that their coeffecients are 1 in the generator chain.
+
 
 // Create an example geometry
 
@@ -137,13 +140,12 @@ HomCut("t10_homcut.msh") = {{69}, {70, 71, 72, 73}};
 // Only find and print the ranks of the relative homology spaces (Betti numbers).
 HomRank {{69},{70, 71, 72, 73}};
 
-// Hide mesh elements
+// Hide mesh elements for instant visualization
 Mesh.Points = 0;
 Mesh.Lines = 0;
 Mesh.Triangles = 0;
 Mesh.Tetrahedra = 0; 
 
-
 // More examples (uncomment):
 //  HomGen("t10_homgen.msh_1") = {{69}, {}}; 
 //  HomGen("t10_homgen.msh_2") = {{69}, {74}}; 
-- 
GitLab