tpl_maxflow_mincost.H

# ifndef TPL_MAXFLOW_MINCOST_H
# define TPL_MAXFLOW_MINCOST_H

# include <Bellman_Ford.H>
# include <tpl_netgraph.H>

namespace Aleph {


    template <typename Arc_Info, typename F_Type = double>
class Net_Cost_Arc : public Net_Arc<Arc_Info, F_Type>
{ 
public:

  typedef F_Type Flow_Type;
  Flow_Type cost;

  Flow_Type flow_cost() const { return this->flow*cost; }

  Net_Cost_Arc() :  cost(0) { /* empty */ }

  Net_Cost_Arc(const Arc_Info & info) 
    :  Net_Arc<Arc_Info, F_Type>(info), cost(0) 
  { /* empty */ }
};

    template <class NodeT, class ArcT>
class Net_Max_Flow_Min_Cost : public Net_Graph<NodeT, ArcT>
{
public:

  typedef Net_Graph<NodeT, ArcT> Net;

  typedef typename Net::Digraph Digraph;

  typedef Net_Max_Flow_Min_Cost<NodeT, ArcT> Net_MFMC;

  typedef ArcT Arc;

  typedef NodeT Node;

  typedef typename Arc::Flow_Type Flow_Type;

  typedef typename Node::Node_Type Node_Type;

  typedef typename Arc::Arc_Type Arc_Type;

  Flow_Type & get_cost(Arc * a) { return a->cost; }

  const Flow_Type & flow_cost(Arc * a) const { return a->flow_cost(); }

  virtual Arc * insert_arc(Node * src_node, Node * tgt_node,
                           const Flow_Type & cap, const Flow_Type & cost)
  {
    Arc * a = Net::insert_arc(src_node, tgt_node, Arc_Type(), cap, 0);
    a->cost = cost;
    return a;
  }

  virtual Arc * insert_arc(Node * src_node, Node * tgt_node)
  {
    Arc * a = Net::insert_arc(src_node, tgt_node, Arc_Type());
    a->cost = 0;
    return a;
  }

  virtual Arc * insert_arc(Node * src_node, Node * tgt_node, 
                           const Arc_Type & arc_info)
  {
    Arc * a = Net::insert_arc(src_node, tgt_node, arc_info);
    a->cost = 0;
    return a;
  }

  Flow_Type compute_flow_cost()
  {
    Flow_Type total = 0;
    for (Arc_Iterator<Net_MFMC> it(*this); it.has_current(); it.next())
      {
        Arc * a = it.get_current();
        if (not a->is_residual)
          total += a->flow_cost();
      }

    return total;
  }
};

    template <class Net> 
class Access_Cost
{
public:

  typedef typename Net::Flow_Type Distance_Type;

  static const Distance_Type Max_Distance;

  static const Distance_Type Zero_Distance;

  typename Net::Flow_Type operator () (typename Net::Arc * a) const
  {
    return a->is_residual ? -((typename Net::Arc *) a->img_arc)->cost : a->cost;
  }
};

  template <class Net> 
const typename Access_Cost<Net>::Distance_Type 
  Access_Cost<Net>::Zero_Distance = 0;


    template <class Net>
  const typename Access_Cost<Net>::Distance_Type 
Access_Cost<Net>::Max_Distance = 
    numeric_limits<typename Access_Cost::Distance_Type>::max();


    template <class Net, template <class> class Max_Flow_Algo>
void max_flow_min_cost_by_cycle_canceling
  (Net & net, bool leave_residual = false)
{    // tipos sinónimos para hacer más simple la lectura
  typedef Aleph::less<typename Access_Cost<Net>::Distance_Type> Cmp;
  typedef Aleph::plus<typename Access_Cost<Net>::Distance_Type> Plus;
  typedef Res_F<Net> Filter;

  Max_Flow_Algo <Net> () (net, true); // obtiene flujo máximo inicial

  Path<Net> cycle(net); // aquí se guardan ciclos negativos
      // eliminar ciclos negativos mientras éstos existan 
  while (Bellman_Ford <Net, Access_Cost<Net>, Filter> (net)(cycle))
    if (increase_flow <Net> (net, cycle) == 0)
      break;

  if (not leave_residual)
    net.unmake_residual_net();
}

