diff --git a/Mesh/meshGFaceRecombine.cpp b/Mesh/meshGFaceRecombine.cpp
index e17b72ff11c39d4dfb96520aa6f167e13cb45f39..341994d4c6e67b1212763c9c056a24d369800642 100644
--- a/Mesh/meshGFaceRecombine.cpp
+++ b/Mesh/meshGFaceRecombine.cpp
@@ -14,7 +14,6 @@
  #define REC2D_DRAW
 
 #include <cmath>
-#include "OpenFile.h"//pour setboundingbox
 #include "meshGFaceRecombine.h"
 #include "MTriangle.h"
 #include "MQuadrangle.h"
@@ -35,7 +34,8 @@ Rec2DData *Rec2DData::_current = NULL;
 double **Rec2DVertex::_qualVSnum = NULL;
 double **Rec2DVertex::_gains = NULL;
 
-int otherParity(int a) {
+int otherParity(int a)
+{
   if (a % 2)
     return a - 1;
   return a + 1;
@@ -43,10 +43,11 @@ int otherParity(int a) {
 
 /**  Recombine2D  **/
 /*******************/
-Recombine2D::Recombine2D(GFace *gf) : _gf(gf)
+Recombine2D::Recombine2D(GFace *gf) : _gf(gf), _strategy(0), _numChange(0)
 {
   if (Recombine2D::_current != NULL) {
     Msg::Warning("[Recombine2D] An instance already in execution");
+    _data = NULL;
     return;
   }
   Recombine2D::_current = this;
@@ -57,7 +58,6 @@ Recombine2D::Recombine2D(GFace *gf) : _gf(gf)
   backgroundMesh::set(_gf);
   _bgm = backgroundMesh::current();
   _data = new Rec2DData();
-  _numChange = 0;
   
   static int po = -1;
   if (++po < 1)
@@ -69,7 +69,7 @@ Recombine2D::Recombine2D(GFace *gf) : _gf(gf)
   
   static int pu = -1;
   if (++pu < 1)
-    Msg::Error("FIXME QualAngle et application action...");
+    Msg::Warning("FIXME Update of Action pointing to edge and vertex (when change)");
   
   // Be able to compute geometrical angle at corners
   std::map<MVertex*, AngleData> mapCornerVert;
@@ -91,16 +91,13 @@ Recombine2D::Recombine2D(GFace *gf) : _gf(gf)
     // triangles
     for (unsigned int i = 0; i < gf->triangles.size(); ++i) {
       MTriangle *t = gf->triangles[i];
-      t->setVolumePositive();
-      static int tia = -1;
-      if (++tia < 1) Msg::Warning("FIXME Sometimes negative volume");
-      Rec2DElement *rel = new Rec2DElement(t);
+      
       Rec2DVertex *rv[3];
       for (int j = 0; j < 3; ++j) {
         MVertex *v = t->getVertex(j);
         if ( (itCorner = mapCornerVert.find(v)) != mapCornerVert.end() ) {
           if (!itCorner->second._rv)
-            itCorner->second._rv = new Rec2DVertex(v, false);
+            itCorner->second._rv = new Rec2DVertex(v);
           rv[j] = itCorner->second._rv;
           itCorner->second._mElements.push_back((MElement*)t);
         }
@@ -111,26 +108,25 @@ Recombine2D::Recombine2D(GFace *gf) : _gf(gf)
           rv[j] = new Rec2DVertex(v);
           mapVert[v] = rv[j];
         }
-        rv[j]->add(rel, false);
       }
+      Rec2DEdge *re[3];
       for (int j = 0; j < 3; ++j) {
-        Rec2DEdge *re;
-        if ( !(re = Rec2DVertex::getCommonEdge(rv[j], rv[(j+1)%3])) )
-          re = new Rec2DEdge(rv[j], rv[(j+1)%3]);
-        rel->add(re);
+        if ( !(re[j] = Rec2DVertex::getCommonEdge(rv[j], rv[(j+1)%3])) )
+          re[j] = new Rec2DEdge(rv[j], rv[(j+1)%3]);
       }
+      
+      new Rec2DElement(t, re, rv);
     }
     // quadrangles
     for (unsigned int i = 0; i < gf->quadrangles.size(); ++i) {
       MQuadrangle *q = gf->quadrangles[i];
-      q->setVolumePositive();
-      Rec2DElement *rel = new Rec2DElement(q);
+      
       Rec2DVertex *rv[4];
       for (int j = 0; j < 4; ++j) {
         MVertex *v = q->getVertex(j);
         if ( (itCorner = mapCornerVert.find(v)) != mapCornerVert.end() ) {
           if (!itCorner->second._rv)
-            itCorner->second._rv = new Rec2DVertex(v, false);
+            itCorner->second._rv = new Rec2DVertex(v);
           rv[j] = itCorner->second._rv;
           itCorner->second._mElements.push_back((MElement*)q);
         }
@@ -141,14 +137,14 @@ Recombine2D::Recombine2D(GFace *gf) : _gf(gf)
           rv[j] = new Rec2DVertex(v);
           mapVert[v] = rv[j];
         }
-        rv[j]->add(rel, false);
       }
+      Rec2DEdge *re[4];
       for (int j = 0; j < 4; ++j) {
-        Rec2DEdge *re;
-        if ( !(re = Rec2DVertex::getCommonEdge(rv[i], rv[(i+1)%3])) )
-          re = new Rec2DEdge(rv[i], rv[(i+1)%3]);
-        rel->add(re);
+        if ( !(re[j] = Rec2DVertex::getCommonEdge(rv[i], rv[(i+1)%3])) )
+          re[j] = new Rec2DEdge(rv[i], rv[(i+1)%3]);
       }
+      
+      new Rec2DElement(q, re, rv);
     }
   }
   // update corner
@@ -162,7 +158,7 @@ Recombine2D::Recombine2D(GFace *gf) : _gf(gf)
     }
   }
   mapCornerVert.clear();
