RhodeCode

  /// the network reliability, especially to test how many links have

  /// to be destroyed in the network to split it to at least two

  /// distinict subnetworks.

  ///

  /// The complexity of the algorithm is \f$ O(nm\log(n)) \f$ but with

  /// \ref FibHeap "Fibonacci heap" it can be decreased to

  /// \f$ O(nm+n^2\log(n)) \f$.  When the edges have unit capacities,

  /// \c BucketHeap can be used which yields \f$ O(nm) \f$ time

  /// complexity.

  ///

  /// \warning The value type of the capacity map should be able to

  /// hold any cut value of the graph, otherwise the result can

  /// overflow.

  /// \note This capacity is supposed to be integer type.

#ifdef DOXYGEN

  template <typename GR, typename CM, typename TR>

#else

  template <typename GR,

            typename CM = typename GR::template EdgeMap<int>,

            typename TR = NagamochiIbarakiDefaultTraits<GR, CM> >

#endif

  class NagamochiIbaraki {

  public:

    typedef TR Traits;

    /// The type of the underlying graph.

    typedef typename Traits::Graph Graph;

    /// The type of the capacity map.

    typedef typename Traits::CapacityMap CapacityMap;

    /// The value type of the capacity map.

    typedef typename Traits::CapacityMap::Value Value;

    /// The heap type used by the algorithm.

    typedef typename Traits::Heap Heap;

    /// The cross reference type used for the heap.

    typedef typename Traits::HeapCrossRef HeapCrossRef;

    ///\name Named template parameters

    ///@{

    struct SetUnitCapacityTraits : public Traits {

      typedef ConstMap<typename Graph::Edge, Const<int, 1> > CapacityMap;

      static CapacityMap *createCapacityMap(const Graph&) {

        return new CapacityMap();

      }

    };

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

    /// the capacity map to a constMap<Edge, int, 1>() instance

    ///

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

    /// the capacity map to a constMap<Edge, int, 1>() instance

    struct SetUnitCapacity

      : public NagamochiIbaraki<Graph, CapacityMap,

                                SetUnitCapacityTraits> {

      typedef NagamochiIbaraki<Graph, CapacityMap,

                               SetUnitCapacityTraits> Create;

    };

    template <class H, class CR>

    struct SetHeapTraits : public Traits {

      typedef CR HeapCrossRef;

      typedef H Heap;

      static HeapCrossRef *createHeapCrossRef(int num) {

        LEMON_ASSERT(false, "HeapCrossRef is not initialized");

        return 0; // ignore warnings

      }

      static Heap *createHeap(HeapCrossRef &) {

        LEMON_ASSERT(false, "Heap is not initialized");

        return 0; // ignore warnings

      }

    };

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

    /// heap and cross reference type

    ///

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

    /// cross reference type. The heap has to maximize the priorities.

    template <class H, class CR = RangeMap<int> >

    struct SetHeap

      : public NagamochiIbaraki<Graph, CapacityMap, SetHeapTraits<H, CR> > {

      typedef NagamochiIbaraki< Graph, CapacityMap, SetHeapTraits<H, CR> >

      Create;

    };

    template <class H, class CR>

    struct SetStandardHeapTraits : public Traits {

      typedef CR HeapCrossRef;

      typedef H Heap;

      static HeapCrossRef *createHeapCrossRef(int size) {

        return new HeapCrossRef(size);

      }

      static Heap *createHeap(HeapCrossRef &crossref) {

        return new Heap(crossref);

      }

    };

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

    /// heap and cross reference type with automatic allocation

    ///

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

    /// cross reference type with automatic allocation. They should

    /// have standard constructor interfaces to be able to

    /// automatically created by the algorithm (i.e. the graph should

    /// be passed to the constructor of the cross reference and the

    /// cross reference should be passed to the constructor of the

    /// heap). However, external heap and cross reference objects

    /// could also be passed to the algorithm using the \ref heap()

    /// function before calling \ref run() or \ref init(). The heap

    /// has to maximize the priorities.

    /// \sa SetHeap

    template <class H, class CR = RangeMap<int> >

    struct SetStandardHeap

      : public NagamochiIbaraki<Graph, CapacityMap,

                                SetStandardHeapTraits<H, CR> > {

      typedef NagamochiIbaraki<Graph, CapacityMap,

                               SetStandardHeapTraits<H, CR> > Create;

    };

    ///@}

  private:

    const Graph &_graph;

    const CapacityMap *_capacity;

    bool _local_capacity; // unit capacity

    struct ArcData {

      typename Graph::Node target;

      int prev, next;

    };

    struct EdgeData {

      Value capacity;

      Value cut;

    };

    struct NodeData {

      int first_arc;

      typename Graph::Node prev, next;

      int curr_arc;

      typename Graph::Node last_rep;

      Value sum;

    };

    typename Graph::template NodeMap<NodeData> *_nodes;

    std::vector<ArcData> _arcs;

    std::vector<EdgeData> _edges;

    typename Graph::Node _first_node;

    int _node_num;

    Value _min_cut;

    HeapCrossRef *_heap_cross_ref;

    bool _local_heap_cross_ref;

    Heap *_heap;

    bool _local_heap;

    typedef typename Graph::template NodeMap<typename Graph::Node> NodeList;

    NodeList *_next_rep;

    typedef typename Graph::template NodeMap<bool> MinCutMap;

    MinCutMap *_cut_map;

    void createStructures() {

      if (!_nodes) {

        _nodes = new (typename Graph::template NodeMap<NodeData>)(_graph);

      }

      if (!_capacity) {

        _local_capacity = true;

        _capacity = Traits::createCapacityMap(_graph);

      }

      if (!_heap_cross_ref) {

        _local_heap_cross_ref = true;

        _heap_cross_ref = Traits::createHeapCrossRef(_graph);

      }

      if (!_heap) {

        _local_heap = true;

        _heap = Traits::createHeap(*_heap_cross_ref);

      }

      if (!_next_rep) {

        _next_rep = new NodeList(_graph);

      }

      if (!_cut_map) {

        _cut_map = new MinCutMap(_graph);

      }

    }

    typedef NagamochiIbaraki Create;

    /// \brief Constructor.

    ///

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

    /// \param capacity The capacity map used by the algorithm.

    NagamochiIbaraki(const Graph& graph, const CapacityMap& capacity)

      : _graph(graph), _capacity(&capacity), _local_capacity(false),

        _nodes(0), _arcs(), _edges(), _min_cut(),

        _heap_cross_ref(0), _local_heap_cross_ref(false),

        _heap(0), _local_heap(false),

        _next_rep(0), _cut_map(0) {}

    /// \brief Constructor.

    ///

    /// This constructor can be used only when the Traits class

    /// defines how can the local capacity map be instantiated.

    /// If the SetUnitCapacity used the algorithm automatically

    /// constructs the capacity map.

    ///

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

    NagamochiIbaraki(const Graph& graph)

      : _graph(graph), _capacity(0), _local_capacity(false),

        _nodes(0), _arcs(), _edges(), _min_cut(),

        _heap_cross_ref(0), _local_heap_cross_ref(false),

        _heap(0), _local_heap(false),

        _next_rep(0), _cut_map(0) {}

    /// \brief Destructor.

    ///

    /// Destructor.

    ~NagamochiIbaraki() {

      if (_local_capacity) delete _capacity;

      if (_nodes) delete _nodes;

      if (_local_heap) delete _heap;

      if (_local_heap_cross_ref) delete _heap_cross_ref;

      if (_next_rep) delete _next_rep;

      if (_cut_map) delete _cut_map;

    }

    /// \brief Sets the heap and the cross reference used by algorithm.

    ///

    /// Sets the heap and the cross reference used by algorithm.

    /// If you don't use this function before calling \ref run(),

    /// it will allocate one. The destuctor deallocates this

    /// automatically allocated heap and cross reference, of course.

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

    NagamochiIbaraki &heap(Heap& hp, HeapCrossRef &cr)

    {

      if (_local_heap_cross_ref) {

        delete _heap_cross_ref;

        _local_heap_cross_ref = false;

      }

      _heap_cross_ref = &cr;

      if (_local_heap) {

        delete _heap;

        _local_heap = false;

      }

      _heap = &hp;

      return *this;

    }

    /// \name Execution control

    /// The simplest way to execute the algorithm is to use

    /// one of the member functions called \c run().

    /// \n

    /// If you need more control on the execution,

    /// first you must call \ref init() and then call the start()

    /// or proper times the processNextPhase() member functions.

    ///@{

    /// \brief Initializes the internal data structures.

    ///

    /// Initializes the internal data structures.

    void init() {

      createStructures();

      int edge_num = countEdges(_graph);

      _edges.resize(edge_num);

      _arcs.resize(2 * edge_num);

      typename Graph::Node prev = INVALID;

      _node_num = 0;

      for (typename Graph::NodeIt n(_graph); n != INVALID; ++n) {

        (*_cut_map)[n] = false;

        (*_next_rep)[n] = INVALID;

        (*_nodes)[n].last_rep = n;

        (*_nodes)[n].first_arc = -1;

        (*_nodes)[n].curr_arc = -1;

        (*_nodes)[n].prev = prev;

        if (prev != INVALID) {

          (*_nodes)[prev].next = n;

        }

        (*_nodes)[n].next = INVALID;

        (*_nodes)[n].sum = 0;

        prev = n;

        ++_node_num;

      }

      _first_node = typename Graph::NodeIt(_graph);

      int index = 0;

      for (typename Graph::NodeIt n(_graph); n != INVALID; ++n) {

        for (typename Graph::OutArcIt a(_graph, n); a != INVALID; ++a) {

          typename Graph::Node m = _graph.target(a);

          if (!(n < m)) continue;

          (*_nodes)[n].sum += (*_capacity)[a];

          (*_nodes)[m].sum += (*_capacity)[a];

          int c = (*_nodes)[m].curr_arc;

          if (c != -1 && _arcs[c ^ 1].target == n) {

            _edges[c >> 1].capacity += (*_capacity)[a];

          } else {

            _edges[index].capacity = (*_capacity)[a];

            _arcs[index << 1].prev = -1;

            if ((*_nodes)[n].first_arc != -1) {

              _arcs[(*_nodes)[n].first_arc].prev = (index << 1);

            }

            _arcs[index << 1].next = (*_nodes)[n].first_arc;

            (*_nodes)[n].first_arc = (index << 1);

            _arcs[index << 1].target = m;

            (*_nodes)[m].curr_arc = (index << 1);

            _arcs[(index << 1) | 1].prev = -1;

            if ((*_nodes)[m].first_arc != -1) {

              _arcs[(*_nodes)[m].first_arc].prev = ((index << 1) | 1);

            }

            _arcs[(index << 1) | 1].next = (*_nodes)[m].first_arc;

            (*_nodes)[m].first_arc = ((index << 1) | 1);

            _arcs[(index << 1) | 1].target = n;

            ++index;

          }

        }

      }

      typename Graph::Node cut_node = INVALID;

      _min_cut = std::numeric_limits<Value>::max();

      for (typename Graph::Node n = _first_node; 

           n != INVALID; n = (*_nodes)[n].next) {

        if ((*_nodes)[n].sum < _min_cut) {

          cut_node = n;

          _min_cut = (*_nodes)[n].sum;

        }

      }

      (*_cut_map)[cut_node] = true;

      if (_min_cut == 0) {

        _first_node = INVALID;

      }

    }

  public:

    /// \brief Processes the next phase

    ///

    /// Processes the next phase in the algorithm. It must be called

    /// at most one less the number of the nodes in the graph.

    ///

    ///\return %True when the algorithm finished.

    bool processNextPhase() {

      if (_first_node == INVALID) return true;

      _heap->clear();

      for (typename Graph::Node n = _first_node; 

           n != INVALID; n = (*_nodes)[n].next) {

        (*_heap_cross_ref)[n] = Heap::PRE_HEAP;

      }

      std::vector<typename Graph::Node> order;

      order.reserve(_node_num);

      int sep = 0;

      Value alpha = 0;

      Value pmc = std::numeric_limits<Value>::max();

      _heap->push(_first_node, static_cast<Value>(0));

      while (!_heap->empty()) {

        typename Graph::Node n = _heap->top();

        Value v = _heap->prio();

        _heap->pop();

        for (int a = (*_nodes)[n].first_arc; a != -1; a = _arcs[a].next) {

          switch (_heap->state(_arcs[a].target)) {

          case Heap::PRE_HEAP: