Index: /issm/trunk/src/c/Bamgx/Mesh2.h
===================================================================
--- /issm/trunk/src/c/Bamgx/Mesh2.h	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/Mesh2.h	(revision 2850)
@@ -791,5 +791,4 @@
 	     catch(...) { this->~Triangles(); throw; } }
   Triangles(Triangles &,Geometry * pGh=0,Triangles* pBTh=0,Int4 nbvxx=0 ); // COPY OPERATEUR
-  //  Triangles(Triangles &){ std::cerr << " BUG call copy opretor of Triangles" << std::endl;MeshError(111);}
   Triangles(const Triangles &,const int *flag,const int *bb); // truncature
 
Index: /issm/trunk/src/c/Bamgx/objects/GeometricalEdge.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/GeometricalEdge.cpp	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/objects/GeometricalEdge.cpp	(revision 2850)
@@ -8,4 +8,7 @@
 #include "../QuadTree.h"
 #include "../SetOfE4.h"
+
+#undef __FUNCT__ 
+#define __FUNCT__ "GeometricalEdge"
 
 using namespace std;
Index: /issm/trunk/src/c/Bamgx/objects/Geometry.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/Geometry.cpp	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/objects/Geometry.cpp	(revision 2850)
@@ -11,4 +11,7 @@
 #include "../../include/macros.h"
 #include "../../toolkits/toolkits.h"
