[Lemon-commits] [lemon_svn] marci: r1026 - hugo/trunk/src/work/marci/lp
Lemon SVN
svn at lemon.cs.elte.hu
Mon Nov 6 20:42:47 CET 2006
Author: marci
Date: Tue Aug 17 15:20:46 2004
New Revision: 1026
Added:
hugo/trunk/src/work/marci/lp/lp_solver_wrapper.h
hugo/trunk/src/work/marci/lp/makefile
hugo/trunk/src/work/marci/lp/max_flow_by_lp.cc
Log:
An experimental LPSolverWrapper class which uses glpk. For a short
demo, max flow problems are solved with it. This demo does not
demonstrates, but the main aims of this class are row and column
generation capabilities, i.e. to be a core for easily
implementable branch-and-cut a column generetion algorithms.
Added: hugo/trunk/src/work/marci/lp/lp_solver_wrapper.h
==============================================================================
--- (empty file)
+++ hugo/trunk/src/work/marci/lp/lp_solver_wrapper.h Tue Aug 17 15:20:46 2004
@@ -0,0 +1,382 @@
+// -*- c++ -*-
+#ifndef HUGO_LP_SOLVER_WRAPPER_H
+#define HUGO_LP_SOLVER_WRAPPER
+
+// #include <stdio.h>
+#include <stdlib.h>
+// #include <stdio>
+//#include <stdlib>
+#include "glpk.h"
+
+#include <iostream>
+#include <vector>
+#include <string>
+#include <list>
+#include <memory>
+#include <utility>
+
+//#include <sage_graph.h>
+//#include <hugo/list_graph.h>
+//#include <hugo/graph_wrapper.h>
+#include <hugo/invalid.h>
+//#include <bfs_dfs.h>
+//#include <stp.h>
+//#include <hugo/max_flow.h>
+//#include <augmenting_flow.h>
+//#include <iter_map.h>
+
+using std::cout;
+using std::cin;
+using std::endl;
+
+namespace hugo {
+
+ /// \brief A partitioned vector with iterable classes.
+ ///
+ /// This class implements a container in which the data is stored in an
+ /// stl vector, the range is partitioned into sets and each set is
+ /// doubly linked in a list.
+ /// That is, each class is iterable by hugo iterators, and any member of
+ /// the vector can bo moved to an other class.
+ template <typename T>
+ class IterablePartition {
+ protected:
+ struct Node {
+ T data;
+ int prev; //invalid az -1
+ int next;
+ };
+ std::vector<Node> nodes;
+ struct Tip {
+ int first;
+ int last;
+ };
+ std::vector<Tip> tips;
+ public:
+ /// The classes are indexed by integers from \c 0 to \c classNum()-1.
+ int classNum() const { return tips.size(); }
+ /// This hugo style iterator iterates through a class.
+ class ClassIt;
+ /// Constructor. The number of classes is to be given which is fixed
+ /// over the life of the container.
+ /// The partition classes are indexed from 0 to class_num-1.
+ IterablePartition(int class_num) {
+ for (int i=0; i<class_num; ++i) {
+ Tip t;
+ t.first=t.last=-1;
+ tips.push_back(t);
+ }
+ }
+ protected:
+ void befuz(ClassIt it, int class_id) {
+ if (tips[class_id].first==-1) {
+ if (tips[class_id].last==-1) {
+ nodes[it.i].prev=nodes[it.i].next=-1;
+ tips[class_id].first=tips[class_id].last=it.i;
+ }
+ } else {
+ nodes[it.i].prev=tips[class_id].last;
+ nodes[it.i].next=-1;
+ nodes[tips[class_id].last].next=it.i;
+ tips[class_id].last=it.i;
+ }
+ }
+ void kifuz(ClassIt it, int class_id) {
+ if (tips[class_id].first==it.i) {
+ if (tips[class_id].last==it.i) {
+ tips[class_id].first=tips[class_id].last=-1;
+ } else {
+ tips[class_id].first=nodes[it.i].next;
+ nodes[nodes[it.i].next].prev=-1;
+ }
+ } else {
+ if (tips[class_id].last==it.i) {
+ tips[class_id].last=nodes[it.i].prev;
+ nodes[nodes[it.i].prev].next=-1;
+ } else {
+ nodes[nodes[it.i].next].prev=nodes[it.i].prev;
+ nodes[nodes[it.i].prev].next=nodes[it.i].next;
+ }
+ }
+ }
+ public:
+ /// A new element with data \c t is pushed into the vector and into class
+ /// \c class_id.
+ ClassIt push_back(const T& t, int class_id) {
+ Node n;
+ n.data=t;
+ nodes.push_back(n);
+ int i=nodes.size()-1;
+ befuz(i, class_id);
+ return i;
+ }
+ /// A member is moved to an other class.
+ void set(ClassIt it, int old_class_id, int new_class_id) {
+ kifuz(it.i, old_class_id);
+ befuz(it.i, new_class_id);
+ }
+ /// Returns the data pointed by \c it.
+ T& operator[](ClassIt it) { return nodes[it.i].data; }
+ /// Returns the data pointed by \c it.
+ const T& operator[](ClassIt it) const { return nodes[it.i].data; }
+ class ClassIt {
+ friend class IterablePartition;
+ protected:
+ int i;
+ public:
+ /// Default constructor.
+ ClassIt() { }
+ /// This constructor constructs an iterator which points
+ /// to the member of th container indexed by the integer _i.
+ ClassIt(const int& _i) : i(_i) { }
+ /// Invalid constructor.
+ ClassIt(const Invalid&) : i(-1) { }
+ };
+ /// First member of class \c class_id.
+ ClassIt& first(ClassIt& it, int class_id) const {
+ it.i=tips[class_id].first;
+ return it;
+ }
+ /// Next member.
+ ClassIt& next(ClassIt& it) const {
+ it.i=nodes[it.i].next;
+ return it;
+ }
+ /// True iff the iterator is valid.
+ bool valid(const ClassIt& it) const { return it.i!=-1; }
+ };
+
+ /// \brief Wrappers for LP solvers
+ ///
+ /// This class implements a hugo wrapper for glpk.
+ /// Later other LP-solvers will be wrapped into hugo.
+ /// The aim of this class is to give a general surface to different
+ /// solvers, i.e. it makes possible to write algorithms using LP's,
+ /// in which the solver can be changed to an other one easily.
+ class LPSolverWrapper {
+ public:
+
+// class Row {
+// protected:
+// int i;
+// public:
+// Row() { }
+// Row(const Invalid&) : i(0) { }
+// Row(const int& _i) : i(_i) { }
+// operator int() const { return i; }
+// };
+// class RowIt : public Row {
+// public:
+// RowIt(const Row& row) : Row(row) { }
+// };
+
+// class Col {
+// protected:
+// int i;
+// public:
+// Col() { }
+// Col(const Invalid&) : i(0) { }
+// Col(const int& _i) : i(_i) { }
+// operator int() const { return i; }
+// };
+// class ColIt : public Col {
+// ColIt(const Col& col) : Col(col) { }
+// };
+
+ public:
+ LPX* lp;
+ typedef IterablePartition<int>::ClassIt RowIt;
+ IterablePartition<int> row_iter_map;
+ typedef IterablePartition<int>::ClassIt ColIt;
+ IterablePartition<int> col_iter_map;
+ //std::vector<int> row_id_to_lp_row_id;
+ //std::vector<int> col_id_to_lp_col_id;
+ const int VALID_ID;
+ const int INVALID_ID;
+
+ public:
+ LPSolverWrapper() : lp(lpx_create_prob()),
+ row_iter_map(2),
+ col_iter_map(2),
+ //row_id_to_lp_row_id(), col_id_to_lp_col_id(),
+ VALID_ID(0), INVALID_ID(1) {
+ lpx_set_int_parm(lp, LPX_K_DUAL, 1);
+ }
+ ~LPSolverWrapper() {
+ lpx_delete_prob(lp);
+ }
+ void setMinimize() {
+ lpx_set_obj_dir(lp, LPX_MIN);
+ }
+ void setMaximize() {
+ lpx_set_obj_dir(lp, LPX_MAX);
+ }
+ ColIt addCol() {
+ int i=lpx_add_cols(lp, 1);
+ ColIt col_it;
+ col_iter_map.first(col_it, INVALID_ID);
+ if (col_iter_map.valid(col_it)) { //van hasznalhato hely
+ col_iter_map.set(col_it, INVALID_ID, VALID_ID);
+ col_iter_map[col_it]=i;
+ //col_id_to_lp_col_id[col_iter_map[col_it]]=i;
+ } else { //a cucc vegere kell inzertalni mert nincs szabad hely
+ //col_id_to_lp_col_id.push_back(i);
+ //int j=col_id_to_lp_col_id.size()-1;
+ col_it=col_iter_map.push_back(i, VALID_ID);
+ }
+// edge_index_map.set(e, i);
+// lpx_set_col_bnds(lp, i, LPX_DB, 0.0, 1.0);
+// lpx_set_obj_coef(lp, i, cost[e]);
+ return col_it;
+ }
+ RowIt addRow() {
+ int i=lpx_add_rows(lp, 1);
+ RowIt row_it;
+ row_iter_map.first(row_it, INVALID_ID);
+ if (row_iter_map.valid(row_it)) { //van hasznalhato hely
+ row_iter_map.set(row_it, INVALID_ID, VALID_ID);
+ row_iter_map[row_it]=i;
+ } else { //a cucc vegere kell inzertalni mert nincs szabad hely
+ row_it=row_iter_map.push_back(i, VALID_ID);
+ }
+ return row_it;
+ }
+ //pair<RowIt, double>-bol kell megadni egy std range-et
+ template <typename Begin, typename End>
+ void setColCoeffs(const ColIt& col_it,
+ Begin begin, End end) {
+ int mem_length=1+lpx_get_num_rows(lp);
+ int* indices = new int[mem_length];
+ double* doubles = new double[mem_length];
+ int length=0;
+ for ( ; begin!=end; ++begin) {
+ ++length;
+ indices[length]=row_iter_map[begin->first];
+ doubles[length]=begin->second;
+ }
+ lpx_set_mat_col(lp, col_iter_map[col_it], length, indices, doubles);
+ delete [] indices;
+ delete [] doubles;
+ }
+ //pair<ColIt, double>-bol kell megadni egy std range-et
+ template <typename Begin, typename End>
+ void setRowCoeffs(const RowIt& row_it,
+ Begin begin, End end) {
+ int mem_length=1+lpx_get_num_cols(lp);
+ int* indices = new int[mem_length];
+ double* doubles = new double[mem_length];
+ int length=0;
+ for ( ; begin!=end; ++begin) {
+ ++length;
+ indices[length]=col_iter_map[begin->first];
+ doubles[length]=begin->second;
+ }
+ lpx_set_mat_row(lp, row_iter_map[row_it], length, indices, doubles);
+ delete [] indices;
+ delete [] doubles;
+ }
+ void eraseCol(const ColIt& col_it) {
+ col_iter_map.set(col_it, VALID_ID, INVALID_ID);
+ int cols[2];
+ cols[1]=col_iter_map[col_it];
+ lpx_del_cols(lp, 1, cols);
+ col_iter_map[col_it]=0; //glpk specifikus
+ ColIt it;
+ for (col_iter_map.first(it, VALID_ID);
+ col_iter_map.valid(it); col_iter_map.next(it)) {
+ if (col_iter_map[it]>cols[1]) --col_iter_map[it];
+ }
+ }
+ void eraseRow(const RowIt& row_it) {
+ row_iter_map.set(row_it, VALID_ID, INVALID_ID);
+ int rows[2];
+ rows[1]=row_iter_map[row_it];
+ lpx_del_rows(lp, 1, rows);
+ row_iter_map[row_it]=0; //glpk specifikus
+ RowIt it;
+ for (row_iter_map.first(it, VALID_ID);
+ row_iter_map.valid(it); row_iter_map.next(it)) {
+ if (row_iter_map[it]>rows[1]) --row_iter_map[it];
+ }
+ }
+ void setColBounds(const ColIt& col_it, int bound_type,
+ double lo, double up) {
+ lpx_set_col_bnds(lp, col_iter_map[col_it], bound_type, lo, up);
+ }
+ void setObjCoef(const ColIt& col_it, double obj_coef) {
+ lpx_set_obj_coef(lp, col_iter_map[col_it], obj_coef);
+ }
+ void setRowBounds(const RowIt& row_it, int bound_type,
+ double lo, double up) {
+ lpx_set_row_bnds(lp, row_iter_map[row_it], bound_type, lo, up);
+ }
+// void setObjCoef(const RowIt& row_it, double obj_coef) {
+// lpx_set_obj_coef(lp, row_iter_map[row_it], obj_coef);
+// }
+ void solveSimplex() { lpx_simplex(lp); }
+ void solvePrimalSimplex() { lpx_simplex(lp); }
+ void solveDualSimplex() { lpx_simplex(lp); }
+ double getPrimal(const ColIt& col_it) {
+ return lpx_get_col_prim(lp, col_iter_map[col_it]);
+ }
+ double getObjVal() { return lpx_get_obj_val(lp); }
+ int rowNum() const { return lpx_get_num_rows(lp); }
+ int colNum() const { return lpx_get_num_cols(lp); }
+ int warmUp() { return lpx_warm_up(lp); }
+ void printWarmUpStatus(int i) {
+ switch (i) {
+ case LPX_E_OK: cout << "LPX_E_OK" << endl; break;
+ case LPX_E_EMPTY: cout << "LPX_E_EMPTY" << endl; break;
+ case LPX_E_BADB: cout << "LPX_E_BADB" << endl; break;
+ case LPX_E_SING: cout << "LPX_E_SING" << endl; break;
+ }
+ }
+ int getPrimalStatus() { return lpx_get_prim_stat(lp); }
+ void printPrimalStatus(int i) {
+ switch (i) {
+ case LPX_P_UNDEF: cout << "LPX_P_UNDEF" << endl; break;
+ case LPX_P_FEAS: cout << "LPX_P_FEAS" << endl; break;
+ case LPX_P_INFEAS: cout << "LPX_P_INFEAS" << endl; break;
+ case LPX_P_NOFEAS: cout << "LPX_P_NOFEAS" << endl; break;
+ }
+ }
+ int getDualStatus() { return lpx_get_dual_stat(lp); }
+ void printDualStatus(int i) {
+ switch (i) {
+ case LPX_D_UNDEF: cout << "LPX_D_UNDEF" << endl; break;
+ case LPX_D_FEAS: cout << "LPX_D_FEAS" << endl; break;
+ case LPX_D_INFEAS: cout << "LPX_D_INFEAS" << endl; break;
+ case LPX_D_NOFEAS: cout << "LPX_D_NOFEAS" << endl; break;
+ }
+ }
+ /// Returns the status of the slack variable assigned to row \c row_it.
+ int getRowStat(const RowIt& row_it) {
+ return lpx_get_row_stat(lp, row_iter_map[row_it]);
+ }
+ void printRowStatus(int i) {
+ switch (i) {
+ case LPX_BS: cout << "LPX_BS" << endl; break;
+ case LPX_NL: cout << "LPX_NL" << endl; break;
+ case LPX_NU: cout << "LPX_NU" << endl; break;
+ case LPX_NF: cout << "LPX_NF" << endl; break;
+ case LPX_NS: cout << "LPX_NS" << endl; break;
+ }
+ }
+ /// Returns the status of the variable assigned to column \c col_it.
+ int getColStat(const ColIt& col_it) {
+ return lpx_get_col_stat(lp, col_iter_map[col_it]);
+ }
+ void printColStatus(int i) {
+ switch (i) {
+ case LPX_BS: cout << "LPX_BS" << endl; break;
+ case LPX_NL: cout << "LPX_NL" << endl; break;
+ case LPX_NU: cout << "LPX_NU" << endl; break;
+ case LPX_NF: cout << "LPX_NF" << endl; break;
+ case LPX_NS: cout << "LPX_NS" << endl; break;
+ }
+ }
+ };
+
+} //namespace hugo
+
+#endif //HUGO_LP_SOLVER_WRAPPER_H
Added: hugo/trunk/src/work/marci/lp/makefile
==============================================================================
--- (empty file)
+++ hugo/trunk/src/work/marci/lp/makefile Tue Aug 17 15:20:46 2004
@@ -0,0 +1,73 @@
+#INCLUDEDIRS ?= -I.. -I. -I./{marci,jacint,alpar,klao,akos}
+GLPKROOT = /usr/local/glpk-4.4
+INCLUDEDIRS ?= -I../../{marci,alpar,klao,akos,athos} -I. -I../../.. -I../.. -I.. -I$(GLPKROOT)/include
+#INCLUDEDIRS ?= -I../.. -I../.. -I../../{marci,jacint,alpar,klao,akos} -I/usr/local/glpk-4.4/include
+CXXFLAGS = -g -O2 -W -Wall $(INCLUDEDIRS) -ansi -pedantic
+LDFLAGS = -L$(GLPKROOT)/lib -lglpk
+
+BINARIES = max_flow_by_lp sample sample2 sample11 sample15
+
+#include ../makefile
+
+# Hat, ez elismerem, hogy nagyon ronda, de mukodik minden altalam
+# ismert rendszeren :-) (Misi)
+ifdef GCCVER
+CXX := g++-$(GCCVER)
+else
+CXX := $(shell type -p g++-3.3 || type -p g++-3.2 || type -p g++-3.0 || type -p g++-3 || echo g++)
+endif
+
+ifdef DEBUG
+CXXFLAGS += -DDEBUG
+endif
+
+CC := $(CXX)
+
+
+all: $(BINARIES)
+
+################
+# Minden binarishoz egy sor, felsorolva, hogy mely object file-okbol
+# all elo.
+# Kiveve ha siman file.cc -> file esetrol van szo, amikor is nem kell
+# irni semmit.
+
+#proba: proba.o seged.o
+
+################
+
+
+# .depend dep depend:
+# -$(CXX) $(CXXFLAGS) -M $(BINARIES:=.cc) > .depend
+
+#makefile: .depend
+#sinclude .depend
+#moert nem megy az eredeti /usr/bin/ld-vel?
+
+# %: %.o
+# $(CXX) -o $@ $< $(LDFLAGS)
+
+# %.o: %.cc
+# $(CXX) $(CXXFLAGS) -c $<
+
+%: %.cc
+ $(CXX) $(CXXFLAGS) -o $@ $< $(LDFLAGS)
+
+sample11prof: sample11prof.o
+ $(CXX) -pg -o sample11prof sample11prof.o -L$(GLPKROOT)/lib -lglpk
+sample11prof.o: sample11.cc
+ $(CXX) -pg $(CXXFLAGS) -c -o sample11prof.o sample11.cc
+
+# sample.o: sample.cc
+# $(CXX) $(CXXFLAGS) -c -o sample.o sample.cc
+
+# sample2: sample2.o
+# $(CXX) -o sample2 sample2.o -L/usr/local/glpk-4.4/lib -lglpk
+# sample2.o: sample2.cc
+# $(CXX) $(CXXFLAGS) -c -o sample2.o sample2.cc
+
+
+clean:
+ $(RM) *.o $(BINARIES) .depend
+
+.PHONY: all clean dep depend
Added: hugo/trunk/src/work/marci/lp/max_flow_by_lp.cc
==============================================================================
--- (empty file)
+++ hugo/trunk/src/work/marci/lp/max_flow_by_lp.cc Tue Aug 17 15:20:46 2004
@@ -0,0 +1,233 @@
+// -*- c++ -*-
+#include <iostream>
+#include <fstream>
+
+#include <sage_graph.h>
+#include <hugo/smart_graph.h>
+#include <hugo/dimacs.h>
+#include <hugo/time_measure.h>
+//#include <graph_wrapper.h>
+#include <hugo/max_flow.h>
+#include <augmenting_flow.h>
+//#include <preflow_res.h>
+#include <for_each_macros.h>
+#include <lp_solver_wrapper.h>
+
+using namespace hugo;
+
+// Use a DIMACS max flow file as stdin.
+// max_flow_demo < dimacs_max_flow_file
+
+template<typename Edge, typename EdgeIndexMap>
+class PrimalMap {
+protected:
+ LPSolverWrapper* lp;
+ EdgeIndexMap* edge_index_map;
+public:
+ PrimalMap(LPSolverWrapper& _lp, EdgeIndexMap& _edge_index_map) :
+ lp(&_lp), edge_index_map(&_edge_index_map) { }
+ double operator[](Edge e) const {
+ return lp->getPrimal((*edge_index_map)[e]);
+ }
+};
+
+int main(int, char **) {
+
+ typedef SageGraph MutableGraph;
+ //typedef SmartGraph Graph;
+ typedef SageGraph Graph;
+ typedef Graph::Node Node;
+ typedef Graph::EdgeIt EdgeIt;
+
+ Graph g;
+
+ Node s, t;
+ Graph::EdgeMap<int> cap(g);
+ //readDimacsMaxFlow(std::cin, g, s, t, cap);
+ readDimacs(std::cin, g, cap, s, t);
+ Timer ts;
+ Graph::EdgeMap<int> flow(g); //0 flow
+ MaxFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >
+ max_flow_test(g, s, t, cap, flow);
+ AugmentingFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >
+ augmenting_flow_test(g, s, t, cap, flow);
+
+ Graph::NodeMap<bool> cut(g);
+
+ {
+ std::cout << "preflow ..." << std::endl;
+ ts.reset();
+ max_flow_test.run();
+ std::cout << "elapsed time: " << ts << std::endl;
+ std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
+ max_flow_test.actMinCut(cut);
+
+ FOR_EACH_LOC(Graph::EdgeIt, e, g) {
+ if (cut[g.tail(e)] && !cut[g.head(e)] && !flow[e]==cap[e])
+ std::cout << "Slackness does not hold!" << std::endl;
+ if (!cut[g.tail(e)] && cut[g.head(e)] && flow[e]>0)
+ std::cout << "Slackness does not hold!" << std::endl;
+ }
+ }
+
+// {
+// std::cout << "preflow ..." << std::endl;
+// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
+// ts.reset();
+// max_flow_test.preflow(MaxFlow<Graph, int, Graph::EdgeMap<int>, Graph::EdgeMap<int> >::GEN_FLOW);
+// std::cout << "elapsed time: " << ts << std::endl;
+// std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
+
+// FOR_EACH_LOC(Graph::EdgeIt, e, g) {
+// if (cut[g.tail(e)] && !cut[g.head(e)] && !flow[e]==cap[e])
+// std::cout << "Slackness does not hold!" << std::endl;
+// if (!cut[g.tail(e)] && cut[g.head(e)] && flow[e]>0)
+// std::cout << "Slackness does not hold!" << std::endl;
+// }
+// }
+
+// {
+// std::cout << "wrapped preflow ..." << std::endl;
+// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
+// ts.reset();
+// pre_flow_res.run();
+// std::cout << "elapsed time: " << ts << std::endl;
+// std::cout << "flow value: "<< pre_flow_test.flowValue() << std::endl;
+// }
+
+ {
+ std::cout << "physical blocking flow augmentation ..." << std::endl;
+ FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
+ ts.reset();
+ int i=0;
+ while (augmenting_flow_test.augmentOnBlockingFlow<MutableGraph>()) { ++i; }
+ std::cout << "elapsed time: " << ts << std::endl;
+ std::cout << "number of augmentation phases: " << i << std::endl;
+ std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
+
+ FOR_EACH_LOC(Graph::EdgeIt, e, g) {
+ if (cut[g.tail(e)] && !cut[g.head(e)] && !flow[e]==cap[e])
+ std::cout << "Slackness does not hold!" << std::endl;
+ if (!cut[g.tail(e)] && cut[g.head(e)] && flow[e]>0)
+ std::cout << "Slackness does not hold!" << std::endl;
+ }
+ }
+
+// {
+// std::cout << "faster physical blocking flow augmentation ..." << std::endl;
+// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
+// ts.reset();
+// int i=0;
+// while (max_flow_test.augmentOnBlockingFlow1<MutableGraph>()) { ++i; }
+// std::cout << "elapsed time: " << ts << std::endl;
+// std::cout << "number of augmentation phases: " << i << std::endl;
+// std::cout << "flow value: "<< max_flow_test.flowValue() << std::endl;
+// }
+
+ {
+ std::cout << "on-the-fly blocking flow augmentation ..." << std::endl;
+ FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
+ ts.reset();
+ int i=0;
+ while (augmenting_flow_test.augmentOnBlockingFlow2()) { ++i; }
+ std::cout << "elapsed time: " << ts << std::endl;
+ std::cout << "number of augmentation phases: " << i << std::endl;
+ std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
+
+ FOR_EACH_LOC(Graph::EdgeIt, e, g) {
+ if (cut[g.tail(e)] && !cut[g.head(e)] && !flow[e]==cap[e])
+ std::cout << "Slackness does not hold!" << std::endl;
+ if (!cut[g.tail(e)] && cut[g.head(e)] && flow[e]>0)
+ std::cout << "Slackness does not hold!" << std::endl;
+ }
+ }
+
+// {
+// std::cout << "on-the-fly shortest path augmentation ..." << std::endl;
+// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
+// ts.reset();
+// int i=0;
+// while (augmenting_flow_test.augmentOnShortestPath()) { ++i; }
+// std::cout << "elapsed time: " << ts << std::endl;
+// std::cout << "number of augmentation phases: " << i << std::endl;
+// std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
+
+// FOR_EACH_LOC(Graph::EdgeIt, e, g) {
+// if (cut[g.tail(e)] && !cut[g.head(e)] && !flow[e]==cap[e])
+// std::cout << "Slackness does not hold!" << std::endl;
+// if (!cut[g.tail(e)] && cut[g.head(e)] && flow[e]>0)
+// std::cout << "Slackness does not hold!" << std::endl;
+// }
+// }
+
+// {
+// std::cout << "on-the-fly shortest path augmentation ..." << std::endl;
+// FOR_EACH_LOC(Graph::EdgeIt, e, g) flow.set(e, 0);
+// ts.reset();
+// int i=0;
+// while (augmenting_flow_test.augmentOnShortestPath2()) { ++i; }
+// std::cout << "elapsed time: " << ts << std::endl;
+// std::cout << "number of augmentation phases: " << i << std::endl;
+// std::cout << "flow value: "<< augmenting_flow_test.flowValue() << std::endl;
+
+// FOR_EACH_LOC(Graph::EdgeIt, e, g) {
+// if (cut[g.tail(e)] && !cut[g.head(e)] && !flow[e]==cap[e])
+// std::cout << "Slackness does not hold!" << std::endl;
+// if (!cut[g.tail(e)] && cut[g.head(e)] && flow[e]>0)
+// std::cout << "Slackness does not hold!" << std::endl;
+// }
+// }
+
+ ts.reset();
+ LPSolverWrapper lp;
+ lp.setMaximize();
+ typedef LPSolverWrapper::ColIt ColIt;
+ typedef LPSolverWrapper::RowIt RowIt;
+ typedef Graph::EdgeMap<ColIt> EdgeIndexMap;
+ EdgeIndexMap edge_index_map(g);
+ PrimalMap<Graph::Edge, EdgeIndexMap> lp_flow(lp, edge_index_map);
+ Graph::EdgeIt e;
+ for (g.first(e); g.valid(e); g.next(e)) {
+ ColIt col_it=lp.addCol();
+ edge_index_map.set(e, col_it);
+ lp.setColBounds(col_it, LPX_DB, 0.0, cap[e]);
+ }
+ Graph::NodeIt n;
+ for (g.first(n); g.valid(n); g.next(n)) {
+ if (n!=s) {
+ //hurokelek miatt
+ Graph::EdgeMap<int> coeffs(g, 0);
+ {
+ Graph::InEdgeIt e;
+ for (g.first(e, n); g.valid(e); g.next(e)) coeffs.set(e, coeffs[e]+1);
+ }
+ {
+ Graph::OutEdgeIt e;
+ for (g.first(e, n); g.valid(e); g.next(e)) coeffs.set(e, coeffs[e]-1);
+ }
+ if (n==t) {
+ Graph::EdgeIt e;
+ //std::vector< std::pair<ColIt, double> > row;
+ for (g.first(e); g.valid(e); g.next(e)) {
+ if (coeffs[e]!=0)
+ lp.setObjCoef(edge_index_map[e], coeffs[e]);
+ }
+ } else {
+ RowIt row_it=lp.addRow();
+ Graph::EdgeIt e;
+ std::vector< std::pair<ColIt, double> > row;
+ for (g.first(e); g.valid(e); g.next(e)) {
+ if (coeffs[e]!=0)
+ row.push_back(std::make_pair(edge_index_map[e], coeffs[e]));
+ }
+ lp.setRowCoeffs(row_it, row.begin(), row.end());
+ lp.setRowBounds(row_it, LPX_FX, 0.0, 0.0);
+ }
+ }
+ }
+ lp.solveSimplex();
+ std::cout << "flow value: "<< lp.getObjVal() << std::endl;
+ std::cout << "elapsed time: " << ts << std::endl;
+
+ return 0;
+}
More information about the Lemon-commits
mailing list