Dear Gmsh developers,

I am having a performance problem with Gmsh 2.5.0 library API where
entity additions (GModel::addVertex()/Line()/...) are painfully slow
(taking hours) with hundreds of thousands (100,000s) of entities. With
a bit of profiling it turned out that the problem was coming from
getMaxElementaryNumber() counting the max elementary number each time
an entity is added, which is obviously of an O(n^2) cost. Hence I
modified the relevant codes as attached in order to store the number
in a array instead. In my case the performance improved to order of 10
secondes which is more than 100x of speedup. The patch is only meant
to "work for me" which means it may not be generic enough for
everyone, but I hope this kind of improvement goes into the next
release of Gmsh.

Takuya

Takuya OSHIMA, Ph.D.
Faculty of Engineering, Niigata University
8050 Ikarashi-Ninocho, Nishi-ku, Niigata, 950-2181, JAPAN

>From 4eb87c65261aa787ea7953cf3559e485b44ba427 Mon Sep 17 00:00:00 2001
From: Takuya OSHIMA <[email protected]>
Date: Fri, 12 Aug 2011 20:02:45 +0900
Subject: [PATCH 4/4] Performance improvement in entity insertions to GModel

---
 Geo/GModel.cpp |   70 ++++++++++++++++++++++++++++++++++++++++++++++++--------
 Geo/GModel.h   |   36 +++++++++++++++++++++++++---
 2 files changed, 92 insertions(+), 14 deletions(-)

diff --git a/Geo/GModel.cpp b/Geo/GModel.cpp
index dffdfb7..36bfa71 100644
--- a/Geo/GModel.cpp
+++ b/Geo/GModel.cpp
@@ -142,6 +142,8 @@ GModel::GModel(std::string name)
 #if defined(HAVE_MESH)
   _fields = new FieldManager();
 #endif
+  for(int i = 0; i < 4; ++i)
+    maxElementaryNumber[i] = 0;
 }
 
 GModel::~GModel()
@@ -230,6 +232,8 @@ void GModel::destroy()
     delete *it;
   vertices.clear();
 
+  _recalculateMaxElementaryNumber(-1);
+
   destroyMeshCaches();
 
   MVertex::resetGlobalNumber();
@@ -336,25 +340,45 @@ GVertex *GModel::getVertexByTag(int n) const
 void GModel::remove(GRegion *r)
 {
   riter it = std::find(firstRegion(), lastRegion(), r);
-  if(it != (riter)regions.end()) regions.erase(it);
+  if(it != (riter)regions.end()){
+    const int tag = (*it)->tag(), dim = (*it)->dim();
+    regions.erase(it);
+    if(maxElementaryNumber[dim] == std::abs(tag))
+      _recalculateMaxElementaryNumber(dim);
+  }
 }
 
 void GModel::remove(GFace *f)
 {
   fiter it = std::find(firstFace(), lastFace(), f);
-  if(it != faces.end()) faces.erase(it);
+  if(it != faces.end()){
+    const int tag = (*it)->tag(), dim = (*it)->dim();
+    faces.erase(it);
+    if(maxElementaryNumber[dim] == std::abs(tag))
+      _recalculateMaxElementaryNumber(dim);
+  }
 }
 
 void GModel::remove(GEdge *e)
 {
   eiter it = std::find(firstEdge(), lastEdge(), e);
-  if(it != edges.end()) edges.erase(it);
+  if(it != edges.end()){
+    const int tag = (*it)->tag(), dim = (*it)->dim();
+    edges.erase(it);
+    if(maxElementaryNumber[dim] == std::abs(tag))
+      _recalculateMaxElementaryNumber(dim);
+  }
 }
 
 void GModel::remove(GVertex *v)
 {
   viter it = std::find(firstVertex(), lastVertex(), v);
-  if(it != vertices.end()) vertices.erase(it);
+  if(it != vertices.end()){
+    const int tag = (*it)->tag(), dim = (*it)->dim();
+    vertices.erase(it);
+    if(maxElementaryNumber[dim] == std::abs(tag))
+      _recalculateMaxElementaryNumber(dim);
+  }
 }
 
 void GModel::snapVertices()
@@ -401,15 +425,38 @@ void GModel::getEntities(std::vector<GEntity*> &entities)
   entities.insert(entities.end(), regions.begin(), regions.end());
 }
 
