/*******************************************************************************
+
+  LEDA 3.5
+
+  _ortho.c
+
+  This file is part of the LEDA research version (LEDA-R) that can be 
+  used free of charge in academic research and teaching. Any commercial
+  use of this software requires a license which is distributed by the
+  LEDA Software GmbH, Postfach 151101, 66041 Saarbruecken, FRG
+  (fax +49 681 31104).
+
+  Copyright (c) 1991-1997  by  Max-Planck-Institut fuer Informatik
+  Im Stadtwald, 66123 Saarbruecken, Germany     
+  All rights reserved.
+ 
*******************************************************************************/
#define EPS ""
#define INFINITY MAXINT
#define REV(e) P.reversal(e)
#define SUCC(e) P.face_cycle_succ(e)
#define PRED(e) P.face_cycle_pred(e)
#define IS_CAGE(f)  (P.get_cage(P.first_face_edge(f)))
#define IS_OUTER(f) (!P.get_inner(P.first_face_edge(f)))
#define NEXT(d) ((direction)((d+5)%4))
#define PREV(d) ((direction)((d+3)%4))
#define OPP(d)  ((direction)((d+6)%4))

#define ERR_EMPTY_GRAPH         "ORTHO: input graph is empty"
#define ERR_NOT_CONNECTED       "ORTHO: input graph is not connected"
#define ERR_NO_PLANAR_MAP       "ORTHO: this is no planar map"
#define ERR_INVALID_NETWORK     "ORTHO: invalid network"
#define ERR_NO_FEASIBLE_FLOW    "ORTHO: no feasible flow"
#define ERR_NO_ORTHO_REP        "ORTHO: orthogonal rep. not valid"
#define ERR_BAD_STRING          "ORTHO: trying to set nonempty string"
#define ERR_BAD_ANGLE           "ORTHO: bad angle"
#define ERR_BAD_DIRECTION       "ORTHO: bad direction"
#define ERR_OPEN_CPLEX_ENV      "CPLEX: failed to open environment"
#define ERR_LOAD_LP             "CPLEX: failed to load LP"
#define ERR_CPLEX_ADD_ROWS      "CPLEX: CPXaddrows failed"
#define ERR_CPLEX_WRITE         "CPLEX: CPXwriteLP failed"
#define ERR_CPLEX_OPT           "CPLEX: failed to optimize LP"
#define ERR_CPLEX_GET_SOLUTION  "CPLEX: failed to obtain solution"
#define ERR_CPLEX_FREE          "CPLEX: CPXfreeprob failed"
#define ERR_CPLEX_CLOSE         "CPLEX: CPXcloseCPLEX failed"

enum direction{north,east,south,west,unexplored};
enum v_type{real,bend,dissection,big};

#include <LEDA/planar_map.h>
#include <LEDA/node_map.h>
#include <LEDA/edge_map.h>
#include <LEDA/face_map.h>
#include <LEDA/stack.h>
#include <LEDA/graph_alg.h>


void longest_paths(const GRAPH<int,int>&G,node_array<int>&l){
node_array<int>INDEG(G,0);node v;edge e;
stack<node>S;
forall_nodes(v,G){
INDEG[v]= indeg(v);
if(INDEG[v]==0)S.push(v);
}
while(!S.empty()){
v= S.pop();
forall_out_edges(e,v){
node w= target(e);
l[w]= Max(l[w],l[v]+G.inf(e));
if(--INDEG[w]==0)S.push(w);
}
}
}


int angle(const char c){
int result= 0;
switch(c){
case'0':result= 90;break;
case'1':result= 270;break;
default:error_handler(1,ERR_BAD_ANGLE);
}
return result;
}




bool Euler(const graph&G){
int n= 0;node v;forall_nodes(v,G)if(outdeg(v)!=0)++n;
int m= G.number_of_edges()/2;
int f= G.number_of_faces();
return(m==n+f-2);
}