+
+#undef __FUNCT__ 
+#define __FUNCT__ "Geometry"
 
 namespace bamg {
@@ -540,5 +543,4 @@
 		const Vertex &v0=e[0],&v1=e[1];
 		V.m = Metric(1.0-s, v0,s, v1);
-#define MXE__LINE  __LINE__+1
 		const int mxe =100;
 		GeometricalEdge *ge[mxe+1];
@@ -566,18 +568,14 @@
 			<<  " sg 1= " << vg1 
 			<< "--------------------------------------------" << endl;
-		while (eg0 != (GeometricalEdge*) vg0  &&  (*eg0)(sens0) != (GeometricalVertex*) vg0)
-		  { 
+		while (eg0 != (GeometricalEdge*) vg0  &&  (*eg0)(sens0) != (GeometricalVertex*) vg0){ 
 			if (bge<=0) {
 				//          int kkk;
-				// if (NbTry) cout <<"Read (int) to Show Sioux window", cin >> kkk ;
-				if(NbTry) 
-				  {
-					cerr << " -- Fatal Error: on the class triangles before call Geometry::ProjectOnCurve" << endl; 
-					cerr << "   The mesh of the  Geometry is to fine: ";
-					cerr << "     1)  a mesh edge  contening more than "<< mxe/2 << " geometrical edges." << endl;
-					cerr << "     2)  code bug : be sure that we call   Triangles::SetVertexFieldOn() before " << endl;
-					cerr << "   To solve the problem do a coarsening of the geometrical mesh " << endl;
-					cerr << " or change the constant value of mxe in " << __FILE__ << " line " << MXE__LINE << "( dangerous way )" << endl;	  
-					MeshError(222);
+				if(NbTry) {
+					printf("Fatal Error: on the class triangles before call Geometry::ProjectOnCurve\n");
+					printf("That bug might come from:\n");
+					printf(" 1)  a mesh edge  contening more than %i geometrical edges\n",mxe/2);
+					printf(" 2)  code bug : be sure that we call   Triangles::SetVertexFieldOn() before\n");
+					printf("To solve the problem do a coarsening of the geometrical mesh or change the constant value of mxe (dangerous)\n");
+					throw ErrorException(__FUNCT__,exprintf("see above"));
 				  }
 				NbTry++;
@@ -603,10 +601,10 @@
 				NbTry++;
 				if (NbTry<2) goto retry;
-				cerr << "   The mesh of the  Geometry is to fine:" ;
-				cerr << "     1)  a mesh edge  contening more than "<< mxe/2 << " geometrical edges." << endl;
-				cerr << "     2)  code bug : be sure that we call   Triangles::SetVertexFieldOn() before " << endl;
-				cerr << "   To solve the problem do a coarsening of the geometrical mesh " << endl;
-				cerr << " or change the constant value of mxe in " << __FILE__ << " line " << MXE__LINE << "( dangerous way )" << endl;	
-				MeshError(223);
+				printf("Fatal Error: on the class triangles before call Geometry::ProjectOnCurve\n");
+				printf("That bug might come from:\n");
+				printf(" 1)  a mesh edge  contening more than %i geometrical edges\n",mxe/2);
+				printf(" 2)  code bug : be sure that we call   Triangles::SetVertexFieldOn() before\n");
+				printf("To solve the problem do a coarsening of the geometrical mesh or change the constant value of mxe (dangerous)\n");
+				throw ErrorException(__FUNCT__,exprintf("see above"));
 			}
 
@@ -702,73 +700,56 @@
 		Int4 * ev = new Int4 [ 2 * nbe ];
 		float  * eangle = new float[ nbe ];
+		double eps = 1e-20;
+		QuadTree quadtree; // to find same vertices
+		Vertex * v0 = vertices; 
+		GeometricalVertex  * v0g = (GeometricalVertex  *) (void *) v0;   
+
+		for (i=0;i<nbv;i++) 
+		 vertices[i].link = vertices +i;
+		for (i=0;i<nbv;i++) 
 		  {
-			double eps = 1e-20;
-			QuadTree quadtree; // to find same vertices
-			Vertex * v0 = vertices; 
-			GeometricalVertex  * v0g = (GeometricalVertex  *) (void *) v0;   
-			int k=0;
-			for (i=0;i<nbv;i++) 
-			 vertices[i].link = vertices +i;
-			for (i=0;i<nbv;i++) 
-			  {
-				vertices[i].i = toI2(vertices[i].r); // set integer coordinate
-				Vertex *v= quadtree.NearestVertex(vertices[i].i.x,vertices[i].i.y); 
-				if( v && Norme1(v->r - vertices[i]) < eps )
-				  { // link v & vertices[i] 
-					// vieille ruse pour recuperer j 
-					GeometricalVertex * vg = (GeometricalVertex  *) (void *) v;
-					int j = vg-v0g;
-					assert( v ==  & (Vertex &) vertices[j]);
-					vertices[i].link = vertices + j;
-					k++;	      
+			vertices[i].i = toI2(vertices[i].r); // set integer coordinate
+			Vertex *v= quadtree.NearestVertex(vertices[i].i.x,vertices[i].i.y); 
+			if( v && Norme1(v->r - vertices[i]) < eps )
+			  { // link v & vertices[i] 
+				// vieille ruse pour recuperer j 
+				GeometricalVertex * vg = (GeometricalVertex  *) (void *) v;
+				int j = vg-v0g;
+				assert( v ==  & (Vertex &) vertices[j]);
+				vertices[i].link = vertices + j;
+				k++;	      
+			  }
+			else  quadtree.Add(vertices[i]); 
+		  }
+		if (k) {
+			printf("number of distinct vertices= %i, over %i\n",nbv - k,nbv);
+			printf("List of duplicate vertices:\n");
+			for (i=0;i<nbv;i++){
+				if (!vertices[i].IsThe()) printf("  %i and %i\n",i,Number(vertices[i].The()));
+			}
+			throw ErrorException(__FUNCT__,exprintf("See above"));
+		}
+
+		//  verification of cracked edge
+		for (i=0;i<nbe;i++){
+			if (edges[i].Cracked() ) {
+				//    verification of crack
+				GeometricalEdge & e1=edges[i];
+				GeometricalEdge & e2=*e1.link;
+				cerr << i << " " << e1[0].The() << " " << e2[0].The() << " " <<  e1[1].The() << " " << e2[1].The() << endl;
+				if ( e1[0].The() == e2[0].The() && e1[1].The() == e2[1].The() )
+				  {
 				  }
-				else  quadtree.Add(vertices[i]); 
-			  }
-			if (k) {
-				cout << " Number of distinte vertices " << nbv - k << " Over " << nbv << endl;
-				//if (verbosity>10) 
-				  {
-					cout << " The duplicate vertex " << endl;
-					for (i=0;i<nbv;i++)
-					 if (!vertices[i].IsThe())
-					  cout << " " << i << " and " << Number(vertices[i].The()) << endl;
-					MeshError(102);
-					//throw(ErrorExec("exit",1));    
-				  }
-			}
-
-			//  verification of cracked edge
-			int err =0;
-			for (i=0;i<nbe;i++)
-			 if (edges[i].Cracked() )
-				{
-				 //    verification of crack
-				 GeometricalEdge & e1=edges[i];
-				 GeometricalEdge & e2=*e1.link;
-				 cerr << i << " " << e1[0].The() << " " << e2[0].The() << " " <<  e1[1].The() << " " << e2[1].The() << endl;
-				 if ( e1[0].The() == e2[0].The() && e1[1].The() == e2[1].The() )
+				else 
+				 if ( e1[0].The() == e2[1].The() && e1[1].The() == e2[0].The() )
 					{
+					 //nothing
 					}
-				 else 
-				  if ( e1[0].The() == e2[1].The() && e1[1].The() == e2[0].The() )
-					 {
-					 }
-				  else
-					 {
-					  err++;
-					  cerr << " Cracked edges with no same vertex " << &e1-edges << " " << &e2 -edges << endl;
-					 }
-				}
-			 else
-				{
-				 //  if (!edges[i][0].IsThe()) err++;
-				 // if (!edges[i][1].IsThe()) err++;
-				}
-			if (err)
-			  {
-				cerr << " Some vertex was not distint and not on cracked edge " << err<< endl;
-				MeshError(222);
-			  }
-		  }
+				 else {
+					 throw ErrorException(__FUNCT__,exprintf("Cracked edges with no same vertex")); 
+				 }
+			}
+		}
+
 		if(verbosity>7) 
 		 for (i=0;i<nbv;i++)
@@ -776,24 +757,22 @@
 			cout << "     The geo vertices  " << i << " is required" << endl;
 
-		for (i=0;i<nbv;i++) 
-		 hv[i]=-1;// empty list
-
-		for (i=0;i<nbe;i++) 
-		  {
+		for (i=0;i<nbv;i++) hv[i]=-1;// empty list
+
+		for (i=0;i<nbe;i++) {
 			R2 v10  =  edges[i].v[1]->r -  edges[i].v[0]->r;
 			Real8 lv10 = Norme2(v10);
 			if(lv10 == 0) {
-				cerr << "The length  of " <<i<< "th Egde is 0 " << endl ;
-				MeshError(1);}
-				eangle[i] = atan2(v10.y,v10.x)  ; // angle in [ -Pi,Pi ]
-				if(verbosity>9) 
-				 cout << "     angle edge " << i <<" " << eangle[i]*180/Pi<< v10<<endl;
-				for (jj=0;jj<2;jj++)
-				  { // generation of list
-					Int4 v =  Number(edges[i].v[jj]);
-					ev[k] = hv[v];
-					hv[v] = k++;
-				  }
-		  }
+				throw ErrorException(__FUNCT__,exprintf("Length of edge %i is 0",i));
+			}
+			eangle[i] = atan2(v10.y,v10.x)  ; // angle in [ -Pi,Pi ]
+			if(verbosity>9) 
+			 cout << "     angle edge " << i <<" " << eangle[i]*180/Pi<< v10<<endl;
+			for (jj=0;jj<2;jj++)
+			  { // generation of list
+				Int4 v =  Number(edges[i].v[jj]);
+				ev[k] = hv[v];
+				hv[v] = k++;
+			  }
+		}
 		// bulle sort on the angle of edge  
 		for (i=0;i<nbv;i++) {
@@ -820,14 +799,4 @@
 			} // end while (exch)
 
-			if (ord >= 1 ) 
-			  { /*
-					 Int4 n = hv[i];
-					 while ( n >=0) 
-					 { Int4 i1 = n/2,j1 = n%2;
-				//float a = 180*(j1 ? OppositeAngle(eangle[i1]): eangle[i1])/Pi;
-				n = ev[n];
-				}
-				*/
-			  } 
 			if(ord == 2) { // angulare test to find a corner 
 				Int4 n1 = hv[i];
@@ -882,14 +851,7 @@
 			 Int4 n1 = ev[k++]; 
 			 Int4 i1 = n1/2 ,j1=n1%2;
-			 if( edges[i1].v[j1] != edges[i].v[jj]) 
-				{ cerr << " Bug Adj edge " << i << " " << jj << 
-				 " et " << i1 << " " << j1 << " k=" << k;
-				 cerr << Number(edges[i].v[jj]) <<" <> " 
-					<< Number(edges[i1].v[j1])  <<endl;
-				 cerr << "edge " << Number(edges[i].v[0]) << " " 
-					<< Number(edges[i].v[1]) << endl; 
-				 //    cerr << "in file " <<filename <<endl;
-				 MeshError(1);
-				}
+			 if( edges[i1].v[j1] != edges[i].v[jj]) {
+				 throw ErrorException(__FUNCT__,exprintf("Bug Adj edge"));
+			 }
 			 edges[i1].Adj[j1] = edges + i;
 			 edges[i1].SensAdj[j1] = jj;
Index: /issm/trunk/src/c/Bamgx/objects/ListofIntersectionTriangles.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/ListofIntersectionTriangles.cpp	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/objects/ListofIntersectionTriangles.cpp	(revision 2850)
@@ -11,4 +11,7 @@
 #include "../../include/macros.h"
 #include "../../toolkits/toolkits.h"
+
+#undef __FUNCT__ 
+#define __FUNCT__ "ListofIntersectionTriangles"
 
 namespace bamg {
@@ -181,15 +184,8 @@
 						j= VerticesOfTriangularEdge[ocut][1];
 					}
-					else  {
-						cerr << " Bug Split Edge " << endl;
-						cerr << " dt[0]= " << dt[0] 
-						  << " dt[1]= " << dt[1] 
-						  << " dt[2]= "<< dt[2] << endl;
-						cerr << i0 << " " << i1 << " " << i2 << endl;
-						cerr << " A = " << A << " B= " << B << endl;
-						cerr << " Triangle t = " <<  *t << endl;
-						cerr << (*t)[0] << (*t)[1] << (*t)[0] << endl;
-						cerr << " nbt = " << nbt << endl;
-						MeshError(100);}}
+					else {
+						throw ErrorException(__FUNCT__,exprintf("Bug Split Edge"));
+					}
+				}
 
 						k = OppositeVertex[ocut];
Index: /issm/trunk/src/c/Bamgx/objects/MatVVP2x2.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/MatVVP2x2.cpp	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/objects/MatVVP2x2.cpp	(revision 2850)
@@ -11,4 +11,7 @@
 #include "../../include/macros.h"
 #include "../../toolkits/toolkits.h"
+
+#undef __FUNCT__ 
+#define __FUNCT__ "MatVVP2x2"
 
 namespace bamg {
Index: /issm/trunk/src/c/Bamgx/objects/MetricAnIso.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/MetricAnIso.cpp	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/objects/MetricAnIso.cpp	(revision 2850)
@@ -8,4 +8,7 @@
 #include "../QuadTree.h"
 #include "../SetOfE4.h"
+
+#undef __FUNCT__ 
+#define __FUNCT__ "MetricAnIso"
 
 using namespace std;
Index: /issm/trunk/src/c/Bamgx/objects/QuadTree.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/QuadTree.cpp	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/objects/QuadTree.cpp	(revision 2850)
@@ -5,4 +5,7 @@
 #include "../Mesh2.h"
 #include "../QuadTree.h"
+
+#undef __FUNCT__ 
+#define __FUNCT__ "QuadTree"
 
 namespace bamg {
Index: /issm/trunk/src/c/Bamgx/objects/SetOfE4.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/SetOfE4.cpp	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/objects/SetOfE4.cpp	(revision 2850)
@@ -1,5 +1,12 @@
 #include <iostream> 
+
+#include "../../shared/shared.h"
+#include "../../include/macros.h"
+#include "../../toolkits/toolkits.h"
 #include "../meshtype.h"
 #include "../SetOfE4.h"
+
+#undef __FUNCT__ 
+#define __FUNCT__ "SetOfEdges4"
 
 using namespace std;
@@ -24,13 +31,12 @@
 	Int4 SetOfEdges4::find(Int4 ii,Int4 jj) { 
 		if (tete == 0 ) {
-			cerr <<"SetOfEdges4::find \nplus de tete de liste\n";
-			MeshError(888);}
-			Int4 n = tete[ Abs( ii ) % nx ];
+			throw ErrorException(__FUNCT__,exprintf("SetOfEdges4::find no more tete de liste (?)"));
+		}
+		Int4 n = tete[ Abs( ii ) % nx ];
 
-			while (n >= 0) 
-			 if (ii == Edges[n].i && jj == Edges[n].j)
-			  return n;
-			 else n = Edges[n].next;
-			return -1; // n'existe pas
+		while (n >= 0) 
+		 if (ii == Edges[n].i && jj == Edges[n].j) return n;
+		 else n = Edges[n].next;
+		return -1; //do not exist, return -1
 	}
 	/*}}}1*/
@@ -38,22 +44,21 @@
 	Int4 SetOfEdges4::add(Int4 ii,Int4 jj) {
 		if (tete == 0 ) {
-			cerr << "SetOfEdges4::add\n plus de tete de liste \n" << endl;
-			MeshError(888);}
+			throw ErrorException(__FUNCT__,exprintf("SetOfEdges4::add no more tete de liste (?)"));
+		}
+		Int4 h;
+		Int4 n = tete[ h = Abs( ii ) % nx ];
+		while (n >= 0) 
+		 if (ii == Edges[n].i && jj == Edges[n].j)
+		  return n;
+		 else n = Edges[n].next;
+		if (nbax <=NbOfEdges ) {
+			throw ErrorException(__FUNCT__,exprintf("SetOfEdges4::add overflow: NbOfEdges=%i > nbax=%i",NbOfEdges,nbax));
+		}
 
-			Int4 h;
-			Int4 n = tete[ h = Abs( ii ) % nx ];
-			while (n >= 0) 
-			 if (ii == Edges[n].i && jj == Edges[n].j)
-			  return n;
-			 else n = Edges[n].next;
-			if (nbax <=NbOfEdges ) {
-				cerr << " SetOfEdges4::add\noverflow de la pile "  << nbax << " " << NbOfEdges << endl;
-				MeshError(888);}
-
-				Edges[NbOfEdges].i=ii;
-				Edges[NbOfEdges].j=jj;
-				Edges[NbOfEdges].next= tete[h];
-				tete[h] = NbOfEdges;
-				return NbOfEdges ++;
+		Edges[NbOfEdges].i=ii;
+		Edges[NbOfEdges].j=jj;
+		Edges[NbOfEdges].next= tete[h];
+		tete[h] = NbOfEdges;
+		return NbOfEdges ++;
 	}
 	/*}}}1*/
Index: /issm/trunk/src/c/Bamgx/objects/Triangle.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/Triangle.cpp	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/objects/Triangle.cpp	(revision 2850)
@@ -11,4 +11,7 @@
 #include "../../include/macros.h"
 #include "../../toolkits/toolkits.h"
+
+#undef __FUNCT__ 
+#define __FUNCT__ "Triangle"
 
 namespace bamg {
Index: /issm/trunk/src/c/Bamgx/objects/Triangles.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/Triangles.cpp	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/objects/Triangles.cpp	(revision 2850)
@@ -11,4 +11,7 @@
 #include "../../include/macros.h"
 #include "../../toolkits/toolkits.h"
+
+#undef __FUNCT__ 
+#define __FUNCT__ "Triangles"
 
 namespace bamg {
@@ -110,6 +113,5 @@
 		  assert(kt==nbt);
 		  if (nbt ==0 && nbv ==0) {
-			  cout << "Error all triangles was remove " << endl;
-			  MeshError(999,this);
+			  throw ErrorException(__FUNCT__,exprintf("All triangles have been removed"));
 		  }
 		  delete [] kk;
@@ -733,7 +735,8 @@
 		// --------------------------
 		long int verbosity=0;
-		if (verbosity>1) 
-		 cout << "  -- construction of the geometry from the 2d mesh " << endl;
-		if (nbt<=0 || nbv <=0 ) { MeshError(101);}
+		if (verbosity>1) printf("   construction of the geometry from the 2d mesh\n");
+		if (nbt<=0 || nbv <=0 ) {
+			throw ErrorException(__FUNCT__,exprintf("nbt or nbv is negative"));
+		}
 
 		// construction of the edges 
@@ -765,43 +768,30 @@
 			edge4->addtrie(Number(edges[i][0]),Number(edges[i][1]));
 		  }
-		if (nbe !=  edge4->nb())
-		  { 
-			cerr << " Some Double edge in the mesh, the number is " << nbe 
-			  << " nbe4=" << edge4->nb()  << endl;
-			MeshError(1002);
+		if (nbe !=  edge4->nb()){ 
+			throw ErrorException(__FUNCT__,exprintf("Some Double edge in the mesh, the number is %i, nbe4=%i",nbe,edge4->nb())); 
 		  }
-		for (i=0;i<nbt;i++)
-		 for  (j=0;j<3;j++)
-			{
-			 // Int4 i0,i1;
-			 Int4 k =edge4->addtrie(Number(triangles[i][VerticesOfTriangularEdge[j][0]]),
-						 Number(triangles[i][VerticesOfTriangularEdge[j][1]]));
-			 Int4 invisible = triangles[i].Hidden(j);
-			 if(st[k]==-1)
-			  st[k]=3*i+j;
-			 else if(st[k]>=0) {
-				 assert( ! triangles[i].TriangleAdj(j) && !triangles[st[k] / 3].TriangleAdj((int) (st[k]%3)));
-
-				 triangles[i].SetAdj2(j,triangles + st[k] / 3,(int) (st[k]%3));
-				 if (invisible)  triangles[i].SetHidden(j);
-				 if (k<nbe) {
-					 triangles[i].SetLocked(j);
-				 }
-				 st[k]=-2-st[k]; }
-			 else {
-				 cerr << " The edge (" 
-					<< Number(triangles[i][VerticesOfTriangularEdge[j][0]])
-					<< " , " 
-					<< Number(triangles[i][VerticesOfTriangularEdge[j][1]])
-					<< " ) is in more than 2 triangles " <<k <<endl;
-				 cerr << " Edge " << j << " Of Triangle " << i << endl;
-				 cerr << " Edge " << (-st[k]+2)%3 << " Of Triangle " << (-st[k]+2)/3  << endl;
-				 cerr << " Edge " << triangles[(-st[k]+2)/3].NuEdgeTriangleAdj((int)((-st[k]+2)%3))
-					<< " Of Triangle " <<  Number(triangles[(-st[k]+2)/3].TriangleAdj((int)((-st[k]+2)%3))) 
-					<< endl;
-				 MeshError(9999);}	
-
-
+		for (i=0;i<nbt;i++){
+			for  (j=0;j<3;j++) {
+				Int4 k =edge4->addtrie(Number(triangles[i][VerticesOfTriangularEdge[j][0]]), Number(triangles[i][VerticesOfTriangularEdge[j][1]]));
+				Int4 invisible = triangles[i].Hidden(j);
+				if(st[k]==-1) st[k]=3*i+j;
+				else if(st[k]>=0) {
+					assert( ! triangles[i].TriangleAdj(j) && !triangles[st[k] / 3].TriangleAdj((int) (st[k]%3)));
+					triangles[i].SetAdj2(j,triangles + st[k] / 3,(int) (st[k]%3));
+					if (invisible)  triangles[i].SetHidden(j);
+					if (k<nbe) {
+						triangles[i].SetLocked(j);
+					}
+					st[k]=-2-st[k]; 
+				}
+				else {
+					printf("The edge (%i,%i) belongs to more than 2 triangles (%i)\n",Number(triangles[i][VerticesOfTriangularEdge[j][0]]),Number(triangles[i][VerticesOfTriangularEdge[j][1]]),k);
+					printf("Edge %i of triangle %i\n",j,i);
+					printf("Edge %i of triangle %i\n",(-st[k]+2)%3,(-st[k]+2)/3);
+					printf("Edge %i of triangle %i\n",triangles[(-st[k]+2)/3].NuEdgeTriangleAdj((int)((-st[k]+2)%3)),Number(triangles[(-st[k]+2)/3].TriangleAdj((int)((-st[k]+2)%3))));
+					throw ErrorException(__FUNCT__,exprintf("An edge belongs to more than 2 triangles"));
+				}	
 			}
+		}
 		Int4 nbedges = edge4->nb(); // the total number of edges 
 		delete edge4;
@@ -809,13 +799,11 @@
 
 		if(verbosity>5) {
-			if (name)
-			 cout << "    On Mesh " << name << endl;
-			cout << "    - The number of Vertices  = " << nbv << endl;
-			cout << "    - The number of Triangles = " << nbt << endl;
-			cout << "    - The number of given edge = " << nbe << endl;
-			cout << "    - The number of all edges = " << nbedges << endl;
-			cout << "    - The Euler number = 1-Nb Of Hole = " << nbt-nbedges+nbv << endl; }
-
-
+			printf("         info on Mesh %s:\n",name);
+			printf("            - number of vertices    = %i \n",nbv); 
+			printf("            - number of triangles   = %i \n",nbt); 
+			printf("            - number of given edges = %i \n",nbe); 
+			printf("            - number of all edges   = %i \n"  ,edge4->nb()); 
+			printf("            - Euler number 1 - nb of holes = %i \n"  ,nbt-edge4->nb()+nbv); 
+		}
 			// check the consistant of edge[].adj and the geometrical required  vertex
 			k=0;
@@ -824,6 +812,5 @@
 
 			for (i=0;i<nbedges;i++)
-			 if (st[i] <-1) // edge internal
-				{ 
+			 if (st[i] <-1) {// edge internal
 				 it =  (-2-st[i])/3;
 				 j  =  (int) ((-2-st[i])%3);
@@ -832,5 +819,5 @@
 				 if (triangles[it].color != tt.color|| i < nbeold) // Modif FH 06122055 // between 2 sub domai
 				  k++;
-				}
+			 }
 			 else if (st[i] >=0) // edge alone 
 			  // if (i >= nbeold) 
@@ -842,11 +829,5 @@
 			k += kk;
 			kk=0;
-			if (k)
-			  {
-
-				//      if (nbe) {
-				//	cerr << k << " boundary edges  are not defined as edges " << endl;
-				//	MeshError(9998);
-				// }
+			if (k) {
 				// construction of the edges 
 				nbe = k;
@@ -1138,5 +1119,5 @@
 				  }
 				else
-				 MeshError(103);
+				 throw ErrorException(__FUNCT__,exprintf("%i should be >=0"));
 			  }
 
@@ -1157,39 +1138,25 @@
 		R2 A=vA,B=vB;
 		Vertex * pvA=&vA, * pvB=&vB;
-		if (vA.vint == IsVertexOnVertex)
-		  {
-			//  cout << " Debut vertex = " << BTh.Number(vA.onbv) ;
+		if (vA.vint == IsVertexOnVertex){
 			pA=vA.onbv;
-		  }
-		else if (vA.vint == IsVertexOnEdge)
-		  {
+		}
+		else if (vA.vint == IsVertexOnEdge){
 			pA=vA.onbe->be;
 			tA=vA.onbe->abcisse;
-			// cout << " Debut edge = " << BTh.Number(vA.onbv) << " " << tA ;
-
-		  }
-		else
-		  {cerr << "ProjectOnCurve On Vertex " << BTh.Number(vA) << " " << endl;
-			cerr << " forget call to SetVertexFieldOnBTh" << endl;
-			MeshError(-1);
-		  } 
-
-		if (vB.vint == IsVertexOnVertex)
-		  {
-			// cout << " Fin vertex = " << BTh.Number(vB.onbv) << endl;
+		}
+		else {
+			throw ErrorException(__FUNCT__,exprintf("ProjectOnCurve On Vertex %i forget call to SetVertexFieldOnBTh",BTh.Number(vA)));
+		} 
+
+		if (vB.vint == IsVertexOnVertex){
 			pB=vB.onbv;
-		  }
-		else if(vB.vint == IsVertexOnEdge)
-		  {
+		}
+		else if(vB.vint == IsVertexOnEdge){
 			pB=vB.onbe->be;
 			tB=vB.onbe->abcisse;
-			// cout << " Fin edge = " << BTh.Number(vB.onbe->be) << " " <<  tB ;
-
-		  }
-		else
-		  {cerr << "ProjectOnCurve On Vertex " << BTh.Number(vB) << " " << endl;
-			cerr << " forget call to SetVertexFieldOnBTh" << endl;
-			MeshError(-1);
-		  } 
+		}
+		else {
+			throw ErrorException(__FUNCT__,exprintf("ProjectOnCurve On Vertex %i forget call to SetVertexFieldOnBTh",BTh.Number(vB)));
+		} 
 		Edge * e = &BhAB;
 		assert( pA && pB && e);
@@ -1233,6 +1200,6 @@
 		Real8 abscisse = -1;
 
-		for (int cas=0;cas<2;cas++)
-		  {// 2 times algo:
+		for (int cas=0;cas<2;cas++){
+			// 2 times algo:
 			//    1 for computing the length l
 			//    2 for find the vertex 
@@ -1289,12 +1256,4 @@
 					assert(thetab>=0 && thetab<=1);
 					BR = VertexOnEdge(&R,eee,thetab);
-
-					//	cout << kkk << " eee = " << BTh.Number(eee) << "  v0=  " 
-					//     << BTh.Number(v0) << " " << te0
-					//      << " v1 = " << BTh.Number(v1) <<  " " << tB  << endl;
-
-					//out << Number(R) << " Opt  = " <<  thetab << " on  " <<  BTh.Number(eee) 
-					//	    << " = " << R << endl;
-
 					return  Gh.ProjectOnCurve(*eee,thetab,R,GR);
 				  }
@@ -1303,8 +1262,6 @@
 
 		  }
-		cerr << " Big Bug" << endl;
-		MeshError(678);
+		throw ErrorException(__FUNCT__,exprintf("Big bug..."));
 		return 0; // just for the compiler 
-
 	}                  
 	/*}}}1*/
@@ -1529,7 +1486,5 @@
 		//  
 
-		for (i=0;i<nbt;i++)
-		  {  
-
+		for (i=0;i<nbt;i++) {
 			Triangle & t = triangles[i];
 			assert(t.link);
@@ -1553,37 +1508,26 @@
 					Vertex &A=vertices[ks];
 					Real8 aa,bb,cc,dd;
-					if ((dd=Area2(v0.r,v1.r,A.r)) >=0)
-					  { // warning PB roundoff error 
+					if ((dd=Area2(v0.r,v1.r,A.r)) >=0){
+						// warning PB roundoff error 
 						if (t.link && ( (aa=Area2( A.r    , t[1].r , t[2].r )) < 0.0 
 										||   (bb=Area2( t[0].r , A.r    , t[2].r )) < 0.0  
-										||   (cc=Area2( t[0].r , t[1].r , A.r    )) < 0.0))
-						 ferr++, cerr << " Error : " <<  ke + nbvold << " not in triangle " 
-							<< i << " In=" << !!t.link
-							<< " " <<  aa  << " " << bb << " " << cc << " " << dd << endl;
-
-					  }
-
-					else
-					  {
+										||   (cc=Area2( t[0].r , t[1].r , A.r    )) < 0.0)){
+							printf("%i not in triangle %i In= %i %i %i %i %i\n",ke + nbvold,i,!!t.link,aa,bb,cc,dd);
+							throw ErrorException(__FUNCT__,exprintf("Number of triangles with P2 interpolation Problem"));
+						}
+					}
+					else {
 						if (tt.link && ( (aa=Area2( A.r     , tt[1].r , tt[2].r )) < 0 
 										||   (bb=Area2( tt[0].r , A.r     , tt[2].r )) < 0 
-										||   (cc=Area2( tt[0].r , tt[1].r , A.r     )) < 0)) 
-						 ferr++, cerr << " Warning : " <<  ke + nbvold << " not in triangle " << ii 
-							<< " In=" << !!tt.link 
-							<< " " <<  aa  << " " << bb << " " << cc << " " << dd << endl;
-
-					  } 
-
+										||   (cc=Area2( tt[0].r , tt[1].r , A.r     )) < 0)){
+							printf("%i not in triangle %i In= %i %i %i %i %i\n",ke + nbvold,ii,!!tt.link,aa,bb,cc,dd);
+							throw ErrorException(__FUNCT__,exprintf("Number of triangles with P2 interpolation Problem"));
+						}
+					} 
 				  }
 			  }
