LEMON/LEMON-official Changeset - r693:e01957e96c67

Repositories » LEMON » LEMON-official

Summary
Changelog
Files
Switch To
- loading...
Options
- Lightweight changelog
- Search
Pull Requests

Changeset - r693:e01957e96c67

« 691:8d289c
« 692:cb8270

842:dc9316 »
698:3adf5e »
697:a8dfe8 »
694:dcba64 »

r693:e01957e96c67

Wed, 29 Apr 2009 20:22:14

[Not Reviewed]

0 comments (0 inline)

alpar (Alpar Juttner)
alpar@cs.elte.hu

Merge

0 15 0

merge default

0 files changed with 665 insertions and 722 deletions:

CMakeLists.txt

Makefile.am

cmake/FindCOIN.cmake

cmake/FindCPLEX.cmake

cmake/FindGLPK.cmake

doc/groups.dox

lemon/Makefile.am

lemon/circulation.h

lemon/core.h

lemon/network_simplex.h

348

520

lemon/preflow.h

test/lp_test.cc

test/min_cost_flow_test.cc

161

105

test/mip_test.cc

tools/dimacs-solver.cc

↑ Collapse diff ↑

CMakeLists.txt

Show inline comments

@@ -14,25 +14,21 @@
 		INCLUDE(FindDoxygen)
 		INCLUDE(FindGhostscript)
 		FIND_PACKAGE(GLPK 4.33)
 		FIND_PACKAGE(CPLEX)
 		FIND_PACKAGE(COIN)
 		ADD_DEFINITIONS(-DHAVE_CONFIG_H)
 		IF(MSVC)
 		  SET(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} /wd4250 /wd4355 /wd4800 /wd4996")
 		# Suppressed warnings:
 		# C4250: 'class1' : inherits 'class2::member' via dominance
 		# C4355: 'this' : used in base member initializer list
 		# C4800: 'type' : forcing value to bool 'true' or 'false' (performance warning)
 		# C4996: 'function': was declared deprecated
 		ENDIF(MSVC)
 		ADD_DEFINITIONS(-DHAVE_CONFIG_H)
 		INCLUDE(CheckTypeSize)
 		CHECK_TYPE_SIZE("long long" LEMON_LONG_LONG)
 		ENABLE_TESTING()
 		ADD_SUBDIRECTORY(lemon)

Makefile.am

Show inline comments

@@ -8,17 +8,18 @@
 		EXTRA_DIST = \
 			AUTHORS \
 			LICENSE \
 			m4/lx_check_cplex.m4 \
 			m4/lx_check_glpk.m4 \
 			m4/lx_check_soplex.m4 \
 			m4/lx_check_clp.m4 \
 			m4/lx_check_cbc.m4 \
 			m4/lx_check_coin.m4 \
 			CMakeLists.txt \
 			cmake/FindGhostscript.cmake \
 			cmake/FindCPLEX.cmake \
 			cmake/FindGLPK.cmake \
 			cmake/FindCOIN.cmake \
 			cmake/version.cmake.in \
 			cmake/version.cmake \
 			cmake/nsis/lemon.ico \
 			cmake/nsis/uninstall.ico
 		pkgconfigdir = $(libdir)/pkgconfig

cmake/FindCOIN.cmake

Show inline comments

 		SET(COIN_ROOT_DIR "" CACHE PATH "COIN root directory")
 		FIND_PATH(COIN_INCLUDE_DIR coin/CoinUtilsConfig.h
 		  PATHS ${COIN_ROOT_DIR}/include)
 		FIND_LIBRARY(COIN_CBC_LIBRARY libCbc
 		  PATHS ${COIN_ROOT_DIR}/lib)
 		FIND_LIBRARY(COIN_CBC_SOLVER_LIBRARY libCbcSolver
 		  PATHS ${COIN_ROOT_DIR}/lib)
 		FIND_LIBRARY(COIN_CGL_LIBRARY libCgl
 		  PATHS ${COIN_ROOT_DIR}/lib)
 		FIND_LIBRARY(COIN_CLP_LIBRARY libClp
 		  PATHS ${COIN_ROOT_DIR}/lib)
 		FIND_LIBRARY(COIN_COIN_UTILS_LIBRARY libCoinUtils
 		  PATHS ${COIN_ROOT_DIR}/lib)
 		FIND_LIBRARY(COIN_OSI_LIBRARY libOsi
 		  PATHS ${COIN_ROOT_DIR}/lib)
 		FIND_LIBRARY(COIN_OSI_CBC_LIBRARY libOsiCbc
 		  PATHS ${COIN_ROOT_DIR}/lib)
 		FIND_LIBRARY(COIN_OSI_CLP_LIBRARY libOsiClp
 		  PATHS ${COIN_ROOT_DIR}/lib)
 		FIND_LIBRARY(COIN_OSI_VOL_LIBRARY libOsiVol
 		  PATHS ${COIN_ROOT_DIR}/lib)
 		FIND_LIBRARY(COIN_VOL_LIBRARY libVol
 		  PATHS ${COIN_ROOT_DIR}/lib)
 		  HINTS ${COIN_ROOT_DIR}/include
+		)
 		FIND_LIBRARY(COIN_CBC_LIBRARY
 		  NAMES Cbc libCbc
 		  HINTS ${COIN_ROOT_DIR}/lib
+		)
 		FIND_LIBRARY(COIN_CBC_SOLVER_LIBRARY
 		  NAMES CbcSolver libCbcSolver
 		  HINTS ${COIN_ROOT_DIR}/lib
+		)
 		FIND_LIBRARY(COIN_CGL_LIBRARY
 		  NAMES Cgl libCgl
 		  HINTS ${COIN_ROOT_DIR}/lib
+		)
 		FIND_LIBRARY(COIN_CLP_LIBRARY
 		  NAMES Clp libClp
 		  HINTS ${COIN_ROOT_DIR}/lib
+		)
 		FIND_LIBRARY(COIN_COIN_UTILS_LIBRARY
 		  NAMES CoinUtils libCoinUtils
 		  HINTS ${COIN_ROOT_DIR}/lib
+		)
 		FIND_LIBRARY(COIN_OSI_LIBRARY
 		  NAMES Osi libOsi
 		  HINTS ${COIN_ROOT_DIR}/lib
+		)
 		FIND_LIBRARY(COIN_OSI_CBC_LIBRARY
 		  NAMES OsiCbc libOsiCbc
 		  HINTS ${COIN_ROOT_DIR}/lib
+		)
 		FIND_LIBRARY(COIN_OSI_CLP_LIBRARY
 		  NAMES OsiClp libOsiClp
 		  HINTS ${COIN_ROOT_DIR}/lib
+		)
 		FIND_LIBRARY(COIN_OSI_VOL_LIBRARY
 		  NAMES OsiVol libOsiVol
 		  HINTS ${COIN_ROOT_DIR}/lib
+		)
 		FIND_LIBRARY(COIN_VOL_LIBRARY
 		  NAMES Vol libVol
 		  HINTS ${COIN_ROOT_DIR}/lib
+		)
 		INCLUDE(FindPackageHandleStandardArgs)
 		FIND_PACKAGE_HANDLE_STANDARD_ARGS(COIN DEFAULT_MSG
 		  COIN_INCLUDE_DIR
 		  COIN_CBC_LIBRARY
 		  COIN_CBC_SOLVER_LIBRARY

cmake/FindCPLEX.cmake

Show inline comments

 		SET(CPLEX_ROOT_DIR "" CACHE PATH "CPLEX root directory")
 		FIND_PATH(CPLEX_INCLUDE_DIR
 		  ilcplex/cplex.h
 		  PATHS "C:/ILOG/CPLEX91/include")
 		  PATHS "C:/ILOG/CPLEX91/include"
 		  PATHS "/opt/ilog/cplex91/include"
 		  HINTS ${CPLEX_ROOT_DIR}/include
+		)
 		FIND_LIBRARY(CPLEX_LIBRARY
 		  NAMES cplex91
 		  PATHS "C:/ILOG/CPLEX91/lib/msvc7/stat_mda")
 		  cplex91
 		  PATHS "C:/ILOG/CPLEX91/lib/msvc7/stat_mda"
 		  PATHS "/opt/ilog/cplex91/bin"
 		  HINTS ${CPLEX_ROOT_DIR}/bin
+		)
 		INCLUDE(FindPackageHandleStandardArgs)
 		FIND_PACKAGE_HANDLE_STANDARD_ARGS(CPLEX DEFAULT_MSG CPLEX_LIBRARY CPLEX_INCLUDE_DIR)
 		FIND_PATH(CPLEX_BIN_DIR
 		  cplex91.dll
-		  PATHS "C:/ILOG/CPLEX91/bin/x86_win32")
 		  PATHS "C:/ILOG/CPLEX91/bin/x86_win32"
+		)
 		IF(CPLEX_FOUND)
 		  SET(CPLEX_INCLUDE_DIRS ${CPLEX_INCLUDE_DIR})
 		  SET(CPLEX_LIBRARIES ${CPLEX_LIBRARY})
 		  IF(CMAKE_SYSTEM_NAME STREQUAL "Linux")
 		    SET(CPLEX_LIBRARIES "${CPLEX_LIBRARIES};m;pthread")
 		  ENDIF(CMAKE_SYSTEM_NAME STREQUAL "Linux")
 		ENDIF(CPLEX_FOUND)
 		MARK_AS_ADVANCED(CPLEX_LIBRARY CPLEX_INCLUDE_DIR CPLEX_BIN_DIR)
 		IF(CPLEX_FOUND)
 		  SET(LEMON_HAVE_LP TRUE)

cmake/FindGLPK.cmake

Show inline comments

 		SET(GLPK_ROOT_DIR "" CACHE PATH "GLPK root directory")
 		SET(GLPK_REGKEY "[HKEY_LOCAL_MACHINE\\SOFTWARE\\GnuWin32\\Glpk;InstallPath]")
 		GET_FILENAME_COMPONENT(GLPK_ROOT_PATH ${GLPK_REGKEY} ABSOLUTE)
 		FIND_PATH(GLPK_INCLUDE_DIR
 		  glpk.h
-		  PATHS ${GLPK_REGKEY}/include)
 		  PATHS ${GLPK_REGKEY}/include
 		  HINTS ${GLPK_ROOT_DIR}/include
+		)
 		FIND_LIBRARY(GLPK_LIBRARY
 		  glpk
 		  PATHS ${GLPK_REGKEY}/lib
 		  HINTS ${GLPK_ROOT_DIR}/lib
+		)
 		FIND_LIBRARY(GLPK_LIBRARY
 		  NAMES glpk
-		  PATHS ${GLPK_REGKEY}/lib)
+		IF(GLPK_INCLUDE_DIR AND GLPK_LIBRARY)
 		  FILE(READ ${GLPK_INCLUDE_DIR}/glpk.h GLPK_GLPK_H)
 		  STRING(REGEX MATCH "define[ ]+GLP_MAJOR_VERSION[ ]+[0-9]+" GLPK_MAJOR_VERSION_LINE "${GLPK_GLPK_H}")
 		  STRING(REGEX REPLACE "define[ ]+GLP_MAJOR_VERSION[ ]+([0-9]+)" "\\1" GLPK_VERSION_MAJOR "${GLPK_MAJOR_VERSION_LINE}")
 		  STRING(REGEX MATCH "define[ ]+GLP_MINOR_VERSION[ ]+[0-9]+" GLPK_MINOR_VERSION_LINE "${GLPK_GLPK_H}")
 		  STRING(REGEX REPLACE "define[ ]+GLP_MINOR_VERSION[ ]+([0-9]+)" "\\1" GLPK_VERSION_MINOR "${GLPK_MINOR_VERSION_LINE}")
 		  SET(GLPK_VERSION_STRING "${GLPK_VERSION_MAJOR}.${GLPK_VERSION_MINOR}")
 		  IF(GLPK_FIND_VERSION)
 		    IF(GLPK_FIND_VERSION_COUNT GREATER 2)
 		      MESSAGE(SEND_ERROR "unexpected version string")
 		    ENDIF(GLPK_FIND_VERSION_COUNT GREATER 2)
 		    MATH(EXPR GLPK_REQUESTED_VERSION "${GLPK_FIND_VERSION_MAJOR}*100 + ${GLPK_FIND_VERSION_MINOR}")
 		    MATH(EXPR GLPK_FOUND_VERSION "${GLPK_VERSION_MAJOR}*100 + ${GLPK_VERSION_MINOR}")
 		    IF(GLPK_FOUND_VERSION LESS GLPK_REQUESTED_VERSION)
 		      SET(GLPK_PROPER_VERSION_FOUND FALSE)
 		    ELSE(GLPK_FOUND_VERSION LESS GLPK_REQUESTED_VERSION)
 		      SET(GLPK_PROPER_VERSION_FOUND TRUE)
 		    ENDIF(GLPK_FOUND_VERSION LESS GLPK_REQUESTED_VERSION)
 		  ELSE(GLPK_FIND_VERSION)
 		    SET(GLPK_PROPER_VERSION_FOUND TRUE)
 		  ENDIF(GLPK_FIND_VERSION)
 		ENDIF(GLPK_INCLUDE_DIR AND GLPK_LIBRARY)
 		INCLUDE(FindPackageHandleStandardArgs)
 		FIND_PACKAGE_HANDLE_STANDARD_ARGS(GLPK DEFAULT_MSG GLPK_LIBRARY GLPK_INCLUDE_DIR)
+		FIND_PACKAGE_HANDLE_STANDARD_ARGS(GLPK DEFAULT_MSG GLPK_LIBRARY GLPK_INCLUDE_DIR GLPK_PROPER_VERSION_FOUND)
 		IF(GLPK_FOUND)
 		  SET(GLPK_INCLUDE_DIRS ${GLPK_INCLUDE_DIR})
 		  SET(GLPK_LIBRARIES ${GLPK_LIBRARY})
 		  SET(GLPK_BIN_DIR ${GLPK_ROOT_PATH}/bin)
 		ENDIF(GLPK_FOUND)

