Changeset 911:2914b6f0fde0 in lemon for lemon/network_simplex.h
 Timestamp:
 02/26/10 14:00:20 (13 years ago)
 Branch:
 default
 Parents:
 909:2c35bef44dd1 (diff), 910:f3bc4e9b5f3a (diff)
Note: this is a merge changeset, the changes displayed below correspond to the merge itself.
Use the (diff) links above to see all the changes relative to each parent.  Phase:
 public
 Files:

 2 edited
Legend:
 Unmodified
 Added
 Removed

lemon/network_simplex.h
r898 r911 165 165 166 166 typedef std::vector<int> IntVector; 167 typedef std::vector<char> CharVector;168 167 typedef std::vector<Value> ValueVector; 169 168 typedef std::vector<Cost> CostVector; 169 typedef std::vector<char> BoolVector; 170 // Note: vector<char> is used instead of vector<bool> for efficiency reasons 170 171 171 172 // State constants for arcs … … 214 215 IntVector _last_succ; 215 216 IntVector _dirty_revs; 216 CharVector _forward;217 CharVector _state;217 BoolVector _forward; 218 BoolVector _state; 218 219 int _root; 219 220 … … 246 247 const IntVector &_target; 247 248 const CostVector &_cost; 248 const CharVector &_state;249 const BoolVector &_state; 249 250 const CostVector &_pi; 250 251 int &_in_arc; … … 267 268 bool findEnteringArc() { 268 269 Cost c; 269 for (int e = _next_arc; e <_search_arc_num; ++e) {270 for (int e = _next_arc; e != _search_arc_num; ++e) { 270 271 c = _state[e] * (_cost[e] + _pi[_source[e]]  _pi[_target[e]]); 271 272 if (c < 0) { … … 275 276 } 276 277 } 277 for (int e = 0; e <_next_arc; ++e) {278 for (int e = 0; e != _next_arc; ++e) { 278 279 c = _state[e] * (_cost[e] + _pi[_source[e]]  _pi[_target[e]]); 279 280 if (c < 0) { … … 298 299 const IntVector &_target; 299 300 const CostVector &_cost; 300 const CharVector &_state;301 const BoolVector &_state; 301 302 const CostVector &_pi; 302 303 int &_in_arc; … … 315 316 bool findEnteringArc() { 316 317 Cost c, min = 0; 317 for (int e = 0; e <_search_arc_num; ++e) {318 for (int e = 0; e != _search_arc_num; ++e) { 318 319 c = _state[e] * (_cost[e] + _pi[_source[e]]  _pi[_target[e]]); 319 320 if (c < min) { … … 337 338 const IntVector &_target; 338 339 const CostVector &_cost; 339 const CharVector &_state;340 const BoolVector &_state; 340 341 const CostVector &_pi; 341 342 int &_in_arc; … … 356 357 { 357 358 // The main parameters of the pivot rule 358 const double BLOCK_SIZE_FACTOR = 0.5;359 const double BLOCK_SIZE_FACTOR = 1.0; 359 360 const int MIN_BLOCK_SIZE = 10; 360 361 … … 369 370 int cnt = _block_size; 370 371 int e; 371 for (e = _next_arc; e <_search_arc_num; ++e) {372 for (e = _next_arc; e != _search_arc_num; ++e) { 372 373 c = _state[e] * (_cost[e] + _pi[_source[e]]  _pi[_target[e]]); 373 374 if (c < min) { … … 380 381 } 381 382 } 382 for (e = 0; e <_next_arc; ++e) {383 for (e = 0; e != _next_arc; ++e) { 383 384 c = _state[e] * (_cost[e] + _pi[_source[e]]  _pi[_target[e]]); 384 385 if (c < min) { … … 410 411 const IntVector &_target; 411 412 const CostVector &_cost; 412 const CharVector &_state;413 const BoolVector &_state; 413 414 const CostVector &_pi; 414 415 int &_in_arc; … … 471 472 min = 0; 472 473 _curr_length = 0; 473 for (e = _next_arc; e <_search_arc_num; ++e) {474 for (e = _next_arc; e != _search_arc_num; ++e) { 474 475 c = _state[e] * (_cost[e] + _pi[_source[e]]  _pi[_target[e]]); 475 476 if (c < 0) { … … 482 483 } 483 484 } 484 for (e = 0; e <_next_arc; ++e) {485 for (e = 0; e != _next_arc; ++e) { 485 486 c = _state[e] * (_cost[e] + _pi[_source[e]]  _pi[_target[e]]); 486 487 if (c < 0) { … … 513 514 const IntVector &_target; 514 515 const CostVector &_cost; 515 const CharVector &_state;516 const BoolVector &_state; 516 517 const CostVector &_pi; 517 518 int &_in_arc; … … 566 567 // Check the current candidate list 567 568 int e; 568 for (int i = 0; i <_curr_length; ++i) {569 for (int i = 0; i != _curr_length; ++i) { 569 570 e = _candidates[i]; 570 571 _cand_cost[e] = _state[e] * … … 579 580 int limit = _head_length; 580 581 581 for (e = _next_arc; e <_search_arc_num; ++e) {582 for (e = _next_arc; e != _search_arc_num; ++e) { 582 583 _cand_cost[e] = _state[e] * 583 584 (_cost[e] + _pi[_source[e]]  _pi[_target[e]]); … … 591 592 } 592 593 } 593 for (e = 0; e <_next_arc; ++e) {594 for (e = 0; e != _next_arc; ++e) { 594 595 _cand_cost[e] = _state[e] * 595 596 (_cost[e] + _pi[_source[e]]  _pi[_target[e]]); … … 1361 1362 1362 1363 // Update _rev_thread using the new _thread values 1363 for (int i = 0; i <int(_dirty_revs.size()); ++i) {1364 for (int i = 0; i != int(_dirty_revs.size()); ++i) { 1364 1365 u = _dirty_revs[i]; 1365 1366 _rev_thread[_thread[u]] = u; … … 1431 1432 _pi[u] += sigma; 1432 1433 } 1434 } 1435 1436 // Heuristic initial pivots 1437 bool initialPivots() { 1438 Value curr, total = 0; 1439 std::vector<Node> supply_nodes, demand_nodes; 1440 for (NodeIt u(_graph); u != INVALID; ++u) { 1441 curr = _supply[_node_id[u]]; 1442 if (curr > 0) { 1443 total += curr; 1444 supply_nodes.push_back(u); 1445 } 1446 else if (curr < 0) { 1447 demand_nodes.push_back(u); 1448 } 1449 } 1450 if (_sum_supply > 0) total = _sum_supply; 1451 if (total <= 0) return true; 1452 1453 IntVector arc_vector; 1454 if (_sum_supply >= 0) { 1455 if (supply_nodes.size() == 1 && demand_nodes.size() == 1) { 1456 // Perform a reverse graph search from the sink to the source 1457 typename GR::template NodeMap<bool> reached(_graph, false); 1458 Node s = supply_nodes[0], t = demand_nodes[0]; 1459 std::vector<Node> stack; 1460 reached[t] = true; 1461 stack.push_back(t); 1462 while (!stack.empty()) { 1463 Node u, v = stack.back(); 1464 stack.pop_back(); 1465 if (v == s) break; 1466 for (InArcIt a(_graph, v); a != INVALID; ++a) { 1467 if (reached[u = _graph.source(a)]) continue; 1468 int j = _arc_id[a]; 1469 if (_cap[j] >= total) { 1470 arc_vector.push_back(j); 1471 reached[u] = true; 1472 stack.push_back(u); 1473 } 1474 } 1475 } 1476 } else { 1477 // Find the min. cost incomming arc for each demand node 1478 for (int i = 0; i != int(demand_nodes.size()); ++i) { 1479 Node v = demand_nodes[i]; 1480 Cost c, min_cost = std::numeric_limits<Cost>::max(); 1481 Arc min_arc = INVALID; 1482 for (InArcIt a(_graph, v); a != INVALID; ++a) { 1483 c = _cost[_arc_id[a]]; 1484 if (c < min_cost) { 1485 min_cost = c; 1486 min_arc = a; 1487 } 1488 } 1489 if (min_arc != INVALID) { 1490 arc_vector.push_back(_arc_id[min_arc]); 1491 } 1492 } 1493 } 1494 } else { 1495 // Find the min. cost outgoing arc for each supply node 1496 for (int i = 0; i != int(supply_nodes.size()); ++i) { 1497 Node u = supply_nodes[i]; 1498 Cost c, min_cost = std::numeric_limits<Cost>::max(); 1499 Arc min_arc = INVALID; 1500 for (OutArcIt a(_graph, u); a != INVALID; ++a) { 1501 c = _cost[_arc_id[a]]; 1502 if (c < min_cost) { 1503 min_cost = c; 1504 min_arc = a; 1505 } 1506 } 1507 if (min_arc != INVALID) { 1508 arc_vector.push_back(_arc_id[min_arc]); 1509 } 1510 } 1511 } 1512 1513 // Perform heuristic initial pivots 1514 for (int i = 0; i != int(arc_vector.size()); ++i) { 1515 in_arc = arc_vector[i]; 1516 if (_state[in_arc] * (_cost[in_arc] + _pi[_source[in_arc]]  1517 _pi[_target[in_arc]]) >= 0) continue; 1518 findJoinNode(); 1519 bool change = findLeavingArc(); 1520 if (delta >= MAX) return false; 1521 changeFlow(change); 1522 if (change) { 1523 updateTreeStructure(); 1524 updatePotential(); 1525 } 1526 } 1527 return true; 1433 1528 } 1434 1529 … … 1455 1550 PivotRuleImpl pivot(*this); 1456 1551 1552 // Perform heuristic initial pivots 1553 if (!initialPivots()) return UNBOUNDED; 1554 1457 1555 // Execute the Network Simplex algorithm 1458 1556 while (pivot.findEnteringArc()) { 
lemon/network_simplex.h
r910 r911 196 196 IntVector _source; 197 197 IntVector _target; 198 bool _arc_mixing; 198 199 199 200 // Node and arc data … … 635 636 NetworkSimplex(const GR& graph, bool arc_mixing = false) : 636 637 _graph(graph), _node_id(graph), _arc_id(graph), 638 _arc_mixing(arc_mixing), 637 639 MAX(std::numeric_limits<Value>::max()), 638 640 INF(std::numeric_limits<Value>::has_infinity ? … … 645 647 "The cost type of NetworkSimplex must be signed"); 646 648 647 // Resize vectors 648 _node_num = countNodes(_graph); 649 _arc_num = countArcs(_graph); 650 int all_node_num = _node_num + 1; 651 int max_arc_num = _arc_num + 2 * _node_num; 652 653 _source.resize(max_arc_num); 654 _target.resize(max_arc_num); 655 656 _lower.resize(_arc_num); 657 _upper.resize(_arc_num); 658 _cap.resize(max_arc_num); 659 _cost.resize(max_arc_num); 660 _supply.resize(all_node_num); 661 _flow.resize(max_arc_num); 662 _pi.resize(all_node_num); 663 664 _parent.resize(all_node_num); 665 _pred.resize(all_node_num); 666 _forward.resize(all_node_num); 667 _thread.resize(all_node_num); 668 _rev_thread.resize(all_node_num); 669 _succ_num.resize(all_node_num); 670 _last_succ.resize(all_node_num); 671 _state.resize(max_arc_num); 672 673 // Copy the graph 674 int i = 0; 675 for (NodeIt n(_graph); n != INVALID; ++n, ++i) { 676 _node_id[n] = i; 677 } 678 if (arc_mixing) { 679 // Store the arcs in a mixed order 680 int k = std::max(int(std::sqrt(double(_arc_num))), 10); 681 int i = 0, j = 0; 682 for (ArcIt a(_graph); a != INVALID; ++a) { 683 _arc_id[a] = i; 684 _source[i] = _node_id[_graph.source(a)]; 685 _target[i] = _node_id[_graph.target(a)]; 686 if ((i += k) >= _arc_num) i = ++j; 687 } 688 } else { 689 // Store the arcs in the original order 690 int i = 0; 691 for (ArcIt a(_graph); a != INVALID; ++a, ++i) { 692 _arc_id[a] = i; 693 _source[i] = _node_id[_graph.source(a)]; 694 _target[i] = _node_id[_graph.target(a)]; 695 } 696 } 697 698 // Reset parameters 649 // Reset data structures 699 650 reset(); 700 651 } … … 844 795 /// \endcode 845 796 /// 846 /// This function can be called more than once. All the parameters847 /// that have been given are kept for the next call, unless848 /// \ref reset() is called, thus only the modified parameters849 /// have to be set again. See \ref reset() for examples.850 /// However, the underlying digraph must not be modified after this851 /// class have been constructed, since it copies and extends the graph.797 /// This function can be called more than once. All the given parameters 798 /// are kept for the next call, unless \ref resetParams() or \ref reset() 799 /// is used, thus only the modified parameters have to be set again. 800 /// If the underlying digraph was also modified after the construction 801 /// of the class (or the last \ref reset() call), then the \ref reset() 802 /// function must be called. 852 803 /// 853 804 /// \param pivot_rule The pivot rule that will be used during the … … 863 814 /// 864 815 /// \see ProblemType, PivotRule 816 /// \see resetParams(), reset() 865 817 ProblemType run(PivotRule pivot_rule = BLOCK_SEARCH) { 866 818 if (!init()) return INFEASIBLE; … … 874 826 /// \ref costMap(), \ref supplyMap(), \ref stSupply(), \ref supplyType(). 875 827 /// 876 /// It is useful for multiple run() calls. If this function is not 877 /// used, all the parameters given before are kept for the next 878 /// \ref run() call. 879 /// However, the underlying digraph must not be modified after this 880 /// class have been constructed, since it copies and extends the graph. 828 /// It is useful for multiple \ref run() calls. Basically, all the given 829 /// parameters are kept for the next \ref run() call, unless 830 /// \ref resetParams() or \ref reset() is used. 831 /// If the underlying digraph was also modified after the construction 832 /// of the class or the last \ref reset() call, then the \ref reset() 833 /// function must be used, otherwise \ref resetParams() is sufficient. 881 834 /// 882 835 /// For example, … … 888 841 /// .supplyMap(sup).run(); 889 842 /// 890 /// // Run again with modified cost map (reset () is not called,843 /// // Run again with modified cost map (resetParams() is not called, 891 844 /// // so only the cost map have to be set again) 892 845 /// cost[e] += 100; 893 846 /// ns.costMap(cost).run(); 894 847 /// 895 /// // Run again from scratch using reset ()848 /// // Run again from scratch using resetParams() 896 849 /// // (the lower bounds will be set to zero on all arcs) 897 /// ns.reset ();850 /// ns.resetParams(); 898 851 /// ns.upperMap(capacity).costMap(cost) 899 852 /// .supplyMap(sup).run(); … … 901 854 /// 902 855 /// \return <tt>(*this)</tt> 903 NetworkSimplex& reset() { 856 /// 857 /// \see reset(), run() 858 NetworkSimplex& resetParams() { 904 859 for (int i = 0; i != _node_num; ++i) { 905 860 _supply[i] = 0; … … 915 870 } 916 871 872 /// \brief Reset the internal data structures and all the parameters 873 /// that have been given before. 874 /// 875 /// This function resets the internal data structures and all the 876 /// paramaters that have been given before using functions \ref lowerMap(), 877 /// \ref upperMap(), \ref costMap(), \ref supplyMap(), \ref stSupply(), 878 /// \ref supplyType(). 879 /// 880 /// It is useful for multiple \ref run() calls. Basically, all the given 881 /// parameters are kept for the next \ref run() call, unless 882 /// \ref resetParams() or \ref reset() is used. 883 /// If the underlying digraph was also modified after the construction 884 /// of the class or the last \ref reset() call, then the \ref reset() 885 /// function must be used, otherwise \ref resetParams() is sufficient. 886 /// 887 /// See \ref resetParams() for examples. 888 /// 889 /// \return <tt>(*this)</tt> 890 /// 891 /// \see resetParams(), run() 892 NetworkSimplex& reset() { 893 // Resize vectors 894 _node_num = countNodes(_graph); 895 _arc_num = countArcs(_graph); 896 int all_node_num = _node_num + 1; 897 int max_arc_num = _arc_num + 2 * _node_num; 898 899 _source.resize(max_arc_num); 900 _target.resize(max_arc_num); 901 902 _lower.resize(_arc_num); 903 _upper.resize(_arc_num); 904 _cap.resize(max_arc_num); 905 _cost.resize(max_arc_num); 906 _supply.resize(all_node_num); 907 _flow.resize(max_arc_num); 908 _pi.resize(all_node_num); 909 910 _parent.resize(all_node_num); 911 _pred.resize(all_node_num); 912 _forward.resize(all_node_num); 913 _thread.resize(all_node_num); 914 _rev_thread.resize(all_node_num); 915 _succ_num.resize(all_node_num); 916 _last_succ.resize(all_node_num); 917 _state.resize(max_arc_num); 918 919 // Copy the graph 920 int i = 0; 921 for (NodeIt n(_graph); n != INVALID; ++n, ++i) { 922 _node_id[n] = i; 923 } 924 if (_arc_mixing) { 925 // Store the arcs in a mixed order 926 int k = std::max(int(std::sqrt(double(_arc_num))), 10); 927 int i = 0, j = 0; 928 for (ArcIt a(_graph); a != INVALID; ++a) { 929 _arc_id[a] = i; 930 _source[i] = _node_id[_graph.source(a)]; 931 _target[i] = _node_id[_graph.target(a)]; 932 if ((i += k) >= _arc_num) i = ++j; 933 } 934 } else { 935 // Store the arcs in the original order 936 int i = 0; 937 for (ArcIt a(_graph); a != INVALID; ++a, ++i) { 938 _arc_id[a] = i; 939 _source[i] = _node_id[_graph.source(a)]; 940 _target[i] = _node_id[_graph.target(a)]; 941 } 942 } 943 944 // Reset parameters 945 resetParams(); 946 return *this; 947 } 948 917 949 /// @} 918 950
Note: See TracChangeset
for help on using the changeset viewer.