-		  }
-		if(ferr)
-		  {
-			cerr << " Number of triangles with P2 interpolation Probleme " << ferr << endl;;
-			MeshError(9);
-		  }
-
-		for (i=0;i<nbt;i++)
-		  {
+		}
+
+		for (i=0;i<nbt;i++){
 			ksplit[i]=1; // no split by default
 			const Triangle & t = triangles[ i];
@@ -1966,8 +1910,4 @@
 						 // angle b12 > angle ba2 => cotg(angle b12) < cotg(angle ba2)
 						 OnSwap =  ((double) cosb12 * (double)  sinba2) <  ((double) cosba2 * (double) sinb12);
-						 //  	 if(CurrentTh) 
-						 //  	   cout << "swap " << CurrentTh->Number(sa) << " " << CurrentTh->Number(sb) << " " ;
-						 //  	 cout <<  cosb12 << " " <<  sinba2 << " "  <<  cosba2 << " " << sinb12 
-						 //  	      << " Onswap = " <<  OnSwap << endl;
 						 break;
 					 }
@@ -2073,11 +2013,7 @@
 		Icoor2 detOld = t->det;
 
-		if ( ( infv <0 ) && (detOld <0) ||  ( infv >=0  ) && (detOld >0) ) 
-		  {
-			cerr << "  infv " << infv << " det = " << detOld << endl;
-			cerr << Number(s) << " "<< Number(s0) << " "  
-			  << Number(s1) << " "  << Number(s2) << endl;
-			MeshError(3);
-		  }
+		if (( infv <0 ) && (detOld <0) ||  ( infv >=0  ) && (detOld >0) ){
+			throw ErrorException(__FUNCT__,exprintf("infv=%g det=%g"));
+		}
 
 		// if det3 do not exist then constuct det3
