RhodeCode

if HAVE_SOPLEX

lemon_libemon_la_SOURCES += lemon/soplex.cc

endif

if HAVE_CLP

lemon_libemon_la_SOURCES += lemon/clp.cc

endif

if HAVE_CBC

lemon_libemon_la_SOURCES += lemon/cbc.cc

endif

lemon_HEADERS += \

	lemon/adaptors.h \

	lemon/arg_parser.h \

	lemon/assert.h \

	lemon/bellman_ford.h \

	lemon/bfs.h \

	lemon/bin_heap.h \

	lemon/binomial_heap.h \

	lemon/bucket_heap.h \

	lemon/capacity_scaling.h \

	lemon/cbc.h \

	lemon/circulation.h \

	lemon/clp.h \

	lemon/color.h \

	lemon/concept_check.h \

	lemon/connectivity.h \

	lemon/core.h \

	lemon/cost_scaling.h \

	lemon/counter.h \

	lemon/cplex.h \

	lemon/cycle_canceling.h \

	lemon/dfs.h \

	lemon/dheap.h \

	lemon/dijkstra.h \

	lemon/dim2.h \

	lemon/dimacs.h \

	lemon/edge_set.h \

	lemon/elevator.h \

	lemon/error.h \

	lemon/euler.h \

	lemon/fib_heap.h \

	lemon/full_graph.h \

	lemon/glpk.h \

	lemon/gomory_hu.h \

	lemon/graph_to_eps.h \

	lemon/grid_graph.h \

	lemon/hypercube_graph.h \

	lemon/kruskal.h \

	lemon/hao_orlin.h \

	lemon/lgf_reader.h \

	lemon/lgf_writer.h \

	lemon/list_graph.h \

	lemon/lp.h \

	lemon/lp_base.h \

	lemon/lp_skeleton.h \

	lemon/maps.h \

	lemon/matching.h \

	lemon/math.h \

	lemon/min_cost_arborescence.h \

	lemon/nauty_reader.h \

	lemon/network_simplex.h \

	lemon/pairing_heap.h \

	lemon/path.h \

	lemon/planarity.h \

	lemon/preflow.h \

	lemon/quad_heap.h \

	lemon/radix_heap.h \

	lemon/radix_sort.h \

	lemon/random.h \

	lemon/smart_graph.h \

	lemon/soplex.h \

	lemon/static_graph.h \

	lemon/suurballe.h \

	lemon/time_measure.h \

	lemon/tolerance.h \

	lemon/unionfind.h \

	lemon/bits/windows.h

bits_HEADERS += \

	lemon/bits/alteration_notifier.h \

	lemon/bits/array_map.h \

	lemon/bits/bezier.h \

	lemon/bits/default_map.h \

	lemon/bits/edge_set_extender.h \

	lemon/bits/enable_if.h \

	lemon/bits/graph_adaptor_extender.h \

	lemon/bits/graph_extender.h \

	lemon/bits/map_extender.h \

	lemon/bits/path_dump.h \

	lemon/bits/solver_bits.h \

	lemon/bits/traits.h \

	lemon/bits/variant.h \

	lemon/bits/vector_map.h

concept_HEADERS += \

/* -*- C++ -*-

 *

 * This file is a part of LEMON, a generic C++ optimization library

 *

 * Copyright (C) 2003-2008

 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport

 * (Egervary Research Group on Combinatorial Optimization, EGRES).

 *

 * Permission to use, modify and distribute this software is granted

 * provided that this copyright notice appears in all copies. For

 * precise terms see the accompanying LICENSE file.

 *

 * This software is provided "AS IS" with no warranty of any kind,

 * express or implied, and with no claim as to its suitability for any

 * purpose.

 *

 */

#ifndef LEMON_CYCLE_CANCELING_H

#define LEMON_CYCLE_CANCELING_H

/// \ingroup min_cost_flow_algs

/// \file

/// \brief Cycle-canceling algorithms for finding a minimum cost flow.

#include <vector>

#include <limits>

#include <lemon/core.h>

#include <lemon/maps.h>

#include <lemon/path.h>

#include <lemon/math.h>

