// -*- c++ -*-

 #include <iostream>

 #include <fstream>

 #include <lemon/smart_graph.h>

 #include <lemon/list_graph.h>

 #include <lemon/dimacs.h>

 #include <lemon/time_measure.h>

 #include <lp_solver_wrapper_3.h>

 using std::cout;

 using std::endl;

 using namespace lemon;

 template<typename Edge, typename EdgeIndexMap> 

 class PrimalMap {

 protected:

   LPGLPK* lp;

   EdgeIndexMap* edge_index_map;

 public:

   PrimalMap(LPGLPK& _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]);

   }

 };

 // Use a DIMACS max flow file as stdin.

 // max_flow_expresion < dimacs_max_flow_file

 int main(int, char **) {

   typedef ListGraph Graph;

   typedef Graph::Node Node;

   typedef Graph::Edge Edge;

   typedef Graph::EdgeIt EdgeIt;

   Graph g;

   Node s, t;

   Graph::EdgeMap<int> cap(g);

   readDimacs(std::cin, g, cap, s, t);

   Timer ts;

   typedef LPGLPK LPSolver;

   LPSolver lp;

   lp.setMaximize();

   typedef LPSolver::ColIt ColIt;

   typedef LPSolver::RowIt RowIt;

   typedef Graph::EdgeMap<ColIt> EdgeIndexMap;

   EdgeIndexMap edge_index_map(g);

   PrimalMap<Edge, EdgeIndexMap> lp_flow(lp, edge_index_map);

   // capacity function

   for (Graph::EdgeIt e(g); e!=INVALID; ++e) {

     ColIt col_it=lp.addCol();

     edge_index_map.set(e, col_it);

     if (cap[e]==0)

       lp.setColBounds(col_it, LPSolver::FIXED, 0, cap[e]);

     else 

       lp.setColBounds(col_it, LPSolver::DOUBLE, 0, cap[e]);

   }

   for (Graph::NodeIt n(g); n!=INVALID; ++n) {

     LPSolver::Expression expr;

     for (Graph::OutEdgeIt e(g, n); e!=INVALID; ++e)

       expr+=edge_index_map[e];

     for (Graph::InEdgeIt e(g, n); e!=INVALID; ++e)

       expr-=edge_index_map[e];

     // cost function

     if (n==s) {

       lp.setObjCoeffs(expr);      

     }

     // flow conservation

     if ((n!=s) && (n!=t)) {

       RowIt row_it=lp.addRow();

       lp.setRowCoeffs(row_it, expr);

       lp.setRowBounds(row_it, LPSolver::FIXED, 0.0, 0.0);

     }

   }

   lp.solveSimplex();

   //std::cout << lp.colNum() << std::endl;

   //std::cout << lp.rowNum() << std::endl;

   //std::cout << "flow value: "<< lp.getObjVal() << std::endl;

   for (Graph::EdgeIt e(g); e!=INVALID; ++e) 

     flow.set(e, lp_flow[e]);

 }