@@ -2114,68 +2050,63 @@
 					 }}
 				 else {
-					 cerr << " bug  " << nbd0 <<endl;
-					 cerr << " Bug double points in " << endl ;
-					 cerr << " s = " << Number(s) << " " <<  s << endl;
-					 cerr << " s0 = "<< Number(s0) << " "  << s0 << endl;
-					 cerr << " s1 = "<< Number(s1) << " "  << s1 << endl;
-					 cerr << " s2 = "<< Number(s2) << " "  << s2 << endl;
-					 MeshError(5,this);}
-
-					 // remove de MarkUnSwap edge
-					 t->SetUnMarkUnSwap(0);     
-					 t->SetUnMarkUnSwap(1);     
-					 t->SetUnMarkUnSwap(2);
-
-					 tt[0]= t;
-					 tt[1]= &triangles[nbt++];
-					 tt[2]= &triangles[nbt++];
-
-					 if (nbt>nbtx) {
-						 cerr << " No enougth triangles " << endl;
-						 MeshError(999,this);
-					 }
-
-					 *tt[1]=   *tt[2]= *t;
-					 // gestion of the link
-					 tt[0]->link=tt[1];
-					 tt[1]->link=tt[2]; 
-
-					 (* tt[0])(OppositeVertex[0])=&s;
-					 (* tt[1])(OppositeVertex[1])=&s;
-					 (* tt[2])(OppositeVertex[2])=&s;
-
-					 tt[0]->det=det3[0];
-					 tt[1]->det=det3[1];
-					 tt[2]->det=det3[2];         
-
-					 //  update adj des triangles externe 
-					 tt[0]->SetAdjAdj(0);
-					 tt[1]->SetAdjAdj(1);
-					 tt[2]->SetAdjAdj(2);
-					 //  update des adj des 3 triangle interne
-					 const int i0 = 0;
-					 const int i1= NextEdge[i0];
-					 const int i2 = PreviousEdge[i0];
-
-					 tt[i0]->SetAdj2(i2,tt[i2],i0);
-					 tt[i1]->SetAdj2(i0,tt[i0],i1);
-					 tt[i2]->SetAdj2(i1,tt[i1],i2);
-
-					 tt[0]->SetTriangleContainingTheVertex();
-					 tt[1]->SetTriangleContainingTheVertex();
-					 tt[2]->SetTriangleContainingTheVertex();
-
-
-					 // swap if the point s is on a edge
-					 if(izerodet>=0) {
-						 //  cout << " the point s is on a edge =>swap " << iedge << " "  << *tt[izerodet] << endl;
-						 int rswap =tt[izerodet]->swap(iedge);
-
-						 if (!rswap) 
-							{
-							 cout << " Pb swap the point s is on a edge =>swap " << iedge << " "  << *tt[izerodet] << endl;
-							}
-						 assert(rswap);
-					 }
+					 printf("bug (%i): Bug double points in\n",nbd0);
+					 throw ErrorException(__FUNCT__,exprintf("See above"));
+				 }
+
+				// remove de MarkUnSwap edge
+				t->SetUnMarkUnSwap(0);     
+				t->SetUnMarkUnSwap(1);     
+				t->SetUnMarkUnSwap(2);
+
+				tt[0]= t;
+				tt[1]= &triangles[nbt++];
+				tt[2]= &triangles[nbt++];
+
+				if (nbt>nbtx) {
+					throw ErrorException(__FUNCT__,exprintf("Not ebough triangles"));
+				}
+
+				*tt[1]=   *tt[2]= *t;
+				// gestion of the link
+				tt[0]->link=tt[1];
+				tt[1]->link=tt[2]; 
+
+				(* tt[0])(OppositeVertex[0])=&s;
+				(* tt[1])(OppositeVertex[1])=&s;
+				(* tt[2])(OppositeVertex[2])=&s;
+
+				tt[0]->det=det3[0];
+				tt[1]->det=det3[1];
+				tt[2]->det=det3[2];         
+
+				//  update adj des triangles externe 
+				tt[0]->SetAdjAdj(0);
+				tt[1]->SetAdjAdj(1);
+				tt[2]->SetAdjAdj(2);
+				//  update des adj des 3 triangle interne
+				const int i0 = 0;
+				const int i1= NextEdge[i0];
+				const int i2 = PreviousEdge[i0];
+
+				tt[i0]->SetAdj2(i2,tt[i2],i0);
+				tt[i1]->SetAdj2(i0,tt[i0],i1);
+				tt[i2]->SetAdj2(i1,tt[i1],i2);
+
+				tt[0]->SetTriangleContainingTheVertex();
+				tt[1]->SetTriangleContainingTheVertex();
+				tt[2]->SetTriangleContainingTheVertex();
+
+
+				// swap if the point s is on a edge
+				if(izerodet>=0) {
+					//  cout << " the point s is on a edge =>swap " << iedge << " "  << *tt[izerodet] << endl;
+					int rswap =tt[izerodet]->swap(iedge);
+
+					if (!rswap) 
+					  {
+						cout << " Pb swap the point s is on a edge =>swap " << iedge << " "  << *tt[izerodet] << endl;
+					  }
+					assert(rswap);
+				}
 	}
 	/*}}}1*/
@@ -2187,6 +2118,5 @@
 		Int4 nbvold=nbv;
 		long int verbosity=2;
-		for (it=0;it<nbt;it++)
-		  {
+		for (it=0;it<nbt;it++){
 			Triangle &t=triangles[it];
 			if (t.link)
@@ -2214,13 +2144,11 @@
 					 }
 			  }
-		  }
+		}
 		ReMakeTriangleContainingTheVertex();    
