Changeset 2635:23570e368afa in lemon0.x
 Timestamp:
 06/01/09 17:35:27 (10 years ago)
 Branch:
 default
 Phase:
 public
 Convert:
 svn:c9d7d8f590d60310b91f818b3a526b0e/lemon/trunk@3520
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

lemon/network_simplex.h
r2634 r2635 242 242 _edge(ns._edge), _source(ns._source), _target(ns._target), 243 243 _cost(ns._cost), _state(ns._state), _pi(ns._pi), 244 _in_edge(ns._in_edge), _edge_num(ns._edge_num ), _next_edge(0)244 _in_edge(ns._in_edge), _edge_num(ns._edge_num + ns._node_num), _next_edge(0) 245 245 { 246 246 // The main parameters of the pivot rule … … 592 592 int _edge_num; 593 593 594 // The depth vector of the spanning tree structure595 IntVector _depth;596 594 // The parent vector of the spanning tree structure 597 595 IntVector _parent; … … 600 598 // The thread vector of the spanning tree structure 601 599 IntVector _thread; 600 601 IntVector _rev_thread; 602 IntVector _succ_num; 603 IntVector _last_succ; 604 605 IntVector _dirty_revs; 606 602 607 // The forward vector of the spanning tree structure 603 608 BoolVector _forward; … … 883 888 _pi.resize(all_node_num, 0); 884 889 885 _depth.resize(all_node_num);886 890 _parent.resize(all_node_num); 887 891 _pred.resize(all_node_num); 892 _forward.resize(all_node_num); 888 893 _thread.resize(all_node_num); 889 _forward.resize(all_node_num); 894 _rev_thread.resize(all_node_num); 895 _succ_num.resize(all_node_num); 896 _last_succ.resize(all_node_num); 890 897 _state.resize(all_edge_num, STATE_LOWER); 891 898 … … 916 923 // Set data for the artificial root node 917 924 _root = _node_num; 918 _depth[_root] = 0;919 925 _parent[_root] = 1; 920 926 _pred[_root] = 1; 921 927 _thread[_root] = 0; 928 _rev_thread[0] = _root; 929 _succ_num[_root] = all_node_num; 930 _last_succ[_root] = _root  1; 922 931 _supply[_root] = 0; 923 932 _pi[_root] = 0; … … 956 965 Capacity max_cap = std::numeric_limits<Capacity>::max(); 957 966 for (int u = 0, e = _edge_num; u != _node_num; ++u, ++e) { 958 _thread[u] = u + 1;959 _depth[u] = 1;960 967 _parent[u] = _root; 961 968 _pred[u] = e; 969 _thread[u] = u + 1; 970 _rev_thread[u + 1] = u; 971 _succ_num[u] = 1; 972 _last_succ[u] = u; 973 _cap[e] = max_cap; 974 _state[e] = STATE_TREE; 962 975 if (_supply[u] >= 0) { 976 _forward[u] = true; 977 _pi[u] = 0; 978 _source[e] = u; 979 _target[e] = _root; 963 980 _flow[e] = _supply[u]; 964 _forward[u] = true; 965 _pi[u] = max_cost; 966 } else { 967 _flow[e] = _supply[u]; 981 _cost[e] = 0; 982 } 983 else { 968 984 _forward[u] = false; 969 985 _pi[u] = max_cost; 970 } 971 _cost[e] = max_cost; 972 _cap[e] = max_cap; 973 _state[e] = STATE_TREE; 986 _source[e] = _root; 987 _target[e] = u; 988 _flow[e] = _supply[u]; 989 _cost[e] = max_cost; 990 } 974 991 } 975 992 … … 981 998 int u = _source[_in_edge]; 982 999 int v = _target[_in_edge]; 983 while (_depth[u] > _depth[v]) u = _parent[u];984 while (_depth[v] > _depth[u]) v = _parent[v];985 1000 while (u != v) { 986 u = _parent[u]; 987 v = _parent[v]; 1001 if (_succ_num[u] < _succ_num[v]) { 1002 u = _parent[u]; 1003 } else { 1004 v = _parent[v]; 1005 } 988 1006 } 989 1007 join = u; … … 1060 1078 } 1061 1079 } 1062 1063 // Update _thread and _parent vectors 1064 void updateThreadParent() { 1065 int u; 1080 1081 // Update the tree structure 1082 void updateTreeStructure() { 1083 int u, w; 1084 int old_rev_thread = _rev_thread[u_out]; 1085 int old_succ_num = _succ_num[u_out]; 1086 int old_last_succ = _last_succ[u_out]; 1066 1087 v_out = _parent[u_out]; 1067 1088 1068 // Handle the case when join and v_out coincide 1069 bool par_first = false; 1070 if (join == v_out) { 1071 for (u = join; u != u_in && u != v_in; u = _thread[u]) ; 1072 if (u == v_in) { 1073 par_first = true; 1074 while (_thread[u] != u_out) u = _thread[u]; 1075 first = u; 1076 } 1077 } 1078 1079 // Find the last successor of u_in (u) and the node after it (right) 1080 // according to the thread index 1081 for (u = u_in; _depth[_thread[u]] > _depth[u_in]; u = _thread[u]) ; 1082 right = _thread[u]; 1083 if (_thread[v_in] == u_out) { 1084 for (last = u; _depth[last] > _depth[u_out]; last = _thread[last]) ; 1085 if (last == u_out) last = _thread[last]; 1086 } 1087 else last = _thread[v_in]; 1088 1089 // Update stem nodes 1089 u = _last_succ[u_in]; // the last successor of u_in 1090 right = _thread[u]; // the node after it 1091 1092 // Handle the case when old_rev_thread equals to v_in 1093 // (it also means that join and v_out coincide) 1094 if (old_rev_thread == v_in) { 1095 last = _thread[_last_succ[u_out]]; 1096 } else { 1097 last = _thread[v_in]; 1098 } 1099 1100 // Update _thread and _parent along the stem nodes (i.e. the nodes 1101 // between u_in and u_out, whose parent have to be changed) 1090 1102 _thread[v_in] = stem = u_in; 1103 _dirty_revs.clear(); 1104 _dirty_revs.push_back(v_in); 1091 1105 par_stem = v_in; 1092 1106 while (stem != u_out) { 1093 _thread[u] = new_stem = _parent[stem]; 1094 1095 // Find the node just before the stem node (u) according to 1096 // the original thread index 1097 for (u = new_stem; _thread[u] != stem; u = _thread[u]) ; 1098 _thread[u] = right; 1099 1100 // Change the parent node of stem and shift stem and par_stem nodes 1101 _parent[stem] = par_stem; 1107 // Insert the next stem node into the thread list 1108 new_stem = _parent[stem]; 1109 _thread[u] = new_stem; 1110 _dirty_revs.push_back(u); 1111 1112 // Remove the subtree of stem from the thread list 1113 w = _rev_thread[stem]; 1114 _thread[w] = right; 1115 _rev_thread[right] = w; 1116 1117 // Change the parent node and shift stem nodes 1118 _parent[stem] = par_stem; 1102 1119 par_stem = stem; 1103 1120 stem = new_stem; 1104 1121 1105 // Find the last successor of stem (u) and the node after it (right)1106 // according to the thread index1107 for (u = stem; _depth[_thread[u]] > _depth[stem]; u = _thread[u]);1122 // Update u and right nodes 1123 u = _last_succ[stem] == _last_succ[par_stem] ? 1124 _rev_thread[par_stem] : _last_succ[stem]; 1108 1125 right = _thread[u]; 1109 1126 } 1110 1127 _parent[u_out] = par_stem; 1128 _last_succ[u_out] = u; 1111 1129 _thread[u] = last; 1112 1113 if (join == v_out && par_first) { 1114 if (first != v_in) _thread[first] = right; 1130 _rev_thread[last] = u; 1131 1132 // Remove the subtree of u_out from the thread list except for 1133 // the case when old_rev_thread equals to v_in 1134 // (it also means that join and v_out coincide) 1135 if (old_rev_thread != v_in) { 1136 _thread[old_rev_thread] = right; 1137 _rev_thread[right] = old_rev_thread; 1138 } 1139 1140 // Update _rev_thread using the new _thread values 1141 for (int i = 0; i < int(_dirty_revs.size()); ++i) { 1142 u = _dirty_revs[i]; 1143 _rev_thread[_thread[u]] = u; 1144 } 1145 1146 // Update _last_succ for the stem nodes from u_out to u_in 1147 for (u = u_out; u != u_in; u = _parent[u]) { 1148 _last_succ[_parent[u]] = _last_succ[u]; 1149 } 1150 1151 // Set limits for updating _last_succ form v_in and v_out 1152 // towards the root 1153 int up_limit_in = 1; 1154 int up_limit_out = 1; 1155 if (_last_succ[join] == v_in) { 1156 up_limit_out = join; 1115 1157 } else { 1116 for (u = v_out; _thread[u] != u_out; u = _thread[u]) ; 1117 _thread[u] = right; 1118 } 1119 } 1120 1121 // Update _pred and _forward vectors 1122 void updatePredEdge() { 1123 int u = u_out, v; 1158 up_limit_in = join; 1159 } 1160 1161 // Update _last_succ from v_in towards the root 1162 for (u = v_in; u != up_limit_in && _last_succ[u] == v_in; 1163 u = _parent[u]) { 1164 _last_succ[u] = _last_succ[u_out]; 1165 } 1166 // Update _last_succ from v_out towards the root 1167 if (join != old_rev_thread && v_in != old_rev_thread) { 1168 for (u = v_out; u != up_limit_out && _last_succ[u] == old_last_succ; 1169 u = _parent[u]) { 1170 _last_succ[u] = old_rev_thread; 1171 } 1172 } else { 1173 for (u = v_out; u != up_limit_out && _last_succ[u] == old_last_succ; 1174 u = _parent[u]) { 1175 _last_succ[u] = _last_succ[u_out]; 1176 } 1177 } 1178 1179 // Update _pred and _forward for the stem nodes from u_out to u_in 1180 u = u_out; 1124 1181 while (u != u_in) { 1125 v= _parent[u];1126 _pred[u] = _pred[ v];1127 _forward[u] = !_forward[ v];1128 u = v;1182 w = _parent[u]; 1183 _pred[u] = _pred[w]; 1184 _forward[u] = !_forward[w]; 1185 u = w; 1129 1186 } 1130 1187 _pred[u_in] = _in_edge; 1131 1188 _forward[u_in] = (u_in == _source[_in_edge]); 1132 } 1133 1134 // Update _depth and _potential vectors 1135 void updateDepthPotential() { 1136 _depth[u_in] = _depth[v_in] + 1; 1189 1190 // Update _succ_num from v_in to join 1191 for (u = v_in; u != join; u = _parent[u]) { 1192 _succ_num[u] += old_succ_num; 1193 } 1194 // Update _succ_num from v_out to join 1195 for (u = v_out; u != join; u = _parent[u]) { 1196 _succ_num[u] = old_succ_num; 1197 } 1198 // Update _succ_num for the stem nodes from u_out to u_in 1199 int tmp = 0; 1200 u = u_out; 1201 while (u != u_in) { 1202 w = _parent[u]; 1203 tmp = _succ_num[u]  _succ_num[w] + tmp; 1204 _succ_num[u] = tmp; 1205 u = w; 1206 } 1207 _succ_num[u_in] = old_succ_num; 1208 } 1209 1210 // Update potentials 1211 void updatePotential() { 1137 1212 Cost sigma = _forward[u_in] ? 1138 1213 _pi[v_in]  _pi[u_in]  _cost[_pred[u_in]] : 1139 1214 _pi[v_in]  _pi[u_in] + _cost[_pred[u_in]]; 1140 _pi[u_in] += sigma; 1141 for(int u = _thread[u_in]; _parent[u] != 1; u = _thread[u]) { 1142 _depth[u] = _depth[_parent[u]] + 1; 1143 if (_depth[u] <= _depth[u_in]) break; 1215 // Update in the lower subtree (which has been moved) 1216 int end = _thread[_last_succ[u_in]]; 1217 for (int u = u_in; u != end; u = _thread[u]) { 1144 1218 _pi[u] += sigma; 1145 1219 } … … 1167 1241 bool start() { 1168 1242 PivotRuleImplementation pivot(*this); 1169 1243 1170 1244 // Execute the network simplex algorithm 1171 1245 while (pivot.findEnteringEdge()) { … … 1174 1248 changeFlow(change); 1175 1249 if (change) { 1176 updateThreadParent(); 1177 updatePredEdge(); 1178 updateDepthPotential(); 1179 } 1180 } 1181 1250 updateTreeStructure(); 1251 updatePotential(); 1252 } 1253 } 1254 1182 1255 // Check if the flow amount equals zero on all the artificial edges 1183 1256 for (int e = _edge_num; e != _edge_num + _node_num; ++e) { … … 1200 1273 (*_potential_result)[_node[i]] = _pi[i]; 1201 1274 } 1202 1275 1203 1276 return true; 1204 1277 }
Note: See TracChangeset
for help on using the changeset viewer.