class ortho_map:public planar_map{
node_map<node>v_in_G;
node_map<v_type>v_T;
edge_map<edge>e_in_G;
edge_map<edge>e_in_P;
edge_map<int>a;
edge_map<string>s;
edge_map<bool>inner;
edge_map<bool>cage;
edge_array<direction>dir;
edge_array<int>length;
node_array<int>x,y;

public:
edge next_level_edge(const node,const direction);
ortho_map(const graph&);
void print();
//bool check();
void init_maps(const edge,const int,const string,const bool,const bool);
node get_orig(const node v)const{return v_in_G[v];}
edge get_orig(const edge e)const{return e_in_G[e];}
direction get_dir(const edge e)const{return dir[e];}
edge get_copy(const edge e)const{return e_in_P[e];}
v_type get_type(const node v)const{return v_T[v];}
void set_inner(const edge e,const bool t){inner[e]= t;}
void set_cage(const edge e,const bool t){cage[e]= t;}
bool get_inner(const edge e)const{return inner[e];}
bool get_cage(const edge e)const{return cage[e];}
void set_type(const node v,const v_type t){v_T[v]= t;}
int get_a(const edge e)const{return a[e];}
void set_a(const edge e,const int angle){a[e]= angle;}
string get_s(const edge e)const{return s[e];}
void set_s(const edge e,const string s_new){s[e]= s_new;}
void set_s(const edge e,const int,const int,const bool);
int get_x(const node v)const{return x[v];}
int get_y(const node v)const{return y[v];}
void set_rev_s(const edge e,const int,const int,const bool);
void set_length(const edge e,const int l_new){length[e]= l_new;}
int get_length(const edge e)const{return length[e];}
edge split_map_edge(edge);
edge new_edge(edge,edge);
edge split_bend_edge(edge);
void PrintEdge(const edge);
void init_rest(){dir.init(*this,unexplored);
length.init(*this,0);
x.init(*this,-1);y.init(*this,-1);}



void ortho_map::assign_directions(edge e,direction d){
while(dir[e]==unexplored){
dir[e]= d;edge r= reversal(e);
if(dir[r]==unexplored){
if(inner[r])assign_directions(r,OPP(d));
else dir[r]= OPP(d);
}
switch(a[e]){
case 90:d= NEXT(d);break;
case 270:d= PREV(d);break;
case 360:d= OPP(d);break;
}
e= face_cycle_succ(e);
}
}

#line 126 "ortho.w"

/*:30*/
#line 58 "ortho_map.w"
;
void determine_position(const node,const int x,const int y,node_array<bool>&);
void norm_positions();
edge succ_corner_edge(const edge);
~ortho_map(){clear();}
LEDA_MEMORY(ortho_map)
};




edge ortho_map::next_level_edge(const node v,const direction d){
edge e;
switch(d){
case east:forall_out_edges(e,v)if(dir[e]==north)return e;break;
case west:forall_in_edges(e,v)if(dir[e]==north)return e;break;
case north:forall_in_edges(e,v)if(dir[e]==east)return e;break;
case south:forall_out_edges(e,v)if(dir[e]==east)return e;break;
default:break;
}
return NULL;
}




ortho_map::ortho_map(const graph&G):planar_map(G){
if(G.number_of_nodes()==0)error_handler(1,ERR_EMPTY_GRAPH);
if(!Is_Connected(G))error_handler(1,ERR_NOT_CONNECTED);
if(!Euler(G))error_handler(1,ERR_NO_PLANAR_MAP);
v_in_G.init(*this);e_in_G.init(*this);
e_in_P.init(G);


node v_P= first_node(),v_G= G.first_node();
while(v_P){
v_in_G[v_P]= v_G;
v_P= succ_node(v_P);
v_G= G.succ_node(v_G);
}
forall_nodes(v_P,*this){
v_G= v_in_G[v_P];
edge e_P= first_adj_edge(v_P),e_G= G.first_adj_edge(v_G);
while(e_P){
e_in_G[e_P]= e_G;e_in_P[e_G]= e_P;
e_G= adj_succ(e_G);e_P= adj_succ(e_P);
}
}


;
a.init(*this,90);s.init(*this,EPS);v_T.init(*this,real);
inner.init(*this,true);cage.init(*this,false);
}





edge ortho_map::split_map_edge(edge e){
edge n= planar_map::split_edge(e);
edge er= reversal(e),nr= reversal(n);
a[er]= a[nr];a[n]= a[e];
s[er]= EPS;s[n]= EPS;
inner[er]= inner[nr];cage[er]= cage[nr];
inner[n]= inner[e];cage[n]= cage[e];
e_in_G[n]= e_in_G[e];
e_in_G[er]= e_in_G[nr];
v_in_G[source(n)]= NULL;

return n;
}