-		if (nbvold!=nbv) 
-		  {
+		if (nbvold!=nbv){
 			Int4  iv = nbvold;
 			Int4 NbSwap = 0;
 			Icoor2 dete[3];  
-			for (Int4 i=nbvold;i<nbv;i++) 
-			  {// for all the new point
+			for (Int4 i=nbvold;i<nbv;i++) {// for all the new point
 				Vertex & vi = vertices[i];
 				vi.i = toI2(vi.r);
@@ -2234,12 +2162,6 @@
 				Triangle *tcvi = FindTriangleContening(vi.i,dete);
 				if (tcvi && !tcvi->link) {
-					cout << i <<  " PB insert point " << Number(vi) << vi << Number(vi) 
-					  << " tcvi = " << tcvi << " " << tcvi->link << endl;
-					cout << (*tcvi)[1] <<  (*tcvi)[2] << endl;
-					tcvi = FindTriangleContening(vi.i,dete);
-					cout << (*tcvi)[1] <<  (*tcvi)[2] << endl;
-					MeshError(1001,this);
+					printf("problem inserting point\n");
 				}
-
 
 				quadtree->Add(vi);
@@ -2249,17 +2171,14 @@
 				iv++;
 				//      }
-			  }
-			if (verbosity>3) 
-			  {
-				cout << "    Nb Of New Point " << iv ;
-				cout << " Nb swap = " << NbSwap << " to  split internal edges with border vertices" ;}
-
+			}
+			if (verbosity>3) {
+				printf("   number of points: %i\n",iv);
+				printf("   number of swap to  split internal edges with border vertices: %i\n",NbSwap);
 				nbv = iv;
-		  }
-		if (NbSplitEdge >  nbv-nbvold)
-		 cout << " Warning not enough vertices  to split all internal edges "  << endl
-			<< "    we lost " << NbSplitEdge - ( nbv-nbvold) << " Edges Sorry " << endl;
-		if (verbosity>2)
-		 cout << "SplitInternalEdgeWithBorderVertices: Number of splited edge " << NbSplitEdge << endl;
+			}
+		}
+		if (NbSplitEdge>nbv-nbvold) printf("WARNING: not enough vertices  to split all internal edges, we lost %i edges...\n",NbSplitEdge - ( nbv-nbvold));
+		if (verbosity>2) printf("SplitInternalEdgeWithBorderVertices: Number of splited edge %i\n",NbSplitEdge);
+
 		return  NbSplitEdge;
 	}
@@ -2315,15 +2234,7 @@
 					// << "   " <<  Number(vi) << " <--> " << Number(vj) <<endl;
 					Triangle *tcvj = FindTriangleContening(vj.i,dete);
-					if (tcvj && !tcvj->link) 
-					  {
-						cerr << i <<  " PB insert point " << Number(vj) << vj << Number(vi) 
-						  << " tcvj = " << tcvj << " " << tcvj->link << endl;
-						cerr << (*tcvj)[1] <<  (*tcvj)[2] << endl;
-						tcvj = FindTriangleContening(vj.i,dete);
-						cout << (*tcvj)[1] <<  (*tcvj)[2] << endl;
-						MeshError(1001,this);
-					  }
-
-
+					if (tcvj && !tcvj->link){
+						throw ErrorException(__FUNCT__,exprintf("problem inserting point"));
+					}
 					quadtree->Add(vj);
 					assert (tcvj && tcvj->det >= 0) ;// internal 
@@ -2334,19 +2245,16 @@
 			  } 
 			if (verbosity>3) {
-				cout << "    Nb Of New Point " << iv << " Nb Of To close Points " << nbv-iv ;
-				cout << " Nb swap = " << NbSwap << " after " ;}
-
+				printf("   number of new points: %i\n",iv);
+				printf("   number of to close (?) points: %i\n",nbv-iv);
+				printf("   number of swap after: %i\n",NbSwap);
+			}
 				nbv = iv;
 		}
 
-		for (i=nbvold;i<nbv;i++) 
-		 NbSwap += vertices[i].Optim(1);  
-		if (verbosity>3) 
-		 cout << " NbSwap = " <<  NbSwap << endl;
-
+		for (i=nbvold;i<nbv;i++) NbSwap += vertices[i].Optim(1);  
+		if (verbosity>3) printf("   NbSwap=%i\n",NbSwap);
 
 		NbTSwap +=  NbSwap ;
 		return nbv-nbvold;
-
 	}
 	/*}}}1*/
@@ -2864,6 +2772,5 @@
 		SetIntCoor();
 		Int4 i;
-		for (i=0;i<nbv;i++) 
-		 ordre[i]= &vertices[i] ;
+		for (i=0;i<nbv;i++) ordre[i]= &vertices[i] ;
 
 		// construction d'un ordre aleatoire 
@@ -2873,115 +2780,111 @@
 		 ordre[is3]= &vertices[k3 = (k3 + PrimeNumber)% nbv];
 
-
-
-
-		for (i=2 ; det( ordre[0]->i, ordre[1]->i, ordre[i]->i ) == 0;) 
-		 if  ( ++i >= nbv) {
-			 cerr << " All the vertices are aline " << endl;
-			 MeshError(998,this); }
-
-			 // echange i et 2 dans ordre afin 
-			 // que les 3 premiers ne soit pas aligne
-			 Exchange( ordre[2], ordre[i]);
-
-			 // on ajoute un point a l'infini pour construire le maillage
-			 // afin d'avoir une definition simple des aretes frontieres
-			 nbt = 2;
-
-
-			 // on construit un maillage trivale forme
-			 // d'une arete et de 2 triangles
-			 // construit avec le 2 aretes orientes et 
-			 Vertex *  v0=ordre[0], *v1=ordre[1];
-
-			 triangles[0](0) = 0; // sommet pour infini 
-			 triangles[0](1) = v0;
-			 triangles[0](2) = v1;
-
-			 triangles[1](0) = 0;// sommet pour infini 
-			 triangles[1](2) = v0;
-			 triangles[1](1) = v1;
-			 const int e0 = OppositeEdge[0];
-			 const int e1 = NextEdge[e0];
-			 const int e2 = PreviousEdge[e0];
-			 triangles[0].SetAdj2(e0, &triangles[1] ,e0);
-			 triangles[0].SetAdj2(e1, &triangles[1] ,e2);
-			 triangles[0].SetAdj2(e2, &triangles[1] ,e1);
-
-			 triangles[0].det = -1;  // faux triangles
-			 triangles[1].det = -1;  // faux triangles
-
-			 triangles[0].SetTriangleContainingTheVertex();
-			 triangles[1].SetTriangleContainingTheVertex();
-
-			 triangles[0].link=&triangles[1];
-			 triangles[1].link=&triangles[0];
-
-			 //  nbtf = 2;
-			 if (  !quadtree )  quadtree = new QuadTree(this,0);
-			 quadtree->Add(*v0);
-			 quadtree->Add(*v1);
-
-			 // on ajoute les sommets un Ò un 
-			 Int4 NbSwap=0;
-
-			 time1=CPUtime();
-
-			 if (verbosity>3) cout << "  -- Begin of insertion process " << endl;
-
-			 for (Int4 icount=2; icount<nbv; icount++) {
-				 Vertex *vi  = ordre[icount];
-				 //    cout << " Insert " << Number(vi) << endl;
-				 Icoor2 dete[3];
-				 Triangle *tcvi = FindTriangleContening(vi->i,dete);
-				 quadtree->Add(*vi); 
-				 Add(*vi,tcvi,dete);
-				 NbSwap += vi->Optim(1,0);
-
-			 }// fin de boucle en icount
-			 time2=CPUtime();
-			 if (verbosity>3) 
-			  cout << "    NbSwap of insertion " <<    NbSwap 
-				 << " NbSwap/Nbv " <<  (float) NbSwap / (float) nbv 
-				 << " NbUnSwap " << NbUnSwap << " Nb UnSwap/Nbv " 
-				 << (float)NbUnSwap /(float) nbv 
-				 <<endl;
-			 NbUnSwap = 0;
-			 // construction d'un ordre aleatoire 
-			 //  const int PrimeNumber= (nbv % 999983) ? 1000003: 999983 ;
+		for (i=2 ; det( ordre[0]->i, ordre[1]->i, ordre[i]->i ) == 0;){
+			if  ( ++i >= nbv) {
+				throw ErrorException(__FUNCT__,exprintf("all the vertices are aligned"));
+			}
+		}
+		// echange i et 2 dans ordre afin 
+		// que les 3 premiers ne soit pas aligne
+		Exchange( ordre[2], ordre[i]);
+
+		// on ajoute un point a l'infini pour construire le maillage
+		// afin d'avoir une definition simple des aretes frontieres
+		nbt = 2;
+
+		// on construit un maillage trivale forme
+		// d'une arete et de 2 triangles
+		// construit avec le 2 aretes orientes et 
+		Vertex *  v0=ordre[0], *v1=ordre[1];
+
+		triangles[0](0) = 0; // sommet pour infini 
+		triangles[0](1) = v0;
+		triangles[0](2) = v1;
+
+		triangles[1](0) = 0;// sommet pour infini 
+		triangles[1](2) = v0;
+		triangles[1](1) = v1;
+		const int e0 = OppositeEdge[0];
+		const int e1 = NextEdge[e0];
+		const int e2 = PreviousEdge[e0];
+		triangles[0].SetAdj2(e0, &triangles[1] ,e0);
+		triangles[0].SetAdj2(e1, &triangles[1] ,e2);
+		triangles[0].SetAdj2(e2, &triangles[1] ,e1);
+
+		triangles[0].det = -1;  // faux triangles
+		triangles[1].det = -1;  // faux triangles
+
+		triangles[0].SetTriangleContainingTheVertex();
+		triangles[1].SetTriangleContainingTheVertex();
+
+		triangles[0].link=&triangles[1];
+		triangles[1].link=&triangles[0];
+
+		//  nbtf = 2;
+		if (  !quadtree )  quadtree = new QuadTree(this,0);
+		quadtree->Add(*v0);
+		quadtree->Add(*v1);
+
+		// on ajoute les sommets un Ò un 
+		Int4 NbSwap=0;
+
+		time1=CPUtime();
+
+		if (verbosity>3) cout << "  -- Begin of insertion process " << endl;
+
+		for (Int4 icount=2; icount<nbv; icount++) {
+			Vertex *vi  = ordre[icount];
+			//    cout << " Insert " << Number(vi) << endl;
+			Icoor2 dete[3];
+			Triangle *tcvi = FindTriangleContening(vi->i,dete);
+			quadtree->Add(*vi); 
+			Add(*vi,tcvi,dete);
+			NbSwap += vi->Optim(1,0);
+
+		}// fin de boucle en icount
+		time2=CPUtime();
+		if (verbosity>3) 
+		 cout << "    NbSwap of insertion " <<    NbSwap 
+			<< " NbSwap/Nbv " <<  (float) NbSwap / (float) nbv 
+			<< " NbUnSwap " << NbUnSwap << " Nb UnSwap/Nbv " 
+			<< (float)NbUnSwap /(float) nbv 
+			<<endl;
+		NbUnSwap = 0;
+		// construction d'un ordre aleatoire 
+		//  const int PrimeNumber= (nbv % 999983) ? 1000003: 999983 ;
 #ifdef NBLOOPOPTIM
 
