Changeset 941:186aa53d2802 in lemon0.x for src/lemon/min_cost_flow.h
 Timestamp:
 10/08/04 15:07:51 (16 years ago)
 Branch:
 default
 Phase:
 public
 Convert:
 svn:c9d7d8f590d60310b91f818b3a526b0e/lemon/trunk@1283
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

src/lemon/min_cost_flow.h
r921 r941 71 71 typedef typename Graph::template EdgeMap<int> EdgeIntMap; 72 72 73 74 73 typedef ResGraphWrapper<const Graph,int,CapacityMap,EdgeIntMap> ResGW; 75 74 typedef typename ResGW::Edge ResGraphEdge; 76 75 76 protected: 77 78 const Graph& g; 79 const LengthMap& length; 80 const CapacityMap& capacity; 81 82 EdgeIntMap flow; 83 typedef typename Graph::template NodeMap<Length> PotentialMap; 84 PotentialMap potential; 85 86 Node s; 87 Node t; 88 89 Length total_length; 90 77 91 class ModLengthMap { 78 92 typedef typename Graph::template NodeMap<Length> NodeMap; 79 const ResGW& G;80 const LengthMap & ol;93 const ResGW& g; 94 const LengthMap &length; 81 95 const NodeMap &pot; 82 96 public : 83 97 typedef typename LengthMap::KeyType KeyType; 84 98 typedef typename LengthMap::ValueType ValueType; 99 100 ModLengthMap(const ResGW& _g, 101 const LengthMap &_length, const NodeMap &_pot) : 102 g(_g), /*rev(_rev),*/ length(_length), pot(_pot) { } 85 103 86 104 ValueType operator[](typename ResGW::Edge e) const { 87 if ( G.forward(e))88 return ol[e](pot[G.head(e)]pot[G.tail(e)]);105 if (g.forward(e)) 106 return length[e](pot[g.head(e)]pot[g.tail(e)]); 89 107 else 90 return  ol[e](pot[G.head(e)]pot[G.tail(e)]);108 return length[e](pot[g.head(e)]pot[g.tail(e)]); 91 109 } 92 110 93 ModLengthMap(const ResGW& _G, 94 const LengthMap &o, const NodeMap &p) : 95 G(_G), /*rev(_rev),*/ ol(o), pot(p){}; 96 };//ModLengthMap 97 98 99 protected: 100 101 //Input 102 const Graph& G; 103 const LengthMap& length; 104 const CapacityMap& capacity; 105 106 107 //auxiliary variables 108 109 //To store the flow 110 EdgeIntMap flow; 111 //To store the potential (dual variables) 112 typedef typename Graph::template NodeMap<Length> PotentialMap; 113 PotentialMap potential; 114 115 116 Length total_length; 117 111 }; //ModLengthMap 112 113 ResGW res_graph; 114 ModLengthMap mod_length; 115 Dijkstra<ResGW, ModLengthMap> dijkstra; 118 116 119 117 public : 120 118 121 /// The constructor of the class. 122 123 ///\param _G The directed graph the algorithm runs on. 124 ///\param _length The length (weight or cost) of the edges. 125 ///\param _cap The capacity of the edges. 126 MinCostFlow(Graph& _G, LengthMap& _length, CapacityMap& _cap) : G(_G), 127 length(_length), capacity(_cap), flow(_G), potential(_G){ } 128 129 130 ///Runs the algorithm. 131 132 ///Runs the algorithm. 133 ///Returns k if there is a flow of value at least k edgedisjoint 134 ///from s to t. 135 ///Otherwise it returns the maximum value of a flow from s to t. 136 /// 137 ///\param s The source node. 138 ///\param t The target node. 139 ///\param k The value of the flow we are looking for. 140 /// 141 ///\todo May be it does make sense to be able to start with a nonzero 142 /// feasible primaldual solution pair as well. 143 int run(Node s, Node t, int k) { 144 145 //Resetting variables from previous runs 146 total_length = 0; 147 148 for (typename Graph::EdgeIt e(G); e!=INVALID; ++e) flow.set(e, 0); 149 150 //Initialize the potential to zero 151 for (typename Graph::NodeIt n(G); n!=INVALID; ++n) potential.set(n, 0); 152 153 154 //We need a residual graph 155 ResGW res_graph(G, capacity, flow); 156 157 158 ModLengthMap mod_length(res_graph, length, potential); 159 160 Dijkstra<ResGW, ModLengthMap> dijkstra(res_graph, mod_length); 161 162 int i; 163 for (i=0; i<k; ++i){ 164 dijkstra.run(s); 165 if (!dijkstra.reached(t)){ 166 //There are no flow of value k from s to t 167 break; 168 }; 119 /*! \brief The constructor of the class. 120 121 \param _g The directed graph the algorithm runs on. 122 \param _length The length (weight or cost) of the edges. 123 \param _cap The capacity of the edges. 124 \param _s Source node. 125 \param _t Target node. 126 */ 127 MinCostFlow(Graph& _g, LengthMap& _length, CapacityMap& _cap, 128 Node _s, Node _t) : 129 g(_g), length(_length), capacity(_cap), flow(_g), potential(_g), 130 s(_s), t(_t), 131 res_graph(g, capacity, flow), 132 mod_length(res_graph, length, potential), 133 dijkstra(res_graph, mod_length) { 134 reset(); 135 } 136 137 /*! Tries to augment the flow between s and t by 1. 138 The return value shows if the augmentation is successful. 139 */ 140 bool augment() { 141 dijkstra.run(s); 142 if (!dijkstra.reached(t)) { 143 144 //Unsuccessful augmentation. 145 return false; 146 } else { 147 148 //We have to change the potential 149 for(typename ResGW::NodeIt n(res_graph); n!=INVALID; ++n) 150 potential[n] += dijkstra.distMap()[n]; 169 151 170 //We have to change the potential171 for(typename ResGW::NodeIt n(res_graph); n!=INVALID; ++n)172 potential[n] += dijkstra.distMap()[n];173 174 175 152 //Augmenting on the sortest path 176 153 Node n=t; … … 187 164 } 188 165 189 166 return true; 190 167 } 191 192 168 } 169 170 /*! \brief Runs the algorithm. 171 172 Runs the algorithm. 173 Returns k if there is a flow of value at least k from s to t. 174 Otherwise it returns the maximum value of a flow from s to t. 175 176 \param k The value of the flow we are looking for. 177 178 \todo May be it does make sense to be able to start with a nonzero 179 feasible primaldual solution pair as well. 180 181 \todo If the actual flow value is bigger than k, then everything is 182 cleared and the algorithm starts from zero flow. Is it a good approach? 183 */ 184 int run(int k) { 185 if (flowValue()>k) reset(); 186 while (flowValue()<k && augment()) { } 187 return flowValue(); 188 } 189 190 /*! \brief The class is reset to zero flow and potential. 191 The class is reset to zero flow and potential. 192 */ 193 void reset() { 194 total_length=0; 195 for (typename Graph::EdgeIt e(g); e!=INVALID; ++e) flow.set(e, 0); 196 for (typename Graph::NodeIt n(g); n!=INVALID; ++n) potential.set(n, 0); 197 } 198 199 /*! Returns the value of the actual flow. 200 */ 201 int flowValue() const { 202 int i=0; 203 for (typename Graph::OutEdgeIt e(g, s); e!=INVALID; ++e) i+=flow[e]; 204 for (typename Graph::InEdgeIt e(g, s); e!=INVALID; ++e) i=flow[e]; 193 205 return i; 194 206 } 195 207 196 197 198 /// Gives back the total weight of the found flow. 199 200 ///This function gives back the total weight of the found flow. 201 ///Assumes that \c run() has been run and nothing changed since then. 208 /// Total weight of the found flow. 209 210 /// This function gives back the total weight of the found flow. 202 211 Length totalLength(){ 203 212 return total_length; … … 207 216 208 217 ///Returns a const reference to the EdgeMap \c flow. 209 ///\pre \ref run() must210 ///be called before using this function.211 218 const EdgeIntMap &getFlow() const { return flow;} 212 219 213 / //Returns a const reference to the NodeMap \c potential (the dual solution).214 215 ///Returns a const reference to the NodeMap \c potential (the dual solution).216 /// \pre \ref run() must be called before using this function.220 /*! \brief Returns a const reference to the NodeMap \c potential (the dual solution). 221 222 Returns a const reference to the NodeMap \c potential (the dual solution). 223 */ 217 224 const PotentialMap &getPotential() const { return potential;} 218 225 219 / // Checking the complementary slackness optimality criteria220 221 ///This function checks, whether the given solution is optimal222 ///If executed after the call of \c run() then it should return with true.223 ///This function only checks optimality, doesn't bother with feasibility.224 ///It is meant for testing purposes.225 ///226 /*! \brief Checking the complementary slackness optimality criteria. 227 228 This function checks, whether the given flow and potential 229 satisfiy the complementary slackness cnditions (i.e. these are optimal). 230 This function only checks optimality, doesn't bother with feasibility. 231 For testing purpose. 232 */ 226 233 bool checkComplementarySlackness(){ 227 234 Length mod_pot; 228 235 Length fl_e; 229 for(typename Graph::EdgeIt e( G); e!=INVALID; ++e) {236 for(typename Graph::EdgeIt e(g); e!=INVALID; ++e) { 230 237 //C^{\Pi}_{i,j} 231 mod_pot = length[e]potential[ G.head(e)]+potential[G.tail(e)];238 mod_pot = length[e]potential[g.head(e)]+potential[g.tail(e)]; 232 239 fl_e = flow[e]; 233 240 if (0<fl_e && fl_e<capacity[e]) { … … 247 254 } 248 255 249 250 256 }; //class MinCostFlow 251 257
Note: See TracChangeset
for help on using the changeset viewer.