diff --git a/Geo/OCCFace.cpp b/Geo/OCCFace.cpp
index afbf418dd0f16d08c484d464bbd45c492e7a1c56..7804b46947fc04ad2ac2ad17a95347078b4ce01f 100644
--- a/Geo/OCCFace.cpp
+++ b/Geo/OCCFace.cpp
@@ -256,7 +256,7 @@ GPoint OCCFace::closestPoint(const SPoint3 &qp,
   GeomAPI_ProjectPointOnSurf proj(pnt, occface, umin, umax, vmin, vmax);
 
   if(!proj.NbPoints()) {
-    Msg::Warning("OCC Project Point on Surface FAIL");
+    Msg::Debug("OCC Project Point on Surface FAIL");
     GPoint gp(0, 0);
     gp.setNoSuccess();
     return gp;
diff --git a/Mesh/meshGFace.cpp b/Mesh/meshGFace.cpp
index cfbed9bb8ab1227ead1612e4201da0c43977dc42..2650d32faf8a177ff6cc51303c3b9d676a4f6c47 100644
--- a/Mesh/meshGFace.cpp
+++ b/Mesh/meshGFace.cpp
@@ -131,8 +131,8 @@ bool pointInsideParametricDomain(std::vector<SPoint2> &bnd, SPoint2 &p,
 
 void trueBoundary(const char *iii, GFace *gf, std::vector<SPoint2> &bnd)
 {
-  //    FILE* view_t = Fopen(iii,"w");
-  //    fprintf(view_t,"View \"True Boundary\"{\n");
+  ///FILE* view_t = Fopen(iii,"w");
+  //      fprintf(view_t,"View \"True Boundary\"{\n");
   std::vector<GEdge *> edg = gf->edges();
 
   std::set<GEdge *> edges(edg.begin(), edg.end());
@@ -149,16 +149,16 @@ void trueBoundary(const char *iii, GFace *gf, std::vector<SPoint2> &bnd)
         double xi = r.low() + (r.high() - r.low()) * t;
         p[k] = ge->reparamOnFace(gf, xi, i);
         if(k > 0) {
-	  //          	  fprintf(view_t,"SL(%g,%g,%g,%g,%g,%g){1,1};\n",p[k-1].x(),p[k-1].y(),0.0,
-	  //          		  p[k].x(),p[k].y(),0.0);
+	  //	            	  fprintf(view_t,"SL(%g,%g,%g,%g,%g,%g){1,1};\n",p[k-1].x(),p[k-1].y(),0.0,
+	  //	            		  p[k].x(),p[k].y(),0.0);
           bnd.push_back(p[k - 1]);
           bnd.push_back(p[k]);
         }
       }
     }
   }
-  //    fprintf(view_t,"};\n");
-  //    fclose(view_t);
+  //      fprintf(view_t,"};\n");
+  //      fclose(view_t);
 }
 
 static void computeElementShapes(GFace *gf, double &worst, double &avg,
@@ -2732,17 +2732,44 @@ void deMeshGFace::operator()(GFace *gf)
   gf->correspondingVertices.clear();
 }
 
-/*
-static double TRIANGLE_VALIDITY(MTriangle *t)
+
+static double TRIANGLE_VALIDITY(GFace *gf, MTriangle *t)
 {
+  SPoint2 p1,p2,p3;
+  reparamMeshVertexOnFace(t->getVertex(0),gf, p1);
+  reparamMeshVertexOnFace(t->getVertex(1),gf, p2);
+  reparamMeshVertexOnFace(t->getVertex(2),gf, p3);
+  SVector3 N1 = gf->normal(p1);
+  SVector3 N2 = gf->normal(p2);
+  SVector3 N3 = gf->normal(p3);
+  SVector3 N = N1 + N2 + N3;
+  N.normalize();
+  SVector3 d1(t->getVertex(1)->x() - t->getVertex(0)->x(),
+	      t->getVertex(1)->y() - t->getVertex(0)->y(),
+	      t->getVertex(1)->z() - t->getVertex(0)->z());
+  SVector3 d2(t->getVertex(2)->x() - t->getVertex(0)->x(),
+	      t->getVertex(2)->y() - t->getVertex(0)->y(),
+	      t->getVertex(2)->z() - t->getVertex(0)->z());
+
+  SVector3 c = crossprod (d1,d2);
+  
+  return dot (N,c);
+   
 }
+
 static bool isMeshValid (GFace *gf){
-  
+  int invalid = 0;
   for (size_t i=0;i<gf->triangles.size(); i++){
+    double v = TRIANGLE_VALIDITY(gf, gf->triangles[i]);
+    if (v < 0)invalid ++;
   }
-
+  if (invalid == 0 || invalid == gf->triangles.size())return true;
+  //  printf("%d %d %d\n",gf->tag(),gf->triangles.size(),invalid);
+  return false;
+  
+  
 }
-*/
+
 
 // for debugging, change value from -1 to -100;
 int debugSurface = -1; //-100;
@@ -2837,6 +2864,15 @@ void meshGFace::operator()(GFace *gf, bool print)
   if(gf->getNumMeshElements() == 0) {
     Msg::Warning("Surface %d consists of no elements", gf->tag());
   }
+
+  if (algoDelaunay2D(gf) && !isMeshValid (gf)){
+    Msg::Warning ("Delaunay based mesher failed on surface %d --> moving to meshadapt",gf->tag());
+    deMeshGFace killer;
+    killer (gf);
+    gf->setMeshingAlgo (1);
+    (*this)(gf,print);
+  }
+
 }
 
 static bool getGFaceNormalFromVert(GFace *gf, MElement *el, SVector3 &nf)