-			 k3 = rand()%nbv ; 
-			 for (int is4=0; is4<nbv; is4++) 
-			  ordre[is4]= &vertices[k3 = (k3 + PrimeNumber)% nbv];
-
-			 double timeloop = time2 ;
-			 for(int Nbloop=0;Nbloop<NBLOOPOPTIM;Nbloop++) 
-				{
-				 double time000 = timeloop;
-				 Int4  NbSwap = 0;
-				 for (int is1=0; is1<nbv; is1++) 
-				  NbSwap += ordre[is1]->Optim(0,0);
-				 timeloop = CPUtime();
-				 if (verbosity>3) 
-				  cout << "    Optim Loop "<<Nbloop<<" NbSwap: " <<  NbSwap 
-					 << " NbSwap/Nbv " 	   <<  (float) NbSwap / (float) nbv 
-					 << " CPU=" << timeloop - time000 << "  s, " 
-					 << " NbUnSwap/Nbv " << (float)NbUnSwap /(float) nbv  
-					 <<  endl;
-				 NbUnSwap = 0;
-				 if(!NbSwap) break;
-				}
-			 ReMakeTriangleContainingTheVertex(); 
-			 // because we break the TriangleContainingTheVertex
+		k3 = rand()%nbv ; 
+		for (int is4=0; is4<nbv; is4++) 
+		 ordre[is4]= &vertices[k3 = (k3 + PrimeNumber)% nbv];
+
+		double timeloop = time2 ;
+		for(int Nbloop=0;Nbloop<NBLOOPOPTIM;Nbloop++) 
+		  {
+			double time000 = timeloop;
+			Int4  NbSwap = 0;
+			for (int is1=0; is1<nbv; is1++) 
+			 NbSwap += ordre[is1]->Optim(0,0);
+			timeloop = CPUtime();
+			if (verbosity>3) 
+			 cout << "    Optim Loop "<<Nbloop<<" NbSwap: " <<  NbSwap 
+				<< " NbSwap/Nbv " 	   <<  (float) NbSwap / (float) nbv 
+				<< " CPU=" << timeloop - time000 << "  s, " 
+				<< " NbUnSwap/Nbv " << (float)NbUnSwap /(float) nbv  
+				<<  endl;
+			NbUnSwap = 0;
+			if(!NbSwap) break;
+		  }
+		ReMakeTriangleContainingTheVertex(); 
+		// because we break the TriangleContainingTheVertex
 #endif
-			 time3=CPUtime();
-			 if (verbosity>4) 
-			  cout << "    init " << time1 - time0 << " initialisation,  " 
-				 << time2 - time1 << "s, insert point  " 
-				 << time3 -time2 << "s, optim " << endl
-				 << "     Init Total Cpu Time = " << time3 - time0 << "s " << endl;
-
-			 CurrentTh=OldCurrentTh;
+		time3=CPUtime();
+		if (verbosity>4) 
+		 cout << "    init " << time1 - time0 << " initialisation,  " 
+			<< time2 - time1 << "s, insert point  " 
+			<< time3 -time2 << "s, optim " << endl
+			<< "     Init Total Cpu Time = " << time3 - time0 << "s " << endl;
+
+		CurrentTh=OldCurrentTh;
 	}
 	/*}}}1*/
@@ -2997,6 +2900,6 @@
 		  k++,cerr << " det T" << t << " = " << 0 << endl;
 		if (k!=0) {
-			cerr << " ther is  " << k << "  triangles of mes = 0 " << endl;
-			MeshError(11,this);}
+			throw ErrorException(__FUNCT__,exprintf("there is %i triangles of mes = 0",k));
+		}
 
 			TriangleAdjacent ta(0,0);
@@ -3026,7 +2929,5 @@
 
 			if (k!=0) {
-				cerr << " they is " << k << " lost edges " << endl;
-				cerr << " The boundary is crossing may be!" << endl;
-				MeshError(10,this);
+				throw ErrorException(__FUNCT__,exprintf("There are %i lost edges, the boundary might be crossing",k));
 			}
 			for (Int4 j=0;j<nbv;j++)
@@ -3118,10 +3019,7 @@
 			  }   
 			it++;} // end while (it<nbt)
-			if (nbt == NbOutT ||  !NbSubDomTot) 
-			  {
-				cout << "\n error : " <<  NbOutT << " " << NbSubDomTot <<" " << nbt << endl;
-				cerr << "Error: The boundary is not close => All triangles are outside " << endl;
-				MeshError(888,this);
-			  }
+			if (nbt == NbOutT ||  !NbSubDomTot) {
+				throw ErrorException(__FUNCT__,exprintf("The boundary is not close: all triangles are outside"));
+			}
 
 			delete [] HeapArete;
@@ -3233,6 +3131,5 @@
 				 mark[it]=triangles[it].link ? -1 : -2;
 				Int4 inew =0;
-				for (Int4 i=0;i<NbSubDomains;i++) 
-				  {
+				for (Int4 i=0;i<NbSubDomains;i++) {
 					GeometricalEdge &eg = *Gh.subdomains[i].edge;
 					subdomains[i].ref = Gh.subdomains[i].ref;
@@ -3258,6 +3155,5 @@
 					TriangleAdjacent  ta(t,EdgesVertexTriangle[v0->vint][0]);// previous edges
 
-					while (1) 
-					  {
+					while (1) {
 						assert( v0 == ta.EdgeVertex(1) );
 						//	 cout << " recherche " << Number( ta.EdgeVertex(0)) << endl;
@@ -3268,46 +3164,38 @@
 							 subdomains[i].head=t=ta;
 							//cout << "      triangle  =" << Number(t) << " = " << (*t)[0].r <<  (*t)[1].r <<  (*t)[2].r << endl;
-							if(t<triangles || t >= triangles+nbt || t->det < 0 
-										|| t->link == 0) // Ajoute aout 200 
+							if(t<triangles || t >= triangles+nbt || t->det < 0 || t->link == 0) {
+								throw ErrorException(__FUNCT__,exprintf("bad definition of SubSomain %i",i));
+							}
+							Int4 it = Number(t);
+							if (mark[it] >=0) {
+								if(verbosity>10){
+									cerr << "     Warning: the sub domain " << i << " ref = " << subdomains[i].ref 
+									  << " is previouly defined with "  <<mark[it] << " ref = " << subdomains[mark[it]].ref
+									  << " skip this def " << endl;
+								}
+								break;
+							}
+							if(i != inew) 
+							 Exchange(subdomains[i],subdomains[inew]);
+							inew++;
+							Triangle *tt=t;
+							Int4 kkk=0;
+							do 
 							  {
-								cerr << " Error in the def of sub domain "<<i
-								  << " form border " << NbSubDomains - i  << "/" << NbSubDomains
-								  << ": Bad sens  " << Gh.Number(eg) <<" "<< sens <<  endl;  
-								err++;
-								break;}
-								Int4 it = Number(t);
-								if (mark[it] >=0) {
-									if(verbosity>10)
-									 cerr << "     Warning: the sub domain " << i << " ref = " << subdomains[i].ref 
-										<< " is previouly defined with "  <<mark[it] << " ref = " << subdomains[mark[it]].ref
-										<< " skip this def " << endl;
-									break;}
-									if(i != inew) 
-									 Exchange(subdomains[i],subdomains[inew]);
-									inew++;
-									Triangle *tt=t;
-									Int4 kkk=0;
-									do 
-									  {
-										kkk++;
-										assert(mark[Number(tt)]<0);
-										mark[Number(tt)]=i;
-										tt=tt->link;
-									  } while (tt!=t);
-									if(verbosity>7)
-									 cout << "     Nb de triangles dans le sous domaine " << i << " de ref " << subdomains[i].ref << " = " << kkk << endl;
-									break;}
-									ta = Previous(Adj(ta));         
-									if(t == (Triangle *) ta) {
-										err++;
-										cerr << " Error in the def of sub domain " << i 
-										  << " edge=" << Gh.Number(eg) << " " << sens << endl;
-										break;}
-										//         cout << " NB of remove subdomain " << NbSubDomTot-NbSubDomains<< endl;
-
-					  }
-
-				  }
-				if (err) MeshError(777,this);
+								kkk++;
+								assert(mark[Number(tt)]<0);
+								mark[Number(tt)]=i;
+								tt=tt->link;
+							  } while (tt!=t);
+							if(verbosity>7)
+							 cout << "     Nb de triangles dans le sous domaine " << i << " de ref " << subdomains[i].ref << " = " << kkk << endl;
+							break;
+						}
+						ta = Previous(Adj(ta));         
+						if(t == (Triangle *) ta) {
+							throw ErrorException(__FUNCT__,exprintf("bad definition of SubSomain %i",i));
+						}
+					}
+				}
 
 				if (inew < NbSubDomains) {
@@ -3549,9 +3437,7 @@
 		VertexOnBThVertex = new VertexOnVertex[NbVerticesOnGeomVertex];
 		//
-		if( NbVerticesOnGeomVertex >= nbvx) 
-		  {
-			cerr << " Too much vertices on geometry " << NbVerticesOnGeomVertex << " >= " << nbvx << endl; 
-			MeshError(1,this);
-		  }
+		if( NbVerticesOnGeomVertex >= nbvx) {
+			throw ErrorException(__FUNCT__,exprintf("too many vertices on geometry: %i >= %i",NbVerticesOnGeomVertex,nbvx));
+		}
 		assert(vertices);
 		for (i=0;i<Gh.nbv;i++)
@@ -3764,10 +3650,7 @@
 								 assert (A1-vertices>=0 && A1-vertices <nbv);
 								 break;}
-								 if (!ee.adj[k1])
-									{cerr << "Error adj edge " << BTh.Number(ee) << ", nbe = "  << nbe 
-									 << " Gh.vertices " << Gh.vertices 
-										<< " k1 = " << k1 << " on=" << *ee[k1].on << endl;
-									 cerr << ee[k1].on->gv-Gh.vertices << endl;
-									}
+								 if (!ee.adj[k1]) {
+									 throw ErrorException(__FUNCT__,exprintf(" adj edge %i, nbe=%i, Gh.vertices=%i",BTh.Number(ee),nbe,Gh.vertices));
+								 }
 								 pe = ee.adj[k1]; // next edge
 								 k0 = pe->Intersection(ee); 
@@ -3823,13 +3706,9 @@
 		// end new code	    
 		// do the allocation