#include <lemon/static_graph.h>

#include <lemon/adaptors.h>

#include <lemon/circulation.h>

#include <lemon/bellman_ford.h>

namespace lemon {

  /// \addtogroup min_cost_flow_algs

  /// @{

  /// \brief Implementation of cycle-canceling algorithms for

  /// finding a \ref min_cost_flow "minimum cost flow".

  ///

  /// \ref CycleCanceling implements three different cycle-canceling

  /// algorithms for finding a \ref min_cost_flow "minimum cost flow"

  /// \ref amo93networkflows, \ref klein67primal,

  /// \ref goldberg89cyclecanceling.

  /// The most efficent one (both theoretically and practically)

  /// is the \ref CANCEL_AND_TIGHTEN "Cancel and Tighten" algorithm,

  /// thus it is the default method.

  /// It is strongly polynomial, but in practice, it is typically much

  /// slower than the scaling algorithms and NetworkSimplex.

  ///

  /// 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 V The number type used for flow amounts, capacity bounds

  /// and supply values in the algorithm. By default, it is \c int.

  /// \tparam C The number type used for costs and potentials in the

  /// algorithm. By default, it is the same as \c V.

  ///

  /// \warning Both number types must be signed and all input data must

  /// be integer.

  /// \warning This algorithm does not support negative costs for such

  /// arcs that have infinite upper bound.

  ///

  /// \note For more information about the three available methods,

  /// see \ref Method.

#ifdef DOXYGEN

  template <typename GR, typename V, typename C>

#else

  template <typename GR, typename V = int, typename C = V>

#endif

  class CycleCanceling

  {

  public:

    /// The type of the digraph

    typedef GR Digraph;

          } else {

            // Search for node disjoint negative cycles

            std::vector<int> state(_res_node_num, 0);

            int id = 0;

            for (int u = 0; u != _res_node_num; ++u) {

              if (state[u] != 0) continue;

              ++id;

              int v = u;

              for (; v != -1 && state[v] == 0; v = pred[v] == INVALID ?

                   -1 : rgr.id(rgr.source(pred[v]))) {

                state[v] = id;

              }

              if (v != -1 && state[v] == id) {

                // A negative cycle is found

                cycle_found = true;

                cycle.clear();

                StaticDigraph::Arc a = pred[v];

                Value d, delta = _res_cap[rgr.id(a)];

                cycle.push_back(rgr.id(a));

                while (rgr.id(rgr.source(a)) != v) {

                  a = pred_map[rgr.source(a)];

                  d = _res_cap[rgr.id(a)];

                  if (d < delta) delta = d;

                  cycle.push_back(rgr.id(a));

                }

                // Augment along the cycle

                for (int i = 0; i < int(cycle.size()); ++i) {

                  int j = cycle[i];

                  _res_cap[j] -= delta;

                  _res_cap[_reverse[j]] += delta;

                }

              }

            }

          }

          // Increase iteration limit if no cycle is found