-  // update neighbouring, create the 'Rec2DTwoTri2Quad'
+  // update boundary, create the 'Rec2DTwoTri2Quad'
   {
     Rec2DData::iter_re it = Rec2DData::firstEdge();
     for (; it != Rec2DData::lastEdge(); ++it) {
@@ -190,7 +186,6 @@ Recombine2D::Recombine2D(GFace *gf) : _gf(gf)
     Rec2DData::iter_rv it = Rec2DData::firstVertex();
     for (; it != Rec2DData::lastVertex(); ++it) {
       Rec2DVertex *rv = *it;
-      rv->initQualAngle();
       if (rv->getOnBoundary()) {
         if (!rv->getParity()) {
           int base = Rec2DData::getNewParity();
@@ -226,10 +221,12 @@ bool Recombine2D::recombine()
 #endif
     
     if (!_remainAllQuad(nextAction)) {
-      nextAction->color(190, 0, 0);
       delete nextAction;
+#ifdef REC2D_DRAW
+      nextAction->color(190, 0, 0);
       while (Cpu()-time < REC2D_WAIT_TIME)
         FlGui::instance()->check();
+#endif
       continue;
     }
     ++_numChange;
@@ -259,33 +256,66 @@ bool Recombine2D::recombine()
   return 1;
 }
 
+bool Recombine2D::recombine(int depth)
+{
+  double bestGlobalQuality;
+  Rec2DNode *root = new Rec2DNode(NULL, NULL, bestGlobalQuality, depth);
+  Rec2DNode *currentNode = root->selectBestNode();
+  
+  int i = 0;
+  while (Rec2DData::hasAction()) {
+    Msg::Info("{%d}", ++i);
+    Msg::Info(" ");
+    currentNode->develop(depth, bestGlobalQuality);
+    currentNode = currentNode->selectBestNode();
+  }
+  currentNode->selectBestNode();
+  
+  Msg::Info("==> %g", Rec2DData::getGlobalQuality());
+}
+
 bool Recombine2D::developTree()
 {
-  double bestGlobalValue;
-  Rec2DNode root(NULL, NULL, bestGlobalValue);
+  double bestGlobalQuality;
+  Rec2DNode root(NULL, NULL, bestGlobalQuality);
   _data->printState();
   
-  Msg::Info("best global value : %g", bestGlobalValue);
+  Msg::Info("best global value : %g", bestGlobalQuality);
   Msg::Info("num end node : %d", Rec2DData::getNumEndNode());
   
   Rec2DData::sortEndNode();
   Rec2DData::drawEndNode(100);
-  SetBoundingBox();
   //_gf->triangles.clear();
 }
 
-void Recombine2D::nextTreeActions(Rec2DAction *ra,
-                                  std::vector<Rec2DAction*> &actions)
+void Recombine2D::nextTreeActions(std::vector<Rec2DAction*> &actions,
+                                  std::vector<Rec2DElement*> &neighbourEl)
 {
-  Rec2DAction *next = Rec2DData::getBestNonHiddenAction();
-  if (!next)
-    return;
-  std::vector<Rec2DElement*> elements;
-  next->getElements(elements);
-  for (int i = 0; i < elements[0]->getNumActions(); ++i) {
-    if (!Rec2DData::isOutOfDate(elements[0]->getAction(i)))
-      actions.push_back(elements[0]->getAction(i));
+  actions.clear();
+  Rec2DAction *ra = NULL;
+  Rec2DElement *rel = NULL;
+  switch (_current->_strategy) {
+    case 0 : // random triangle of random action
+      //ra = Rec2DData::getRandomAction();
+      //rel = ra->getRandomTriangle();
+      //break;
+      
+    case 1 : // random triangle of best action
+      //ra = Rec2DData::getBestAction();
+      //rel = ra->getRandomTriangle();
+      //break;
+      
+    case 2 : // random neighbour triangle
+      
+      
+    case 3 : // triangle of best neighbour action
+    default :
+      ra = Rec2DData::getBestAction();
+      if (!ra) return;
+      rel = ra->getRandomElement();
+      break;
   }
+  rel->getActions(actions);
 }
 
 double Recombine2D::_geomAngle(MVertex *v,
@@ -593,11 +623,16 @@ void Rec2DData::remove(Rec2DAction *ra)
 {
   std::list<Rec2DAction*>::iterator it = _current->_actions.begin();
   while (it != _current->_actions.end()) {
-    if (*it == ra)
+    if (*it == ra) {
       it = _current->_actions.erase(it);
+      return;
+    }
     else
       ++it;
   }
+  static int a = -1;
+  if (++a < 1)
+  Msg::Error("[Rec2DData] removing too much action");
 }
 
 void Rec2DData::printState()
@@ -607,7 +642,7 @@ void Rec2DData::printState()
   Msg::Info("numEdge %d (%d), valEdge %g => %g", _numEdge, _edges.size(), (double)_valEdge, (double)_valEdge/_numEdge);
   Msg::Info("numVert %d (%d), valVert %g => %g", _numVert, _vertices.size(), (double)_valVert, (double)_valVert/_numVert);
   Msg::Info("Element (%d)", _elements.size());
-  Msg::Info("global Value %g", Rec2DData::getGlobalValue());
+  Msg::Info("global Value %g", Rec2DData::getGlobalQuality());
   Msg::Info("num action %d", _actions.size());
   std::map<int, std::vector<Rec2DVertex*> >::iterator it = _parities.begin();
   for (; it != _parities.end(); ++it) {
@@ -640,6 +675,21 @@ int Rec2DData::getNewParity()
   return (it->first/2)*2 + 2;
 }
 
+Rec2DDataChange* Rec2DData::getNewDataChange()
+{
+  _current->_changes.push_back(new Rec2DDataChange());
+  return _current->_changes.back();
+}
+
+bool Rec2DData::revertDataChange(Rec2DDataChange *rdc)
+{
+  if (_current->_changes.back() != rdc)
+    return false;
+  _current->_changes.pop_back();
+  rdc->revert();
+  delete rdc;
+}
+
 void Rec2DData::removeParity(Rec2DVertex *rv, int p)
 {
   std::map<int, std::vector<Rec2DVertex*> >::iterator it;
@@ -825,13 +875,13 @@ void Rec2DData::revertAssumedParities()
   }
 }
 
-double Rec2DData::getGlobalValue()
+double Rec2DData::getGlobalQuality()
 {
   double a = (_current->_valVert) / (_current->_numVert);
   return a * (_current->_valEdge) / (_current->_numEdge);
 }
 
-double Rec2DData::getGlobalValue(int numEdge, double valEdge,
+double Rec2DData::getGlobalQuality(int numEdge, double valEdge,
                                    int numVert, double valVert)
 {
   double a = (_current->_valVert + valVert) / (_current->_numVert + numVert);
@@ -845,39 +895,13 @@ Rec2DAction* Rec2DData::getBestAction()
   return *std::max_element(_current->_actions.begin(),
                            _current->_actions.end(), lessRec2DAction());
 }
-
-Rec2DAction* Rec2DData::getBestNonHiddenAction()
-{
-  _current->_actions.sort(gterRec2DAction());
-  std::list<Rec2DAction*>::iterator it = _current->_actions.begin();
-  while (it != _current->_actions.end() && Rec2DData::isOutOfDate(*it)) ++it;
-  if (it == _current->_actions.end())
-    return NULL;
-  return *it;
-}
-
-bool Rec2DData::isOutOfDate(Rec2DAction *ra)
-{
-  std::vector<Rec2DElement*> elements;
-  ra->getElements(elements);
-  bool b = false;
-  for (unsigned int i = 0; i < elements.size(); ++i) {
-    if (_current->_hiddenElements.find(elements[i]) !=
-        _current->_hiddenElements.end()               ) {
-      b = true;
-      break;
-    }
-  }
-  return b;
-}
-
 void Rec2DData::drawEndNode(int num)
 {
   double dx = .0, dy = .0;
   for (unsigned int i = 0; i < num && i < _current->_endNodes.size(); ++i) {
     std::map<int, std::vector<double> > data;
     Rec2DNode *currentNode = _current->_endNodes[i];
-    Msg::Info("%d -> %g", i+1, currentNode->getGlobVal());
+    Msg::Info("%d -> %g", i+1, currentNode->getGlobQual());
     int k = 0;
     if ( !((i+1) % ((int)std::sqrt(num)+1)) ) {
       dx = .0;
@@ -886,8 +910,8 @@ void Rec2DData::drawEndNode(int num)
     else
       dx += 1.2;
     while (currentNode && currentNode->getAction()) {
-      //Msg::Info("%g", currentNode->getGlobVal());
-      //data[currentNode->getNum()].push_back(currentNode->getGlobVal());
+      //Msg::Info("%g", currentNode->getGlobQual());
+      //data[currentNode->getNum()].push_back(currentNode->getGlobQual());
       data[currentNode->getAction()->getNum(dx, dy)].push_back(++k);
       currentNode = currentNode->getFather();
     }
@@ -902,6 +926,55 @@ void Rec2DData::sortEndNode()
             moreRec2DNode()             );
 }
 
+/**  Rec2DDataChange  **/
+/***********************/
+void Rec2DDataChange::hide(Rec2DEdge *re)
+{
+  _hiddenEdge.push_back(re);
+  re->hide();
+}
+
+void Rec2DDataChange::hide(Rec2DElement *rel)
+{
+  _hiddenElement.push_back(rel);
+  rel->hide();
+}
+
+void Rec2DDataChange::hide(std::vector<Rec2DAction*> actions)
+{
+  _hiddenAction.insert(_hiddenAction.end(), actions.begin(), actions.end());
+  for (unsigned int i = 0; i < actions.size(); ++i) {
+    actions[i]->hide();
+  }
+}
+
+void Rec2DDataChange::append(Rec2DElement *rel)
+{
+  _newElement.push_back(rel);
+}
+
+void Rec2DDataChange::revert()
+{
+  for (unsigned int i = 0; i < _newAction.size(); ++i)
+    delete _newAction[i];
+  for (unsigned int i = 0; i < _newElement.size(); ++i)
+    delete _newElement[i];
+  for (unsigned int i = 0; i < _newEdge.size(); ++i)
+    delete _newEdge[i];
+  for (unsigned int i = 0; i < _newVertex.size(); ++i)
+    delete _newVertex[i];
+  
+  for (unsigned int i = 0; i < _hiddenVertex.size(); ++i)
+    _hiddenVertex[i]->reveal();
+  for (unsigned int i = 0; i < _hiddenEdge.size(); ++i)
+    _hiddenEdge[i]->reveal();
+  for (unsigned int i = 0; i < _hiddenElement.size(); ++i)
+    _hiddenElement[i]->reveal();
+  for (unsigned int i = 0; i < _hiddenAction.size(); ++i)
+    _hiddenAction[i]->reveal();
+}
+
+
 /**  Rec2DAction  **/
 /*******************/
 bool lessRec2DAction::operator()(Rec2DAction *ra1, Rec2DAction *ra2) const
@@ -915,7 +988,7 @@ bool gterRec2DAction::operator()(Rec2DAction *ra1, Rec2DAction *ra2) const
 }
 
 Rec2DAction::Rec2DAction()
-: _lastUpdate(Recombine2D::getNumChange()-1), _globValIfExecuted(.0)
+: _lastUpdate(Recombine2D::getNumChange()-1), _globQualIfExecuted(.0)
 {
 }
 
@@ -927,8 +1000,8 @@ bool Rec2DAction::operator<(Rec2DAction &other)
 double Rec2DAction::getReward()
 {
   if (_lastUpdate < Recombine2D::getNumChange())
-    _computeGlobVal();
-  return _globValIfExecuted/* - Rec2DData::getGlobalValue()*/;
+    _computeGlobQual();
+  return _globQualIfExecuted/* - Rec2DData::getGlobalQuality()*/;
 }
 
 
@@ -971,8 +1044,10 @@ Rec2DTwoTri2Quad::Rec2DTwoTri2Quad(Rec2DElement *el0, Rec2DElement *el1)
   Rec2DData::add(this);
 }
 
-Rec2DTwoTri2Quad::~Rec2DTwoTri2Quad()
+void Rec2DTwoTri2Quad::hide()
 {
+  static int a = -1;
+  if (++a < 1) Msg::Warning("FIXME changer le syst�me d'execution action");
   if (_triangles[0])
     _triangles[0]->remove(this);
   if (_triangles[1])
@@ -980,18 +1055,30 @@ Rec2DTwoTri2Quad::~Rec2DTwoTri2Quad()
   Rec2DData::remove(this);
 }
 
-void Rec2DTwoTri2Quad::_computeGlobVal()
+void Rec2DTwoTri2Quad::reveal()
+{
+  static int a = -1;
+  if (++a < 1) Msg::Warning("FIXME regarder si action sans son triangle");
+  
+  if (_triangles[0])
+    _triangles[0]->add(this);
+  if (_triangles[1])
+    _triangles[1]->add(this);
+  Rec2DData::add(this);
+}
+
+void Rec2DTwoTri2Quad::_computeGlobQual()
 {
   double valEdge = -REC2D_EDGE_BASE * _edges[4]->getQual();
   for (int i = 0; i < 4; ++i)
     valEdge += REC2D_EDGE_QUAD * _edges[i]->getQual();
   
-  double valVert = _vertices[0]->getGain(-1) + _vertices[1]->getGain(-1);
+  double valVert;
   valVert += _vertices[0]->getGainMerge(_triangles[0], _triangles[1]);
   valVert += _vertices[1]->getGainMerge(_triangles[0], _triangles[1]);
   
-  _globValIfExecuted =
-    Rec2DData::getGlobalValue(4*REC2D_EDGE_QUAD - REC2D_EDGE_BASE,
+  _globQualIfExecuted =
+    Rec2DData::getGlobalQuality(4*REC2D_EDGE_QUAD - REC2D_EDGE_BASE,
                               valEdge, 0, valVert                 );
   _lastUpdate = Recombine2D::getNumChange();
 }
@@ -1057,25 +1144,23 @@ void Rec2DTwoTri2Quad::apply(std::vector<Rec2DVertex*> &newPar)
   _triangles[1] = NULL;
   delete _edges[4];
   
-  /*new Rec2DCollapse(*/new Rec2DElement(_edges)/*)*/;
-}
-
-void Rec2DTwoTri2Quad::choose(Rec2DElement *&rel)
-{
-  _triangles[0]->hide();
-  _triangles[1]->hide();
-  _edges[4]->hide();
+  new Rec2DElement((MQuadrangle*)NULL, _edges);
   
-  rel = new Rec2DElement(_edges, true);
+  Recombine2D::incNumChange();
 }
 
-void Rec2DTwoTri2Quad::unChoose(Rec2DElement *rel)
-{
-  delete rel;
+void Rec2DTwoTri2Quad::apply(Rec2DDataChange *rdc)
+{  
+  std::vector<Rec2DAction*> actions;
+  _triangles[0]->getUniqueActions(actions);
+  _triangles[1]->getUniqueActions(actions);
+  rdc->hide(actions);
+  rdc->hide(_triangles[0]);
+  rdc->hide(_triangles[1]);
+  rdc->hide(_edges[4]);
+  rdc->append(new Rec2DElement((MQuadrangle*)NULL, _edges));
   
-  _edges[4]->reveal();
-  _triangles[0]->reveal();
-  _triangles[1]->reveal();
+  Recombine2D::incNumChange();
 }
 
 bool Rec2DTwoTri2Quad::isObsolete()
@@ -1162,6 +1247,19 @@ void Rec2DTwoTri2Quad::getElements(std::vector<Rec2DElement*> &elem)
   elem.push_back(_triangles[1]);
 }
 
+void Rec2DTwoTri2Quad::getNeighbourElements(std::vector<Rec2DElement*> &elem)
+{
+  elem.clear();
+  _triangles[0]->getMoreNeighbours(elem);
+  _triangles[1]->getMoreNeighbours(elem);
+  for (unsigned int i = 0; i < elem.size(); ++i) {
+    if (elem[i] == _triangles[0] || elem[i] == _triangles[0]) {
+      elem[i] = elem.back();
+      elem.pop_back();
+    }
+  }
+}
+
 int Rec2DTwoTri2Quad::getNum(double shiftx, double shifty)
 {
   MQuadrangle *quad;
@@ -1189,6 +1287,11 @@ int Rec2DTwoTri2Quad::getNum(double shiftx, double shifty)
   return quad->getNum();
 }
 
+Rec2DElement* Rec2DTwoTri2Quad::getRandomElement()
+{
+  return _triangles[rand() % 2];
+}
+
 
 /**  Rec2DEdge  **/
 /*****************/
@@ -1196,34 +1299,22 @@ Rec2DEdge::Rec2DEdge(Rec2DVertex *rv0, Rec2DVertex *rv1)
 : _rv0(rv0), _rv1(rv1), _boundary(-1), _lastUpdate(-2),
   _weight(REC2D_EDGE_BASE+2*REC2D_EDGE_QUAD)
 {
-  _rv0->add(this);
-  _rv1->add(this);
-  Rec2DData::add(this);
-  _computeQual();
-  Rec2DData::addEdge(_weight, _weight*getQual());
-}
-
-Rec2DEdge::~Rec2DEdge()
-{
-  _rv0->remove(this);
-  _rv1->remove(this);
-  Rec2DData::remove(this);
-  Rec2DData::addEdge(-_weight, -_weight*getQual());
+  reveal();
 }
 
 void Rec2DEdge::hide()
 {
-  //_rv0->remove(this);
-  //_rv1->remove(this);
-  //Rec2DData::remove(this);
+  _rv0->rmv(this);
+  _rv1->rmv(this);
+  Rec2DData::remove(this);
   Rec2DData::addEdge(-_weight, -getVal());
 }
 
 void Rec2DEdge::reveal()
 {
-  //_rv0->add(this);
-  //_rv1->add(this);
-  //Rec2DData::remove(this);
+  _rv0->add(this);
+  _rv1->add(this);
+  Rec2DData::add(this);
   Rec2DData::addEdge(_weight, getVal());
 }
 
@@ -1270,7 +1361,7 @@ void Rec2DEdge::_addWeight(int w)
   Rec2DData::addEdge(w, w*getQual());
 }
 
-Rec2DElement* Rec2DEdge::getSingleElement(Rec2DEdge *re)
+Rec2DElement* Rec2DEdge::getUniqueElement(Rec2DEdge *re)
 {
   std::vector<Rec2DElement*> elem;
   Rec2DVertex::getCommonElements(re->getVertex(0), re->getVertex(1), elem);
@@ -1344,49 +1435,62 @@ Rec2DVertex* Rec2DEdge::getOtherVertex(Rec2DVertex *rv) const
 
 /**  Rec2DVertex  **/
 /*******************/
-Rec2DVertex::Rec2DVertex(MVertex *v, bool toSave)
+Rec2DVertex::Rec2DVertex(MVertex *v)
 : _v(v), _angle(4*M_PI), _onWhat(1), _parity(0), _assumedParity(0),
-  _lastMove(-1), _isMesh(toSave), _sumQualAngle(.0)
+  _lastMove(Recombine2D::getNumChange()), _sumQualAngle(.0)
 {
   reparamMeshVertexOnFace(_v, Recombine2D::getGFace(), _param);
-  if (_isMesh) {
-    Rec2DData::add(this);
-    Rec2DData::addVert(1, getQual());
-  }
+  Rec2DData::add(this);
+#ifdef REC2D_DRAW
   if (_v)
     _v->setIndex(_parity);
+#endif
 }
 
 Rec2DVertex::Rec2DVertex(Rec2DVertex *rv, double ang)
-: _v(rv->_v), _angle(ang), _onWhat(-1), _parity(0), _assumedParity(0),
-  _lastMove(-1), _isMesh(true), _sumQualAngle(.0)
+: _v(rv->_v), _angle(ang), _onWhat(-1), _parity(rv->_parity),
+  _assumedParity(rv->_assumedParity), _lastMove(rv->_lastMove),
+  _sumQualAngle(rv->_sumQualAngle), _edges(rv->_edges),
+  _elements(rv->_elements), _param(_param)
 {
-  _edges = rv->_edges;
-  _elements = rv->_elements;
-  _param = rv->_param;
+  static int a = -1;
+  if (++a < -1) Msg::Warning("FIXME Edges really necessary ?");
+  
   for (int i = 0; i < _edges.size(); ++i) {
     _edges[i]->swap(rv, this);
   }
   Rec2DData::add(this);
-  Rec2DData::addVert(1, getQual());
+  if (!_elements.empty())
+    Rec2DData::addVert(2, getQual());
   delete rv;
+#ifdef REC2D_DRAW
   if (_v)
     _v->setIndex(_parity);
+#endif
 }
 
-Rec2DVertex::~Rec2DVertex()
+void Rec2DVertex::hide()
 {
   if (_parity)
     Rec2DData::removeParity(this, _parity);
   if (_assumedParity)
     Rec2DData::removeAssumedParity(this, _assumedParity);
-  if (_isMesh) {
-    Rec2DData::remove(this);
-    if  (_elements.empty())
-      Rec2DData::addVert(-2, -getQual());
-    else
-      Rec2DData::addVert(-2, -getQual() - _sumQualAngle/_elements.size());
-  }
+  
+  Rec2DData::remove(this);
+  if  (!_elements.empty())
+    Rec2DData::addVert(-2, -getQual());
+}
+
+void Rec2DVertex::reveal()
+{
+  if (_parity)
+    Rec2DData::addParity(this, _parity);
+  if (_assumedParity)
+    Rec2DData::addAssumedParity(this, _assumedParity);
+  
+  Rec2DData::add(this);
+  if  (!_elements.empty())
+    Rec2DData::addVert(2, getQual());
 }
 
 void Rec2DVertex::initStaticTable()
@@ -1416,17 +1520,6 @@ void Rec2DVertex::initStaticTable()
   }
 }
 
-void Rec2DVertex::initQualAngle()
-{
-  if (_sumQualAngle != .0) {
-    Msg::Error("[Rec2DVertex] Mmmh... Sorry, I won't do that !");
-    return;
-  }
-  for (unsigned int i = 0; i < _elements.size(); ++i)
-    _sumQualAngle += _angle2Qual(_elements[i]->getAngle(this));
-  Rec2DData::addVert(1, _sumQualAngle/_elements.size());
-}
-
 Rec2DEdge* Rec2DVertex::getCommonEdge(Rec2DVertex *rv0, Rec2DVertex *rv1)
 {
   for (unsigned int i = 0; i < rv0->_edges.size(); ++i) {
@@ -1498,13 +1591,9 @@ bool Rec2DVertex::_recursiveBoundParity(Rec2DVertex *prevRV, int p0, int p1)
 void Rec2DVertex::setOnBoundary()
 {
   if (_onWhat > 0) {
-    if (_isMesh) {
-      Rec2DData::addValVert(-getQual());
-      _onWhat = 0;
-      Rec2DData::addValVert(getQual());
-    }
-    else
-      _onWhat = 0;
+    Rec2DData::addValVert(-getQual());
+    _onWhat = 0;
+    Rec2DData::addValVert(getQual());
   }
 }
 
@@ -1605,17 +1694,19 @@ void Rec2DVertex::getTriangles(std::set<Rec2DElement*> &tri)
   }
 }
 
-double Rec2DVertex::getQual(int numEl) const
+double Rec2DVertex::getQualDegree(int numEl) const
 {
   int nEl = numEl > -1 ? numEl : _elements.size();
   if (_onWhat > -1)
     return _qualVSnum[_onWhat][nEl];
-  if (nEl == 0)
+  if (nEl == 0) {
+    Msg::Error("[Rec2DVertex] I don't want this anymore !");
     return -10.;
+  }
   return 1. - fabs(2./M_PI * _angle/nEl - 1.);
 }
 
-double Rec2DVertex::getGain(int n) const
+double Rec2DVertex::getGainDegree(int n) const
 {
   if (!n)
     return .0;
@@ -1637,11 +1728,15 @@ double Rec2DVertex::getGain(int n) const
     }
   }
   
-  if (_elements.size() == 0)
+  if (_elements.size() == 0) {
+    Msg::Error("[Rec2DVertex] I don't want this anymore !");
     return 11 - fabs(2./M_PI * _angle/(_elements.size() + n) - 1.);
+  }
   
-  if (_elements.size() + n == 0)
+  if (_elements.size() + n == 0) {
+    Msg::Error("[Rec2DVertex] I don't want this anymore !");
     return fabs(2./M_PI * _angle/_elements.size() - 1.) - 11;
+  }
   
   return fabs(2./M_PI * _angle/_elements.size() - 1.)
          - fabs(2./M_PI * _angle/(_elements.size() + n) - 1.);
@@ -1649,13 +1744,13 @@ double Rec2DVertex::getGain(int n) const
 
 double Rec2DVertex::getGainMerge(Rec2DElement *rel1, Rec2DElement *rel2)
 {
-  double diffQualAngle;
-  diffQualAngle = _angle2Qual(rel1->getAngle(this) + rel2->getAngle(this));
-  diffQualAngle -= _angle2Qual(rel1->getAngle(this));
-  diffQualAngle -= _angle2Qual(rel2->getAngle(this));
+  double qualAngle = _sumQualAngle;
+  qualAngle += _angle2Qual(rel1->getAngle(this) + rel2->getAngle(this));
+  qualAngle -= _angle2Qual(rel1->getAngle(this));
+  qualAngle -= _angle2Qual(rel2->getAngle(this));
   
-  return (_sumQualAngle + diffQualAngle) / (_elements.size() - 1)
-         - _sumQualAngle / _elements.size();
+  return qualAngle / (_elements.size()-1) - getQualAngle();
+         + getGainDegree(-1);
 }
 
 void Rec2DVertex::add(Rec2DEdge *re)
@@ -1678,7 +1773,7 @@ bool Rec2DVertex::has(Rec2DEdge *re) const
   return false;
 }
 
-void Rec2DVertex::remove(Rec2DEdge *re)
+void Rec2DVertex::rmv(Rec2DEdge *re)
 {
   unsigned int i = 0;
   while (i < _edges.size()) {
@@ -1692,7 +1787,7 @@ void Rec2DVertex::remove(Rec2DEdge *re)
   Msg::Error("[Rec2DVertex] Didn't removed edge, didn't have it");
 }
 
-void Rec2DVertex::add(Rec2DElement *rel, bool compQualAngle)
+void Rec2DVertex::add(Rec2DElement *rel)
 {
   for (unsigned int i = 0; i < _elements.size(); ++i) {
     if (_elements[i] == rel) {
@@ -1700,17 +1795,12 @@ void Rec2DVertex::add(Rec2DElement *rel, bool compQualAngle)
       return;
     }
   }
-  if (_isMesh) {
-    if (compQualAngle && !_elements.empty())
-      Rec2DData::addValVert(getGain(1) - _sumQualAngle/_elements.size());
-    else
-      Rec2DData::addValVert(getGain(1));
-    if (compQualAngle)
-      _sumQualAngle += _angle2Qual(rel->getAngle(this));
-  }
+  if (!_elements.empty())
+    Rec2DData::addVert(-2, -getQual());
+  
   _elements.push_back(rel);
-  if (_isMesh && compQualAngle)
-    Rec2DData::addValVert(_sumQualAngle/_elements.size());
+  _sumQualAngle += _angle2Qual(rel->getAngle(this));
+  Rec2DData::addVert(2, getQual());
 }
 
 bool Rec2DVertex::has(Rec2DElement *rel) const
@@ -1722,18 +1812,18 @@ bool Rec2DVertex::has(Rec2DElement *rel) const
   return false;
 }
 
-void Rec2DVertex::remove(Rec2DElement *rel)
+void Rec2DVertex::rmv(Rec2DElement *rel)
 {
   unsigned int i = 0;
   while (i < _elements.size()) {
     if (_elements[i] == rel) {
-      if (_isMesh)
-        Rec2DData::addValVert(getGain(-1) - _sumQualAngle/_elements.size());
+      Rec2DData::addVert(-2, -getQual());
       _sumQualAngle -= _angle2Qual(_elements[i]->getAngle(this));
       _elements[i] = _elements.back();
       _elements.pop_back();
-      if (_isMesh && !_elements.empty())
-        Rec2DData::addValVert(_sumQualAngle/_elements.size());
+      
+      if (!_elements.empty())
+        Rec2DData::addVert(2, getQual());
       return;
     }
     ++i;
@@ -1744,82 +1834,85 @@ void Rec2DVertex::remove(Rec2DElement *rel)
 
 /**  Rec2DElement  **/
 /********************/
-Rec2DElement::Rec2DElement(MTriangle *t) : _mEl((MElement *)t), _numEdge(3)
+Rec2DElement::Rec2DElement(MTriangle *t, Rec2DEdge **re, Rec2DVertex **rv)
+: _mEl((MElement *)t), _numEdge(3)
 {
-  for (int i = 0; i < 4; ++i) {
-    _edges[i] = NULL;
-    _elements[i] = NULL;
+  _edges[3] = NULL;
+  _elements[3] = NULL;
+  
+  for (int i = 0; i < 3; ++i) {
+    _edges[i] = re[i];
+    _edges[i]->addHasTri();
   }
-  Rec2DData::add(this);
-}
-
-Rec2DElement::Rec2DElement(MQuadrangle *q) : _mEl((MElement *)q), _numEdge(4)
-{
-  for (int i = 0; i < 4; ++i) {
-    _edges[i] = NULL;
-    _elements[i] = NULL;
+  
+  for (int i = 0; i < 3; ++i) {
+    _elements[i] = Rec2DEdge::getUniqueElement(_edges[i]);
+    if (_elements[i])
+      _elements[i]->addNeighbour(_edges[i], this);
   }
+  
+  if (rv) {
+    for (int i = 0; i < 3; ++i)
+      rv[i]->add(this);
+  }
+  else {
+    std::vector<Rec2DVertex*> vert;
+    getVertices(vert);
+    for (int i = 0; i < 3; ++i)
+      vert[i]->add(this);
+  }
+  
   Rec2DData::add(this);
 }
 
-Rec2DElement::Rec2DElement(Rec2DEdge **edges, bool tree) : _mEl(NULL), _numEdge(4)
+Rec2DElement::Rec2DElement(MQuadrangle *q, Rec2DEdge **re, Rec2DVertex **rv)
+: _mEl((MElement *)q), _numEdge(4)
 {
   for (int i = 0; i < 4; ++i) {
-    _edges[i] = edges[i];
+    _edges[i] = re[i];
     _edges[i]->addHasQuad();
-    if (tree) {
-      _elements[i] = NULL;
-    }
-    else {
-      _elements[i] = Rec2DEdge::getSingleElement(edges[i]);
-      if (_elements[i])
-        _elements[i]->addNeighbour(_edges[i], this);
-    }
   }
-  std::vector<Rec2DVertex*> vertices(4);
-  getVertices(vertices);
-  for (unsigned int i = 0; i < 4; ++i) {
-    vertices[i]->add(this);
+  
+  for (int i = 0; i < 4; ++i) {
+    _elements[i] = Rec2DEdge::getUniqueElement(_edges[i]);
+    if (_elements[i])
+      _elements[i]->addNeighbour(_edges[i], this);
+  }
+  
+  if (rv) {
+    for (int i = 0; i < 4; ++i)
+      rv[i]->add(this);
   }
+  else {
+    std::vector<Rec2DVertex*> vert;
+    getVertices(vert);
+    for (int i = 0; i < 4; ++i)
+      vert[i]->add(this);
+  }
+  
   Rec2DData::add(this);
 }
 
-Rec2DElement::~Rec2DElement()
+void Rec2DElement::hide()
 {
   if (_actions.size())
-    Msg::Error("[Rec2DElement] I didn't want to be deleted :'(");
+    Msg::Error("[Rec2DElement] I didn't want to be hidden :'(");
   for (int i = 0; i < _numEdge; ++i) {
     if (_numEdge == 3)
       _edges[i]->remHasTri();
     else
       _edges[i]->remHasQuad();
     if (_elements[i])
-      _elements[i]->removeNeighbour(_edges[i], this);
+      _elements[i]->rmvNeighbour(_edges[i], this);
   }
   std::vector<Rec2DVertex*> vertices(_numEdge);
   getVertices(vertices);
   for (unsigned int i = 0; i < _numEdge; ++i) {
-    vertices[i]->remove(this);
+    vertices[i]->rmv(this);
   }
   Rec2DData::remove(this);
 }
 
-void Rec2DElement::hide()
-{
-  for (int i = 0; i < _numEdge; ++i) {
-    if (_numEdge == 3)
-      _edges[i]->remHasTri();
-    else
-      _edges[i]->remHasQuad();
-  }
-  std::vector<Rec2DVertex*> vertices(_numEdge);
-  getVertices(vertices);
-  for (unsigned int i = 0; i < _numEdge; ++i) {
-    vertices[i]->remove(this);
-  }
-  Rec2DData::addHidden(this);
-}
-
 void Rec2DElement::reveal()
 {
   for (int i = 0; i < _numEdge; ++i) {
@@ -1827,13 +1920,15 @@ void Rec2DElement::reveal()
       _edges[i]->addHasTri();
     else
       _edges[i]->addHasQuad();
+    if (_elements[i])
+      _elements[i]->addNeighbour(_edges[i], this);
   }
   std::vector<Rec2DVertex*> vertices(_numEdge);
   getVertices(vertices);
   for (unsigned int i = 0; i < _numEdge; ++i) {
     vertices[i]->add(this);
   }
-  Rec2DData::remHidden(this);
+  Rec2DData::add(this);
 }
 
 void Rec2DElement::add(Rec2DEdge *re)
@@ -1906,7 +2001,7 @@ void Rec2DElement::addNeighbour(Rec2DEdge *re, Rec2DElement *rel)
 }
 
 
-void Rec2DElement::removeNeighbour(Rec2DEdge *re, Rec2DElement *rel)
+void Rec2DElement::rmvNeighbour(Rec2DEdge *re, Rec2DElement *rel)
 {
   for (int i = 0; i < _numEdge; ++i) {
     if (_edges[i] == re) {
@@ -1999,6 +2094,14 @@ void Rec2DElement::getVertices(std::vector<Rec2DVertex*> &verts) const
   }
 }
 
+void Rec2DElement::getMoreNeighbours(std::vector<Rec2DElement*> &elem) const
+{
+  for (int i = 0; i < _numEdge; ++i) {
+    if (_elements[i])
+      elem.push_back(_elements[i]);
+  }
+}
+
 void Rec2DElement::getUniqueActions(std::vector<Rec2DAction*> &vectRA) const
 {
   for (unsigned int i = 0; i < _actions.size(); ++i) {
@@ -2093,45 +2196,152 @@ bool moreRec2DNode::operator()(Rec2DNode *rn1, Rec2DNode *rn2) const
   return rn1->getNumTri() < rn2->getNumTri();
 }
 
-Rec2DNode::Rec2DNode(Rec2DNode *father, Rec2DAction *ra, double &bestEndGlobVal)
-: _father(father), _ra(ra), _globalValue(.0), _bestEndGlobVal(.0)
+Rec2DNode::Rec2DNode(Rec2DNode *father, Rec2DAction *ra,
+                     double &bestEndGlobQual, int depth )
+: _father(father), _ra(ra), _globalQuality(.0), _bestEndGlobQual(.0),
+  _dataChange(NULL)
 {
-  _son[0] = NULL;
-  _son[1] = NULL;
-  _son[2] = NULL;
-  Rec2DElement *newElement;
+  for (int i = 0; i < REC2D_NUM_SON; ++i)
+    _son[i] = NULL;
+  
+  std::vector<Rec2DElement*> neighbours;
   if (_ra) {
-    _ra->choose(newElement);
-    _remainingTri = father->getNumTri() - _ra->getNumElement();
+    _ra->getNeighbourElements(neighbours);
+    _dataChange = Rec2DData::getNewDataChange();
+    _ra->apply(_dataChange);
   }
+  
+  if (father)
+    _remainingTri = father->getNumTri() - _ra->getNumElement();
   else
     _remainingTri = Rec2DData::getNumElement();
+  _globalQuality = Rec2DData::getGlobalQuality();
   
+  if (depth != 0) {
+    std::vector<Rec2DAction*> actions;
+    Recombine2D::nextTreeActions(actions, neighbours);
+    
+    if (actions.empty()) {
+      _bestEndGlobQual = _globalQuality;
+      if (depth < 0)
+        Rec2DData::addEndNode(this);
+    }
+    else {
+      for (int i = 0; i < actions.size(); ++i) {
+        double bestSonEGQ;
+        _son[i] = new Rec2DNode(this, actions[i], bestSonEGQ, depth-1);
+        if (bestSonEGQ > _bestEndGlobQual) {
+          _bestEndGlobQual = bestSonEGQ;
+          Rec2DNode *tmp = _son[0];
+          _son[0] = _son[i];
+          _son[i] = tmp;
+        }
+      }
+    }
+  }
+  else
+    _bestEndGlobQual = _globalQuality;
   
-  _globalValue = Rec2DData::getGlobalValue();
+  bestEndGlobQual = _bestEndGlobQual;
   
-  std::vector<Rec2DAction*> actions;
-  Recombine2D::nextTreeActions(_ra, actions);
+  if (_dataChange) {
+    Rec2DData::revertDataChange(_dataChange);
+    _dataChange = NULL;
+  }
+}
+
+Rec2DNode::~Rec2DNode()
+{
+  for (int i = 0; i < REC2D_NUM_SON; ++i)
+    delete _son[i];
+}
+
+Rec2DNode* Rec2DNode::selectBestNode()
+{
+  for (int i = 1; i < REC2D_NUM_SON; ++i) {
+    delete _son[i];
+    _son[i] = NULL;
+  }
   
-  if (actions.empty()) {
-    Rec2DData::addEndNode(this);
-    _bestEndGlobVal = _globalValue;
+  if (_ra) {
+    _dataChange = Rec2DData::getNewDataChange();
+    _ra->apply(_dataChange);
   }
-  else
-    for (int i = 0; i < actions.size(); ++i) {
-      double bestSonEndGlobVal;
-      _son[i] = new Rec2DNode(this, actions[i], bestSonEndGlobVal);
-      _bestEndGlobVal = std::max(_bestEndGlobVal, bestSonEndGlobVal);
+  
+  return _son[0];
+}
+
+void Rec2DNode::develop(int depth, double &bestEndGlobQual)
+{
+  _bestEndGlobQual = .0;
+  if (_son[0] && depth > 0) {
+    int i = 0;
+    double bestSonEGQ;
+    while (_son[i]) {
+      _son[i]->develop(depth-1, bestSonEGQ);
+      if (bestSonEGQ > _bestEndGlobQual) {
+        _bestEndGlobQual = bestSonEGQ;
+        Rec2DNode *tmp = _son[0];
+        _son[0] = _son[i];
+        _son[i] = tmp;
+      }
+      ++i;
+    }
+  }
+  else if (depth > 0) {
+    std::vector<Rec2DElement*> neighbours;
+    if (_ra) {
+      _ra->getNeighbourElements(neighbours);
+      _dataChange = Rec2DData::getNewDataChange();
+      _ra->apply(_dataChange);
     }
     
-  bestEndGlobVal = _bestEndGlobVal;
+    std::vector<Rec2DAction*> actions;
+    Recombine2D::nextTreeActions(actions, neighbours);
+    
+    if (!actions.empty()) {
+      for (int i = 0; i < actions.size(); ++i) {
+        double bestSonEGQ;
+        _son[i] = new Rec2DNode(this, actions[i], bestSonEGQ, depth-1);
+        if (bestSonEGQ > _bestEndGlobQual) {
+          _bestEndGlobQual = bestSonEGQ;
+          if (i > 0) {
+            Rec2DNode *tmp = _son[0];
+            _son[0] = _son[i];
+            _son[i] = tmp;
+          }
+        }
+      }
+    }
+    else
+      _bestEndGlobQual = _globalQuality;
+  }
+  else
+    Msg::Error("[Rec2DNode] That should not happen :'(");
+  
+  bestEndGlobQual = _bestEndGlobQual;
+}
+
+bool Rec2DNode::hasAction()
+{
+  Msg::Warning("FIXME Need definition");
+  /*Rec2DElement *newElement;
+  if (_ra)
+    _ra->choose(newElement);
+    
+  std::vector<Rec2DAction*> actions;
+  Recombine2D::nextTreeActions(_ra, actions);
+  
+  bool b = !actions.empty();
   
   if (_ra)
     _ra->unChoose(newElement);
+  
+  return b;*/
 }
 
 bool Rec2DNode::operator<(Rec2DNode &other)
 {
-  return _globalValue < other._globalValue;
+  return _globalQuality < other._globalQuality;
 }
 
diff --git a/Mesh/meshGFaceRecombine.h b/Mesh/meshGFaceRecombine.h
index eac015caa806c0792ecf8f7b6c9c4debbea2ab53..5905d2c2479dead9df8ebc4c54712a98e269bcee 100644
--- a/Mesh/meshGFaceRecombine.h
+++ b/Mesh/meshGFaceRecombine.h
@@ -13,12 +13,13 @@
 #define REC2D_EDGE_BASE 2
 #define REC2D_EDGE_QUAD 1
 #define REC2D_ALIGNMENT .5
+#define REC2D_NUM_SON 3
 
 #include "GFace.h"
 #include "BackgroundMesh.h"
 //#include "GModel.h"
 //#include "MEdge.h"
-#include "MQuadrangle.h"
+//#include "MQuadrangle.h"
 
 class Rec2DNode;
 class Rec2DVertex;
@@ -26,6 +27,7 @@ class Rec2DEdge;
 class Rec2DElement;
 class Rec2DAction;
 class Rec2DData;
+class Rec2DDataChange;
 struct lessRec2DAction {
   bool operator()(Rec2DAction*, Rec2DAction*) const;
 };
@@ -43,6 +45,7 @@ struct moreRec2DNode {
   bool operator()(Rec2DNode*, Rec2DNode*) const;
 };
 
+//
 class Recombine2D {
   private :
     GFace *_gf;
@@ -51,17 +54,23 @@ class Recombine2D {
     static Recombine2D *_current;
     
     int _numChange;
+    int _strategy;
     
   public :
     Recombine2D(GFace*);
     ~Recombine2D();
     
     bool recombine();
+    bool recombine(int depth);
     bool developTree();
-    static void nextTreeActions(Rec2DAction*, std::vector<Rec2DAction*>&);
+    static void nextTreeActions(std::vector<Rec2DAction*>&,
+                                std::vector<Rec2DElement*> &neighbours);
+    
+    inline void setStrategy(int s) {_strategy = s;}
     
     static inline GFace* getGFace() {return _current->_gf;}
     static inline int getNumChange() {return _current->_numChange;}
+    static inline void incNumChange() {++_current->_numChange;}
     static inline backgroundMesh* bgm() {return _current->_bgm;}
     static void add(MQuadrangle *q) {_current->_gf->quadrangles.push_back(q);}
     
@@ -81,10 +90,11 @@ class Rec2DData {
     
     std::set<Rec2DEdge*> _edges;
     std::set<Rec2DVertex*> _vertices;
-    std::set<Rec2DElement*> _elements, _hiddenElements;
+    std::set<Rec2DElement*> _elements;
     
     std::list<Rec2DAction*> _actions;
     std::vector<Rec2DNode*> _endNodes;
+    std::vector<Rec2DDataChange*> _changes;
     
     std::map<int, std::vector<Rec2DVertex*> > _parities;
     std::map<int, std::vector<Rec2DVertex*> > _assumedParities;
@@ -100,11 +110,13 @@ class Rec2DData {
     std::vector<MQuadrangle*> _quad;
 #endif
     
-    static int getNumEndNode() {return _current->_endNodes.size();}
-    static int getNumElement() {return _current->_elements.size();}
+    static inline int getNumEndNode() {return _current->_endNodes.size();}
+    static inline int getNumElement() {return _current->_elements.size();}
+    static Rec2DDataChange* getNewDataChange();
+    static bool revertDataChange(Rec2DDataChange*);
     
-    static double getGlobalValue();
-    static double getGlobalValue(int numEdge, double valEdge,
+    static double getGlobalQuality();
+    static double getGlobalQuality(int numEdge, double valEdge,
                                  int numVert, double valVert );
     static inline void addVert(int num, double val) {
       _current->_numVert += num;
@@ -123,7 +135,7 @@ class Rec2DData {
     static inline int getNumVert() {return _current->_numVert;}
     static inline double getValVert() {return (double)_current->_valVert;}
     static Rec2DAction* getBestAction();
-    static Rec2DAction* getBestNonHiddenAction();
+    static inline bool hasAction() {return !_current->_actions.empty();}
     
     typedef std::set<Rec2DEdge*>::iterator iter_re;
     typedef std::set<Rec2DVertex*>::iterator iter_rv;
@@ -153,13 +165,6 @@ class Rec2DData {
     }
     static void sortEndNode();
     static inline void drawEndNode(int num);
-    static inline void addHidden(Rec2DElement *rel) {
-      _current->_hiddenElements.insert(rel);
-    }
-    static inline void remHidden(Rec2DElement *rel) {
-      _current->_hiddenElements.erase(rel);
-    }
-    static bool isOutOfDate(Rec2DAction*);
     
     static int getNewParity();
     static void removeParity(Rec2DVertex*, int);
@@ -178,32 +183,53 @@ class Rec2DData {
     static void revertAssumedParities();
 };
 
+class Rec2DDataChange {
+  private :
+    std::vector<Rec2DEdge*> _hiddenEdge, _newEdge;
+    std::vector<Rec2DVertex*> _hiddenVertex, _newVertex;
+    std::vector<Rec2DElement*> _hiddenElement, _newElement;
+    std::vector<Rec2DAction*> _hiddenAction, _newAction;
+    std::vector<std::pair<Rec2DVertex*, SPoint2> > _oldCoordinate;
+    
+  public :
+    void hide(Rec2DEdge*);
+    void hide(Rec2DElement*);
+    void hide(std::vector<Rec2DAction*>);
+    
+    void append(Rec2DElement*);
+    
+    void revert();
+};
+
 class Rec2DAction {
   protected :
-    double _globValIfExecuted;
+    double _globQualIfExecuted;
     int _lastUpdate;
     
   public :
     Rec2DAction();
-    virtual inline ~Rec2DAction() {Rec2DData::remove(this);}
+    virtual ~Rec2DAction() {}
+    virtual void hide() = 0;
+    virtual void reveal() = 0;
     
     bool operator<(Rec2DAction&);
     virtual double getReward();
     virtual void color(int, int, int) = 0;
     virtual void apply(std::vector<Rec2DVertex*> &newPar) = 0;
+    virtual void apply(Rec2DDataChange*) = 0;
     virtual bool isObsolete() = 0;
     virtual bool isAssumedObsolete() = 0;
     virtual void getAssumedParities(int*) = 0;
     virtual bool whatWouldYouDo(std::map<Rec2DVertex*, std::vector<int> >&) = 0;
     virtual Rec2DVertex* getVertex(int) = 0;
-    virtual void choose(Rec2DElement*&) = 0;
-    virtual void unChoose(Rec2DElement*) = 0;
     virtual int getNumElement() = 0;
     virtual void getElements(std::vector<Rec2DElement*>&) = 0;
+    virtual void getNeighbourElements(std::vector<Rec2DElement*>&) = 0;
     virtual int getNum(double shiftx, double shifty) = 0;
+    virtual Rec2DElement* getRandomElement() = 0;
     
   private :
-    virtual void _computeGlobVal() = 0;
+    virtual void _computeGlobQual() = 0;
 };
 
 class Rec2DTwoTri2Quad : public Rec2DAction {
@@ -214,24 +240,27 @@ class Rec2DTwoTri2Quad : public Rec2DAction {
     
   public :
     Rec2DTwoTri2Quad(Rec2DElement*, Rec2DElement*);
-    ~Rec2DTwoTri2Quad();
+    ~Rec2DTwoTri2Quad() {hide();}
+    virtual void hide();
+    virtual void reveal();
     
     virtual void color(int, int, int);
     virtual void apply(std::vector<Rec2DVertex*> &newPar);
+    virtual void apply(Rec2DDataChange*);
     virtual bool isObsolete();
     virtual bool isAssumedObsolete();
     static bool isObsolete(int*);
     virtual void getAssumedParities(int*);
     virtual bool whatWouldYouDo(std::map<Rec2DVertex*, std::vector<int> >&);
     virtual inline Rec2DVertex* getVertex(int i) {return _vertices[i];} //-
-    virtual void choose(Rec2DElement*&);
-    virtual void unChoose(Rec2DElement*);
     virtual inline int getNumElement() {return 2;}
     virtual void getElements(std::vector<Rec2DElement*>&);
+    virtual void getNeighbourElements(std::vector<Rec2DElement*>&);
     virtual int getNum(double shiftx, double shifty);
+    virtual Rec2DElement* getRandomElement();
     
   private :
-    virtual void _computeGlobVal();
+    virtual void _computeGlobQual();
 };
 
 class Rec2DEdge {
@@ -243,7 +272,7 @@ class Rec2DEdge {
     
   public :
     Rec2DEdge(Rec2DVertex*, Rec2DVertex*);
-    ~Rec2DEdge();
+    ~Rec2DEdge() {hide();}
     void hide();
     void reveal();
     
@@ -258,7 +287,7 @@ class Rec2DEdge {
     
     inline Rec2DVertex* getVertex(int i) const {if (i) return _rv1; return _rv0;}
     Rec2DVertex* getOtherVertex(Rec2DVertex*) const;
-    static Rec2DElement* getSingleElement(Rec2DEdge*);
+    static Rec2DElement* getUniqueElement(Rec2DEdge*);
     
     void swap(Rec2DVertex *oldRV, Rec2DVertex *newRV);
     
@@ -287,20 +316,21 @@ class Rec2DVertex {
     int _lastMove, _onWhat; // _onWhat={-1:corner,0:edge,1:face}
     int _parity, _assumedParity;
     SPoint2 _param;
-    bool _isMesh;
     
     static double **_qualVSnum;
     static double **_gains;
     
   public :
-    Rec2DVertex(MVertex*, bool toSave = true);
+    Rec2DVertex(MVertex*);
     Rec2DVertex(Rec2DVertex*, double angle);
-    ~Rec2DVertex();
+    ~Rec2DVertex() {hide();}
+    void hide();
+    void reveal();
     
-    double getQual(int numEl = -1) const;
-    double getGain(int) const;
-    void initQualAngle();
-    inline double getQualAngle() {return _sumQualAngle/_elements.size();}
+    inline double getQual() const {return getQualDegree() + getQualAngle();}
+    inline double getQualAngle() const {return _sumQualAngle/_elements.size();}
+    double getQualDegree(int numEl = -1) const;
+    double getGainDegree(int) const;
     double getGainMerge(Rec2DElement*, Rec2DElement*);
     
     inline void setOnBoundary();
@@ -330,11 +360,11 @@ class Rec2DVertex {
     
     void add(Rec2DEdge*);
     bool has(Rec2DEdge*) const;
-    void remove(Rec2DEdge*);
+    void rmv(Rec2DEdge*);
     
-    void add(Rec2DElement*, bool b = true);
+    void add(Rec2DElement*);
     bool has(Rec2DElement*) const;
-    void remove(Rec2DElement*);
+    void rmv(Rec2DElement*);
     
     static void initStaticTable();
     static Rec2DEdge* getCommonEdge(Rec2DVertex*, Rec2DVertex*);
@@ -355,10 +385,9 @@ class Rec2DElement {
     std::vector<Rec2DAction*> _actions;
     
   public :
-    Rec2DElement(MTriangle*);
-    Rec2DElement(MQuadrangle*);
-    Rec2DElement(Rec2DEdge**, bool tree = false);
-    ~Rec2DElement();
+    Rec2DElement(MTriangle*, Rec2DEdge**, Rec2DVertex **rv = NULL);
+    Rec2DElement(MQuadrangle*, Rec2DEdge**, Rec2DVertex **rv = NULL);
+    ~Rec2DElement() {hide();}
     void hide();
     void reveal();
     
@@ -370,7 +399,7 @@ class Rec2DElement {
     void add(Rec2DAction*);
     void remove(Rec2DAction*);
     void addNeighbour(Rec2DEdge*, Rec2DElement*);
-    void removeNeighbour(Rec2DEdge*, Rec2DElement*);
+    void rmvNeighbour(Rec2DEdge*, Rec2DElement*);
     
     inline MElement* getMElement() const {return _mEl;}
 #ifdef REC2D_DRAW
@@ -397,10 +426,12 @@ class Rec2DElement {
     
     inline int getNumActions() const {return _actions.size();}
     inline Rec2DAction* getAction(int i) const {return _actions[i];}
+    inline void getActions(std::vector<Rec2DAction*> &v) const {v = _actions;};
     void getUniqueActions(std::vector<Rec2DAction*>&) const;
     void getAssumedParities(int*) const;
     void getMoreEdges(std::vector<Rec2DEdge*>&) const;
     void getVertices(std::vector<Rec2DVertex*>&) const;
+    void getMoreNeighbours(std::vector<Rec2DElement*>&) const;
     Rec2DVertex* getOtherVertex(Rec2DVertex*, Rec2DVertex*) const;
     static Rec2DEdge* getCommonEdge(Rec2DElement*, Rec2DElement*);
     
@@ -411,20 +442,31 @@ class Rec2DElement {
 class Rec2DNode {
   private :
     Rec2DNode *_father;
-    Rec2DNode *_son[3];
+    Rec2DNode *_son[REC2D_NUM_SON];
     Rec2DAction *_ra;
-    double _globalValue, _bestEndGlobVal;
+    double _globalQuality, _bestEndGlobQual;
     int _remainingTri;
     
+    Rec2DDataChange *_dataChange;
+    
   public :
-    Rec2DNode(Rec2DNode *father, Rec2DAction*, double &bestEndGlobVal);
+    Rec2DNode(Rec2DNode *father, Rec2DAction*,
+              double &bestEndGlobQual, int depth = -1);
+    ~Rec2DNode();
+    
+    Rec2DNode* selectBestNode();
+    void recoverSequence();
+    void rmvSon(Rec2DNode*);
+    void develop(int depth, double &bestEndGlobQual);
+    inline bool hasSon() {return _son[0];}
+    bool hasAction();
     
     bool operator<(Rec2DNode&);
     inline Rec2DNode* getFather() {return _father;}
-    //inline int getNum() {return _ra->getNum();}
     inline Rec2DAction* getAction() {return _ra;}
-    inline double getGlobVal() {return _globalValue;}
+    inline double getGlobQual() {return _globalQuality;}
     inline int getNumTri() {return _remainingTri;}
 };
 
+
 #endif
\ No newline at end of file