Modified matching procedure

Common/OpenFile.cpp

@@ -362,10 +362,9 @@ int MergeFile(std::string fileName, bool warnIfMissing)
         GModel* tmp2 = GModel::current();
         GModel* tmp = new GModel();
         tmp->readMSH(fileName);
-        if(GeomMeshMatcher::instance()->match(tmp2, tmp))
-          fileName = "out.msh";
+        GeomMeshMatcher::instance()->match(tmp2, tmp);
         delete tmp;
-      }
+      } else
 	status = GModel::current()->readMSH(fileName);
 #if defined(HAVE_POST)
       if(status > 1) status = PView::readMSH(fileName);

Geo/GeomMeshMatcher.cpp

@@ -21,7 +21,11 @@
 #include "MQuadrangle.h"
 #include "MVertex.h"
 #include "MLine.h"
+#include "MPyramid.h"
+#include "MPrism.h"
 #include "MPoint.h"
+#include "MHexahedron.h"
+#include "MTetrahedron.h"
 
 GeomMeshMatcher *GeomMeshMatcher::_gmm_instance = 0;
 
@@ -591,12 +595,93 @@ int GeomMeshMatcher::forceTomatch(GModel *geom, GModel *mesh, const double TOL)
   return 0;
 }
 
+static void copy_vertices (GVertex *to, GVertex *from, std::map<MVertex*,MVertex*> &_mesh_to_geom){
+  if (from) {
+    for (int i=0;i<from->mesh_vertices.size();i++){
+      MVertex *v_from = from->mesh_vertices[i];
+      MVertex *v_to = new MVertex (v_from->x(),v_from->y(),v_from->z(), to);
+      to->mesh_vertices.push_back(v_to);
+      _mesh_to_geom[v_from] = v_to;
+    }
+  }
+}
+static void copy_vertices (GRegion *to, GRegion *from, std::map<MVertex*,MVertex*> &_mesh_to_geom){
+  for (int i=0;i<from->mesh_vertices.size();i++){
+    MVertex *v_from = from->mesh_vertices[i];
+    MVertex *v_to = new MVertex (v_from->x(),v_from->y(),v_from->z(), to);
+    to->mesh_vertices.push_back(v_to);
+    _mesh_to_geom[v_from] = v_to;
+  }
+}
+static void copy_vertices (GEdge *to, GEdge *from, std::map<MVertex*,MVertex*> &_mesh_to_geom){
+  for (int i=0;i<from->mesh_vertices.size();i++){
+    MVertex *v_from = from->mesh_vertices[i];
+    double t = to->parFromPoint ( SPoint3(v_from->x(),v_from->y(),v_from->z() ) );
+    MEdgeVertex *v_to = new MEdgeVertex (v_from->x(),v_from->y(),v_from->z(), to, t );
+    to->mesh_vertices.push_back(v_to);
+    _mesh_to_geom[v_from] = v_to;
+  }
+}
+static void copy_vertices (GFace *to, GFace *from, std::map<MVertex*,MVertex*> &_mesh_to_geom){
+  for (int i=0;i<from->mesh_vertices.size();i++){
+    MVertex *v_from = from->mesh_vertices[i];
+    SPoint2 uv = to->parFromPoint ( SPoint3(v_from->x(),v_from->y(),v_from->z() ) );
+    MFaceVertex *v_to = new MFaceVertex (v_from->x(),v_from->y(),v_from->z(), to, uv.x(),uv.y() );
+    to->mesh_vertices.push_back(v_to);
+    _mesh_to_geom[v_from] = v_to;
+  }
+}
+
+template <class ELEMENT>
+static void copy_elements (std::vector<ELEMENT*> &to, std::vector<ELEMENT*> &from, std::map<MVertex*,MVertex*> &_mesh_to_geom){
+  MElementFactory toto;
+  for (int i=0;i < to.size();i++){
+    ELEMENT *e = from[i];
+    std::vector<MVertex*> nodes;
+    for(int j=0;j<e->getNumVertices();j++)
+      nodes.push_back(_mesh_to_geom[e->getVertex(j)]);
+    to.push_back( (ELEMENT*)(toto.create(e->getType(), nodes) ));
+  }
+}
+
+
+void copy_vertices (GModel *geom, GModel *mesh, std::map<MVertex*,MVertex*> &_mesh_to_geom){
+  // copy all elements
+  for (GModel::viter vit = geom->firstVertex() ; vit != geom->lastVertex(); ++vit)
+    copy_vertices(*vit,mesh->getVertexByTag((*vit)->tag()),_mesh_to_geom);
+  for (GModel::eiter eit = geom->firstEdge() ; eit != geom->lastEdge(); ++eit)
+    copy_vertices(*eit,mesh->getEdgeByTag((*eit)->tag()),_mesh_to_geom);
+  for (GModel::fiter fit = geom->firstFace() ; fit != geom->lastFace(); ++fit)
+    copy_vertices(*fit,mesh->getFaceByTag((*fit)->tag()),_mesh_to_geom);
+  for (GModel::riter rit = geom->firstRegion() ; rit != geom->lastRegion(); ++rit)
+    copy_vertices(*rit,mesh->getRegionByTag((*rit)->tag()),_mesh_to_geom);
+}
+void copy_elements (GModel *geom, GModel *mesh, std::map<MVertex*,MVertex*> &_mesh_to_geom){
+  // copy all elements
+  for (GModel::viter vit = geom->firstVertex() ; vit != geom->lastVertex(); ++vit)
+    copy_elements<MPoint>((*vit)->points,mesh->getVertexByTag((*vit)->tag())->points,_mesh_to_geom);
+  for (GModel::eiter eit = geom->firstEdge() ; eit != geom->lastEdge(); ++eit)
+    copy_elements<MLine>((*eit)->lines,mesh->getEdgeByTag((*eit)->tag())->lines,_mesh_to_geom);
+  for (GModel::fiter fit = geom->firstFace() ; fit != geom->lastFace(); ++fit) {
+    copy_elements<MTriangle>((*fit)->triangles,mesh->getFaceByTag((*fit)->tag())->triangles,_mesh_to_geom);
+    copy_elements<MQuadrangle>((*fit)->quadrangles,mesh->getFaceByTag((*fit)->tag())->quadrangles,_mesh_to_geom);
+  }
+  for (GModel::riter rit = geom->firstRegion() ; rit != geom->lastRegion(); ++rit) {
+    copy_elements<MTetrahedron>((*rit)->tetrahedra,mesh->getRegionByTag((*rit)->tag())->tetrahedra,_mesh_to_geom);
+    copy_elements<MHexahedron>((*rit)->hexahedra,mesh->getRegionByTag((*rit)->tag())->hexahedra,_mesh_to_geom);
+    copy_elements<MPrism>((*rit)->prisms,mesh->getRegionByTag((*rit)->tag())->prisms,_mesh_to_geom);
+    copy_elements<MPyramid>((*rit)->pyramids,mesh->getRegionByTag((*rit)->tag())->pyramids,_mesh_to_geom);
+  }
+}
+
+
 int GeomMeshMatcher::match(GModel *geom, GModel *mesh)
 {
   mesh->createTopologyFromMesh();
   bool ok = true;
   // This will match VERTICES
   std::vector<Pair<GVertex*, GVertex*> > *coresp_v = matchVertices(geom, mesh,ok);
+  ok = true;
   if (ok) {
     // This will match EDGES
     std::vector<Pair<GEdge*, GEdge*> > *coresp_e = matchEdges(geom, mesh, coresp_v,ok);
@@ -608,7 +693,10 @@ int GeomMeshMatcher::match(GModel *geom, GModel *mesh)
         //std::vector<Pair<GRegion*, GRegion*> >* coresp_r =
         matchRegions(geom, mesh, coresp_f,ok);
         if (ok) {
-          mesh->writeMSH("out.msh",2.0,false,true);
+          //mesh->writeMSH("out.msh",2.0,false,true);
+	  std::map<MVertex*,MVertex*> _mesh_to_geom;
+	  copy_vertices(geom, mesh, _mesh_to_geom);
+	  copy_elements(geom, mesh, _mesh_to_geom);
           return 1;
         } else {
           Msg::Error("Could not match every region !");