doc/groups.dox

Show inline comments

@@ -349,23 +349,23 @@
 		circulations.
 		The \e minimum \e cost \e flow \e problem is to find a feasible flow of
 		minimum total cost from a set of supply nodes to a set of demand nodes
 		in a network with capacity constraints (lower and upper bounds)
 		and arc costs.
 		Formally, let \f$G=(V,A)\f$ be a digraph,
 		\f$lower, upper: A\rightarrow\mathbf{Z}^+_0\f$ denote the lower and
 		Formally, let \f$G=(V,A)\f$ be a digraph, \f$lower: A\rightarrow\mathbf{Z}\f$,
 		\f$upper: A\rightarrow\mathbf{Z}\cup\{+\infty\}\f$ denote the lower and
 		upper bounds for the flow values on the arcs, for which
 		\f$0 \leq lower(uv) \leq upper(uv)\f$ holds for all \f$uv\in A\f$.
 		\f$cost: A\rightarrow\mathbf{Z}^+_0\f$ denotes the cost per unit flow
-		on the arcs, and \f$sup: V\rightarrow\mathbf{Z}\f$ denotes the
+		\f$lower(uv) \leq upper(uv)\f$ must hold for all \f$uv\in A\f$,
 		\f$cost: A\rightarrow\mathbf{Z}\f$ denotes the cost per unit flow
 		on the arcs and \f$sup: V\rightarrow\mathbf{Z}\f$ denotes the
 		signed supply values of the nodes.
 		If \f$sup(u)>0\f$, then \f$u\f$ is a supply node with \f$sup(u)\f$
 		supply, if \f$sup(u)<0\f$, then \f$u\f$ is a demand node with
 		\f$-sup(u)\f$ demand.
-		A minimum cost flow is an \f$f: A\rightarrow\mathbf{Z}^+_0\f$ solution
 		A minimum cost flow is an \f$f: A\rightarrow\mathbf{Z}\f$ solution
 		of the following optimization problem.
 		\f[ \min\sum_{uv\in A} f(uv) \cdot cost(uv) \f]
 		\f[ \sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) \geq
 		    sup(u) \quad \forall u\in V \f]
 		\f[ lower(uv) \leq f(uv) \leq upper(uv) \quad \forall uv\in A \f]
@@ -401,30 +401,30 @@
 		However if the sum of the supply values is zero, then these two problems
 		are equivalent. So if you need the equality form, you have to ensure this
 		additional contraint for the algorithms.
 		The dual solution of the minimum cost flow problem is represented by node
 		potentials \f$\pi: V\rightarrow\mathbf{Z}\f$.
-		An \f$f: A\rightarrow\mathbf{Z}^+_0\f$ feasible solution of the problem
 		An \f$f: A\rightarrow\mathbf{Z}\f$ feasible solution of the problem
 		is optimal if and only if for some \f$\pi: V\rightarrow\mathbf{Z}\f$
 		node potentials the following \e complementary \e slackness optimality
 		conditions hold.
 		 - For all \f$uv\in A\f$ arcs:
 		   - if \f$cost^\pi(uv)>0\f$, then \f$f(uv)=lower(uv)\f$;
 		   - if \f$lower(uv)<f(uv)<upper(uv)\f$, then \f$cost^\pi(uv)=0\f$;
 		   - if \f$cost^\pi(uv)<0\f$, then \f$f(uv)=upper(uv)\f$.
 		 - For all \f$u\in V\f$:
+		 - For all \f$u\in V\f$ nodes:
 		   - if \f$\sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) \neq sup(u)\f$,
 		     then \f$\pi(u)=0\f$.
 		Here \f$cost^\pi(uv)\f$ denotes the \e reduced \e cost of the arc
-		\f$uv\in A\f$ with respect to the node potentials \f$\pi\f$, i.e.
+		\f$uv\in A\f$ with respect to the potential function \f$\pi\f$, i.e.
 		\f[ cost^\pi(uv) = cost(uv) + \pi(u) - \pi(v).\f]
 		All algorithms provide dual solution (node potentials) as well
+		All algorithms provide dual solution (node potentials) as well,
 		if an optimal flow is found.
 		LEMON contains several algorithms for solving minimum cost flow problems.
 		 - \ref NetworkSimplex Primal Network Simplex algorithm with various
 		   pivot strategies.
 		 - \ref CostScaling Push-Relabel and Augment-Relabel algorithms based on

lemon/Makefile.am

Show inline comments

@@ -12,13 +12,14 @@
 			lemon/color.cc \
 			lemon/lp_base.cc \
 			lemon/lp_skeleton.cc \
 			lemon/random.cc \
 			lemon/bits/windows.cc
+		nodist_lemon_HEADERS = lemon/config.h
 		lemon_libemon_la_CXXFLAGS = \
 			$(AM_CXXFLAGS) \
 			$(GLPK_CFLAGS) \
 			$(CPLEX_CFLAGS) \
 			$(SOPLEX_CXXFLAGS) \
 			$(CLP_CXXFLAGS) \
@@ -54,17 +55,17 @@
 		lemon_HEADERS += \
 			lemon/adaptors.h \
 			lemon/arg_parser.h \
 			lemon/assert.h \
 			lemon/bfs.h \
 			lemon/bin_heap.h \
 			lemon/cbc.h \
 			lemon/circulation.h \
 			lemon/clp.h \
 			lemon/color.h \
 			lemon/concept_check.h \
 			lemon/config.h \
 			lemon/connectivity.h \
 			lemon/counter.h \
 			lemon/core.h \
 			lemon/cplex.h \
 			lemon/dfs.h \
 			lemon/dijkstra.h \

lemon/circulation.h

Show inline comments

@@ -61,21 +61,21 @@
 		    ///
 		    /// The type of the map that stores the signed supply values of the
 		    /// nodes.
 		    /// It must conform to the \ref concepts::ReadMap "ReadMap" concept.
 		    typedef SM SupplyMap;
 		    /// \brief The type of the flow values.
 		    typedef typename SupplyMap::Value Flow;
 		    /// \brief The type of the flow and supply values.
 		    typedef typename SupplyMap::Value Value;
 		    /// \brief The type of the map that stores the flow values.
 		    ///
 		    /// The type of the map that stores the flow values.
 		    /// It must conform to the \ref concepts::ReadWriteMap "ReadWriteMap"
 		    /// concept.
-		    typedef typename Digraph::template ArcMap<Flow> FlowMap;
+		    typedef typename Digraph::template ArcMap<Value> FlowMap;
 		    /// \brief Instantiates a FlowMap.
 		    ///
 		    /// This function instantiates a \ref FlowMap.
 		    /// \param digraph The digraph for which we would like to define
 		    /// the flow map.
@@ -101,13 +101,13 @@
 		      return new Elevator(digraph, max_level);
+		    }
 		    /// \brief The tolerance used by the algorithm
 		    ///
 		    /// The tolerance used by the algorithm to handle inexact computation.
-		    typedef lemon::Tolerance<Flow> Tolerance;
+		    typedef lemon::Tolerance<Value> Tolerance;
 		  };
 		  /**
 		     \brief Push-relabel algorithm for the network circulation problem.
@@ -184,14 +184,14 @@
 		  public:
 		    ///The \ref CirculationDefaultTraits "traits class" of the algorithm.
 		    typedef TR Traits;
 		    ///The type of the digraph the algorithm runs on.
 		    typedef typename Traits::Digraph Digraph;
 		    ///The type of the flow values.
 		    typedef typename Traits::Flow Flow;
 		    ///The type of the flow and supply values.
 		    typedef typename Traits::Value Value;
 		    ///The type of the lower bound map.
 		    typedef typename Traits::LowerMap LowerMap;
 		    ///The type of the upper bound (capacity) map.
 		    typedef typename Traits::UpperMap UpperMap;
 		    ///The type of the supply map.
@@ -218,13 +218,13 @@
 		    FlowMap *_flow;
 		    bool _local_flow;
 		    Elevator* _level;
 		    bool _local_level;
-		    typedef typename Digraph::template NodeMap<Flow> ExcessMap;
+		    typedef typename Digraph::template NodeMap<Value> ExcessMap;
 		    ExcessMap* _excess;
 		    Tolerance _tol;
 		    int _el;
 		  public:
@@ -527,13 +527,13 @@
 		          (*_excess)[_g.source(e)] -= (*_up)[e];
 		        } else if (_tol.less(-(*_excess)[_g.target(e)], (*_lo)[e])) {
 		          _flow->set(e, (*_lo)[e]);
 		          (*_excess)[_g.target(e)] += (*_lo)[e];
 		          (*_excess)[_g.source(e)] -= (*_lo)[e];
 		        } else {
-		          Flow fc = -(*_excess)[_g.target(e)];
+		          Value fc = -(*_excess)[_g.target(e)];
 		          _flow->set(e, fc);
 		          (*_excess)[_g.target(e)] = 0;
 		          (*_excess)[_g.source(e)] -= fc;
+		        }
+		      }
@@ -560,17 +560,17 @@
 		      Node act;
 		      Node bact=INVALID;
 		      Node last_activated=INVALID;
 		      while((act=_level->highestActive())!=INVALID) {
 		        int actlevel=(*_level)[act];
 		        int mlevel=_node_num;
-		        Flow exc=(*_excess)[act];
+		        Value exc=(*_excess)[act];
 		        for(OutArcIt e(_g,act);e!=INVALID; ++e) {
 		          Node v = _g.target(e);
-		          Flow fc=(*_up)[e]-(*_flow)[e];
+		          Value fc=(*_up)[e]-(*_flow)[e];
 		          if(!_tol.positive(fc)) continue;
 		          if((*_level)[v]<actlevel) {
 		            if(!_tol.less(fc, exc)) {
 		              _flow->set(e, (*_flow)[e] + exc);
 		              (*_excess)[v] += exc;
 		              if(!_level->active(v) && _tol.positive((*_excess)[v]))
@@ -588,13 +588,13 @@
+		            }
+		          }
 		          else if((*_level)[v]<mlevel) mlevel=(*_level)[v];
+		        }
 		        for(InArcIt e(_g,act);e!=INVALID; ++e) {
 		          Node v = _g.source(e);
-		          Flow fc=(*_flow)[e]-(*_lo)[e];
+		          Value fc=(*_flow)[e]-(*_lo)[e];
 		          if(!_tol.positive(fc)) continue;
 		          if((*_level)[v]<actlevel) {
 		            if(!_tol.less(fc, exc)) {
 		              _flow->set(e, (*_flow)[e] - exc);
 		              (*_excess)[v] += exc;
 		              if(!_level->active(v) && _tol.positive((*_excess)[v]))
@@ -658,19 +658,19 @@
 		    /// these functions.\n
 		    /// Either \ref run() or \ref start() should be called before
 		    /// using them.
 		    ///@{
 		    /// \brief Returns the flow on the given arc.
+		    /// \brief Returns the flow value on the given arc.
 		    ///
 		    /// Returns the flow on the given arc.
+		    /// Returns the flow value on the given arc.
 		    ///
 		    /// \pre Either \ref run() or \ref init() must be called before
 		    /// using this function.
-		    Flow flow(const Arc& arc) const {
+		    Value flow(const Arc& arc) const {
 		      return (*_flow)[arc];
+		    }
 		    /// \brief Returns a const reference to the flow map.
 		    ///
 		    /// Returns a const reference to the arc map storing the found flow.
@@ -747,13 +747,13 @@
 		    ///
 		    bool checkFlow() const {
 		      for(ArcIt e(_g);e!=INVALID;++e)
 		        if((*_flow)[e]<(*_lo)[e]||(*_flow)[e]>(*_up)[e]) return false;
 		      for(NodeIt n(_g);n!=INVALID;++n)
+		        {
-		          Flow dif=-(*_supply)[n];
+		          Value dif=-(*_supply)[n];
 		          for(InArcIt e(_g,n);e!=INVALID;++e) dif-=(*_flow)[e];
 		          for(OutArcIt e(_g,n);e!=INVALID;++e) dif+=(*_flow)[e];
 		          if(_tol.negative(dif)) return false;
+		        }
 		      return true;
+		    }
@@ -762,16 +762,16 @@
 		    ///Check whether or not the last execution provides a barrier.
 		    ///\sa barrier()
 		    ///\sa barrierMap()
 		    bool checkBarrier() const
+		    {
 		      Flow delta=0;
 		      Flow inf_cap = std::numeric_limits<Flow>::has_infinity ?
 		        std::numeric_limits<Flow>::infinity() :
 		        std::numeric_limits<Flow>::max();
 		      Value delta=0;
 		      Value inf_cap = std::numeric_limits<Value>::has_infinity ?
 		        std::numeric_limits<Value>::infinity() :
 		        std::numeric_limits<Value>::max();
 		      for(NodeIt n(_g);n!=INVALID;++n)
 		        if(barrier(n))
 		          delta-=(*_supply)[n];
 		      for(ArcIt e(_g);e!=INVALID;++e)
+		        {
 		          Node s=_g.source(e);

lemon/core.h

Show inline comments

@@ -19,13 +19,13 @@
 		#ifndef LEMON_CORE_H
 		#define LEMON_CORE_H
 		#include <vector>
 		#include <algorithm>
-		#include <lemon/core.h>
+		#include <lemon/config.h>
 		#include <lemon/bits/enable_if.h>
 		#include <lemon/bits/traits.h>
 		#include <lemon/assert.h>
 		///\file
 		///\brief LEMON core utilities.

lemon/network_simplex.h

Show inline comments

@@ -27,15 +27,12 @@
 		#include <vector>
 		#include <limits>
 		#include <algorithm>
 		#include <lemon/core.h>
 		#include <lemon/math.h>
 		#include <lemon/maps.h>
 		#include <lemon/circulation.h>
 		#include <lemon/adaptors.h>
 		namespace lemon {
 		  /// \addtogroup min_cost_flow
 		  /// @{
@@ -47,55 +44,118 @@
 		  /// This algorithm is a specialized version of the linear programming
 		  /// simplex method directly for the minimum cost flow problem.
 		  /// It is one of the most efficient solution methods.
 		  ///
 		  /// In general this class is the fastest implementation available
 		  /// in LEMON for the minimum cost flow problem.
 		  /// Moreover it supports both direction of the supply/demand inequality
 		  /// constraints. For more information see \ref ProblemType.
 		  /// Moreover it supports both directions of the supply/demand inequality
 		  /// constraints. For more information see \ref SupplyType.
 		  ///
 		  /// Most of the parameters of the problem (except for the digraph)
 		  /// can be given using separate functions, and the algorithm can be
 		  /// executed using the \ref run() function. If some parameters are not
 		  /// specified, then default values will be used.
 		  ///
 		  /// \tparam GR The digraph type the algorithm runs on.
-		  /// \tparam F The value type used for flow amounts, capacity bounds
+		  /// \tparam V The value type used for flow amounts, capacity bounds
 		  /// and supply values in the algorithm. By default it is \c int.
 		  /// \tparam C The value type used for costs and potentials in the
-		  /// algorithm. By default it is the same as \c F.
+		  /// algorithm. By default it is the same as \c V.
 		  ///
 		  /// \warning Both value types must be signed and all input data must
 		  /// be integer.
 		  ///
 		  /// \note %NetworkSimplex provides five different pivot rule
 		  /// implementations, from which the most efficient one is used
 		  /// by default. For more information see \ref PivotRule.
-		  template <typename GR, typename F = int, typename C = F>
+		  template <typename GR, typename V = int, typename C = V>
 		  class NetworkSimplex
+		  {
 		  public:
 		    /// The flow type of the algorithm
 		    typedef F Flow;
-		    /// The cost type of the algorithm
+		    /// The type of the flow amounts, capacity bounds and supply values
 		    typedef V Value;
 		    /// The type of the arc costs
 		    typedef C Cost;
 		#ifdef DOXYGEN
 		    /// The type of the flow map
 		    typedef GR::ArcMap<Flow> FlowMap;
 		    /// The type of the potential map
 		    typedef GR::NodeMap<Cost> PotentialMap;
 		#else
 		    /// The type of the flow map
 		    typedef typename GR::template ArcMap<Flow> FlowMap;
 		    /// The type of the potential map
 		    typedef typename GR::template NodeMap<Cost> PotentialMap;
 		#endif
 		  public:
-		    /// \brief Enum type for selecting the pivot rule.
+		    /// \brief Problem type constants for the \c run() function.
 		    ///
-		    /// Enum type for selecting the pivot rule for the \ref run()
+		    /// Enum type containing the problem type constants that can be
 		    /// returned by the \ref run() function of the algorithm.
 		    enum ProblemType {
 		      /// The problem has no feasible solution (flow).
 		      INFEASIBLE,
 		      /// The problem has optimal solution (i.e. it is feasible and
 		      /// bounded), and the algorithm has found optimal flow and node
 		      /// potentials (primal and dual solutions).
 		      OPTIMAL,
 		      /// The objective function of the problem is unbounded, i.e.
 		      /// there is a directed cycle having negative total cost and
 		      /// infinite upper bound.
 		      UNBOUNDED
 		    };
 		    /// \brief Constants for selecting the type of the supply constraints.
 		    ///
 		    /// Enum type containing constants for selecting the supply type,
 		    /// i.e. the direction of the inequalities in the supply/demand
 		    /// constraints of the \ref min_cost_flow "minimum cost flow problem".
 		    ///
 		    /// The default supply type is \c GEQ, since this form is supported
 		    /// by other minimum cost flow algorithms and the \ref Circulation
 		    /// algorithm, as well.
 		    /// The \c LEQ problem type can be selected using the \ref supplyType()
 		    /// function.
 		    ///
 		    /// Note that the equality form is a special case of both supply types.
 		    enum SupplyType {
 		      /// This option means that there are <em>"greater or equal"</em>
 		      /// supply/demand constraints in the definition, i.e. the exact
 		      /// formulation of the problem is the following.
 		      /**
 		          \f[ \min\sum_{uv\in A} f(uv) \cdot cost(uv) \f]
 		          \f[ \sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) \geq
 		              sup(u) \quad \forall u\in V \f]
 		          \f[ lower(uv) \leq f(uv) \leq upper(uv) \quad \forall uv\in A \f]
 		      */
 		      /// It means that the total demand must be greater or equal to the
 		      /// total supply (i.e. \f$\sum_{u\in V} sup(u)\f$ must be zero or
 		      /// negative) and all the supplies have to be carried out from
 		      /// the supply nodes, but there could be demands that are not
 		      /// satisfied.
 		      GEQ,
 		      /// It is just an alias for the \c GEQ option.
 		      CARRY_SUPPLIES = GEQ,
 		      /// This option means that there are <em>"less or equal"</em>
 		      /// supply/demand constraints in the definition, i.e. the exact
 		      /// formulation of the problem is the following.
 		      /**
 		          \f[ \min\sum_{uv\in A} f(uv) \cdot cost(uv) \f]
 		          \f[ \sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) \leq
 		              sup(u) \quad \forall u\in V \f]
 		          \f[ lower(uv) \leq f(uv) \leq upper(uv) \quad \forall uv\in A \f]
 		      */
 		      /// It means that the total demand must be less or equal to the
 		      /// total supply (i.e. \f$\sum_{u\in V} sup(u)\f$ must be zero or
 		      /// positive) and all the demands have to be satisfied, but there
 		      /// could be supplies that are not carried out from the supply
 		      /// nodes.
 		      LEQ,
 		      /// It is just an alias for the \c LEQ option.
 		      SATISFY_DEMANDS = LEQ
 		    };
 		    /// \brief Constants for selecting the pivot rule.
 		    ///
 		    /// Enum type containing constants for selecting the pivot rule for
 		    /// the \ref run() function.
 		    ///
 		    /// \ref NetworkSimplex provides five different pivot rule
 		    /// implementations that significantly affect the running time
 		    /// of the algorithm.
 		    /// By default \ref BLOCK_SEARCH "Block Search" is used, which
 		    /// proved to be the most efficient and the most robust on various
 		    /// test inputs according to our benchmark tests.
