gravatar
kpeter (Peter Kovacs)
kpeter@inf.elte.hu
Minor improvements in CostScaling (#417)
0 1 0
default
1 file changed with 49 insertions and 37 deletions:
↑ Collapse diff ↑
Ignore white space 12 line context
... ...
@@ -572,24 +572,31 @@
572 572
        _scost[_reverse[j]] = 0;
573 573
      }
574 574
      _have_lower = false;
575 575
      return *this;
576 576
    }
577 577

	
578
    /// \brief Reset all the parameters that have been given before.
578
    /// \brief Reset the internal data structures and all the parameters
579
    /// that have been given before.
579 580
    ///
580
    /// This function resets all the paramaters that have been given
581
    /// before using functions \ref lowerMap(), \ref upperMap(),
582
    /// \ref costMap(), \ref supplyMap(), \ref stSupply().
581
    /// This function resets the internal data structures and all the
582
    /// paramaters that have been given before using functions \ref lowerMap(),
583
    /// \ref upperMap(), \ref costMap(), \ref supplyMap(), \ref stSupply().
583 584
    ///
584
    /// It is useful for multiple run() calls. If this function is not
585
    /// used, all the parameters given before are kept for the next
586
    /// \ref run() call.
587
    /// However, the underlying digraph must not be modified after this
588
    /// class have been constructed, since it copies and extends the graph.
585
    /// It is useful for multiple \ref run() calls. By default, all the given
586
    /// parameters are kept for the next \ref run() call, unless
587
    /// \ref resetParams() or \ref reset() is used.
588
    /// If the underlying digraph was also modified after the construction
589
    /// of the class or the last \ref reset() call, then the \ref reset()
590
    /// function must be used, otherwise \ref resetParams() is sufficient.
591
    ///
592
    /// See \ref resetParams() for examples.
593
    ///
589 594
    /// \return <tt>(*this)</tt>
595
    ///
596
    /// \see resetParams(), run()
590 597
    CostScaling& reset() {
591 598
      // Resize vectors
592 599
      _node_num = countNodes(_graph);
593 600
      _arc_num = countArcs(_graph);
594 601
      _res_node_num = _node_num + 1;
595 602
      _res_arc_num = 2 * (_arc_num + _node_num);
... ...
@@ -887,37 +894,35 @@
887 894
          _res_cap[ra] = 0;
888 895
          _cost[a] = 0;
889 896
          _cost[ra] = 0;
890 897
        }
891 898
      }
892 899

	
893
      return OPTIMAL;
894
    }
895

	
896
    // Execute the algorithm and transform the results
897
    void start(Method method) {
898
      // Maximum path length for partial augment
899
      const int MAX_PATH_LENGTH = 4;
900

	
901 900
      // Initialize data structures for buckets
902 901
      _max_rank = _alpha * _res_node_num;
903 902
      _buckets.resize(_max_rank);
904 903
      _bucket_next.resize(_res_node_num + 1);
905 904
      _bucket_prev.resize(_res_node_num + 1);
906 905
      _rank.resize(_res_node_num + 1);
907 906

	
908
      // Execute the algorithm
907
      return OPTIMAL;
908
    }
909

	
910
    // Execute the algorithm and transform the results
