#include<scil/scil.h>
#include<scil/constraints/atour.h>
#include<scil/constraints/path.h>
//#include<scil/constraints/disjunctive.h>
#include<LEDA/tuple.h>
#include<LEDA/stream.h>
#include<LEDA/string.h>
#include<fstream.h>
#include<LEDA/graph_alg.h>
#include<LEDA/p_queue.h>

#define MULT 2

using namespace SCIL;
using namespace LEDA;

typedef two_tuple<node, node> prec;

LEDA::string read_next_line(istream& is) {
  LEDA::string s="#";
  while((!is.eof()) && ((s.head(1)=="#") || (s.length()==0) || (s.head(1)==" "))) {
    s.read_line(is);
  };
  if(is.eof()) s="";
  return s;
};

class GATOUR : public sym_constraint {
  GRAPH<int,int>& G;
  var_map<edge>& VM;
  row_map<node>& NM;
  int n;
  list<node>* L;
  list<cons_obj*> CL;
  node sour,targ;

public: 
  int addition(int c) {
    int d=MULT*10;
    if(c<d) return 10;
    d*=MULT*10;
    if(c<d) return 8;
    d*=MULT*10;
    if(c<d) return 6;
    d*=MULT*10;
    if(c<d) return 4;
    return 2;
  };

  GATOUR(GRAPH<int,int>& G_, int n_, var_map<edge>& VM_, row_map<node>& NM_) : G(G_), VM(VM_), NM(NM_) {
    n=n_;

    L=new list<node>[n];

    node v;
    forall_nodes(v,G) {
      if(G[v]==0) {
        sour=v;
      } else {
        if(G[v]==n*n) targ=v;
        else
          L[G[v]/n].append(v);
      }
    };
  };

  void MYDFS(graph& H, node s, node_array<bool>& R,edge_array<int>& E, edge_array<int>& F) {
    R[s]=true;
    edge e;
    forall_out_edges(e,s) if((F[e]!=E[e]) && (R[target(e)]==false))
      MYDFS(H,target(e),R,E,F);
    forall_in_edges(e,s) if((F[e]!=0) && (R[source(e)]==false))
      MYDFS(H,source(e),R,E,F);
  };

  status feasible(solution& s) {
    return feasible_solution;

    //return infeasible_solution;
  };

  void CDFS(subproblem& s, GRAPH<int,int>& G, node u, int* v, int nv, int n, int c, int& m, list<edge>& L1, list<edge>& L2) {
    edge e;
    forall_out_edges(e, u) {
      if((target(e)==targ) && (nv==n)) {
        if(c+G[e]<m) {
          m=c+G[e];
          L2.clear();
          edge f;
          forall(f,L1) L2.append(f);
          L2.append(e);
        }
      } else {
        if((s.value(VM[e])>0.001) && (v[G[target(e)]/n]==0)) {
          v[G[target(e)]/n]=1;
          L1.append(e);
          CDFS(s, G, target(e), v, nv+1, n, c+G[e], m, L1, L2);
          L1.pop_back();
          v[G[target(e)]/n]=0;
        }
      }
    };
  };

