examples/concurrent_priority_queue/shortpath/shortpath.cpp

   1 /*
   2     Copyright (c) 2005-2019 Intel Corporation
   3
   4     Licensed under the Apache License, Version 2.0 (the "License");
   5     you may not use this file except in compliance with the License.
   6     You may obtain a copy of the License at
   7
   8         http://www.apache.org/licenses/LICENSE-2.0
   9
  10     Unless required by applicable law or agreed to in writing, software
  11     distributed under the License is distributed on an "AS IS" BASIS,
  12     WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  13     See the License for the specific language governing permissions and
  14     limitations under the License.
  15 */
  16
  17 #include <cstdio>
  18 #include <vector>
  19 #include <atomic>
  20 #include <math.h>
  21
  22 #include "tbb/tick_count.h"
  23 #include "tbb/task_group.h"
  24 #include "tbb/concurrent_priority_queue.h"
  25 #include "tbb/spin_mutex.h"
  26 #include "tbb/parallel_for.h"
  27 #include "tbb/blocked_range.h"
  28 #include "tbb/global_control.h"
  29 #include "../../common/utility/utility.h"
  30 #include "../../common/utility/fast_random.h"
  31 #include "../../common/utility/get_default_num_threads.h"
  32
  33 #if defined(_MSC_VER) && defined(_Wp64)
  34     // Workaround for overzealous compiler warnings in /Wp64 mode
  35     #pragma warning (disable: 4267)
  36 #endif /* _MSC_VER && _Wp64 */
  37
  38 using namespace std;
  39 using namespace tbb;
  40
  41 struct point {
  42     double x, y;
  43     point() {}
  44     point(double _x, double _y) : x(_x), y(_y) {}
  45     point(const point& p) : x(p.x), y(p.y) {}
  46 };
  47
  48 double get_distance(const point& p1, const point& p2) {
  49     double xdiff=p1.x-p2.x, ydiff=p1.y-p2.y;
  50     return sqrt(xdiff*xdiff + ydiff*ydiff);
  51 }
  52
  53 // generates random points on 2D plane within a box of maxsize width & height
  54 point generate_random_point(utility::FastRandom& mr) {
  55     const size_t maxsize=500;
  56     double x = (double)(mr.get() % maxsize);
  57     double y = (double)(mr.get() % maxsize);
  58     return point(x,y);
  59 }
  60
  61 // weighted toss makes closer nodes (in the point vector) heavily connected
  62 bool die_toss(size_t a, size_t b, utility::FastRandom& mr) {
  63     int node_diff = std::abs((int)(a-b));
  64     // near nodes
  65     if (node_diff < 16) return true;
  66     // mid nodes
  67     if (node_diff < 64) return ((int)mr.get() % 8 == 0);
  68     // far nodes
  69     if (node_diff < 512) return ((int)mr.get() % 16 == 0);
  70     return false;
  71 }
  72
  73 typedef vector<point> point_set;
  74 typedef size_t vertex_id;
  75 typedef std::pair<vertex_id,double> vertex_rec;
  76 typedef vector<vector<vertex_id> > edge_set;
  77
  78 bool verbose = false;                       // prints bin details and other diagnostics to screen
  79 bool silent = false;                        // suppress all output except for time
  80 size_t N = 1000;                            // number of vertices
  81 size_t src = 0;                             // start of path
  82 size_t dst = N-1;                           // end of path
  83 double INF=100000.0;                        // infinity
  84 size_t grainsize = 16;                      // number of vertices per task on average
  85 size_t max_spawn;                           // max tasks to spawn
  86 std::atomic<size_t> num_spawn;              // number of active tasks
  87
  88 point_set vertices;                         // vertices
  89 edge_set edges;                             // edges
  90 vector<vertex_id> predecessor;              // for recreating path from src to dst
  91
  92 vector<double> f_distance;                  // estimated distances at particular vertex
  93 vector<double> g_distance;                  // current shortest distances from src vertex
  94 spin_mutex    *locks;                       // a lock for each vertex
  95 task_group *sp_group;                       // task group for tasks executing sub-problems
  96
  97 class compare_f {
  98 public:
  99     bool operator()(const vertex_rec& u, const vertex_rec& v) const {
 100         return u.second>v.second;
 101     }
 102 };
 103
 104 concurrent_priority_queue<vertex_rec, compare_f> open_set; // tentative vertices
 105
 106 void shortpath_helper();
 107
 108 #if !__TBB_CPP11_LAMBDAS_PRESENT
 109 class shortpath_helper_functor {
 110 public:
 111     shortpath_helper_functor() {};
 112     void operator() () const { shortpath_helper(); }
 113 };
 114 #endif
 115
 116 void shortpath() {
 117     sp_group = new task_group;
 118     g_distance[src] = 0.0; // src's distance from src is zero
 119     f_distance[src] = get_distance(vertices[src], vertices[dst]); // estimate distance from src to dst
 120     open_set.push(make_pair(src,f_distance[src])); // push src into open_set
 121 #if __TBB_CPP11_LAMBDAS_PRESENT
 122     sp_group->run([](){ shortpath_helper(); });
 123 #else
 124     sp_group->run( shortpath_helper_functor() );
 125 #endif
 126     sp_group->wait();
 127     delete sp_group;
 128 }
 129
 130 void shortpath_helper() {
 131     vertex_rec u_rec;
 132     while (open_set.try_pop(u_rec)) {
 133         vertex_id u = u_rec.first;
 134         if (u==dst) continue;
 135         double f = u_rec.second;
 136         double old_g_u = 0.0;
 137         {
 138             spin_mutex::scoped_lock l(locks[u]);
 139             if (f > f_distance[u]) continue; // prune search space
 140             old_g_u = g_distance[u];
 141         }
 142         for (size_t i=0; i<edges[u].size(); ++i) {
 143             vertex_id v = edges[u][i];
 144             double new_g_v = old_g_u + get_distance(vertices[u], vertices[v]);
 145             double new_f_v = 0.0;
 146             // the push flag lets us move some work out of the critical section below
 147             bool push = false;
 148             {
 149                 spin_mutex::scoped_lock l(locks[v]);
 150                 if (new_g_v < g_distance[v]) {
 151                     predecessor[v] = u;
 152                     g_distance[v] = new_g_v;
 153                     new_f_v = f_distance[v] = g_distance[v] + get_distance(vertices[v], vertices[dst]);
 154                     push = true;
 155                 }
 156             }
 157             if (push) {
 158                 open_set.push(make_pair(v,new_f_v));
 159                 size_t n_spawn = ++num_spawn;
 160                 if (n_spawn < max_spawn) {
 161 #if __TBB_CPP11_LAMBDAS_PRESENT
 162                     sp_group->run([]{ shortpath_helper(); });
 163 #else
 164                     sp_group->run( shortpath_helper_functor() );
 165 #endif
 166                 }
 167                 else --num_spawn;
 168             }
 169         }
 170     }
 171     --num_spawn;
 172 }
 173
 174 void make_path(vertex_id src, vertex_id dst, vector<vertex_id>& path) {
 175     vertex_id at = predecessor[dst];
 176     if (at == N) path.push_back(src);
 177     else if (at == src) { path.push_back(src); path.push_back(dst); }
 178     else { make_path(src, at, path); path.push_back(dst); }
 179 }
 180
 181 void print_path() {
 182     vector<vertex_id> path;
 183     double path_length=0.0;
 184     make_path(src, dst, path);
 185     if (verbose) printf("\n      ");
 186     for (size_t i=0; i<path.size(); ++i) {
 187         if (path[i] != dst) {
 188             double seg_length = get_distance(vertices[path[i]], vertices[path[i+1]]);
 189             if (verbose) printf("%6.1f       ", seg_length);
 190             path_length += seg_length;
 191         }
 192         else if (verbose) printf("\n");
 193     }
 194     if (verbose) {
 195         for (size_t i=0; i<path.size(); ++i) {
 196             if (path[i] != dst) printf("(%4d)------>", (int)path[i]);
 197             else printf("(%4d)\n", (int)path[i]);
 198         }
 199     }
 200     if (verbose) printf("Total distance = %5.1f\n", path_length);
 201     else if (!silent) printf(" %5.1f\n", path_length);
 202 }
 203
 204 #if !__TBB_CPP11_LAMBDAS_PRESENT
 205 class gen_vertices {
 206 public:
 207     gen_vertices() {}
 208     void operator() (blocked_range<size_t>& r) const {
 209         utility::FastRandom my_random((unsigned int)r.begin());
 210         for (size_t i=r.begin(); i!=r.end(); ++i) {
 211             vertices[i] = generate_random_point(my_random);
 212         }
 213     }
 214 };
 215
 216 class gen_edges {
 217 public:
 218     gen_edges() {}
 219     void operator() (blocked_range<size_t>& r) const {
 220         utility::FastRandom my_random((unsigned int)r.begin());
 221         for (size_t i=r.begin(); i!=r.end(); ++i) {
 222             for (size_t j=0; j<i; ++j) {
 223                 if (die_toss(i, j, my_random))
 224                     edges[i].push_back(j);
 225             }
 226         }
 227     }
 228 };
 229
 230 class reset_vertices {
 231 public:
 232     reset_vertices() {}
 233     void operator() (blocked_range<size_t>& r) const {
 234         for (size_t i=r.begin(); i!=r.end(); ++i) {
 235             f_distance[i] = g_distance[i] = INF;
 236             predecessor[i] = N;
 237         }
 238     }
 239 };
 240 #endif
 241
 242 void InitializeGraph() {
 243     global_control c(tbb::global_control::max_allowed_parallelism, utility::get_default_num_threads());
 244     vertices.resize(N);
 245     edges.resize(N);
 246     predecessor.resize(N);
 247     g_distance.resize(N);
 248     f_distance.resize(N);
 249     locks = new spin_mutex[N];
 250     if (verbose) printf("Generating vertices...\n");
 251 #if __TBB_CPP11_LAMBDAS_PRESENT
 252     parallel_for(blocked_range<size_t>(0,N,64),
 253                  [&](blocked_range<size_t>& r) {
 254                      utility::FastRandom my_random(r.begin());
 255                      for (size_t i=r.begin(); i!=r.end(); ++i) {
 256                          vertices[i] = generate_random_point(my_random);
 257                      }
 258                  }, simple_partitioner());
 259 #else
 260     parallel_for(blocked_range<size_t>(0,N,64), gen_vertices(), simple_partitioner());
 261 #endif
 262     if (verbose) printf("Generating edges...\n");
 263 #if __TBB_CPP11_LAMBDAS_PRESENT
 264     parallel_for(blocked_range<size_t>(0,N,64),
 265                  [&](blocked_range<size_t>& r) {
 266                      utility::FastRandom my_random(r.begin());
 267                      for (size_t i=r.begin(); i!=r.end(); ++i) {
 268                          for (size_t j=0; j<i; ++j) {
 269                              if (die_toss(i, j, my_random))
 270                                  edges[i].push_back(j);
 271                          }
 272                      }
 273                  }, simple_partitioner());
 274 #else
 275     parallel_for(blocked_range<size_t>(0,N,64), gen_edges(), simple_partitioner());
 276 #endif
 277     for (size_t i=0; i<N; ++i) {
 278         for (size_t j=0; j<edges[i].size(); ++j) {
 279             vertex_id k = edges[i][j];
 280             edges[k].push_back(i);
 281         }
 282     }
 283     if (verbose) printf("Done.\n");
 284 }
 285
 286 void ReleaseGraph() {
 287     delete []locks;
 288 }
 289
 290 void ResetGraph() {
 291     global_control c(tbb::global_control::max_allowed_parallelism, utility::get_default_num_threads());
 292 #if __TBB_CPP11_LAMBDAS_PRESENT
 293     parallel_for(blocked_range<size_t>(0,N),
 294                  [&](blocked_range<size_t>& r) {
 295                      for (size_t i=r.begin(); i!=r.end(); ++i) {
 296                          f_distance[i] = g_distance[i] = INF;
 297                          predecessor[i] = N;
 298                      }
 299                  });
 300 #else
 301     parallel_for(blocked_range<size_t>(0,N), reset_vertices());
 302 #endif
 303 }
 304
 305 int main(int argc, char *argv[]) {
 306     try {
 307         utility::thread_number_range threads(utility::get_default_num_threads);
 308         utility::parse_cli_arguments(argc, argv,
 309                                      utility::cli_argument_pack()
 310                                      //"-h" option for displaying help is present implicitly
 311                                      .positional_arg(threads,"#threads",utility::thread_number_range_desc)
 312                                      .arg(verbose,"verbose","   print diagnostic output to screen")
 313                                      .arg(silent,"silent","    limits output to timing info; overrides verbose")
 314                                      .arg(N,"N","         number of vertices")
 315                                      .arg(src,"start","      start of path")
 316                                      .arg(dst,"end","        end of path")
 317         );
 318         if (silent) verbose = false;  // make silent override verbose
 319         else
 320             printf("shortpath will run with %d vertices to find shortest path between vertices"
 321                    " %d and %d using %d:%d threads.\n",
 322                    (int)N, (int)src, (int)dst, (int)threads.first, (int)threads.last);
 323
 324         if (dst >= N) {
 325             if (verbose)
 326                 printf("end value %d is invalid for %d vertices; correcting to %d\n", (int)dst, (int)N, (int)N-1);
 327             dst = N-1;
 328         }
 329
 330         num_spawn = 0;
 331         max_spawn = N/grainsize;
 332         tick_count t0, t1;
 333         InitializeGraph();
 334         for (int n_thr=threads.first; n_thr<=threads.last; n_thr=threads.step(n_thr)) {
 335             ResetGraph();
 336             global_control c(tbb::global_control::max_allowed_parallelism, n_thr);
 337             t0 = tick_count::now();
 338             shortpath();
 339             t1 = tick_count::now();
 340             if (!silent) {
 341                 if (predecessor[dst] != N) {
 342                     printf("%d threads: [%6.6f] The shortest path from vertex %d to vertex %d is:",
 343                            (int)n_thr, (t1-t0).seconds(), (int)src, (int)dst);
 344                     print_path();
 345                 }
 346                 else {
 347                     printf("%d threads: [%6.6f] There is no path from vertex %d to vertex %d\n",
 348                            (int)n_thr, (t1-t0).seconds(), (int)src, (int)dst);
 349                 }
 350             } else
 351                 utility::report_elapsed_time((t1-t0).seconds());
 352         }
 353         ReleaseGraph();
 354         return 0;
 355     } catch(std::exception& e) {
 356         cerr<<"error occurred. error text is :\"" <<e.what()<<"\"\n";
 357         return 1;
 358     }
 359 }