911
    void start(Method method) {
912
      const int MAX_PARTIAL_PATH_LENGTH = 4;
913

	
909 914
      switch (method) {
910 915
        case PUSH:
911 916
          startPush();
912 917
          break;
913 918
        case AUGMENT:
914 919
          startAugment(_res_node_num - 1);
915 920
          break;
916 921
        case PARTIAL_AUGMENT:
917
          startAugment(MAX_PATH_LENGTH);
922
          startAugment(MAX_PARTIAL_PATH_LENGTH);
918 923
          break;
919 924
      }
920 925

	
921 926
      // Compute node potentials for the original costs
922 927
      _arc_vec.clear();
923 928
      _cost_vec.clear();
... ...
@@ -948,19 +953,21 @@
948 953
    void initPhase() {
949 954
      // Saturate arcs not satisfying the optimality condition
950 955
      for (int u = 0; u != _res_node_num; ++u) {
951 956
        int last_out = _first_out[u+1];
952 957
        LargeCost pi_u = _pi[u];
953 958
        for (int a = _first_out[u]; a != last_out; ++a) {
954
          int v = _target[a];
955
          if (_res_cap[a] > 0 && _cost[a] + pi_u - _pi[v] < 0) {
956
            Value delta = _res_cap[a];
957
            _excess[u] -= delta;
958
            _excess[v] += delta;
959
            _res_cap[a] = 0;
960
            _res_cap[_reverse[a]] += delta;
959
          Value delta = _res_cap[a];
960
          if (delta > 0) {
961
            int v = _target[a];
962
            if (_cost[a] + pi_u - _pi[v] < 0) {
963
              _excess[u] -= delta;
964
              _excess[v] += delta;
965
              _res_cap[a] = 0;
966
              _res_cap[_reverse[a]] += delta;
967
            }
961 968
          }
962 969
        }
963 970
      }
964 971

	
965 972
      // Find active nodes (i.e. nodes with positive excess)
966 973
      for (int u = 0; u != _res_node_num; ++u) {
... ...
@@ -998,31 +1005,33 @@
998 1005
      }
999 1006
      return done;
1000 1007
    }
1001 1008

	
1002 1009
    // Global potential update heuristic
1003 1010
    void globalUpdate() {
1004
      int bucket_end = _root + 1;
1011
      const int bucket_end = _root + 1;
1005 1012

	
1006 1013
      // Initialize buckets
1007 1014
      for (int r = 0; r != _max_rank; ++r) {
1008 1015
        _buckets[r] = bucket_end;
1009 1016
      }
1010 1017
      Value total_excess = 0;
1018
      int b0 = bucket_end;
1011 1019
      for (int i = 0; i != _res_node_num; ++i) {
1012 1020
        if (_excess[i] < 0) {
1013 1021
          _rank[i] = 0;
1014
          _bucket_next[i] = _buckets[0];
1015
          _bucket_prev[_buckets[0]] = i;
1016
          _buckets[0] = i;
1022
          _bucket_next[i] = b0;
1023
          _bucket_prev[b0] = i;
1024
          b0 = i;
1017 1025
        } else {
1018 1026
          total_excess += _excess[i];
1019 1027
          _rank[i] = _max_rank;
1020 1028
        }
1021 1029
      }
1022 1030
      if (total_excess == 0) return;
1031
      _buckets[0] = b0;
1023 1032

	
1024 1033
      // Search the buckets
1025 1034
      int r = 0;
1026 1035
      for ( ; r != _max_rank; ++r) {
1027 1036
        while (_buckets[r] != bucket_end) {
1028 1037
          // Remove the first node from the current bucket
... ...
@@ -1038,33 +1047,36 @@
1038 1047
              int v = _source[ra];
1039 1048
              int old_rank_v = _rank[v];
1040 1049
              if (r < old_rank_v) {
1041 1050
                // Compute the new rank of v
1042 1051
                LargeCost nrc = (_cost[ra] + _pi[v] - pi_u) / _epsilon;
1043 1052
                int new_rank_v = old_rank_v;
1044
                if (nrc < LargeCost(_max_rank))
1045
                  new_rank_v = r + 1 + int(nrc);
1053
                if (nrc < LargeCost(_max_rank)) {
1054
                  new_rank_v = r + 1 + static_cast<int>(nrc);
1055
                }
1046 1056

	
1047 1057
                // Change the rank of v
1048 1058
                if (new_rank_v < old_rank_v) {
1049 1059
                  _rank[v] = new_rank_v;
1050 1060
                  _next_out[v] = _first_out[v];
1051 1061

	
1052 1062
                  // Remove v from its old bucket
1053 1063
                  if (old_rank_v < _max_rank) {
1054 1064
                    if (_buckets[old_rank_v] == v) {
1055 1065
                      _buckets[old_rank_v] = _bucket_next[v];
1056 1066
                    } else {
1057
                      _bucket_next[_bucket_prev[v]] = _bucket_next[v];
1058
                      _bucket_prev[_bucket_next[v]] = _bucket_prev[v];
1067
                      int pv = _bucket_prev[v], nv = _bucket_next[v];
1068
                      _bucket_next[pv] = nv;
1069
                      _bucket_prev[nv] = pv;
1059 1070
                    }
1060 1071
                  }
1061 1072

	
1062
                  // Insert v to its new bucket
1063
                  _bucket_next[v] = _buckets[new_rank_v];
1064
                  _bucket_prev[_buckets[new_rank_v]] = v;
1073
                  // Insert v into its new bucket
1074
                  int nv = _buckets[new_rank_v];
1075
                  _bucket_next[v] = nv;
1076
                  _bucket_prev[nv] = v;
1065 1077
                  _buckets[new_rank_v] = v;
1066 1078
                }
1067 1079
              }
1068 1080
            }
1069 1081
          }
1070 1082

	
0 comments (0 inline)