          if (!cycle_found) {

            length_bound = static_cast<int>(length_bound * BF_LIMIT_FACTOR);

          }

        }

      }

    }

    // Execute the "Minimum Mean Cycle Canceling" method

    void startMinMeanCycleCanceling() {

      typedef SimplePath<StaticDigraph> SPath;

      typedef typename SPath::ArcIt SPathArcIt;

        ::template SetPath<SPath>::Create MMC;

      SPath cycle;

      MMC mmc(_sgr, _cost_map);

      mmc.cycle(cycle);

      buildResidualNetwork();

        // Find the cycle

        mmc.findCycle();

        // Compute delta value

        Value delta = INF;

        for (SPathArcIt a(cycle); a != INVALID; ++a) {

          Value d = _res_cap[_id_vec[_sgr.id(a)]];

          if (d < delta) delta = d;

        }

        // Augment along the cycle

        for (SPathArcIt a(cycle); a != INVALID; ++a) {

          int j = _id_vec[_sgr.id(a)];

          _res_cap[j] -= delta;

          _res_cap[_reverse[j]] += delta;

        }

        // Rebuild the residual network        

        buildResidualNetwork();

      }

    }

    // Execute the "Cancel And Tighten" method

    void startCancelAndTighten() {

      // Constants for the min mean cycle computations

      const double LIMIT_FACTOR = 1.0;

      const int MIN_LIMIT = 5;

      // Contruct auxiliary data vectors

      DoubleVector pi(_res_node_num, 0.0);

      IntVector level(_res_node_num);

      BoolVector reached(_res_node_num);

      BoolVector processed(_res_node_num);

      IntVector pred_node(_res_node_num);

      IntVector pred_arc(_res_node_num);

      std::vector<int> stack(_res_node_num);

      std::vector<int> proc_vector(_res_node_num);

      // Initialize epsilon

      double epsilon = 0;

      for (int a = 0; a != _res_arc_num; ++a) {

        if (_res_cap[a] > 0 && -_cost[a] > epsilon)

          epsilon = -_cost[a];

      }

      // Start phases

      Tolerance<double> tol;

      tol.epsilon(1e-6);

        if (--iter > 0) {

          for (int u = 0; u != _res_node_num; ++u) {

            level[u] = 0;

          }

          for (int i = proc_head; i > 0; --i) {

            int u = proc_vector[i];

            double p = pi[u];

            int l = level[u] + 1;

            int last_out = _first_out[u+1];

            for (int a = _first_out[u]; a != last_out; ++a) {

              int v = _target[a];

              if (_res_cap[a] > 0 && tol.negative(_cost[a] + p - pi[v]) &&

                  l > level[v]) level[v] = l;

            }

          }

          // Modify potentials

          double q = std::numeric_limits<double>::max();

          for (int u = 0; u != _res_node_num; ++u) {

            int lu = level[u];

            double p, pu = pi[u];

            int last_out = _first_out[u+1];

            for (int a = _first_out[u]; a != last_out; ++a) {

              if (_res_cap[a] == 0) continue;

              int v = _target[a];

              int ld = lu - level[v];

              if (ld > 0) {

                p = (_cost[a] + pu - pi[v] + epsilon) / (ld + 1);

                if (p < q) q = p;

              }

            }

          }

          for (int u = 0; u != _res_node_num; ++u) {

            pi[u] -= q * level[u];

          }

          // Modify epsilon

          epsilon = 0;

          for (int u = 0; u != _res_node_num; ++u) {

            double curr, pu = pi[u];

            int last_out = _first_out[u+1];

            for (int a = _first_out[u]; a != last_out; ++a) {

              if (_res_cap[a] == 0) continue;

              curr = _cost[a] + pu - pi[_target[a]];

              if (-curr > epsilon) epsilon = -curr;

            }

          }

        } else {

          typedef typename BellmanFord<StaticDigraph, CostArcMap>

            ::template SetDistMap<CostNodeMap>::Create BF;

          // Set epsilon to the minimum cycle mean

          buildResidualNetwork();

          MMC mmc(_sgr, _cost_map);

          epsilon = -mmc.cycleMean();

          // Compute feasible potentials for the current epsilon

          for (int i = 0; i != int(_cost_vec.size()); ++i) {

            _cost_vec[i] = cycle_size * _cost_vec[i] - cycle_cost;

          }

          BF bf(_sgr, _cost_map);

          bf.distMap(_pi_map);

          bf.init(0);

          bf.start();

          for (int u = 0; u != _res_node_num; ++u) {

            pi[u] = static_cast<double>(_pi[u]) / cycle_size;

          }

          iter = limit;

        }

      }

    }

  }; //class CycleCanceling

  ///@}

} //namespace lemon

#endif //LEMON_CYCLE_CANCELING_H

/* -*- C++ -*-

 *

 * This file is a part of LEMON, a generic C++ optimization library

 *

 * Copyright (C) 2003-2008

 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport

 * (Egervary Research Group on Combinatorial Optimization, EGRES).

 *

 * Permission to use, modify and distribute this software is granted

 * provided that this copyright notice appears in all copies. For

 * precise terms see the accompanying LICENSE file.

 *

 * This software is provided "AS IS" with no warranty of any kind,

 * express or implied, and with no claim as to its suitability for any

 * purpose.

 *

 */

