mwb_matching.t

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

/*******************************************************************************
+
+  LEDA 4.5  
+
+
+  mwb_matching.t
+
+
+  Copyright (c) 1995-2004
+  by Algorithmic Solutions Software GmbH
+  All rights reserved.
+ 
*******************************************************************************/

// $Revision: 1.4 $  $Date: 2004/02/06 11:20:23 $


#if !defined(LEDA_ROOT_INCL_ID)
#define LEDA_ROOT_INCL_ID 450466
#include <LEDA/PREAMBLE.h>
#endif

#include <LEDA/mwb_matching.h>
#include <LEDA/node_pq.h>
#include <LEDA/stack.h>
#include <LEDA/std/assert.h>

LEDA_BEGIN_NAMESPACE


// mod. 12/2000 by GS
bool PRUNE = true;


////////////// MAX WEIGHT BIPARTITE MATCHING

template <class NT>
__temp_func_inline
list<edge> MAX_WEIGHT_BIPARTITE_MATCHING_T(graph& G,
                               const edge_array<NT>& c)
{ node_array<NT> pot(G);
  list<edge> M = MAX_WEIGHT_BIPARTITE_MATCHING_T(G,c,pot);
#ifdef TEST
  assert(CHECK_MWBM_T(G,c,M,pot));
#endif
  return M;
}


template <class NT>
__temp_func_inline
list<edge> MAX_WEIGHT_BIPARTITE_MATCHING_T(graph& G,                                    
                      const list<node>& A, const list<node>& B,const edge_array<NT>& c)
{ node_array<NT> pot(G);
  list<edge> M = MAX_WEIGHT_BIPARTITE_MATCHING_T(G,A,B,c,pot);
#ifdef TEST
  assert(CHECK_MWBM_T(G,c,M,pot));
#endif
  return M;
}


template <class NT>
__temp_func_inline
list<edge> MAX_WEIGHT_BIPARTITE_MATCHING_T(graph& G,                      
          const edge_array<NT>& c,node_array<NT>& pot)
{
  list<node> A,B;
  node v; edge e;

  if ( !Is_Bipartite(G,A,B) )
    error_handler(1,"MAX_WEIGHT_BIPARTITE_MATCHING: \
G is not bipartite");

  edge_array<int> edge_number(G); int i = 0;
  forall_nodes(v,G) 
    forall_adj_edges(e,v) edge_number[e] = i++;
  
  list<edge> edges_out_of_B;
  forall(v,B) 
  { list<edge> outedges = G.adj_edges(v);
    edges_out_of_B.conc(outedges);
  }
  forall(e,edges_out_of_B) G.rev_edge(e);

  list<edge> result = MAX_WEIGHT_BIPARTITE_MATCHING_T(G,A,B,c,pot);
  
  forall(e,edges_out_of_B) G.rev_edge(e);

  G.sort_edges(edge_number);
#ifdef TEST
  assert(CHECK_MWBM_T(G,c,result,pot));
#endif

  return result;

}



////////////// MAX ASSIGNMENT

template <class NT>
__temp_func_inline
list<edge> MAX_WEIGHT_ASSIGNMENT_T(graph& G,
                               const edge_array<NT>& c)
{ node_array<NT> pot(G);
  list<edge> M = MAX_WEIGHT_ASSIGNMENT_T(G,c,pot);
  if ( M.length() == 0 ) return M;
#ifdef TEST
  assert(CHECK_MAX_WEIGHT_ASSIGNMENT_T(G,c,M,pot));
#endif
  return M;
}


template <class NT>
__temp_func_inline
list<edge> MAX_WEIGHT_ASSIGNMENT_T(graph& G,                                    
                      const list<node>& A, const list<node>& B,const edge_array<NT>& c)
{ node_array<NT> pot(G);
  list<edge> M = MAX_WEIGHT_ASSIGNMENT_T(G,A,B,c,pot);
  if ( M.length() == 0 ) return M;
#ifdef TEST
  assert(CHECK_MAX_WEIGHT_ASSIGNMENT_T(G,c,M,pot));
#endif
  return M;
}


