max_flow0.t

A Library of Efficient Data Types and Algorithms,封装了常用的ADT及其相关算法的软件包

void phase1(const graph& G, node s, node t,
            const cap_array& cap,
            flow_array& flow,
            excess_array& excess,
            succ_array& succ)
{ 
  int n = G.number_of_nodes();
  int m = G.number_of_edges();

  int heuristic  = (h > 0) ? int(h*m) : int(-h*n);
  int limit_heur = heuristic;

  int relabel_count = 0;
  int push_count = 0;
  int insp_count = 0;

  edge e;
  forall_out_edges(e,s) 
  { node u = target(G,e,s);
    NT c = cap[e];
    flow[e] = c;
    excess[s] -= c;
    excess[u] += c;
  }

  dist_array dist;
  dist.use_node_data(G);

  node_level_queue<graph,succ_array>  U(G,succ);
  int* count = new int[n];

  compute_dist0(G,s,t,excess,succ,dist,U,count);

  excess[t] += 1;  // prevents t from beeing inserted (again)

  node v;
  while ((v = U.del_max()) != t)
  {
    NT  ev = excess[v]; 
    int dv = dist[v]; 

    int  dmin = n;    
    edge emin = 0;
    NT   rmin = 0;

    edge e;
    forall_out_edges(e,v)
    { insp_count++;
      NT fe = flow[e];
      NT rc = cap[e] - fe;
      if (rc == 0) continue;
      node w = target(G,e,v);
      int dw = dist[w];
      if (dw == dv-1)
      { push_count++;
        NT ew = excess[w];
        if (ew == 0) U.insert_non_max(w,dw);
        if (ev <= rc) 
        { flow[e] = fe + ev; 
          excess[w] = ew + ev;
          excess[v] = 0; 
          goto NEXT_NODE;
         }
        flow[e] = fe + rc;
        excess[w] = ew + rc;
        ev -= rc;
        continue;
       }
      if (v != w && dw < dmin) { dmin = dw; emin = e; rmin = rc; }
    }

    forall_in_edges(e,v)
    { insp_count++;
      NT fe = flow[e];
      if (fe == 0) continue;
      node w = source(G,e,v);
      int dw = dist[w];
      if (dw == dv-1)
      { push_count++;
        NT ew = excess[w];
        if (ew == 0) U.insert_non_max(w,dw);
        if (ev <= fe) 
        { flow[e] = fe - ev; 
          excess[w] = ew + ev;
          excess[v] = 0; 
          goto NEXT_NODE;
         }
        flow[e] = 0;
        excess[w] = ew + fe;
        ev -= fe;
        continue;
       }
      if (v != w && dw < dmin) { dmin = dw; emin = e; rmin = -fe; }
    } 

    excess[v] = ev;

    // remaining excess at v
    // relabel vertex v (i.e. update dist[v]) because all
    // admissible edges in the residual graph have been saturated 

    if (--count[dv] == 0) // gap
    { handle_gap(G,v,dist,succ,count);
      goto NEXT_NODE;
     }

    relabel_count++;

    dv = dmin + 1;
    dist[v] = dv;

    if (dv >= n) goto NEXT_NODE;

    count[dv]++;

    if (ev <= rmin || ev <= -rmin)
    { node w = G.opposite(emin,v);
      if (excess[w] == 0) U.insert(w,dmin);
      excess[w]  += ev;
      flow[emin] += (rmin > 0) ? ev : -ev;
      excess[v]  = 0;
      goto NEXT_NODE;
     }


    if ((h > 0 && insp_count    >= limit_heur) || 
        (h < 0 && relabel_count >= limit_heur))
    { 
      U.clear();
      compute_dist1(G,s,t,cap,flow,excess,succ,dist,U,count);
      limit_heur += heuristic;
      goto NEXT_NODE;
     }

   U.insert_max(v,dv);

NEXT_NODE:;
  }

  excess[t] -= 1;

  num_pushes[1] += push_count;
  num_relabels[1] += relabel_count;
  num_inspections[1] += insp_count;

  delete[] count;

  fval = excess[t];
}



