lemon/dynamic_tree.h
author kpeter
Thu, 28 Feb 2008 02:54:27 +0000
changeset 2581 054566ac0934
parent 2519 a7376f7ed899
permissions -rw-r--r--
Query improvements in the min cost flow algorithms.

- External flow and potential maps can be used.
- New query functions: flow() and potential().
- CycleCanceling also provides dual solution (node potentials).
- Doc improvements.
     1 /* -*- C++ -*-
     2  *
     3  * This file is a part of LEMON, a generic C++ optimization library
     4  *
     5  * Copyright (C) 2003-2008
     6  * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
     7  * (Egervary Research Group on Combinatorial Optimization, EGRES).
     8  *
     9  * Permission to use, modify and distribute this software is granted
    10  * provided that this copyright notice appears in all copies. For
    11  * precise terms see the accompanying LICENSE file.
    12  *
    13  * This software is provided "AS IS" with no warranty of any kind,
    14  * express or implied, and with no claim as to its suitability for any
    15  * purpose.
    16  *
    17  */
    18 
    19 #ifndef LEMON_DYNAMIC_TREE_H
    20 #define LEMON_DYNAMIC_TREE_H
    21 
    22 /// \ingroup auxdata
    23 /// \file
    24 /// \brief The dynamic tree data structure of Sleator and Tarjan.
    25 ///
    26 
    27 #include <vector>
    28 #include <limits>
    29 #include <lemon/tolerance.h>
    30 
    31 namespace lemon {
    32 
    33   /// \ingroup auxdata
    34   ///
    35   /// \brief The dynamic tree data structure of Sleator and Tarjan.
    36   ///
    37   /// This class provides an implementation of the dynamic tree data
    38   /// structure for maintaining a set of node-disjoint rooted
    39   /// trees. Each item has an associated value, and the item with
    40   /// minimum value can be find in \f$O(\log(n)\f$ on the path from a
    41   /// node to the its root, and the items on such path can be
    42   /// increased or decreased globally with a certain value in the same
    43   /// running time. We regard a tree edge as directed toward the root,
    44   /// that is from child to parent. Its structure can be modified by
    45   /// two basic operations: \e link(v,w) adds an edge between a root v
    46   /// and a node w in a different component; \e cut(v) removes the
    47   /// edge between v and its parent.
    48   /// 
    49   /// \param _Value The value type of the items.  
    50   /// \param _ItemIntMap Converts item type of node to integer.
    51   /// \param _Tolerance The tolerance class to handle computation
    52   /// problems.
    53   /// \param _enableSize If true then the data structre manatain the
    54   /// size of each tree. The feature is used in \ref GoldbergTarjan
    55   /// algorithm. The default value is true.
    56   ///
    57   /// \author Hamori Tamas
    58 #ifdef DOXYGEN
    59   template <typename _Value, typename _ItemIntMap, 
    60 	    typename _Tolerance, bool _enableSize>
    61 #else
    62   template <typename _Value, typename _ItemIntMap, 
    63 	    typename _Tolerance = lemon::Tolerance<_Value>,
    64 	    bool _enableSize = true>
    65 #endif
    66   class DynamicTree {
    67   public:
    68     /// \brief The integer map on the items.
    69     typedef _ItemIntMap ItemIntMap;
    70     /// \brief The item type of nodes.
    71     typedef typename ItemIntMap::Key Item;
    72     /// \brief The value type of the algorithms.
    73     typedef _Value Value;
    74     /// \brief The tolerance used by the algorithm.
    75     typedef _Tolerance Tolerance;
    76 
    77   private:
    78   
    79     class ItemData;
    80 
    81     std::vector<ItemData> _data;
    82     ItemIntMap &_iim;
    83     Value _max_value;
    84     Tolerance _tolerance;
    85 
    86   public:
    87     /// \brief The constructor of the class.
    88     ///
    89     /// \param iim The integer map on the items. 
    90     /// \param tolerance Tolerance class.
    91     DynamicTree(ItemIntMap &iim, const Tolerance& tolerance = Tolerance())
    92       : _iim(iim), _max_value(std::numeric_limits<Value>::max()), 
    93 	_tolerance(tolerance) {}
    94   
    95     ~DynamicTree() {}
    96 
    97     /// \brief Clears the data structure
    98     ///
    99     /// Clears the data structure
   100     void clear() {
   101       _data.clear();
   102     }
   103 
   104     /// \brief Sets the tolerance used by algorithm.
   105     ///
   106     /// Sets the tolerance used by algorithm.
   107     void tolerance(const Tolerance& tolerance) const {
   108       _tolerance = tolerance;
   109       return *this;
   110     } 
   111   
   112     /// \brief Returns the tolerance used by algorithm.
   113     ///
   114     /// Returns the tolerance used by algorithm.
   115     const Tolerance& tolerance() const {
   116       return tolerance;
   117     } 
   118   
   119     /// \brief Create a new tree containing a single node with cost zero.
   120     void makeTree(const Item &item) {
   121       _data[makePath(item)].successor = -1;
   122     }
   123     
   124     /// \brief Return the root of the tree containing node with itemtype
   125     /// \e item.
   126     Item findRoot(const Item &item) {
   127       return _data[findTail(expose(_iim[item]))].id;
   128     }
   129     
   130     /// \brief Return the the value of nodes in the tree containing
   131     /// node with itemtype \e item.
   132     int findSize(const Item &item) {
   133       return _data[expose(_iim[item])].size;
   134     }
   135     
   136     /// \brief Return the minimum cost containing node.
   137     /// 
   138     /// Return into \e d the minimum cost on the tree path from \e item
   139     /// to findRoot(item).  Return the last item (closest to its root)
   140     /// on this path of the minimum cost.
   141     Item findCost(const Item &item, Value& d){
   142       return _data[findPathCost(expose(_iim[item]),d)].id;
   143     }
   144     
   145     /// \brief Add \e x value to the cost of every node on the path from
   146     /// \e item to findRoot(item).
   147     void addCost(const Item &item, Value x) {
   148       addPathCost(expose(_iim[item]), x);
   149     }
   150     
   151     /// \brief Combine the trees containing nodes \e item1 and \e item2
   152     /// by adding an edge from \e item1 \e item2.
   153     /// 
   154     /// This operation assumes that \e item1 is root and \e item2 is in
   155     /// a different tree.
   156     void link(const Item &item1, const Item &item2){
   157       int v = _iim[item1];
   158       int w = _iim[item2];
   159       int p = expose(w);
   160       join(-1, expose(v), p);
   161       _data[v].successor = -1;
   162       _data[v].size += _data[p].size;
   163 
   164     }    
   165     
   166     /// \brief Divide the tree containing node \e item into two trees by
   167     /// deleting the edge out of \e item.
   168     /// 
   169     /// This operation assumes that \e item is not a tree root.
   170     void cut(const Item &item) {
   171       int v = _iim[item];
   172       int p, q;
   173       expose(v);
   174       split(p, v, q);
   175       if (p != -1) {
   176 	_data[p].successor = v;
   177       }
   178       _data[v].size -= _data[q].size;
   179       if (q != -1) {
   180 	_data[q].successor = _data[v].successor;
   181       }
   182       _data[v].successor = -1;
   183     }
   184 
   185     ///\brief 
   186     Item parent(const Item &item){
   187       return _data[_iim[item].p].id;
   188     }
   189 
   190     ///\brief Return the upper bound of the costs.
   191     Value maxValue() const {
   192       return _max_value;
   193     }
   194     
   195   private:
   196 
   197     int makePath(const Item &item) {
   198       _iim.set(item, _data.size());
   199       ItemData v(item);
   200       _data.push_back(v);
   201       return _iim[item];
   202     }
   203 
   204     int findPath(int v) {
   205       splay(v);
   206       return v;
   207     }
   208     
   209     int findPathCost(int p, Value &d) {
   210       while ((_data[p].right != -1 && 
   211 	      !_tolerance.less(0, _data[_data[p].right].dmin)) || 
   212 	     (_data[p].left != -1 && _tolerance.less(0, _data[p].dcost))) {
   213 	if (_data[p].right != -1 && 
   214 	    !_tolerance.less(0, _data[_data[p].right].dmin)) {
   215 	  p = _data[p].right;
   216 	} else if (_data[p].left != -1 && 
   217 		   !_tolerance.less(0, _data[_data[p].left].dmin)) {
   218 	  p = _data[p].left;
   219 	}
   220       }
   221       splay(p);
   222       d = _data[p].dmin;
   223       return p; 
   224     }
   225 
   226     int findTail(int p){
   227       while (_data[p].right != -1) {
   228 	p = _data[p].right;
   229       }
   230       splay(p);
   231       return p;
   232     }
   233     
   234     void addPathCost(int p, Value x) {
   235       if (!_tolerance.less(x, _max_value)) {
   236 	_data[p].dmin = x;
   237 	_data[p].dcost = x;
   238       } else if (!_tolerance.less(-x, _max_value)) {
   239 	_data[p].dmin = 0;
   240 	_data[p].dcost = 0;
   241       } else {
   242 	_data[p].dmin += x;
   243       }
   244     }
   245 
   246     void join(int p, int v, int q) {
   247       Value min = _max_value;
   248       Value pmin = _max_value;
   249       Value vmin = _data[v].dmin;
   250       Value qmin = _max_value;
   251       if (p != -1){
   252 	pmin = _data[p].dmin;
   253       }
   254       if (q != -1){
   255 	qmin = _data[q].dmin;
   256       }
   257         
   258       if (_tolerance.less(vmin, qmin)) {
   259 	if (_tolerance.less(vmin,pmin)) {
   260 	  min = vmin;
   261 	} else {
   262 	  min = pmin;
   263 	}
   264       } else if (_tolerance.less(qmin,pmin)) {
   265 	min = qmin;
   266       } else {
   267 	min = pmin;
   268       }
   269 
   270       if (p != -1){
   271 	_data[p].parent = v;
   272 	_data[p].dmin -= min;
   273       }
   274       if (q!=-1){
   275 	_data[q].parent = v;
   276 	if (_tolerance.less(_data[q].dmin,_max_value)) {
   277 	  _data[q].dmin -= min;
   278 	}
   279       }
   280       _data[v].left = p;
   281       _data[v].right = q;
   282       if (_tolerance.less(min,_max_value)) {
   283 	_data[v].dcost = _data[v].dmin - min;
   284       }
   285       _data[v].dmin = min;
   286     }
   287 
   288     void split(int &p, int v, int &q){
   289       splay(v);
   290       p = -1;
   291       if (_data[v].left != -1){
   292 	p = _data[v].left;
   293 	_data[p].dmin += _data[v].dmin;
   294 	_data[p].parent = -1;
   295 	_data[v].left = -1;
   296       }
   297       q = -1;
   298       if (_data[v].right != -1) {
   299 	q=_data[v].right;
   300 	if (_tolerance.less(_data[q].dmin, _max_value)) {
   301 	  _data[q].dmin += _data[v].dmin;
   302 	}
   303 	_data[q].parent = -1;
   304 	_data[v].right = -1;
   305       } 
   306       if (_tolerance.less(_data[v].dcost, _max_value)) {
   307 	_data[v].dmin += _data[v].dcost;
   308 	_data[v].dcost = 0;
   309       } else {
   310 	_data[v].dmin = _data[v].dcost;
   311       }
   312     }
   313  
   314     int expose(int v) {
   315       int p, q, r, w;
   316       p = -1;
   317       while (v != -1) {
   318 	w = _data[findPath(v)].successor;
   319 	split(q, v, r);
   320 	if (q != -1) {
   321 	  _data[q].successor = v;
   322 	}
   323 	join(p, v, r);
   324 	p = v;
   325 	v = w;
   326       }
   327       _data[p].successor = -1;
   328       return p;
   329     }
   330 
   331     void splay(int v) {
   332       while (_data[v].parent != -1) {
   333 	if (v == _data[_data[v].parent].left) {
   334 	  if (_data[_data[v].parent].parent == -1) {
   335 	    zig(v);
   336 	  } else {
   337 	    if (_data[v].parent == _data[_data[_data[v].parent].parent].left) {
   338 	      zig(_data[v].parent);
   339 	      zig(v);
   340 	    } else {
   341 	      zig(v);
   342 	      zag(v);
   343 	    }
   344 	  }
   345 	} else {
   346 	  if (_data[_data[v].parent].parent == -1) {
   347 	    zag(v);
   348 	  } else {
   349 	    if (_data[v].parent == _data[_data[_data[v].parent].parent].left) {
   350 	      zag(v);
   351 	      zig(v);
   352 	    } else {
   353 	      zag(_data[v].parent);
   354 	      zag(v);
   355 	    }
   356 	  }
   357 	}
   358       }
   359     }
   360 
   361 
   362     void zig(int v) {
   363       Value min = _data[_data[v].parent].dmin;
   364       int a = _data[v].parent;
   365         
   366       Value aa = _data[a].dcost;
   367       if (_tolerance.less(aa, _max_value)) { 
   368 	aa += min;
   369       }
   370 
   371 
   372       int b = v;
   373       Value ab = min + _data[b].dmin;
   374       Value bb = _data[b].dcost;
   375       if (_tolerance.less(bb, _max_value)) { 
   376 	bb += ab;
   377       }
   378 
   379       int c = -1;
   380       Value cc = _max_value;
   381       if (_data[a].right != -1) {
   382 	c = _data[a].right;
   383 	cc = _data[c].dmin;
   384 	if (_tolerance.less(cc, _max_value)) {
   385 	  cc += min;
   386 	}
   387       }
   388 
   389       int d = -1;
   390       Value dd = _max_value;
   391       if (_data[v].left != -1){
   392 	d = _data[v].left;
   393 	dd = ab + _data[d].dmin;
   394       }
   395 
   396       int e = -1;
   397       Value ee = _max_value;
   398       if (_data[v].right != -1) {
   399 	e = _data[v].right;
   400 	ee = ab + _data[e].dmin;
   401       }
   402 
   403       Value min2;
   404       if (_tolerance.less(0, _data[b].dmin) || 
   405 	  (e != -1 && !_tolerance.less(0, _data[e].dmin))) {
   406 	min2 = min;
   407       } else {
   408 	if (_tolerance.less(aa, cc)) {
   409 	  if (_tolerance.less(aa, ee)) {
   410 	    min2 = aa;
   411 	  } else {
   412 	    min2 = ee;
   413 	  }
   414 	} else if (_tolerance.less(cc, ee)) {
   415 	  min2 = cc;
   416 	} else {
   417 	  min2 = ee;
   418 	}
   419       }
   420         
   421       _data[a].dcost = aa;
   422       if (_tolerance.less(aa, _max_value)) { 
   423 	_data[a].dcost -= min2;
   424       }
   425       _data[a].dmin = min2;
   426       if (_tolerance.less(min2,_max_value)) { 
   427 	_data[a].dmin -= min; 
   428       }
   429       _data[a].size -= _data[b].size;
   430       _data[b].dcost = bb;
   431       if (_tolerance.less(bb, _max_value)) { 
   432 	_data[b].dcost -= min;
   433       }
   434       _data[b].dmin = min;
   435       _data[b].size += _data[a].size;
   436       if (c != -1) {
   437 	_data[c].dmin = cc;
   438 	if (_tolerance.less(cc, _max_value)) {
   439 	  _data[c].dmin -= min2;
   440 	}
   441       }
   442       if (d != -1) {
   443 	_data[d].dmin = dd - min;
   444 	_data[a].size += _data[d].size;
   445 	_data[b].size -= _data[d].size;
   446       }
   447       if (e != -1) {
   448 	_data[e].dmin = ee - min2;
   449       }
   450         
   451       int w = _data[v].parent;
   452       _data[v].successor = _data[w].successor;
   453       _data[w].successor = -1;
   454       _data[v].parent = _data[w].parent;
   455       _data[w].parent = v;
   456       _data[w].left = _data[v].right;
   457       _data[v].right = w;
   458       if (_data[v].parent != -1){
   459 	if (_data[_data[v].parent].right == w) {
   460 	  _data[_data[v].parent].right = v;
   461 	} else {
   462 	  _data[_data[v].parent].left = v;
   463 	}
   464       }
   465       if (_data[w].left != -1){
   466 	_data[_data[w].left].parent = w;
   467       }
   468     }
   469 
   470 
   471     void zag(int v) {
   472 
   473       Value min = _data[_data[v].parent].dmin;
   474 
   475       int a = _data[v].parent;
   476       Value aa = _data[a].dcost;
   477       if (_tolerance.less(aa, _max_value)) { 
   478 	aa += min;
   479       }
   480         
   481       int b = v;
   482       Value ab = min + _data[b].dmin;
   483       Value bb = _data[b].dcost;
   484       if (_tolerance.less(bb, _max_value)) {
   485 	bb += ab;
   486       }
   487 
   488       int c = -1;
   489       Value cc = _max_value;
   490       if (_data[a].left != -1){
   491 	c = _data[a].left;
   492 	cc = min + _data[c].dmin;
   493       }
   494 
   495       int d = -1;
   496       Value dd = _max_value;
   497       if (_data[v].right!=-1) {
   498 	d = _data[v].right;
   499 	dd = _data[d].dmin;
   500 	if (_tolerance.less(dd, _max_value)) {
   501 	  dd += ab;
   502 	}
   503       }
   504 
   505       int e = -1;
   506       Value ee = _max_value;
   507       if (_data[v].left != -1){
   508 	e = _data[v].left;
   509 	ee = ab + _data[e].dmin;
   510       }
   511 
   512       Value min2;
   513       if (_tolerance.less(0, _data[b].dmin) || 
   514 	  (e != -1 && !_tolerance.less(0, _data[e].dmin))) {
   515 	min2 = min;
   516       } else {
   517 	if (_tolerance.less(aa, cc)) {
   518 	  if (_tolerance.less(aa, ee)) {
   519 	    min2 = aa;
   520 	  } else {
   521 	    min2 = ee;
   522 	  }
   523 	} else if (_tolerance.less(cc, ee)) {
   524 	  min2 = cc;
   525 	} else {
   526 	  min2 = ee;
   527 	}
   528       }
   529       _data[a].dcost = aa;
   530       if (_tolerance.less(aa, _max_value)) { 
   531 	_data[a].dcost -= min2;
   532       }
   533       _data[a].dmin = min2;
   534       if (_tolerance.less(min2, _max_value)) {
   535 	_data[a].dmin -= min;
   536       }
   537       _data[a].size -= _data[b].size;
   538       _data[b].dcost = bb;
   539       if (_tolerance.less(bb, _max_value)) { 
   540 	_data[b].dcost -= min;
   541       }
   542       _data[b].dmin = min;
   543       _data[b].size += _data[a].size;
   544       if (c != -1) {
   545 	_data[c].dmin = cc - min2;
   546       }
   547       if (d != -1) {
   548 	_data[d].dmin = dd;
   549 	_data[a].size += _data[d].size;
   550 	_data[b].size -= _data[d].size;
   551 	if (_tolerance.less(dd, _max_value)) {
   552 	  _data[d].dmin -= min;
   553 	}
   554       }
   555       if (e != -1) {
   556 	_data[e].dmin = ee - min2;
   557       }
   558         
   559       int w = _data[v].parent;
   560       _data[v].successor = _data[w].successor;
   561       _data[w].successor = -1;
   562       _data[v].parent = _data[w].parent;
   563       _data[w].parent = v;
   564       _data[w].right = _data[v].left;
   565       _data[v].left = w;
   566       if (_data[v].parent != -1){
   567 	if (_data[_data[v].parent].left == w) {
   568 	  _data[_data[v].parent].left = v;
   569 	} else {
   570 	  _data[_data[v].parent].right = v;
   571 	}
   572       }
   573       if (_data[w].right != -1){
   574 	_data[_data[w].right].parent = w;
   575       }
   576     }
   577 
   578   private:
   579 
   580     class ItemData {
   581     public:
   582       Item id;
   583       int size;
   584       int successor;
   585       int parent;
   586       int left;
   587       int right;
   588       Value dmin;
   589       Value dcost;
   590         
   591     public:
   592       ItemData(const Item &item)
   593 	: id(item), size(1), successor(), parent(-1), 
   594 	  left(-1), right(-1), dmin(0), dcost(0) {}
   595     };
   596      
   597   };
   598 
   599   template <typename _Value, typename _ItemIntMap, typename _Tolerance>
   600   class DynamicTree<_Value, _ItemIntMap, _Tolerance, false> {
   601   public:
   602     typedef _ItemIntMap ItemIntMap;
   603     typedef typename ItemIntMap::Key Item;
   604     typedef _Value Value;
   605     typedef _Tolerance Tolerance;
   606 
   607   private:
   608   
   609     class ItemData;
   610 
   611     std::vector<ItemData> _data;
   612     ItemIntMap &_iim;
   613     Value _max_value;
   614     Tolerance _tolerance;
   615 
   616   public:
   617     DynamicTree(ItemIntMap &iim, const Tolerance& tolerance = Tolerance())
   618       : _iim(iim), _max_value(std::numeric_limits<Value>::max()), 
   619 	_tolerance(tolerance) {}
   620   
   621     ~DynamicTree() {}
   622 
   623     void clear() {
   624       _data.clear();
   625     }
   626 
   627     void tolerance(const Tolerance& tolerance) const {
   628       _tolerance = tolerance;
   629       return *this;
   630     } 
   631   
   632     const Tolerance& tolerance() const {
   633       return tolerance;
   634     } 
   635   
   636     void makeTree(const Item &item) {
   637       _data[makePath(item)].successor = -1;
   638     }
   639     
   640     Item findRoot(const Item &item) {
   641       return _data[findTail(expose(_iim[item]))].id;
   642     }
   643     
   644     Item findCost(const Item &item, Value& d){
   645       return _data[findPathCost(expose(_iim[item]),d)].id;
   646     }
   647     
   648     void addCost(const Item &item, Value x){
   649       addPathCost(expose(_iim[item]), x);
   650     }
   651     
   652     void link(const Item &item1, const Item &item2){
   653       int v = _iim[item1];
   654       int w = _iim[item2];
   655       int p = expose(w);
   656       join(-1, expose(v), p);
   657       _data[v].successor = -1;
   658     }    
   659     
   660     void cut(const Item &item) {
   661       int v = _iim[item];
   662       int p, q;
   663       expose(v);
   664       split(p, v, q);
   665       if (p != -1) {
   666 	_data[p].successor = v;
   667       }
   668       if (q != -1) {
   669 	_data[q].successor = _data[v].successor;
   670       }
   671       _data[v].successor = -1;
   672     }
   673 
   674     Item parent(const Item &item){
   675       return _data[_iim[item].p].id;
   676     }
   677 
   678     Value maxValue() const {
   679       return _max_value;
   680     }
   681     
   682   private:
   683 
   684     int makePath(const Item &item) {
   685       _iim.set(item, _data.size());
   686       ItemData v(item);
   687       _data.push_back(v);
   688       return _iim[item];
   689     }
   690 
   691     int findPath(int v) {
   692       splay(v);
   693       return v;
   694     }
   695     
   696     int findPathCost(int p, Value &d) {
   697       while ((_data[p].right != -1 && 
   698 	      !_tolerance.less(0, _data[_data[p].right].dmin)) || 
   699 	     (_data[p].left != -1 && _tolerance.less(0, _data[p].dcost))) {
   700 	if (_data[p].right != -1 && 
   701 	    !_tolerance.less(0, _data[_data[p].right].dmin)) {
   702 	  p = _data[p].right;
   703 	} else if (_data[p].left != -1 && 
   704 		   !_tolerance.less(0, _data[_data[p].left].dmin)){
   705 	  p = _data[p].left;
   706 	}
   707       }
   708       splay(p);
   709       d = _data[p].dmin;
   710       return p; 
   711     }
   712 
   713     int findTail(int p) {
   714       while (_data[p].right != -1) {
   715 	p = _data[p].right;
   716       }
   717       splay(p);
   718       return p;
   719     }
   720     
   721     void addPathCost(int p, Value x) {
   722       if (!_tolerance.less(x, _max_value)) {
   723 	_data[p].dmin = x;_data[p].dcost = x;
   724       } else if (!_tolerance.less(-x, _max_value)) {
   725 	_data[p].dmin = 0;
   726 	_data[p].dcost = 0;
   727       } else {
   728 	_data[p].dmin += x;
   729       }
   730     }
   731 
   732     void join(int p, int v, int q) {
   733       Value min = _max_value;
   734       Value pmin = _max_value;
   735       Value vmin = _data[v].dmin;
   736       Value qmin = _max_value;
   737       if (p != -1){
   738 	pmin = _data[p].dmin;
   739       }
   740       if (q != -1){
   741 	qmin = _data[q].dmin;
   742       }
   743         
   744       if (_tolerance.less(vmin, qmin)) {
   745 	if (_tolerance.less(vmin,pmin)) {
   746 	  min = vmin;
   747 	} else {
   748 	  min = pmin;
   749 	}
   750       } else if (_tolerance.less(qmin,pmin)) {
   751 	min = qmin;
   752       } else {
   753 	min = pmin;
   754       }
   755 
   756       if (p != -1){
   757 	_data[p].parent = v;
   758 	_data[p].dmin -= min;
   759       }
   760       if (q != -1){
   761 	_data[q].parent = v;
   762 	if (_tolerance.less(_data[q].dmin,_max_value)) {
   763 	  _data[q].dmin -= min;
   764 	}
   765       }
   766       _data[v].left = p;
   767       _data[v].right = q;
   768       if (_tolerance.less(min, _max_value)) {
   769 	_data[v].dcost = _data[v].dmin - min;
   770       }
   771       _data[v].dmin = min;
   772     }
   773 
   774     void split(int &p, int v, int &q){
   775       splay(v);
   776       p = -1;
   777       if (_data[v].left != -1){
   778 	p = _data[v].left;
   779 	_data[p].dmin += _data[v].dmin;
   780 	_data[p].parent = -1;
   781 	_data[v].left = -1;
   782       }
   783       q = -1;
   784       if (_data[v].right != -1) {
   785 	q=_data[v].right;
   786 	if (_tolerance.less(_data[q].dmin, _max_value)) {
   787 	  _data[q].dmin += _data[v].dmin;
   788 	}
   789 	_data[q].parent = -1;
   790 	_data[v].right = -1;
   791       } 
   792       if (_tolerance.less(_data[v].dcost, _max_value)) {
   793 	_data[v].dmin += _data[v].dcost;
   794 	_data[v].dcost = 0;
   795       } else {
   796 	_data[v].dmin = _data[v].dcost;
   797       }
   798     }
   799  
   800     int expose(int v) {
   801       int p, q, r, w;
   802       p = -1;
   803       while (v != -1) {
   804 	w = _data[findPath(v)].successor;
   805 	split(q, v, r);
   806 	if (q != -1) {
   807 	  _data[q].successor = v;
   808 	}
   809 	join(p, v, r);
   810 	p = v;
   811 	v = w;
   812       }
   813       _data[p].successor = -1;
   814       return p;
   815     }
   816 
   817     void splay(int v) {
   818       while (_data[v].parent != -1) {
   819 	if (v == _data[_data[v].parent].left) {
   820 	  if (_data[_data[v].parent].parent == -1) {
   821 	    zig(v);
   822 	  } else {
   823 	    if (_data[v].parent == _data[_data[_data[v].parent].parent].left) {
   824 	      zig(_data[v].parent);
   825 	      zig(v);
   826 	    } else {
   827 	      zig(v);
   828 	      zag(v);
   829 	    }
   830 	  }
   831 	} else {
   832 	  if (_data[_data[v].parent].parent == -1) {
   833 	    zag(v);
   834 	  } else {
   835 	    if (_data[v].parent == _data[_data[_data[v].parent].parent].left) {
   836 	      zag(v);
   837 	      zig(v);
   838 	    } else {
   839 	      zag(_data[v].parent);
   840 	      zag(v);
   841 	    }
   842 	  }
   843 	}
   844       }
   845     }
   846 
   847 
   848     void zig(int v) {
   849       Value min = _data[_data[v].parent].dmin;
   850       int a = _data[v].parent;
   851         
   852       Value aa = _data[a].dcost;
   853       if (_tolerance.less(aa, _max_value)) { 
   854 	aa+= min;
   855       }
   856 
   857 
   858       int b = v;
   859       Value ab = min + _data[b].dmin;
   860       Value bb = _data[b].dcost;
   861       if (_tolerance.less(bb, _max_value)) { 
   862 	bb+= ab;
   863       }
   864 
   865       int c = -1;
   866       Value cc = _max_value;
   867       if (_data[a].right != -1) {
   868 	c = _data[a].right;
   869 	cc = _data[c].dmin;
   870 	if (_tolerance.less(cc, _max_value)) {
   871 	  cc+=min;
   872 	}
   873       }
   874 
   875       int d = -1;
   876       Value dd = _max_value;
   877       if (_data[v].left != -1){
   878 	d = _data[v].left;
   879 	dd = ab + _data[d].dmin;
   880       }
   881 
   882       int e = -1;
   883       Value ee = _max_value;
   884       if (_data[v].right != -1) {
   885 	e = _data[v].right;
   886 	ee = ab + _data[e].dmin;
   887       }
   888 
   889       Value min2;
   890       if (_tolerance.less(0, _data[b].dmin) || 
   891 	  (e != -1 && !_tolerance.less(0, _data[e].dmin))) {
   892 	min2 = min;
   893       } else {
   894 	if (_tolerance.less(aa, cc)) {
   895 	  if (_tolerance.less(aa, ee)) {
   896 	    min2 = aa;
   897 	  } else {
   898 	    min2 = ee;
   899 	  }
   900 	} else if (_tolerance.less(cc, ee)) {
   901 	  min2 = cc;
   902 	} else {
   903 	  min2 = ee;
   904 	}
   905       }
   906         
   907       _data[a].dcost = aa;
   908       if (_tolerance.less(aa, _max_value)) { 
   909 	_data[a].dcost -= min2;
   910       }
   911       _data[a].dmin = min2;
   912       if (_tolerance.less(min2,_max_value)) { 
   913 	_data[a].dmin -= min; 
   914       }
   915       _data[b].dcost = bb;
   916       if (_tolerance.less(bb, _max_value)) { 
   917 	_data[b].dcost -= min;
   918       }
   919       _data[b].dmin = min;
   920       if (c != -1) {
   921 	_data[c].dmin = cc;
   922 	if (_tolerance.less(cc, _max_value)) {
   923 	  _data[c].dmin -= min2;
   924 	}
   925       }
   926       if (d != -1) {
   927 	_data[d].dmin = dd - min;
   928       }
   929       if (e != -1) {
   930 	_data[e].dmin = ee - min2;
   931       }
   932         
   933       int w = _data[v].parent;
   934       _data[v].successor = _data[w].successor;
   935       _data[w].successor = -1;
   936       _data[v].parent = _data[w].parent;
   937       _data[w].parent = v;
   938       _data[w].left = _data[v].right;
   939       _data[v].right = w;
   940       if (_data[v].parent != -1){
   941 	if (_data[_data[v].parent].right == w) {
   942 	  _data[_data[v].parent].right = v;
   943 	} else {
   944 	  _data[_data[v].parent].left = v;
   945 	}
   946       }
   947       if (_data[w].left != -1){
   948 	_data[_data[w].left].parent = w;
   949       }
   950     }
   951 
   952 
   953     void zag(int v) {
   954 
   955       Value min = _data[_data[v].parent].dmin;
   956 
   957       int a = _data[v].parent;
   958       Value aa = _data[a].dcost;
   959       if (_tolerance.less(aa, _max_value)) { 
   960 	aa += min;
   961       }
   962         
   963       int b = v;
   964       Value ab = min + _data[b].dmin;
   965       Value bb = _data[b].dcost;
   966       if (_tolerance.less(bb, _max_value)) {
   967 	bb += ab;
   968       }
   969 
   970       int c = -1;
   971       Value cc = _max_value;
   972       if (_data[a].left != -1){
   973 	c = _data[a].left;
   974 	cc = min + _data[c].dmin;
   975       }
   976 
   977       int d = -1;
   978       Value dd = _max_value;
   979       if (_data[v].right!=-1) {
   980 	d = _data[v].right;
   981 	dd = _data[d].dmin;
   982 	if (_tolerance.less(dd, _max_value)) {
   983 	  dd += ab;
   984 	}
   985       }
   986 
   987       int e = -1;
   988       Value ee = _max_value;
   989       if (_data[v].left != -1){
   990 	e = _data[v].left;
   991 	ee = ab + _data[e].dmin;
   992       }
   993 
   994       Value min2;
   995       if (_tolerance.less(0, _data[b].dmin) || 
   996 	  (e != -1 && !_tolerance.less(0, _data[e].dmin))) {
   997 	min2 = min;
   998       } else {
   999 	if (_tolerance.less(aa, cc)) {
  1000 	  if (_tolerance.less(aa, ee)) {
  1001 	    min2 = aa;
  1002 	  } else {
  1003 	    min2 = ee;
  1004 	  }
  1005 	} else if (_tolerance.less(cc, ee)) {
  1006 	  min2 = cc;
  1007 	} else {
  1008 	  min2 = ee;
  1009 	}
  1010       }
  1011       _data[a].dcost = aa;
  1012       if (_tolerance.less(aa, _max_value)) { 
  1013 	_data[a].dcost -= min2;
  1014       }
  1015       _data[a].dmin = min2;
  1016       if (_tolerance.less(min2, _max_value)) {
  1017 	_data[a].dmin -= min;
  1018       }
  1019       _data[b].dcost = bb;
  1020       if (_tolerance.less(bb, _max_value)) { 
  1021 	_data[b].dcost -= min;
  1022       }
  1023       _data[b].dmin = min;
  1024       if (c != -1) {
  1025 	_data[c].dmin = cc - min2;
  1026       }
  1027       if (d != -1) {
  1028 	_data[d].dmin = dd;
  1029 	if (_tolerance.less(dd, _max_value)) {
  1030 	  _data[d].dmin -= min;
  1031 	}
  1032       }
  1033       if (e != -1) {
  1034 	_data[e].dmin = ee - min2;
  1035       }
  1036         
  1037       int w = _data[v].parent;
  1038       _data[v].successor = _data[w].successor;
  1039       _data[w].successor = -1;
  1040       _data[v].parent = _data[w].parent;
  1041       _data[w].parent = v;
  1042       _data[w].right = _data[v].left;
  1043       _data[v].left = w;
  1044       if (_data[v].parent != -1){
  1045 	if (_data[_data[v].parent].left == w) {
  1046 	  _data[_data[v].parent].left = v;
  1047 	} else {
  1048 	  _data[_data[v].parent].right = v;
  1049 	}
  1050       }
  1051       if (_data[w].right != -1){
  1052 	_data[_data[w].right].parent = w;
  1053       }
  1054     }
  1055 
  1056   private:
  1057 
  1058     class ItemData {
  1059     public:
  1060       Item id;
  1061       int successor;
  1062       int parent;
  1063       int left;
  1064       int right;
  1065       Value dmin;
  1066       Value dcost;
  1067         
  1068     public:
  1069       ItemData(const Item &item)
  1070 	: id(item), successor(), parent(-1), 
  1071 	  left(-1), right(-1), dmin(0), dcost(0) {}
  1072     };
  1073      
  1074   };
  1075 
  1076 }
  1077 
  1078 #endif