diff --git a/Geo/GEdge.cpp b/Geo/GEdge.cpp
index 940d6e59ca1f1e87d712d37027ccaefb9968c605..76ae91e04c3ff00e5beb9bbcf214c0fd95c435b7 100644
--- a/Geo/GEdge.cpp
+++ b/Geo/GEdge.cpp
@@ -172,7 +172,7 @@ void GEdge::setVisibility(char val, bool recursive)
 std::string GEdge::getAdditionalInfoString()
 {
   std::ostringstream sstream;
-  if(v0 && v1) sstream << "{" << v0->tag() << "," << v1->tag() << "}";
+  if(v0 && v1) sstream << "{" << v0->tag() << " " << v1->tag() << "}";
   return sstream.str();
 }
 
diff --git a/Geo/GFace.cpp b/Geo/GFace.cpp
index 6dcb922de7a53e7c1a14cb37d49277a81e922375..d90e9b09335054a844011be0c6d2c018929566c3 100644
--- a/Geo/GFace.cpp
+++ b/Geo/GFace.cpp
@@ -614,8 +614,6 @@ void GFace::computeMeshSizeFieldAccuracy(double &avg,double &max_e, double &min_
     min_e = std::min(min_e, lone);
   }
 #endif
-  //printf("Emi efficiency tau (%g) =%g nE=%d nT=%d \n", avg, 100*exp(avg/(double)nE), nE, getNumMeshElements());
-
 }
 
 double GFace::curvatureDiv(const SPoint2 &param) const
@@ -902,7 +900,8 @@ class data_wrapper{
 };
 
 // Callback function for BFGS