template <class NT>
__temp_func_inline
list<edge> MAX_WEIGHT_ASSIGNMENT_T(graph& G,                      
          const edge_array<NT>& c,node_array<NT>& pot)
{
  list<node> A,B;
  node v; edge e;

  if ( !Is_Bipartite(G,A,B) )
    error_handler(1,"MAX_WEIGHT_ASSIGNMENT: \
G is not bipartite");

  edge_array<int> edge_number(G); int i = 0;
  forall_nodes(v,G) 
    forall_adj_edges(e,v) edge_number[e] = i++;
  
  list<edge> edges_out_of_B;
  forall(v,B) 
  { list<edge> outedges = G.adj_edges(v);
    edges_out_of_B.conc(outedges);
  }
  forall(e,edges_out_of_B) G.rev_edge(e);

  list<edge> result = MAX_WEIGHT_ASSIGNMENT_T(G,A,B,c,pot); 

  forall(e,edges_out_of_B) G.rev_edge(e);

  G.sort_edges(edge_number);

  if ( result.length() == 0 ) return result;
#ifdef TEST
  assert(CHECK_MAX_WEIGHT_ASSIGNMENT_T(G,c,result,pot));
#endif
  return result;

}

////////////// MIN ASSIGNMENT

template <class NT>
__temp_func_inline
list<edge> MIN_WEIGHT_ASSIGNMENT_T(graph& G,
                               const edge_array<NT>& c)
{ node_array<NT> pot(G);
  list<edge> M = MIN_WEIGHT_ASSIGNMENT_T(G,c,pot);
  if (M.length() == 0) return M;
#ifdef TEST
  assert(CHECK_MIN_WEIGHT_ASSIGNMENT_T(G,c,result,pot));
#endif
  return M;
}


template <class NT>
__temp_func_inline
list<edge> MIN_WEIGHT_ASSIGNMENT_T(graph& G,                                    
                      const list<node>& A, const list<node>& B,const edge_array<NT>& c)
{ node_array<NT> pot(G);
  list<edge> M = MIN_WEIGHT_ASSIGNMENT_T(G,A,B,c,pot);
  if (M.length() == 0) return M;
#ifdef TEST
  assert(CHECK_MIN_WEIGHT_ASSIGNMENT_T(G,c,result,pot));
#endif
  return M;
}


template <class NT>
__temp_func_inline
list<edge> MIN_WEIGHT_ASSIGNMENT_T(graph& G,                      
          const edge_array<NT>& c,node_array<NT>& pot)
{
  list<node> A,B;
  node v; edge e;

  if ( !Is_Bipartite(G,A,B) )
    error_handler(1,"MIN_WEIGHT_ASSIGNMENT: \
G is not bipartite");

  edge_array<int> edge_number(G); int i = 0;
  forall_nodes(v,G) 
    forall_adj_edges(e,v) edge_number[e] = i++;
  
  list<edge> edges_out_of_B;
  forall(v,B) 
  { list<edge> outedges = G.adj_edges(v);
    edges_out_of_B.conc(outedges);
  }
  forall(e,edges_out_of_B) G.rev_edge(e);

  list<edge> result = MIN_WEIGHT_ASSIGNMENT_T(G,A,B,c,pot); 

  forall(e,edges_out_of_B) G.rev_edge(e);

  G.sort_edges(edge_number);

  if (result.length() == 0) return result;
#ifdef TEST
  assert(CHECK_MIN_WEIGHT_ASSIGNMENT_T(G,c,result,pot));
#endif
  return result;

}

LEDA_END_NAMESPACE
#include <LEDA/misc.h>
LEDA_BEGIN_NAMESPACE

template <class NT>
__temp_func_inline
bool CHECK_MWBM_T(const graph& G, const edge_array<NT>& c,
                  const list<edge>& M, const node_array<NT>& pot)
{ node v; edge e;

  string s = "\nCHECK_MWBM_T: ";
  bool res = true;

  // M is a matching
  node_array<int> deg_in_M(G,0);
  forall(e,M) 
  { deg_in_M[G.source(e)]++;
    deg_in_M[G.target(e)]++;
  }
  forall_nodes(v,G) 
    res = res && leda_assert(deg_in_M[v] <= 1,s + "M is not a matching");

  // node potentials are non-negative
  forall_nodes(v,G) 
    res = res && leda_assert(pot[v] >= 0,s + "negative node potential");;

  // edges have non-negative reduced cost
  forall_edges(e,G)
  { node v = G.source(e); node w = G.target(e); 
    res = res && leda_assert(c[e] <= pot[v] + pot[w],s + "negative reduced cost");
  }

  // edges in M have reduced cost equal to zero
  forall(e,M)
  { node v = G.source(e); node w = G.target(e); 
    res = res && leda_assert(c[e] == pot[v] + pot[w],s + "non-tight matching edge");
  }

  // free nodes have potential equal to zero
  forall_nodes(v,G) 
    res = res && leda_assert((deg_in_M[v] != 0 || pot[v] == 0),
                        s + "free node with non-zero potential");

  return res;
}


/*
template <class NodeArray, class EdgeArray>
__temp_func_inline
bool CHECK_ASSIGNMENT_T(const graph& G, const EdgeArray& c,
        const list<edge>& M, const NodeArray& pot, bool max_cost)
*/

