Ticket #431: 23a4fa3a7e62.patch

lemon/circulation.h

@@ -236 +236 @@
     typedef typename Digraph::template NodeMap<Value> ExcessMap;
     ExcessMap* _excess;
 
+    typename Digraph::template NodeMap<Arc> _next_arc;
+  
     Tolerance _tol;
     int _el;
 
@@ -344 +346 @@
                 const UpperMap &upper, const SupplyMap &supply)
       : _g(graph), _lo(&lower), _up(&upper), _supply(&supply),
         _flow(NULL), _local_flow(false), _level(NULL), _local_level(false),
-        _excess(NULL) {}
+        _excess(NULL), _next_arc(graph) {}
 
     /// Destructor.
     ~Circulation() {
@@ -574 +576 @@
       Node act;
       Node bact=INVALID;
       Node last_activated=INVALID;
+
+      for(NodeIt n(_g);n!=INVALID;++n)
+        if((_next_arc[n]=OutArcIt(_g,n))==INVALID)
+          _next_arc[n]=InArcIt(_g,n);
+      
       while((act=_level->highestActive())!=INVALID) {
         int actlevel=(*_level)[act];
-        int mlevel=_node_num;
         Value exc=(*_excess)[act];
 
-        for(OutArcIt e(_g,act);e!=INVALID; ++e) {
+        Arc next_a = _next_arc[act];
+        if(next_a == INVALID) goto next_l;
+        if (_g.source(next_a)!=act) 
+          goto next_in;
+        for(OutArcIt e(_g,next_a);e!=INVALID; ++e) {
           Node v = _g.target(e);
           Value fc=(*_up)[e]-(*_flow)[e];
           if(!_tol.positive(fc)) continue;
@@ -591 +601 @@
                 _level->activate(v);
               (*_excess)[act] = 0;
               _level->deactivate(act);
+              _next_arc[act]=e;
               goto next_l;
             }
             else {
@@ -601 +612 @@
               exc-=fc;
             }
           }
-          else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
         }
-        for(InArcIt e(_g,act);e!=INVALID; ++e) {
+        next_a = InArcIt(_g,act);
+      next_in:
+        for(InArcIt e(_g,next_a);e!=INVALID; ++e) {
           Node v = _g.source(e);
           Value fc=(*_flow)[e]-(*_lo)[e];
           if(!_tol.positive(fc)) continue;
@@ -615 +627 @@
                 _level->activate(v);
               (*_excess)[act] = 0;
               _level->deactivate(act);
+              _next_arc[act]=e;
               goto next_l;
             }
             else {
@@ -625 +638 @@
               exc-=fc;
             }
           }
-          else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
         }
 
+        if((_next_arc[act]=OutArcIt(_g,act))==INVALID)
+          _next_arc[act]=InArcIt(_g,act);
+
         (*_excess)[act] = exc;
         if(!_tol.positive(exc)) _level->deactivate(act);
-        else if(mlevel==_node_num) {
-          _level->liftHighestActiveToTop();
-          _el = _node_num;
-          return false;
-        }
         else {
-          _level->liftHighestActive(mlevel+1);
-          if(_level->onLevel(actlevel)==0) {
-            _el = actlevel;
+          int mlevel=_node_num;
+          Node rn;
+          for(OutArcIt e(_g,act);e!=INVALID; ++e)
+            if(_tol.positive((*_up)[e]-(*_flow)[e]) &&
+               (*_level)[rn = _g.runningNode(e)]<mlevel)
+              mlevel=(*_level)[rn];
+          for(InArcIt e(_g,act);e!=INVALID; ++e)
+            if(_tol.positive((*_flow)[e]-(*_lo)[e]) &&
+               (*_level)[rn = _g.runningNode(e)]<mlevel)
+              mlevel=(*_level)[rn];
+
+          if(mlevel==_node_num) {
+            _level->liftHighestActiveToTop();
+            _el = _node_num;
             return false;
           }
+          else {
+            _level->liftHighestActive(mlevel+1);
+            if(_level->onLevel(actlevel)==0) {
+              _el = actlevel;
+              return false;
+            }
+          }
         }
       next_l:
         ;