-		if(step==0) 
-		  {
+		if(step==0){
 			//if(!NbOfNewPoints) break;// nothing ????? bug 
-			if(nbv+NbOfNewPoints > nbvx) 
-			  {
-				cerr << " Too much vertices on geometry " << nbv+NbOfNewPoints  << " >= " << nbvx << endl;
-				MeshError(3,this);
-			  }
-			//cout << " NbOfNewEdge" << NbOfNewEdge << " NbOfNewPoints " << NbOfNewPoints << endl;
+			if(nbv+NbOfNewPoints > nbvx) {
+				throw ErrorException(__FUNCT__,exprintf("too many vertices on geometry: %i >= %i",nbv+NbOfNewPoints,nbvx));
+			}
 			edges = new Edge[NbOfNewEdge];
 			nbex = NbOfNewEdge;
@@ -3846,5 +3725,4 @@
 
 		delete [] bcurve;
-
 
 		Insert();
@@ -3878,9 +3756,7 @@
 		VerticesOnGeomVertex = new VertexOnGeom[NbVerticesOnGeomVertex];  
 		//
-		if( NbVerticesOnGeomVertex >= nbvx) 
-		  {
-			cerr << " Too much vertices on geometry " << NbVerticesOnGeomVertex << " >= " << nbvx << endl;
-			MeshError(1,this);
-		  }
+		if( NbVerticesOnGeomVertex >= nbvx) {
+			throw ErrorException(__FUNCT__,exprintf("too many vertices on geometry: %i >= %i",NbVerticesOnGeomVertex,nbvx));
+		}
 		for (i=0;i<Gh.nbv;i++)
 		 if (Gh[i].Required()&& Gh[i].IsThe()  ) {//Gh  vertices Required
@@ -4157,5 +4033,5 @@
 			  //	 assert( e[i]);
 		  }
-		if(kk) MeshError(997,this);
+		if(kk) throw ErrorException(__FUNCT__,exprintf("See above"));
 
 		return e;
@@ -4188,6 +4064,5 @@
 		// computation of the det 
 		int Nberr=0;
-		for (i=0;i<nbt;i++)
-		  {
+		for (i=0;i<nbt;i++) {
 			Vertex & v0 = triangles[i][0];
 			Vertex & v1 = triangles[i][1];
@@ -4196,25 +4071,16 @@
 			  {
 				triangles[i].det= det(v0,v1,v2);
-				if (triangles[i].det <=0 && Nberr++ <10)
-				  {
+				if (triangles[i].det <=0 && Nberr++ <10){
 					if(Nberr==1)
-					 if (strfrom)
-					  cerr << "+++ Fatal Error " << strfrom << "(SetInCoor)  Error :  area of Triangle < 0 " << endl; 
-					 else 
-					  cerr << "+++  Fatal Error Triangle (in SetInCoor) area of Triangle < 0" << endl;
-					cerr << " Triangle " << i << "  det  (I2) = " << triangles[i].det ;
-					cerr << " (R2) " << Det(v1.r-v0.r,v2.r-v0.r);
-					cerr << "; The 3  vertices " << endl;
-					cerr << Number(v0) << " "  << Number(v1) << " " 
-					  << Number(v2) << " : " ;
-					cerr << v0.r << v1.r << v2.r << " ; ";
-					cerr << v0.i << v1.i << v2.i << endl;
-				  }
+					 if (strfrom){
+						 throw ErrorException(__FUNCT__,exprintf("Fatal error %s (SetInCoor) : area of Triangle %i < 0",strfrom,i));
+					 }
+					 else{
+						 throw ErrorException(__FUNCT__,exprintf("Fatal error (SetInCoor) : area of Triangle %i < 0",i));
+					 }
+				}
 			  }
-			else
-			 triangles[i].det= -1; // Null triangle; 
-		  }
-		if (Nberr) MeshError(899,this);
-
+			else triangles[i].det= -1; // Null triangle; 
+		}
 	}
 	/*}}}1*/
@@ -4338,7 +4204,7 @@
 			for (i=2 ; det( ordre[0]->i, ordre[1]->i, ordre[i]->i ) == 0;) 
 			 if  ( ++i >= nbvb) {
-				 cerr << "FillHoleInMesh: All the vertices are aline " << nbvb << endl;
-				 MeshError(998,this); }
-				 Exchange( ordre[2], ordre[i]);
+				 throw ErrorException(__FUNCT__,exprintf("FillHoleInMesh: All the vertices are aligned"));
+			 }
+			Exchange( ordre[2], ordre[i]);
 
 				 Vertex *  v0=ordre[0], *v1=ordre[1];
@@ -4402,6 +4268,5 @@
 				 }
 				 if(nbloss) {
-					 cerr << " we loss some  " << nbloss << " "  << " edges other " << knbe << endl;
-					 MeshError(1100,this);
+					 throw ErrorException(__FUNCT__,exprintf("we lost(?) %i edges other %i",nbloss,knbe));
 				 }
 
@@ -4490,6 +4355,5 @@
 						 */
 				 if (k) {
-					 cerr << "Error Nb of triangles edge alone = " << k << endl;
-					 MeshError(9997,this);
+					 throw ErrorException(__FUNCT__,exprintf("number of triangles edges alone = %i",k));
 				 }
 				 FindSubDomain();
@@ -4882,9 +4746,6 @@
 		SetVertexFieldOn();
 
-
-		if (vlast >= vend)
-		  {  
-			cerr << " Not enougth vertices to crack the mesh we need " << nbv << " vertices " << endl;
-			MeshError(555,this);
+		if (vlast >= vend) {  
+			throw ErrorException(__FUNCT__,exprintf("Not enougth vertices: to crack the mesh we need %i vertices",nbv));
 		  }
 		cout << "  NbCrackedVertices " <<  NbCrackedVertices << endl;
@@ -4905,9 +4766,8 @@
 
 			if (! a || !a->t ) {
-				if (a) 
-				  {cerr << " Attention PB TriangleConteningTheVertex  vertex number=" << Number(a) << endl;
-					cerr  << "We forget a call to ReMakeTriangleContainingTheVertex" << endl;}
-					cerr << " Pb with " << B << toR2(B) << endl;
-					MeshError(7777);
+				if (a) {
+					printf("TriangleConteningTheVertex vertex number %i, another call to ReMakeTriangleContainingTheVertex was required\n", Number(a));
+				}
+				throw ErrorException(__FUNCT__,exprintf("problem in Triangles::FindTriangleContening"));
 			}
 			assert(a>= vertices && a < vertices+nbv);
@@ -6025,5 +5885,5 @@
 	} 
 	/*}}}1*/
