// -*- c++ -*- #ifndef LEMON_MIN_COST_GEN_FLOW_H #define LEMON_MIN_COST_GEN_FLOW_H #include //#include #include #include //#include //#include //#include #include //#include //#include #include <../merge_node_graph_wrapper.h> #include namespace lemon { template class PrimalMap { protected: LPSolverWrapper* lp; EdgeIndexMap* edge_index_map; public: PrimalMap(LPSolverWrapper& _lp, EdgeIndexMap& _edge_index_map) : lp(&_lp), edge_index_map(&_edge_index_map) { } double operator[](Edge e) const { return lp->getPrimal((*edge_index_map)[e]); } }; // excess: rho-delta template , typename LCapMap=typename Graph::template EdgeMap, typename CapMap=typename Graph::template EdgeMap, typename FlowMap=typename Graph::template EdgeMap, typename CostMap=typename Graph::template EdgeMap > class MinCostGenFlow { protected: const Graph& g; const Excess& excess; const LCapMap& lcapacity; const CapMap& capacity; FlowMap& flow; const CostMap& cost; public: MinCostGenFlow(const Graph& _g, const Excess& _excess, const LCapMap& _lcapacity, const CapMap& _capacity, FlowMap& _flow, const CostMap& _cost) : g(_g), excess(_excess), lcapacity(_lcapacity), capacity(_capacity), flow(_flow), cost(_cost) { } bool feasible() { // std::cout << "making new vertices..." << std::endl; typedef ListGraph Graph2; Graph2 g2; typedef MergeEdgeGraphWrapper GW; // std::cout << "merging..." << std::endl; GW gw(g, g2); typename GW::Node s(INVALID, g2.addNode(), true); typename GW::Node t(INVALID, g2.addNode(), true); typedef SmartGraph Graph3; // std::cout << "making extender graph..." << std::endl; typedef NewEdgeSetGraphWrapper2 GWW; // { // checkConcept(); // } GWW gww(gw); typedef AugmentingGraphWrapper GWWW; GWWW gwww(gw, gww); // std::cout << "making new edges..." << std::endl; typename GWWW::template EdgeMap translated_cap(gwww); for (typename GW::EdgeIt e(gw); e!=INVALID; ++e) translated_cap.set(typename GWWW::Edge(e,INVALID,false), capacity[e]-lcapacity[e]); Num expected=0; // std::cout << "making new edges 2..." << std::endl; for (typename Graph::NodeIt n(g); n!=INVALID; ++n) { Num a=0; for (typename Graph::InEdgeIt e(g, n); e!=INVALID; ++e) a+=lcapacity[e]; for (typename Graph::OutEdgeIt e(g, n); e!=INVALID; ++e) a-=lcapacity[e]; if (excess[n]>a) { typename GWW::Edge e= gww.addEdge(typename GW::Node(n,INVALID,false), t); translated_cap.set(typename GWWW::Edge(INVALID, e, true), excess[n]-a); } if (excess[n] translated_flow(gwww, 0); Preflow preflow(gwww, s, t, translated_cap, translated_flow); preflow.run(); // std::cout << "fv: " << preflow.flowValue() << std::endl; // std::cout << "expected: " << expected << std::endl; for (typename Graph::EdgeIt e(g); e!=INVALID; ++e) { typename GW::Edge ew(e, INVALID, false); typename GWWW::Edge ewww(ew, INVALID, false); flow.set(e, translated_flow[ewww]+lcapacity[e]); } return (expected>=preflow.flowValue()); } void run() { LPSolverWrapper lp; lp.setMinimize(); typedef LPSolverWrapper::ColIt ColIt; typedef LPSolverWrapper::RowIt RowIt; typedef typename Graph::template EdgeMap EdgeIndexMap; EdgeIndexMap edge_index_map(g); PrimalMap lp_flow(lp, edge_index_map); for (typename Graph::EdgeIt e(g); e!=INVALID; ++e) { ColIt col_it=lp.addCol(); edge_index_map.set(e, col_it); if (lcapacity[e]==capacity[e]) lp.setColBounds(col_it, LPX_FX, lcapacity[e], capacity[e]); else lp.setColBounds(col_it, LPX_DB, lcapacity[e], capacity[e]); lp.setObjCoef(col_it, cost[e]); } for (typename Graph::NodeIt n(g); n!=INVALID; ++n) { typename Graph::template EdgeMap coeffs(g, 0); for (typename Graph::InEdgeIt e(g, n); e!=INVALID; ++e) coeffs.set(e, coeffs[e]+1); for (typename Graph::OutEdgeIt e(g, n); e!=INVALID; ++e) coeffs.set(e, coeffs[e]-1); RowIt row_it=lp.addRow(); typename std::vector< std::pair > row; //std::cout << "node:" <