COIN-OR::LEMON - Graph Library

Changes in / [845:a5fc1e1e5039:844:a6eb9698c321] in lemon-main


Ignore:
Files:
11 edited

Legend:

Unmodified
Added
Removed
  • demo/arg_parser_demo.cc

    r842 r440  
    6666    .other("...");
    6767
    68   // Throw an exception when problems occurs. The default behavior is to
    69   // exit(1) on these cases, but this makes Valgrind falsely warn
    70   // about memory leaks.
    71   ap.throwOnProblems();
    72  
    7368  // Perform the parsing process
    7469  // (in case of any error it terminates the program)
    75   // The try {} construct is necessary only if the ap.trowOnProblems()
    76   // setting is in use.
    77   try {
    78     ap.parse();
    79   } catch (ArgParserException &) { return 1; }
     70  ap.parse();
    8071
    8172  // Check each option if it has been given and print its value
  • lemon/arg_parser.cc

    r842 r440  
    2121namespace lemon {
    2222
    23   void ArgParser::_terminate(ArgParserException::Reason reason) const
    24   {
    25     if(_exit_on_problems)
    26       exit(1);
    27     else throw(ArgParserException(reason));
    28   }
    29  
    30  
    3123  void ArgParser::_showHelp(void *p)
    3224  {
    3325    (static_cast<ArgParser*>(p))->showHelp();
    34     (static_cast<ArgParser*>(p))->_terminate(ArgParserException::HELP);
     26    exit(1);
    3527  }
    3628
    3729  ArgParser::ArgParser(int argc, const char * const *argv)
    38     :_argc(argc), _argv(argv), _command_name(argv[0]),
    39     _exit_on_problems(true) {
     30    :_argc(argc), _argv(argv), _command_name(argv[0]) {
    4031    funcOption("-help","Print a short help message",_showHelp,this);
    4132    synonym("help","-help");
     
    352343        i!=_others_help.end();++i) showHelp(i);
    353344    for(Opts::const_iterator i=_opts.begin();i!=_opts.end();++i) showHelp(i);
    354     _terminate(ArgParserException::HELP);
     345    exit(1);
    355346  }
    356347
     
    361352    std::cerr << "\nType '" << _command_name <<
    362353      " --help' to obtain a short summary on the usage.\n\n";
    363     _terminate(ArgParserException::UNKNOWN_OPT);
     354    exit(1);
    364355  }
    365356
     
    424415      std::cerr << "\nType '" << _command_name <<
    425416        " --help' to obtain a short summary on the usage.\n\n";
    426       _terminate(ArgParserException::INVALID_OPT);
     417      exit(1);
    427418    }
    428419  }
  • lemon/arg_parser.h

    r842 r440  
    3434
    3535namespace lemon {
    36 
    37   ///Exception used by ArgParser
    38   class ArgParserException : public Exception {
    39   public:
    40     enum Reason {
    41       HELP,         /// <tt>--help</tt> option was given
    42       UNKNOWN_OPT,  /// Unknown option was given
    43       INVALID_OPT   /// Invalid combination of options
    44     };
    45    
    46   private:
    47     Reason _reason;
    48    
    49   public:
    50     ///Constructor
    51     ArgParserException(Reason r) throw() : _reason(r) {}
    52     ///Virtual destructor
    53     virtual ~ArgParserException() throw() {}
    54     ///A short description of the exception
    55     virtual const char* what() const throw() {
    56       switch(_reason)
    57         {
    58         case HELP:
    59           return "lemon::ArgParseException: ask for help";
    60           break;
    61         case UNKNOWN_OPT:
    62           return "lemon::ArgParseException: unknown option";
    63           break;
    64         case INVALID_OPT:
    65           return "lemon::ArgParseException: invalid combination of options";
    66           break;
    67         }
    68       return "";
    69     }
    70     ///Return the reason for the failure
    71     Reason reason() const {return _reason; }
    72   };
    73 
    7436
    7537  ///Command line arguments parser
     
    142104    std::string _command_name;
    143105
    144    
     106
    145107  private:
    146108    //Bind a function to an option.
     
    154116                    const std::string &help,
    155117                    void (*func)(void *),void *data);
    156 
    157     bool _exit_on_problems;
    158    
    159     void _terminate(ArgParserException::Reason reason) const;
    160118
    161119  public:
     
    423381    const std::vector<std::string> &files() const { return _file_args; }
    424382
    425     ///Throw instead of exit in case of problems
    426     void throwOnProblems()
    427     {
    428       _exit_on_problems=false;
    429     }
    430383  };
    431384}
  • lemon/capacity_scaling.h

    r840 r831  
    140140
    141141    typedef std::vector<int> IntVector;
     142    typedef std::vector<char> BoolVector;
    142143    typedef std::vector<Value> ValueVector;
    143144    typedef std::vector<Cost> CostVector;
    144     typedef std::vector<char> BoolVector;
    145     // Note: vector<char> is used instead of vector<bool> for efficiency reasons
    146145
    147146  private:
     
    800799      if (_factor > 1) {
    801800        // With scaling
    802         Value max_sup = 0, max_dem = 0, max_cap = 0;
    803         for (int i = 0; i != _root; ++i) {
     801        Value max_sup = 0, max_dem = 0;
     802        for (int i = 0; i != _node_num; ++i) {
    804803          Value ex = _excess[i];
    805804          if ( ex > max_sup) max_sup =  ex;
    806805          if (-ex > max_dem) max_dem = -ex;
    807           int last_out = _first_out[i+1] - 1;
    808           for (int j = _first_out[i]; j != last_out; ++j) {
    809             if (_res_cap[j] > max_cap) max_cap = _res_cap[j];
    810           }
     806        }
     807        Value max_cap = 0;
     808        for (int j = 0; j != _res_arc_num; ++j) {
     809          if (_res_cap[j] > max_cap) max_cap = _res_cap[j];
    811810        }
    812811        max_sup = std::min(std::min(max_sup, max_dem), max_cap);
  • lemon/cost_scaling.h

    r840 r831  
    202202
    203203    typedef std::vector<int> IntVector;
     204    typedef std::vector<char> BoolVector;
    204205    typedef std::vector<Value> ValueVector;
    205206    typedef std::vector<Cost> CostVector;
    206207    typedef std::vector<LargeCost> LargeCostVector;
    207     typedef std::vector<char> BoolVector;
    208     // Note: vector<char> is used instead of vector<bool> for efficiency reasons
    209208
    210209  private:
     
    250249    bool _have_lower;
    251250    Value _sum_supply;
    252     int _sup_node_num;
    253251
    254252    // Data structures for storing the digraph
     
    279277    int _alpha;
    280278
    281     IntVector _buckets;
    282     IntVector _bucket_next;
    283     IntVector _bucket_prev;
    284     IntVector _rank;
    285     int _max_rank;
    286  
    287279    // Data for a StaticDigraph structure
    288280    typedef std::pair<int, int> IntPair;
     
    837829      }
    838830
    839       _sup_node_num = 0;
    840       for (NodeIt n(_graph); n != INVALID; ++n) {
    841         if (sup[n] > 0) ++_sup_node_num;
    842       }
    843 
    844831      // Find a feasible flow using Circulation
    845832      Circulation<Digraph, ConstMap<Arc, Value>, ValueArcMap, ValueNodeMap>
     
    876863        for (int a = _first_out[_root]; a != _res_arc_num; ++a) {
    877864          int ra = _reverse[a];
    878           _res_cap[a] = 0;
     865          _res_cap[a] = 1;
    879866          _res_cap[ra] = 0;
    880867          _cost[a] = 0;
     
    890877      // Maximum path length for partial augment
    891878      const int MAX_PATH_LENGTH = 4;
    892 
    893       // Initialize data structures for buckets     
    894       _max_rank = _alpha * _res_node_num;
    895       _buckets.resize(_max_rank);
    896       _bucket_next.resize(_res_node_num + 1);
    897       _bucket_prev.resize(_res_node_num + 1);
    898       _rank.resize(_res_node_num + 1);
    899  
     879     
    900880      // Execute the algorithm
    901881      switch (method) {
     
    936916      }
    937917    }
    938    
    939     // Initialize a cost scaling phase
    940     void initPhase() {
    941       // Saturate arcs not satisfying the optimality condition
    942       for (int u = 0; u != _res_node_num; ++u) {
    943         int last_out = _first_out[u+1];
    944         LargeCost pi_u = _pi[u];
    945         for (int a = _first_out[u]; a != last_out; ++a) {
    946           int v = _target[a];
    947           if (_res_cap[a] > 0 && _cost[a] + pi_u - _pi[v] < 0) {
     918
     919    /// Execute the algorithm performing augment and relabel operations
     920    void startAugment(int max_length = std::numeric_limits<int>::max()) {
     921      // Paramters for heuristics
     922      const int BF_HEURISTIC_EPSILON_BOUND = 1000;
     923      const int BF_HEURISTIC_BOUND_FACTOR  = 3;
     924
     925      // Perform cost scaling phases
     926      IntVector pred_arc(_res_node_num);
     927      std::vector<int> path_nodes;
     928      for ( ; _epsilon >= 1; _epsilon = _epsilon < _alpha && _epsilon > 1 ?
     929                                        1 : _epsilon / _alpha )
     930      {
     931        // "Early Termination" heuristic: use Bellman-Ford algorithm
     932        // to check if the current flow is optimal
     933        if (_epsilon <= BF_HEURISTIC_EPSILON_BOUND) {
     934          _arc_vec.clear();
     935          _cost_vec.clear();
     936          for (int j = 0; j != _res_arc_num; ++j) {
     937            if (_res_cap[j] > 0) {
     938              _arc_vec.push_back(IntPair(_source[j], _target[j]));
     939              _cost_vec.push_back(_cost[j] + 1);
     940            }
     941          }
     942          _sgr.build(_res_node_num, _arc_vec.begin(), _arc_vec.end());
     943
     944          BellmanFord<StaticDigraph, LargeCostArcMap> bf(_sgr, _cost_map);
     945          bf.init(0);
     946          bool done = false;
     947          int K = int(BF_HEURISTIC_BOUND_FACTOR * sqrt(_res_node_num));
     948          for (int i = 0; i < K && !done; ++i)
     949            done = bf.processNextWeakRound();
     950          if (done) break;
     951        }
     952
     953        // Saturate arcs not satisfying the optimality condition
     954        for (int a = 0; a != _res_arc_num; ++a) {
     955          if (_res_cap[a] > 0 &&
     956              _cost[a] + _pi[_source[a]] - _pi[_target[a]] < 0) {
    948957            Value delta = _res_cap[a];
    949             _excess[u] -= delta;
    950             _excess[v] += delta;
     958            _excess[_source[a]] -= delta;
     959            _excess[_target[a]] += delta;
    951960            _res_cap[a] = 0;
    952961            _res_cap[_reverse[a]] += delta;
    953962          }
    954963        }
    955       }
    956      
    957       // Find active nodes (i.e. nodes with positive excess)
    958       for (int u = 0; u != _res_node_num; ++u) {
    959         if (_excess[u] > 0) _active_nodes.push_back(u);
    960       }
    961 
    962       // Initialize the next arcs
    963       for (int u = 0; u != _res_node_num; ++u) {
    964         _next_out[u] = _first_out[u];
    965       }
    966     }
    967    
    968     // Early termination heuristic
    969     bool earlyTermination() {
    970       const double EARLY_TERM_FACTOR = 3.0;
    971 
    972       // Build a static residual graph
    973       _arc_vec.clear();
    974       _cost_vec.clear();
    975       for (int j = 0; j != _res_arc_num; ++j) {
    976         if (_res_cap[j] > 0) {
    977           _arc_vec.push_back(IntPair(_source[j], _target[j]));
    978           _cost_vec.push_back(_cost[j] + 1);
    979         }
    980       }
    981       _sgr.build(_res_node_num, _arc_vec.begin(), _arc_vec.end());
    982 
    983       // Run Bellman-Ford algorithm to check if the current flow is optimal
    984       BellmanFord<StaticDigraph, LargeCostArcMap> bf(_sgr, _cost_map);
    985       bf.init(0);
    986       bool done = false;
    987       int K = int(EARLY_TERM_FACTOR * std::sqrt(double(_res_node_num)));
    988       for (int i = 0; i < K && !done; ++i) {
    989         done = bf.processNextWeakRound();
    990       }
    991       return done;
    992     }
    993 
    994     // Global potential update heuristic
    995     void globalUpdate() {
    996       int bucket_end = _root + 1;
    997    
    998       // Initialize buckets
    999       for (int r = 0; r != _max_rank; ++r) {
    1000         _buckets[r] = bucket_end;
    1001       }
    1002       Value total_excess = 0;
    1003       for (int i = 0; i != _res_node_num; ++i) {
    1004         if (_excess[i] < 0) {
    1005           _rank[i] = 0;
    1006           _bucket_next[i] = _buckets[0];
    1007           _bucket_prev[_buckets[0]] = i;
    1008           _buckets[0] = i;
    1009         } else {
    1010           total_excess += _excess[i];
    1011           _rank[i] = _max_rank;
    1012         }
    1013       }
    1014       if (total_excess == 0) return;
    1015 
    1016       // Search the buckets
    1017       int r = 0;
    1018       for ( ; r != _max_rank; ++r) {
    1019         while (_buckets[r] != bucket_end) {
    1020           // Remove the first node from the current bucket
    1021           int u = _buckets[r];
    1022           _buckets[r] = _bucket_next[u];
    1023          
    1024           // Search the incomming arcs of u
    1025           LargeCost pi_u = _pi[u];
    1026           int last_out = _first_out[u+1];
    1027           for (int a = _first_out[u]; a != last_out; ++a) {
    1028             int ra = _reverse[a];
    1029             if (_res_cap[ra] > 0) {
    1030               int v = _source[ra];
    1031               int old_rank_v = _rank[v];
    1032               if (r < old_rank_v) {
    1033                 // Compute the new rank of v
    1034                 LargeCost nrc = (_cost[ra] + _pi[v] - pi_u) / _epsilon;
    1035                 int new_rank_v = old_rank_v;
    1036                 if (nrc < LargeCost(_max_rank))
    1037                   new_rank_v = r + 1 + int(nrc);
    1038                  
    1039                 // Change the rank of v
    1040                 if (new_rank_v < old_rank_v) {
    1041                   _rank[v] = new_rank_v;
    1042                   _next_out[v] = _first_out[v];
    1043                  
    1044                   // Remove v from its old bucket
    1045                   if (old_rank_v < _max_rank) {
    1046                     if (_buckets[old_rank_v] == v) {
    1047                       _buckets[old_rank_v] = _bucket_next[v];
    1048                     } else {
    1049                       _bucket_next[_bucket_prev[v]] = _bucket_next[v];
    1050                       _bucket_prev[_bucket_next[v]] = _bucket_prev[v];
    1051                     }
    1052                   }
    1053                  
    1054                   // Insert v to its new bucket
    1055                   _bucket_next[v] = _buckets[new_rank_v];
    1056                   _bucket_prev[_buckets[new_rank_v]] = v;
    1057                   _buckets[new_rank_v] = v;
    1058                 }
    1059               }
    1060             }
    1061           }
    1062 
    1063           // Finish search if there are no more active nodes
    1064           if (_excess[u] > 0) {
    1065             total_excess -= _excess[u];
    1066             if (total_excess <= 0) break;
    1067           }
    1068         }
    1069         if (total_excess <= 0) break;
    1070       }
    1071      
    1072       // Relabel nodes
    1073       for (int u = 0; u != _res_node_num; ++u) {
    1074         int k = std::min(_rank[u], r);
    1075         if (k > 0) {
    1076           _pi[u] -= _epsilon * k;
     964       
     965        // Find active nodes (i.e. nodes with positive excess)
     966        for (int u = 0; u != _res_node_num; ++u) {
     967          if (_excess[u] > 0) _active_nodes.push_back(u);
     968        }
     969
     970        // Initialize the next arcs
     971        for (int u = 0; u != _res_node_num; ++u) {
    1077972          _next_out[u] = _first_out[u];
    1078973        }
    1079       }
    1080     }
    1081 
    1082     /// Execute the algorithm performing augment and relabel operations
    1083     void startAugment(int max_length = std::numeric_limits<int>::max()) {
    1084       // Paramters for heuristics
    1085       const int EARLY_TERM_EPSILON_LIMIT = 1000;
    1086       const double GLOBAL_UPDATE_FACTOR = 3.0;
    1087 
    1088       const int global_update_freq = int(GLOBAL_UPDATE_FACTOR *
    1089         (_res_node_num + _sup_node_num * _sup_node_num));
    1090       int next_update_limit = global_update_freq;
    1091      
    1092       int relabel_cnt = 0;
    1093      
    1094       // Perform cost scaling phases
    1095       std::vector<int> path;
    1096       for ( ; _epsilon >= 1; _epsilon = _epsilon < _alpha && _epsilon > 1 ?
    1097                                         1 : _epsilon / _alpha )
    1098       {
    1099         // Early termination heuristic
    1100         if (_epsilon <= EARLY_TERM_EPSILON_LIMIT) {
    1101           if (earlyTermination()) break;
    1102         }
    1103        
    1104         // Initialize current phase
    1105         initPhase();
    1106        
     974
    1107975        // Perform partial augment and relabel operations
    1108976        while (true) {
     
    1114982          if (_active_nodes.size() == 0) break;
    1115983          int start = _active_nodes.front();
     984          path_nodes.clear();
     985          path_nodes.push_back(start);
    1116986
    1117987          // Find an augmenting path from the start node
    1118           path.clear();
    1119988          int tip = start;
    1120           while (_excess[tip] >= 0 && int(path.size()) < max_length) {
     989          while (_excess[tip] >= 0 &&
     990                 int(path_nodes.size()) <= max_length) {
    1121991            int u;
    1122             LargeCost min_red_cost, rc, pi_tip = _pi[tip];
    1123             int last_out = _first_out[tip+1];
     992            LargeCost min_red_cost, rc;
     993            int last_out = _sum_supply < 0 ?
     994              _first_out[tip+1] : _first_out[tip+1] - 1;
    1124995            for (int a = _next_out[tip]; a != last_out; ++a) {
    1125               u = _target[a];
    1126               if (_res_cap[a] > 0 && _cost[a] + pi_tip - _pi[u] < 0) {
    1127                 path.push_back(a);
     996              if (_res_cap[a] > 0 &&
     997                  _cost[a] + _pi[_source[a]] - _pi[_target[a]] < 0) {
     998                u = _target[a];
     999                pred_arc[u] = a;
    11281000                _next_out[tip] = a;
    11291001                tip = u;
     1002                path_nodes.push_back(tip);
    11301003                goto next_step;
    11311004              }
     
    11331006
    11341007            // Relabel tip node
    1135             min_red_cost = std::numeric_limits<LargeCost>::max();
    1136             if (tip != start) {
    1137               int ra = _reverse[path.back()];
    1138               min_red_cost = _cost[ra] + pi_tip - _pi[_target[ra]];
    1139             }
     1008            min_red_cost = std::numeric_limits<LargeCost>::max() / 2;
    11401009            for (int a = _first_out[tip]; a != last_out; ++a) {
    1141               rc = _cost[a] + pi_tip - _pi[_target[a]];
     1010              rc = _cost[a] + _pi[_source[a]] - _pi[_target[a]];
    11421011              if (_res_cap[a] > 0 && rc < min_red_cost) {
    11431012                min_red_cost = rc;
     
    11451014            }
    11461015            _pi[tip] -= min_red_cost + _epsilon;
     1016
     1017            // Reset the next arc of tip
    11471018            _next_out[tip] = _first_out[tip];
    1148             ++relabel_cnt;
    11491019
    11501020            // Step back
    11511021            if (tip != start) {
    1152               tip = _source[path.back()];
    1153               path.pop_back();
     1022              path_nodes.pop_back();
     1023              tip = path_nodes.back();
    11541024            }
    11551025
     
    11591029          // Augment along the found path (as much flow as possible)
    11601030          Value delta;
    1161           int pa, u, v = start;
    1162           for (int i = 0; i != int(path.size()); ++i) {
    1163             pa = path[i];
     1031          int u, v = path_nodes.front(), pa;
     1032          for (int i = 1; i < int(path_nodes.size()); ++i) {
    11641033            u = v;
    1165             v = _target[pa];
     1034            v = path_nodes[i];
     1035            pa = pred_arc[v];
    11661036            delta = std::min(_res_cap[pa], _excess[u]);
    11671037            _res_cap[pa] -= delta;
     
    11721042              _active_nodes.push_back(v);
    11731043          }
    1174 
    1175           // Global update heuristic
    1176           if (relabel_cnt >= next_update_limit) {
    1177             globalUpdate();
    1178             next_update_limit += global_update_freq;
    1179           }
    11801044        }
    11811045      }
     
    11851049    void startPush() {
    11861050      // Paramters for heuristics
    1187       const int EARLY_TERM_EPSILON_LIMIT = 1000;
    1188       const double GLOBAL_UPDATE_FACTOR = 2.0;
    1189 
    1190       const int global_update_freq = int(GLOBAL_UPDATE_FACTOR *
    1191         (_res_node_num + _sup_node_num * _sup_node_num));
    1192       int next_update_limit = global_update_freq;
    1193 
    1194       int relabel_cnt = 0;
    1195      
     1051      const int BF_HEURISTIC_EPSILON_BOUND = 1000;
     1052      const int BF_HEURISTIC_BOUND_FACTOR  = 3;
     1053
    11961054      // Perform cost scaling phases
    11971055      BoolVector hyper(_res_node_num, false);
    1198       LargeCostVector hyper_cost(_res_node_num);
    11991056      for ( ; _epsilon >= 1; _epsilon = _epsilon < _alpha && _epsilon > 1 ?
    12001057                                        1 : _epsilon / _alpha )
    12011058      {
    1202         // Early termination heuristic
    1203         if (_epsilon <= EARLY_TERM_EPSILON_LIMIT) {
    1204           if (earlyTermination()) break;
    1205         }
    1206        
    1207         // Initialize current phase
    1208         initPhase();
     1059        // "Early Termination" heuristic: use Bellman-Ford algorithm
     1060        // to check if the current flow is optimal
     1061        if (_epsilon <= BF_HEURISTIC_EPSILON_BOUND) {
     1062          _arc_vec.clear();
     1063          _cost_vec.clear();
     1064          for (int j = 0; j != _res_arc_num; ++j) {
     1065            if (_res_cap[j] > 0) {
     1066              _arc_vec.push_back(IntPair(_source[j], _target[j]));
     1067              _cost_vec.push_back(_cost[j] + 1);
     1068            }
     1069          }
     1070          _sgr.build(_res_node_num, _arc_vec.begin(), _arc_vec.end());
     1071
     1072          BellmanFord<StaticDigraph, LargeCostArcMap> bf(_sgr, _cost_map);
     1073          bf.init(0);
     1074          bool done = false;
     1075          int K = int(BF_HEURISTIC_BOUND_FACTOR * sqrt(_res_node_num));
     1076          for (int i = 0; i < K && !done; ++i)
     1077            done = bf.processNextWeakRound();
     1078          if (done) break;
     1079        }
     1080
     1081        // Saturate arcs not satisfying the optimality condition
     1082        for (int a = 0; a != _res_arc_num; ++a) {
     1083          if (_res_cap[a] > 0 &&
     1084              _cost[a] + _pi[_source[a]] - _pi[_target[a]] < 0) {
     1085            Value delta = _res_cap[a];
     1086            _excess[_source[a]] -= delta;
     1087            _excess[_target[a]] += delta;
     1088            _res_cap[a] = 0;
     1089            _res_cap[_reverse[a]] += delta;
     1090          }
     1091        }
     1092
     1093        // Find active nodes (i.e. nodes with positive excess)
     1094        for (int u = 0; u != _res_node_num; ++u) {
     1095          if (_excess[u] > 0) _active_nodes.push_back(u);
     1096        }
     1097
     1098        // Initialize the next arcs
     1099        for (int u = 0; u != _res_node_num; ++u) {
     1100          _next_out[u] = _first_out[u];
     1101        }
    12091102
    12101103        // Perform push and relabel operations
    12111104        while (_active_nodes.size() > 0) {
    1212           LargeCost min_red_cost, rc, pi_n;
     1105          LargeCost min_red_cost, rc;
    12131106          Value delta;
    12141107          int n, t, a, last_out = _res_arc_num;
    12151108
     1109          // Select an active node (FIFO selection)
    12161110        next_node:
    1217           // Select an active node (FIFO selection)
    12181111          n = _active_nodes.front();
    1219           last_out = _first_out[n+1];
    1220           pi_n = _pi[n];
    1221          
     1112          last_out = _sum_supply < 0 ?
     1113            _first_out[n+1] : _first_out[n+1] - 1;
     1114
    12221115          // Perform push operations if there are admissible arcs
    12231116          if (_excess[n] > 0) {
    12241117            for (a = _next_out[n]; a != last_out; ++a) {
    12251118              if (_res_cap[a] > 0 &&
    1226                   _cost[a] + pi_n - _pi[_target[a]] < 0) {
     1119                  _cost[a] + _pi[_source[a]] - _pi[_target[a]] < 0) {
    12271120                delta = std::min(_res_cap[a], _excess[n]);
    12281121                t = _target[a];
     
    12301123                // Push-look-ahead heuristic
    12311124                Value ahead = -_excess[t];
    1232                 int last_out_t = _first_out[t+1];
    1233                 LargeCost pi_t = _pi[t];
     1125                int last_out_t = _sum_supply < 0 ?
     1126                  _first_out[t+1] : _first_out[t+1] - 1;
    12341127                for (int ta = _next_out[t]; ta != last_out_t; ++ta) {
    12351128                  if (_res_cap[ta] > 0 &&
    1236                       _cost[ta] + pi_t - _pi[_target[ta]] < 0)
     1129                      _cost[ta] + _pi[_source[ta]] - _pi[_target[ta]] < 0)
    12371130                    ahead += _res_cap[ta];
    12381131                  if (ahead >= delta) break;
     
    12411134
    12421135                // Push flow along the arc
    1243                 if (ahead < delta && !hyper[t]) {
     1136                if (ahead < delta) {
    12441137                  _res_cap[a] -= ahead;
    12451138                  _res_cap[_reverse[a]] += ahead;
     
    12481141                  _active_nodes.push_front(t);
    12491142                  hyper[t] = true;
    1250                   hyper_cost[t] = _cost[a] + pi_n - pi_t;
    12511143                  _next_out[n] = a;
    12521144                  goto next_node;
     
    12711163          // Relabel the node if it is still active (or hyper)
    12721164          if (_excess[n] > 0 || hyper[n]) {
    1273              min_red_cost = hyper[n] ? -hyper_cost[n] :
    1274                std::numeric_limits<LargeCost>::max();
     1165            min_red_cost = std::numeric_limits<LargeCost>::max() / 2;
    12751166            for (int a = _first_out[n]; a != last_out; ++a) {
    1276               rc = _cost[a] + pi_n - _pi[_target[a]];
     1167              rc = _cost[a] + _pi[_source[a]] - _pi[_target[a]];
    12771168              if (_res_cap[a] > 0 && rc < min_red_cost) {
    12781169                min_red_cost = rc;
     
    12801171            }
    12811172            _pi[n] -= min_red_cost + _epsilon;
     1173            hyper[n] = false;
     1174
     1175            // Reset the next arc
    12821176            _next_out[n] = _first_out[n];
    1283             hyper[n] = false;
    1284             ++relabel_cnt;
    12851177          }
    12861178       
     
    12921184            _active_nodes.pop_front();
    12931185          }
    1294          
    1295           // Global update heuristic
    1296           if (relabel_cnt >= next_update_limit) {
    1297             globalUpdate();
    1298             for (int u = 0; u != _res_node_num; ++u)
    1299               hyper[u] = false;
    1300             next_update_limit += global_update_freq;
    1301           }
    13021186        }
    13031187      }
  • lemon/cycle_canceling.h

    r840 r830  
    145145   
    146146    typedef std::vector<int> IntVector;
     147    typedef std::vector<char> CharVector;
    147148    typedef std::vector<double> DoubleVector;
    148149    typedef std::vector<Value> ValueVector;
    149150    typedef std::vector<Cost> CostVector;
    150     typedef std::vector<char> BoolVector;
    151     // Note: vector<char> is used instead of vector<bool> for efficiency reasons
    152151
    153152  private:
     
    200199    IntArcMap _arc_idb;
    201200    IntVector _first_out;
    202     BoolVector _forward;
     201    CharVector _forward;
    203202    IntVector _source;
    204203    IntVector _target;
     
    964963      DoubleVector pi(_res_node_num, 0.0);
    965964      IntVector level(_res_node_num);
    966       BoolVector reached(_res_node_num);
    967       BoolVector processed(_res_node_num);
     965      CharVector reached(_res_node_num);
     966      CharVector processed(_res_node_num);
    968967      IntVector pred_node(_res_node_num);
    969968      IntVector pred_arc(_res_node_num);
  • lemon/graph_to_eps.h

    r838 r786  
    685685#else
    686686      os << bits::getWinFormattedDate();
    687       os << std::endl;
    688687#endif
    689688    }
     689    os << std::endl;
    690690
    691691    if (_autoArcWidthScale) {
  • lemon/hartmann_orlin.h

    r841 r825  
    406406    /// \pre \ref run() or \ref findMinMean() must be called before
    407407    /// using this function.
    408     Value cycleLength() const {
    409       return static_cast<Value>(_best_length);
     408    LargeValue cycleLength() const {
     409      return _best_length;
    410410    }
    411411
  • lemon/howard.h

    r841 r825  
    385385    /// \pre \ref run() or \ref findMinMean() must be called before
    386386    /// using this function.
    387     Value cycleLength() const {
    388       return static_cast<Value>(_best_length);
     387    LargeValue cycleLength() const {
     388      return _best_length;
    389389    }
    390390
  • lemon/karp.h

    r841 r825  
    393393    /// \pre \ref run() or \ref findMinMean() must be called before
    394394    /// using this function.
    395     Value cycleLength() const {
    396       return static_cast<Value>(_cycle_length);
     395    LargeValue cycleLength() const {
     396      return _cycle_length;
    397397    }
    398398
  • lemon/network_simplex.h

    r840 r830  
    165165
    166166    typedef std::vector<int> IntVector;
     167    typedef std::vector<char> CharVector;
    167168    typedef std::vector<Value> ValueVector;
    168169    typedef std::vector<Cost> CostVector;
    169     typedef std::vector<char> BoolVector;
    170     // Note: vector<char> is used instead of vector<bool> for efficiency reasons
    171170
    172171    // State constants for arcs
     
    215214    IntVector _last_succ;
    216215    IntVector _dirty_revs;
    217     BoolVector _forward;
    218     BoolVector _state;
     216    CharVector _forward;
     217    CharVector _state;
    219218    int _root;
    220219
     
    247246      const IntVector  &_target;
    248247      const CostVector &_cost;
    249       const BoolVector &_state;
     248      const CharVector &_state;
    250249      const CostVector &_pi;
    251250      int &_in_arc;
     
    268267      bool findEnteringArc() {
    269268        Cost c;
    270         for (int e = _next_arc; e != _search_arc_num; ++e) {
     269        for (int e = _next_arc; e < _search_arc_num; ++e) {
    271270          c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
    272271          if (c < 0) {
     
    276275          }
    277276        }
    278         for (int e = 0; e != _next_arc; ++e) {
     277        for (int e = 0; e < _next_arc; ++e) {
    279278          c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
    280279          if (c < 0) {
     
    299298      const IntVector  &_target;
    300299      const CostVector &_cost;
    301       const BoolVector &_state;
     300      const CharVector &_state;
    302301      const CostVector &_pi;
    303302      int &_in_arc;
     
    316315      bool findEnteringArc() {
    317316        Cost c, min = 0;
    318         for (int e = 0; e != _search_arc_num; ++e) {
     317        for (int e = 0; e < _search_arc_num; ++e) {
    319318          c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
    320319          if (c < min) {
     
    338337      const IntVector  &_target;
    339338      const CostVector &_cost;
    340       const BoolVector &_state;
     339      const CharVector &_state;
    341340      const CostVector &_pi;
    342341      int &_in_arc;
     
    357356      {
    358357        // The main parameters of the pivot rule
    359         const double BLOCK_SIZE_FACTOR = 1.0;
     358        const double BLOCK_SIZE_FACTOR = 0.5;
    360359        const int MIN_BLOCK_SIZE = 10;
    361360
     
    370369        int cnt = _block_size;
    371370        int e;
    372         for (e = _next_arc; e != _search_arc_num; ++e) {
     371        for (e = _next_arc; e < _search_arc_num; ++e) {
    373372          c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
    374373          if (c < min) {
     
    381380          }
    382381        }
    383         for (e = 0; e != _next_arc; ++e) {
     382        for (e = 0; e < _next_arc; ++e) {
    384383          c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
    385384          if (c < min) {
     
    411410      const IntVector  &_target;
    412411      const CostVector &_cost;
    413       const BoolVector &_state;
     412      const CharVector &_state;
    414413      const CostVector &_pi;
    415414      int &_in_arc;
     
    472471        min = 0;
    473472        _curr_length = 0;
    474         for (e = _next_arc; e != _search_arc_num; ++e) {
     473        for (e = _next_arc; e < _search_arc_num; ++e) {
    475474          c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
    476475          if (c < 0) {
     
    483482          }
    484483        }
    485         for (e = 0; e != _next_arc; ++e) {
     484        for (e = 0; e < _next_arc; ++e) {
    486485          c = _state[e] * (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
    487486          if (c < 0) {
     
    514513      const IntVector  &_target;
    515514      const CostVector &_cost;
    516       const BoolVector &_state;
     515      const CharVector &_state;
    517516      const CostVector &_pi;
    518517      int &_in_arc;
     
    567566        // Check the current candidate list
    568567        int e;
    569         for (int i = 0; i != _curr_length; ++i) {
     568        for (int i = 0; i < _curr_length; ++i) {
    570569          e = _candidates[i];
    571570          _cand_cost[e] = _state[e] *
     
    580579        int limit = _head_length;
    581580
    582         for (e = _next_arc; e != _search_arc_num; ++e) {
     581        for (e = _next_arc; e < _search_arc_num; ++e) {
    583582          _cand_cost[e] = _state[e] *
    584583            (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
     
    592591          }
    593592        }
    594         for (e = 0; e != _next_arc; ++e) {
     593        for (e = 0; e < _next_arc; ++e) {
    595594          _cand_cost[e] = _state[e] *
    596595            (_cost[e] + _pi[_source[e]] - _pi[_target[e]]);
     
    13621361
    13631362      // Update _rev_thread using the new _thread values
    1364       for (int i = 0; i != int(_dirty_revs.size()); ++i) {
     1363      for (int i = 0; i < int(_dirty_revs.size()); ++i) {
    13651364        u = _dirty_revs[i];
    13661365        _rev_thread[_thread[u]] = u;
     
    14321431        _pi[u] += sigma;
    14331432      }
    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;
    15281433    }
    15291434
     
    15501455      PivotRuleImpl pivot(*this);
    15511456
    1552       // Perform heuristic initial pivots
    1553       if (!initialPivots()) return UNBOUNDED;
    1554 
    15551457      // Execute the Network Simplex algorithm
    15561458      while (pivot.findEnteringArc()) {
Note: See TracChangeset for help on using the changeset viewer.