template<class cap_array, class flow_array>
void phase2(const graph& G, node s, node t,
            const cap_array& cap,
            flow_array& flow,
            excess_array& excess,
            succ_array& succ)
{ 
  int n = G.number_of_nodes();
  int m = G.number_of_edges();
  
  int heuristic  = (h > 0) ? int(h*m) : int(-h*n);
  int limit_heur = heuristic;

  int relabel_count = 0;
  int push_count = 0;
  int insp_count = 0;

  node_level_queue<graph,succ_array>  U(G,succ);

  dist_array dist;
  dist.use_node_data(G);

  compute_dist2(G,s,t,flow,excess,succ,dist,U);

  node v;
  while ((v = U.del_max()) != t)
  {
    NT  ev = excess[v]; 
    int dv = dist[v]; 
    int dmin = n;

    edge e;
    forall_in_edges(e,v)
    { insp_count++;
      NT fe = flow[e];
      if (fe == 0) continue;
      node w = source(G,e,v);
      int dw = dist[w];
      if (dw == dv-1)
       { push_count++;
         NT ew = excess[w];
         if (ew == 0) U.insert_non_max(w,dw);
         if (ev <= fe) 
         { flow[e] = fe - ev; 
           excess[w] = ew + ev;
           ev = 0; 
           break;
          }
         else
         { flow[e] = 0;
           excess[w] = ew + fe;
           ev -= fe;
          }
         continue;
        }
      if (dw < dmin) dmin = dw;
    } 
    
    excess[v] = ev;
    if (ev == 0) continue;

    if ((h > 0 && insp_count    >= limit_heur) || 
        (h < 0 && relabel_count >= limit_heur))
      { 
        U.clear();
        compute_dist2(G,s,t,flow,excess,succ,dist,U);
        limit_heur += heuristic;
       }
    else
      { relabel_count++;
        dv = dmin+1;
        dist[v] = dv;
        U.insert_max(v,dv);
       }
  }

  num_pushes[2] += push_count;
  num_relabels[2] += relabel_count;
  num_inspections[2] += insp_count;

}






template<class cap_array, class flow_array>
NT run(const graph& G, node s, node t, 
       const cap_array& cap, 
       flow_array& flow)
{ 
  if (s == t) error_handler(1,"MAXFLOW: source == sink");

  float T = used_time();

  reset_counters();

  excess_array excess;
  excess.use_node_data(G,0);

  succ_array succ;
  succ.use_node_data(G);

  flow.init(G,0);

  phase1(G,s,t,cap,flow,excess,succ);
  cputime[1] = used_time(T);

  phase2(G,s,t,cap,flow,excess,succ);
  cputime[2] = used_time(T);

  num_pushes[0]      = num_pushes[1]      + num_pushes[2];
  num_relabels[0]    = num_relabels[1]    + num_relabels[2];
  num_updates[0]     = num_updates[1]     + num_updates[2];
  num_gap_nodes[0]   = num_gap_nodes[1]   + num_gap_nodes[2];
  num_inspections[0] = num_inspections[1] + num_inspections[2];

  cputime[0] = cputime[1] + cputime[2];

  return fval;
}


template<class cap_array, class flow_array>
NT operator()(const graph& G, node s, node t, 
              const cap_array& cap, flow_array& flow) 
{ 
  return run(G,s,t,cap,flow); 
 }




void statistics(ostream& out)
{ 
  out << endl;

  for(int i=1; i<=2; i++)
  { out << 
    string("%8d pushes %7d relabels %3d updates %9d insp %5d gaps  %.2f s",
            pushes(i), relabels(i), updates(i), inspections(i), gap_nodes(i),
            cputime[i]) << endl;
  }

  cout << string("update: %.2f   gaps: %.2f",update_time,gap_time) << endl;

  out << endl;
}