/// \ingroup min_mean_cycle

///

/// \file

/// \brief Hartmann-Orlin's algorithm for finding a minimum mean cycle.

#include <vector>

#include <limits>

#include <lemon/core.h>

#include <lemon/path.h>

#include <lemon/tolerance.h>

#include <lemon/connectivity.h>

namespace lemon {

  ///

  /// \tparam GR The type of the digraph.

  /// It must conform to the \ref concepts::Rea_data "Rea_data" concept.

#ifdef DOXYGEN

#else

#endif

  {

    /// The type of the digraph

    typedef GR Digraph;

    ///

    /// otherwise it is \c double.

    /// The tolerance type used for internal computations

    /// \brief The path type of the found cycles

    ///

    /// The path type of the found cycles.

    /// It must conform to the \ref lemon::concepts::Path "Path" concept

    /// and it must have an \c addFront() function.

    typedef lemon::Path<Digraph> Path;

  };

  {

    typedef GR Digraph;

#ifdef LEMON_HAVE_LONG_LONG

#else

#endif

    typedef lemon::Path<Digraph> Path;

  };

  /// \addtogroup min_mean_cycle

  /// @{

  /// \brief Implementation of the Hartmann-Orlin algorithm for finding

  /// a minimum mean cycle.

  ///

  /// This class implements the Hartmann-Orlin algorithm for finding

  /// \ref amo93networkflows, \ref dasdan98minmeancycle.

  /// It is an improved version of \ref Karp "Karp"'s original algorithm,

  /// it applies an efficient early termination scheme.

  /// It runs in time O(ne) and uses space O(n<sup>2</sup>+e).

  ///

  /// \tparam GR The type of the digraph the algorithm runs on.

  /// map type is \ref concepts::Digraph::ArcMap "GR::ArcMap<int>".

  /// \tparam TR The traits class that defines various types used by the

  /// In most cases, this parameter should not be set directly,

  /// consider to use the named template parameters instead.

#ifdef DOXYGEN

#else

  template < typename GR,

#endif

  {

  public:

    /// The type of the digraph

    typedef typename TR::Digraph Digraph;

    ///

    /// otherwise it is \c double.

    /// The tolerance type

    typedef typename TR::Tolerance Tolerance;

    /// \brief The path type of the found cycles

    ///

    /// The path type of the found cycles.

    /// it is \ref lemon::Path "Path<Digraph>".

    typedef typename TR::Path Path;

    typedef TR Traits;

  private:

    TEMPLATE_DIGRAPH_TYPEDEFS(Digraph);

    // Data sturcture for path data

    struct PathData

    {

      Arc pred;

        dist(d), pred(p) {}

    };

    typedef typename Digraph::template NodeMap<std::vector<PathData> >

      PathDataNodeMap;

  private:

    // The digraph the algorithm runs on

    const Digraph &_gr;

    // Data for storing the strongly connected components

    int _comp_num;

    typename Digraph::template NodeMap<int> _comp;

    std::vector<std::vector<Node> > _comp_nodes;

    std::vector<Node>* _nodes;

    typename Digraph::template NodeMap<std::vector<Arc> > _out_arcs;

    // Data for the found cycles

    bool _curr_found, _best_found;

    int _curr_size, _best_size;

    Node _curr_node, _best_node;

    int _curr_level, _best_level;

    Path *_cycle_path;

    bool _local_path;

    // Node map for storing path data

    PathDataNodeMap _data;

    // The processed nodes in the last round

    std::vector<Node> _process;

    Tolerance _tolerance;

    // Infinite constant

  public:

    /// \name Named Template Parameters

    /// @{

    template <typename T>

      typedef lemon::Tolerance<T> Tolerance;

    };

    /// \brief \ref named-templ-param "Named parameter" for setting

    ///

    /// type. It is used for internal computations in the algorithm.