    template <class Net, template <class> class Cost> static 
typename Net::Flow_Type search_max_arc_cost(Net & net)
{
  Arc_Iterator<Net> it(net);
  typename Net::Flow_Type max_cost = Cost <Net> () (it.get_current());
  for (it.next(); it.has_current(); it.next())
    {
      typename Net::Flow_Type curr_cost = Cost <Net> () (it.get_current());
      if (curr_cost > max_cost)
        max_cost = curr_cost;
    }

  return max_cost;
}

    template <class Net> static 
typename Net::Flow_Type compute_max_possible_flow(Net & net)
{
  typename Net::Node * source = net.get_source();
  typename Net::Node * sink   = net.get_sink();

  return std::min(net.get_out_cap(source), net.get_in_cap(sink));
}

    template <class Net, template <class> class Cost> static
typename Net::Arc * create_dummy_arc(Net & net)
{
  const typename Net::Flow_Type max_cost = 
    net.get_num_nodes() * search_max_arc_cost <Net, Cost> (net);

  net.make_residual_net();

  typename Net::Node * src = net.get_source();
  typename Net::Node * tgt = net.get_sink();
  const typename Net::Flow_Type max_flow = 
    std::min(net.get_out_cap(src), net.get_in_cap(tgt));
  typename Net::Net::Digraph * digraph = &net;
  typename Net::Arc * a = digraph->insert_arc(src, tgt); 
  net.get_cookie() = a;

  try 
    {
      typename Net::Arc * img = digraph->insert_arc(tgt, src); 

      digraph->disconnect_arc(a); // no supra-arco, dejarlo por coste

      a->is_residual = false;
      a->img_arc     = img;
      a->cap         = max_flow;
      a->cost        = max_cost;
      a->flow        = 0;

      img->is_residual = true;
      img->img_arc     = a;
      img->cap         = max_flow;
      img->cost        = max_cost;
      img->flow        = 0;

      return img;

    }
  catch (...)
    {
      digraph->remove_arc(a);
      throw;
    }
}

    template <class Net> static
void destroy_dummy_arc(Net & net)
{
  typename Net::Arc * a = (typename Net::Arc *) net.get_cookie();
  typename Net::Arc * img = (typename Net::Arc *) a->img_arc;
  typename Net::Net::Digraph * digraph = &net;

  digraph->connect_arc(a);
  digraph->remove_arc(img);
  digraph->remove_arc(a);

  net.unmake_residual_net();
}

    template <class Net>
void max_flow_min_cost_by_cycle_canceling(Net & net)
{
  typedef Aleph::less<typename Access_Cost<Net>::Distance_Type> Cmp;
  typedef Aleph::plus<typename Access_Cost<Net>::Distance_Type> Plus;
  typedef Res_F<Net> Res_Filter;

  create_dummy_arc <Net, Access_Cost> (net);

  Path<Net> cycle(net); // aquí se guardan ciclos negativos
  while (Bellman_Ford_Negative_Cycle 
         <Net, Access_Cost<Net>, Cmp,Plus, Res_Filter> () (net, cycle))
    if (increase_flow <Net> (net, cycle) == 0)

      {
        ((typename Net::Arc *) net.get_cookie())->cost = 0;
        break;
      }

  destroy_dummy_arc(net);
}

    template <class Net>
class Max_Flow_Min_Cost
{
public:

  void by_cycle_canceling_and_dummy_arc(Net & net) const
  {
    max_flow_min_cost_by_cycle_canceling(net);
  }

      template <template <class> class Max_Flow_Algo>
  void by_cycle_canceling(Net & net, bool leave_residual = false)
  {
    max_flow_min_cost_by_cycle_canceling<Net, Max_Flow_Algo>
      (net, leave_residual);
  }
};


} // end namespace Aleph

# endif // TPL_MAXFLOW_MINCOST_H