template <class NT>
__temp_func_inline
bool CHECK_ASSIGNMENT_T(const graph& G, const edge_array<NT>& c,
        const list<edge>& M, const node_array<NT>& pot, int max_cost)
{ node v; edge e;
  string s = "CHECK_ASSIGNMENT_T";
  if ( max_cost ) s += "(max): "; else s += "(min): ";
  bool res = true;

  node_array<int> deg_in_M(G,0);
  forall(e,M) 
  { deg_in_M[G.source(e)]++;
    deg_in_M[G.target(e)]++;
  }
  forall_nodes(v,G) res = res && leda_assert(deg_in_M[v] <= 1,s + "M is not a matching");
  forall_edges(e,G)
  { node v = G.source(e); node w = G.target(e); 
    if (max_cost)
      { res = res && leda_assert(c[e] <= pot[v] + pot[w],s + "illegal reduced cost"); }
    else
      { res = res && leda_assert(c[e] >= pot[v] + pot[w],s + "illegal reduced cost"); }
     
  }
  forall(e,M)
  { node v = G.source(e); node w = G.target(e); 
    res = res && leda_assert(c[e] == pot[v] + pot[w],s + "non_tight edge in M");
  }
  return res;
}



template <class NT>
__temp_func_inline
bool CHECK_MIN_WEIGHT_ASSIGNMENT_T(const graph& G, const edge_array<NT>& c,
        const list<edge>& M, const node_array<NT>& pot)
{ return CHECK_ASSIGNMENT_T(G, c, M, pot, 0); }

template <class NT>
__temp_func_inline
bool CHECK_MAX_WEIGHT_ASSIGNMENT_T(const graph& G, const edge_array<NT>& c,
        const list<edge>& M, const node_array<NT>& pot)
{ return CHECK_ASSIGNMENT_T(G, c, M, pot, 1); }






inline void augment_path_to(graph& G, node v, const node_array<edge>& pred)
{ edge e = pred[v];
  while (e)
  { G.rev_edge(e);
    e = pred[G.target(e)]; // not source (!!!)
  }
}


template <class NT>
inline void augment(graph& G, node a, const edge_array<NT>& c,
                node_array<NT>& pot, node_array<bool>& free, 
                node_array<edge>& pred, node_array<NT>& dist, 
                node_pq<NT>& PQ)
{ 
  dist[a] = 0;

  node best_node_in_A = a; 
  NT   minA           = pot[a];
  NT   Delta;

  // mod. 12/2000 by GS
  NT   upper_bound    = minA;

  stack<node> RA;  RA.push(a);
  stack<node> RB;

  node a1 = a; edge e;
  
  forall_adj_edges(e,a1)
  { node b = G.target(e);
    NT db = dist[a1] + (pot[a1] + pot[b] - c[e]); 

    // mod. 12/2000 by GS
    if (PRUNE && db >= upper_bound) continue;
    if (free[b]) upper_bound = db;

    if ( pred[b] == nil ) 
    { dist[b] = db; pred[b] = e; RB.push(b);  
      PQ.insert(b,db); 
    }
    else
    if ( db < dist[b] )
    { dist[b] = db; pred[b] = e;
      PQ.decrease_p(b,db); 
    }
  }



  while ( true )
  {     
    node b;
    NT db = 0;
    if (PQ.empty()) b = nil;
    else { b = PQ.del_min(); db = dist[b]; }
       

    
    if ( b == nil || db >= minA ) 
    { Delta = minA;
      
      augment_path_to(G,best_node_in_A,pred);  

      free[a] = false; free[best_node_in_A] = true; // order is important
      break;
 
    }
    else
    { if ( free[b] )
      { Delta = db;
        
        augment_path_to(G,b,pred);
        free[a] = free[b] = false;
        break;
 
      }
      else
      { 
        e = G.first_adj_edge(b); 
        node a1 = G.target(e);
        pred[a1] = e; RA.push(a1);
        dist[a1] = db; 
              
        if (db + pot[a1] < minA)
        { best_node_in_A = a1;
          minA = db + pot[a1];

          // mod. 12/2000 by GS
	  upper_bound = leda_min(upper_bound, minA);
        }

        
        forall_adj_edges(e,a1)
        { node b = G.target(e);
          NT db = dist[a1] + (pot[a1] + pot[b] - c[e]); 

          // mod. 12/2000 by GS
          if (PRUNE && db >= upper_bound) continue;
          if (free[b]) upper_bound = db;

          if ( pred[b] == nil ) 
          { dist[b] = db; pred[b] = e; RB.push(b);  
            PQ.insert(b,db); 
          }
          else
          if ( db < dist[b] )