@@ -128,77 +188,21 @@
 		      /// It is a modified version of the Candidate List method.
 		      /// It keeps only the several best eligible arcs from the former
 		      /// candidate list and extends this list in every iteration.
 		      ALTERING_LIST
 		    };
 		    /// \brief Enum type for selecting the problem type.
 		    ///
 		    /// Enum type for selecting the problem type, i.e. the direction of
 		    /// the inequalities in the supply/demand constraints of the
 		    /// \ref min_cost_flow "minimum cost flow problem".
 		    ///
 		    /// The default problem type is \c GEQ, since this form is supported
 		    /// by other minimum cost flow algorithms and the \ref Circulation
 		    /// algorithm as well.
 		    /// The \c LEQ problem type can be selected using the \ref problemType()
 		    /// function.
 		    ///
 		    /// Note that the equality form is a special case of both problem type.
 		    enum ProblemType {
 		      /// This option means that there are "<em>greater or equal</em>"
 		      /// constraints in the defintion, i.e. the exact formulation of the
 		      /// problem is the following.
 		      /**
 		          \f[ \min\sum_{uv\in A} f(uv) \cdot cost(uv) \f]
 		          \f[ \sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) \geq
 		              sup(u) \quad \forall u\in V \f]
 		          \f[ lower(uv) \leq f(uv) \leq upper(uv) \quad \forall uv\in A \f]
 		      */
 		      /// It means that the total demand must be greater or equal to the
 		      /// total supply (i.e. \f$\sum_{u\in V} sup(u)\f$ must be zero or
 		      /// negative) and all the supplies have to be carried out from
 		      /// the supply nodes, but there could be demands that are not
 		      /// satisfied.
 		      GEQ,
 		      /// It is just an alias for the \c GEQ option.
 		      CARRY_SUPPLIES = GEQ,
 		      /// This option means that there are "<em>less or equal</em>"
 		      /// constraints in the defintion, i.e. the exact formulation of the
 		      /// problem is the following.
 		      /**
 		          \f[ \min\sum_{uv\in A} f(uv) \cdot cost(uv) \f]
 		          \f[ \sum_{uv\in A} f(uv) - \sum_{vu\in A} f(vu) \leq
 		              sup(u) \quad \forall u\in V \f]
 		          \f[ lower(uv) \leq f(uv) \leq upper(uv) \quad \forall uv\in A \f]
 		      */
 		      /// It means that the total demand must be less or equal to the
 		      /// total supply (i.e. \f$\sum_{u\in V} sup(u)\f$ must be zero or
 		      /// positive) and all the demands have to be satisfied, but there
 		      /// could be supplies that are not carried out from the supply
 		      /// nodes.
 		      LEQ,
 		      /// It is just an alias for the \c LEQ option.
 		      SATISFY_DEMANDS = LEQ
 		    };
 		  private:
 		    TEMPLATE_DIGRAPH_TYPEDEFS(GR);
 		    typedef typename GR::template ArcMap<Flow> FlowArcMap;
 		    typedef typename GR::template ArcMap<Cost> CostArcMap;
 		    typedef typename GR::template NodeMap<Flow> FlowNodeMap;
 		    typedef std::vector<Arc> ArcVector;
 		    typedef std::vector<Node> NodeVector;
 		    typedef std::vector<int> IntVector;
 		    typedef std::vector<bool> BoolVector;