  status separate(subproblem& s) {
    int ccount=0;
    {
      bool found_solution;
      int v[n];
      edge e,e1,f;
      forall_out_edges(f,sour) if(s.value(VM[f])!=0) {
        found_solution=true;
        solution SOL;
        for(int i=0; i<n; i++) v[i]=0;
        int vi=0;
        int d=0;
        node w=sour;
        while(vi<n-1) {
          vi++;
          e1=nil;
          if(w!=sour) {
            forall_out_edges(e,w) if(target(e)!=targ) 
              if((v[G[target(e)]%n]==0) && ((e1==nil) || (1000*s.value(VM[e])-G[e]>1000*s.value(VM[e1])-G[e1]))) 
                e1=e;
          } else e1=f;
          //assert(e1!=nil);
	  if(e1==nil) { found_solution=false; break; }
          if(vi==1) v[G[target(e1)]/n]=1;
          v[G[target(e1)]%n]=1;
          w=target(e1);
          SOL.set_value(VM[e1],1);
          d+=G[e1];
        }
        forall_out_edges(e,w) if(target(e)==targ) {
          SOL.set_value(VM[e],1);
        }
        //cout<<d<<" sol\n";
        if(found_solution) s.save_solution(SOL);
      }
    }
/*
    {
      int v[n];
      for(int i=0; i<n; i++) v[i]=0;
      int d=10000000;
      list<edge> L1, L2;
      CDFS(s,G,sour,v,0,n,0,d, L1, L2);
      cout<<d<<endl;
      if(d<10000000) {
        solution SOL;
        edge e;
        forall(e,L2) {
          cout<<G[e]<<" ";
          SOL.set_value(VM[e],1);
        };
        cout<<endl;
        s.save_solution(SOL);
      };
    };
*/

    return no_constraint_found;

    status sta=no_constraint_found;

    cons_obj* co;
    forall(co,CL) if(co->violated(s)) {
      sta=constraint_found;
      s.add_basic_constraint(co);
      ccount++;
    }
    if(sta==constraint_found) {
      cout<<"POOL SEPARATION\n";
      return sta;
    };

    node v,u;
    edge e,f;
    node_array<node> GtoH(G);
    graph H;
    edge_array<int> C(H,2*G.number_of_edges()+4*n*n,0);
    var_map<edge> HM;
    node HE[n];
    forall_nodes(v,G) GtoH[v]=H.new_node();
    for(int i=0; i<n; i++) HE[i]=H.new_node();
    forall_edges(e,G) if((s.value(VM[e])>0.001) && (target(e)!=targ)) {
      f=H.new_edge(GtoH[source(e)], GtoH[target(e)]);
      HM[f]=VM[e];
      C[f]=G[e];
      f=H.new_edge(GtoH[target(e)], GtoH[source(e)]);
      HM[f]=VM[e];
      C[f]=G[e];
    }
    forall_in_edges(e,targ) {
      f=H.new_edge(GtoH[source(e)],HE[G[source(e)]%n]);
      HM[f]=VM[e];
      C[f]=G[e];
      f=H.new_edge(HE[G[source(e)]%n],GtoH[targ]);
      HM[f]=VM[e];
      C[f]=G[e];
      f=H.new_edge(HE[G[source(e)]%n],GtoH[source(e)]);
      HM[f]=VM[e];
      C[f]=G[e];
      f=H.new_edge(GtoH[targ],HE[G[source(e)]%n]);
      HM[f]=VM[e];
      C[f]=G[e];
    }

    p_queue<int,cons_obj*> PQ;

    edge_array<int> E(H, H.number_of_edges()+2*n*n,0);
    edge_array<int> F(H, H.number_of_edges()+2*n*n,0);   
    list<edge> EL, FL;
    node_array<bool> R(H);

    forall_edges(e, H) {
      E[e]=leda_max((int) (10000.0*s.value(HM[e]))+addition(C[e]),0);
    }
/*
    forall_nodes(u,G) if(s.value(NM[u])>0.6) {
      forall_nodes(v,G) if((v>u) && (s.value(NM[v])>0.6)) {
        cout<<"test\n";
        int d=MAX_FLOW(H,GtoH[u],GtoH[v],E,F);
        if(d<20000*(s.value(NM[u])+s.value(NM[v])-1)-100) {
          cout<<"cut\n";
          forall_nodes(v, H) R[v]=false;
          MYDFS(H, GtoH[u], R, E, F);
          row r;
          forall_edges(e,G) 
            if(R[GtoH[source(e)]]!=R[GtoH[target(e)]]) r+=VM[e];
          s.add_basic_constraint(r>=NM[u]+NM[v]);
          ccount++;
          sta=constraint_found;
        };
      };
    };
*/

    edge ts=H.new_edge(GtoH[targ],GtoH[sour]);
    E[ts]=20000;
    edge ts1=H.new_edge(GtoH[sour],GtoH[targ]);
    E[ts1]=20000;
    for(int i=0; i<n-1; i++) {
      forall(v, L[i]) {
        EL.append(H.new_edge(HE[i], GtoH[v]));
      }
      forall(e, EL) E[e]=20000;

      for(int j=i+1; j<n; j++) {
        forall(v, L[j]) {
          FL.append(H.new_edge(GtoH[v],HE[j]));
        }
        forall(e, FL) E[e]=20000;

        int d=MAX_FLOW(H,HE[i],HE[j],E,F);

        if(d<=18000) {
          //cout<<d<<" found constraint\n";

          forall_nodes(v, H) R[v]=false;
          MYDFS(H, HE[i], R, E, F);
          row r;
          forall_edges(e,G) 
            if((R[GtoH[source(e)]]) && (!R[GtoH[target(e)]])) {
              r+=VM[e];
            }
          s.add_basic_constraint(r>=1);
          ccount++;
          //if((20000-d)/r.size()>0) PQ.insert(-(20000-d)/r.size(),r>=2);
          //else cout<<d<<" "<<r.size()<<" denied\n";
          sta=constraint_found;
        };

        forall(e,FL) H.del_edge(e);
        FL.clear();
      }
      forall(e,EL) H.del_edge(e);
      EL.clear();
    };     
    if(sta==constraint_found) cout<<"YES\n"; else cout<<"NO\n";
/*
    h_array<var, int> M;
    cons_obj* c;
    row_entry re;
    bool b;
    int i=0;
    while(!PQ.empty()) {
      c=PQ[PQ.find_min()];
      if(!c->violated(s)) cout<<"ERROR\n";
      //cout<<c->slack(s)<<endl;
      PQ.del_min();
      row r;
      c->non_zero_entries(r);
      //cout<<r.size()<<endl;
      b=true;
      //forall(re,r) if(M[re.Var]>=3) b=false;
      if(i<1500) {
        //forall(re,r) M[re.Var]++;
        s.add_basic_constraint(c);
        //CL.append(new cons_obj(*c));
        count++;
      } else delete c;
    }
*/
    cout<<ccount<<" constraints separated\n";

    return sta;
  };
};


