3 * This file is a part of LEMON, a generic C++ optimization library
5 * Copyright (C) 2003-2007
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
9 * Permission to use, modify and distribute this software is granted
10 * provided that this copyright notice appears in all copies. For
11 * precise terms see the accompanying LICENSE file.
13 * This software is provided "AS IS" with no warranty of any kind,
14 * express or implied, and with no claim as to its suitability for any
19 #include <lemon/list_graph.h>
20 #include <lemon/graph_utils.h>
21 #include <lemon/random.h>
22 #include <lemon/dim2.h>
23 #include <lemon/bfs.h>
24 #include <lemon/counter.h>
25 #include <lemon/suurballe.h>
26 #include <lemon/graph_to_eps.h>
27 #include <lemon/graph_writer.h>
28 #include <lemon/arg_parser.h>
29 #include <lemon/euler.h>
32 #include <lemon/unionfind.h>
33 #include <lemon/time_measure.h>
35 using namespace lemon;
37 typedef dim2::Point<double> Point;
39 UGRAPH_TYPEDEFS(ListUGraph);
48 std::vector<Node> nodes;
49 ListUGraph::NodeMap<Point> coords(g);
54 for(UEdgeIt e(g);e!=INVALID;++e)
55 tlen+=sqrt((coords[g.source(e)]-coords[g.target(e)]).normSquare());
61 const double EPSILON=1e-8;
62 bool tsp_improve(Node u, Node v)
64 double luv=std::sqrt((coords[v]-coords[u]).normSquare());
69 for(IncEdgeIt e(g,v2);(n=g.runningNode(e))==u2;++e);
72 if(luv+std::sqrt((coords[v2]-coords[u2]).normSquare())-EPSILON>
73 std::sqrt((coords[u]-coords[u2]).normSquare())+
74 std::sqrt((coords[v]-coords[v2]).normSquare()))
76 g.erase(findUEdge(g,u,v));
77 g.erase(findUEdge(g,u2,v2));
87 bool tsp_improve(Node u)
89 for(IncEdgeIt e(g,u);e!=INVALID;++e)
90 if(tsp_improve(u,g.runningNode(e))) return true;
99 for(NodeIt n(g);n!=INVALID;++n)
100 if(tsp_improve(n)) b=true;
106 for(int i=0;i<N;i++) g.addEdge(nodes[i],nodes[(i+1)%N]);
115 Line(Point _a,Point _b) :a(_a),b(_b) {}
116 Line(Node _a,Node _b) : a(coords[_a]),b(coords[_b]) {}
117 Line(const Edge &e) : a(coords[g.source(e)]),b(coords[g.target(e)]) {}
118 Line(const UEdge &e) : a(coords[g.source(e)]),b(coords[g.target(e)]) {}
121 inline std::ostream& operator<<(std::ostream &os, const Line &l)
123 os << l.a << "->" << l.b;
127 bool cross(Line a, Line b)
129 Point ao=rot90(a.b-a.a);
130 Point bo=rot90(b.b-b.a);
131 return (ao*(b.a-a.a))*(ao*(b.b-a.a))<0 &&
132 (bo*(a.a-b.a))*(bo*(a.b-b.a))<0;
142 bool pedgeLess(Pedge a,Pedge b)
147 std::vector<UEdge> edges;
151 Counter cnt("Number of edges added: ");
152 std::vector<Pedge> pedges;
153 for(NodeIt n(g);n!=INVALID;++n)
154 for(NodeIt m=++(NodeIt(n));m!=INVALID;++m)
159 p.len=(coords[m]-coords[n]).normSquare();
162 std::sort(pedges.begin(),pedges.end(),pedgeLess);
163 for(std::vector<Pedge>::iterator pi=pedges.begin();pi!=pedges.end();++pi)
165 Line li(pi->a,pi->b);
167 for(;e!=INVALID && !cross(e,li);++e) ;
170 ne=g.addEdge(pi->a,pi->b);
179 Counter cnt("Number of edges removed: ");
180 Bfs<ListUGraph> bfs(g);
181 for(std::vector<UEdge>::reverse_iterator ei=edges.rbegin();
182 ei!=edges.rend();++ei)
184 Node a=g.source(*ei);
185 Node b=g.target(*ei);
188 if(bfs.predEdge(b)==INVALID || bfs.dist(b)>d)
196 Counter cnt("Number of edges removed: ");
197 for(std::vector<UEdge>::reverse_iterator ei=edges.rbegin();
198 ei!=edges.rend();++ei)
200 Node a=g.source(*ei);
201 Node b=g.target(*ei);
203 ConstMap<Edge,int> cegy(1);
204 Suurballe<ListUGraph,ConstMap<Edge,int> > sur(g,cegy,a,b);
206 if(k<2 || sur.totalLength()>d)
209 // else std::cout << "Remove edge " << g.id(a) << "-" << g.id(b) << '\n';
213 void sparseTriangle(int d)
215 Counter cnt("Number of edges added: ");
216 std::vector<Pedge> pedges;
217 for(NodeIt n(g);n!=INVALID;++n)
218 for(NodeIt m=++(NodeIt(n));m!=INVALID;++m)
223 p.len=(coords[m]-coords[n]).normSquare();
226 std::sort(pedges.begin(),pedges.end(),pedgeLess);
227 for(std::vector<Pedge>::iterator pi=pedges.begin();pi!=pedges.end();++pi)
229 Line li(pi->a,pi->b);
231 for(;e!=INVALID && !cross(e,li);++e) ;
234 ConstMap<Edge,int> cegy(1);
235 Suurballe<ListUGraph,ConstMap<Edge,int> >
236 sur(g,cegy,pi->a,pi->b);
238 if(k<2 || sur.totalLength()>d)
240 ne=g.addEdge(pi->a,pi->b);
251 std::vector<Pedge> pedges;
253 std::cout << T.realTime() << "s: Setting up the edges...\n";
254 for(NodeIt n(g);n!=INVALID;++n)
256 for(NodeIt m=++(NodeIt(n));m!=INVALID;++m)
261 p.len=(coords[m]-coords[n]).normSquare();
264 if(progress && en>=pr*double(N)/100)
266 std::cout << pr << "% \r" << std::flush;
270 std::cout << T.realTime() << "s: Sorting the edges...\n";
271 std::sort(pedges.begin(),pedges.end(),pedgeLess);
272 std::cout << T.realTime() << "s: Creating the tree...\n";
273 ListUGraph::NodeMap<int> comp(g);
274 UnionFind<ListUGraph::NodeMap<int> > uf(comp);
275 for (NodeIt it(g); it != INVALID; ++it)
277 for(std::vector<Pedge>::iterator pi=pedges.begin();pi!=pedges.end();++pi)
279 if ( uf.join(pi->a,pi->b) ) {
280 g.addEdge(pi->a,pi->b);
286 std::cout << T.realTime() << "s: Done\n";
289 Node common(UEdge e, UEdge f)
291 return (g.source(e)==g.source(f)||g.source(e)==g.target(f))?
292 g.source(e):g.target(e);
297 std::cout << "Find a tree..." << std::endl;
301 std::cout << "Total edge length (tree) : " << totalLen() << std::endl;
303 std::cout << "Make it Euler..." << std::endl;
306 std::vector<Node> leafs;
307 for(NodeIt n(g);n!=INVALID;++n)
308 if(countIncEdges(g,n)%2==1) leafs.push_back(n);
309 for(unsigned int i=0;i<leafs.size();i+=2)
310 g.addEdge(leafs[i],leafs[i+1]);
313 for(NodeIt n(g);n!=INVALID;++n)
314 if(countIncEdges(g,n)%2)
315 std::cout << "GEBASZ!!!" << std::endl;
317 std::cout << "Number of edges : " << countUEdges(g) << std::endl;
319 // for(NodeIt n(g);n!=INVALID;++n)
320 // if(countIncEdges(g,n)>2)
322 // std::cout << std::endl;
324 std::cout << "Total edge length (euler) : " << totalLen() << std::endl;
326 ListUGraph::UEdgeMap<UEdge> enext(g);
328 UEulerIt<ListUGraph> e(g);
331 // std::cout << "Tour edge: " << g.id(UEdge(e)) << std::endl;
332 for(++e;e!=INVALID;++e)
334 // std::cout << "Tour edge: " << g.id(UEdge(e)) << std::endl;
341 std::cout << "Creating a tour from that..." << std::endl;
343 int nnum = countNodes(g);
344 int ednum = countUEdges(g);
346 for(UEdge p=UEdgeIt(g);ednum>nnum;p=enext[p])
348 // std::cout << "Checking edge " << g.id(p) << std::endl;
352 Node n1=g.oppositeNode(n2,e);
353 Node n3=g.oppositeNode(n2,f);
354 if(countIncEdges(g,n2)>2)
356 // std::cout << "Remove an Edge" << std::endl;
360 UEdge ne=g.addEdge(n1,n3);
367 std::cout << "Total edge length (tour) : " << totalLen() << std::endl;
371 std::cout << "Total edge length (2-opt tour) : " << totalLen() << std::endl;
375 int main(int argc,const char **argv)
377 ArgParser ap(argc,argv);
380 bool disc_d, square_d, gauss_d;
381 // bool tsp_a,two_a,tree_a;
386 std::string ndist("disc");
387 ap.refOption("n", "Number of nodes (default is 100)", N)
388 .intOption("g", "Girth parameter (default is 10)", 10)
389 .refOption("cities", "Number of cities (default is 1)", num_of_cities)
390 .refOption("area", "Full relative area of the cities (default is 1)", area)
391 .refOption("disc", "Nodes are evenly distributed on a unit disc (default)",disc_d)
392 .optionGroup("dist", "disc")
393 .refOption("square", "Nodes are evenly distributed on a unit square", square_d)
394 .optionGroup("dist", "square")
396 "Nodes are located according to a two-dim gauss distribution",
398 .optionGroup("dist", "gauss")
399 // .mandatoryGroup("dist")
400 .onlyOneGroup("dist")
401 .boolOption("eps", "Also generate .eps output (prefix.eps)")
402 .boolOption("dir", "Directed graph is generated (each edges are replaced by two directed ones)")
403 .boolOption("2con", "Create a two connected planar graph")
404 .optionGroup("alg","2con")
405 .boolOption("tree", "Create a min. cost spanning tree")
406 .optionGroup("alg","tree")
407 .boolOption("tsp", "Create a TSP tour")
408 .optionGroup("alg","tsp")
409 .boolOption("tsp2", "Create a TSP tour (tree based)")
410 .optionGroup("alg","tsp2")
412 .other("[prefix]","Prefix of the output files. Default is 'lgf-gen-out'")
416 switch(ap.files().size())
419 prefix="lgf-gen-out";
422 prefix=ap.files()[0];
425 std::cerr << "\nAt most one prefix can be given\n\n";
430 std::vector<double> sizes;
431 std::vector<double> cum_sizes;
432 for(int s=0;s<num_of_cities;s++)
434 // sum_sizes+=rnd.exponential();
438 cum_sizes.push_back(sum_sizes);
441 for(int s=0;s<num_of_cities;s++)
443 Point center=(num_of_cities==1?Point(0,0):rnd.disc());
445 for(;i<N*(cum_sizes[s]/sum_sizes);i++) {
448 coords[n]=center+rnd.gauss2()*area*
449 std::sqrt(sizes[s]/sum_sizes);
452 for(;i<N*(cum_sizes[s]/sum_sizes);i++) {
455 coords[n]=center+Point(rnd()*2-1,rnd()*2-1)*area*
456 std::sqrt(sizes[s]/sum_sizes);
458 else if(disc_d || true)
459 for(;i<N*(cum_sizes[s]/sum_sizes);i++) {
462 coords[n]=center+rnd.disc()*area*
463 std::sqrt(sizes[s]/sum_sizes);
469 std::cout << "#2-opt improvements: " << tsp_impr_num << std::endl;
473 std::cout << "#2-opt improvements: " << tsp_impr_num << std::endl;
475 else if(ap["2con"]) {
476 std::cout << "Make triangles\n";
478 sparseTriangle(ap["g"]);
479 std::cout << "Make it sparser\n";
482 else if(ap["tree"]) {
487 std::cout << "Number of nodes : " << countNodes(g) << std::endl;
488 std::cout << "Number of edges : " << countUEdges(g) << std::endl;
490 for(UEdgeIt e(g);e!=INVALID;++e)
491 tlen+=sqrt((coords[g.source(e)]-coords[g.target(e)]).normSquare());
492 std::cout << "Total edge length : " << tlen << std::endl;
494 graphToEps(g,prefix+".eps").
495 scale(600).nodeScale(.2).edgeWidthScale(.001).preScale(false).
496 coords(coords).run();
499 GraphWriter<ListUGraph>(prefix+".lgf",g).
500 writeNodeMap("coordinates_x",scaleMap(xMap(coords),600)).
501 writeNodeMap("coordinates_y",scaleMap(yMap(coords),600)).
503 else UGraphWriter<ListUGraph>(prefix+".lgf",g).
504 writeNodeMap("coordinates_x",scaleMap(xMap(coords),600)).
505 writeNodeMap("coordinates_y",scaleMap(yMap(coords),600)).