COIN-OR::LEMON - Graph Library

source: lemon-0.x/tools/lgf-gen.cc @ 2402:da8eb8f4ea41

Last change on this file since 2402:da8eb8f4ea41 was 2402:da8eb8f4ea41, checked in by Alpar Juttner, 17 years ago

An improved version of ArgParser?: You don't need to give an explicit storage
for each option.
TODO: Documentation must be updated

File size: 9.6 KB
Line 
1/* -*- C++ -*-
2 *
3 * This file is a part of LEMON, a generic C++ optimization library
4 *
5 * Copyright (C) 2003-2007
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
8 *
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.
12 *
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
15 * purpose.
16 *
17 */
18
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 <cmath>
30#include <algorithm>
31#include <lemon/unionfind.h>
32#include <lemon/time_measure.h>
33
34using namespace lemon;
35
36typedef dim2::Point<double> Point;
37
38UGRAPH_TYPEDEFS(ListUGraph);
39
40bool progress=true;
41
42int N;
43// int girth;
44
45ListUGraph g;
46
47std::vector<Node> nodes;
48ListUGraph::NodeMap<Point> coords(g);
49
50int tsp_impr_num=0;
51
52const double EPSILON=1e-8;
53bool tsp_improve(Node u, Node v)
54{
55  double luv=std::sqrt((coords[v]-coords[u]).normSquare());
56  Node u2=u;
57  Node v2=v;
58  do {
59    Node n;
60    for(IncEdgeIt e(g,v2);(n=g.runningNode(e))==u2;++e);
61    u2=v2;
62    v2=n;
63    if(luv+std::sqrt((coords[v2]-coords[u2]).normSquare())-EPSILON>
64       std::sqrt((coords[u]-coords[u2]).normSquare())+
65       std::sqrt((coords[v]-coords[v2]).normSquare()))
66      {
67        g.erase(findUEdge(g,u,v));
68        g.erase(findUEdge(g,u2,v2));
69        g.addEdge(u2,u);
70        g.addEdge(v,v2);
71        tsp_impr_num++;
72        return true;
73      }
74  } while(v2!=u);
75  return false;
76}
77
78bool tsp_improve(Node u)
79{
80  for(IncEdgeIt e(g,u);e!=INVALID;++e)
81    if(tsp_improve(u,g.runningNode(e))) return true;
82  return false;
83}
84
85void tsp_improve()
86{
87  bool b;
88  do {
89    b=false;
90    for(NodeIt n(g);n!=INVALID;++n)
91      if(tsp_improve(n)) b=true;
92  } while(b);
93}
94
95void tsp()
96{
97  for(int i=0;i<N;i++) g.addEdge(nodes[i],nodes[(i+1)%N]);
98  tsp_improve();
99}
100
101class Line
102{
103public:
104  Point a;
105  Point b;
106  Line(Point _a,Point _b) :a(_a),b(_b) {}
107  Line(Node _a,Node _b) : a(coords[_a]),b(coords[_b]) {}
108  Line(const Edge &e) : a(coords[g.source(e)]),b(coords[g.target(e)]) {}
109  Line(const UEdge &e) : a(coords[g.source(e)]),b(coords[g.target(e)]) {}
110};
111 
112inline std::ostream& operator<<(std::ostream &os, const Line &l)
113{
114  os << l.a << "->" << l.b;
115  return os;
116}
117
118bool cross(Line a, Line b)
119{
120  Point ao=rot90(a.b-a.a);
121  Point bo=rot90(b.b-b.a);
122  return (ao*(b.a-a.a))*(ao*(b.b-a.a))<0 &&
123    (bo*(a.a-b.a))*(bo*(a.b-b.a))<0;
124}
125
126struct Pedge
127{
128  Node a;
129  Node b;
130  double len;
131};
132
133bool pedgeLess(Pedge a,Pedge b)
134{
135  return a.len<b.len;
136}
137
138std::vector<UEdge> edges;
139
140void triangle()
141{
142  Counter cnt("Number of edges added: ");
143  std::vector<Pedge> pedges;
144  for(NodeIt n(g);n!=INVALID;++n)
145    for(NodeIt m=++(NodeIt(n));m!=INVALID;++m)
146      {
147        Pedge p;
148        p.a=n;
149        p.b=m;
150        p.len=(coords[m]-coords[n]).normSquare();
151        pedges.push_back(p);
152      }
153  std::sort(pedges.begin(),pedges.end(),pedgeLess);
154  for(std::vector<Pedge>::iterator pi=pedges.begin();pi!=pedges.end();++pi)
155    {
156      Line li(pi->a,pi->b);
157      UEdgeIt e(g);
158      for(;e!=INVALID && !cross(e,li);++e) ;
159      UEdge ne;
160      if(e==INVALID) {
161        ne=g.addEdge(pi->a,pi->b);
162        edges.push_back(ne);
163        cnt++;
164      }
165    }
166}
167
168void sparse(int d)
169{
170  Counter cnt("Number of edges removed: ");
171  Bfs<ListUGraph> bfs(g);
172  for(std::vector<UEdge>::reverse_iterator ei=edges.rbegin();
173      ei!=edges.rend();++ei)
174    {
175      Node a=g.source(*ei);
176      Node b=g.target(*ei);
177      g.erase(*ei);
178      bfs.run(a,b);
179      if(bfs.predEdge(b)==INVALID || bfs.dist(b)>d)
180        g.addEdge(a,b);
181      else cnt++;
182    }
183}
184
185void sparse2(int d)
186{
187  Counter cnt("Number of edges removed: ");
188  for(std::vector<UEdge>::reverse_iterator ei=edges.rbegin();
189      ei!=edges.rend();++ei)
190    {
191      Node a=g.source(*ei);
192      Node b=g.target(*ei);
193      g.erase(*ei);
194      ConstMap<Edge,int> cegy(1);
195      Suurballe<ListUGraph,ConstMap<Edge,int> > sur(g,cegy,a,b);
196      int k=sur.run(2);
197      if(k<2 || sur.totalLength()>d)
198        g.addEdge(a,b);
199      else cnt++;
200//       else std::cout << "Remove edge " << g.id(a) << "-" << g.id(b) << '\n';
201    }
202}
203
204void sparseTriangle(int d)
205{
206  Counter cnt("Number of edges added: ");
207  std::vector<Pedge> pedges;
208  for(NodeIt n(g);n!=INVALID;++n)
209    for(NodeIt m=++(NodeIt(n));m!=INVALID;++m)
210      {
211        Pedge p;
212        p.a=n;
213        p.b=m;
214        p.len=(coords[m]-coords[n]).normSquare();
215        pedges.push_back(p);
216      }
217  std::sort(pedges.begin(),pedges.end(),pedgeLess);
218  for(std::vector<Pedge>::iterator pi=pedges.begin();pi!=pedges.end();++pi)
219    {
220      Line li(pi->a,pi->b);
221      UEdgeIt e(g);
222      for(;e!=INVALID && !cross(e,li);++e) ;
223      UEdge ne;
224      if(e==INVALID) {
225        ConstMap<Edge,int> cegy(1);
226        Suurballe<ListUGraph,ConstMap<Edge,int> >
227          sur(g,cegy,pi->a,pi->b);
228        int k=sur.run(2);
229        if(k<2 || sur.totalLength()>d)
230          {
231            ne=g.addEdge(pi->a,pi->b);
232            edges.push_back(ne);
233            cnt++;
234          }
235      }
236    }
237}
238
239void minTree() {
240  int en=0;
241  int pr=0;
242  std::vector<Pedge> pedges;
243  Timer T;
244  std::cout << T.realTime() << "s: Setting up the edges...\n";
245  for(NodeIt n(g);n!=INVALID;++n)
246    {
247      for(NodeIt m=++(NodeIt(n));m!=INVALID;++m)
248        {
249          Pedge p;
250          p.a=n;
251          p.b=m;
252          p.len=(coords[m]-coords[n]).normSquare();
253          pedges.push_back(p);
254        }
255      en++;
256      if(progress && en>=pr*double(N)/100)
257        {
258          std::cout << pr << "%  \r" << std::flush;
259          pr++;
260        }
261    }
262  std::cout << T.realTime() << "s: Sorting the edges...\n";
263  std::sort(pedges.begin(),pedges.end(),pedgeLess);
264  std::cout << T.realTime() << "s: Creating the tree...\n";
265  ListUGraph::NodeMap<int> comp(g);
266  UnionFind<ListUGraph::NodeMap<int> > uf(comp);
267  for (NodeIt it(g); it != INVALID; ++it)
268    uf.insert(it); 
269  for(std::vector<Pedge>::iterator pi=pedges.begin();pi!=pedges.end();++pi)
270    {
271      if ( uf.join(pi->a,pi->b) ) {
272        g.addEdge(pi->a,pi->b);
273        if(en>=N-1)
274          break;
275      }
276    }
277  std::cout << T.realTime() << "s: Done\n";
278}
279
280
281
282int main(int argc,char **argv)
283{
284  ArgParser ap(argc,argv);
285
286//   bool eps;
287  bool disc_d, square_d, gauss_d;
288//   bool tsp_a,two_a,tree_a;
289  int num_of_cities=1;
290  double area=1;
291  N=100;
292//   girth=10;
293  std::string ndist("disc");
294  ap.refOption("n", "Number of nodes (default is 100)", N)
295    .intOption("g", "Girth parameter (default is 10)", 10)
296    .refOption("cities", "Number of cities (default is 1)", num_of_cities)
297    .refOption("area", "Full relative area of the cities (default is 1)", area)
298    .refOption("disc", "Nodes are evenly distributed on a unit disc (default)",disc_d)
299    .optionGroup("dist", "disc")
300    .refOption("square", "Nodes are evenly distributed on a unit square", square_d)
301    .optionGroup("dist", "square")
302    .refOption("gauss",
303            "Nodes are located according to a two-dim gauss distribution",
304            gauss_d)
305    .optionGroup("dist", "gauss")
306//     .mandatoryGroup("dist")
307    .onlyOneGroup("dist")
308    .boolOption("eps", "Also generate .eps output (prefix.eps)")
309    .boolOption("2con", "Create a two connected planar graph")
310    .optionGroup("alg","2con")
311    .boolOption("tree", "Create a min. cost spanning tree")
312    .optionGroup("alg","tree")
313    .boolOption("tsp", "Create a TSP tour")
314    .optionGroup("alg","tsp")
315    .onlyOneGroup("alg")
316    .other("[prefix]","Prefix of the output files. Default is 'lgf-gen-out'")
317    .run();
318 
319  std::string prefix;
320  switch(ap.files().size())
321    {
322    case 0:
323      prefix="lgf-gen-out";
324      break;
325    case 1:
326      prefix=ap.files()[0];
327      break;
328    default:
329      std::cerr << "\nAt most one prefix can be given\n\n";
330      exit(1);
331    }
332 
333  double sum_sizes=0;
334  std::vector<double> sizes;
335  std::vector<double> cum_sizes;
336  for(int s=0;s<num_of_cities;s++)
337    {
338      //        sum_sizes+=rnd.exponential();
339      double d=rnd();
340      sum_sizes+=d;
341      sizes.push_back(d);
342      cum_sizes.push_back(sum_sizes);
343    }
344  int i=0;
345  for(int s=0;s<num_of_cities;s++)
346    {
347      Point center=(num_of_cities==1?Point(0,0):rnd.disc());
348      if(gauss_d)
349        for(;i<N*(cum_sizes[s]/sum_sizes);i++) {
350          Node n=g.addNode();
351          nodes.push_back(n);
352          coords[n]=center+rnd.gauss2()*area*
353            std::sqrt(sizes[s]/sum_sizes);
354        }
355      else if(square_d)
356        for(;i<N*(cum_sizes[s]/sum_sizes);i++) {
357          Node n=g.addNode();
358          nodes.push_back(n);
359          coords[n]=center+Point(rnd()*2-1,rnd()*2-1)*area*
360            std::sqrt(sizes[s]/sum_sizes);
361        }
362      else if(disc_d || true)
363        for(;i<N*(cum_sizes[s]/sum_sizes);i++) {
364          Node n=g.addNode();
365          nodes.push_back(n);
366          coords[n]=center+rnd.disc()*area*
367            std::sqrt(sizes[s]/sum_sizes);
368        }
369    }
370 
371  if(ap["tsp"]) {
372    tsp();
373    std::cout << "#2-opt improvements: " << tsp_impr_num << std::endl;
374  }
375  else if(ap["2con"]) {
376    std::cout << "Make triangles\n";
377    //   triangle();
378    sparseTriangle(ap["g"]);
379    std::cout << "Make it sparser\n";
380    sparse2(ap["g"]);
381  }
382  else if(ap["tree"]) {
383    minTree();
384  }
385 
386
387  std::cout << "Number of nodes    : " << countNodes(g) << std::endl;
388  std::cout << "Number of edges    : " << countUEdges(g) << std::endl;
389  double tlen=0;
390  for(UEdgeIt e(g);e!=INVALID;++e)
391    tlen+=sqrt((coords[g.source(e)]-coords[g.target(e)]).normSquare());
392  std::cout << "Total edge length  : " << tlen << std::endl;
393  if(ap["eps"])
394    graphToEps(g,prefix+".eps").
395      scale(600).nodeScale(.2).edgeWidthScale(.001).preScale(false).
396      coords(coords).run();
397
398  UGraphWriter<ListUGraph>(prefix+".lgf",g).
399    writeNodeMap("coordinates_x",scaleMap(xMap(coords),600)).
400    writeNodeMap("coordinates_y",scaleMap(yMap(coords),600)).
401    run();
402}
403
Note: See TracBrowser for help on using the repository browser.