template<class cap_array, class flow_array>
bool check(const graph& G, node s, node t, 
           const cap_array& cap, const flow_array& flow,
           string& msg)
{
  edge e;
  forall_edges(e,G) 
  { if (flow[e] < 0 || flow[e] > cap[e]) 
    { msg = string("illegal flow value for edge %d",G.index(e));
      return false;
     }
   }
  
  excess_array excess(G);

  node v;
  forall_nodes(v,G) excess[v] = 0;

  forall_nodes(v,G) 
  { edge e;
    forall_out_edges(e,v) 
    { node w = target(G,e,v);
      excess[v] -= flow[e]; 
      excess[w] += flow[e];
     }
  }

  forall_nodes(v,G) 
  { if (v == s  || v == t || excess[v] == 0) continue;
    msg = "node with non-zero excess";
    return false;
  }

/*
  if (fval != excess[t])
  { msg = "fval != excess[t]";
    return false;
   }
*/
  
  node_array<bool,graph> reached(G,false);
  reached[s] = true;

  list<node> Q;
  Q.append(s); 

  while ( !Q.empty() )
  { node v = Q.pop(); 
    forall_out_edges(e,v) 
    { node w = target(G,e,v);
      if ( flow[e] < cap[e] && !reached[w] ) 
      { reached[w] = true; 
        Q.append(w); 
       }
    }
    forall_in_edges(e,v) 
    { node w = source(G,e,v);
      if ( flow[e] > 0 && !reached[w] ) 
      { reached[w] = true; 
        Q.append(w); 
       }
    }
  }

  if (reached[t])
  { msg = "t is reachable in G_f";
    return false;
   }

  forall_nodes(v,G)
  { if (v != s && v != t && excess[v] != 0)
    { msg = string("node %d has non-zero excess",G.index(v));
      return false;
     }
   }

  return true;
}


template <class cap_array, class flow_array>
void print(const graph& G, node s, node t, 
           const cap_array& cap, const flow_array& flow)
{
  cout << "s = " << G.index(s) << endl;
  cout << "t = " << G.index(t) << endl;
  cout << endl;

  node v;
  forall_nodes(v,G)
  { cout << string("%2d :",G.index(v));
    edge e;
    forall_out_edges(e,v)
      cout << string(" --%d/%d-->%d",flow[e],cap[e],G.index(target(G,e,v)));
    cout << endl;
   }
  cout << endl;
}
  
 
};


template<class NT>
NT MAX_FLOW_T(const graph& G, node s, node t, const edge_array<NT>& cap, 
                                              edge_array<NT>& flow, 
                                              int& num_pushes, 
                                              int& num_edge_inspections,
                                              int& num_relabels,
                                              int& num_global_relabels,
                                              int& num_gaps,
                                              float h)
{ max_flow<NT,graph> mf;
  mf.set_heuristic(h);
  NT f = mf.run(G,s,t,cap,flow);
  num_pushes = mf.pushes(0);
  num_relabels = mf.relabels(0);
  num_global_relabels = mf.updates(0);
  num_edge_inspections = mf.inspections(0);
  num_gaps = mf.gap_nodes(0);
  return f;
}



template<class NT>
NT MAX_FLOW_T(const graph& G, node s, node t, const edge_array<NT>& cap, 
                                              edge_array<NT>& flow)
{ max_flow<NT,graph> mf;
  return mf.run(G,s,t,cap,flow);
}


template <class NT>
bool CHECK_MAX_FLOW_T(const graph& G, node s, node t,
                      const edge_array<NT>& cap, const edge_array<NT>& flow)
{ max_flow<NT,graph> mf;
  string msg;
  bool ok = mf.check(G,s,t,cap,flow,msg);
  if (!ok) error_handler(1,string("CHECK_MAX_FLOW_T: ") + msg + ".");
  return ok;
}
  


LEDA_END_NAMESPACE

#if LEDA_ROOT_INCL_ID == 450450
#undef LEDA_ROOT_INCL_ID
#include <LEDA/POSTAMBLE.h>
#endif