-void bfgs_callback(const alglib::real_1d_array& x, double& func, alglib::real_1d_array& grad, void* ptr)
+void bfgs_callback(const alglib::real_1d_array& x, double& func,
+                   alglib::real_1d_array& grad, void* ptr)
 {
   data_wrapper* w = static_cast<data_wrapper*>(ptr);
   SPoint3 p = w->get_point();
@@ -950,9 +949,11 @@ GPoint GFace::closestPoint(const SPoint3 &queryPoint, const double initialGuess[
   //  fprintf(F,"View \" \" {\n");
   //  fprintf(F,"SP(%g,%g,%g) {%g};\n",queryPoint.x(),queryPoint.y(),queryPoint.z(),0.0);
   double initial_guesses = 10.0;
-  for(double u = uu.low(); u<=uu.high()+1.e-5; u+=(uu.high()-uu.low())/initial_guesses) {
+  for(double u = uu.low(); u <= uu.high() + 1.e-5;
+      u += (uu.high() - uu.low()) / initial_guesses) {
     //    printf("%f\n", u);
-    for(double v = vv.low(); v<=vv.high()+1.e-5; v+=(vv.high()-vv.low())/initial_guesses) {
+    for(double v = vv.low(); v <= vv.high() + 1.e-5;
+        v += (vv.high() - vv.low()) / initial_guesses) {
       GPoint pnt = point(u, v);
       SPoint3 spnt(pnt.x(), pnt.y(), pnt.z());
       double dist = queryPoint.distance(spnt);
@@ -993,19 +994,19 @@ GPoint GFace::closestPoint(const SPoint3 &queryPoint, const double initialGuess[
   double epsx = 0.0;
   alglib::ae_int_t maxits = 500;
 
-  minlbfgssetcond(state,epsg,epsf,epsx,maxits);
+  minlbfgssetcond(state, epsg, epsf, epsx, maxits);
 
   // Solving the problem
   data_wrapper w;
   w.set_point(queryPoint);
   w.set_face(this);
 
-  minlbfgsoptimize(state,bfgs_callback,NULL,&w);
+  minlbfgsoptimize(state, bfgs_callback, NULL, &w);
 
   // Getting the results
   alglib::minlbfgsreport rep;
 
-  minlbfgsresults(state,x,rep);
+  minlbfgsresults(state, x, rep);
 
   GPoint pntF = point(x[0], x[1]);
   if (rep.terminationtype != 4){
@@ -1286,7 +1287,6 @@ bool GFace::fillPointCloud(double maxDist,
 			   std::vector<SPoint2> *uvpoints,
                            std::vector<SVector3> *normals)
 {
-
   if(!points) return false;
 
   if (buildSTLTriangulation()){
diff --git a/Geo/GPoint.cpp b/Geo/GPoint.cpp
deleted file mode 100644
index 4774096e7f55f5c0366254fcb0983d6bfaca5cea..0000000000000000000000000000000000000000
--- a/Geo/GPoint.cpp
+++ /dev/null
@@ -1,5 +0,0 @@
-// Gmsh - Copyright (C) 1997-2013 C. Geuzaine, J.-F. Remacle
-//
-// See the LICENSE.txt file for license information. Please report all
-// bugs and problems to the public mailing list <gmsh@geuz.org>.
-
diff --git a/Geo/GRegion.cpp b/Geo/GRegion.cpp
index 8331eb8fbdb0687eec2307b30b12c82f187a7616..c77b8edfee009f883407fb95b6d1db5181f52668 100644
--- a/Geo/GRegion.cpp
+++ b/Geo/GRegion.cpp
@@ -209,13 +209,10 @@ void GRegion::setVisibility(char val, bool recursive)
 std::string GRegion::getAdditionalInfoString()
 {
   std::ostringstream sstream;
-  if(l_faces.size() > 20){
-    sstream << "{" << l_faces.front()->tag() << ",...," << l_faces.back()->tag() << "}";
-  }
-  else if(l_faces.size()){
+  if(l_faces.size()){
     sstream << "{";
     for(std::list<GFace*>::iterator it = l_faces.begin(); it != l_faces.end(); ++it){
-      if(it != l_faces.begin()) sstream << ",";
+      if(it != l_faces.begin()) sstream << " ";
       sstream << (*it)->tag();
     }
     sstream << "}";
@@ -288,15 +285,15 @@ bool GRegion::edgeConnected(GRegion *r) const
   return false;
 }
 
-void GRegion::replaceFaces (std::list<GFace*> &new_faces)
+void GRegion::replaceFaces(std::list<GFace*> &new_faces)
 {
-  replaceFacesInternal (new_faces);
+  replaceFacesInternal(new_faces);
   if (l_faces.size() != new_faces.size()){
-    Msg::Fatal("impossible to replace faces in region %d (%d vs %d)",tag(),
-	       l_faces.size(),new_faces.size());
+    Msg::Error("Impossible to replace faces in region %d (%d vs %d)",
+               tag(), l_faces.size(), new_faces.size());
   }
 
-  std::list<GFace*>::iterator it  = l_faces.begin();
+  std::list<GFace*>::iterator it = l_faces.begin();
   std::list<GFace*>::iterator it2 = new_faces.begin();
   std::list<int>::iterator it3 = l_dirs.begin();
   std::list<int> newdirs;
@@ -312,18 +309,17 @@ void GRegion::replaceFaces (std::list<GFace*> &new_faces)
   l_dirs = newdirs;
 }
 
-double GRegion::computeSolidProperties (std::vector<double> cg,
-					std::vector<double> inertia)
+double GRegion::computeSolidProperties(std::vector<double> cg,
+                                       std::vector<double> inertia)
 {
   std::list<GFace*>::iterator it = l_faces.begin();
-  std::list<int>::iterator itdir =  l_dirs.begin();
+  std::list<int>::iterator itdir = l_dirs.begin();
   double volumex = 0;
   double volumey = 0;
   double volumez = 0;
   double surface = 0;
   cg[0] = cg[1] = cg[2] = 0.0;
   for ( ; it != l_faces.end(); ++it,++itdir){
-    //printf("face %d dir %d %d elements\n",(*it)->tag(),*itdir,(int)(*it)->triangles.size());
     for (unsigned int i = 0; i < (*it)->triangles.size(); ++i){
       MTriangle *e = (*it)->triangles[i];
       int npt;
@@ -332,65 +328,65 @@ double GRegion::computeSolidProperties (std::vector<double> cg,
       for (int j=0;j<npt;j++){
 	SPoint3 pt;
 	// compute x,y,z of the integration point
-	e->pnt(pts[j].pt[0],pts[j].pt[1],pts[j].pt[2],pt);
+	e->pnt(pts[j].pt[0], pts[j].pt[1], pts[j].pt[2], pt);
 	double jac[3][3];
 	// compute normal
-	double detJ = e->getJacobian(pts[j].pt[0],pts[j].pt[1],pts[j].pt[2],jac);
-	SVector3 n (jac[2][0],jac[2][1],jac[2][2]);
+	double detJ = e->getJacobian(pts[j].pt[0], pts[j].pt[1], pts[j].pt[2], jac);
+	SVector3 n(jac[2][0], jac[2][1], jac[2][2]);
 	n.normalize();
 	n *= (double)*itdir;
-	surface += detJ*pts[j].weight;
-	volumex += detJ*n.x()*pt.x()*pts[j].weight;
-	volumey += detJ*n.y()*pt.y()*pts[j].weight;
-	volumez += detJ*n.z()*pt.z()*pts[j].weight;
-	cg[0]  += detJ*n.x()*(pt.x()*pt.x())*pts[j].weight*0.5;
-	cg[1]  += detJ*n.y()*(pt.y()*pt.y())*pts[j].weight*0.5;
-	cg[2]  += detJ*n.z()*(pt.z()*pt.z())*pts[j].weight*0.5;
+	surface += detJ* pts[j].weight;
+	volumex += detJ * n.x() * pt.x() * pts[j].weight;
+	volumey += detJ * n.y() * pt.y() * pts[j].weight;
+	volumez += detJ * n.z() * pt.z() * pts[j].weight;
+	cg[0] += detJ * n.x() * (pt.x() * pt.x()) * pts[j].weight * 0.5;
+	cg[1] += detJ * n.y() * (pt.y() * pt.y()) * pts[j].weight * 0.5;
+	cg[2] += detJ * n.z() * (pt.z() * pt.z()) * pts[j].weight * 0.5;
       }
     }
   }
 
-  printf("%g -- %g %g %g\n",surface,volumex,volumey,volumez);
+  printf("%g -- %g %g %g\n", surface, volumex, volumey, volumez);
 
   double volume = volumex;
 
-  cg[0]/=volume;
-  cg[1]/=volume;
-  cg[2]/=volume;
+  cg[0] /= volume;
+  cg[1] /= volume;
+  cg[2] /= volume;
 
   it = l_faces.begin();
-  itdir =  l_dirs.begin();
-  inertia[0] =
-    inertia[1] =
-    inertia[2] =
-    inertia[3] =
-    inertia[4] =
-    inertia[5] = 0.0;
+  itdir = l_dirs.begin();
+  inertia[0] = inertia[1] = inertia[2] = inertia[3] = inertia[4] = inertia[5] = 0.0;
 
   for ( ; it != l_faces.end(); ++it,++itdir){
     for (unsigned int i = 0; i < (*it)->getNumMeshElements(); ++i){
       MElement *e = (*it)->getMeshElement(i);
       int npt;
       IntPt *pts;
-      e->getIntegrationPoints (2*(e->getPolynomialOrder()-1)+3, &npt, &pts);
+      e->getIntegrationPoints(2 * (e->getPolynomialOrder() - 1) + 3, &npt, &pts);
       for (int j = 0; j < npt; j++){
 	SPoint3 pt;
 	// compute x,y,z of the integration point
-	e->pnt(pts[j].pt[0],pts[j].pt[1],pts[j].pt[2],pt);
+	e->pnt(pts[j].pt[0], pts[j].pt[1], pts[j].pt[2], pt);
 	double jac[3][3];
 	// compute normal
-	double detJ = e->getJacobian(pts[j].pt[0],pts[j].pt[1],pts[j].pt[2],jac);
-	SVector3 n (jac[2][0],jac[2][1],jac[2][2]);
+	double detJ = e->getJacobian(pts[j].pt[0], pts[j].pt[1], pts[j].pt[2], jac);
+	SVector3 n(jac[2][0], jac[2][1], jac[2][2]);
 	n *= (double)*itdir;
-	inertia[0]  += pts[j].weight*detJ*n.x()*(pt.x()-cg[0])*(pt.x()-cg[0])*(pt.x()-cg[0])/3.0;
-	inertia[1]  += pts[j].weight*detJ*n.y()*(pt.y()-cg[1])*(pt.y()-cg[1])*(pt.y()-cg[1])/3.0;
-	inertia[2]  += pts[j].weight*detJ*n.z()*(pt.z()-cg[2])*(pt.z()-cg[2])*(pt.z()-cg[2])/3.0;
-	inertia[3]  += pts[j].weight*detJ*n.x()*(pt.y()-cg[1])*(pt.x()-cg[0])*(pt.x()-cg[0])/3.0;
-	inertia[4]  += pts[j].weight*detJ*n.x()*(pt.z()-cg[2])*(pt.x()-cg[0])*(pt.x()-cg[0])/3.0;
-	inertia[5]  += pts[j].weight*detJ*n.y()*(pt.z()-cg[2])*(pt.y()-cg[1])*(pt.y()-cg[1])/3.0;
+	inertia[0] += pts[j].weight * detJ * n.x() *
+          (pt.x() - cg[0]) * (pt.x() - cg[0]) * (pt.x() - cg[0]) / 3.0;
+	inertia[1] += pts[j].weight * detJ * n.y() *
+          (pt.y() - cg[1]) * (pt.y() - cg[1]) * (pt.y() - cg[1]) / 3.0;
+	inertia[2] += pts[j].weight * detJ * n.z() *
+          (pt.z() - cg[2]) * (pt.z() - cg[2]) * (pt.z() - cg[2]) / 3.0;
+	inertia[3] += pts[j].weight * detJ * n.x() *
+          (pt.y() - cg[1]) * (pt.x() - cg[0]) * (pt.x() - cg[0]) / 3.0;
+	inertia[4] += pts[j].weight * detJ * n.x() *
+          (pt.z() - cg[2]) * (pt.x() - cg[0]) * (pt.x() - cg[0]) / 3.0;
+	inertia[5] += pts[j].weight * detJ * n.y() *
+          (pt.z() - cg[2]) * (pt.y() - cg[1]) * (pt.y() - cg[1]) / 3.0;
       }
     }
   }
   return volume;
 }
-
diff --git a/Geo/gmshEdge.cpp b/Geo/gmshEdge.cpp
index 69b29c8e1780bdf202086ad241ae437d10482e38..394d2d419c99b67d543c11d5f9e3c377bd7e297a 100644
--- a/Geo/gmshEdge.cpp
+++ b/Geo/gmshEdge.cpp
@@ -77,7 +77,7 @@ std::string gmshEdge::getAdditionalInfoString()
     std::ostringstream sstream;
     sstream << "{";
     for(int i = 0; i < List_Nbr(c->Control_Points); i++){
-      if(i) sstream << ",";
+      if(i) sstream << " ";
       Vertex *v;
       List_Read(c->Control_Points, i, &v);
       sstream << v->Num;
diff --git a/Mesh/meshGEdge.cpp b/Mesh/meshGEdge.cpp
index 8ef3ca968eb61eb675c4524515a8e1ba5b9b82f1..21e4d1241029e7a964013318ff847932d570cb6a 100644
--- a/Mesh/meshGEdge.cpp
+++ b/Mesh/meshGEdge.cpp
@@ -372,8 +372,8 @@ void meshGEdge::operator() (GEdge *ge)
     N = ge->meshAttributes.nbPointsTransfinite;
   }
   else{
-    if (CTX::instance()->mesh.algo2d == ALGO_2D_BAMG 
-	//	|| CTX::instance()->mesh.algo2d == ALGO_2D_PACK_PRLGRMS 
+    if (CTX::instance()->mesh.algo2d == ALGO_2D_BAMG
+	//	|| CTX::instance()->mesh.algo2d == ALGO_2D_PACK_PRLGRMS
 	|| blf){
       a = Integration(ge, t_begin, t_end, F_Lc_aniso, Points,
                       CTX::instance()->mesh.lcIntegrationPrecision);
@@ -393,7 +393,7 @@ void meshGEdge::operator() (GEdge *ge)
     N = std::max(ge->minimumMeshSegments() + 1, (int)(a + 1.99));
   }
 
-  // force odd number of points for if blossom is used for recombination
+  // force odd number of points if blossom is used for recombination
   if(ge->meshAttributes.Method != MESH_TRANSFINITE &&
      CTX::instance()->mesh.algoRecombine == 1 && N % 2 == 0){
     if(/* 1 ||*/ CTX::instance()->mesh.recombineAll){