COIN-OR::LEMON - Graph Library

Changeset 659:c5984e925384 in lemon-0.x


Ignore:
Timestamp:
05/25/04 14:31:18 (15 years ago)
Author:
athos
Branch:
default
Phase:
public
Convert:
svn:c9d7d8f5-90d6-0310-b91f-818b3a526b0e/lemon/trunk@862
Message:

Almost ready.

Location:
src/work/athos
Files:
1 added
1 edited

Legend:

Unmodified
Added
Removed
  • src/work/athos/mincostflow.h

    r657 r659  
    6060      //typedef typename ResGraphType::template NodeMap<Length> NodeMap;
    6161      typedef typename Graph::template NodeMap<Length> NodeMap;
    62       const ResGraphType& G;
     62      const ResGraphType& res_graph;
    6363      //      const EdgeIntMap& rev;
    6464      const LengthMap &ol;
     
    6969       
    7070      ValueType operator[](typename ResGraphType::Edge e) const {     
    71         if (G.forward(e))
    72           return  ol[e]-(pot[G.head(e)]-pot[G.tail(e)]);   
     71        if (res_graph.forward(e))
     72          return  ol[e]-(pot[res_graph.head(e)]-pot[res_graph.tail(e)]);   
    7373        else
    74           return -ol[e]-(pot[G.head(e)]-pot[G.tail(e)]);   
     74          return -ol[e]-(pot[res_graph.head(e)]-pot[res_graph.tail(e)]);   
    7575      }     
    7676       
    77       ModLengthMap(const ResGraphType& _G,
     77      ModLengthMap(const ResGraphType& _res_graph,
    7878                   const LengthMap &o,  const NodeMap &p) :
    79         G(_G), /*rev(_rev),*/ ol(o), pot(p){};
     79        res_graph(_res_graph), /*rev(_rev),*/ ol(o), pot(p){};
    8080    };//ModLengthMap
    8181
     
    8484   
    8585    //Input
    86     const Graph& G;
     86    const Graph& graph;
    8787    const LengthMap& length;
    8888    const SupplyDemandMap& supply_demand;//supply or demand of nodes
     
    105105
    106106
    107     MinCostFlow(Graph& _G, LengthMap& _length, SupplyDemandMap& _supply_demand) : G(_G),
    108       length(_length), supply_demand(_supply_demand), flow(_G), potential(_G){ }
     107    MinCostFlow(Graph& _graph, LengthMap& _length, SupplyDemandMap& _supply_demand) : graph(_graph),
     108      length(_length), supply_demand(_supply_demand), flow(_graph), potential(_graph){ }
    109109
    110110   
     
    115115    ///\todo May be it does make sense to be able to start with a nonzero
    116116    /// feasible primal-dual solution pair as well.
    117     int run() {
     117    void run() {
    118118
    119119      //Resetting variables from previous runs
     
    125125
    126126      //A heap for the excess nodes
    127       HeapMap excess_nodes_map(G,-1);
     127      HeapMap excess_nodes_map(graph,-1);
    128128      HeapType excess_nodes(excess_nodes_map);
    129129
    130130      //A heap for the deficit nodes
    131       HeapMap deficit_nodes_map(G,-1);
     131      HeapMap deficit_nodes_map(graph,-1);
    132132      HeapType deficit_nodes(deficit_nodes_map);
    133133
    134134      //A container to store nonabundant arcs
    135135      list<Edge> nonabundant_arcs;
    136      
    137       FOR_EACH_LOC(typename Graph::EdgeIt, e, G){
     136
     137       
     138      FOR_EACH_LOC(typename Graph::EdgeIt, e, graph){
    138139        flow.set(e,0);
    139140        nonabundant_arcs.push_back(e);
     
    146147
    147148      //A union-find structure to store the abundant components
    148       UFE::MapType abund_comp_map(G);
     149      UFE::MapType abund_comp_map(graph);
    149150      UFE abundant_components(abund_comp_map);
    150151
    151152
    152153
    153       FOR_EACH_LOC(typename Graph::NodeIt, n, G){
     154      FOR_EACH_LOC(typename Graph::NodeIt, n, graph){
    154155        excess_deficit.set(n,supply_demand[n]);
    155156        //A supply node
     
    175176      SupplyDemand max_excess = delta;
    176177     
     178      //ConstMap<Edge,SupplyDemand> ConstEdgeMap;
     179
    177180      //We need a residual graph which is uncapacitated
    178       ResGraphType res_graph(G, flow);
    179 
    180 
     181      ResGraphType res_graph(graph, flow);
     182     
     183      //An EdgeMap to tell which arcs are abundant
     184      template typename Graph::EdgeMap<bool> abundant_arcs(graph);
     185
     186      //Let's construct the sugraph consisting only of the abundant edges
     187      typedef ConstMap< typename Graph::Node, bool > ConstNodeMap;
     188      ConstNodeMap const_true_map(true);
     189      typedef SubGraphWrapper< Graph, ConstNodeMap,
     190         template typename Graph::EdgeMap<bool> >
     191        AbundantGraph;
     192      AbundantGraph abundant_graph(graph, const_true_map, abundant_arcs );
     193     
     194      //Let's construct the residual graph for the abundant graph
     195      typedef ResGraphWrapper<const AbundantGraph,int,CapacityMap,EdgeIntMap>
     196        ResAbGraph;
     197      //Again uncapacitated
     198      ResAbGraph res_ab_graph(abundant_graph, flow);
     199     
     200      //We need things for the bfs
     201      typename ResAbGraph::NodeMap<bool> bfs_reached(res_ab_graph);
     202      typename ResAbGraph::NodeMap<typename ResAbGraph::Edge>
     203        bfs_pred(res_ab_graph);
     204      NullMap<typename ResAbGraph::Node, int> bfs_dist_dummy(res_ab_graph);
     205      //We want to run bfs-es (more) on this graph 'res_ab_graph'
     206      Bfs < ResAbGraph ,
     207        typename ResAbGraph::NodeMap<bool>,
     208        typename ResAbGraph::NodeMap<typename ResAbGraph::Edge>,
     209        NullMap<typename ResAbGraph::Node, int> >
     210        bfs(res_ab_graph, bfs_reached, bfs_pred, bfs_dist_dummy);
    181211     
    182212      ModLengthMap mod_length(res_graph, length, potential);
     
    206236              //Merge
    207237              if (a != b){
    208                 //Find path and augment
    209                 //!!!
    210                 //Find path and augment
    211238                abundant_components.join(a,b);
     239                //We want to push the smaller
     240                //Which has greater absolut value excess/deficit
     241                Node root=(abs(excess_deficit[a])>abs(excess_deficit[b]))?a:b;
     242                //Which is the other
     243                Node non_root = ( a == root ) ? b : a ;
     244                abundant_components.makeRep(root);
     245                SupplyDemand qty_to_augment = abs(excess_deficit[non_root]);
     246                //Push the positive value
     247                if (excess_deficit[non_root] < 0)
     248                  swap(root, non_root);
     249                //If the non_root node has excess/deficit at all
     250                if (qty_to_augment>0){
     251                  //Find path and augment
     252                  bfs.run(non_root);
     253                  //root should be reached
     254                 
     255                  //Augmenting on the found path
     256                  Node n=root;
     257                  ResGraphEdge e;
     258                  while (n!=non_root){
     259                    e = bfs_pred(n);
     260                    n = res_graph.tail(e);
     261                    res_graph.augment(e,qty_to_augment);
     262                  }
     263         
     264                  //We know that non_root had positive excess
     265                  excess_nodes[non_root] -= qty_to_augment;
     266                  //But what about root node
     267                  //It might have been positive and so became larger
     268                  if (excess_deficit[root]>0){
     269                    excess_nodes[root] += qty_to_augment;
     270                  }
     271                  else{
     272                    //Or negative but not turned into positive
     273                    deficit_nodes[root] -= qty_to_augment;
     274                  }
     275
     276                  //Update the excess_deficit map
     277                  excess_deficit[non_root] -= qty_to_augment;
     278                  excess_deficit[root] += qty_to_augment;
     279
     280                 
     281                }
    212282              }
    213283              //What happens to i?
    214               nonabundant_arcs.erase(i);
     284              //Marci and Zsolt says I shouldn't do such things
     285              nonabundant_arcs.erase(i++);
     286              abundant_arcs[i] = true;
    215287            }
    216288            else
     
    223295        SupplyDemand max_excess = excess_nodes[s];
    224296        Node t = deficit_nodes.top();
    225         if (max_excess < dificit_nodes[t]){
    226           max_excess = dificit_nodes[t];
    227         }
    228 
    229 
    230         while(max_excess > 0){
    231 
     297        if (max_excess < deficit_nodes[t]){
     298          max_excess = deficit_nodes[t];
     299        }
     300
     301
     302        while(max_excess > (n-1)*delta/n){
     303         
    232304         
    233305          //s es t valasztasa
    234 
     306         
    235307          //Dijkstra part       
    236308          dijkstra.run(s);
    237 
     309         
    238310          /*We know from theory that t can be reached
    239311          if (!dijkstra.reached(t)){
     
    264336            */
    265337          }
     338         
     339          //Update the excess_deficit map
     340          excess_deficit[s] -= delta;
     341          excess_deficit[t] += delta;
     342         
    266343
    267344          //Update the excess_nodes heap
     
    286363          }
    287364          //Dijkstra part ends here
     365         
     366          //Choose s and t again
     367          s = excess_nodes.top();
     368          max_excess = excess_nodes[s];
     369          t = deficit_nodes.top();
     370          if (max_excess < deficit_nodes[t]){
     371            max_excess = deficit_nodes[t];
     372          }
     373
    288374        }
    289375
     
    298384        //Update the max_excess
    299385        max_excess = 0;
    300         FOR_EACH_LOC(typename Graph::NodeIt, n, G){
     386        FOR_EACH_LOC(typename Graph::NodeIt, n, graph){
    301387          if (max_excess < excess_deficit[n]){
    302388            max_excess = excess_deficit[n];
     
    337423      Length mod_pot;
    338424      Length fl_e;
    339       FOR_EACH_LOC(typename Graph::EdgeIt, e, G){
     425      FOR_EACH_LOC(typename Graph::EdgeIt, e, graph){
    340426        //C^{\Pi}_{i,j}
    341         mod_pot = length[e]-potential[G.head(e)]+potential[G.tail(e)];
     427        mod_pot = length[e]-potential[graph.head(e)]+potential[graph.tail(e)];
    342428        fl_e = flow[e];
    343429        //      std::cout << fl_e << std::endl;
Note: See TracChangeset for help on using the changeset viewer.