[Lemon-commits] [lemon_svn] deba: r2664 - in hugo/trunk: lemon test
Lemon SVN
svn at lemon.cs.elte.hu
Mon Nov 6 20:54:19 CET 2006
Author: deba
Date: Fri Mar 31 13:10:44 2006
New Revision: 2664
Added:
hugo/trunk/test/arborescence_test.cc
Modified:
hugo/trunk/lemon/min_cost_arborescence.h
hugo/trunk/test/Makefile.am
Log:
Bugfix in the minimum cost arborescence algorithm
Dual solution computation and interface for algorithm
Optimality test on random graph
Modified: hugo/trunk/lemon/min_cost_arborescence.h
==============================================================================
--- hugo/trunk/lemon/min_cost_arborescence.h (original)
+++ hugo/trunk/lemon/min_cost_arborescence.h Fri Mar 31 13:10:44 2006
@@ -26,6 +26,7 @@
#include <vector>
#include <lemon/list_graph.h>
+#include <lemon/bin_heap.h>
namespace lemon {
@@ -56,11 +57,9 @@
/// in the arborescence.
///
/// The type of the map that stores which edges are in the arborescence.
- /// It must meet the \ref concept::ReadWriteMap "ReadWriteMap" concept.
- /// Initially it will be setted to false on each edge. The algorithm
- /// may set each value one time to true and maybe after it to false again.
- /// Therefore you cannot use maps like BackInserteBoolMap with this
- /// algorithm.
+ /// It must meet the \ref concept::WritedMap "WriteMap" concept.
+ /// Initially it will be setted to false on each edge. After it
+ /// will set all arborescence edges once.
typedef typename Graph::template EdgeMap<bool> ArborescenceMap;
/// \brief Instantiates a ArborescenceMap.
@@ -72,6 +71,20 @@
return new ArborescenceMap(_graph);
}
+ /// \brief The type of the PredMap
+ ///
+ /// The type of the PredMap. It is a node map with an edge value type.
+ typedef typename Graph::template NodeMap<typename Graph::Edge> PredMap;
+
+ /// \brief Instantiates a PredMap.
+ ///
+ /// This function instantiates a \ref PredMap.
+ /// \param _graph is the graph, to which we would like to define the
+ /// PredMap.
+ static PredMap *createPredMap(const Graph &_graph){
+ return new PredMap(_graph);
+ }
+
};
/// \ingroup spantree
@@ -86,6 +99,9 @@
/// given sources and spans all the nodes which are reachable from the
/// sources. The time complexity of the algorithm is O(n^2 + e).
///
+ /// The algorithm provides also an optimal dual solution to arborescence
+ /// that way the optimality of the algorithm can be proofed easily.
+ ///
/// \param _Graph The graph type the algorithm runs on. The default value
/// is \ref ListGraph. The value of _Graph is not used directly by
/// MinCostArborescence, it is only passed to
@@ -135,6 +151,8 @@
typedef typename Traits::CostMap CostMap;
///The type of the costs of the edges.
typedef typename Traits::Value Value;
+ ///The type of the predecessor map.
+ typedef typename Traits::PredMap PredMap;
///The type of the map that stores which edges are in the arborescence.
typedef typename Traits::ArborescenceMap ArborescenceMap;
@@ -157,14 +175,20 @@
};
- const Graph* graph;
- const CostMap* cost;
+ const Graph *graph;
+ const CostMap *cost;
- ArborescenceMap* _arborescence_map;
- bool local_arborescence_map;
+ PredMap *_pred;
+ bool local_pred;
- typedef typename Graph::template NodeMap<int> LevelMap;
- LevelMap *_level;
+ ArborescenceMap *_arborescence;
+ bool local_arborescence;
+
+ typedef typename Graph::template EdgeMap<int> EdgeOrder;
+ EdgeOrder *_edge_order;
+
+ typedef typename Graph::template NodeMap<int> NodeOrder;
+ NodeOrder *_node_order;
typedef typename Graph::template NodeMap<CostEdge> CostEdgeMap;
CostEdgeMap *_cost_edges;
@@ -179,7 +203,30 @@
std::vector<StackLevel> level_stack;
std::vector<Node> queue;
- int node_counter;
+ typedef std::vector<typename Graph::Node> DualNodeList;
+
+ DualNodeList _dual_node_list;
+
+ struct DualVariable {
+ int begin, end;
+ Value value;
+
+ DualVariable(int _begin, int _end, Value _value)
+ : begin(_begin), end(_end), value(_value) {}
+
+ };
+
+ typedef std::vector<DualVariable> DualVariables;
+
+ DualVariables _dual_variables;
+
+ typedef typename Graph::template NodeMap<int> HeapCrossRef;
+
+ HeapCrossRef *_heap_cross_ref;
+
+ typedef BinHeap<Node, int, HeapCrossRef> Heap;
+
+ Heap *_heap;
public:
@@ -208,6 +255,26 @@
typedef MinCostArborescence<Graph, CostMap,
DefArborescenceMapTraits<T> > Create;
};
+
+ template <class T>
+ struct DefPredMapTraits : public Traits {
+ typedef T PredMap;
+ static PredMap *createPredMap(const Graph &)
+ {
+ throw UninitializedParameter();
+ }
+ };
+
+ /// \brief \ref named-templ-param "Named parameter" for
+ /// setting PredMap type
+ ///
+ /// \ref named-templ-param "Named parameter" for setting
+ /// PredMap type
+ template <class T>
+ struct DefPredMap
+ : public MinCostArborescence<Graph, CostMap, DefPredMapTraits<T> > {
+ typedef MinCostArborescence<Graph, CostMap, DefPredMapTraits<T> > Create;
+ };
/// @}
@@ -216,9 +283,10 @@
/// \param _graph The graph the algorithm will run on.
/// \param _cost The cost map used by the algorithm.
MinCostArborescence(const Graph& _graph, const CostMap& _cost)
- : graph(&_graph), cost(&_cost),
- _arborescence_map(0), local_arborescence_map(false),
- _level(0), _cost_edges(0) {}
+ : graph(&_graph), cost(&_cost), _pred(0), local_pred(false),
+ _arborescence(0), local_arborescence(false),
+ _edge_order(0), _node_order(0), _cost_edges(0),
+ _heap_cross_ref(0), _heap(0) {}
/// \brief Destructor.
~MinCostArborescence() {
@@ -230,7 +298,24 @@
/// Sets the arborescence map.
/// \return \c (*this)
MinCostArborescence& arborescenceMap(ArborescenceMap& m) {
- _arborescence_map = &m;
+ if (local_arborescence) {
+ delete _arborescence;
+ }
+ local_arborescence = false;
+ _arborescence = &m;
+ return *this;
+ }
+
+ /// \brief Sets the arborescence map.
+ ///
+ /// Sets the arborescence map.
+ /// \return \c (*this)
+ MinCostArborescence& predMap(PredMap& m) {
+ if (local_pred) {
+ delete _pred;
+ }
+ local_pred = false;
+ _pred = &m;
return *this;
}
@@ -246,7 +331,7 @@
///
/// Returns a reference to the arborescence map.
const ArborescenceMap& arborescenceMap() const {
- return *_arborescence_map;
+ return *_arborescence;
}
/// \brief Returns true if the edge is in the arborescence.
@@ -254,23 +339,141 @@
/// Returns true if the edge is in the arborescence.
/// \param edge The edge of the graph.
/// \pre \ref run() must be called before using this function.
- bool arborescenceEdge(Edge edge) const {
- return (*_arborescence_map)[edge];
+ bool arborescence(Edge edge) const {
+ return (*_pred)[graph->target(edge)] == edge;
+ }
+
+ /// \brief Returns a reference to the pred map.
+ ///
+ /// Returns a reference to the pred map.
+ const PredMap& predMap() const {
+ return *_pred;
+ }
+
+ /// \brief Returns the predecessor edge of the given node.
+ ///
+ /// Returns the predecessor edge of the given node.
+ bool pred(Node node) const {
+ return (*_pred)[node];
}
/// \brief Returns the cost of the arborescence.
///
/// Returns the cost of the arborescence.
- Value arborescenceCost() const {
+ Value arborescenceValue() const {
Value sum = 0;
for (EdgeIt it(*graph); it != INVALID; ++it) {
- if (arborescenceEdge(it)) {
+ if (arborescence(it)) {
sum += (*cost)[it];
}
}
return sum;
}
+ /// \brief Indicates that a node is reachable from the sources.
+ ///
+ /// Indicates that a node is reachable from the sources.
+ bool reached(Node node) const {
+ return (*_node_order)[node] != -3;
+ }
+
+ /// \brief Indicates that a node is processed.
+ ///
+ /// Indicates that a node is processed. The arborescence path exists
+ /// from the source to the given node.
+ bool processed(Node node) const {
+ return (*_node_order)[node] == -1;
+ }
+
+ /// \brief Returns the number of the dual variables in basis.
+ ///
+ /// Returns the number of the dual variables in basis.
+ int dualSize() const {
+ return _dual_variables.size();
+ }
+
+ /// \brief Returns the value of the dual solution.
+ ///
+ /// Returns the value of the dual solution. It should be
+ /// equal to the arborescence value.
+ Value dualValue() const {
+ Value sum = 0;
+ for (int i = 0; i < (int)_dual_variables.size(); ++i) {
+ sum += _dual_variables[i].value;
+ }
+ return sum;
+ }
+
+ /// \brief Returns the number of the nodes in the dual variable.
+ ///
+ /// Returns the number of the nodes in the dual variable.
+ int dualSize(int k) const {
+ return _dual_variables[k].end - _dual_variables[k].begin;
+ }
+
+ /// \brief Returns the value of the dual variable.
+ ///
+ /// Returns the the value of the dual variable.
+ const Value& dualValue(int k) const {
+ return _dual_variables[k].value;
+ }
+
+ /// \brief Lemon iterator for get a dual variable.
+ ///
+ /// Lemon iterator for get a dual variable. This class provides
+ /// a common style lemon iterator which gives back a subset of
+ /// the nodes.
+ class DualIt {
+ public:
+
+ /// \brief Constructor.
+ ///
+ /// Constructor for get the nodeset of the variable.
+ DualIt(const MinCostArborescence& algorithm, int variable)
+ : _algorithm(&algorithm), _variable(variable)
+ {
+ _index = _algorithm->_dual_variables[_variable].begin;
+ }
+
+ /// \brief Invalid constructor.
+ ///
+ /// Invalid constructor.
+ DualIt(Invalid) : _algorithm(0) {}
+
+ /// \brief Conversion to node.
+ ///
+ /// Conversion to node.
+ operator Node() const {
+ return _algorithm ? _algorithm->_dual_node_list[_index] : INVALID;
+ }
+
+ /// \brief Increment operator.
+ ///
+ /// Increment operator.
+ DualIt& operator++() {
+ ++_index;
+ if (_algorithm->_dual_variables[_variable].end == _index) {
+ _algorithm = 0;
+ }
+ return *this;
+ }
+
+ bool operator==(const DualIt& it) const {
+ return (Node)(*this) == (Node)it;
+ }
+ bool operator!=(const DualIt& it) const {
+ return (Node)(*this) != (Node)it;
+ }
+ bool operator<(const DualIt& it) const {
+ return (Node)(*this) < (Node)it;
+ }
+
+ private:
+ const MinCostArborescence* _algorithm;
+ int _variable;
+ int _index;
+ };
+
/// @}
/// \name Execution control
@@ -279,7 +482,7 @@
/// If you need more control on the execution,
/// first you must call \ref init(), then you can add several
/// source nodes with \ref addSource().
- /// Finally \ref start() will perform the actual path
+ /// Finally \ref start() will perform the arborescence
/// computation.
///@{
@@ -290,13 +493,18 @@
///
void init() {
initStructures();
+ _heap->clear();
for (NodeIt it(*graph); it != INVALID; ++it) {
(*_cost_edges)[it].edge = INVALID;
- (*_level)[it] = -3;
+ _node_order->set(it, -3);
+ _heap_cross_ref->set(it, Heap::PRE_HEAP);
}
for (EdgeIt it(*graph); it != INVALID; ++it) {
- _arborescence_map->set(it, false);
+ _arborescence->set(it, false);
+ _edge_order->set(it, -1);
}
+ _dual_node_list.clear();
+ _dual_variables.clear();
}
/// \brief Adds a new source node.
@@ -309,14 +517,15 @@
Node node = nodes.back();
nodes.pop_back();
for (OutEdgeIt it(*graph, node); it != INVALID; ++it) {
- if ((*_level)[graph->target(it)] == -3) {
- (*_level)[graph->target(it)] = -2;
- nodes.push_back(graph->target(it));
- queue.push_back(graph->target(it));
+ Node target = graph->target(it);
+ if ((*_node_order)[target] == -3) {
+ (*_node_order)[target] = -2;
+ nodes.push_back(target);
+ queue.push_back(target);
}
}
}
- (*_level)[source] = -1;
+ (*_node_order)[source] = -1;
}
/// \brief Processes the next node in the priority queue.
@@ -327,19 +536,20 @@
///
/// \warning The queue must not be empty!
Node processNextNode() {
- node_counter = 0;
Node node = queue.back();
queue.pop_back();
- if ((*_level)[node] == -2) {
+ if ((*_node_order)[node] == -2) {
Edge edge = prepare(node);
- while ((*_level)[graph->source(edge)] != -1) {
- if ((*_level)[graph->source(edge)] >= 0) {
- edge = contract(bottom((*_level)[graph->source(edge)]));
+ Node source = graph->source(edge);
+ while ((*_node_order)[source] != -1) {
+ if ((*_node_order)[source] >= 0) {
+ edge = contract(source);
} else {
- edge = prepare(graph->source(edge));
+ edge = prepare(source);
}
+ source = graph->source(edge);
}
- finalize(graph->target(edge));
+ finalize(edge);
level_stack.clear();
}
return node;
@@ -400,46 +610,74 @@
protected:
void initStructures() {
- if (!_arborescence_map) {
- local_arborescence_map = true;
- _arborescence_map = Traits::createArborescenceMap(*graph);
+ if (!_pred) {
+ local_pred = true;
+ _pred = Traits::createPredMap(*graph);
+ }
+ if (!_arborescence) {
+ local_arborescence = true;
+ _arborescence = Traits::createArborescenceMap(*graph);
}
- if (!_level) {
- _level = new LevelMap(*graph);
+ if (!_edge_order) {
+ _edge_order = new EdgeOrder(*graph);
+ }
+ if (!_node_order) {
+ _node_order = new NodeOrder(*graph);
}
if (!_cost_edges) {
_cost_edges = new CostEdgeMap(*graph);
}
+ if (!_heap_cross_ref) {
+ _heap_cross_ref = new HeapCrossRef(*graph, -1);
+ }
+ if (!_heap) {
+ _heap = new Heap(*_heap_cross_ref);
+ }
}
void destroyStructures() {
- if (_level) {
- delete _level;
+ if (local_arborescence) {
+ delete _arborescence;
+ }
+ if (local_pred) {
+ delete _pred;
+ }
+ if (!_edge_order) {
+ delete _edge_order;
+ }
+ if (_node_order) {
+ delete _node_order;
}
if (!_cost_edges) {
delete _cost_edges;
}
- if (local_arborescence_map) {
- delete _arborescence_map;
+ if (!_heap) {
+ delete _heap;
+ }
+ if (!_heap_cross_ref) {
+ delete _heap_cross_ref;
}
}
Edge prepare(Node node) {
std::vector<Node> nodes;
- (*_level)[node] = node_counter;
+ (*_node_order)[node] = _dual_node_list.size();
+ StackLevel level;
+ level.node_level = _dual_node_list.size();
+ _dual_node_list.push_back(node);
for (InEdgeIt it(*graph, node); it != INVALID; ++it) {
Edge edge = it;
+ Node source = graph->source(edge);
Value value = (*cost)[it];
- if (graph->source(edge) == node ||
- (*_level)[graph->source(edge)] == -3) continue;
- if ((*_cost_edges)[graph->source(edge)].edge == INVALID) {
- (*_cost_edges)[graph->source(edge)].edge = edge;
- (*_cost_edges)[graph->source(edge)].value = value;
- nodes.push_back(graph->source(edge));
+ if (source == node || (*_node_order)[source] == -3) continue;
+ if ((*_cost_edges)[source].edge == INVALID) {
+ (*_cost_edges)[source].edge = edge;
+ (*_cost_edges)[source].value = value;
+ nodes.push_back(source);
} else {
- if ((*_cost_edges)[graph->source(edge)].value > value) {
- (*_cost_edges)[graph->source(edge)].edge = edge;
- (*_cost_edges)[graph->source(edge)].value = value;
+ if ((*_cost_edges)[source].value > value) {
+ (*_cost_edges)[source].edge = edge;
+ (*_cost_edges)[source].value = value;
}
}
}
@@ -449,35 +687,37 @@
minimum = (*_cost_edges)[nodes[i]];
}
}
- StackLevel level;
- level.node_level = node_counter;
+ _edge_order->set(minimum.edge, _dual_variables.size());
+ DualVariable var(_dual_node_list.size() - 1,
+ _dual_node_list.size(), minimum.value);
+ _dual_variables.push_back(var);
for (int i = 0; i < (int)nodes.size(); ++i) {
(*_cost_edges)[nodes[i]].value -= minimum.value;
level.edges.push_back((*_cost_edges)[nodes[i]]);
(*_cost_edges)[nodes[i]].edge = INVALID;
}
level_stack.push_back(level);
- ++node_counter;
- _arborescence_map->set(minimum.edge, true);
return minimum.edge;
}
- Edge contract(int node_bottom) {
+ Edge contract(Node node) {
+ int node_bottom = bottom(node);
std::vector<Node> nodes;
while (!level_stack.empty() &&
level_stack.back().node_level >= node_bottom) {
for (int i = 0; i < (int)level_stack.back().edges.size(); ++i) {
Edge edge = level_stack.back().edges[i].edge;
+ Node source = graph->source(edge);
Value value = level_stack.back().edges[i].value;
- if ((*_level)[graph->source(edge)] >= node_bottom) continue;
- if ((*_cost_edges)[graph->source(edge)].edge == INVALID) {
- (*_cost_edges)[graph->source(edge)].edge = edge;
- (*_cost_edges)[graph->source(edge)].value = value;
- nodes.push_back(graph->source(edge));
+ if ((*_node_order)[source] >= node_bottom) continue;
+ if ((*_cost_edges)[source].edge == INVALID) {
+ (*_cost_edges)[source].edge = edge;
+ (*_cost_edges)[source].value = value;
+ nodes.push_back(source);
} else {
- if ((*_cost_edges)[graph->source(edge)].value > value) {
- (*_cost_edges)[graph->source(edge)].edge = edge;
- (*_cost_edges)[graph->source(edge)].value = value;
+ if ((*_cost_edges)[source].value > value) {
+ (*_cost_edges)[source].edge = edge;
+ (*_cost_edges)[source].value = value;
}
}
}
@@ -489,6 +729,9 @@
minimum = (*_cost_edges)[nodes[i]];
}
}
+ _edge_order->set(minimum.edge, _dual_variables.size());
+ DualVariable var(node_bottom, _dual_node_list.size(), minimum.value);
+ _dual_variables.push_back(var);
StackLevel level;
level.node_level = node_bottom;
for (int i = 0; i < (int)nodes.size(); ++i) {
@@ -497,34 +740,45 @@
(*_cost_edges)[nodes[i]].edge = INVALID;
}
level_stack.push_back(level);
- _arborescence_map->set(minimum.edge, true);
return minimum.edge;
}
- int bottom(int level) {
+ int bottom(Node node) {
int k = level_stack.size() - 1;
- while (level_stack[k].node_level > level) {
+ while (level_stack[k].node_level > (*_node_order)[node]) {
--k;
}
return level_stack[k].node_level;
}
- void finalize(Node source) {
- std::vector<Node> nodes;
- nodes.push_back(source);
- while (!nodes.empty()) {
- Node node = nodes.back();
- nodes.pop_back();
- for (OutEdgeIt it(*graph, node); it != INVALID; ++it) {
- if ((*_level)[graph->target(it)] >= 0 && (*_arborescence_map)[it]) {
- (*_level)[graph->target(it)] = -1;
- nodes.push_back(graph->target(it));
- } else {
- _arborescence_map->set(it, false);
+ void finalize(Edge edge) {
+ Node node = graph->target(edge);
+ _heap->push(node, (*_edge_order)[edge]);
+ _pred->set(node, edge);
+ while (!_heap->empty()) {
+ Node source = _heap->top();
+ _heap->pop();
+ _node_order->set(source, -1);
+ for (OutEdgeIt it(*graph, source); it != INVALID; ++it) {
+ if ((*_edge_order)[it] < 0) continue;
+ Node target = graph->target(it);
+ switch(_heap->state(target)) {
+ case Heap::PRE_HEAP:
+ _heap->push(target, (*_edge_order)[it]);
+ _pred->set(target, it);
+ break;
+ case Heap::IN_HEAP:
+ if ((*_edge_order)[it] < (*_heap)[target]) {
+ _heap->decrease(target, (*_edge_order)[it]);
+ _pred->set(target, it);
+ }
+ break;
+ case Heap::POST_HEAP:
+ break;
}
}
+ _arborescence->set((*_pred)[source], true);
}
- (*_level)[source] = -1;
}
};
@@ -551,11 +805,9 @@
::Create mca(graph, cost);
mca.arborescenceMap(arborescence);
mca.run(source);
- return mca.arborescenceCost();
+ return mca.arborescenceValue();
}
}
#endif
-
-// Hilbert - Huang
Modified: hugo/trunk/test/Makefile.am
==============================================================================
--- hugo/trunk/test/Makefile.am (original)
+++ hugo/trunk/test/Makefile.am Fri Mar 31 13:10:44 2006
@@ -13,6 +13,7 @@
check_PROGRAMS = \
all_pairs_shortest_path_test \
+ arborescence_test \
bfs_test \
counter_test \
dfs_test \
@@ -52,6 +53,7 @@
XFAIL_TESTS = test_tools_fail$(EXEEXT)
all_pairs_shortest_path_test_SOURCES = all_pairs_shortest_path_test.cc
+arborescence_test_SOURCES = arborescence_test.cc
bfs_test_SOURCES = bfs_test.cc
counter_test_SOURCES = counter_test.cc
dfs_test_SOURCES = dfs_test.cc
Added: hugo/trunk/test/arborescence_test.cc
==============================================================================
--- (empty file)
+++ hugo/trunk/test/arborescence_test.cc Fri Mar 31 13:10:44 2006
@@ -0,0 +1,104 @@
+#include <iostream>
+#include <set>
+#include <vector>
+#include <iterator>
+
+#include <cmath>
+#include <cstdlib>
+
+#include <lemon/smart_graph.h>
+#include <lemon/min_cost_arborescence.h>
+
+#include <lemon/graph_utils.h>
+#include <lemon/time_measure.h>
+
+#include <lemon/tolerance.h>
+
+using namespace lemon;
+using namespace std;
+
+int main(int argc, const char *argv[]) {
+ srand(time(0));
+ typedef SmartGraph Graph;
+ GRAPH_TYPEDEFS(Graph);
+
+ typedef Graph::EdgeMap<double> CostMap;
+
+ const int n = argc > 1 ? atoi(argv[1]) : 100;
+ const int e = argc > 2 ? atoi(argv[2]) : (int)(n * log(n));
+
+ Graph graph;
+ CostMap cost(graph);
+ vector<Node> nodes;
+
+ for (int i = 0; i < n; ++i) {
+ nodes.push_back(graph.addNode());
+ }
+
+ for (int i = 0; i < e; ++i) {
+ int s = (int)(n * (double)rand() / (RAND_MAX + 1.0));
+ int t = (int)(n * (double)rand() / (RAND_MAX + 1.0));
+ double c = rand() / (1.0 + RAND_MAX) * 100.0 + 20.0;
+ Edge edge = graph.addEdge(nodes[s], nodes[t]);
+ cost[edge] = c;
+ }
+
+ Node source = nodes[(int)(n * (double)rand() / (RAND_MAX + 1.0))];
+
+ MinCostArborescence<Graph, CostMap> mca(graph, cost);
+ mca.run(source);
+
+ vector<pair<double, set<Node> > > dualSolution(mca.dualSize());
+
+ for (int i = 0; i < mca.dualSize(); ++i) {
+ dualSolution[i].first = mca.dualValue(i);
+ for (MinCostArborescence<Graph, CostMap>::DualIt it(mca, i);
+ it != INVALID; ++it) {
+ dualSolution[i].second.insert(it);
+ }
+ }
+
+ Tolerance<double> tolerance;
+
+ for (EdgeIt it(graph); it != INVALID; ++it) {
+ if (mca.reached(graph.source(it))) {
+ double sum = 0.0;
+ for (int i = 0; i < (int)dualSolution.size(); ++i) {
+ if (dualSolution[i].second.find(graph.target(it))
+ != dualSolution[i].second.end() &&
+ dualSolution[i].second.find(graph.source(it))
+ == dualSolution[i].second.end()) {
+ sum += dualSolution[i].first;
+ }
+ }
+ if (mca.arborescence(it)) {
+ LEMON_ASSERT(!tolerance.less(sum, cost[it]), "INVALID DUAL");
+ }
+ LEMON_ASSERT(!tolerance.less(cost[it], sum), "INVALID DUAL");
+ }
+ }
+
+
+ LEMON_ASSERT(!tolerance.different(mca.dualValue(), mca.arborescenceValue()),
+ "INVALID DUAL");
+
+
+ LEMON_ASSERT(mca.reached(source), "INVALID ARBORESCENCE");
+ for (EdgeIt it(graph); it != INVALID; ++it) {
+ LEMON_ASSERT(!mca.reached(graph.source(it)) ||
+ mca.reached(graph.target(it)), "INVALID ARBORESCENCE");
+ }
+
+ for (NodeIt it(graph); it != INVALID; ++it) {
+ if (!mca.reached(it)) continue;
+ int cnt = 0;
+ for (InEdgeIt jt(graph, it); jt != INVALID; ++jt) {
+ if (mca.arborescence(jt)) {
+ ++cnt;
+ }
+ }
+ LEMON_ASSERT((it == source ? cnt == 0 : cnt == 1), "INVALID ARBORESCENCE");
+ }
+
+ return 0;
+}
More information about the Lemon-commits
mailing list