From 688ffb561392ba0bb53788a6ed6e86fed4006ec1 Mon Sep 17 00:00:00 2001
From: Emilie Marchandise <emilie.marchandise@uclouvain.be>
Date: Fri, 7 Oct 2011 13:17:59 +0000
Subject: [PATCH] Solved bug in GfaceCompound with lsys delete

---
 Geo/GFaceCompound.cpp | 59 +++++++++++++++++++++++++------------------
 Geo/GFaceCompound.h   |  9 +++----
 Geo/GModelIO_Geo.cpp  |  2 +-
 Mesh/meshGFace.cpp    |  4 +--
 4 files changed, 41 insertions(+), 33 deletions(-)

diff --git a/Geo/GFaceCompound.cpp b/Geo/GFaceCompound.cpp
index dcfef1ad22..e62438a7cc 100644
--- a/Geo/GFaceCompound.cpp
+++ b/Geo/GFaceCompound.cpp
@@ -919,29 +919,14 @@ void GFaceCompound::computeALoop(std::set<GEdge*> &_unique, std::list<GEdge*> &l
 }
 
 GFaceCompound::GFaceCompound(GModel *m, int tag, std::list<GFace*> &compound,
-			     std::list<GEdge*> &U0, 
-                             linearSystem<double> *lsys, typeOfMapping mpg, 
+			     std::list<GEdge*> &U0, typeOfMapping mpg, 
                              int allowPartition)
   : GFace(m, tag), _compound(compound),  oct(0), _U0(U0), 
-    _lsys(lsys),_mapping(mpg), _allowPartition(allowPartition)
+    _mapping(mpg), _allowPartition(allowPartition)
 {
   ONE = new simpleFunction<double>(1.0);
   MONE = new simpleFunction<double>(-1.0);
 
-  if (!_lsys) {
-#if defined(HAVE_PETSC) && !defined(HAVE_TAUCS)
-    _lsys = new linearSystemPETSc<double>;
-#elif defined(HAVE_GMM) && !defined(HAVE_TAUCS)
-    linearSystemGmm<double> *_lsysb = new linearSystemGmm<double>;
-    _lsysb->setGmres(1);
-    _lsys = _lsysb;
-#elif defined(HAVE_TAUCS) 
-    _lsys = new linearSystemCSRTaucs<double>;
-#else
-    _lsys = new linearSystemFull<double>;
-#endif
- }
-
   for(std::list<GFace*>::iterator it = _compound.begin(); it != _compound.end(); ++it){
     //EMI FIX
     //if ((*it)->tag() == 3) _mapping = CONFORMAL;
@@ -965,7 +950,6 @@ GFaceCompound::~GFaceCompound()
     Octree_Delete(oct);
     delete [] _gfct;
   }
-  if (_lsys)delete _lsys;
   delete ONE;
   delete MONE;
 }
@@ -1057,6 +1041,20 @@ SPoint2 GFaceCompound::getCoordinates(MVertex *v) const
 
 void GFaceCompound::parametrize(iterationStep step, typeOfMapping tom) const
 {  
+
+linearSystem<double> *_lsys = 0;
+#if defined(HAVE_PETSC) && !defined(HAVE_TAUCS)
+  _lsys = new linearSystemPETSc<double>;
+#elif defined(HAVE_GMM) && !defined(HAVE_TAUCS)
+  linearSystemGmm<double> *_lsysb = new linearSystemGmm<double>;
+  _lsysb->setGmres(1);
+  _lsys = _lsysb;
+#elif defined(HAVE_TAUCS) 
+  _lsys = new linearSystemCSRTaucs<double>;
+#else
+  _lsys = new linearSystemFull<double>;
+#endif
+
   dofManager<double> myAssembler(_lsys);
 
   if(_type == UNITCIRCLE){
@@ -1133,11 +1131,10 @@ void GFaceCompound::parametrize(iterationStep step, typeOfMapping tom) const
     }
   }
 
-  _lsys->clear();
+  delete _lsys;
 
 }
 
-
 bool GFaceCompound::parametrize_conformal_spectral() const
 {
 #if !defined(HAVE_PETSC) && !defined(HAVE_SLEPC)
@@ -1247,9 +1244,9 @@ bool GFaceCompound::parametrize_conformal_spectral() const
       k = k+2;
     }
     
-    lsysA->clear();
-    lsysB->clear();
-    
+    delete lsysA;
+    delete lsysB;
+
     return checkOverlap();
     
   }