edge ortho_map::split_bend_edge(edge e){
string s_e= s[e];int a_e= a[e];
bool inner_e= inner[e];bool cage_e= cage[e];
edge er= reversal(e);
string s_er= s[er];int a_er= a[er];
bool inner_er= inner[er];bool cage_er= cage[er];

edge n= split_map_edge(e);
er= reversal(e);edge nr= reversal(n);

a[e]= angle(s_e[s_e.length()-1]);a[er]= a_er;
a[n]= a_e;a[nr]= angle(s_er[0]);
s[e]= s_e.head(s_e.length()-1);s[er]= s_er.tail(s_er.length()-1);
s[n]= EPS;s[nr]= EPS;
inner[e]= inner[n]= inner_e;
cage[e]= cage[n]= cage_e;
inner[er]= inner[nr]= inner_er;
cage[er]= cage[nr]= cage_er;
v_T[source(n)]= bend;

return n;
}

edge ortho_map::new_edge(edge e1,edge e2){
edge n= planar_map::new_edge(e1,e2);
e_in_G[n]= NULL;e_in_G[reversal(n)]= NULL;
return n;
}


void ortho_map::print(){
face f;edge e;
forall_faces(f,(*this)){
if(!inner[first_face_edge(f)])cout<<"outer ";
if(cage[first_face_edge(f)])cout<<"cage ";
cout<<"face: "<<endl;
forall_face_edges(e,f){PrintEdge(e);newline;}
}
}

/*
bool ortho_map::check(){
bool result= true;return result;
face f;edge e;node v;
forall_faces(f,(*this)){
int rotation= 0;
forall_face_edges(e,f)rotation+= zeroes(s[e])-ones(s[e])+2-a[e]/90;
if(!inner[first_face_edge(f)])result= result&&(rotation==-4);
else result= result&&(rotation==4);
}
forall_nodes(v,(*this))if(outdeg(v)!=0){
int angle_sum= 0;
forall_in_edges(e,v)angle_sum+= a[e];
result= result&&(angle_sum==360);
}
return result;
}
*/

void ortho_map::set_s(const edge e,const int flow,const int flow_rev,
const bool bridge){
if(s[e]!=EPS)error_handler(1,ERR_BAD_STRING);
for(int i= 0;i<flow;i++)s[e]+= "0";
if(!bridge){
 for(int i= 0;i<flow_rev;i++)s[e]+= "1";
}
}

void ortho_map::set_rev_s(const edge e,const int flow,const int flow_rev,
const bool bridge){
if(s[e]!=EPS)error_handler(1,ERR_BAD_STRING);
for(int i= 0;i<flow;i++)s[e]+= "1";
if(!bridge){
for(int i= 0;i<flow_rev;i++)s[e]= string('0')+s[e];
}
}

void ortho_map::PrintEdge(const edge e){
planar_map::print_edge(e);
cout<<"\t("<<s[e]<<" "<<a[e]<<" ["<<inner[e]<<cage[e]<<"])";
}

void ortho_map::init_maps(const edge e,const int a_new,const string s_new,
const bool i_new,const bool c_new){
a[e]= a_new;s[e]= s_new;
inner[e]= i_new;cage[e]= c_new;
}


void ortho_map::determine_position(const node v,const int x_new,
const int y_new,node_array<bool>&seen){
if(seen[v])return;
edge e;x[v]= x_new;y[v]= y_new;
seen[v]= true;
forall_out_edges(e,v){
if(dir[e]==north)
determine_position(target(e),x_new+length[e],y_new,seen);
if(dir[e]==west)
determine_position(target(e),x_new,y_new+length[e],seen);
}
forall_in_edges(e,v){
if(dir[e]==north)
determine_position(source(e),x_new-length[e],y_new,seen);
if(dir[e]==west)
determine_position(source(e),x_new,y_new-length[e],seen);
}
}

void ortho_map::norm_positions(){
int xmin= 0,ymin= 0;node v;
forall_nodes(v,(*this))if(v_T[v]!=big){
xmin= Min(xmin,x[v]);
ymin= Min(ymin,y[v]);
}
forall_nodes(v,(*this)){
x[v]-= xmin;
y[v]-= ymin;
}
}

edge ortho_map::succ_corner_edge(const edge e){
edge e_c;
for(e_c= face_cycle_succ(e);a[e_c]==180;e_c= face_cycle_succ(e_c));
return e_c;
}