int main(int argc, char** argv) {

  // READ_GRAPH

  ifstream is(argv[1]);

  GRAPH<int,int> G;

  LEDA::string_istream si(read_next_line(is));
  int n,i,j,k,l;

  si>>n;
  node N[n][n];
  LEDA::string Names[n];

  cout<<n<<" nodes in the graph\n";

  for(int i=0; i<n; i++) Names[i]=read_next_line(is);

  node s=G.new_node(0);
  node t=G.new_node(n*n);

  for(int i=0; i<n; i++)
    for(int j=0; j<n; j++)
      if(i!=j) { 
        N[i][j]=G.new_node(n*i+j);
        G.new_edge(N[i][j],t,0);
      }

  while(!is.eof()) {
    string_istream si(read_next_line(is));
    si>>i;
    if(!si.eof()) {
      si>>j; si>>k; 
      leda_string ls;
      char ch;
      ls.read_line(si);
      l=0;
      bool neg=false;
      for(int in=0; in<ls.length(); in++) {
        ch=ls[in];
        if(ch=='-') neg=true;
        if(ch==',') l*=MULT;
        if(ch=='0') l=10*l  ;
        if(ch=='1') l=10*l+1;
        if(ch=='2') l=10*l+2;
        if(ch=='3') l=10*l+3;
        if(ch=='4') l=10*l+4;
        if(ch=='5') l=10*l+5;
        if(ch=='6') l=10*l+6;
        if(ch=='7') l=10*l+7;
        if(ch=='8') l=10*l+8;
        if(ch=='9') l=10*l+9;
      }
      if(neg) l=-l;
      if(i!=-1) 
        G.new_edge(N[i][j],N[j][k],l);
      else
        G.new_edge(s,N[j][k],l);
    }
  };

  // END_READ_GRAPH

  row_map<node> EM;

  ILP_Problem IP(Optsense_Min);

  edge e,f;
  var_map<edge> VM;
  row_map<edge> VM1;
  forall_edges(e, G) {
    VM[e]=IP.add_variable(G[e],0,1, Vartype_Integer);
  };

  graph H;
  node HN[n+1];
  row sr[n];
  row tr[n];
  row_map<edge> VM2;
  for(int i=0; i<=n; i++) HN[i]=H.new_node();
  for(int i=0; i<n; i++) {
    for(int j=0; j<n; j++) if(i!=j) {
      row r;
      forall_in_edges(e, N[i][j]) {
        r+=VM[e];
        if(source(e)==s) sr[i]+=VM[e];
      }
      f=H.new_edge(HN[i],HN[j]);
      //VM1[f]=IP.add_variable(0,0,1, Vartype_Integer);
      VM2[f]=r;
      EM[N[i][j]]=VM2[f];
      //IP.add_basic_constraint(r==VM1[f]);
    }
  }
  forall_in_edges(e,t) { 
    tr[G[source(e)]%n]+=VM[e];
  }
  for(int i=0; i<n; i++) {
    f=H.new_edge(HN[i],HN[n]);
    //VM1[f]=IP.add_variable(0,0,1,Vartype_Integer);
    VM2[f]=tr[i];
    //IP.add_basic_constraint(tr[i]==VM1[f]);
    f=H.new_edge(HN[n],HN[i]);
    //VM1[f]=IP.add_variable(0,0,1,Vartype_Integer);
    VM2[f]=sr[i];
    //IP.add_basic_constraint(sr[i]==VM1[f]);
  };

  cout<<H.number_of_nodes()<<" nodes in H\n";

  IP.add_sym_constraint(new ATOUR(H, VM2));

  IP.add_sym_constraint(new GATOUR(G, n, VM, EM));

  node u,v;
  forall_nodes(u,G) {
    row r;
    if(u!=s) forall_in_edges(e,u) r-=VM[e];
    if(u!=t) forall_out_edges(e,u) r+=VM[e];
    if(u==s) IP.add_basic_constraint(r==1);
    else {
      if(u==t) IP.add_basic_constraint(r==-1);
      else IP.add_basic_constraint(r==0);
    };
  };
/*
  for(int i=0; i<n; i++) {
    row r1;
    for(int j=0; j<n; j++) if(j!=i)
      forall_out_edges(e, N[i][j])
        r1+=VM[e];
    forall_in_edges(e,t) if(G[source(e)]%n==i) r1+=VM[e];
    IP.add_basic_constraint(r1==1);
  }
*/
  IP.optimize();

  int d=0;
  forall_edges(e, G) if(IP.get_solution(VM[e])>=0.5) {
    cout<<"("<<G[source(e)]/n<<","<<G[source(e)]%n<<")-("
             <<G[target(e)]/n<<","<<G[target(e)]%n<<"):"<<G[e]<<" edge\n";
    d+=G[e];
  } else G.del_edge(e);

  u=s;
  do {
    u=G.opposite(u,G.first_adj_edge(u));
    if(u!=t) cout<<G[u]/n<<" "<<Names[G[u]/n]<<endl;
  } while(u!=t);
  cout<<"\n\n";

  node_array<int> NN(H);
  for(int i=0; i<=n; i++) NN[HN[i]]=i;
  forall_edges(e,H) {
    double d=IP.get_solution(VM2[e]);
    //assert(d<=1.001);
    if(d>=0.5) cout<<d<<" "<<NN[source(e)]<<" "<<NN[target(e)]<<" edge\n";
  };

  cerr<<argv[1]<<" "<<d<<endl;

  exit(0);

  return (int) (d+0.5);
};