-
+	/*FUNCTION AGoodNumberPrimeWith{{{1*/
 	Int4 AGoodNumberPrimeWith(Int4 n){
 		const Int4 BigPrimeNumber[] ={ 567890359L,
@@ -6042,10 +5902,12 @@
 			return pi; 
 	}
-
+	/*}}}1*/
+	/*FUNCTION MeshError{{{1*/
 	void MeshError(int Err,Triangles *Th){ 
 		cerr << " Fatal error in the meshgenerator " << Err << endl ;
 		exit(1); 
 	}
-
+	/*}}}1*/
+	/*FUNCTION  ostream& operator{{{1*/
 	ostream& operator <<(ostream& f, const  Triangle & ta) {
 		if(CurrentTh)
@@ -6067,162 +5929,156 @@
 			<< "{" << ta.at[2] << " " << ta.aa[2] << "} " 
 			<< "]" ;
-		return f;}
-
-
-		int SwapForForcingEdge(Vertex   *  & pva ,Vertex  * &   pvb ,
-					TriangleAdjacent & tt1,Icoor2 & dets1, Icoor2 & detsa,Icoor2 & detsb, int & NbSwap)
-		  { // l'arete ta coupe l'arete pva pvb
-			// de cas apres le swap sa coupe toujours
-			// on cherche l'arete suivante 
-			// on suppose que detsa >0 et detsb <0
-			// attention la routine echange pva et pvb 
-
-			if(tt1.Locked()) return 0; // frontiere croise 
-
-			TriangleAdjacent tt2 = Adj(tt1);
-			Triangle *t1=tt1,*t2=tt2;// les 2 triangles adjacent
-			Int1 a1=tt1,a2=tt2;// les 2 numero de l arete dans les 2 triangles
-			assert ( a1 >= 0 && a1 < 3 );
-
-			Vertex & sa= (* t1)[VerticesOfTriangularEdge[a1][0]];
-			Vertex & s1= (*t1)[OppositeVertex[a1]];
-			Vertex & s2= (*t2)[OppositeVertex[a2]];
-
-
-			Icoor2 dets2 = det(*pva,*pvb,s2);
-			Icoor2 det1=t1->det , det2=t2->det ;
-			Icoor2 detT = det1+det2;
-			assert((det1>0 ) && (det2 > 0));
-			assert ( (detsa < 0) && (detsb >0) ); // [a,b] cut infinite line va,bb
-			Icoor2 ndet1 = bamg::det(s1,sa,s2);
-			Icoor2 ndet2 = detT - ndet1;
-
-			int ToSwap =0; //pas de swap
-			if ((ndet1 >0) && (ndet2 >0)) 
-			  { // on peut swaper  
-				if ((dets1 <=0 && dets2 <=0) || (dets2 >=0 && detsb >=0))
-				 ToSwap =1; 
-				else // swap alleatoire 
-				 if (BinaryRand()) 
-				  ToSwap =2; 
-			  }
-			if (ToSwap) NbSwap++,
-			 bamg::swap(t1,a1,t2,a2,&s1,&s2,ndet1,ndet2);
-
-			int ret=1;
-
-			if (dets2 < 0) {// haut
-				dets1 = ToSwap ? dets1 : detsa ;
-				detsa = dets2; 
-				tt1 =  Previous(tt2) ;}
-			else if (dets2 > 0){// bas 
-				dets1 = ToSwap ? dets1 : detsb ;
-				detsb = dets2;
-				//xxxx tt1 = ToSwap ? tt1 : Next(tt2);
-				if(!ToSwap) tt1 =  Next(tt2);
+		return f;
+	}
+	/*}}}1*/
+	/*FUNCTION SwapForForcingEdge{{{1*/
+	int SwapForForcingEdge(Vertex   *  & pva ,Vertex  * &   pvb ,TriangleAdjacent & tt1,Icoor2 & dets1, Icoor2 & detsa,Icoor2 & detsb, int & NbSwap) {
+		// l'arete ta coupe l'arete pva pvb
+		// de cas apres le swap sa coupe toujours
+		// on cherche l'arete suivante 
+		// on suppose que detsa >0 et detsb <0
+		// attention la routine echange pva et pvb 
+
+		if(tt1.Locked()) return 0; // frontiere croise 
+
+		TriangleAdjacent tt2 = Adj(tt1);
+		Triangle *t1=tt1,*t2=tt2;// les 2 triangles adjacent
+		Int1 a1=tt1,a2=tt2;// les 2 numero de l arete dans les 2 triangles
+		assert ( a1 >= 0 && a1 < 3 );
+
+		Vertex & sa= (* t1)[VerticesOfTriangularEdge[a1][0]];
+		Vertex & s1= (*t1)[OppositeVertex[a1]];
+		Vertex & s2= (*t2)[OppositeVertex[a2]];
+
+
+		Icoor2 dets2 = det(*pva,*pvb,s2);
+		Icoor2 det1=t1->det , det2=t2->det ;
+		Icoor2 detT = det1+det2;
+		assert((det1>0 ) && (det2 > 0));
+		assert ( (detsa < 0) && (detsb >0) ); // [a,b] cut infinite line va,bb
+		Icoor2 ndet1 = bamg::det(s1,sa,s2);
+		Icoor2 ndet2 = detT - ndet1;
+
+		int ToSwap =0; //pas de swap
+		if ((ndet1 >0) && (ndet2 >0)) 
+		  { // on peut swaper  
+			if ((dets1 <=0 && dets2 <=0) || (dets2 >=0 && detsb >=0))
+			 ToSwap =1; 
+			else // swap alleatoire 
+			 if (BinaryRand()) 
+			  ToSwap =2; 
+		  }
+		if (ToSwap) NbSwap++,
+		 bamg::swap(t1,a1,t2,a2,&s1,&s2,ndet1,ndet2);
+
+		int ret=1;
+
+		if (dets2 < 0) {// haut
+			dets1 = ToSwap ? dets1 : detsa ;
+			detsa = dets2; 
+			tt1 =  Previous(tt2) ;}
+		else if (dets2 > 0){// bas 
+			dets1 = ToSwap ? dets1 : detsb ;
+			detsb = dets2;
+			//xxxx tt1 = ToSwap ? tt1 : Next(tt2);
+			if(!ToSwap) tt1 =  Next(tt2);
+		}
+		else { // changement de sens 
+			ret = -1;
+			Exchange(pva,pvb);
+			Exchange(detsa,detsb);
+			Exchange(dets1,dets2);
+			Exchange(tt1,tt2);
+			dets1=-dets1;
+			dets2=-dets2;
+			detsa=-detsa;
+			detsb=-detsb;
+
+			if (ToSwap) 
+			 if (dets2 < 0) {// haut
+				 dets1 = (ToSwap ? dets1 : detsa) ;
+				 detsa = dets2; 
+				 tt1 =  Previous(tt2) ;}
+			 else if (dets2 > 0){// bas 
+				 dets1 = (ToSwap ? dets1 : detsb) ;
+				 detsb =  dets2;
+				 if(!ToSwap) tt1 =  Next(tt2);
+			 }
+			 else {// on a fin ???
+				 tt1 = Next(tt2);
+				 ret =0;}
+
+		}
+		return ret;
+	}
+	/*}}}1*/
+	/*FUNCTION ForceEdge{{{1*/
+
+	int ForceEdge(Vertex &a, Vertex & b,TriangleAdjacent & taret)  { 
+		int NbSwap =0;
+		assert(a.t && b.t); // the 2 vertex is in a mesh 
+		int k=0;
+		taret=TriangleAdjacent(0,0); // erreur 
+
+		TriangleAdjacent tta(a.t,EdgesVertexTriangle[a.vint][0]);
+		Vertex   *v1, *v2 = tta.EdgeVertex(0),*vbegin =v2;
+		// we turn around a in the  direct sens  
+
+		Icoor2 det2 = v2 ? det(*v2,a,b): -1 , det1;
+		if(v2) // normal case 
+		 det2 = det(*v2,a,b);
+		else { // no chance infini vertex try the next
+			tta= Previous(Adj(tta));
+			v2 = tta.EdgeVertex(0);
+			vbegin =v2;
+			assert(v2);
+			det2 = det(*v2,a,b);
+			//   cout << " No Change try the next" << endl;
+		}
+
+		while (v2 != &b) {
+			TriangleAdjacent tc = Previous(Adj(tta));    
+			v1 = v2; 
+			v2 = tc.EdgeVertex(0);
+			det1 = det2;
+			det2 =  v2 ? det(*v2,a,b): det2; 
+
+			if((det1 < 0) && (det2 >0)) { 
+				// try to force the edge 
+				Vertex * va = &a, *vb = &b;
+				tc = Previous(tc);
+				assert ( v1 && v2);
+				Icoor2 detss = 0,l=0,ks;
+				// cout << "Real ForcingEdge " << *va << *vb << detss << endl;
+				while ((ks=SwapForForcingEdge(  va,  vb, tc, detss, det1,det2,NbSwap)))
+				 if(l++ > 10000000) {
+					 throw ErrorException(__FUNCT__,exprintf("Loop in forcing Egde, nb de swap=%i, nb of try swap (%i)too big",NbSwap,l));
+				 }
+				Vertex *aa = tc.EdgeVertex(0), *bb = tc.EdgeVertex(1);
+				if (( aa == &a ) && (bb == &b) ||  (bb ==  &a ) && (aa == &b)) {
+					tc.SetLock();
+					a.Optim(1,0);
+					b.Optim(1,0);
+					taret = tc;
+					return NbSwap;
+				}
+				else 
+				  {
+					taret = tc;
+					return -2; // error  boundary is crossing
+				  }
 			}
-			else { // changement de sens 
-				ret = -1;
-				Exchange(pva,pvb);
-				Exchange(detsa,detsb);
-				Exchange(dets1,dets2);
-				Exchange(tt1,tt2);
-				dets1=-dets1;
-				dets2=-dets2;
-				detsa=-detsa;
-				detsb=-detsb;
-
-				if (ToSwap) 
-				 if (dets2 < 0) {// haut
-					 dets1 = (ToSwap ? dets1 : detsa) ;
-					 detsa = dets2; 
-					 tt1 =  Previous(tt2) ;}
-				 else if (dets2 > 0){// bas 
-					 dets1 = (ToSwap ? dets1 : detsb) ;
-					 detsb =  dets2;
-					 if(!ToSwap) tt1 =  Next(tt2);
-				 }
-				 else {// on a fin ???
-					 tt1 = Next(tt2);
-					 ret =0;}
-
-			}
-			return ret;
-		  }
-
-		int ForceEdge(Vertex &a, Vertex & b,TriangleAdjacent & taret)  
-		  { 
-			int NbSwap =0;
-			assert(a.t && b.t); // the 2 vertex is in a mesh 
-			int k=0;
-			taret=TriangleAdjacent(0,0); // erreur 
-
-			TriangleAdjacent tta(a.t,EdgesVertexTriangle[a.vint][0]);
-			Vertex   *v1, *v2 = tta.EdgeVertex(0),*vbegin =v2;
-			// we turn around a in the  direct sens  
-
-			Icoor2 det2 = v2 ? det(*v2,a,b): -1 , det1;
-			if(v2) // normal case 
-			 det2 = det(*v2,a,b);
-			else { // no chance infini vertex try the next
-				tta= Previous(Adj(tta));
-				v2 = tta.EdgeVertex(0);
-				vbegin =v2;
-				assert(v2);
-				det2 = det(*v2,a,b);
-				//   cout << " No Change try the next" << endl;
-			}
-
-			while (v2 != &b) {
-				TriangleAdjacent tc = Previous(Adj(tta));    
-				v1 = v2; 
-				v2 = tc.EdgeVertex(0);
-				det1 = det2;
-				det2 =  v2 ? det(*v2,a,b): det2; 
-
-				if((det1 < 0) && (det2 >0)) { 
-					// try to force the edge 
-					Vertex * va = &a, *vb = &b;
-					tc = Previous(tc);
-					assert ( v1 && v2);
-					Icoor2 detss = 0,l=0,ks;
-					// cout << "Real ForcingEdge " << *va << *vb << detss << endl;
-					while ((ks=SwapForForcingEdge(  va,  vb, tc, detss, det1,det2,NbSwap)))
-					 if(l++ > 10000000) {
-						 cerr << " Loop in forcing Egde AB" 
-							<<"\n vertex A " << a 
-							<<"\n vertex B " <<  b 
-							<<"\n nb de swap " << NbSwap 
-							<<"\n nb of try  swap too big = " <<  l << " gearter than " <<  1000000 << endl;
-
-						 if ( CurrentTh ) 
-						  cerr << " vertex number " << CurrentTh->Number(a) << " " <<  CurrentTh->Number(b) << endl;
-						 MeshError(990);
-					 }
-					Vertex *aa = tc.EdgeVertex(0), *bb = tc.EdgeVertex(1);
-					if (( aa == &a ) && (bb == &b) ||  (bb ==  &a ) && (aa == &b)) {
-						tc.SetLock();
-						a.Optim(1,0);
-						b.Optim(1,0);
-						taret = tc;
-						return NbSwap;
-					}
-					else 
-					  {
-						taret = tc;
-						return -2; // error  boundary is crossing
-					  }
-				}
-				tta = tc;
-				assert(k++<2000);
-				if ( vbegin == v2 ) return -1;// error 
-			}
-
-			tta.SetLock();
-			taret=tta;
-			a.Optim(1,0);
-			b.Optim(1,0);
-			return NbSwap; 
-		  }
+			tta = tc;
+			assert(k++<2000);
+			if ( vbegin == v2 ) return -1;// error 
+		}
+
+		tta.SetLock();
+		taret=tta;
+		a.Optim(1,0);
+		b.Optim(1,0);
+		return NbSwap; 
+	}
+	/*}}}1*/
 
 }
Index: /issm/trunk/src/c/Bamgx/objects/Vertex.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/Vertex.cpp	(revision 2849)
+++ /issm/trunk/src/c/Bamgx/objects/Vertex.cpp	(revision 2850)
@@ -11,4 +11,7 @@
 #include "../../include/macros.h"
 #include "../../toolkits/toolkits.h"
+
+#undef __FUNCT__ 
+#define __FUNCT__ "Vertex"
 
 namespace bamg {