/*
bool Euler(const graph&);
int angle(const char);
int zeroes(const string&);
int ones(const string&);
void longest_paths(const GRAPH<int,int>&,node_array<int>&);
*/


/*
#include "common.h"
#include "ortho_map.h"
*/

#include <LEDA/stack.h>
#include <LEDA/array.h>
#include <LEDA/set.h>
#include <LEDA/integer_matrix.h>
#include <LEDA/graph_alg.h>

#ifdef CPLEX
extern"C"
#include<cplex.h>
#endif




typedef list<int>  intlist;
typedef list<node> nodelist;

int ORTHO_EMBEDDING(const graph&G,
                    node_array<int>&x_pos,
                    node_array<int>&y_pos,
                    edge_array<intlist>&x_bends,
                    edge_array<intlist>&y_bends, bool  )
{
ortho_map P(G);


node_map<edge>corr_cage_edge(P);
edge_array<nodelist>b_nodes(G);
edge_array<node>b_nodes_first(G,NULL);
edge_array<node>b_nodes_last(G,NULL);
list<node>all_nodes= P.all_nodes();
node v;edge e;face f;
forall(v,all_nodes)if(outdeg(v)>4){
P.set_type(v,big);
int d= outdeg(v),i= 0;
array<edge>out(d);edge e_cage;


forall_out_edges(e,v){
out[i++]= e;edge e_orig= P.get_orig(e);
edge e_split= P.split_map_edge(e);
node c_i= source(e_split);
P.set_type(c_i,bend);
b_nodes_first[e_orig]= c_i;
b_nodes_last[G.reversal(e_orig)]= c_i;
}





for(i= 0;i<d;i++){
e_cage= P.new_edge(SUCC(out[i]),REV(out[(i+1)%d]));
P.init_maps(e_cage,90,EPS,true,true);
edge r_cage= P.reversal(e_cage);
P.init_maps(r_cage,90,EPS,true,false);
}



corr_cage_edge[v]= e_cage;
forall_out_edges(e,v)P.join_faces(e);



}
//cout<<"cages created"<<endl;





//if(check)if(!Euler(P))error_handler(1,ERR_NO_PLANAR_MAP);


P.compute_faces();
face f_0;int f_0_deg= 0;
forall_faces(f,P)
if(!IS_CAGE(f)&&P.size(f)>f_0_deg){
f_0= f;f_0_deg= P.size(f_0);
}
forall_face_edges(e,f_0)P.set_inner(e,false);



;


graph N;
node s,t;
edge_array<int>cap,cost,l;
int z= 0;



list<node>V_hat;
list<node>F;
node_map<node>NtoV(N);
node_map<face>NtoF(N);
face_map<node>FtoN(P);


;


s= N.new_node();
t= N.new_node();
forall_nodes(v,P)if(outdeg(v)>0&&outdeg(v)<=3){
node n= N.new_node();
NtoV[n]= v;
V_hat.append(n);
}
forall_faces(f,P){
node n= N.new_node();
NtoF[n]= f;FtoN[f]= n;
F.append(n);
}


;


list<edge>A_s_v,A_s_f,A_v,A_v_cage,A_f,A_f_t;


forall_faces(f,P){
int size= P.size(f);
if(!IS_OUTER(f)&&size<=3)
A_s_f.append(N.new_edge(s,FtoN[f]));
if(IS_OUTER(f)||size>=5)A_f_t.append(N.new_edge(FtoN[f],t));
}
forall(v,V_hat)A_s_v.append(N.new_edge(s,v));


;


forall(v,V_hat){
set<face>F_adj;
forall_adj_faces(f,NtoV[v])F_adj.insert(f);
forall(f,F_adj)
if(IS_CAGE(f))A_v_cage.append(N.new_edge(v,FtoN[f]));
else A_v.append(N.new_edge(v,FtoN[f]));
}


;


#define MAX_BENDS_PER_EDGE INFINITY
edge_array<bool>marked(P,false);
edge_map<edge>partner(N),a_in_P(N);edge a1,a2;
forall_faces(f,P)
forall_face_edges(e,f)if(!marked[e]){
face g= P.face_of(REV(e));
if(f==g){
a1= N.new_edge(FtoN[f],FtoN[f]);
A_f.append(a1);
partner[a1]= a1;a_in_P[a1]= e;
}else{
a1= N.new_edge(FtoN[f],FtoN[g]);
a2= N.new_edge(FtoN[g],FtoN[f]);
A_f.append(a1);A_f.append(a2);
partner[a1]= a2;partner[a2]= a1;
a_in_P[a1]= e;a_in_P[a2]= e;
edge e_loop;
for(e_loop= e;
P.face_of(REV(e_loop))==g&&!marked[e_loop];
e_loop= SUCC(e_loop))
marked[e_loop]= marked[REV(e_loop)]= true;
for(e_loop= PRED(e);
P.face_of(REV(e_loop))==g&&!marked[e_loop];
e_loop= PRED(e_loop))
marked[e_loop]= marked[REV(e_loop)]= true;
}
}


;


;


cap.init(N,0);cost.init(N,0);int z_s= 0;edge a;
l.init(N,0);
forall(a,A_s_f)cap[a]= 4-P.size(NtoF[target(a)]);
forall(a,A_s_v)cap[a]= 4-outdeg(NtoV[target(a)]);
forall(a,A_v)cap[a]= 3;
forall(a,A_v_cage){cap[a]= 3;l[a]= 1;}
forall(a,A_f){
face g= NtoF[target(a)];face f= NtoF[source(a)];
cap[a]= IS_CAGE(g)?0:MAX_BENDS_PER_EDGE;
cost[a]= IS_CAGE(f)?0:1;
}
forall(a,A_f_t){
face f= NtoF[source(a)];
if(IS_OUTER(f))cap[a]= P.size(f)+4;
else cap[a]= P.size(f)-4;
}

forall_out_edges(a,s)z_s+= cap[a];
forall_in_edges(a,t)z+= cap[a];

A_v.conc(A_v_cage);

if(z!=z_s)error_handler(1,ERR_INVALID_NETWORK);


;
//cout<<"network constructed"<<endl;






edge_array<int>flow(N);
node_array<int>supply(N,0);
supply[s]= z;supply[t]= -z;
bool feasible= MIN_COST_FLOW(N,l,cap,cost,supply,flow);
if(!feasible)error_handler(1,ERR_NO_FEASIBLE_FLOW);
//cout<<"min-cost flow computed"<<endl;







forall(a,A_v){
node v= NtoV[source(a)];face f= NtoF[target(a)];
edge e_in;bool flag= false;
forall_in_edges(e_in,v)
if(P.face_of(e_in)==f&&!flag){
P.set_a(e_in,(flow[a]+1)*90);
flag= true;
}
}






marked.init(N,false);int no_of_bends= 0;
forall(a,A_f)if(!marked[a]){
edge a_rev= partner[a];
marked[a]= true;marked[a_rev]= true;
e= a_in_P[a];
bool bridge= (a==a_rev);
P.set_s(e,flow[a],flow[a_rev],bridge);
P.set_rev_s(REV(e),flow[a],flow[a_rev],bridge);
no_of_bends+= (flow[a]+flow[a_rev]);
}





marked.init(P,false);
list<edge>all_edges= P.all_edges();
forall(e,all_edges)if(!marked[e]){
marked[e]= true;marked[REV(e)]= true;
edge e_orig= P.get_orig(e);
while(P.get_s(e)!=EPS){
edge e_split= P.split_bend_edge(e);
if(e_orig){
b_nodes[e_orig].push(source(e_split));
b_nodes[G.reversal(e_orig)].append(source(e_split));
}
}
if(e_orig){
if(b_nodes_first[e_orig]){
b_nodes[e_orig].push(b_nodes_first[e_orig]);
b_nodes[G.reversal(e_orig)].append(b_nodes_first[e_orig]);
}
if(b_nodes_last[e_orig]){
b_nodes[e_orig].append(b_nodes_last[e_orig]);
b_nodes[G.reversal(e_orig)].push(b_nodes_last[e_orig]);
}
}
}



//cout<<"orthogonal representation constructed"<<endl;
//if(check)if(!P.check())error_handler(1,ERR_NO_ORTHO_REP);
//cout<<"no. of bends in 4-graph: "<<no_of_bends<<endl;






list<face>all_faces;


forall_faces(f,P)if(IS_OUTER(f))all_faces.push(f);
else{
if(IS_CAGE(f)){
int pred_angle= P.get_a(PRED(P.first_face_edge(f)));
forall_face_edges(e,f){
int act_angle= P.get_a(e);
if(pred_angle==90&&act_angle==90){
edge e_split= P.split_map_edge(e);
P.set_a(e,180);P.set_a(REV(e_split),180);
P.set_type(source(e_split),dissection);
}
pred_angle= P.get_a(e);
}
}
all_faces.append(f);
}



forall(f,all_faces){

stack<edge>S;
int size= 0;
forall_face_edges(e,f)size+= abs(2-P.get_a(e)/90);
e= P.succ_corner_edge(P.first_face_edge(f));int state= 0;
while(size>4){
if(state==0&&!S.empty()){
e= P.succ_corner_edge(S.top());
state= 1;
}
int angle= P.get_a(e);
switch(angle){
case 90:if(state==2)state= 3;
if(state==1)state= 2;
break;
case 270:
case 360:if(state==2&&IS_OUTER(f))state= 4;
else{
S.push(e);state= 1;
if(angle==360)S.push(e);
}
break;
}
if(state==3){

edge e1= S.pop();
edge e2= P.succ_corner_edge(e1);
edge e3= P.succ_corner_edge(e2);
edge e4= P.succ_corner_edge(e3);int a4= P.get_a(e4);
edge e5= P.split_map_edge(e4);
P.set_type(source(e5),dissection);
P.set_a(e1,P.get_a(e1)-90);
P.set_a(e5,a4);P.set_a(REV(e5),180);
P.set_a(e4,90);
edge e6= P.new_edge(P.face_cycle_succ(e1),e5);
P.init_maps(e6,90,EPS,P.get_inner(e1),false);
P.init_maps(REV(e6),90,EPS,true,false);
P.set_inner(e2,true);P.set_inner(e3,true);P.set_inner(e4,true);
size-= 2;state= 0;e= e6;



}
if(state==4){

edge e1= S.pop();
edge e2= P.succ_corner_edge(e1);
edge e3= P.succ_corner_edge(e2);
edge e4= P.new_edge(P.face_cycle_succ(e1),P.face_cycle_succ(e3));
P.init_maps(e4,P.get_a(e3)-90,"1",false,false);
P.init_maps(REV(e4),90,"0",true,false);
edge e5= P.split_bend_edge(e4);
P.set_a(e1,P.get_a(e1)-90);
P.set_inner(e2,true);
P.init_maps(e3,90,EPS,true,false);
P.set_type(source(e5),dissection);
S.push(e1);
size-= 2;state= 0;



}
e= P.succ_corner_edge(e);
}


}
//if(check&&!P.check())error_handler(1,ERR_NO_ORTHO_REP);


P.init_rest();
P.assign_directions(P.first_edge(),north);









#ifdef CPLEX





int n= 0;int m= 0;
node_array<int>v_num(P);edge_array<int>e_num(P,-1);
forall_nodes(v,P)v_num[v]= n++;
forall_edges(e,P)if(P.get_dir(e)==north||P.get_dir(e)==west)
e_num[e]= m++;
int basic_rows= 2*m;
int rows= basic_rows;
int cols= 2*n;
int basic_nonzeroes= 2*basic_rows;
int nonzeroes= basic_nonzeroes;
double*obj= new double[cols];
double*rhs= new double[rows];
char*sense= new char[rows];
int*matbeg= new int[cols];
int*matcnt= new int[cols];
int*matind= new int[nonzeroes];
double*matval= new double[nonzeroes];
double*lb= new double[cols];
double*ub= new double[cols];

CPXENVptr env= NULL;
CPXLPptr lp= NULL;

int status;







for(int j= 0;j<cols;j++){
obj[j]= 0.0;
matbeg[j]= 0;matcnt[j]= 0;
lb[j]= 0.0;ub[j]= INFBOUND;
}
for(int k= 0;k<basic_nonzeroes;k++){
matind[k]= 0;matval[k]= 0.0;
}





#define X(v) (v_num[v])
#define Y(v) (n + v_num[v])
forall_edges(e,P){
node s= source(e),t= target(e);
switch(P.get_dir(e)){
case north:
obj[v_num[s]]-= 1.0;
obj[v_num[t]]+= 1.0;
break;
case west:
obj[v_num[s]+n]-= 1.0;
obj[v_num[t]+n]+= 1.0;
break;
default:break;
}
}





int nonzero_cnt= 0;
forall_nodes(v,P){
int act_col= X(v);
matbeg[act_col]= nonzero_cnt;
forall_inout_edges(e,v){
bool out_edge= (v==source(e));int act_row;
switch(P.get_dir(e)){
case north:
act_row= 2*e_num[e];
matcnt[act_col]++;
if(out_edge)matval[nonzero_cnt]= 1.0;
else matval[nonzero_cnt]= -1.0;
matind[nonzero_cnt]= act_row;
sense[act_row]= 'L';rhs[act_row]= -1.0;
++nonzero_cnt;
break;
case west:
act_row= 2*e_num[e]+1;
matcnt[act_col]++;
if(out_edge)matval[nonzero_cnt]= 1.0;
else matval[nonzero_cnt]= -1.0;
matind[nonzero_cnt]= act_row;
sense[act_row]= 'E';rhs[act_row]= 0.0;
++nonzero_cnt;
break;
default:break;
}
}
}
forall_nodes(v,P){
int act_col= Y(v);
matbeg[act_col]= nonzero_cnt;
forall_inout_edges(e,v){
bool out_edge= (v==source(e));int act_row;
switch(P.get_dir(e)){
case west:
act_row= 2*e_num[e];
matcnt[act_col]++;
if(out_edge)matval[nonzero_cnt]= 1.0;
else matval[nonzero_cnt]= -1.0;
matind[nonzero_cnt]= act_row;
sense[act_row]= 'L';rhs[act_row]= -1.0;
++nonzero_cnt;
break;
case north:
act_row= 2*e_num[e]+1;
matcnt[act_col]++;
if(out_edge)matval[nonzero_cnt]= 1.0;
else matval[nonzero_cnt]= -1.0;
matind[nonzero_cnt]= act_row;
sense[act_row]= 'E';rhs[act_row]= 0.0;
++nonzero_cnt;
break;
default:break;
}
}
}
env= CPXopenCPLEX(&status);
if(env==NULL)error_handler(1,ERR_OPEN_CPLEX_ENV);

lp= CPXloadlp(env,"find_coords",cols,basic_rows,CPX_MIN,obj,rhs,
sense,matbeg,matcnt,matind,
matval,lb,ub,NULL,cols,rows,nonzeroes);
if(lp==NULL)error_handler(1,ERR_LOAD_LP);





int solstat;
double objval;
double*sol= new double[cols];
double*pi= new double[rows];
double*slack= new double[rows];
double*dj= new double[cols];

status= CPXoptimize(env,lp);
if(status){
CPXcloseCPLEX(&env);
error_handler(1,ERR_CPLEX_OPT);
}

status= CPXsolution(env,lp,&solstat,&objval,sol,pi,slack,dj);
if(status)error_handler(1,ERR_CPLEX_GET_SOLUTION);
int max_x= 0,max_y= 0;
forall_nodes(v,P)if(P.get_type(v)==real){
max_x= Max(max_x,(x_pos[P.get_orig(v)]= (int)sol[X(v)]));
max_y= Max(max_y,(y_pos[P.get_orig(v)]= (int)sol[Y(v)]));
}

if(env)status= CPXfreeprob(env,&lp);
if(status)error_handler(1,ERR_CPLEX_FREE);
status= CPXcloseCPLEX(&env);
if(status)error_handler(1,ERR_CPLEX_CLOSE);







#else




graph N_h,N_v;
face_map<node>FtoN_h(P),FtoN_v(P);
edge_map<edge>EtoA(P);


node s_h= N_h.new_node(),t_h= N_h.new_node(),s_v= N_v.new_node(),
t_v= N_v.new_node();

forall_faces(f,P){
if(!IS_OUTER(f)){
FtoN_h[f]= N_h.new_node();
FtoN_v[f]= N_v.new_node();
}
}





forall_faces(f,P)forall_face_edges(e,f){
face g= P.face_of(REV(e));
bool out_f= IS_OUTER(f),out_g= IS_OUTER(g);
switch(P.get_dir(e)){
case north:
if(out_f&&!out_g)EtoA[e]= N_h.new_edge(FtoN_h[g],t_h);
if(out_g&&!out_f)EtoA[e]= N_h.new_edge(s_h,FtoN_h[f]);
if(out_f&&out_g)EtoA[e]= N_h.new_edge(s_h,t_h);
if(!out_f&&!out_g)EtoA[e]= 
N_h.new_edge(FtoN_h[g],FtoN_h[f]);
break;
case west:
if(out_f&&!out_g)EtoA[e]= N_v.new_edge(FtoN_v[g],t_v);
if(out_g&&!out_f)EtoA[e]= N_v.new_edge(s_v,FtoN_v[f]);
if(out_f&&out_g)EtoA[e]= N_v.new_edge(s_v,t_v);
if(!out_f&&!out_g)EtoA[e]= 
N_v.new_edge(FtoN_v[g],FtoN_v[f]);
break;
default:break;
}
}
edge e_h= N_h.new_edge(s_h,t_h);
edge e_v= N_v.new_edge(s_v,t_v);




edge_array<int>cap_h(N_h,INFINITY),cap_v(N_v,INFINITY);
edge_array<int>cost_h(N_h,1),cost_v(N_v,1);
edge_array<int>l_h(N_h,1),l_v(N_v,1);
edge_array<int>flow_h(N_h,0),flow_v(N_v,0);
node_array<int>b_h(N_h,0),b_v(N_v,0);
forall_nodes(v,N_h)
b_h[s_h]+= abs(outdeg(v)-indeg(v));
b_h[s_h]= b_h[s_h]/2+1;
b_h[t_h]= -b_h[s_h];
forall_nodes(v,N_v)
b_v[s_v]+= abs(outdeg(v)-indeg(v));
b_v[s_v]= b_v[s_v]/2+1;
b_v[t_v]= -b_v[s_v];
cost_h[e_h]= 0;cost_v[e_v]= 0;
l_h[e_h]= 0;l_v[e_v]= 0;





feasible= MIN_COST_FLOW(N_h,l_h,cap_h,cost_h,b_h,flow_h);
if(!feasible)error_handler(1,ERR_NO_FEASIBLE_FLOW);
feasible= MIN_COST_FLOW(N_v,l_v,cap_v,cost_v,b_v,flow_v);
if(!feasible)error_handler(1,ERR_NO_FEASIBLE_FLOW);





forall_edges(e,P){
switch(P.get_dir(e)){
case north:
P.set_length(e,flow_h[EtoA[e]]);
P.set_length(REV(e),flow_h[EtoA[e]]);
break;
case west:
P.set_length(e,flow_v[EtoA[e]]);
P.set_length(REV(e),flow_v[EtoA[e]]);
break;
default:break;
}
}





for(v= P.choose_node();P.get_type(v)!=real;v= P.succ_node(v));
x_pos.init(G);y_pos.init(G);node_array<bool>seen(P,false);
P.determine_position(v,0,0,seen);
P.norm_positions();
int max_x= 0,max_y= 0;
forall_nodes(v,P)if(P.get_type(v)==real){
int x = P.get_x(v);
int y = P.get_y(v);
x_pos[P.get_orig(v)]=x;
y_pos[P.get_orig(v)]=y;
max_x= Max(max_x,x);
max_y= Max(max_y,y);
}



//cout<<"coordinates computed by network method"<<endl;



#endif



forall_nodes(v,P)if(P.get_type(v)==big){
int x_big,y_big;
list<int>lx,ly;
face f_cage= P.face_of(corr_cage_edge[v]);
forall_face_edges(e,f_cage){
if(P.get_type(source(e))!=dissection){
int x= P.get_x(source(e)),y= P.get_y(source(e));
lx.append(x);ly.append(y);
}
}
lx.sort();ly.sort();list_item it;
x_big= (lx.head()+lx.tail())/2;
y_big= (ly.head()+ly.tail())/2;
forall_items(it,lx)if(lx[it]!=lx.head()&&lx[it]!=lx.tail()
&&lx[it]==lx[lx.cyclic_succ(it)])x_big= lx[it];
forall_items(it,ly)if(ly[it]!=ly.head()&&ly[it]!=ly.tail()
&&ly[it]==ly[ly.cyclic_succ(it)])y_big= ly[it];
x_pos[P.get_orig(v)]= x_big;
y_pos[P.get_orig(v)]= y_big;
}

x_bends.init(G);y_bends.init(G);
forall_edges(e,G)
forall(v,b_nodes[e]){
x_bends[e].append(P.get_x(v));
y_bends[e].append(P.get_y(v));
}

//cout<<"cages resolved"<<endl<<endl;








return Max(max_x,max_y);
}