-		    typedef std::vector<Flow> FlowVector;
+		    typedef std::vector<Value> ValueVector;
 		    typedef std::vector<Cost> CostVector;
 		    // State constants for arcs
 		    enum ArcStateEnum {
 		      STATE_UPPER = -1,
 		      STATE_TREE  =  0,
@@ -210,38 +214,29 @@
 		    // Data related to the underlying digraph
 		    const GR &_graph;
 		    int _node_num;
 		    int _arc_num;
 		    // Parameters of the problem
 		    FlowArcMap *_plower;
 		    FlowArcMap *_pupper;
 		    CostArcMap *_pcost;
 		    FlowNodeMap *_psupply;
 		    bool _pstsup;
 		    Node _psource, _ptarget;
 		    Flow _pstflow;
 		    ProblemType _ptype;
 		    // Result maps
 		    FlowMap *_flow_map;
 		    PotentialMap *_potential_map;
 		    bool _local_flow;
 		    bool _local_potential;
 		    bool _have_lower;
 		    SupplyType _stype;
 		    Value _sum_supply;
 		    // Data structures for storing the digraph
 		    IntNodeMap _node_id;
-		    ArcVector _arc_ref;
+		    IntArcMap _arc_id;
 		    IntVector _source;
 		    IntVector _target;
 		    // Node and arc data
-		    FlowVector _cap;
+		    ValueVector _lower;
 		    ValueVector _upper;
 		    ValueVector _cap;
 		    CostVector _cost;
 		    FlowVector _supply;
 		    FlowVector _flow;
 		    ValueVector _supply;
 		    ValueVector _flow;
 		    CostVector _pi;
 		    // Data for storing the spanning tree structure
 		    IntVector _parent;
 		    IntVector _pred;
 		    IntVector _thread;
@@ -254,13 +249,22 @@
 		    int _root;
 		    // Temporary data used in the current pivot iteration
 		    int in_arc, join, u_in, v_in, u_out, v_out;
 		    int first, second, right, last;
 		    int stem, par_stem, new_stem;
-		    Flow delta;
+		    Value delta;
 		  public:
 		    /// \brief Constant for infinite upper bounds (capacities).
 		    ///
 		    /// Constant for infinite upper bounds (capacities).
 		    /// It is \c std::numeric_limits<Value>::infinity() if available,
 		    /// \c std::numeric_limits<Value>::max() otherwise.
 		    const Value INF;
 		  private:
 		    // Implementation of the First Eligible pivot rule
 		    class FirstEligiblePivotRule
+		    {
@@ -656,323 +660,290 @@
 		    /// \brief Constructor.
 		    ///
 		    /// The constructor of the class.
 		    ///
 		    /// \param graph The digraph the algorithm runs on.
 		    NetworkSimplex(const GR& graph) :
 		      _graph(graph),
 		      _plower(NULL), _pupper(NULL), _pcost(NULL),
 		      _psupply(NULL), _pstsup(false), _ptype(GEQ),
 		      _flow_map(NULL), _potential_map(NULL),
 		      _local_flow(false), _local_potential(false),
 		      _node_id(graph)
 		      _graph(graph), _node_id(graph), _arc_id(graph),
 		      INF(std::numeric_limits<Value>::has_infinity ?
 		          std::numeric_limits<Value>::infinity() :
 		          std::numeric_limits<Value>::max())
+		    {
 		      LEMON_ASSERT(std::numeric_limits<Flow>::is_integer &&
 		                   std::numeric_limits<Flow>::is_signed,
 		        "The flow type of NetworkSimplex must be signed integer");
 		      LEMON_ASSERT(std::numeric_limits<Cost>::is_integer &&
 		                   std::numeric_limits<Cost>::is_signed,
 		        "The cost type of NetworkSimplex must be signed integer");
+		    }
+		      // Check the value types
 		      LEMON_ASSERT(std::numeric_limits<Value>::is_signed,
 		        "The flow type of NetworkSimplex must be signed");
 		      LEMON_ASSERT(std::numeric_limits<Cost>::is_signed,
 		        "The cost type of NetworkSimplex must be signed");
 		      // Resize vectors
 		      _node_num = countNodes(_graph);
 		      _arc_num = countArcs(_graph);
 		      int all_node_num = _node_num + 1;
 		      int all_arc_num = _arc_num + _node_num;
 		    /// Destructor.
 		    ~NetworkSimplex() {
 		      if (_local_flow) delete _flow_map;
 		      if (_local_potential) delete _potential_map;
 		      _source.resize(all_arc_num);
 		      _target.resize(all_arc_num);
 		      _lower.resize(all_arc_num);
 		      _upper.resize(all_arc_num);
 		      _cap.resize(all_arc_num);
 		      _cost.resize(all_arc_num);
 		      _supply.resize(all_node_num);
 		      _flow.resize(all_arc_num);
 		      _pi.resize(all_node_num);
 		      _parent.resize(all_node_num);
 		      _pred.resize(all_node_num);
 		      _forward.resize(all_node_num);
 		      _thread.resize(all_node_num);
 		      _rev_thread.resize(all_node_num);
 		      _succ_num.resize(all_node_num);
 		      _last_succ.resize(all_node_num);
 		      _state.resize(all_arc_num);
 		      // Copy the graph (store the arcs in a mixed order)
 		      int i = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n, ++i) {
 		        _node_id[n] = i;
+		      }
 		      int k = std::max(int(std::sqrt(double(_arc_num))), 10);
 		      i = 0;
 		      for (ArcIt a(_graph); a != INVALID; ++a) {
 		        _arc_id[a] = i;
 		        _source[i] = _node_id[_graph.source(a)];
 		        _target[i] = _node_id[_graph.target(a)];
 		        if ((i += k) >= _arc_num) i = (i % k) + 1;
+		      }
 		      // Initialize maps
 		      for (int i = 0; i != _node_num; ++i) {
 		        _supply[i] = 0;
+		      }
 		      for (int i = 0; i != _arc_num; ++i) {
 		        _lower[i] = 0;
 		        _upper[i] = INF;
 		        _cost[i] = 1;
+		      }
 		      _have_lower = false;
 		      _stype = GEQ;
+		    }
 		    /// \name Parameters
 		    /// The parameters of the algorithm can be specified using these
 		    /// functions.
 		    /// @{
 		    /// \brief Set the lower bounds on the arcs.
 		    ///
 		    /// This function sets the lower bounds on the arcs.
 		    /// If neither this function nor \ref boundMaps() is used before
 		    /// calling \ref run(), the lower bounds will be set to zero
-		    /// on all arcs.
+		    /// If it is not used before calling \ref run(), the lower bounds
 		    /// will be set to zero on all arcs.
 		    ///
 		    /// \param map An arc map storing the lower bounds.
-		    /// Its \c Value type must be convertible to the \c Flow type
+		    /// Its \c Value type must be convertible to the \c Value type
 		    /// of the algorithm.
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    template <typename LOWER>
 		    NetworkSimplex& lowerMap(const LOWER& map) {
 		      delete _plower;
 		      _plower = new FlowArcMap(_graph);
 		    template <typename LowerMap>
 		    NetworkSimplex& lowerMap(const LowerMap& map) {
 		      _have_lower = true;
 		      for (ArcIt a(_graph); a != INVALID; ++a) {
-		        (*_plower)[a] = map[a];
+		        _lower[_arc_id[a]] = map[a];
+		      }
 		      return *this;
+		    }
 		    /// \brief Set the upper bounds (capacities) on the arcs.
 		    ///
 		    /// This function sets the upper bounds (capacities) on the arcs.
 		    /// If none of the functions \ref upperMap(), \ref capacityMap()
 		    /// and \ref boundMaps() is used before calling \ref run(),
 		    /// the upper bounds (capacities) will be set to
 		    /// \c std::numeric_limits<Flow>::max() on all arcs.
 		    /// If it is not used before calling \ref run(), the upper bounds
 		    /// will be set to \ref INF on all arcs (i.e. the flow value will be
 		    /// unbounded from above on each arc).
 		    ///
 		    /// \param map An arc map storing the upper bounds.
-		    /// Its \c Value type must be convertible to the \c Flow type
+		    /// Its \c Value type must be convertible to the \c Value type
 		    /// of the algorithm.
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    template<typename UPPER>
 		    NetworkSimplex& upperMap(const UPPER& map) {
 		      delete _pupper;
 		      _pupper = new FlowArcMap(_graph);
 		    template<typename UpperMap>
 		    NetworkSimplex& upperMap(const UpperMap& map) {
 		      for (ArcIt a(_graph); a != INVALID; ++a) {
-		        (*_pupper)[a] = map[a];
+		        _upper[_arc_id[a]] = map[a];
+		      }
 		      return *this;
+		    }
 		    /// \brief Set the upper bounds (capacities) on the arcs.
 		    ///
 		    /// This function sets the upper bounds (capacities) on the arcs.
 		    /// It is just an alias for \ref upperMap().
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    template<typename CAP>
 		    NetworkSimplex& capacityMap(const CAP& map) {
 		      return upperMap(map);
+		    }
 		    /// \brief Set the lower and upper bounds on the arcs.
 		    ///
 		    /// This function sets the lower and upper bounds on the arcs.
 		    /// If neither this function nor \ref lowerMap() is used before
 		    /// calling \ref run(), the lower bounds will be set to zero
 		    /// on all arcs.
 		    /// If none of the functions \ref upperMap(), \ref capacityMap()
 		    /// and \ref boundMaps() is used before calling \ref run(),
 		    /// the upper bounds (capacities) will be set to
 		    /// \c std::numeric_limits<Flow>::max() on all arcs.
 		    ///
 		    /// \param lower An arc map storing the lower bounds.
 		    /// \param upper An arc map storing the upper bounds.
 		    ///
 		    /// The \c Value type of the maps must be convertible to the
 		    /// \c Flow type of the algorithm.
 		    ///
 		    /// \note This function is just a shortcut of calling \ref lowerMap()
 		    /// and \ref upperMap() separately.
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    template <typename LOWER, typename UPPER>
 		    NetworkSimplex& boundMaps(const LOWER& lower, const UPPER& upper) {
 		      return lowerMap(lower).upperMap(upper);
+		    }
 		    /// \brief Set the costs of the arcs.
 		    ///
 		    /// This function sets the costs of the arcs.
 		    /// If it is not used before calling \ref run(), the costs
 		    /// will be set to \c 1 on all arcs.
 		    ///
 		    /// \param map An arc map storing the costs.
 		    /// Its \c Value type must be convertible to the \c Cost type
 		    /// of the algorithm.
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    template<typename COST>
 		    NetworkSimplex& costMap(const COST& map) {
 		      delete _pcost;
 		      _pcost = new CostArcMap(_graph);
 		    template<typename CostMap>
 		    NetworkSimplex& costMap(const CostMap& map) {
 		      for (ArcIt a(_graph); a != INVALID; ++a) {
-		        (*_pcost)[a] = map[a];
+		        _cost[_arc_id[a]] = map[a];
+		      }
 		      return *this;
+		    }
 		    /// \brief Set the supply values of the nodes.
 		    ///
 		    /// This function sets the supply values of the nodes.
 		    /// If neither this function nor \ref stSupply() is used before
 		    /// calling \ref run(), the supply of each node will be set to zero.
 		    /// (It makes sense only if non-zero lower bounds are given.)
 		    ///
 		    /// \param map A node map storing the supply values.
-		    /// Its \c Value type must be convertible to the \c Flow type
+		    /// Its \c Value type must be convertible to the \c Value type
 		    /// of the algorithm.
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    template<typename SUP>
 		    NetworkSimplex& supplyMap(const SUP& map) {
 		      delete _psupply;
 		      _pstsup = false;
-		      _psupply = new FlowNodeMap(_graph);
+		    template<typename SupplyMap>
 		    NetworkSimplex& supplyMap(const SupplyMap& map) {
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
-		        (*_psupply)[n] = map[n];
+		        _supply[_node_id[n]] = map[n];
+		      }
 		      return *this;
+		    }
 		    /// \brief Set single source and target nodes and a supply value.
 		    ///
 		    /// This function sets a single source node and a single target node
 		    /// and the required flow value.
 		    /// If neither this function nor \ref supplyMap() is used before
 		    /// calling \ref run(), the supply of each node will be set to zero.
 		    /// (It makes sense only if non-zero lower bounds are given.)
 		    ///
 		    /// Using this function has the same effect as using \ref supplyMap()
 		    /// with such a map in which \c k is assigned to \c s, \c -k is
 		    /// assigned to \c t and all other nodes have zero supply value.
 		    ///
 		    /// \param s The source node.
 		    /// \param t The target node.
 		    /// \param k The required amount of flow from node \c s to node \c t
 		    /// (i.e. the supply of \c s and the demand of \c t).
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    NetworkSimplex& stSupply(const Node& s, const Node& t, Flow k) {
 		      delete _psupply;
 		      _psupply = NULL;
 		      _pstsup = true;
 		      _psource = s;
 		      _ptarget = t;
-		      _pstflow = k;
+		    NetworkSimplex& stSupply(const Node& s, const Node& t, Value k) {
 		      for (int i = 0; i != _node_num; ++i) {
 		        _supply[i] = 0;
+		      }
 		      _supply[_node_id[s]] =  k;
 		      _supply[_node_id[t]] = -k;
 		      return *this;
+		    }
-		    /// \brief Set the problem type.
+		    /// \brief Set the type of the supply constraints.
 		    ///
 		    /// This function sets the problem type for the algorithm.
 		    /// If it is not used before calling \ref run(), the \ref GEQ problem
 		    /// This function sets the type of the supply/demand constraints.
 		    /// If it is not used before calling \ref run(), the \ref GEQ supply
 		    /// type will be used.
 		    ///
-		    /// For more information see \ref ProblemType.
+		    /// For more information see \ref SupplyType.
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    NetworkSimplex& problemType(ProblemType problem_type) {
 		      _ptype = problem_type;
 		    NetworkSimplex& supplyType(SupplyType supply_type) {
 		      _stype = supply_type;
 		      return *this;
+		    }
 		    /// \brief Set the flow map.
 		    ///
 		    /// This function sets the flow map.
 		    /// If it is not used before calling \ref run(), an instance will
 		    /// be allocated automatically. The destructor deallocates this
 		    /// automatically allocated map, of course.
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    NetworkSimplex& flowMap(FlowMap& map) {
 		      if (_local_flow) {
 		        delete _flow_map;
 		        _local_flow = false;
+		      }
 		      _flow_map = &map;
 		      return *this;
+		    }
 		    /// \brief Set the potential map.
 		    ///
 		    /// This function sets the potential map, which is used for storing
 		    /// the dual solution.
 		    /// If it is not used before calling \ref run(), an instance will
 		    /// be allocated automatically. The destructor deallocates this
 		    /// automatically allocated map, of course.
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    NetworkSimplex& potentialMap(PotentialMap& map) {
 		      if (_local_potential) {
 		        delete _potential_map;
 		        _local_potential = false;
+		      }
 		      _potential_map = &map;
 		      return *this;
+		    }
 		    /// @}
 		    /// \name Execution Control
 		    /// The algorithm can be executed using \ref run().
 		    /// @{
 		    /// \brief Run the algorithm.
 		    ///
 		    /// This function runs the algorithm.
 		    /// The paramters can be specified using functions \ref lowerMap(),
 		    /// \ref upperMap(), \ref capacityMap(), \ref boundMaps(),
 		    /// \ref costMap(), \ref supplyMap(), \ref stSupply(),
-		    /// \ref problemType(), \ref flowMap() and \ref potentialMap().
+		    /// \ref upperMap(), \ref costMap(), \ref supplyMap(), \ref stSupply(),
 		    /// \ref supplyType().
 		    /// For example,
 		    /// \code
 		    ///   NetworkSimplex<ListDigraph> ns(graph);
-		    ///   ns.boundMaps(lower, upper).costMap(cost)
+		    ///   ns.lowerMap(lower).upperMap(upper).costMap(cost)
 		    ///     .supplyMap(sup).run();
 		    /// \endcode
 		    ///
 		    /// This function can be called more than once. All the parameters
 		    /// that have been given are kept for the next call, unless
 		    /// \ref reset() is called, thus only the modified parameters
 		    /// have to be set again. See \ref reset() for examples.
 		    /// However the underlying digraph must not be modified after this
 		    /// class have been constructed, since it copies and extends the graph.
 		    ///
 		    /// \param pivot_rule The pivot rule that will be used during the
 		    /// algorithm. For more information see \ref PivotRule.
 		    ///
 		    /// \return \c true if a feasible flow can be found.
 		    bool run(PivotRule pivot_rule = BLOCK_SEARCH) {
-		      return init() && start(pivot_rule);
+		    /// \return \c INFEASIBLE if no feasible flow exists,
 		    /// \n \c OPTIMAL if the problem has optimal solution
 		    /// (i.e. it is feasible and bounded), and the algorithm has found
 		    /// optimal flow and node potentials (primal and dual solutions),
 		    /// \n \c UNBOUNDED if the objective function of the problem is
 		    /// unbounded, i.e. there is a directed cycle having negative total
 		    /// cost and infinite upper bound.
 		    ///
 		    /// \see ProblemType, PivotRule
 		    ProblemType run(PivotRule pivot_rule = BLOCK_SEARCH) {
 		      if (!init()) return INFEASIBLE;
 		      return start(pivot_rule);
+		    }
 		    /// \brief Reset all the parameters that have been given before.
 		    ///
 		    /// This function resets all the paramaters that have been given
 		    /// before using functions \ref lowerMap(), \ref upperMap(),
 		    /// \ref capacityMap(), \ref boundMaps(), \ref costMap(),
 		    /// \ref supplyMap(), \ref stSupply(), \ref problemType(),
-		    /// \ref flowMap() and \ref potentialMap().
+		    /// \ref costMap(), \ref supplyMap(), \ref stSupply(), \ref supplyType().
 		    ///
 		    /// It is useful for multiple run() calls. If this function is not
 		    /// used, all the parameters given before are kept for the next
 		    /// \ref run() call.
 		    /// However the underlying digraph must not be modified after this
 		    /// class have been constructed, since it copies and extends the graph.
 		    ///
 		    /// For example,
 		    /// \code
 		    ///   NetworkSimplex<ListDigraph> ns(graph);
 		    ///
 		    ///   // First run
-		    ///   ns.lowerMap(lower).capacityMap(cap).costMap(cost)
+		    ///   ns.lowerMap(lower).upperMap(upper).costMap(cost)
 		    ///     .supplyMap(sup).run();
 		    ///
 		    ///   // Run again with modified cost map (reset() is not called,
 		    ///   // so only the cost map have to be set again)
 		    ///   cost[e] += 100;
 		    ///   ns.costMap(cost).run();
 		    ///
 		    ///   // Run again from scratch using reset()
 		    ///   // (the lower bounds will be set to zero on all arcs)
 		    ///   ns.reset();
-		    ///   ns.capacityMap(cap).costMap(cost)
+		    ///   ns.upperMap(capacity).costMap(cost)
 		    ///     .supplyMap(sup).run();
 		    /// \endcode
 		    ///
 		    /// \return <tt>(*this)</tt>
 		    NetworkSimplex& reset() {
 		      delete _plower;
 		      delete _pupper;
 		      delete _pcost;
 		      delete _psupply;
 		      _plower = NULL;
 		      _pupper = NULL;
 		      _pcost = NULL;
 		      _psupply = NULL;
 		      _pstsup = false;
 		      _ptype = GEQ;
 		      if (_local_flow) delete _flow_map;
 		      if (_local_potential) delete _potential_map;
 		      _flow_map = NULL;
 		      _potential_map = NULL;
 		      _local_flow = false;
 		      _local_potential = false;
+		      for (int i = 0; i != _node_num; ++i) {
 		        _supply[i] = 0;
+		      }
 		      for (int i = 0; i != _arc_num; ++i) {
 		        _lower[i] = 0;
 		        _upper[i] = INF;
 		        _cost[i] = 1;
+		      }
 		      _have_lower = false;
 		      _stype = GEQ;
 		      return *this;
+		    }
 		    /// @}
 		    /// \name Query Functions
@@ -982,33 +953,30 @@
 		    /// @{
 		    /// \brief Return the total cost of the found flow.
 		    ///
 		    /// This function returns the total cost of the found flow.
-		    /// The complexity of the function is O(e).
+		    /// Its complexity is O(e).
 		    ///
 		    /// \note The return type of the function can be specified as a
 		    /// template parameter. For example,
 		    /// \code
 		    ///   ns.totalCost<double>();
 		    /// \endcode
 		    /// It is useful if the total cost cannot be stored in the \c Cost
 		    /// type of the algorithm, which is the default return type of the
 		    /// function.
 		    ///
 		    /// \pre \ref run() must be called before using this function.
 		    template <typename Num>
 		    Num totalCost() const {
 		      Num c = 0;
 		      if (_pcost) {
 		        for (ArcIt e(_graph); e != INVALID; ++e)
 		          c += (*_flow_map)[e] * (*_pcost)[e];
 		      } else {
 		        for (ArcIt e(_graph); e != INVALID; ++e)
-		          c += (*_flow_map)[e];
+		    template <typename Number>
 		    Number totalCost() const {
 		      Number c = 0;
 		      for (ArcIt a(_graph); a != INVALID; ++a) {
 		        int i = _arc_id[a];
 		        c += Number(_flow[i]) * Number(_cost[i]);
+		      }
 		      return c;
+		    }
 		#ifndef DOXYGEN
 		    Cost totalCost() const {
@@ -1018,293 +986,137 @@
 		    /// \brief Return the flow on the given arc.
 		    ///
 		    /// This function returns the flow on the given arc.
 		    ///
 		    /// \pre \ref run() must be called before using this function.
 		    Flow flow(const Arc& a) const {
 		      return (*_flow_map)[a];
 		    Value flow(const Arc& a) const {
 		      return _flow[_arc_id[a]];
+		    }
-		    /// \brief Return a const reference to the flow map.
+		    /// \brief Return the flow map (the primal solution).
 		    ///
 		    /// This function returns a const reference to an arc map storing
 		    /// the found flow.
 		    /// This function copies the flow value on each arc into the given
 		    /// map. The \c Value type of the algorithm must be convertible to
 		    /// the \c Value type of the map.
 		    ///
 		    /// \pre \ref run() must be called before using this function.
 		    const FlowMap& flowMap() const {
 		      return *_flow_map;
 		    template <typename FlowMap>
 		    void flowMap(FlowMap &map) const {
 		      for (ArcIt a(_graph); a != INVALID; ++a) {
 		        map.set(a, _flow[_arc_id[a]]);
+		      }
+		    }
 		    /// \brief Return the potential (dual value) of the given node.
 		    ///
 		    /// This function returns the potential (dual value) of the
 		    /// given node.
 		    ///
 		    /// \pre \ref run() must be called before using this function.
 		    Cost potential(const Node& n) const {
-		      return (*_potential_map)[n];
+		      return _pi[_node_id[n]];
+		    }
 		    /// \brief Return a const reference to the potential map
 		    /// (the dual solution).
 		    /// \brief Return the potential map (the dual solution).
 		    ///
 		    /// This function returns a const reference to a node map storing
 		    /// the found potentials, which form the dual solution of the
-		    /// \ref min_cost_flow "minimum cost flow" problem.
+		    /// This function copies the potential (dual value) of each node
 		    /// into the given map.
 		    /// The \c Cost type of the algorithm must be convertible to the
 		    /// \c Value type of the map.
 		    ///
 		    /// \pre \ref run() must be called before using this function.
 		    const PotentialMap& potentialMap() const {
 		      return *_potential_map;
 		    template <typename PotentialMap>
 		    void potentialMap(PotentialMap &map) const {
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        map.set(n, _pi[_node_id[n]]);
+		      }
+		    }
 		    /// @}
 		  private:
 		    // Initialize internal data structures
 		    bool init() {
 		      // Initialize result maps
 		      if (!_flow_map) {
 		        _flow_map = new FlowMap(_graph);
 		        _local_flow = true;
 		      if (_node_num == 0) return false;
 		      // Check the sum of supply values
 		      _sum_supply = 0;
 		      for (int i = 0; i != _node_num; ++i) {
 		        _sum_supply += _supply[i];
+		      }
 		      if (!_potential_map) {
 		        _potential_map = new PotentialMap(_graph);
-		        _local_potential = true;
+		      if ( !((_stype == GEQ && _sum_supply <= 0) ||
 		             (_stype == LEQ && _sum_supply >= 0)) ) return false;
 		      // Remove non-zero lower bounds
 		      if (_have_lower) {
 		        for (int i = 0; i != _arc_num; ++i) {
 		          Value c = _lower[i];
 		          if (c >= 0) {
 		            _cap[i] = _upper[i] < INF ? _upper[i] - c : INF;
 		          } else {
 		            _cap[i] = _upper[i] < INF + c ? _upper[i] - c : INF;
+		          }
 		          _supply[_source[i]] -= c;
 		          _supply[_target[i]] += c;
+		        }
 		      } else {
 		        for (int i = 0; i != _arc_num; ++i) {
 		          _cap[i] = _upper[i];
+		        }
+		      }
 		      // Initialize vectors
 		      _node_num = countNodes(_graph);
 		      _arc_num = countArcs(_graph);
 		      int all_node_num = _node_num + 1;
 		      int all_arc_num = _arc_num + _node_num;
 		      if (_node_num == 0) return false;
 		      _arc_ref.resize(_arc_num);
 		      _source.resize(all_arc_num);
 		      _target.resize(all_arc_num);
 		      _cap.resize(all_arc_num);
 		      _cost.resize(all_arc_num);
 		      _supply.resize(all_node_num);
 		      _flow.resize(all_arc_num);
 		      _pi.resize(all_node_num);
 		      _parent.resize(all_node_num);
 		      _pred.resize(all_node_num);
 		      _forward.resize(all_node_num);
 		      _thread.resize(all_node_num);
 		      _rev_thread.resize(all_node_num);
 		      _succ_num.resize(all_node_num);
 		      _last_succ.resize(all_node_num);
 		      _state.resize(all_arc_num);
 		      // Initialize node related data
 		      bool valid_supply = true;
 		      Flow sum_supply = 0;
 		      if (!_pstsup && !_psupply) {
 		        _pstsup = true;
 		        _psource = _ptarget = NodeIt(_graph);
 		        _pstflow = 0;
+		      }
 		      if (_psupply) {
 		        int i = 0;
 		        for (NodeIt n(_graph); n != INVALID; ++n, ++i) {
 		          _node_id[n] = i;
 		          _supply[i] = (*_psupply)[n];
 		          sum_supply += _supply[i];
 		      // Initialize artifical cost
 		      Cost ART_COST;
 		      if (std::numeric_limits<Cost>::is_exact) {
 		        ART_COST = std::numeric_limits<Cost>::max() / 4 + 1;
 		      } else {
 		        ART_COST = std::numeric_limits<Cost>::min();
 		        for (int i = 0; i != _arc_num; ++i) {
 		          if (_cost[i] > ART_COST) ART_COST = _cost[i];
+		        }
 		        valid_supply = (_ptype == GEQ && sum_supply <= 0) ||
 		                       (_ptype == LEQ && sum_supply >= 0);
 		      } else {
 		        int i = 0;
 		        for (NodeIt n(_graph); n != INVALID; ++n, ++i) {
 		          _node_id[n] = i;
 		          _supply[i] = 0;
+		        }
 		        _supply[_node_id[_psource]] =  _pstflow;
 		        _supply[_node_id[_ptarget]] = -_pstflow;
+		      }
 		      if (!valid_supply) return false;
 		      // Infinite capacity value
 		      Flow inf_cap =
 		        std::numeric_limits<Flow>::has_infinity ?
 		        std::numeric_limits<Flow>::infinity() :
 		        std::numeric_limits<Flow>::max();
 		      // Initialize artifical cost
 		      Cost art_cost;
 		      if (std::numeric_limits<Cost>::is_exact) {
 		        art_cost = std::numeric_limits<Cost>::max() / 4 + 1;
 		      } else {
 		        art_cost = std::numeric_limits<Cost>::min();
 		        for (int i = 0; i != _arc_num; ++i) {
 		          if (_cost[i] > art_cost) art_cost = _cost[i];
+		        }
-		        art_cost = (art_cost + 1) * _node_num;
+		        ART_COST = (ART_COST + 1) * _node_num;
+		      }
 		      // Run Circulation to check if a feasible solution exists
 		      typedef ConstMap<Arc, Flow> ConstArcMap;
 		      ConstArcMap zero_arc_map(0), inf_arc_map(inf_cap);
 		      FlowNodeMap *csup = NULL;
 		      bool local_csup = false;
 		      if (_psupply) {
 		        csup = _psupply;
 		      } else {
 		        csup = new FlowNodeMap(_graph, 0);
 		        (*csup)[_psource] =  _pstflow;
 		        (*csup)[_ptarget] = -_pstflow;
 		        local_csup = true;
 		      // Initialize arc maps
 		      for (int i = 0; i != _arc_num; ++i) {
 		        _flow[i] = 0;
 		        _state[i] = STATE_LOWER;
+		      }
 		      bool circ_result = false;
 		      if (_ptype == GEQ || (_ptype == LEQ && sum_supply == 0)) {
 		        // GEQ problem type
 		        if (_plower) {
 		          if (_pupper) {
 		            Circulation<GR, FlowArcMap, FlowArcMap, FlowNodeMap>
 		              circ(_graph, *_plower, *_pupper, *csup);
 		            circ_result = circ.run();
 		          } else {
 		            Circulation<GR, FlowArcMap, ConstArcMap, FlowNodeMap>
 		              circ(_graph, *_plower, inf_arc_map, *csup);
 		            circ_result = circ.run();
+		          }
 		        } else {
 		          if (_pupper) {
 		            Circulation<GR, ConstArcMap, FlowArcMap, FlowNodeMap>
 		              circ(_graph, zero_arc_map, *_pupper, *csup);
 		            circ_result = circ.run();
 		          } else {
 		            Circulation<GR, ConstArcMap, ConstArcMap, FlowNodeMap>
 		              circ(_graph, zero_arc_map, inf_arc_map, *csup);
 		            circ_result = circ.run();
+		          }
+		        }
 		      } else {
 		        // LEQ problem type
 		        typedef ReverseDigraph<const GR> RevGraph;
 		        typedef NegMap<FlowNodeMap> NegNodeMap;
 		        RevGraph rgraph(_graph);
 		        NegNodeMap neg_csup(*csup);
 		        if (_plower) {
 		          if (_pupper) {
 		            Circulation<RevGraph, FlowArcMap, FlowArcMap, NegNodeMap>
 		              circ(rgraph, *_plower, *_pupper, neg_csup);
 		            circ_result = circ.run();
 		          } else {
 		            Circulation<RevGraph, FlowArcMap, ConstArcMap, NegNodeMap>
 		              circ(rgraph, *_plower, inf_arc_map, neg_csup);
 		            circ_result = circ.run();
+		          }
 		        } else {
 		          if (_pupper) {
 		            Circulation<RevGraph, ConstArcMap, FlowArcMap, NegNodeMap>
 		              circ(rgraph, zero_arc_map, *_pupper, neg_csup);
 		            circ_result = circ.run();
 		          } else {
 		            Circulation<RevGraph, ConstArcMap, ConstArcMap, NegNodeMap>
 		              circ(rgraph, zero_arc_map, inf_arc_map, neg_csup);
 		            circ_result = circ.run();
+		          }
+		        }
+		      }
 		      if (local_csup) delete csup;
 		      if (!circ_result) return false;
 		      // Set data for the artificial root node
 		      _root = _node_num;
 		      _parent[_root] = -1;
 		      _pred[_root] = -1;
 		      _thread[_root] = 0;
 		      _rev_thread[0] = _root;
-		      _succ_num[_root] = all_node_num;
+		      _succ_num[_root] = _node_num + 1;
 		      _last_succ[_root] = _root - 1;
 		      _supply[_root] = -sum_supply;
 		      if (sum_supply < 0) {
 		        _pi[_root] = -art_cost;
 		      } else {
 		        _pi[_root] = art_cost;
+		      }
 		      // Store the arcs in a mixed order
 		      int k = std::max(int(std::sqrt(double(_arc_num))), 10);
 		      int i = 0;
 		      for (ArcIt e(_graph); e != INVALID; ++e) {
 		        _arc_ref[i] = e;
 		        if ((i += k) >= _arc_num) i = (i % k) + 1;
+		      }
 		      // Initialize arc maps
 		      if (_pupper && _pcost) {
 		        for (int i = 0; i != _arc_num; ++i) {
 		          Arc e = _arc_ref[i];
 		          _source[i] = _node_id[_graph.source(e)];
 		          _target[i] = _node_id[_graph.target(e)];
 		          _cap[i] = (*_pupper)[e];
 		          _cost[i] = (*_pcost)[e];
 		          _flow[i] = 0;
 		          _state[i] = STATE_LOWER;
+		        }
 		      } else {
 		        for (int i = 0; i != _arc_num; ++i) {
 		          Arc e = _arc_ref[i];
 		          _source[i] = _node_id[_graph.source(e)];
 		          _target[i] = _node_id[_graph.target(e)];
 		          _flow[i] = 0;
 		          _state[i] = STATE_LOWER;
+		        }
 		        if (_pupper) {
 		          for (int i = 0; i != _arc_num; ++i)
 		            _cap[i] = (*_pupper)[_arc_ref[i]];
 		        } else {
 		          for (int i = 0; i != _arc_num; ++i)
 		            _cap[i] = inf_cap;
+		        }
 		        if (_pcost) {
 		          for (int i = 0; i != _arc_num; ++i)
 		            _cost[i] = (*_pcost)[_arc_ref[i]];
 		        } else {
 		          for (int i = 0; i != _arc_num; ++i)
 		            _cost[i] = 1;
+		        }
+		      }
 		      // Remove non-zero lower bounds
 		      if (_plower) {
 		        for (int i = 0; i != _arc_num; ++i) {
 		          Flow c = (*_plower)[_arc_ref[i]];
 		          if (c != 0) {
 		            _cap[i] -= c;
 		            _supply[_source[i]] -= c;
 		            _supply[_target[i]] += c;
+		          }
+		        }
+		      }
+		      _supply[_root] = -_sum_supply;
 		      _pi[_root] = _sum_supply < 0 ? -ART_COST : ART_COST;
 		      // Add artificial arcs and initialize the spanning tree data structure
 		      for (int u = 0, e = _arc_num; u != _node_num; ++u, ++e) {
 		        _parent[u] = _root;
 		        _pred[u] = e;
 		        _thread[u] = u + 1;
 		        _rev_thread[u + 1] = u;
 		        _succ_num[u] = 1;
 		        _last_succ[u] = u;
 		        _parent[u] = _root;
 		        _pred[u] = e;
 		        _cost[e] = art_cost;
 		        _cap[e] = inf_cap;
 		        _cost[e] = ART_COST;
 		        _cap[e] = INF;
 		        _state[e] = STATE_TREE;
-		        if (_supply[u] > 0 || (_supply[u] == 0 && sum_supply <= 0)) {
+		        if (_supply[u] > 0 || (_supply[u] == 0 && _sum_supply <= 0)) {
 		          _flow[e] = _supply[u];
 		          _forward[u] = true;
-		          _pi[u] = -art_cost + _pi[_root];
+		          _pi[u] = -ART_COST + _pi[_root];
 		        } else {
 		          _flow[e] = -_supply[u];
 		          _forward[u] = false;
-		          _pi[u] = art_cost + _pi[_root];
+		          _pi[u] = ART_COST + _pi[_root];
+		        }
+		      }
 		      return true;
+		    }
@@ -1333,29 +1145,31 @@
 		      } else {
 		        first  = _target[in_arc];
 		        second = _source[in_arc];
+		      }
 		      delta = _cap[in_arc];
 		      int result = 0;
-		      Flow d;
+		      Value d;
 		      int e;
 		      // Search the cycle along the path form the first node to the root
 		      for (int u = first; u != join; u = _parent[u]) {
 		        e = _pred[u];
-		        d = _forward[u] ? _flow[e] : _cap[e] - _flow[e];
 		        d = _forward[u] ?
 		          _flow[e] : (_cap[e] == INF ? INF : _cap[e] - _flow[e]);
 		        if (d < delta) {
 		          delta = d;
 		          u_out = u;
 		          result = 1;
+		        }
+		      }
 		      // Search the cycle along the path form the second node to the root
 		      for (int u = second; u != join; u = _parent[u]) {
 		        e = _pred[u];
-		        d = _forward[u] ? _cap[e] - _flow[e] : _flow[e];
 		        d = _forward[u] ?
 		          (_cap[e] == INF ? INF : _cap[e] - _flow[e]) : _flow[e];
 		        if (d <= delta) {
 		          delta = d;
 		          u_out = u;
 		          result = 2;
+		        }
+		      }
@@ -1371,13 +1185,13 @@
+		    }
 		    // Change _flow and _state vectors
 		    void changeFlow(bool change) {
 		      // Augment along the cycle
 		      if (delta > 0) {
-		        Flow val = _state[in_arc] * delta;
+		        Value val = _state[in_arc] * delta;
 		        _flow[in_arc] += val;
 		        for (int u = _source[in_arc]; u != join; u = _parent[u]) {
 		          _flow[_pred[u]] += _forward[u] ? -val : val;
+		        }
 		        for (int u = _target[in_arc]; u != join; u = _parent[u]) {
 		          _flow[_pred[u]] += _forward[u] ? val : -val;
@@ -1523,13 +1337,13 @@
 		      for (int u = u_in; u != end; u = _thread[u]) {
 		        _pi[u] += sigma;
+		      }
+		    }
 		    // Execute the algorithm
-		    bool start(PivotRule pivot_rule) {
+		    ProblemType start(PivotRule pivot_rule) {
 		      // Select the pivot rule implementation
 		      switch (pivot_rule) {
 		        case FIRST_ELIGIBLE:
 		          return start<FirstEligiblePivotRule>();
 		        case BEST_ELIGIBLE:
 		          return start<BestEligiblePivotRule>();
@@ -1537,47 +1351,61 @@
 		          return start<BlockSearchPivotRule>();
 		        case CANDIDATE_LIST:
 		          return start<CandidateListPivotRule>();
 		        case ALTERING_LIST:
 		          return start<AlteringListPivotRule>();
+		      }
-		      return false;
+		      return INFEASIBLE; // avoid warning
+		    }
 		    template <typename PivotRuleImpl>
-		    bool start() {
+		    ProblemType start() {
 		      PivotRuleImpl pivot(*this);
 		      // Execute the Network Simplex algorithm
 		      while (pivot.findEnteringArc()) {
 		        findJoinNode();
 		        bool change = findLeavingArc();
 		        if (delta >= INF) return UNBOUNDED;
 		        changeFlow(change);
 		        if (change) {
 		          updateTreeStructure();
 		          updatePotential();
+		        }
+		      }
 		      // Copy flow values to _flow_map
 		      if (_plower) {
 		        for (int i = 0; i != _arc_num; ++i) {
 		          Arc e = _arc_ref[i];
 		          _flow_map->set(e, (*_plower)[e] + _flow[i]);
+		        }
 		      } else {
 		        for (int i = 0; i != _arc_num; ++i) {
 		          _flow_map->set(_arc_ref[i], _flow[i]);
 		      // Check feasibility
 		      if (_sum_supply < 0) {
 		        for (int u = 0, e = _arc_num; u != _node_num; ++u, ++e) {
 		          if (_supply[u] >= 0 && _flow[e] != 0) return INFEASIBLE;
+		        }
+		      }
 		      // Copy potential values to _potential_map
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
-		        _potential_map->set(n, _pi[_node_id[n]]);
+		      else if (_sum_supply > 0) {
 		        for (int u = 0, e = _arc_num; u != _node_num; ++u, ++e) {
 		          if (_supply[u] <= 0 && _flow[e] != 0) return INFEASIBLE;
+		        }
+		      }
 		      else {
 		        for (int u = 0, e = _arc_num; u != _node_num; ++u, ++e) {
 		          if (_flow[e] != 0) return INFEASIBLE;
+		        }
+		      }
-		      return true;
+		      // Transform the solution and the supply map to the original form
 		      if (_have_lower) {
 		        for (int i = 0; i != _arc_num; ++i) {
 		          Value c = _lower[i];
 		          if (c != 0) {
 		            _flow[i] += c;
 		            _supply[_source[i]] += c;
 		            _supply[_target[i]] -= c;
+		          }
+		        }
+		      }
 		      return OPTIMAL;
+		    }
 		  }; //class NetworkSimplex
 		  ///@}

lemon/preflow.h

Show inline comments

@@ -43,19 +43,19 @@
 		    ///
 		    /// The type of the map that stores the arc capacities.
 		    /// It must meet the \ref concepts::ReadMap "ReadMap" concept.
 		    typedef CAP CapacityMap;
 		    /// \brief The type of the flow values.
-		    typedef typename CapacityMap::Value Flow;
+		    typedef typename CapacityMap::Value Value;
 		    /// \brief The type of the map that stores the flow values.
 		    ///
 		    /// The type of the map that stores the flow values.
 		    /// It must meet the \ref concepts::ReadWriteMap "ReadWriteMap" concept.
-		    typedef typename Digraph::template ArcMap<Flow> FlowMap;
+		    typedef typename Digraph::template ArcMap<Value> FlowMap;
 		    /// \brief Instantiates a FlowMap.
 		    ///
 		    /// This function instantiates a \ref FlowMap.
 		    /// \param digraph The digraph for which we would like to define
 		    /// the flow map.
@@ -81,13 +81,13 @@
 		      return new Elevator(digraph, max_level);
+		    }
 		    /// \brief The tolerance used by the algorithm
 		    ///
 		    /// The tolerance used by the algorithm to handle inexact computation.
-		    typedef lemon::Tolerance<Flow> Tolerance;
+		    typedef lemon::Tolerance<Value> Tolerance;
 		  };
 		  /// \ingroup max_flow
 		  ///
@@ -122,13 +122,13 @@
 		    typedef TR Traits;
 		    ///The type of the digraph the algorithm runs on.
 		    typedef typename Traits::Digraph Digraph;
 		    ///The type of the capacity map.
 		    typedef typename Traits::CapacityMap CapacityMap;
 		    ///The type of the flow values.
-		    typedef typename Traits::Flow Flow;
+		    typedef typename Traits::Value Value;
 		    ///The type of the flow map.
 		    typedef typename Traits::FlowMap FlowMap;
 		    ///The type of the elevator.
 		    typedef typename Traits::Elevator Elevator;
 		    ///The type of the tolerance.
@@ -148,13 +148,13 @@
 		    FlowMap* _flow;
 		    bool _local_flow;
 		    Elevator* _level;
 		    bool _local_level;
-		    typedef typename Digraph::template NodeMap<Flow> ExcessMap;
+		    typedef typename Digraph::template NodeMap<Value> ExcessMap;
 		    ExcessMap* _excess;
 		    Tolerance _tolerance;
 		    bool _phase;
@@ -467,13 +467,13 @@
 		      for (ArcIt e(_graph); e != INVALID; ++e) {
 		        _flow->set(e, flowMap[e]);
+		      }
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
-		        Flow excess = 0;
+		        Value excess = 0;
 		        for (InArcIt e(_graph, n); e != INVALID; ++e) {
 		          excess += (*_flow)[e];
+		        }
 		        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
 		          excess -= (*_flow)[e];
+		        }
@@ -516,25 +516,25 @@
+		        }
 		        queue.swap(nqueue);
+		      }
 		      _level->initFinish();
 		      for (OutArcIt e(_graph, _source); e != INVALID; ++e) {
-		        Flow rem = (*_capacity)[e] - (*_flow)[e];
+		        Value rem = (*_capacity)[e] - (*_flow)[e];
 		        if (_tolerance.positive(rem)) {
 		          Node u = _graph.target(e);
 		          if ((*_level)[u] == _level->maxLevel()) continue;
 		          _flow->set(e, (*_capacity)[e]);
 		          (*_excess)[u] += rem;
 		          if (u != _target && !_level->active(u)) {
 		            _level->activate(u);
+		          }
+		        }
+		      }
 		      for (InArcIt e(_graph, _source); e != INVALID; ++e) {
-		        Flow rem = (*_flow)[e];
+		        Value rem = (*_flow)[e];
 		        if (_tolerance.positive(rem)) {
 		          Node v = _graph.source(e);
 		          if ((*_level)[v] == _level->maxLevel()) continue;
 		          _flow->set(e, 0);
 		          (*_excess)[v] += rem;
 		          if (v != _target && !_level->active(v)) {
@@ -561,17 +561,17 @@
 		      Node n = _level->highestActive();
 		      int level = _level->highestActiveLevel();
 		      while (n != INVALID) {
 		        int num = _node_num;
 		        while (num > 0 && n != INVALID) {
-		          Flow excess = (*_excess)[n];
+		          Value excess = (*_excess)[n];
 		          int new_level = _level->maxLevel();
 		          for (OutArcIt e(_graph, n); e != INVALID; ++e) {
-		            Flow rem = (*_capacity)[e] - (*_flow)[e];
+		            Value rem = (*_capacity)[e] - (*_flow)[e];
 		            if (!_tolerance.positive(rem)) continue;
 		            Node v = _graph.target(e);
 		            if ((*_level)[v] < level) {
 		              if (!_level->active(v) && v != _target) {
 		                _level->activate(v);
+		              }
@@ -588,13 +588,13 @@
 		            } else if (new_level > (*_level)[v]) {
 		              new_level = (*_level)[v];
+		            }
+		          }
 		          for (InArcIt e(_graph, n); e != INVALID; ++e) {
-		            Flow rem = (*_flow)[e];
+		            Value rem = (*_flow)[e];
 		            if (!_tolerance.positive(rem)) continue;
 		            Node v = _graph.source(e);
 		            if ((*_level)[v] < level) {
 		              if (!_level->active(v) && v != _target) {
 		                _level->activate(v);
+		              }
@@ -634,17 +634,17 @@
 		          level = _level->highestActiveLevel();
 		          --num;
+		        }
 		        num = _node_num * 20;
 		        while (num > 0 && n != INVALID) {
-		          Flow excess = (*_excess)[n];
+		          Value excess = (*_excess)[n];
 		          int new_level = _level->maxLevel();
 		          for (OutArcIt e(_graph, n); e != INVALID; ++e) {
-		            Flow rem = (*_capacity)[e] - (*_flow)[e];
+		            Value rem = (*_capacity)[e] - (*_flow)[e];
 		            if (!_tolerance.positive(rem)) continue;
 		            Node v = _graph.target(e);
 		            if ((*_level)[v] < level) {
 		              if (!_level->active(v) && v != _target) {
 		                _level->activate(v);
+		              }
@@ -661,13 +661,13 @@
 		            } else if (new_level > (*_level)[v]) {
 		              new_level = (*_level)[v];
+		            }
+		          }
 		          for (InArcIt e(_graph, n); e != INVALID; ++e) {
-		            Flow rem = (*_flow)[e];
+		            Value rem = (*_flow)[e];
 		            if (!_tolerance.positive(rem)) continue;
 		            Node v = _graph.source(e);
 		            if ((*_level)[v] < level) {
 		              if (!_level->active(v) && v != _target) {
 		                _level->activate(v);
+		              }
@@ -775,18 +775,18 @@
 		          _level->activate(n);
+		        }
+		      }
 		      Node n;
 		      while ((n = _level->highestActive()) != INVALID) {
-		        Flow excess = (*_excess)[n];
+		        Value excess = (*_excess)[n];
 		        int level = _level->highestActiveLevel();
 		        int new_level = _level->maxLevel();
 		        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
-		          Flow rem = (*_capacity)[e] - (*_flow)[e];
+		          Value rem = (*_capacity)[e] - (*_flow)[e];
 		          if (!_tolerance.positive(rem)) continue;
 		          Node v = _graph.target(e);
 		          if ((*_level)[v] < level) {
 		            if (!_level->active(v) && v != _source) {
 		              _level->activate(v);
+		            }
@@ -803,13 +803,13 @@
 		          } else if (new_level > (*_level)[v]) {
 		            new_level = (*_level)[v];
+		          }
+		        }
 		        for (InArcIt e(_graph, n); e != INVALID; ++e) {
-		          Flow rem = (*_flow)[e];
+		          Value rem = (*_flow)[e];
 		          if (!_tolerance.positive(rem)) continue;
 		          Node v = _graph.source(e);
 		          if ((*_level)[v] < level) {
 		            if (!_level->active(v) && v != _source) {
 		              _level->activate(v);
+		            }
@@ -894,24 +894,24 @@
 		    /// Returns the value of the maximum flow by returning the excess
 		    /// of the target node. This value equals to the value of
 		    /// the maximum flow already after the first phase of the algorithm.
 		    ///
 		    /// \pre Either \ref run() or \ref init() must be called before
 		    /// using this function.
-		    Flow flowValue() const {
+		    Value flowValue() const {
 		      return (*_excess)[_target];
+		    }
 		    /// \brief Returns the flow on the given arc.
+		    /// \brief Returns the flow value on the given arc.
 		    ///
 		    /// Returns the flow on the given arc. This method can
+		    /// Returns the flow value on the given arc. This method can
 		    /// be called after the second phase of the algorithm.
 		    ///
 		    /// \pre Either \ref run() or \ref init() must be called before
 		    /// using this function.
-		    Flow flow(const Arc& arc) const {
+		    Value flow(const Arc& arc) const {
 		      return (*_flow)[arc];
+		    }
 		    /// \brief Returns a const reference to the flow map.
 		    ///
 		    /// Returns a const reference to the arc map storing the found flow.

test/lp_test.cc

Show inline comments

@@ -18,15 +18,13 @@
 		#include <sstream>
 		#include <lemon/lp_skeleton.h>
 		#include "test_tools.h"
 		#include <lemon/tolerance.h>
 		#ifdef HAVE_CONFIG_H
 		#include <lemon/config.h>
 		#endif
 		#ifdef LEMON_HAVE_GLPK
 		#include <lemon/glpk.h>
 		#endif
 		#ifdef LEMON_HAVE_CPLEX

test/min_cost_flow_test.cc

Show inline comments

@@ -15,12 +15,13 @@
 		 * purpose.
+		 *
 		 */
 		#include <iostream>
 		#include <fstream>
 		#include <limits>
 		#include <lemon/list_graph.h>
 		#include <lemon/lgf_reader.h>
 		#include <lemon/network_simplex.h>
@@ -30,131 +31,124 @@
 		#include "test_tools.h"
 		using namespace lemon;
 		char test_lgf[] =
 		  "@nodes\n"
 		  "label  sup1 sup2 sup3 sup4 sup5\n"
 		  "    1    20   27    0   20   30\n"
 		  "    2    -4    0    0   -8   -3\n"
 		  "    3     0    0    0    0    0\n"
 		  "    4     0    0    0    0    0\n"
 		  "    5     9    0    0    6   11\n"
 		  "    6    -6    0    0   -5   -6\n"
 		  "    7     0    0    0    0    0\n"
 		  "    8     0    0    0    0    3\n"
 		  "    9     3    0    0    0    0\n"
 		  "   10    -2    0    0   -7   -2\n"
 		  "   11     0    0    0  -10    0\n"
 		  "   12   -20  -27    0  -30  -20\n"
 		  "\n"
 		  "label  sup1 sup2 sup3 sup4 sup5 sup6\n"
 		  "    1    20   27    0   30   20   30\n"
 		  "    2    -4    0    0    0   -8   -3\n"
 		  "    3     0    0    0    0    0    0\n"
 		  "    4     0    0    0    0    0    0\n"
 		  "    5     9    0    0    0    6   11\n"
 		  "    6    -6    0    0    0   -5   -6\n"
 		  "    7     0    0    0    0    0    0\n"
 		  "    8     0    0    0    0    0    3\n"
 		  "    9     3    0    0    0    0    0\n"
 		  "   10    -2    0    0    0   -7   -2\n"
 		  "   11     0    0    0    0  -10    0\n"
 		  "   12   -20  -27    0  -30  -30  -20\n"
 		  "\n"
 		  "@arcs\n"
 		  "       cost  cap low1 low2\n"
 		  " 1  2    70   11    0    8\n"
 		  " 1  3   150    3    0    1\n"
 		  " 1  4    80   15    0    2\n"
 		  " 2  8    80   12    0    0\n"
 		  " 3  5   140    5    0    3\n"
 		  " 4  6    60   10    0    1\n"
 		  " 4  7    80    2    0    0\n"
 		  " 4  8   110    3    0    0\n"
 		  " 5  7    60   14    0    0\n"
 		  " 5 11   120   12    0    0\n"
 		  " 6  3     0    3    0    0\n"
 		  " 6  9   140    4    0    0\n"
 		  " 6 10    90    8    0    0\n"
 		  " 7  1    30    5    0    0\n"
 		  " 8 12    60   16    0    4\n"
 		  " 9 12    50    6    0    0\n"
 		  "10 12    70   13    0    5\n"
 		  "10  2   100    7    0    0\n"
 		  "10  7    60   10    0    0\n"
 		  "11 10    20   14    0    6\n"
 		  "12 11    30   10    0    0\n"
 		  "       cost  cap low1 low2 low3\n"
 		  " 1  2    70   11    0    8    8\n"
 		  " 1  3   150    3    0    1    0\n"
 		  " 1  4    80   15    0    2    2\n"
 		  " 2  8    80   12    0    0    0\n"
 		  " 3  5   140    5    0    3    1\n"
 		  " 4  6    60   10    0    1    0\n"
 		  " 4  7    80    2    0    0    0\n"
 		  " 4  8   110    3    0    0    0\n"
 		  " 5  7    60   14    0    0    0\n"
 		  " 5 11   120   12    0    0    0\n"
 		  " 6  3     0    3    0    0    0\n"
 		  " 6  9   140    4    0    0    0\n"
 		  " 6 10    90    8    0    0    0\n"
 		  " 7  1    30    5    0    0   -5\n"
 		  " 8 12    60   16    0    4    3\n"
 		  " 9 12    50    6    0    0    0\n"
 		  "10 12    70   13    0    5    2\n"
 		  "10  2   100    7    0    0    0\n"
 		  "10  7    60   10    0    0   -3\n"
 		  "11 10    20   14    0    6  -20\n"
 		  "12 11    30   10    0    0  -10\n"
 		  "\n"
 		  "@attributes\n"
 		  "source 1\n"
 		  "target 12\n";
-		enum ProblemType {
+		enum SupplyType {
 		  EQ,
 		  GEQ,
 		  LEQ
 		};
 		// Check the interface of an MCF algorithm
-		template <typename GR, typename Flow, typename Cost>
+		template <typename GR, typename Value, typename Cost>
 		class McfClassConcept
+		{
 		public:
 		  template <typename MCF>
 		  struct Constraints {
 		    void constraints() {
 		      checkConcept<concepts::Digraph, GR>();
 		      MCF mcf(g);
 		      const MCF& const_mcf = mcf;
 		      b = mcf.reset()
 		             .lowerMap(lower)
 		             .upperMap(upper)
 		             .capacityMap(upper)
 		             .boundMaps(lower, upper)
 		             .costMap(cost)
 		             .supplyMap(sup)
 		             .stSupply(n, n, k)
 		             .flowMap(flow)
 		             .potentialMap(pot)
 		             .run();
 		      const MCF& const_mcf = mcf;
 		      const typename MCF::FlowMap &fm = const_mcf.flowMap();
 		      const typename MCF::PotentialMap &pm = const_mcf.potentialMap();
 		      v = const_mcf.totalCost();
-		      double x = const_mcf.template totalCost<double>();
+		      c = const_mcf.totalCost();
 		      x = const_mcf.template totalCost<double>();
 		      v = const_mcf.flow(a);
 		      v = const_mcf.potential(n);
 		      ignore_unused_variable_warning(fm);
 		      ignore_unused_variable_warning(pm);
-		      ignore_unused_variable_warning(x);
+		      c = const_mcf.potential(n);
 		      const_mcf.flowMap(fm);
 		      const_mcf.potentialMap(pm);
+		    }
 		    typedef typename GR::Node Node;
 		    typedef typename GR::Arc Arc;
 		    typedef concepts::ReadMap<Node, Flow> NM;
 		    typedef concepts::ReadMap<Arc, Flow> FAM;
 		    typedef concepts::ReadMap<Node, Value> NM;
 		    typedef concepts::ReadMap<Arc, Value> VAM;
 		    typedef concepts::ReadMap<Arc, Cost> CAM;
 		    typedef concepts::WriteMap<Arc, Value> FlowMap;
 		    typedef concepts::WriteMap<Node, Cost> PotMap;
 		    const GR &g;
 		    const FAM &lower;
 		    const FAM &upper;
 		    const VAM &lower;
 		    const VAM &upper;
 		    const CAM &cost;
 		    const NM &sup;
 		    const Node &n;
 		    const Arc &a;
 		    const Flow &k;
 		    Flow v;
 		    const Value &k;
 		    FlowMap fm;
 		    PotMap pm;
 		    bool b;
 		    typename MCF::FlowMap &flow;
-		    typename MCF::PotentialMap &pot;
+		    double x;
 		    typename MCF::Value v;
 		    typename MCF::Cost c;
 		  };
 		};
 		// Check the feasibility of the given flow (primal soluiton)
 		template < typename GR, typename LM, typename UM,
 		           typename SM, typename FM >
 		bool checkFlow( const GR& gr, const LM& lower, const UM& upper,
 		                const SM& supply, const FM& flow,
-		                ProblemType type = EQ )
+		                SupplyType type = EQ )
+		{
 		  TEMPLATE_DIGRAPH_TYPEDEFS(GR);
 		  for (ArcIt e(gr); e != INVALID; ++e) {
 		    if (flow[e] < lower[e] || flow[e] > upper[e]) return false;
+		  }
@@ -205,130 +199,192 @@
 		  return opt;
+		}
 		// Run a minimum cost flow algorithm and check the results
 		template < typename MCF, typename GR,
 		           typename LM, typename UM,
 		           typename CM, typename SM >
 		void checkMcf( const MCF& mcf, bool mcf_result,
 		           typename CM, typename SM,
 		           typename PT >
 		void checkMcf( const MCF& mcf, PT mcf_result,
 		               const GR& gr, const LM& lower, const UM& upper,
 		               const CM& cost, const SM& supply,
-		               bool result, typename CM::Value total,
+		               PT result, bool optimal, typename CM::Value total,
 		               const std::string &test_id = "",
-		               ProblemType type = EQ )
+		               SupplyType type = EQ )
+		{
 		  check(mcf_result == result, "Wrong result " + test_id);
 		  if (result) {
 		    check(checkFlow(gr, lower, upper, supply, mcf.flowMap(), type),
 		  if (optimal) {
 		    typename GR::template ArcMap<typename SM::Value> flow(gr);
 		    typename GR::template NodeMap<typename CM::Value> pi(gr);
 		    mcf.flowMap(flow);
 		    mcf.potentialMap(pi);
 		    check(checkFlow(gr, lower, upper, supply, flow, type),
 		          "The flow is not feasible " + test_id);
 		    check(mcf.totalCost() == total, "The flow is not optimal " + test_id);
 		    check(checkPotential(gr, lower, upper, cost, supply, mcf.flowMap(),
 		                         mcf.potentialMap()),
 		    check(checkPotential(gr, lower, upper, cost, supply, flow, pi),
 		          "Wrong potentials " + test_id);
+		  }
+		}
 		int main()
+		{
 		  // Check the interfaces
+		  {
 		    typedef int Flow;
 		    typedef int Cost;
 		    typedef concepts::Digraph GR;
 		    checkConcept< McfClassConcept<GR, Flow, Cost>,
 		                  NetworkSimplex<GR, Flow, Cost> >();
 		    checkConcept< McfClassConcept<GR, int, int>,
 		                  NetworkSimplex<GR> >();
 		    checkConcept< McfClassConcept<GR, double, double>,
 		                  NetworkSimplex<GR, double> >();
 		    checkConcept< McfClassConcept<GR, int, double>,
 		                  NetworkSimplex<GR, int, double> >();
+		  }
 		  // Run various MCF tests
 		  typedef ListDigraph Digraph;
 		  DIGRAPH_TYPEDEFS(ListDigraph);
 		  // Read the test digraph
 		  Digraph gr;
 		  Digraph::ArcMap<int> c(gr), l1(gr), l2(gr), u(gr);
 		  Digraph::NodeMap<int> s1(gr), s2(gr), s3(gr), s4(gr), s5(gr);
 		  Digraph::ArcMap<int> c(gr), l1(gr), l2(gr), l3(gr), u(gr);
 		  Digraph::NodeMap<int> s1(gr), s2(gr), s3(gr), s4(gr), s5(gr), s6(gr);
 		  ConstMap<Arc, int> cc(1), cu(std::numeric_limits<int>::max());
 		  Node v, w;
 		  std::istringstream input(test_lgf);
 		  DigraphReader<Digraph>(gr, input)
 		    .arcMap("cost", c)
 		    .arcMap("cap", u)
 		    .arcMap("low1", l1)
 		    .arcMap("low2", l2)
 		    .arcMap("low3", l3)
 		    .nodeMap("sup1", s1)
 		    .nodeMap("sup2", s2)
 		    .nodeMap("sup3", s3)
 		    .nodeMap("sup4", s4)
 		    .nodeMap("sup5", s5)
 		    .nodeMap("sup6", s6)
 		    .node("source", v)
 		    .node("target", w)
 		    .run();
 		  // Build a test digraph for testing negative costs
 		  Digraph ngr;
 		  Node n1 = ngr.addNode();
 		  Node n2 = ngr.addNode();
 		  Node n3 = ngr.addNode();
 		  Node n4 = ngr.addNode();
 		  Node n5 = ngr.addNode();
 		  Node n6 = ngr.addNode();
 		  Node n7 = ngr.addNode();
 		  Arc a1 = ngr.addArc(n1, n2);
 		  Arc a2 = ngr.addArc(n1, n3);
 		  Arc a3 = ngr.addArc(n2, n4);
 		  Arc a4 = ngr.addArc(n3, n4);
 		  Arc a5 = ngr.addArc(n3, n2);
 		  Arc a6 = ngr.addArc(n5, n3);
 		  Arc a7 = ngr.addArc(n5, n6);
 		  Arc a8 = ngr.addArc(n6, n7);
 		  Arc a9 = ngr.addArc(n7, n5);
 		  Digraph::ArcMap<int> nc(ngr), nl1(ngr, 0), nl2(ngr, 0);
 		  ConstMap<Arc, int> nu1(std::numeric_limits<int>::max()), nu2(5000);
 		  Digraph::NodeMap<int> ns(ngr, 0);
 		  nl2[a7] =  1000;
 		  nl2[a8] = -1000;
 		  ns[n1] =  100;
 		  ns[n4] = -100;
 		  nc[a1] =  100;
 		  nc[a2] =   30;
 		  nc[a3] =   20;
 		  nc[a4] =   80;
 		  nc[a5] =   50;
 		  nc[a6] =   10;
 		  nc[a7] =   80;
 		  nc[a8] =   30;
 		  nc[a9] = -120;
 		  // A. Test NetworkSimplex with the default pivot rule
+		  {
 		    NetworkSimplex<Digraph> mcf(gr);
 		    // Check the equality form
 		    mcf.upperMap(u).costMap(c);
 		    checkMcf(mcf, mcf.supplyMap(s1).run(),
 		             gr, l1, u, c, s1, true,  5240, "#A1");
+		             gr, l1, u, c, s1, mcf.OPTIMAL, true,   5240, "#A1");
 		    checkMcf(mcf, mcf.stSupply(v, w, 27).run(),
 		             gr, l1, u, c, s2, true,  7620, "#A2");
+		             gr, l1, u, c, s2, mcf.OPTIMAL, true,   7620, "#A2");
 		    mcf.lowerMap(l2);
 		    checkMcf(mcf, mcf.supplyMap(s1).run(),
 		             gr, l2, u, c, s1, true,  5970, "#A3");
+		             gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#A3");
 		    checkMcf(mcf, mcf.stSupply(v, w, 27).run(),
 		             gr, l2, u, c, s2, true,  8010, "#A4");
+		             gr, l2, u, c, s2, mcf.OPTIMAL, true,   8010, "#A4");
 		    mcf.reset();
 		    checkMcf(mcf, mcf.supplyMap(s1).run(),
 		             gr, l1, cu, cc, s1, true,  74, "#A5");
+		             gr, l1, cu, cc, s1, mcf.OPTIMAL, true,   74, "#A5");
 		    checkMcf(mcf, mcf.lowerMap(l2).stSupply(v, w, 27).run(),
 		             gr, l2, cu, cc, s2, true,  94, "#A6");
+		             gr, l2, cu, cc, s2, mcf.OPTIMAL, true,   94, "#A6");
 		    mcf.reset();
 		    checkMcf(mcf, mcf.run(),
 		             gr, l1, cu, cc, s3, true,   0, "#A7");
 		    checkMcf(mcf, mcf.boundMaps(l2, u).run(),
-		             gr, l2, u, cc, s3, false,   0, "#A8");
+		             gr, l1, cu, cc, s3, mcf.OPTIMAL, true,    0, "#A7");
 		    checkMcf(mcf, mcf.lowerMap(l2).upperMap(u).run(),
 		             gr, l2, u, cc, s3, mcf.INFEASIBLE, false, 0, "#A8");
 		    mcf.reset().lowerMap(l3).upperMap(u).costMap(c).supplyMap(s4);
 		    checkMcf(mcf, mcf.run(),
 		             gr, l3, u, c, s4, mcf.OPTIMAL, true,   6360, "#A9");
 		    // Check the GEQ form
-		    mcf.reset().upperMap(u).costMap(c).supplyMap(s4);
+		    mcf.reset().upperMap(u).costMap(c).supplyMap(s5);
 		    checkMcf(mcf, mcf.run(),
 		             gr, l1, u, c, s4, true,  3530, "#A9", GEQ);
 		    mcf.problemType(mcf.GEQ);
 		             gr, l1, u, c, s5, mcf.OPTIMAL, true,   3530, "#A10", GEQ);
 		    mcf.supplyType(mcf.GEQ);
 		    checkMcf(mcf, mcf.lowerMap(l2).run(),
 		             gr, l2, u, c, s4, true,  4540, "#A10", GEQ);
 		    mcf.problemType(mcf.CARRY_SUPPLIES).supplyMap(s5);
 		             gr, l2, u, c, s5, mcf.OPTIMAL, true,   4540, "#A11", GEQ);
 		    mcf.supplyType(mcf.CARRY_SUPPLIES).supplyMap(s6);
 		    checkMcf(mcf, mcf.run(),
-		             gr, l2, u, c, s5, false,    0, "#A11", GEQ);
+		             gr, l2, u, c, s6, mcf.INFEASIBLE, false,  0, "#A12", GEQ);
 		    // Check the LEQ form
 		    mcf.reset().problemType(mcf.LEQ);
 		    mcf.upperMap(u).costMap(c).supplyMap(s5);
 		    mcf.reset().supplyType(mcf.LEQ);
 		    mcf.upperMap(u).costMap(c).supplyMap(s6);
 		    checkMcf(mcf, mcf.run(),
-		             gr, l1, u, c, s5, true,  5080, "#A12", LEQ);
+		             gr, l1, u, c, s6, mcf.OPTIMAL, true,   5080, "#A13", LEQ);
 		    checkMcf(mcf, mcf.lowerMap(l2).run(),
 		             gr, l2, u, c, s5, true,  5930, "#A13", LEQ);
 		    mcf.problemType(mcf.SATISFY_DEMANDS).supplyMap(s4);
 		             gr, l2, u, c, s6, mcf.OPTIMAL, true,   5930, "#A14", LEQ);
 		    mcf.supplyType(mcf.SATISFY_DEMANDS).supplyMap(s5);
 		    checkMcf(mcf, mcf.run(),
-		             gr, l2, u, c, s4, false,    0, "#A14", LEQ);
+		             gr, l2, u, c, s5, mcf.INFEASIBLE, false,  0, "#A15", LEQ);
 		    // Check negative costs
 		    NetworkSimplex<Digraph> nmcf(ngr);
 		    nmcf.lowerMap(nl1).costMap(nc).supplyMap(ns);
 		    checkMcf(nmcf, nmcf.run(),
 		      ngr, nl1, nu1, nc, ns, nmcf.UNBOUNDED, false,    0, "#A16");
 		    checkMcf(nmcf, nmcf.upperMap(nu2).run(),
 		      ngr, nl1, nu2, nc, ns, nmcf.OPTIMAL, true,  -40000, "#A17");
 		    nmcf.reset().lowerMap(nl2).costMap(nc).supplyMap(ns);
 		    checkMcf(nmcf, nmcf.run(),
 		      ngr, nl2, nu1, nc, ns, nmcf.UNBOUNDED, false,    0, "#A18");
+		  }
 		  // B. Test NetworkSimplex with each pivot rule
+		  {
 		    NetworkSimplex<Digraph> mcf(gr);
-		    mcf.supplyMap(s1).costMap(c).capacityMap(u).lowerMap(l2);
+		    mcf.supplyMap(s1).costMap(c).upperMap(u).lowerMap(l2);
 		    checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::FIRST_ELIGIBLE),
 		             gr, l2, u, c, s1, true,  5970, "#B1");
+		             gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#B1");
 		    checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::BEST_ELIGIBLE),
 		             gr, l2, u, c, s1, true,  5970, "#B2");
+		             gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#B2");
 		    checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::BLOCK_SEARCH),
 		             gr, l2, u, c, s1, true,  5970, "#B3");
+		             gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#B3");
 		    checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::CANDIDATE_LIST),
 		             gr, l2, u, c, s1, true,  5970, "#B4");
+		             gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#B4");
 		    checkMcf(mcf, mcf.run(NetworkSimplex<Digraph>::ALTERING_LIST),
 		             gr, l2, u, c, s1, true,  5970, "#B5");
+		             gr, l2, u, c, s1, mcf.OPTIMAL, true,   5970, "#B5");
+		  }
 		  return 0;
+		}

test/mip_test.cc

Show inline comments

@@ -15,15 +15,13 @@
 		 * purpose.
+		 *
 		 */
 		#include "test_tools.h"
 		#ifdef HAVE_CONFIG_H
 		#include <lemon/config.h>
 		#endif
 		#ifdef LEMON_HAVE_CPLEX
 		#include <lemon/cplex.h>
 		#endif
 		#ifdef LEMON_HAVE_GLPK

tools/dimacs-solver.cc

Show inline comments

@@ -116,14 +116,14 @@
 		      std::cerr << "LEQ supply contraints are used for NetworkSimplex\n\n";
+		  }
 		  if (report) std::cerr << "Read the file: " << ti << '\n';
 		  ti.restart();
 		  NetworkSimplex<Digraph, Value> ns(g);
 		  ns.lowerMap(lower).capacityMap(cap).costMap(cost).supplyMap(sup);
 		  if (sum_sup > 0) ns.problemType(ns.LEQ);
 		  ns.lowerMap(lower).upperMap(cap).costMap(cost).supplyMap(sup);
 		  if (sum_sup > 0) ns.supplyType(ns.LEQ);
 		  if (report) std::cerr << "Setup NetworkSimplex class: " << ti << '\n';
 		  ti.restart();
 		  bool res = ns.run();
 		  if (report) {
 		    std::cerr << "Run NetworkSimplex: " << ti << "\n\n";
 		    std::cerr << "Feasible flow: " << (res ? "found" : "not found") << '\n';

0 comments (0 inline)

RhodeCode

Login to your account