[Lemon-commits] [lemon_svn] jacint: r178 - hugo/trunk/src/work/jacint
Lemon SVN
svn at lemon.cs.elte.hu
Mon Nov 6 20:37:57 CET 2006
Author: jacint
Date: Thu Feb 26 12:38:51 2004
New Revision: 178
Added:
hugo/trunk/src/work/jacint/dimacs_jgraph.hh
hugo/trunk/src/work/jacint/j_graph.h
hugo/trunk/src/work/jacint/preflow_jgraph.cc
hugo/trunk/src/work/jacint/preflow_jgraph.h
Modified:
hugo/trunk/src/work/jacint/makefile
Log:
Alpar SmartGraph-janak atirasa
Added: hugo/trunk/src/work/jacint/dimacs_jgraph.hh
==============================================================================
--- (empty file)
+++ hugo/trunk/src/work/jacint/dimacs_jgraph.hh Thu Feb 26 12:38:51 2004
@@ -0,0 +1,61 @@
+#ifndef DIMACS_HH
+#define DIMACS_HH
+
+#include <iostream>
+#include <string>
+#include <vector>
+
+namespace hugo {
+
+ template<typename Graph, typename CapacityMap>
+ void readDimacsMaxFlow(std::istream& is, Graph &G, typename Graph::TrivNodeIt &s, typename Graph::TrivNodeIt &t, CapacityMap& capacity) {
+ G.clear();
+ int cap;
+ char d;
+ std::string problem;
+ char c;
+ int i, j;
+ std::string str;
+ int n, m;
+ std::vector<typename Graph::TrivNodeIt> nodes;
+ while (is>>c) {
+ switch (c) {
+ case 'c': //comment
+ getline(is, str);
+ break;
+ case 't': //type
+ getline(is, str);
+ break;
+ case 'p': //problem definition
+ is >> problem >> n >> m;
+ getline(is, str);
+ nodes.resize(n+1);
+ for (int k=1; k<=n; ++k) nodes[k]=G.addNode();
+ break;
+ case 'n': //node definition
+ if (problem=="sp") { //shortest path problem
+ is >> i;
+ getline(is, str);
+ s=nodes[i];
+ }
+ if (problem=="max") { //max flow problem
+ is >> i >> d;
+ getline(is, str);
+ if (d=='s') s=nodes[i];
+ if (d=='t') t=nodes[i];
+ }
+ break;
+ case 'a':
+ is >> i >> j >> cap;
+ getline(is, str);
+ typename Graph::TrivEdgeIt e=G.addEdge(nodes[i], nodes[j]);
+ capacity.update();
+ capacity.set(e, cap);
+ break;
+ }
+ }
+ }
+
+} //namespace marci
+
+#endif //DIMACS_HH
Added: hugo/trunk/src/work/jacint/j_graph.h
==============================================================================
--- (empty file)
+++ hugo/trunk/src/work/jacint/j_graph.h Thu Feb 26 12:38:51 2004
@@ -0,0 +1,363 @@
+// -*- mode:C++ -*-
+
+#ifndef J_GRAPH_H
+#define J_GRAPH_H
+
+#include <iostream>
+#include <vector>
+
+namespace hugo {
+
+ class JGraph {
+
+ struct NodeT
+ {
+ int in, out;
+ NodeT() : in(), out() {}
+ };
+
+ struct EdgeT
+ {
+ int head, tail, in, out;
+ };
+
+
+ std::vector<NodeT> nodes;
+ std::vector<EdgeT> edges;
+
+
+ /* template <typename Key> class DynMapBase
+ {
+ protected:
+ SmartGraph* G;
+ public:
+ virtual void add(const Key k) = NULL;
+ virtual void erase(const Key k) = NULL;
+ DynMapBase(SmartGraph &_G) : G(&_G) {}
+ virtual ~DynMapBase() {}
+ friend class SmartGraph;
+ };
+
+ template <typename T> class DynEdgeMap;
+ template <typename T> class DynEdgeMap;
+ */
+
+ public:
+
+ class TrivNodeIt;
+ class TrivEdgeIt;
+ class NodeIt;
+ class EdgeIt;
+ class OutEdgeIt;
+ class InEdgeIt;
+
+ /* protected:
+ std::vector<DynMapBase<NodeIt> * > dyn_node_maps;
+ std::vector<DynMapBase<EdgeIt> * > dyn_edge_maps;
+ */
+
+ template <typename T> class NodeMap;
+ template <typename T> class EdgeMap;
+
+ public:
+
+ JGraph() : nodes(1), edges(1) {
+ edges[0].head=0;
+ edges[0].tail=0;
+ edges[0].in=0; //numEdges (numNodes is maintained in nodes[0].in)
+ edges[0].out=0;
+ }
+
+
+ /* ~SmartGraph()
+ {
+ for(std::vector<DynMapBase<NodeIt> * >::iterator i=dyn_node_maps.begin();
+ i!=dyn_node_maps.end(); ++i) (**i).G=NULL;
+ for(std::vector<DynMapBase<EdgeIt> * >::iterator i=dyn_edge_maps.begin();
+ i!=dyn_edge_maps.end(); ++i) (**i).G=NULL;
+ }*/
+
+ int numNodes() const { return nodes[0].in; }
+ int numEdges() const { return edges[0].in; }
+
+ TrivNodeIt tail(TrivEdgeIt e) const { return edges[e.n].tail; }
+ TrivNodeIt head(TrivEdgeIt e) const { return edges[e.n].head; }
+
+ /* EachNodeIt& getFirst(EachNodeIt& v) const {
+ v=EachNodeIt(*this); return v; }
+ EachEdgeIt& getFirst(EachEdgeIt& e) const {
+ e=EachEdgeIt(*this); return e; }
+ OutEdgeIt& getFirst(OutEdgeIt& e, const NodeIt v) const {
+ e=OutEdgeIt(*this,v); return e; }
+ InEdgeIt& getFirst(InEdgeIt& e, const NodeIt v) const {
+ e=InEdgeIt(*this,v); return e; }
+ */
+
+ NodeIt firstNode() const {
+ return NodeIt( numNodes() );
+ }
+ EdgeIt firstEdge() const {
+ return EdgeIt( numEdges() );
+ }
+ InEdgeIt firstIn(TrivNodeIt v) const {
+ return InEdgeIt(nodes[v.n].in);
+ }
+ OutEdgeIt firstOut(TrivNodeIt v) const {
+ return OutEdgeIt(nodes[v.n].out);
+ }
+
+
+
+ void next(NodeIt& it) const {--it.n;}
+ void next(OutEdgeIt& it) const { it.n=edges[it.n].out; }
+ void next(InEdgeIt& it) const { it.n=edges[it.n].in; }
+ void next(EdgeIt& it) const {--it.n;}
+
+ int id(TrivNodeIt v) const { return v.n; }
+ int id(TrivEdgeIt e) const { return e.n; }
+
+ TrivNodeIt addNode() {
+ TrivNodeIt n;
+ n.n=++nodes[0].in;
+ nodes.push_back(NodeT()); //FIXME
+
+ /* for(std::vector<DynMapBase<NodeIt> * >::iterator i=dyn_node_maps.begin();
+ i!=dyn_node_maps.end(); ++i) (**i).add(n.n);*/
+
+ return n;
+ }
+
+
+ TrivEdgeIt addEdge(TrivNodeIt u, TrivNodeIt v) {
+ if ( u && v ) {
+ TrivEdgeIt e;
+ e.n=++edges[0].in;
+ edges.push_back(EdgeT()); //FIXME: Hmmm...
+ edges[e.n].tail=u.n; edges[e.n].head=v.n;
+ edges[e.n].out=nodes[u.n].out;
+ edges[e.n].in=nodes[v.n].in;
+ nodes[u.n].out=nodes[v.n].in=e.n;
+ return e;
+ }
+ else return TrivEdgeIt();
+
+ /* for(std::vector<DynMapBase<EdgeIt> * >::iterator i=dyn_edge_maps.begin();
+ i!=dyn_edge_maps.end(); ++i) (**i).add(e.n);*/
+
+
+ }
+
+
+ void clear() {
+ nodes.resize(1);
+ edges.resize(1); edges[0].in=0;
+ }
+
+
+ class TrivNodeIt {
+ friend class JGraph;
+ template <typename T> friend class NodeMap;
+
+ friend class TrivEdgeIt;
+ friend class OutEdgeIt;
+ friend class InEdgeIt;
+
+ protected:
+ int n;
+ public:
+ TrivNodeIt() : n() {}
+ TrivNodeIt(int number) {n=number;}
+ bool operator==(const TrivNodeIt i) const {return n==i.n;}
+ bool operator!=(const TrivNodeIt i) const {return n!=i.n;}
+
+ operator bool() const { return n!=0; }
+
+ };
+
+
+ class NodeIt : public TrivNodeIt {
+ friend class JGraph;
+ public:
+ NodeIt() : TrivNodeIt() { }
+ NodeIt(int i) : TrivNodeIt(i) { }
+ operator bool() const { return n!=0; }
+ };
+
+
+ class TrivEdgeIt {
+ friend class JGraph;
+ template <typename T> friend class EdgeMap;
+
+ friend class TrivNodeIt;
+ friend class NodeIt;
+ protected:
+ int n;
+ public:
+ TrivEdgeIt() : n() { }
+ TrivEdgeIt(int number) {n=number;}
+ bool operator==(const TrivEdgeIt i) const {return n==i.n;}
+ bool operator!=(const TrivEdgeIt i) const {return n!=i.n;}
+ operator bool() const { return n!=0; }
+
+ };
+
+
+ class EdgeIt : public TrivEdgeIt {
+ friend class JGraph;
+ public:
+ EdgeIt() : TrivEdgeIt() { }
+ EdgeIt(int i) : TrivEdgeIt(i) { }
+ operator bool() const { return n!=0; }
+ };
+
+
+ class OutEdgeIt : public TrivEdgeIt {
+ friend class JGraph;
+ public:
+ OutEdgeIt() : TrivEdgeIt() { }
+ OutEdgeIt(int i) : TrivEdgeIt(i) { }
+ };
+
+
+ class InEdgeIt : public TrivEdgeIt {
+ friend class JGraph;
+ public:
+ InEdgeIt() : TrivEdgeIt() { }
+ InEdgeIt(int i) : TrivEdgeIt(i) { }
+ };
+
+
+ // Map types
+
+ template <typename T>
+ class NodeMap {
+ const JGraph& G;
+ std::vector<T> container;
+ public:
+ NodeMap(const JGraph& _G) : G(_G), container( G.numNodes()+1 ) { }
+ NodeMap(const JGraph& _G, T a) :
+ G(_G), container(G.numNodes()+1, a) { }
+ void set(TrivNodeIt n, T a) { container[n.n]=a; }
+ T get(TrivNodeIt n) const { return container[n.n]; }
+ /*T& operator[](NodeIt n) { return container[n.n]; }
+ const T& operator[](NodeIt n) const { return container[n.n]; }*/
+ void update() { container.resize(G.numNodes()+1); }
+ void update(T a) { container.resize(G.numNodes()+1, a); }
+ };
+
+ template <typename T>
+ class EdgeMap {
+ const JGraph& G;
+ std::vector<T> container;
+ public:
+ EdgeMap(const JGraph& _G) : G(_G), container(G.numEdges()+1) { }
+ EdgeMap(const JGraph& _G, T a) :
+ G(_G), container(G.numEdges()+1, a) { }
+ void set(TrivEdgeIt e, T a) { container[e.n]=a; }
+ T get(TrivEdgeIt e) const { return container[e.n]; }
+ /*T& operator[](EdgeIt e) { return container[e.n]; }
+ const T& operator[](EdgeIt e) const { return container[e.n]; } */
+ void update() { container.resize(G.numEdges()+1); }
+ void update(T a) { container.resize(G.numEdges()+1, a); }
+ };
+
+ /*template <typename T> class DynNodeMap : public DynMapBase<NodeIt>
+ {
+ std::vector<T> container;
+
+ public:
+ typedef T ValueType;
+ typedef NodeIt KeyType;
+
+ DynNodeMap(JGraph &_G) :
+ DynMapBase<NodeIt>(_G), container(_G.maxNodeId())
+ {
+ //FIXME: What if there are empty Id's?
+ //FIXME: Can I use 'this' in a constructor?
+ G->dyn_node_maps.push_back(this);
+ }
+ ~DynNodeMap()
+ {
+ if(G) {
+ std::vector<DynMapBase<NodeIt>* >::iterator i;
+ for(i=G->dyn_node_maps.begin();
+ i!=G->dyn_node_maps.end() && *i!=this; ++i) ;
+ //if(*i==this) G->dyn_node_maps.erase(i); //FIXME: Way too slow...
+ //A better way to do that: (Is this really important?)
+ if(*i==this) {
+ *i=G->dyn_node_maps.back();
+ G->dyn_node_maps.pop_back();
+ }
+ }
+ }
+
+ void add(const NodeIt k)
+ {
+ if(k.n>=container.size()) container.resize(k.n+1);
+ }
+ void erase(const NodeIt k)
+ {
+ //FIXME: Please implement me.
+ }
+
+ void set(NodeIt n, T a) { container[n.n]=a; }
+ T get(NodeIt n) const { return container[n.n]; }
+ T& operator[](NodeIt n) { return container[n.n]; }
+ const T& operator[](NodeIt n) const { return container[n.n]; }
+
+ void update() {} //Useless for DynMaps
+ void update(T a) {} //Useless for DynMaps
+ };
+
+ template <typename T> class DynEdgeMap : public DynMapBase<EdgeIt>
+ {
+ std::vector<T> container;
+
+ public:
+ typedef T ValueType;
+ typedef EdgeIt KeyType;
+
+ DynEdgeMap(JGraph &_G) :
+ DynMapBase<EdgeIt>(_G), container(_G.maxEdgeId())
+ {
+ //FIXME: What if there are empty Id's?
+ //FIXME: Can I use 'this' in a constructor?
+ G->dyn_edge_maps.push_back(this);
+ }
+ ~DynEdgeMap()
+ {
+ if(G) {
+ std::vector<DynMapBase<EdgeIt>* >::iterator i;
+ for(i=G->dyn_edge_maps.begin();
+ i!=G->dyn_edge_maps.end() && *i!=this; ++i) ;
+ //if(*i==this) G->dyn_edge_maps.erase(i); //Way too slow...
+ //A better way to do that: (Is this really important?)
+ if(*i==this) {
+ *i=G->dyn_edge_maps.back();
+ G->dyn_edge_maps.pop_back();
+ }
+ }
+ }
+
+ void add(const EdgeIt k)
+ {
+ if(k.n>=int(container.size())) container.resize(k.n+1);
+ }
+ void erase(const EdgeIt k)
+ {
+ //FIXME: Please implement me.
+ }
+
+ void set(EdgeIt n, T a) { container[n.n]=a; }
+ T get(EdgeIt n) const { return container[n.n]; }
+ T& operator[](EdgeIt n) { return container[n.n]; }
+ const T& operator[](EdgeIt n) const { return container[n.n]; }
+
+ void update() {} //Useless for DynMaps
+ void update(T a) {} //Useless for DynMaps
+ };*/
+
+ };
+} //namespace hugo
+
+#endif //J_GRAPH_H
Modified: hugo/trunk/src/work/jacint/makefile
==============================================================================
--- hugo/trunk/src/work/jacint/makefile (original)
+++ hugo/trunk/src/work/jacint/makefile Thu Feb 26 12:38:51 2004
@@ -3,36 +3,44 @@
CXXFLAGS = -W -Wall -ansi -pedantic
LEDAROOT = /ledasrc/LEDA-4.1
-BINARIES = preflow preflow_param preflow_hl0 preflow_hl1 preflow_hl2 preflow_hl3 preflow_hl4 preflow_max_flow
+BINARIES = preflow_jgraph
all: $(BINARIES)
makefile: .depend
sinclude .depend
+
+
+preflow_jgraph:
+ $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../jacint -o preflow_jgraph preflow_jgraph.cc
+
+preflowalpar:
+ $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../alpar -I../jacint -o preflowalpar preflowalpar.cc
+
preflow_max_flow:
- $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../marci -I../jacint -o preflow_max_flow preflow_max_flow.cc
+ $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../jacint -o preflow_max_flow preflow_max_flow.cc
preflow:
- $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../marci -I../jacint -o preflow preflow.cc
+ $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../jacint -o preflow preflow.cc
preflow_hl0:
- $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../marci -I../jacint -o preflow_hl0 preflow_hl0.cc
+ $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../jacint -o preflow_hl0 preflow_hl0.cc
preflow_hl1:
- $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../marci -I../jacint -o preflow_hl1 preflow_hl1.cc
+ $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../jacint -o preflow_hl1 preflow_hl1.cc
preflow_hl2:
- $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../marci -I../jacint -o preflow_hl2 preflow_hl2.cc
+ $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../jacint -o preflow_hl2 preflow_hl2.cc
preflow_hl3:
- $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../marci -I../jacint -o preflow_hl3 preflow_hl3.cc
+ $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../jacint -o preflow_hl3 preflow_hl3.cc
preflow_hl4:
- $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../marci -I../jacint -o preflow_hl4 preflow_hl4.cc
+ $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../jacint -o preflow_hl4 preflow_hl4.cc
preflow_param:
- $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../marci -I../jacint -o preflow_param preflow_param.cc
+ $(CXX3) $(CXXFLAGS) -O3 -I. -I.. -I../jacint -o preflow_param preflow_param.cc
clean:
$(RM) *.o $(BINARIES) .depend
Added: hugo/trunk/src/work/jacint/preflow_jgraph.cc
==============================================================================
--- (empty file)
+++ hugo/trunk/src/work/jacint/preflow_jgraph.cc Thu Feb 26 12:38:51 2004
@@ -0,0 +1,72 @@
+#include <iostream>
+#include <fstream>
+
+#include <j_graph.h>
+#include <dimacs_jgraph.hh>
+#include <preflow_jgraph.h>
+#include <time_measure.h>
+
+using namespace hugo;
+
+// Use a DIMACS max flow file as stdin.
+// read_dimacs_demo < dimacs_max_flow_file
+int main(int, char **) {
+ typedef JGraph::TrivNodeIt TrivNodeIt;
+ typedef JGraph::EdgeIt EdgeIt;
+
+ JGraph G;
+ TrivNodeIt s, t;
+ JGraph::EdgeMap<int> cap(G);
+ readDimacsMaxFlow(std::cin, G, s, t, cap);
+
+ std::cout << "preflow demo jacint ..." << std::endl;
+
+ double mintime=1000000;
+
+ for ( int i=1; i!=11; ++i ) {
+ double pre_time=currTime();
+ preflow<JGraph, int> max_flow_test(G, s, t, cap);
+ double post_time=currTime();
+ if ( mintime > post_time-pre_time ) mintime = post_time-pre_time;
+ }
+
+ double pre_time=currTime();
+ preflow<JGraph, int> max_flow_test(G, s, t, cap);
+ double post_time=currTime();
+
+ JGraph::NodeMap<bool> cut(G);
+ max_flow_test.minCut(cut);
+ int min_cut_value=0;
+ for(EdgeIt e=G.firstEdge(); e; G.next(e)) {
+ if (cut.get(G.tail(e)) && !cut.get(G.head(e))) min_cut_value+=cap.get(e);
+ }
+
+ JGraph::NodeMap<bool> cut1(G);
+ max_flow_test.minMinCut(cut1);
+ int min_min_cut_value=0;
+ for(EdgeIt e=G.firstEdge(); e; G.next(e)) {
+ if (cut.get(G.tail(e)) && !cut.get(G.head(e)))
+ min_min_cut_value+=cap.get(e);
+ }
+
+ JGraph::NodeMap<bool> cut2(G);
+ max_flow_test.maxMinCut(cut2);
+ int max_min_cut_value=0;
+ for(EdgeIt e=G.firstEdge(); e; G.next(e)) {
+ if (cut2.get(G.tail(e)) && !cut2.get(G.head(e)))
+ max_min_cut_value+=cap.get(e);
+ }
+
+ std::cout << "min time of 10 runs: " << mintime << " sec"<< std::endl;
+ std::cout << "phase 0: " << max_flow_test.time-pre_time
+ << " sec"<< std::endl;
+ std::cout << "phase 1: " << post_time-max_flow_test.time
+ << " sec"<< std::endl;
+ std::cout << "flow value: "<< max_flow_test.maxFlow() << std::endl;
+ std::cout << "min cut value: "<< min_cut_value << std::endl;
+ std::cout << "min min cut value: "<< min_min_cut_value << std::endl;
+ std::cout << "max min cut value: "<< max_min_cut_value <<
+ std::endl<< std::endl;
+
+ return 0;
+}
Added: hugo/trunk/src/work/jacint/preflow_jgraph.h
==============================================================================
--- (empty file)
+++ hugo/trunk/src/work/jacint/preflow_jgraph.h Thu Feb 26 12:38:51 2004
@@ -0,0 +1,536 @@
+// -*- C++ -*-
+/*
+preflow.h with 'j_graph interface'
+by jacint.
+Heuristics:
+ 2 phase
+ gap
+ list 'level_list' on the nodes on level i implemented by hand
+ stack 'active' on the active nodes on level i implemented by hand
+ runs heuristic 'highest label' for H1*n relabels
+ runs heuristic 'bound decrease' for H0*n relabels, starts with 'highest label'
+
+Parameters H0 and H1 are initialized to 20 and 10.
+
+The best preflow I could ever write.
+
+The constructor runs the algorithm.
+
+Members:
+
+T maxFlow() : returns the value of a maximum flow
+
+T flowOnEdge(EdgeIt e) : for a fixed maximum flow x it returns x(e)
+
+FlowMap Flow() : returns the fixed maximum flow x
+
+void minMinCut(CutMap& M) : sets M to the characteristic vector of the
+ minimum min cut. M should be a map of bools initialized to false.
+
+void maxMinCut(CutMap& M) : sets M to the characteristic vector of the
+ maximum min cut. M should be a map of bools initialized to false.
+
+void minCut(CutMap& M) : sets M to the characteristic vector of
+ a min cut. M should be a map of bools initialized to false.
+
+*/
+
+#ifndef PREFLOW_H
+#define PREFLOW_H
+
+#define H0 20
+#define H1 1
+
+#include <vector>
+#include <queue>
+
+#include<iostream>
+
+#include <time_measure.h>
+
+namespace hugo {
+
+ template <typename Graph, typename T,
+ typename FlowMap=typename Graph::EdgeMap<T>,
+ typename CapMap=typename Graph::EdgeMap<T> >
+ class preflow {
+
+ typedef typename Graph::TrivNodeIt NodeIt;
+ typedef typename Graph::TrivEdgeIt EdgeIt;
+ typedef typename Graph::NodeIt EachNodeIt;
+ typedef typename Graph::OutEdgeIt OutEdgeIt;
+ typedef typename Graph::InEdgeIt InEdgeIt;
+
+ Graph& G;
+ NodeIt s;
+ NodeIt t;
+ FlowMap flow;
+ CapMap& capacity;
+ T value;
+
+ public:
+ double time;
+ preflow(Graph& _G, NodeIt _s, NodeIt _t, CapMap& _capacity ) :
+ G(_G), s(_s), t(_t), flow(_G, 0), capacity(_capacity)
+ {
+
+ bool phase=0; //phase 0 is the 1st phase, phase 1 is the 2nd
+ int n=G.numNodes();
+ int heur0=(int)(H0*n); //time while running 'bound decrease'
+ int heur1=(int)(H1*n); //time while running 'highest label'
+ int heur=heur1; //starting time interval (#of relabels)
+ bool what_heur=1;
+ /*
+ what_heur is 0 in case 'bound decrease'
+ and 1 in case 'highest label'
+ */
+ bool end=false;
+ /*
+ Needed for 'bound decrease', 'true'
+ means no active nodes are above bound b.
+ */
+ int relabel=0;
+ int k=n-2; //bound on the highest level under n containing a node
+ int b=k; //bound on the highest level under n of an active node
+
+ typename Graph::NodeMap<int> level(G,n);
+ typename Graph::NodeMap<T> excess(G);
+
+ std::vector<NodeIt> active(n);
+ typename Graph::NodeMap<NodeIt> next(G);
+ //Stack of the active nodes in level i < n.
+ //We use it in both phases.
+
+ typename Graph::NodeMap<NodeIt> left(G);
+ typename Graph::NodeMap<NodeIt> right(G);
+ std::vector<NodeIt> level_list(n);
+ /*
+ List of the nodes in level i<n.
+ */
+
+ /*Reverse_bfs from t, to find the starting level.*/
+ level.set(t,0);
+ std::queue<NodeIt> bfs_queue;
+ bfs_queue.push(t);
+
+ while (!bfs_queue.empty()) {
+
+ NodeIt v=bfs_queue.front();
+ bfs_queue.pop();
+ int l=level.get(v)+1;
+
+ for(InEdgeIt e=G.firstIn(v); e; G.next(e)) {
+ NodeIt w=G.tail(e);
+ if ( level.get(w) == n && w != s ) {
+ bfs_queue.push(w);
+ NodeIt first=level_list[l];
+ if ( first ) left.set(first,w);
+ right.set(w,first);
+ level_list[l]=w;
+ level.set(w, l);
+ }
+ }
+ }
+
+ level.set(s,n);
+
+
+ /* Starting flow. It is everywhere 0 at the moment. */
+ for(OutEdgeIt e=G.firstOut(s); e; G.next(e))
+ {
+ T c=capacity.get(e);
+ if ( c == 0 ) continue;
+ NodeIt w=G.head(e);
+ if ( level.get(w) < n ) {
+ if ( excess.get(w) == 0 && w!=t ) {
+ next.set(w,active[level.get(w)]);
+ active[level.get(w)]=w;
+ }
+ flow.set(e, c);
+ excess.set(w, excess.get(w)+c);
+ }
+ }
+
+ /*
+ End of preprocessing
+ */
+
+
+
+ /*
+ Push/relabel on the highest level active nodes.
+ */
+ while ( true ) {
+
+ if ( b == 0 ) {
+ if ( phase ) break;
+
+ if ( !what_heur && !end && k > 0 ) {
+ b=k;
+ end=true;
+ } else {
+ phase=1;
+ time=currTime();
+ level.set(s,0);
+ std::queue<NodeIt> bfs_queue;
+ bfs_queue.push(s);
+
+ while (!bfs_queue.empty()) {
+
+ NodeIt v=bfs_queue.front();
+ bfs_queue.pop();
+ int l=level.get(v)+1;
+
+ for(InEdgeIt e=G.firstIn(v); e; G.next(e)) {
+ if ( capacity.get(e) == flow.get(e) ) continue;
+ NodeIt u=G.tail(e);
+ if ( level.get(u) >= n ) {
+ bfs_queue.push(u);
+ level.set(u, l);
+ if ( excess.get(u) > 0 ) {
+ next.set(u,active[l]);
+ active[l]=u;
+ }
+ }
+ }
+
+ for(OutEdgeIt e=G.firstOut(v); e; G.next(e)) {
+ if ( 0 == flow.get(e) ) continue;
+ NodeIt u=G.head(e);
+ if ( level.get(u) >= n ) {
+ bfs_queue.push(u);
+ level.set(u, l);
+ if ( excess.get(u) > 0 ) {
+ next.set(u,active[l]);
+ active[l]=u;
+ }
+ }
+ }
+ }
+ b=n-2;
+ }
+
+ }
+
+
+ if ( !active[b] ) --b;
+ else {
+ end=false;
+
+ NodeIt w=active[b];
+ active[b]=next.get(w);
+ int lev=level.get(w);
+ T exc=excess.get(w);
+ int newlevel=n; //bound on the next level of w
+
+ for(OutEdgeIt e=G.firstOut(w); e; G.next(e)) {
+
+ if ( flow.get(e) == capacity.get(e) ) continue;
+ NodeIt v=G.head(e);
+ //e=wv
+
+ if( lev > level.get(v) ) {
+ /*Push is allowed now*/
+
+ if ( excess.get(v)==0 && v!=t && v!=s ) {
+ int lev_v=level.get(v);
+ next.set(v,active[lev_v]);
+ active[lev_v]=v;
+ }
+
+ T cap=capacity.get(e);
+ T flo=flow.get(e);
+ T remcap=cap-flo;
+
+ if ( remcap >= exc ) {
+ /*A nonsaturating push.*/
+
+ flow.set(e, flo+exc);
+ excess.set(v, excess.get(v)+exc);
+ exc=0;
+ break;
+
+ } else {
+ /*A saturating push.*/
+
+ flow.set(e, cap);
+ excess.set(v, excess.get(v)+remcap);
+ exc-=remcap;
+ }
+ } else if ( newlevel > level.get(v) ){
+ newlevel = level.get(v);
+ }
+
+ } //for out edges wv
+
+
+ if ( exc > 0 ) {
+ for( InEdgeIt e=G.firstIn(w); e; G.next(e)) {
+
+ if( flow.get(e) == 0 ) continue;
+ NodeIt v=G.tail(e);
+ //e=vw
+
+ if( lev > level.get(v) ) {
+ /*Push is allowed now*/
+
+ if ( excess.get(v)==0 && v!=t && v!=s ) {
+ int lev_v=level.get(v);
+ next.set(v,active[lev_v]);
+ active[lev_v]=v;
+ }
+
+ T flo=flow.get(e);
+
+ if ( flo >= exc ) {
+ /*A nonsaturating push.*/
+
+ flow.set(e, flo-exc);
+ excess.set(v, excess.get(v)+exc);
+ exc=0;
+ break;
+ } else {
+/*A saturating push.*/
+
+ excess.set(v, excess.get(v)+flo);
+ exc-=flo;
+ flow.set(e,0);
+ }
+ } else if ( newlevel > level.get(v) ) {
+ newlevel = level.get(v);
+ }
+ } //for in edges vw
+
+ } // if w still has excess after the out edge for cycle
+
+ excess.set(w, exc);
+
+ /*
+ Relabel
+ */
+
+
+ if ( exc > 0 ) {
+ //now 'lev' is the old level of w
+
+ if ( phase ) {
+ level.set(w,++newlevel);
+ next.set(w,active[newlevel]);
+ active[newlevel]=w;
+ b=newlevel;
+ } else {
+ //unlacing starts
+ NodeIt right_n=right.get(w);
+ NodeIt left_n=left.get(w);
+
+ if ( right_n ) {
+ if ( left_n ) {
+ right.set(left_n, right_n);
+ left.set(right_n, left_n);
+ } else {
+ level_list[lev]=right_n;
+ left.set(right_n, NodeIt());
+ }
+ } else {
+ if ( left_n ) {
+ right.set(left_n, NodeIt());
+ } else {
+ level_list[lev]=NodeIt();
+
+ }
+ }
+ //unlacing ends
+
+ //gapping starts
+ if ( !level_list[lev] ) {
+
+ for (int i=lev; i!=k ; ) {
+ NodeIt v=level_list[++i];
+ while ( v ) {
+ level.set(v,n);
+ v=right.get(v);
+ }
+ level_list[i]=NodeIt();
+ if ( !what_heur ) active[i]=NodeIt();
+ }
+
+ level.set(w,n);
+ b=lev-1;
+ k=b;
+ //gapping ends
+ } else {
+
+ if ( newlevel == n ) level.set(w,n);
+ else {
+ level.set(w,++newlevel);
+ next.set(w,active[newlevel]);
+ active[newlevel]=w;
+ if ( what_heur ) b=newlevel;
+ if ( k < newlevel ) ++k;
+ NodeIt first=level_list[newlevel];
+ if ( first ) left.set(first,w);
+ right.set(w,first);
+ left.set(w,NodeIt());
+ level_list[newlevel]=w;
+ }
+ }
+
+
+ ++relabel;
+ if ( relabel >= heur ) {
+ relabel=0;
+ if ( what_heur ) {
+ what_heur=0;
+ heur=heur0;
+ end=false;
+ } else {
+ what_heur=1;
+ heur=heur1;
+ b=k;
+ }
+ }
+ } //phase 0
+
+
+ } // if ( exc > 0 )
+
+
+ } // if stack[b] is nonempty
+
+ } // while(true)
+
+
+ value = excess.get(t);
+ /*Max flow value.*/
+
+ } //void run()
+
+
+
+
+
+ /*
+ Returns the maximum value of a flow.
+ */
+
+ T maxFlow() {
+ return value;
+ }
+
+
+
+ /*
+ For the maximum flow x found by the algorithm,
+ it returns the flow value on edge e, i.e. x(e).
+ */
+
+ T flowOnEdge(EdgeIt e) {
+ return flow.get(e);
+ }
+
+
+
+ FlowMap Flow() {
+ return flow;
+ }
+
+
+
+ void Flow(FlowMap& _flow ) {
+ for(EachNodeIt v=G.firstNode() ; v; G.next(v))
+ _flow.set(v,flow.get(v));
+ }
+
+
+
+ /*
+ Returns the minimum min cut, by a bfs from s in the residual graph.
+ */
+
+ template<typename _CutMap>
+ void minMinCut(_CutMap& M) {
+
+ std::queue<NodeIt> queue;
+
+ M.set(s,true);
+ queue.push(s);
+
+ while (!queue.empty()) {
+ NodeIt w=queue.front();
+ queue.pop();
+
+ for(OutEdgeIt e=G.firstOut(w) ; e; G.next(e)) {
+ NodeIt v=G.head(e);
+ if (!M.get(v) && flow.get(e) < capacity.get(e) ) {
+ queue.push(v);
+ M.set(v, true);
+ }
+ }
+
+ for(InEdgeIt e=G.firstIn(w) ; e; G.next(e)) {
+ NodeIt v=G.tail(e);
+ if (!M.get(v) && flow.get(e) > 0 ) {
+ queue.push(v);
+ M.set(v, true);
+ }
+ }
+ }
+ }
+
+
+
+ /*
+ Returns the maximum min cut, by a reverse bfs
+ from t in the residual graph.
+ */
+
+ template<typename _CutMap>
+ void maxMinCut(_CutMap& M) {
+
+ std::queue<NodeIt> queue;
+
+ M.set(t,true);
+ queue.push(t);
+
+ while (!queue.empty()) {
+ NodeIt w=queue.front();
+ queue.pop();
+
+ for(InEdgeIt e=G.firstIn(w) ; e; G.next(e)) {
+ NodeIt v=G.tail(e);
+ if (!M.get(v) && flow.get(e) < capacity.get(e) ) {
+ queue.push(v);
+ M.set(v, true);
+ }
+ }
+
+ for(OutEdgeIt e=G.firstOut(w) ; e; G.next(e)) {
+ NodeIt v=G.head(e);
+ if (!M.get(v) && flow.get(e) > 0 ) {
+ queue.push(v);
+ M.set(v, true);
+ }
+ }
+ }
+
+ for(EachNodeIt v=G.firstNode() ; v; G.next(v)) {
+ M.set(v, !M.get(v));
+ }
+
+ }
+
+
+
+ template<typename CutMap>
+ void minCut(CutMap& M) {
+ minMinCut(M);
+ }
+
+
+ };
+}//namespace marci
+#endif
+
+
+
+
More information about the Lemon-commits
mailing list