-int GModel::getMaxElementaryNumber(int dim)
+void GModel::_recalculateMaxElementaryNumber(int dim)
 {
   std::vector<GEntity*> entities;
   getEntities(entities);
-  int num = 0;
-  for(unsigned int i = 0; i < entities.size(); i++)
-    if(dim < 0 || entities[i]->dim() == dim)
-      num = std::max(num, std::abs(entities[i]->tag()));
-  return num;
+  if(dim >= 0){
+    int num = 0;
+    for(unsigned int i = 0; i < entities.size(); i++)
+      if(entities[i]->dim() == dim)
+	num = std::max(num, std::abs(entities[i]->tag()));
+    maxElementaryNumber[dim] = num;
+  }else{
+    for(int i = 0; i < 4; ++i)
+      maxElementaryNumber[i] = 0;
+    for(unsigned int i = 0; i < entities.size(); i++){
+      const int dim = entities[i]->dim();
+      if(dim >= 0)
+	maxElementaryNumber[dim] = std::max(maxElementaryNumber[dim],
+					    std::abs(entities[i]->tag()));
+    }
+  }
+}
+
+int GModel::getMaxElementaryNumber(int dim)
+{
+  if(dim >= 0)
+    return maxElementaryNumber[dim];
+  else{
+    int num = 0;
+    for(int i = 0; i < 4; ++i)
+      num = std::max(num, maxElementaryNumber[i]);
+    return num;
+  }
 }
 
 bool GModel::noPhysicalGroups()
@@ -1992,6 +2039,9 @@ void GModel::glue(double eps)
 void GModel::insertRegion(GRegion *r)
 {
   regions.insert(r);
+  const int dim = r->dim();
+  maxElementaryNumber[dim] = std::max(maxElementaryNumber[dim],
+				      std::abs(r->tag()));
 }
 
 GEdge *getNewModelEdge(GFace *gf1, GFace *gf2, 
diff --git a/Geo/GModel.h b/Geo/GModel.h
index e1d5ea5..07f9527 100644
--- a/Geo/GModel.h
+++ b/Geo/GModel.h
@@ -102,6 +102,9 @@ class GModel
   void _storePhysicalTagsInEntities(int dim,
                                     std::map<int, std::map<int, std::string> > &map);
 
+  // recalculate max elementary number
+  void _recalculateMaxElementaryNumber(int dim);
+
   // entity that is currently being meshed (used for error reporting)
   GEntity *_currentMeshEntity;
 
@@ -116,6 +119,7 @@ class GModel
   std::set<GFace*, GEntityLessThan> faces;
   std::set<GEdge*, GEntityLessThan> edges;
   std::set<GVertex*, GEntityLessThan> vertices;
+  int maxElementaryNumber[4];
 
   void insertRegion (GRegion*);
 
@@ -214,10 +218,34 @@ class GModel
   std::vector<GVertex*> bindingsGetVertices();
 
   // add/remove an entity in the model
-  void add(GRegion *r) { regions.insert(r); }
-  void add(GFace *f) { faces.insert(f); }
-  void add(GEdge *e) { edges.insert(e); }
-  void add(GVertex *v) { vertices.insert(v); }
+  void add(GRegion *r)
+  {
+    regions.insert(r);
+    const int dim = r->dim();
+    maxElementaryNumber[dim] = std::max(maxElementaryNumber[dim],
+					std::abs(r->tag()));
+  }
+  void add(GFace *f)
+  {
+    faces.insert(f);
+    const int dim = f->dim();
+    maxElementaryNumber[dim] = std::max(maxElementaryNumber[dim],
+					std::abs(f->tag()));
+  }
+  void add(GEdge *e)
+  {
+    edges.insert(e);
+    const int dim = e->dim();
+    maxElementaryNumber[dim] = std::max(maxElementaryNumber[dim],
+					std::abs(e->tag()));
+  }
+  void add(GVertex *v)
+  {
+    vertices.insert(v);
+    const int dim = v->dim();
+    maxElementaryNumber[dim] = std::max(maxElementaryNumber[dim],
+					std::abs(v->tag()));
+  }
   void remove(GRegion *r);
   void remove(GFace *f);
   void remove(GEdge *e);
-- 
1.7.6

_______________________________________________
gmsh mailing list
[email protected]
http://www.geuz.org/mailman/listinfo/gmsh

Reply via email to