@@ -1259,6 +1256,20 @@ bool GFaceCompound::parametrize_conformal_spectral() const
 
 bool GFaceCompound::parametrize_conformal() const
 {
+
+  linearSystem<double> *_lsys = 0;
+#if defined(HAVE_PETSC) && !defined(HAVE_TAUCS)
+  _lsys = new linearSystemPETSc<double>;
+#elif defined(HAVE_GMM) && !defined(HAVE_TAUCS)
+  linearSystemGmm<double> *_lsysb = new linearSystemGmm<double>;
+  _lsysb->setGmres(1);
+  _lsys = _lsysb;
+#elif defined(HAVE_TAUCS) 
+  _lsys = new linearSystemCSRTaucs<double>;
+#else
+  _lsys = new linearSystemFull<double>;
+#endif
+
   dofManager<double> myAssembler(_lsys);
 
   MVertex *v1  = _ordered[0];
@@ -1324,7 +1335,7 @@ bool GFaceCompound::parametrize_conformal() const
     coordinates[v] = SPoint3(value1,value2,0.0);
   }
 
-  _lsys->clear();
+  delete _lsys; 
 
   //check for overlapping triangles
   return checkOverlap();
diff --git a/Geo/GFaceCompound.h b/Geo/GFaceCompound.h
index 0ff98f451f..0d37acb896 100644
--- a/Geo/GFaceCompound.h
+++ b/Geo/GFaceCompound.h
@@ -98,8 +98,7 @@ class GFaceCompound : public GFace {
                    double &_u, double &_v) const;
   virtual double locCurvature(MTriangle *t, double u, double v) const;
   void printStuff(int iNewton=0) const;
-  bool trivial() const ;
-  linearSystem <double> *_lsys;
+  bool trivial() const;
   double getSizeH() const;
   double getSizeBB(const std::list<GEdge* > &elist) const;
   SBoundingBox3d boundEdges(const std::list<GEdge* > &elist) const;
@@ -109,8 +108,7 @@ class GFaceCompound : public GFace {
 
  public: 
   GFaceCompound(GModel *m, int tag, std::list<GFace*> &compound,
-		std::list<GEdge*> &U0, linearSystem<double>* lsys =0,
-                typeOfMapping typ = HARMONIC, int allowPartition=1);
+		std::list<GEdge*> &U0, typeOfMapping typ = HARMONIC, int allowPartition=1);
   virtual ~GFaceCompound();
   Range<double> parBounds(int i) const 
   { return trivial() ? (*(_compound.begin()))->parBounds(i) : Range<double>(-1, 1); }
@@ -150,8 +148,7 @@ class GFaceCompound : public GFace {
   typedef enum {HARMONIC=1,CONFORMAL=2, CONVEXCOMBINATION=3, MULTISCALE=4, RBF=5} typeOfMapping;
   typedef enum {UNITCIRCLE, SQUARE} typeOfIsomorphism;
  GFaceCompound(GModel *m, int tag, std::list<GFace*> &compound,
-	       std::list<GEdge*> &U0, linearSystem<double>* lsys =0,
-                typeOfMapping typ = HARMONIC, int allowPartition=1)
+	       std::list<GEdge*> &U0, typeOfMapping typ = HARMONIC, int allowPartition=1)
     : GFace(m, tag)
   {
     Msg::Error("Gmsh has to be compiled with solver support to use GFaceCompounds");
diff --git a/Geo/GModelIO_Geo.cpp b/Geo/GModelIO_Geo.cpp
index 8c84e6d6f4..6afe76ae17 100644
--- a/Geo/GModelIO_Geo.cpp
+++ b/Geo/GModelIO_Geo.cpp
@@ -191,7 +191,7 @@ int GModel::importGEOInternals()
 	if (param == 1) typ =  GFaceCompound::CONFORMAL;
 	if (param == 2) typ =  GFaceCompound::RBF;
 	int algo = CTX::instance()->mesh.remeshAlgo;
-        f = new GFaceCompound(this, s->Num, comp, U0, 0, typ, algo);
+        f = new GFaceCompound(this, s->Num, comp, U0, typ, algo);
 
 	f->meshAttributes.recombine = s->Recombine;
 	f->meshAttributes.recombineAngle = s->RecombineAngle;
diff --git a/Mesh/meshGFace.cpp b/Mesh/meshGFace.cpp
index 76a6c52ca7..393818c585 100644
--- a/Mesh/meshGFace.cpp
+++ b/Mesh/meshGFace.cpp
@@ -1514,7 +1514,7 @@ void deMeshGFace::operator() (GFace *gf)
 }
 
 // for debugging, change value from -1 to -100;
-int debugSurface = -100 ; //-1; 
+int debugSurface = -1 ; //-1; 
 
 void meshGFace::operator() (GFace *gf)
 {
@@ -1712,7 +1712,7 @@ void partitionAndRemesh(GFaceCompound *gf)
     Msg::Info("Parametrize Compound Surface (%d) = %d discrete face",
               num_gfc, pf->tag());
     
-    GFaceCompound *gfc = new GFaceCompound(gf->model(), num_gfc, f_compound, U0, 0,
+    GFaceCompound *gfc = new GFaceCompound(gf->model(), num_gfc, f_compound, U0,
                                             gf->getTypeOfMapping());
     gfc->meshAttributes.recombine = gf->meshAttributes.recombine;
     gf->model()->add(gfc);
-- 
GitLab