    template <typename T>

    };

    template <typename T>

    struct SetPathTraits : public Traits {

      typedef T Path;

    };

    /// \brief \ref named-templ-param "Named parameter" for setting

    /// \c %Path type.

    ///

    /// \ref named-templ-param "Named parameter" for setting the \c %Path

    /// type of the found cycles.

    /// It must conform to the \ref lemon::concepts::Path "Path" concept

    /// and it must have an \c addFront() function.

    template <typename T>

    struct SetPath

    };

    /// @}

  protected:

  public:

    /// \brief Constructor.

    ///

    /// The constructor of the class.

    ///

    /// \param digraph The digraph the algorithm runs on.

      _cycle_path(NULL), _local_path(false), _data(digraph),

    {}

    /// Destructor.

      if (_local_path) delete _cycle_path;

    }

    /// \brief Set the path structure for storing the found cycle.

    ///

    /// This function sets an external path structure for storing the

    /// found cycle.

    ///

    /// If you don't call this function before calling \ref run() or

    /// structure. The destuctor deallocates this automatically

    /// allocated object, of course.

    ///

    /// \note The algorithm calls only the \ref lemon::Path::addFront()

    /// "addFront()" function of the given path structure.

    ///

    /// \return <tt>(*this)</tt>

      if (_local_path) {

        delete _cycle_path;

        _local_path = false;

      }

      _cycle_path = &path;

      return *this;

    }

    /// \brief Set the tolerance used by the algorithm.

    ///

    /// This function sets the tolerance object used by the algorithm.

    ///

    /// \return <tt>(*this)</tt>

      _tolerance = tolerance;

      return *this;

    }

    /// \brief Return a const reference to the tolerance.

    ///

    /// This function returns a const reference to the tolerance object

    /// used by the algorithm.

    const Tolerance& tolerance() const {

      return _tolerance;

    }

    /// \name Execution control

    /// The simplest way to execute the algorithm is to call the \ref run()

    /// function.\n

    /// @{

    /// \brief Run the algorithm.

    ///

    /// This function runs the algorithm.

    /// It can be called more than once (e.g. if the underlying digraph

    ///

    /// \return \c true if a directed cycle exists in the digraph.

    ///

    /// \note <tt>mmc.run()</tt> is just a shortcut of the following code.

    /// \code

    /// \endcode

    bool run() {

    }

    /// \brief Find the minimum cycle mean.

    ///

    /// cycles in the digraph.

    ///

    /// \return \c true if a directed cycle exists in the digraph.

      // Initialization and find strongly connected components

      init();

      findComponents();

      // Find the minimum cycle mean in the components

      for (int comp = 0; comp < _comp_num; ++comp) {

        if (!initComponent(comp)) continue;

        processRounds();

        // Update the best cycle (global minimum mean cycle)

        if ( _curr_found && (!_best_found || 

          _best_found = true;

          _best_size = _curr_size;

          _best_node = _curr_node;

          _best_level = _curr_level;

        }

      }

      return _best_found;

    }

    /// \brief Find a minimum mean directed cycle.

    ///

    ///

    /// \return \c true if a directed cycle exists in the digraph.

    ///

    bool findCycle() {

      if (!_best_found) return false;

      IntNodeMap reached(_gr, -1);

      int r = _best_level + 1;

      Node u = _best_node;

      while (reached[u] < 0) {

        reached[u] = --r;

        u = _gr.source(_data[u][r].pred);

      }

      r = reached[u];

      Arc e = _data[u][r].pred;

      _cycle_path->addFront(e);

      _best_size = 1;

      Node v;

      while ((v = _gr.source(e)) != u) {

        e = _data[v][--r].pred;

        _cycle_path->addFront(e);

        ++_best_size;

      }

      return true;

    }

    /// @}

    /// \name Query Functions

    /// The results of the algorithm can be obtained using these

    /// functions.\n

    /// The algorithm should be executed before using them.

    /// @{

    ///

    ///

    /// using this function.

    }

    /// \brief Return the number of arcs on the found cycle.

    ///

    /// This function returns the number of arcs on the found cycle.

    ///

    /// using this function.

      return _best_size;

    }

    ///

    ///

    /// \note <tt>alg.cycleMean()</tt> is just a shortcut of the

    /// following code.

    /// \code

    /// \endcode

    ///

    /// using this function.

    double cycleMean() const {

    }

    /// \brief Return the found cycle.

    ///

    /// This function returns a const reference to the path structure

    /// storing the found cycle.

    ///

    /// \pre \ref run() or \ref findCycle() must be called before using

    /// this function.

    const Path& cycle() const {

      return *_cycle_path;

    }

    ///@}

  private:

    // Initialization

    void init() {

      if (!_cycle_path) {

        _local_path = true;

        _cycle_path = new Path;

      }

      _cycle_path->clear();

      _best_found = false;

      _best_size = 1;

      _cycle_path->clear();

      for (NodeIt u(_gr); u != INVALID; ++u)

        _data[u].clear();

    }

    // Find strongly connected components and initialize _comp_nodes

    // and _out_arcs

    void findComponents() {

      _comp_num = stronglyConnectedComponents(_gr, _comp);

      _comp_nodes.resize(_comp_num);

      if (_comp_num == 1) {

        _comp_nodes[0].clear();

        for (NodeIt n(_gr); n != INVALID; ++n) {

          _comp_nodes[0].push_back(n);

          _out_arcs[n].clear();

          for (OutArcIt a(_gr, n); a != INVALID; ++a) {

            _out_arcs[n].push_back(a);

          }

        }

      } else {

        for (int i = 0; i < _comp_num; ++i)

          _comp_nodes[i].clear();

        for (NodeIt n(_gr); n != INVALID; ++n) {

          int k = _comp[n];

          _comp_nodes[k].push_back(n);

          _out_arcs[n].clear();

          for (OutArcIt a(_gr, n); a != INVALID; ++a) {

            if (_comp[_gr.target(a)] == k) _out_arcs[n].push_back(a);

          }

        }

      }

    }

    // Initialize path data for the current component

    bool initComponent(int comp) {

      _nodes = &(_comp_nodes[comp]);

      int n = _nodes->size();

      if (n < 1 || (n == 1 && _out_arcs[(*_nodes)[0]].size() == 0)) {

        return false;

      }      

      for (int i = 0; i < n; ++i) {

        _data[(*_nodes)[i]].resize(n + 1, PathData(INF));

      }

      return true;

    }

    // Process all rounds of computing path data for the current component.

    // node to node v containing exactly k arcs.

    void processRounds() {

      Node start = (*_nodes)[0];

      _data[start][0] = PathData(0);

      _process.clear();

      _process.push_back(start);

      int k, n = _nodes->size();

      int next_check = 4;

      bool terminate = false;

      for (k = 1; k <= n && int(_process.size()) < n && !terminate; ++k) {

        processNextBuildRound(k);

        if (k == next_check || k == n) {

          terminate = checkTermination(k);

          next_check = next_check * 3 / 2;

        }

      }

      for ( ; k <= n && !terminate; ++k) {

        processNextFullRound(k);

        if (k == next_check || k == n) {

          terminate = checkTermination(k);

          next_check = next_check * 3 / 2;

        }

      }

    }

    // Process one round and rebuild _process

    void processNextBuildRound(int k) {

      std::vector<Node> next;

      Node u, v;

      Arc e;

      for (int i = 0; i < int(_process.size()); ++i) {

        u = _process[i];

        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {

          e = _out_arcs[u][j];

          v = _gr.target(e);

          if (_tolerance.less(d, _data[v][k].dist)) {

            if (_data[v][k].dist == INF) next.push_back(v);

            _data[v][k] = PathData(d, e);

          }

        }

      }

      _process.swap(next);

    }

    // Process one round using _nodes instead of _process

    void processNextFullRound(int k) {

      Node u, v;

      Arc e;

      for (int i = 0; i < int(_nodes->size()); ++i) {

        u = (*_nodes)[i];

        for (int j = 0; j < int(_out_arcs[u].size()); ++j) {

          e = _out_arcs[u][j];

          v = _gr.target(e);

          if (_tolerance.less(d, _data[v][k].dist)) {

            _data[v][k] = PathData(d, e);

          }

        }

      }

    }

    // Check early termination

    bool checkTermination(int k) {

      typedef std::pair<int, int> Pair;

      typename GR::template NodeMap<Pair> level(_gr, Pair(-1, 0));

      int n = _nodes->size();