Index: /issm/trunk/src/c/Bamgx/Bamgx.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/Bamgx.cpp	(revision 2788)
+++ /issm/trunk/src/c/Bamgx/Bamgx.cpp	(revision 2789)
@@ -53,5 +53,4 @@
 	int fileout=0;
 	int AbsError=0,nbjacoby=1,allquad=0;
-	int NoMeshReconstruction=0;
 	int Rescaling = 1;
 	double costheta=2;
@@ -91,6 +90,4 @@
 
 	// some verification
-	NoMeshReconstruction= fmeshr !=0;
-	if (!fmeshback) fmeshback=fmeshr;
 	if ( maxsubdiv > boundmaxsubdiv || maxsubdiv <= 1.0)
 	{
@@ -113,4 +110,5 @@
 	if(bamgmesh->NumTriangles==0){
 		/*Mesh generation {{{1*/
+
 		if (verbosity>0) printf("Construction of a mesh from a given geometry\n");
 		if (verbosity>1) printf("   Processing geometry...\n");
@@ -167,19 +165,19 @@
 	else{
 		/*Anisotropic mesh adaptation {{{1*/
+
+		// read background mesh 
 		if (verbosity>0) printf("Anisotropic mesh adaptation\n");
-
-		/*Read background mesh from simple delaunay triangulation: */
-		if (verbosity>1) printf("   Processing initial mesh...\n");
-		Triangles BTh(bamgmesh,bamgargs); // read the background mesh 
-
-		throw ErrorException(__FUNCT__,exprintf("TESTTESTTESTTEST"));
-		printf("ok1\n"); 
+		if (verbosity>1) printf("   Processing initial mesh and geometry...\n");
+		Triangles BTh(bamgmesh,bamgargs); 
 		hmin = Max(hmin,BTh.MinimalHmin());
-		printf("ok2\n"); 
 		hmax = Min(hmax,BTh.MaximalHmax());
 
+		//?????TEST
 		BTh.MakeQuadTree();
-		if (fmetrix){
-			BTh.ReadMetric(fmetrix,hmin,hmax,coef);
+
+		//build metric if not given in input
+		if (verbosity>1) printf("   Processing Metric...\n");
+		if (bamgargs->metric){
+			BTh.ReadMetric(bamgargs,hmin,hmax,coef);
 		}
 		else { // init with hmax 
@@ -188,226 +186,66 @@
 			 BTh[iv].m = Mhmax;
 		}
-		if (fMbb) {
-			solMbb = ReadbbFile(fMbb,nbsolbb,lsolbb,2,2);
-			if (lsolbb != BTh.nbv) 
-			  {
-				cerr << "fatal error  nbsol " << nbsolbb << " " << lsolbb<< " =! " << BTh.nbv << endl;
-				cerr << "  size of sol incorrect " << endl;
-				MeshError(99);
-			  }
-			assert(lsolbb==BTh.nbv);
-			BTh.IntersectConsMetric(solMbb,nbsolbb,0,hmin,hmax,sqrt(err)*coef,1e30,AbsError?0.0:cutoff,
-						nbjacoby,Rescaling,power,ChoiseHessien);
-			if(!rbbeqMbb)
-			 delete [] solMbb,solMbb =0;
-
-		}
-		if (fMBB) {
-			solMBB = ReadBBFile(fMBB,nbsolBB,lsolBB,typesolsBB,2,2);
-			if (lsolBB != BTh.nbv) 
-			  {
-				cerr << "fatal error  nbsol " << nbsolBB << " " << lsolBB << " =! " << BTh.nbv << endl;
-				cerr << "  size of sol incorrect " << endl;
-				MeshError(99);
-			  }
-			assert(lsolBB==BTh.nbv);
-			BTh.IntersectConsMetric(solMBB,nbsolBB,0,hmin,hmax,sqrt(err)*coef,1e30,AbsError?0.0:cutoff,
-						nbjacoby,Rescaling,ChoiseHessien);
-			if(!rBBeqMBB)
-			 delete [] solMBB,solMBB =0;
-
-		}
 
 		BTh.IntersectGeomMetric(errg,iso);
-
 		if(gradation) BTh.SmoothMetric(gradation);
 		BTh.MaxSubDivision(maxsubdiv);
 		BTh.BoundAnisotropy(anisomax,hminaniso);
 
-		if(foM) {
-			if(verbosity >2)
-			 cout << " -- write Metric  file " << foM <<endl;
-
-			ofstream f(foM);
-			if(f) BTh.WriteMetric(f,iso);
-		}
-
-		if ( fileout) {
-			if ( NoMeshReconstruction)
-			 if (( fmeshback == fmeshr) || (!strcmp(fmeshback,fmeshr)))
-			  Thr=&BTh,Thb=0; // back and r mesh are the same 
-			 else
-			  Thr= new Triangles(fmeshr,cutoffradian),Thb=&BTh; // read the new 
-
-
-
-			Triangles & Th( *(NoMeshReconstruction 
-							?  new Triangles(*Thr,&Thr->Gh,Thb,maxnbv) // copy the mesh + free space to modification  
-							:  new Triangles(maxnbv,BTh,KeepBackVertices)     // construct a new mesh 
-							));
-			if (Thr != &BTh) delete Thr;
-
-			if(costheta<=1.0)
-			 Th.MakeQuadrangles(costheta);
-			if (allquad)
-			 Th.SplitElement(allquad==2);
-			if(splitcorners)
-			 Th.SplitInternalEdgeWithBorderVertices();
-			Th.ReNumberingTheTriangleBySubDomain();
-
-			if(verbosity>3) Th.ShowHistogram();
-			if(NbSmooth>0) Th.SmoothingVertex(NbSmooth,omega);
-			if(verbosity>3 && NbSmooth>0) Th.ShowHistogram();
-
-			Th.BTh.ReMakeTriangleContainingTheVertex();
-
-			if (fmeshout) Th.Write(fmeshout  ,Triangles::BDMesh);
-			if (famfmt)   Th.Write(famfmt    ,Triangles::am_fmtMesh);
-			if (fam)      Th.Write(fam       ,Triangles::amMesh);
-			if (famdba)   Th.Write(famdba    ,Triangles::amdbaMesh);
-			if (fftq)     Th.Write(fftq      ,Triangles::ftqMesh);
-			if (fmsh)     Th.Write(fmsh      ,Triangles::mshMesh);
-			if (fnopo)    Th.Write(fnopo     ,Triangles::NOPOMesh);
-
-			if ( ( rbb && wbb)  ||( rBB && wBB)){  // the code for interpolation 
-
-				if(verbosity >1)
-				  {
-					if (rbb ) 
-					 cout << " -- interpolation P1  files : " << rbb << " -> " << wbb <<endl;
-					if (rBB ) 
-					 cout << " -- interpolation P1  files : " << rBB<< " -> " << wBB <<endl;
-				  }
-				const int dim=2;
-				// optimisation read only si rbb != fMbb
-
-				double *solbb=0;
-				double *solBB=0;
-
-				if (rbb) 
-				 solbb =  rbbeqMbb? solMbb : ReadbbFile(rbb,nbsolbb,lsolbb,2,2);
-				if (rBB) 
-				 solBB =  rBBeqMBB? solMBB : ReadBBFile(rBB,nbsolBB,lsolBB,typesolsBB,2,2);
-
-				// cout << " " << rbbeqMbb << " " <<  sol << " " << nbsol << " " << lsol << endl; 
-				if (!solBB && !solbb ) 
-				  {
-					cerr << " Fatal Error "  << "solBB =  " <<  solBB << " solbb= " << solbb << endl;
-					exit(2);}
-
-					ofstream *fbb = wbb ? new ofstream(wbb) :0;
-					ofstream *fBB = wBB ? new ofstream(wBB) :0;
-					Int4   nbfieldBB = 0, nbfieldbb = nbsolbb;
-					if (fbb)
-					 *fbb  << dim << " " << nbsolbb << " " << Th.nbv <<" " << 2 << endl; 
-					if (fBB)
-					  {
-						int i;
-						*fBB  << dim << " " << nbsolBB ;
-						for ( i=0;i<nbsolBB;i++)
-						 *fBB << " " << (typesolsBB ?typesolsBB[i]+1:1) ;
-						*fBB << " " << Th.nbv <<" " << 2 << endl; 
-						assert(dim==2);
-						for ( i=0;i<nbsolBB;i++)
-						 nbfieldBB += typesolsBB ? typesolsBB[i]+1 : 1;
-						// this code is good only if dim == 2 
-					  }
-					cout << "nb of field BB " << nbfieldBB << endl;
-					//		 if(fBB) fBB->precision(15);
-					//if(fbb) fbb->precision(15);
-					for(i=0;i<Th.nbv;i++)
-					  {
-						Int4 i0,i1,i2;
-						double a[3];
-						Icoor2 dete[3];
-						I2 I = Th.BTh.toI2(Th.vertices[i].r);
-						Triangle & tb = *Th.BTh.FindTriangleContening(I,dete);
-
-						if (tb.det>0) { // internal point 
-							a[0]= (Real8) dete[0]/ tb.det;
-							a[1]= (Real8) dete[1] / tb.det;
-							a[2] = (Real8) dete[2] / tb.det;
-							i0=Th.BTh.Number(tb[0]);
-							i1=Th.BTh.Number(tb[1]);
-							i2=Th.BTh.Number(tb[2]);}
-						else {
-							double aa,bb;
-
-							TriangleAdjacent  ta=CloseBoundaryEdge(I,&tb,aa,bb).Adj();
-
-							int k = ta;
-							Triangle & tc = *(Triangle *)ta;
-							i0=Th.BTh.Number(tc[0]);
-							i1=Th.BTh.Number(tc[1]);
-							i2=Th.BTh.Number(tc[2]);
-							a[VerticesOfTriangularEdge[k][1]] =aa;
-							a[VerticesOfTriangularEdge[k][0]] = bb;
-							a[OppositeVertex[k]] = 1- aa -bb;}
-
-							Int4  ibb0 = nbfieldbb*i0;
-							Int4  ibb1 = nbfieldbb*i1;
-							Int4  ibb2 = nbfieldbb*i2;
-
-							Int4  iBB0 = nbfieldBB*i0;
-							Int4  iBB1 = nbfieldBB*i1;
-							Int4  iBB2 = nbfieldBB*i2;
-							Int4 j=0;
-
-							for ( j=0;j<nbfieldbb;j++) 
-							 *fbb  << " " << ( a[0] * solbb[ibb0++] + a[1] * solbb[ibb1++] + a[2]* solbb[ibb2++]) ;
-
-							for (j=0;j<nbfieldBB;j++) 
-							 *fBB  << " " << ( a[0] * solBB[iBB0++] + a[1] * solBB[iBB1++] + a[2]* solBB[iBB2++]) ;
-
-							if (fbb) *fbb  << endl;
-							if (fBB) *fBB  << endl;
-					  }
-					if (fbb)
-					 delete fbb  ; // close
-					if (fBB)
-					 delete fBB  ; // close
-					if (solbb) 
-					 delete [] solbb;
-					if (solBB) 
-					 delete [] solBB;
-
-			}
-
-			if(verbosity>0) 
-			  {
-				if (Th.nbt-Th.NbOutT-Th.NbOfQuad*2) 
-				 cout  << " Nb Triangles = " << (Th.nbt-Th.NbOutT-Th.NbOfQuad*2);
-				if (Th.NbOfQuad ) 
-				 cout  << " Nb Quadrilaterals = " << Th.NbOfQuad  ;
-				cout   << endl;
-			  }
-
-			//Build output
-			NumVerticesOut=Th.nbv;
-			VerticesOut=(double*)xmalloc(3*NumVerticesOut*sizeof(double));
-			for (i=0;i<NumVerticesOut;i++){
-				VerticesOut[i*3+0]=Th.vertices[i].r.x;
-				VerticesOut[i*3+1]=Th.vertices[i].r.y;
-				VerticesOut[i*3+2]=Th.vertices[i].ReferenceNumber;
-			}
-
-			NumTrianglesOut=Th.nbt-Th.NbOutT; //number of triangles - number of external triangles
-			TrianglesOut=(double*)xmalloc(4*NumTrianglesOut*sizeof(double));
-			num=0;
-			for(i=0;i<Th.nbt;i++){ 
-				Triangle &t=Th.triangles[i];
-				if (t.link){
-					//write the element only if it is part of the mesh (no boundary element)
-					TrianglesOut[4*num+0]=Th.Number(t[0])+1;
-					TrianglesOut[4*num+1]=Th.Number(t[1])+1;
-					TrianglesOut[4*num+2]=Th.Number(t[2])+1;
-					TrianglesOut[4*num+3]=t.color;
-					num=num+1;
-				}
-			}
-
-			/*clean up*/
-			delete &Th;
-		}
+		//Build new mesh Thr
+		Thr=&BTh,Thb=0;
+		Triangles & Th( *(0 ?  new Triangles(*Thr,&Thr->Gh,Thb,maxnbv) :  new Triangles(maxnbv,BTh,KeepBackVertices)));
+		if (Thr != &BTh) delete Thr;
+
+		if(costheta<=1.0){
+			Th.MakeQuadrangles(costheta);
+		}
+		if (allquad){
+			Th.SplitElement(allquad==2);
+		}
+		if(splitcorners){
+			Th.SplitInternalEdgeWithBorderVertices();
+		}
+		Th.ReNumberingTheTriangleBySubDomain();
+
+		if(NbSmooth>0) Th.SmoothingVertex(NbSmooth,omega);
+
+		Th.BTh.ReMakeTriangleContainingTheVertex();
+
+		//info
+		if(verbosity>0) {
+			if (Th.nbt-Th.NbOutT-Th.NbOfQuad*2){
+				printf("   new number of triangles = %i\n",(Th.nbt-Th.NbOutT-Th.NbOfQuad*2));
+			}
+			if (Th.NbOfQuad ){
+				printf("   new number of quads = %i\n",Th.NbOfQuad);
+			}
+		}
+
+		//Build output
+		NumVerticesOut=Th.nbv;
+		VerticesOut=(double*)xmalloc(3*NumVerticesOut*sizeof(double));
+		for (i=0;i<NumVerticesOut;i++){
+			VerticesOut[i*3+0]=Th.vertices[i].r.x;
+			VerticesOut[i*3+1]=Th.vertices[i].r.y;
+			VerticesOut[i*3+2]=Th.vertices[i].ReferenceNumber;
+		}
+
+		NumTrianglesOut=Th.nbt-Th.NbOutT; //number of triangles - number of external triangles
+		TrianglesOut=(double*)xmalloc(4*NumTrianglesOut*sizeof(double));
+		num=0;
+		for(i=0;i<Th.nbt;i++){ 
+			Triangle &t=Th.triangles[i];
+			if (t.link){
+				//write the element only if it is part of the mesh (no boundary element)
+				TrianglesOut[4*num+0]=Th.Number(t[0])+1;
+				TrianglesOut[4*num+1]=Th.Number(t[1])+1;
+				TrianglesOut[4*num+2]=Th.Number(t[2])+1;
+				TrianglesOut[4*num+3]=t.color;
+				num=num+1;
+			}
+		}
+
+		/*clean up*/
+		//delete &Th;
 		/*}}}*/
 	}
Index: /issm/trunk/src/c/Bamgx/Mesh2.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/Mesh2.cpp	(revision 2788)
+++ /issm/trunk/src/c/Bamgx/Mesh2.cpp	(revision 2789)
@@ -2629,5 +2629,5 @@
 		void Triangles::FindSubDomain(int OutSide=0)
 		{
-			long int verbosity=2;
+			long int verbosity=20;
 			//#define DRAWING1
 
@@ -4032,36 +4032,30 @@
 		}
 
-		Triangles::~Triangles() 
-		{
+		Triangles::~Triangles() {
 			long int verbosity=2;
-			assert(NbRef<=0);
+			//assert(NbRef<=0);
 			if (CurrentTh == this) CurrentTh=0;
-			if(verbosity>10)
-				cout << " ~Triangles "<< this  <<" "<< identity << endl;
-			if(vertices)  delete [] vertices;
+			//if(vertices)  delete [] vertices; //TEST  crash if not commented
 			if(edges)     delete [] edges;
 			if(triangles) delete [] triangles;
 			if(quadtree)  delete  quadtree;
-			if(ordre)     delete [] ordre;
+			//if(ordre)     delete [] ordre; //TEST  crash if not commented
 			if( subdomains) delete []  subdomains;
 			if (VerticesOnGeomEdge) delete [] VerticesOnGeomEdge;
 			if (VerticesOnGeomVertex) delete [] VerticesOnGeomVertex;
-			if (name) delete [] name;
+			//if (name) delete [] name; //TEST crash if not commented
 			if (identity) delete [] identity;
 			if (VertexOnBThVertex) delete [] VertexOnBThVertex;
 			if (VertexOnBThEdge) delete [] VertexOnBThEdge;
 
-			if (&Gh) 
-			{
+			if (&Gh) {
 				if (Gh.NbRef>0) Gh.NbRef--;
 				else if (Gh.NbRef==0) delete &Gh;
 			}
-			if (&BTh && (&BTh != this))
-			{
+			if (&BTh && (&BTh != this)) {
 				if (BTh.NbRef>0) BTh.NbRef--;
 				else if (BTh.NbRef==0) delete &BTh;
 			}
 			PreInit(0); // set all to zero 
-
 		}
 
@@ -4628,5 +4622,5 @@
 void Triangles::MakeQuadTree()
 {  
-	long int verbosity=2;
+	long int verbosity=9;
 	if(verbosity>8)
 		cout << "      MakeQuadTree" << endl;
Index: /issm/trunk/src/c/Bamgx/Mesh2.h
===================================================================
--- /issm/trunk/src/c/Bamgx/Mesh2.h	(revision 2788)
+++ /issm/trunk/src/c/Bamgx/Mesh2.h	(revision 2789)
@@ -906,5 +906,5 @@
   void Read_msh(MeshIstream &);
 
-  void ReadMetric(const char * fmetrix,const Real8 hmin,const Real8 hmax,const Real8 coef);
+  void ReadMetric(BamgArgs* bamgargs,const Real8 hmin,const Real8 hmax,const Real8 coef);
   void IntersectConsMetric(const double * s,const Int4 nbsol,const int * typsols,
 			   const  Real8 hmin,const Real8 hmax, const Real8 coef,
Index: /issm/trunk/src/c/Bamgx/MeshRead.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/MeshRead.cpp	(revision 2788)
+++ /issm/trunk/src/c/Bamgx/MeshRead.cpp	(revision 2789)
@@ -1149,5 +1149,4 @@
   SetIntCoor();
   FillHoleInMesh();
-  printf("ok0\n");
 }
 
Index: /issm/trunk/src/c/Bamgx/Metric.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/Metric.cpp	(revision 2788)
+++ /issm/trunk/src/c/Bamgx/Metric.cpp	(revision 2789)
@@ -301,5 +301,5 @@
 void  Triangles::BoundAnisotropy(Real8 anisomax,Real8 hminaniso)
 {
-	long int verbosity=0;
+	long int verbosity=10;
 
   double lminaniso = 1/ (Max(hminaniso*hminaniso,1e-100));
@@ -805,51 +805,35 @@
 
 
-void Triangles::ReadMetric(const char * fmetrix,const Real8 hmin1=1.0e-30,const Real8 hmax1=1.0e30,const Real8 coef=1)
+void Triangles::ReadMetric(BamgArgs* bamgargs,const Real8 hmin1=1.0e-30,const Real8 hmax1=1.0e30,const Real8 coef=1)
 {
-	long int verbosity=0;
-
-  Real8 hmin = Max(hmin1,MinimalHmin());
-  Real8 hmax = Min(hmax1,MaximalHmax());
-  MeshIstream f_metrix(fmetrix);
-  Int4 k,j;
-  f_metrix >>  k >> j ;
-  if(verbosity>1)
-    cout << " metrix: open " << fmetrix 
-	 << ", le coef = " << coef
-	 << ", hmin = " << hmin 
-	 << ", hmax = " << hmax 
-	 << (  (j == 1)? " Iso " : " AnIso " )<< endl;
-  
-  if (k != nbv || !(j == 1 || j == 3)) 
-    {
-      cerr << " Error Pb metrix " << k << " <> " 
-	   <<  nbv << " or  1 or 3 <> " << j << endl;
-      MeshError(1002);
-    }
-  
-  cout << " j = " << j << endl;
-  //  Int4 nberr = 0;
-  for (Int4 iv=0;iv<nbv;iv++)
-    {
-      Real8 h;
-      if (j == 1) 
-	{
-	f_metrix >>  h ;
-	vertices[iv].m=Metric(Max(hmin,Min(hmax, h*coef)));
+	int  i,j;
+
+	if(bamgargs->verbose>3) printf("      processing metric\n");
+
+	Real8 hmin = Max(hmin1,MinimalHmin());
+	Real8 hmax = Min(hmax1,MaximalHmax());
+
+	//for now we only use j==3
+	j=3;
+
+	for (i=0;i<nbv;i++){
+		Real8 h;
+		if (j == 1){
+			h=bamgargs->metric[i];
+			vertices[i].m=Metric(Max(hmin,Min(hmax, h*coef)));
+		}
+		else if (j==3){
+			Real8 a,b,c;	     
+			a=bamgargs->metric[i*3+0];
+			b=bamgargs->metric[i*3+1];
+			c=bamgargs->metric[i*3+2];
+			MetricAnIso M(a,b,c);
+			MatVVP2x2 Vp(M/coef);
+
+			Vp.Maxh(hmin);
+			Vp.Minh(hmax);
+			vertices[i].m = Vp;
+		}
 	}
-      else if (j==3) 
-	{
-	  Real8 a,b,c;	     
-	  f_metrix >>  a >> b >> c  ;
-	  MetricAnIso M(a,b,c);
-	  MatVVP2x2 Vp(M/coef);
-	  
-	  Vp.Maxh(hmin);
-	  Vp.Minh(hmax);
-	  vertices[iv].m = Vp;
-	  
-	}
-    }
- 
 }
 
@@ -876,5 +860,5 @@
 void  Triangles::MaxSubDivision(Real8 maxsubdiv)
 {
-	long int verbosity=0;
+	long int verbosity=10;
 
 const  Real8 maxsubdiv2 = maxsubdiv*maxsubdiv;
@@ -958,5 +942,5 @@
 void Triangles::SmoothMetric(Real8 raisonmax) 
 { 
-	long int verbosity=0;
+	long int verbosity=10;
 
   if(raisonmax<1.1) return;
Index: /issm/trunk/src/m/classes/public/bamg.m
===================================================================
--- /issm/trunk/src/m/classes/public/bamg.m	(revision 2788)
+++ /issm/trunk/src/m/classes/public/bamg.m	(revision 2789)
@@ -60,5 +60,5 @@
 	bamg_mesh.Vertices=[md.x md.y ones(md.numberofgrids,1)];
 	bamg_mesh.NumTriangles=md.numberofelements;
-	bamg_mesh.Triangles=[md.elements ones(md.numberofgrids,1)];
+	bamg_mesh.Triangles=[md.elements ones(md.numberofelements,1)];
 end
 %}}}
@@ -69,5 +69,5 @@
 bamg_options.hmin=getfieldvalue(options,'hmin',10^-100);
 bamg_options.hmax=getfieldvalue(options,'hmax',10^100);
-bamg_options.gradation=getfieldvalue(options,'gradation',1.2);
+bamg_options.gradation=getfieldvalue(options,'gradation',1.5);
 bamg_options.cutoff=getfieldvalue(options,'cutoff',10^-5);
 bamg_options.verbose=getfieldvalue(options,'verbose',1);
Index: /issm/trunk/src/mex/Bamg/Bamg.cpp
===================================================================
--- /issm/trunk/src/mex/Bamg/Bamg.cpp	(revision 2788)
+++ /issm/trunk/src/mex/Bamg/Bamg.cpp	(revision 2789)
@@ -75,7 +75,7 @@
 	FetchData(&VerticesMesh,NULL,NULL,mxGetField(BAMGMESH,0,"Vertices"));
 	bamgmesh.Vertices=VerticesMesh;
-	FetchData(&NumTrianglesMesh,mxGetField(BAMGMESH,0,"NumVertices"));
+	FetchData(&NumTrianglesMesh,mxGetField(BAMGMESH,0,"NumTriangles"));
 	bamgmesh.NumTriangles=NumTrianglesMesh;
-	FetchData(&TrianglesMesh,NULL,NULL,mxGetField(BAMGMESH,0,"Vertices"));
+	FetchData(&TrianglesMesh,NULL,NULL,mxGetField(BAMGMESH,0,"Triangles"));
 	bamgmesh.Triangles=TrianglesMesh;
 	bamgmesh.hVertices=NULL;
@@ -123,6 +123,6 @@
 	xfree((void**)&hVerticesGeom);
 	xfree((void**)&SubDomainGeom);
-	xfree((void**)&TrianglesMesh);
-	xfree((void**)&VerticesMesh);
+	//xfree((void**)&TrianglesMesh);
+	//xfree((void**)&VerticesMesh);
 
 	/*end module: */