File indexing completed on 2024-12-22 04:10:11
0001 // SPDX-FileCopyrightText: 2006 Stephan Diederich 0002 // 0003 // This code may be used under either of the following two licences: 0004 // 0005 // SPDX-License-Identifier: MIT OR BSL-1.0 0006 0007 #ifndef BOOST_BOYKOV_KOLMOGOROV_MAX_FLOW_HPP 0008 #define BOOST_BOYKOV_KOLMOGOROV_MAX_FLOW_HPP 0009 0010 #include <boost/config.hpp> 0011 #include <boost/assert.hpp> 0012 #include <vector> 0013 #include <list> 0014 #include <utility> 0015 #include <iosfwd> 0016 #include <algorithm> // for std::min and std::max 0017 0018 #include <boost/pending/queue.hpp> 0019 #include <boost/limits.hpp> 0020 #include <boost/property_map/property_map.hpp> 0021 #include <boost/none_t.hpp> 0022 #include <boost/graph/graph_concepts.hpp> 0023 #include <boost/graph/named_function_params.hpp> 0024 #include <boost/graph/lookup_edge.hpp> 0025 #include <boost/concept/assert.hpp> 0026 0027 // The algorithm implemented here is described in: 0028 // 0029 // Boykov, Y., Kolmogorov, V. "An Experimental Comparison of Min-Cut/Max-Flow 0030 // Algorithms for Energy Minimization in Vision", In IEEE Transactions on 0031 // Pattern Analysis and Machine Intelligence, vol. 26, no. 9, pp. 1124-1137, 0032 // Sep 2004. 0033 // 0034 // For further reading, also see: 0035 // 0036 // Kolmogorov, V. "Graph Based Algorithms for Scene Reconstruction from Two or 0037 // More Views". PhD thesis, Cornell University, Sep 2003. 0038 0039 namespace boost { 0040 0041 namespace detail { 0042 0043 template <class Graph, 0044 class EdgeCapacityMap, 0045 class ResidualCapacityEdgeMap, 0046 class ReverseEdgeMap, 0047 class PredecessorMap, 0048 class ColorMap, 0049 class DistanceMap, 0050 class IndexMap> 0051 class bk_max_flow { 0052 typedef typename property_traits<EdgeCapacityMap>::value_type tEdgeVal; 0053 typedef graph_traits<Graph> tGraphTraits; 0054 typedef typename tGraphTraits::vertex_iterator vertex_iterator; 0055 typedef typename tGraphTraits::vertex_descriptor vertex_descriptor; 0056 typedef typename tGraphTraits::edge_descriptor edge_descriptor; 0057 typedef typename tGraphTraits::edge_iterator edge_iterator; 0058 typedef typename tGraphTraits::out_edge_iterator out_edge_iterator; 0059 typedef boost::queue<vertex_descriptor> tQueue; //queue of vertices, used in adoption-stage 0060 typedef typename property_traits<ColorMap>::value_type tColorValue; 0061 typedef color_traits<tColorValue> tColorTraits; 0062 typedef typename property_traits<DistanceMap>::value_type tDistanceVal; 0063 0064 public: 0065 bk_max_flow(Graph& g, 0066 EdgeCapacityMap cap, 0067 ResidualCapacityEdgeMap res, 0068 ReverseEdgeMap rev, 0069 PredecessorMap pre, 0070 ColorMap color, 0071 DistanceMap dist, 0072 IndexMap idx, 0073 vertex_descriptor src, 0074 vertex_descriptor sink): 0075 m_g(g), 0076 m_index_map(idx), 0077 m_cap_map(cap), 0078 m_res_cap_map(res), 0079 m_rev_edge_map(rev), 0080 m_pre_map(pre), 0081 m_tree_map(color), 0082 m_dist_map(dist), 0083 m_source(src), 0084 m_sink(sink), 0085 m_active_nodes(), 0086 m_in_active_list_vec(num_vertices(g), false), 0087 m_in_active_list_map(make_iterator_property_map(m_in_active_list_vec.begin(), m_index_map)), 0088 m_has_parent_vec(num_vertices(g), false), 0089 m_has_parent_map(make_iterator_property_map(m_has_parent_vec.begin(), m_index_map)), 0090 m_time_vec(num_vertices(g), 0), 0091 m_time_map(make_iterator_property_map(m_time_vec.begin(), m_index_map)), 0092 m_flow(0), 0093 m_time(1), 0094 m_last_grow_vertex(graph_traits<Graph>::null_vertex()){ 0095 // initialize the color-map with gray-values 0096 vertex_iterator vi, v_end; 0097 for(boost::tie(vi, v_end) = vertices(m_g); vi != v_end; ++vi){ 0098 set_tree(*vi, tColorTraits::gray()); 0099 } 0100 // Initialize flow to zero which means initializing 0101 // the residual capacity equal to the capacity 0102 edge_iterator ei, e_end; 0103 for(boost::tie(ei, e_end) = edges(m_g); ei != e_end; ++ei) { 0104 put(m_res_cap_map, *ei, get(m_cap_map, *ei)); 0105 BOOST_ASSERT(get(m_rev_edge_map, get(m_rev_edge_map, *ei)) == *ei); //check if the reverse edge map is build up properly 0106 } 0107 //init the search trees with the two terminals 0108 set_tree(m_source, tColorTraits::black()); 0109 set_tree(m_sink, tColorTraits::white()); 0110 put(m_time_map, m_source, 1); 0111 put(m_time_map, m_sink, 1); 0112 } 0113 0114 tEdgeVal max_flow(){ 0115 //augment direct paths from SOURCE->SINK and SOURCE->VERTEX->SINK 0116 augment_direct_paths(); 0117 //start the main-loop 0118 while(true){ 0119 bool path_found; 0120 edge_descriptor connecting_edge; 0121 boost::tie(connecting_edge, path_found) = grow(); //find a path from source to sink 0122 if(!path_found){ 0123 //we're finished, no more paths were found 0124 break; 0125 } 0126 ++m_time; 0127 augment(connecting_edge); //augment that path 0128 adopt(); //rebuild search tree structure 0129 } 0130 return m_flow; 0131 } 0132 0133 // the complete class is protected, as we want access to members in 0134 // derived test-class (see test/boykov_kolmogorov_max_flow_test.cpp) 0135 protected: 0136 void augment_direct_paths(){ 0137 // in a first step, we augment all direct paths from source->NODE->sink 0138 // and additionally paths from source->sink. This improves especially 0139 // graphcuts for segmentation, as most of the nodes have source/sink 0140 // connects but shouldn't have an impact on other maxflow problems 0141 // (this is done in grow() anyway) 0142 out_edge_iterator ei, e_end; 0143 for(boost::tie(ei, e_end) = out_edges(m_source, m_g); ei != e_end; ++ei){ 0144 edge_descriptor from_source = *ei; 0145 vertex_descriptor current_node = target(from_source, m_g); 0146 if(current_node == m_sink){ 0147 tEdgeVal cap = get(m_res_cap_map, from_source); 0148 put(m_res_cap_map, from_source, 0); 0149 m_flow += cap; 0150 continue; 0151 } 0152 edge_descriptor to_sink; 0153 bool is_there; 0154 boost::tie(to_sink, is_there) = lookup_edge(current_node, m_sink, m_g); 0155 if(is_there){ 0156 tEdgeVal cap_from_source = get(m_res_cap_map, from_source); 0157 tEdgeVal cap_to_sink = get(m_res_cap_map, to_sink); 0158 if(cap_from_source > cap_to_sink){ 0159 set_tree(current_node, tColorTraits::black()); 0160 add_active_node(current_node); 0161 set_edge_to_parent(current_node, from_source); 0162 put(m_dist_map, current_node, 1); 0163 put(m_time_map, current_node, 1); 0164 // add stuff to flow and update residuals. we don't need to 0165 // update reverse_edges, as incoming/outgoing edges to/from 0166 // source/sink don't count for max-flow 0167 put(m_res_cap_map, from_source, get(m_res_cap_map, from_source) - cap_to_sink); 0168 put(m_res_cap_map, to_sink, 0); 0169 m_flow += cap_to_sink; 0170 } else if(cap_to_sink > 0){ 0171 set_tree(current_node, tColorTraits::white()); 0172 add_active_node(current_node); 0173 set_edge_to_parent(current_node, to_sink); 0174 put(m_dist_map, current_node, 1); 0175 put(m_time_map, current_node, 1); 0176 // add stuff to flow and update residuals. we don't need to update 0177 // reverse_edges, as incoming/outgoing edges to/from source/sink 0178 // don't count for max-flow 0179 put(m_res_cap_map, to_sink, get(m_res_cap_map, to_sink) - cap_from_source); 0180 put(m_res_cap_map, from_source, 0); 0181 m_flow += cap_from_source; 0182 } 0183 } else if(get(m_res_cap_map, from_source)){ 0184 // there is no sink connect, so we can't augment this path, but to 0185 // avoid adding m_source to the active nodes, we just activate this 0186 // node and set the appropriate things 0187 set_tree(current_node, tColorTraits::black()); 0188 set_edge_to_parent(current_node, from_source); 0189 put(m_dist_map, current_node, 1); 0190 put(m_time_map, current_node, 1); 0191 add_active_node(current_node); 0192 } 0193 } 0194 for(boost::tie(ei, e_end) = out_edges(m_sink, m_g); ei != e_end; ++ei){ 0195 edge_descriptor to_sink = get(m_rev_edge_map, *ei); 0196 vertex_descriptor current_node = source(to_sink, m_g); 0197 if(get(m_res_cap_map, to_sink)){ 0198 set_tree(current_node, tColorTraits::white()); 0199 set_edge_to_parent(current_node, to_sink); 0200 put(m_dist_map, current_node, 1); 0201 put(m_time_map, current_node, 1); 0202 add_active_node(current_node); 0203 } 0204 } 0205 } 0206 0207 /** 0208 * Returns a pair of an edge and a boolean. if the bool is true, the 0209 * edge is a connection of a found path from s->t , read "the link" and 0210 * source(returnVal, m_g) is the end of the path found in the source-tree 0211 * target(returnVal, m_g) is the beginning of the path found in the sink-tree 0212 */ 0213 std::pair<edge_descriptor, bool> grow(){ 0214 BOOST_ASSERT(m_orphans.empty()); 0215 vertex_descriptor current_node; 0216 while((current_node = get_next_active_node()) != graph_traits<Graph>::null_vertex()){ //if there is one 0217 BOOST_ASSERT(get_tree(current_node) != tColorTraits::gray() && 0218 (has_parent(current_node) || 0219 current_node == m_source || 0220 current_node == m_sink)); 0221 0222 if(get_tree(current_node) == tColorTraits::black()){ 0223 //source tree growing 0224 out_edge_iterator ei, e_end; 0225 if(current_node != m_last_grow_vertex){ 0226 m_last_grow_vertex = current_node; 0227 boost::tie(m_last_grow_edge_it, m_last_grow_edge_end) = out_edges(current_node, m_g); 0228 } 0229 for(; m_last_grow_edge_it != m_last_grow_edge_end; ++m_last_grow_edge_it) { 0230 edge_descriptor out_edge = *m_last_grow_edge_it; 0231 if(get(m_res_cap_map, out_edge) > 0){ //check if we have capacity left on this edge 0232 vertex_descriptor other_node = target(out_edge, m_g); 0233 if(get_tree(other_node) == tColorTraits::gray()){ //it's a free node 0234 set_tree(other_node, tColorTraits::black()); //acquire other node to our search tree 0235 set_edge_to_parent(other_node, out_edge); //set us as parent 0236 put(m_dist_map, other_node, get(m_dist_map, current_node) + 1); //and update the distance-heuristic 0237 put(m_time_map, other_node, get(m_time_map, current_node)); 0238 add_active_node(other_node); 0239 } else if(get_tree(other_node) == tColorTraits::black()) { 0240 // we do this to get shorter paths. check if we are nearer to 0241 // the source as its parent is 0242 if(is_closer_to_terminal(current_node, other_node)){ 0243 set_edge_to_parent(other_node, out_edge); 0244 put(m_dist_map, other_node, get(m_dist_map, current_node) + 1); 0245 put(m_time_map, other_node, get(m_time_map, current_node)); 0246 } 0247 } else{ 0248 BOOST_ASSERT(get_tree(other_node)==tColorTraits::white()); 0249 //kewl, found a path from one to the other search tree, return 0250 // the connecting edge in src->sink dir 0251 return std::make_pair(out_edge, true); 0252 } 0253 } 0254 } //for all out-edges 0255 } //source-tree-growing 0256 else{ 0257 BOOST_ASSERT(get_tree(current_node) == tColorTraits::white()); 0258 out_edge_iterator ei, e_end; 0259 if(current_node != m_last_grow_vertex){ 0260 m_last_grow_vertex = current_node; 0261 boost::tie(m_last_grow_edge_it, m_last_grow_edge_end) = out_edges(current_node, m_g); 0262 } 0263 for(; m_last_grow_edge_it != m_last_grow_edge_end; ++m_last_grow_edge_it){ 0264 edge_descriptor in_edge = get(m_rev_edge_map, *m_last_grow_edge_it); 0265 if(get(m_res_cap_map, in_edge) > 0){ //check if there is capacity left 0266 vertex_descriptor other_node = source(in_edge, m_g); 0267 if(get_tree(other_node) == tColorTraits::gray()){ //it's a free node 0268 set_tree(other_node, tColorTraits::white()); //acquire that node to our search tree 0269 set_edge_to_parent(other_node, in_edge); //set us as parent 0270 add_active_node(other_node); //activate that node 0271 put(m_dist_map, other_node, get(m_dist_map, current_node) + 1); //set its distance 0272 put(m_time_map, other_node, get(m_time_map, current_node));//and time 0273 } else if(get_tree(other_node) == tColorTraits::white()){ 0274 if(is_closer_to_terminal(current_node, other_node)){ 0275 //we are closer to the sink than its parent is, so we "adopt" him 0276 set_edge_to_parent(other_node, in_edge); 0277 put(m_dist_map, other_node, get(m_dist_map, current_node) + 1); 0278 put(m_time_map, other_node, get(m_time_map, current_node)); 0279 } 0280 } else{ 0281 BOOST_ASSERT(get_tree(other_node)==tColorTraits::black()); 0282 //kewl, found a path from one to the other search tree, 0283 // return the connecting edge in src->sink dir 0284 return std::make_pair(in_edge, true); 0285 } 0286 } 0287 } //for all out-edges 0288 } //sink-tree growing 0289 0290 //all edges of that node are processed, and no more paths were found. 0291 // remove if from the front of the active queue 0292 finish_node(current_node); 0293 } //while active_nodes not empty 0294 0295 //no active nodes anymore and no path found, we're done 0296 return std::make_pair(edge_descriptor(), false); 0297 } 0298 0299 /** 0300 * augments path from s->t and updates residual graph 0301 * source(e, m_g) is the end of the path found in the source-tree 0302 * target(e, m_g) is the beginning of the path found in the sink-tree 0303 * this phase generates orphans on saturated edges, if the attached verts are 0304 * from different search-trees orphans are ordered in distance to 0305 * sink/source. first the farthest from the source are front_inserted into 0306 * the orphans list, and after that the sink-tree-orphans are 0307 * front_inserted. when going to adoption stage the orphans are popped_front, 0308 * and so we process the nearest verts to the terminals first 0309 */ 0310 void augment(edge_descriptor e) { 0311 BOOST_ASSERT(get_tree(target(e, m_g)) == tColorTraits::white()); 0312 BOOST_ASSERT(get_tree(source(e, m_g)) == tColorTraits::black()); 0313 BOOST_ASSERT(m_orphans.empty()); 0314 0315 const tEdgeVal bottleneck = find_bottleneck(e); 0316 //now we push the found flow through the path 0317 //for each edge we saturate we have to look for the verts that belong to that edge, one of them becomes an orphans 0318 //now process the connecting edge 0319 put(m_res_cap_map, e, get(m_res_cap_map, e) - bottleneck); 0320 BOOST_ASSERT(get(m_res_cap_map, e) >= 0); 0321 put(m_res_cap_map, get(m_rev_edge_map, e), get(m_res_cap_map, get(m_rev_edge_map, e)) + bottleneck); 0322 0323 //now we follow the path back to the source 0324 vertex_descriptor current_node = source(e, m_g); 0325 while(current_node != m_source){ 0326 edge_descriptor pred = get_edge_to_parent(current_node); 0327 const tEdgeVal new_pred_cap = get(m_res_cap_map, pred) - bottleneck; 0328 put(m_res_cap_map, pred, new_pred_cap); 0329 BOOST_ASSERT(get(m_res_cap_map, pred) >= 0); 0330 0331 const edge_descriptor pred_rev = get(m_rev_edge_map, pred); 0332 put(m_res_cap_map, pred_rev, get(m_res_cap_map, pred_rev) + bottleneck); 0333 0334 if(new_pred_cap == 0){ 0335 set_no_parent(current_node); 0336 m_orphans.push_front(current_node); 0337 } 0338 current_node = source(pred, m_g); 0339 } 0340 //then go forward in the sink-tree 0341 current_node = target(e, m_g); 0342 while(current_node != m_sink){ 0343 edge_descriptor pred = get_edge_to_parent(current_node); 0344 const tEdgeVal new_pred_cap = get(m_res_cap_map, pred) - bottleneck; 0345 put(m_res_cap_map, pred, new_pred_cap); 0346 BOOST_ASSERT(get(m_res_cap_map, pred) >= 0); 0347 0348 const edge_descriptor pred_rev = get(m_rev_edge_map, pred); 0349 put(m_res_cap_map, pred_rev, get(m_res_cap_map, pred_rev) + bottleneck); 0350 0351 if(new_pred_cap == 0){ 0352 set_no_parent(current_node); 0353 m_orphans.push_front(current_node); 0354 } 0355 current_node = target(pred, m_g); 0356 } 0357 //and add it to the max-flow 0358 m_flow += bottleneck; 0359 } 0360 0361 /** 0362 * returns the bottleneck of a s->t path (end_of_path is last vertex in 0363 * source-tree, begin_of_path is first vertex in sink-tree) 0364 */ 0365 inline tEdgeVal find_bottleneck(edge_descriptor e){ 0366 BOOST_USING_STD_MIN(); 0367 tEdgeVal minimum_cap = get(m_res_cap_map, e); 0368 vertex_descriptor current_node = source(e, m_g); 0369 //first go back in the source tree 0370 while(current_node != m_source){ 0371 edge_descriptor pred = get_edge_to_parent(current_node); 0372 minimum_cap = min BOOST_PREVENT_MACRO_SUBSTITUTION(minimum_cap, get(m_res_cap_map, pred)); 0373 current_node = source(pred, m_g); 0374 } 0375 //then go forward in the sink-tree 0376 current_node = target(e, m_g); 0377 while(current_node != m_sink){ 0378 edge_descriptor pred = get_edge_to_parent(current_node); 0379 minimum_cap = min BOOST_PREVENT_MACRO_SUBSTITUTION(minimum_cap, get(m_res_cap_map, pred)); 0380 current_node = target(pred, m_g); 0381 } 0382 return minimum_cap; 0383 } 0384 0385 /** 0386 * rebuild search trees 0387 * empty the queue of orphans, and find new parents for them or just drop 0388 * them from the search trees 0389 */ 0390 void adopt(){ 0391 while(!m_orphans.empty() || !m_child_orphans.empty()){ 0392 vertex_descriptor current_node; 0393 if(m_child_orphans.empty()){ 0394 //get the next orphan from the main-queue and remove it 0395 current_node = m_orphans.front(); 0396 m_orphans.pop_front(); 0397 } else{ 0398 current_node = m_child_orphans.front(); 0399 m_child_orphans.pop(); 0400 } 0401 if(get_tree(current_node) == tColorTraits::black()){ 0402 //we're in the source-tree 0403 tDistanceVal min_distance = (std::numeric_limits<tDistanceVal>::max)(); 0404 edge_descriptor new_parent_edge; 0405 out_edge_iterator ei, e_end; 0406 for(boost::tie(ei, e_end) = out_edges(current_node, m_g); ei != e_end; ++ei){ 0407 const edge_descriptor in_edge = get(m_rev_edge_map, *ei); 0408 BOOST_ASSERT(target(in_edge, m_g) == current_node); //we should be the target of this edge 0409 if(get(m_res_cap_map, in_edge) > 0){ 0410 vertex_descriptor other_node = source(in_edge, m_g); 0411 if(get_tree(other_node) == tColorTraits::black() && has_source_connect(other_node)){ 0412 if(get(m_dist_map, other_node) < min_distance){ 0413 min_distance = get(m_dist_map, other_node); 0414 new_parent_edge = in_edge; 0415 } 0416 } 0417 } 0418 } 0419 if(min_distance != (std::numeric_limits<tDistanceVal>::max)()){ 0420 set_edge_to_parent(current_node, new_parent_edge); 0421 put(m_dist_map, current_node, min_distance + 1); 0422 put(m_time_map, current_node, m_time); 0423 } else{ 0424 put(m_time_map, current_node, 0); 0425 for(boost::tie(ei, e_end) = out_edges(current_node, m_g); ei != e_end; ++ei){ 0426 edge_descriptor in_edge = get(m_rev_edge_map, *ei); 0427 vertex_descriptor other_node = source(in_edge, m_g); 0428 if(get_tree(other_node) == tColorTraits::black() && other_node != m_source){ 0429 if(get(m_res_cap_map, in_edge) > 0){ 0430 add_active_node(other_node); 0431 } 0432 if(has_parent(other_node) && source(get_edge_to_parent(other_node), m_g) == current_node){ 0433 //we are the parent of that node 0434 //it has to find a new parent, too 0435 set_no_parent(other_node); 0436 m_child_orphans.push(other_node); 0437 } 0438 } 0439 } 0440 set_tree(current_node, tColorTraits::gray()); 0441 } //no parent found 0442 } //source-tree-adoption 0443 else{ 0444 //now we should be in the sink-tree, check that... 0445 BOOST_ASSERT(get_tree(current_node) == tColorTraits::white()); 0446 out_edge_iterator ei, e_end; 0447 edge_descriptor new_parent_edge; 0448 tDistanceVal min_distance = (std::numeric_limits<tDistanceVal>::max)(); 0449 for(boost::tie(ei, e_end) = out_edges(current_node, m_g); ei != e_end; ++ei){ 0450 const edge_descriptor out_edge = *ei; 0451 if(get(m_res_cap_map, out_edge) > 0){ 0452 const vertex_descriptor other_node = target(out_edge, m_g); 0453 if(get_tree(other_node) == tColorTraits::white() && has_sink_connect(other_node)) 0454 if(get(m_dist_map, other_node) < min_distance){ 0455 min_distance = get(m_dist_map, other_node); 0456 new_parent_edge = out_edge; 0457 } 0458 } 0459 } 0460 if(min_distance != (std::numeric_limits<tDistanceVal>::max)()){ 0461 set_edge_to_parent(current_node, new_parent_edge); 0462 put(m_dist_map, current_node, min_distance + 1); 0463 put(m_time_map, current_node, m_time); 0464 } else{ 0465 put(m_time_map, current_node, 0); 0466 for(boost::tie(ei, e_end) = out_edges(current_node, m_g); ei != e_end; ++ei){ 0467 const edge_descriptor out_edge = *ei; 0468 const vertex_descriptor other_node = target(out_edge, m_g); 0469 if(get_tree(other_node) == tColorTraits::white() && other_node != m_sink){ 0470 if(get(m_res_cap_map, out_edge) > 0){ 0471 add_active_node(other_node); 0472 } 0473 if(has_parent(other_node) && target(get_edge_to_parent(other_node), m_g) == current_node){ 0474 //we were it's parent, so it has to find a new one, too 0475 set_no_parent(other_node); 0476 m_child_orphans.push(other_node); 0477 } 0478 } 0479 } 0480 set_tree(current_node, tColorTraits::gray()); 0481 } //no parent found 0482 } //sink-tree adoption 0483 } //while !orphans.empty() 0484 } //adopt 0485 0486 /** 0487 * return next active vertex if there is one, otherwise a null_vertex 0488 */ 0489 inline vertex_descriptor get_next_active_node(){ 0490 while(true){ 0491 if(m_active_nodes.empty()) 0492 return graph_traits<Graph>::null_vertex(); 0493 vertex_descriptor v = m_active_nodes.front(); 0494 0495 //if it has no parent, this node can't be active (if its not source or sink) 0496 if(!has_parent(v) && v != m_source && v != m_sink){ 0497 m_active_nodes.pop(); 0498 put(m_in_active_list_map, v, false); 0499 } else{ 0500 BOOST_ASSERT(get_tree(v) == tColorTraits::black() || get_tree(v) == tColorTraits::white()); 0501 return v; 0502 } 0503 } 0504 } 0505 0506 /** 0507 * adds v as an active vertex, but only if its not in the list already 0508 */ 0509 inline void add_active_node(vertex_descriptor v){ 0510 BOOST_ASSERT(get_tree(v) != tColorTraits::gray()); 0511 if(get(m_in_active_list_map, v)){ 0512 if (m_last_grow_vertex == v) { 0513 m_last_grow_vertex = graph_traits<Graph>::null_vertex(); 0514 } 0515 return; 0516 } else{ 0517 put(m_in_active_list_map, v, true); 0518 m_active_nodes.push(v); 0519 } 0520 } 0521 0522 /** 0523 * finish_node removes a node from the front of the active queue (its called in grow phase, if no more paths can be found using this node) 0524 */ 0525 inline void finish_node(vertex_descriptor v){ 0526 BOOST_ASSERT(m_active_nodes.front() == v); 0527 m_active_nodes.pop(); 0528 put(m_in_active_list_map, v, false); 0529 m_last_grow_vertex = graph_traits<Graph>::null_vertex(); 0530 } 0531 0532 /** 0533 * removes a vertex from the queue of active nodes (actually this does nothing, 0534 * but checks if this node has no parent edge, as this is the criteria for 0535 * being no more active) 0536 */ 0537 inline void remove_active_node(vertex_descriptor v){ 0538 (void)v; // disable unused parameter warning 0539 BOOST_ASSERT(!has_parent(v)); 0540 } 0541 0542 /** 0543 * returns the search tree of v; tColorValue::black() for source tree, 0544 * white() for sink tree, gray() for no tree 0545 */ 0546 inline tColorValue get_tree(vertex_descriptor v) const { 0547 return get(m_tree_map, v); 0548 } 0549 0550 /** 0551 * sets search tree of v; tColorValue::black() for source tree, white() 0552 * for sink tree, gray() for no tree 0553 */ 0554 inline void set_tree(vertex_descriptor v, tColorValue t){ 0555 put(m_tree_map, v, t); 0556 } 0557 0558 /** 0559 * returns edge to parent vertex of v; 0560 */ 0561 inline edge_descriptor get_edge_to_parent(vertex_descriptor v) const{ 0562 return get(m_pre_map, v); 0563 } 0564 0565 /** 0566 * returns true if the edge stored in m_pre_map[v] is a valid entry 0567 */ 0568 inline bool has_parent(vertex_descriptor v) const{ 0569 return get(m_has_parent_map, v); 0570 } 0571 0572 /** 0573 * sets edge to parent vertex of v; 0574 */ 0575 inline void set_edge_to_parent(vertex_descriptor v, edge_descriptor f_edge_to_parent){ 0576 BOOST_ASSERT(get(m_res_cap_map, f_edge_to_parent) > 0); 0577 put(m_pre_map, v, f_edge_to_parent); 0578 put(m_has_parent_map, v, true); 0579 } 0580 0581 /** 0582 * removes the edge to parent of v (this is done by invalidating the 0583 * entry an additional map) 0584 */ 0585 inline void set_no_parent(vertex_descriptor v){ 0586 put(m_has_parent_map, v, false); 0587 } 0588 0589 /** 0590 * checks if vertex v has a connect to the sink-vertex (@p m_sink) 0591 * @param v the vertex which is checked 0592 * @return true if a path to the sink was found, false if not 0593 */ 0594 inline bool has_sink_connect(vertex_descriptor v){ 0595 tDistanceVal current_distance = 0; 0596 vertex_descriptor current_vertex = v; 0597 while(true){ 0598 if(get(m_time_map, current_vertex) == m_time){ 0599 //we found a node which was already checked this round. use it for distance calculations 0600 current_distance += get(m_dist_map, current_vertex); 0601 break; 0602 } 0603 if(current_vertex == m_sink){ 0604 put(m_time_map, m_sink, m_time); 0605 break; 0606 } 0607 if(has_parent(current_vertex)){ 0608 //it has a parent, so get it 0609 current_vertex = target(get_edge_to_parent(current_vertex), m_g); 0610 ++current_distance; 0611 } else{ 0612 //no path found 0613 return false; 0614 } 0615 } 0616 current_vertex=v; 0617 while(get(m_time_map, current_vertex) != m_time){ 0618 put(m_dist_map, current_vertex, current_distance); 0619 --current_distance; 0620 put(m_time_map, current_vertex, m_time); 0621 current_vertex = target(get_edge_to_parent(current_vertex), m_g); 0622 } 0623 return true; 0624 } 0625 0626 /** 0627 * checks if vertex v has a connect to the source-vertex (@p m_source) 0628 * @param v the vertex which is checked 0629 * @return true if a path to the source was found, false if not 0630 */ 0631 inline bool has_source_connect(vertex_descriptor v){ 0632 tDistanceVal current_distance = 0; 0633 vertex_descriptor current_vertex = v; 0634 while(true){ 0635 if(get(m_time_map, current_vertex) == m_time){ 0636 //we found a node which was already checked this round. use it for distance calculations 0637 current_distance += get(m_dist_map, current_vertex); 0638 break; 0639 } 0640 if(current_vertex == m_source){ 0641 put(m_time_map, m_source, m_time); 0642 break; 0643 } 0644 if(has_parent(current_vertex)){ 0645 //it has a parent, so get it 0646 current_vertex = source(get_edge_to_parent(current_vertex), m_g); 0647 ++current_distance; 0648 } else{ 0649 //no path found 0650 return false; 0651 } 0652 } 0653 current_vertex=v; 0654 while(get(m_time_map, current_vertex) != m_time){ 0655 put(m_dist_map, current_vertex, current_distance); 0656 --current_distance; 0657 put(m_time_map, current_vertex, m_time); 0658 current_vertex = source(get_edge_to_parent(current_vertex), m_g); 0659 } 0660 return true; 0661 } 0662 0663 /** 0664 * returns true, if p is closer to a terminal than q 0665 */ 0666 inline bool is_closer_to_terminal(vertex_descriptor p, vertex_descriptor q){ 0667 //checks the timestamps first, to build no cycles, and after that the real distance 0668 return (get(m_time_map, q) <= get(m_time_map, p) && 0669 get(m_dist_map, q) > get(m_dist_map, p)+1); 0670 } 0671 0672 //////// 0673 // member vars 0674 //////// 0675 Graph& m_g; 0676 IndexMap m_index_map; 0677 EdgeCapacityMap m_cap_map; 0678 ResidualCapacityEdgeMap m_res_cap_map; 0679 ReverseEdgeMap m_rev_edge_map; 0680 PredecessorMap m_pre_map; //stores paths found in the growth stage 0681 ColorMap m_tree_map; //maps each vertex into one of the two search tree or none (gray()) 0682 DistanceMap m_dist_map; //stores distance to source/sink nodes 0683 vertex_descriptor m_source; 0684 vertex_descriptor m_sink; 0685 0686 tQueue m_active_nodes; 0687 std::vector<bool> m_in_active_list_vec; 0688 iterator_property_map<std::vector<bool>::iterator, IndexMap> m_in_active_list_map; 0689 0690 std::list<vertex_descriptor> m_orphans; 0691 tQueue m_child_orphans; // we use a second queue for child orphans, as they are FIFO processed 0692 0693 std::vector<bool> m_has_parent_vec; 0694 iterator_property_map<std::vector<bool>::iterator, IndexMap> m_has_parent_map; 0695 0696 std::vector<long> m_time_vec; //timestamp of each node, used for sink/source-path calculations 0697 iterator_property_map<std::vector<long>::iterator, IndexMap> m_time_map; 0698 tEdgeVal m_flow; 0699 long m_time; 0700 vertex_descriptor m_last_grow_vertex; 0701 out_edge_iterator m_last_grow_edge_it; 0702 out_edge_iterator m_last_grow_edge_end; 0703 }; 0704 0705 } //namespace boost::detail 0706 0707 /** 0708 * non-named-parameter version, given everything 0709 * this is the catch all version 0710 */ 0711 template<class Graph, 0712 class CapacityEdgeMap, 0713 class ResidualCapacityEdgeMap, 0714 class ReverseEdgeMap, class PredecessorMap, 0715 class ColorMap, 0716 class DistanceMap, 0717 class IndexMap> 0718 typename property_traits<CapacityEdgeMap>::value_type 0719 boykov_kolmogorov_max_flow(Graph& g, 0720 CapacityEdgeMap cap, 0721 ResidualCapacityEdgeMap res_cap, 0722 ReverseEdgeMap rev_map, 0723 PredecessorMap pre_map, 0724 ColorMap color, 0725 DistanceMap dist, 0726 IndexMap idx, 0727 typename graph_traits<Graph>::vertex_descriptor src, 0728 typename graph_traits<Graph>::vertex_descriptor sink) 0729 { 0730 typedef typename graph_traits<Graph>::vertex_descriptor vertex_descriptor; 0731 typedef typename graph_traits<Graph>::edge_descriptor edge_descriptor; 0732 0733 //as this method is the last one before we instantiate the solver, we do the concept checks here 0734 BOOST_CONCEPT_ASSERT(( VertexListGraphConcept<Graph> )); //to have vertices(), num_vertices(), 0735 BOOST_CONCEPT_ASSERT(( EdgeListGraphConcept<Graph> )); //to have edges() 0736 BOOST_CONCEPT_ASSERT(( IncidenceGraphConcept<Graph> )); //to have source(), target() and out_edges() 0737 BOOST_CONCEPT_ASSERT(( ReadablePropertyMapConcept<CapacityEdgeMap, edge_descriptor> )); //read flow-values from edges 0738 BOOST_CONCEPT_ASSERT(( ReadWritePropertyMapConcept<ResidualCapacityEdgeMap, edge_descriptor> )); //write flow-values to residuals 0739 BOOST_CONCEPT_ASSERT(( ReadablePropertyMapConcept<ReverseEdgeMap, edge_descriptor> )); //read out reverse edges 0740 BOOST_CONCEPT_ASSERT(( ReadWritePropertyMapConcept<PredecessorMap, vertex_descriptor> )); //store predecessor there 0741 BOOST_CONCEPT_ASSERT(( ReadWritePropertyMapConcept<ColorMap, vertex_descriptor> )); //write corresponding tree 0742 BOOST_CONCEPT_ASSERT(( ReadWritePropertyMapConcept<DistanceMap, vertex_descriptor> )); //write distance to source/sink 0743 BOOST_CONCEPT_ASSERT(( ReadablePropertyMapConcept<IndexMap, vertex_descriptor> )); //get index 0...|V|-1 0744 BOOST_ASSERT(num_vertices(g) >= 2 && src != sink); 0745 0746 detail::bk_max_flow< 0747 Graph, CapacityEdgeMap, ResidualCapacityEdgeMap, ReverseEdgeMap, 0748 PredecessorMap, ColorMap, DistanceMap, IndexMap 0749 > algo(g, cap, res_cap, rev_map, pre_map, color, dist, idx, src, sink); 0750 0751 return algo.max_flow(); 0752 } 0753 0754 /** 0755 * non-named-parameter version, given capacity, residual_capacity, 0756 * reverse_edges, and an index map. 0757 */ 0758 template<class Graph, 0759 class CapacityEdgeMap, 0760 class ResidualCapacityEdgeMap, 0761 class ReverseEdgeMap, 0762 class IndexMap> 0763 typename property_traits<CapacityEdgeMap>::value_type 0764 boykov_kolmogorov_max_flow(Graph& g, 0765 CapacityEdgeMap cap, 0766 ResidualCapacityEdgeMap res_cap, 0767 ReverseEdgeMap rev, 0768 IndexMap idx, 0769 typename graph_traits<Graph>::vertex_descriptor src, 0770 typename graph_traits<Graph>::vertex_descriptor sink) 0771 { 0772 typename graph_traits<Graph>::vertices_size_type n_verts = num_vertices(g); 0773 std::vector<typename graph_traits<Graph>::edge_descriptor> predecessor_vec(n_verts); 0774 std::vector<default_color_type> color_vec(n_verts); 0775 std::vector<typename graph_traits<Graph>::vertices_size_type> distance_vec(n_verts); 0776 return 0777 boykov_kolmogorov_max_flow( 0778 g, cap, res_cap, rev, 0779 make_iterator_property_map(predecessor_vec.begin(), idx), 0780 make_iterator_property_map(color_vec.begin(), idx), 0781 make_iterator_property_map(distance_vec.begin(), idx), 0782 idx, src, sink); 0783 } 0784 0785 /** 0786 * non-named-parameter version, some given: capacity, residual_capacity, 0787 * reverse_edges, color_map and an index map. Use this if you are interested in 0788 * the minimum cut, as the color map provides that info. 0789 */ 0790 template<class Graph, 0791 class CapacityEdgeMap, 0792 class ResidualCapacityEdgeMap, 0793 class ReverseEdgeMap, 0794 class ColorMap, 0795 class IndexMap> 0796 typename property_traits<CapacityEdgeMap>::value_type 0797 boykov_kolmogorov_max_flow(Graph& g, 0798 CapacityEdgeMap cap, 0799 ResidualCapacityEdgeMap res_cap, 0800 ReverseEdgeMap rev, 0801 ColorMap color, 0802 IndexMap idx, 0803 typename graph_traits<Graph>::vertex_descriptor src, 0804 typename graph_traits<Graph>::vertex_descriptor sink) 0805 { 0806 typename graph_traits<Graph>::vertices_size_type n_verts = num_vertices(g); 0807 std::vector<typename graph_traits<Graph>::edge_descriptor> predecessor_vec(n_verts); 0808 std::vector<typename graph_traits<Graph>::vertices_size_type> distance_vec(n_verts); 0809 return 0810 boykov_kolmogorov_max_flow( 0811 g, cap, res_cap, rev, 0812 make_iterator_property_map(predecessor_vec.begin(), idx), 0813 color, 0814 make_iterator_property_map(distance_vec.begin(), idx), 0815 idx, src, sink); 0816 } 0817 0818 /** 0819 * named-parameter version, some given 0820 */ 0821 template<class Graph, class P, class T, class R> 0822 typename property_traits<typename property_map<Graph, edge_capacity_t>::const_type>::value_type 0823 boykov_kolmogorov_max_flow(Graph& g, 0824 typename graph_traits<Graph>::vertex_descriptor src, 0825 typename graph_traits<Graph>::vertex_descriptor sink, 0826 const bgl_named_params<P, T, R>& params) 0827 { 0828 return 0829 boykov_kolmogorov_max_flow( 0830 g, 0831 choose_const_pmap(get_param(params, edge_capacity), g, edge_capacity), 0832 choose_pmap(get_param(params, edge_residual_capacity), g, edge_residual_capacity), 0833 choose_const_pmap(get_param(params, edge_reverse), g, edge_reverse), 0834 choose_pmap(get_param(params, vertex_predecessor), g, vertex_predecessor), 0835 choose_pmap(get_param(params, vertex_color), g, vertex_color), 0836 choose_pmap(get_param(params, vertex_distance), g, vertex_distance), 0837 choose_const_pmap(get_param(params, vertex_index), g, vertex_index), 0838 src, sink); 0839 } 0840 0841 /** 0842 * named-parameter version, none given 0843 */ 0844 template<class Graph> 0845 typename property_traits<typename property_map<Graph, edge_capacity_t>::const_type>::value_type 0846 boykov_kolmogorov_max_flow(Graph& g, 0847 typename graph_traits<Graph>::vertex_descriptor src, 0848 typename graph_traits<Graph>::vertex_descriptor sink) 0849 { 0850 bgl_named_params<int, buffer_param_t> params(0); // bogus empty param 0851 return boykov_kolmogorov_max_flow(g, src, sink, params); 0852 } 0853 0854 } // namespace boost 0855 0856 #endif // BOOST_BOYKOV_KOLMOGOROV_MAX_FLOW_HPP 0857