ceph/src/boost/boost/graph/fruchterman_reingold.hpp

   1 // Copyright 2004, 2005 The Trustees of Indiana University.
   2
   3 // Distributed under the Boost Software License, Version 1.0.
   4 // (See accompanying file LICENSE_1_0.txt or copy at
   5 // http://www.boost.org/LICENSE_1_0.txt)
   6
   7 //  Authors: Douglas Gregor
   8 //           Andrew Lumsdaine
   9 #ifndef BOOST_GRAPH_FRUCHTERMAN_REINGOLD_FORCE_DIRECTED_LAYOUT_HPP
  10 #define BOOST_GRAPH_FRUCHTERMAN_REINGOLD_FORCE_DIRECTED_LAYOUT_HPP
  11
  12 #include <boost/config/no_tr1/cmath.hpp>
  13 #include <boost/graph/graph_traits.hpp>
  14 #include <boost/graph/named_function_params.hpp>
  15 #include <boost/graph/iteration_macros.hpp>
  16 #include <boost/graph/topology.hpp> // For topology concepts
  17 #include <boost/graph/detail/mpi_include.hpp>
  18 #include <vector>
  19 #include <list>
  20 #include <algorithm> // for std::min and std::max
  21 #include <numeric> // for std::accumulate
  22 #include <cmath> // for std::sqrt and std::fabs
  23 #include <functional>
  24
  25 namespace boost
  26 {
  27
  28 struct square_distance_attractive_force
  29 {
  30     template < typename Graph, typename T >
  31     T operator()(typename graph_traits< Graph >::edge_descriptor, T k, T d,
  32         const Graph&) const
  33     {
  34         return d * d / k;
  35     }
  36 };
  37
  38 struct square_distance_repulsive_force
  39 {
  40     template < typename Graph, typename T >
  41     T operator()(typename graph_traits< Graph >::vertex_descriptor,
  42         typename graph_traits< Graph >::vertex_descriptor, T k, T d,
  43         const Graph&) const
  44     {
  45         return k * k / d;
  46     }
  47 };
  48
  49 template < typename T > struct linear_cooling
  50 {
  51     typedef T result_type;
  52
  53     linear_cooling(std::size_t iterations)
  54     : temp(T(iterations) / T(10)), step(0.1)
  55     {
  56     }
  57
  58     linear_cooling(std::size_t iterations, T temp)
  59     : temp(temp), step(temp / T(iterations))
  60     {
  61     }
  62
  63     T operator()()
  64     {
  65         T old_temp = temp;
  66         temp -= step;
  67         if (temp < T(0))
  68             temp = T(0);
  69         return old_temp;
  70     }
  71
  72 private:
  73     T temp;
  74     T step;
  75 };
  76
  77 struct all_force_pairs
  78 {
  79     template < typename Graph, typename ApplyForce >
  80     void operator()(const Graph& g, ApplyForce apply_force)
  81     {
  82         typedef typename graph_traits< Graph >::vertex_iterator vertex_iterator;
  83         vertex_iterator v, end;
  84         for (boost::tie(v, end) = vertices(g); v != end; ++v)
  85         {
  86             vertex_iterator u = v;
  87             for (++u; u != end; ++u)
  88             {
  89                 apply_force(*u, *v);
  90                 apply_force(*v, *u);
  91             }
  92         }
  93     }
  94 };
  95
  96 template < typename Topology, typename PositionMap > struct grid_force_pairs
  97 {
  98     typedef typename property_traits< PositionMap >::value_type Point;
  99     BOOST_STATIC_ASSERT(Point::dimensions == 2);
 100     typedef typename Topology::point_difference_type point_difference_type;
 101
 102     template < typename Graph >
 103     explicit grid_force_pairs(
 104         const Topology& topology, PositionMap position, const Graph& g)
 105     : topology(topology), position(position)
 106     {
 107         two_k = 2. * this->topology.volume(this->topology.extent())
 108             / std::sqrt((double)num_vertices(g));
 109     }
 110
 111     template < typename Graph, typename ApplyForce >
 112     void operator()(const Graph& g, ApplyForce apply_force)
 113     {
 114         typedef typename graph_traits< Graph >::vertex_iterator vertex_iterator;
 115         typedef
 116             typename graph_traits< Graph >::vertex_descriptor vertex_descriptor;
 117         typedef std::list< vertex_descriptor > bucket_t;
 118         typedef std::vector< bucket_t > buckets_t;
 119
 120         std::size_t columns = std::size_t(topology.extent()[0] / two_k + 1.);
 121         std::size_t rows = std::size_t(topology.extent()[1] / two_k + 1.);
 122         buckets_t buckets(rows * columns);
 123         vertex_iterator v, v_end;
 124         for (boost::tie(v, v_end) = vertices(g); v != v_end; ++v)
 125         {
 126             std::size_t column = std::size_t(
 127                 (get(position, *v)[0] + topology.extent()[0] / 2) / two_k);
 128             std::size_t row = std::size_t(
 129                 (get(position, *v)[1] + topology.extent()[1] / 2) / two_k);
 130
 131             if (column >= columns)
 132                 column = columns - 1;
 133             if (row >= rows)
 134                 row = rows - 1;
 135             buckets[row * columns + column].push_back(*v);
 136         }
 137
 138         for (std::size_t row = 0; row < rows; ++row)
 139             for (std::size_t column = 0; column < columns; ++column)
 140             {
 141                 bucket_t& bucket = buckets[row * columns + column];
 142                 typedef typename bucket_t::iterator bucket_iterator;
 143                 for (bucket_iterator u = bucket.begin(); u != bucket.end(); ++u)
 144                 {
 145                     // Repulse vertices in this bucket
 146                     bucket_iterator v = u;
 147                     for (++v; v != bucket.end(); ++v)
 148                     {
 149                         apply_force(*u, *v);
 150                         apply_force(*v, *u);
 151                     }
 152
 153                     std::size_t adj_start_row = row == 0 ? 0 : row - 1;
 154                     std::size_t adj_end_row = row == rows - 1 ? row : row + 1;
 155                     std::size_t adj_start_column = column == 0 ? 0 : column - 1;
 156                     std::size_t adj_end_column
 157                         = column == columns - 1 ? column : column + 1;
 158                     for (std::size_t other_row = adj_start_row;
 159                          other_row <= adj_end_row; ++other_row)
 160                         for (std::size_t other_column = adj_start_column;
 161                              other_column <= adj_end_column; ++other_column)
 162                             if (other_row != row || other_column != column)
 163                             {
 164                                 // Repulse vertices in this bucket
 165                                 bucket_t& other_bucket
 166                                     = buckets[other_row * columns
 167                                         + other_column];
 168                                 for (v = other_bucket.begin();
 169                                      v != other_bucket.end(); ++v)
 170                                 {
 171                                     double dist = topology.distance(
 172                                         get(position, *u), get(position, *v));
 173                                     if (dist < two_k)
 174                                         apply_force(*u, *v);
 175                                 }
 176                             }
 177                 }
 178             }
 179     }
 180
 181 private:
 182     const Topology& topology;
 183     PositionMap position;
 184     double two_k;
 185 };
 186
 187 template < typename PositionMap, typename Topology, typename Graph >
 188 inline grid_force_pairs< Topology, PositionMap > make_grid_force_pairs(
 189     const Topology& topology, const PositionMap& position, const Graph& g)
 190 {
 191     return grid_force_pairs< Topology, PositionMap >(topology, position, g);
 192 }
 193
 194 template < typename Graph, typename PositionMap, typename Topology >
 195 void scale_graph(const Graph& g, PositionMap position, const Topology& topology,
 196     typename Topology::point_type upper_left,
 197     typename Topology::point_type lower_right)
 198 {
 199     if (num_vertices(g) == 0)
 200         return;
 201
 202     typedef typename Topology::point_type Point;
 203     typedef typename Topology::point_difference_type point_difference_type;
 204
 205     // Find min/max ranges
 206     Point min_point = get(position, *vertices(g).first), max_point = min_point;
 207     BGL_FORALL_VERTICES_T(v, g, Graph)
 208     {
 209         min_point = topology.pointwise_min(min_point, get(position, v));
 210         max_point = topology.pointwise_max(max_point, get(position, v));
 211     }
 212
 213     Point old_origin = topology.move_position_toward(min_point, 0.5, max_point);
 214     Point new_origin
 215         = topology.move_position_toward(upper_left, 0.5, lower_right);
 216     point_difference_type old_size = topology.difference(max_point, min_point);
 217     point_difference_type new_size
 218         = topology.difference(lower_right, upper_left);
 219
 220     // Scale to bounding box provided
 221     BGL_FORALL_VERTICES_T(v, g, Graph)
 222     {
 223         point_difference_type relative_loc
 224             = topology.difference(get(position, v), old_origin);
 225         relative_loc = (relative_loc / old_size) * new_size;
 226         put(position, v, topology.adjust(new_origin, relative_loc));
 227     }
 228 }
 229
 230 namespace detail
 231 {
 232     template < typename Topology, typename PropMap, typename Vertex >
 233     void maybe_jitter_point(const Topology& topology, const PropMap& pm,
 234         Vertex v, const typename Topology::point_type& p2)
 235     {
 236         double too_close = topology.norm(topology.extent()) / 10000.;
 237         if (topology.distance(get(pm, v), p2) < too_close)
 238         {
 239             put(pm, v,
 240                 topology.move_position_toward(
 241                     get(pm, v), 1. / 200, topology.random_point()));
 242         }
 243     }
 244
 245     template < typename Topology, typename PositionMap,
 246         typename DisplacementMap, typename RepulsiveForce, typename Graph >
 247     struct fr_apply_force
 248     {
 249         typedef
 250             typename graph_traits< Graph >::vertex_descriptor vertex_descriptor;
 251         typedef typename Topology::point_type Point;
 252         typedef typename Topology::point_difference_type PointDiff;
 253
 254         fr_apply_force(const Topology& topology, const PositionMap& position,
 255             const DisplacementMap& displacement, RepulsiveForce repulsive_force,
 256             double k, const Graph& g)
 257         : topology(topology)
 258         , position(position)
 259         , displacement(displacement)
 260         , repulsive_force(repulsive_force)
 261         , k(k)
 262         , g(g)
 263         {
 264         }
 265
 266         void operator()(vertex_descriptor u, vertex_descriptor v)
 267         {
 268             if (u != v)
 269             {
 270                 // When the vertices land on top of each other, move the
 271                 // first vertex away from the boundaries.
 272                 maybe_jitter_point(topology, position, u, get(position, v));
 273
 274                 double dist
 275                     = topology.distance(get(position, u), get(position, v));
 276                 typename Topology::point_difference_type dispv
 277                     = get(displacement, v);
 278                 if (dist == 0.)
 279                 {
 280                     for (std::size_t i = 0; i < Point::dimensions; ++i)
 281                     {
 282                         dispv[i] += 0.01;
 283                     }
 284                 }
 285                 else
 286                 {
 287                     double fr = repulsive_force(u, v, k, dist, g);
 288                     dispv += (fr / dist)
 289                         * topology.difference(
 290                             get(position, v), get(position, u));
 291                 }
 292                 put(displacement, v, dispv);
 293             }
 294         }
 295
 296     private:
 297         const Topology& topology;
 298         PositionMap position;
 299         DisplacementMap displacement;
 300         RepulsiveForce repulsive_force;
 301         double k;
 302         const Graph& g;
 303     };
 304
 305 } // end namespace detail
 306
 307 template < typename Topology, typename Graph, typename PositionMap,
 308     typename AttractiveForce, typename RepulsiveForce, typename ForcePairs,
 309     typename Cooling, typename DisplacementMap >
 310 void fruchterman_reingold_force_directed_layout(const Graph& g,
 311     PositionMap position, const Topology& topology,
 312     AttractiveForce attractive_force, RepulsiveForce repulsive_force,
 313     ForcePairs force_pairs, Cooling cool, DisplacementMap displacement)
 314 {
 315     typedef typename graph_traits< Graph >::vertex_iterator vertex_iterator;
 316     typedef typename graph_traits< Graph >::vertex_descriptor vertex_descriptor;
 317     typedef typename graph_traits< Graph >::edge_iterator edge_iterator;
 318
 319     double volume = topology.volume(topology.extent());
 320
 321     // assume positions are initialized randomly
 322     double k = pow(volume / num_vertices(g),
 323         1. / (double)(Topology::point_difference_type::dimensions));
 324
 325     detail::fr_apply_force< Topology, PositionMap, DisplacementMap,
 326         RepulsiveForce, Graph >
 327         apply_force(topology, position, displacement, repulsive_force, k, g);
 328
 329     do
 330     {
 331         // Calculate repulsive forces
 332         vertex_iterator v, v_end;
 333         for (boost::tie(v, v_end) = vertices(g); v != v_end; ++v)
 334             put(displacement, *v, typename Topology::point_difference_type());
 335         force_pairs(g, apply_force);
 336
 337         // Calculate attractive forces
 338         edge_iterator e, e_end;
 339         for (boost::tie(e, e_end) = edges(g); e != e_end; ++e)
 340         {
 341             vertex_descriptor v = source(*e, g);
 342             vertex_descriptor u = target(*e, g);
 343
 344             // When the vertices land on top of each other, move the
 345             // first vertex away from the boundaries.
 346             ::boost::detail::maybe_jitter_point(
 347                 topology, position, u, get(position, v));
 348
 349             typename Topology::point_difference_type delta
 350                 = topology.difference(get(position, v), get(position, u));
 351             double dist = topology.distance(get(position, u), get(position, v));
 352             double fa = attractive_force(*e, k, dist, g);
 353
 354             put(displacement, v, get(displacement, v) - (fa / dist) * delta);
 355             put(displacement, u, get(displacement, u) + (fa / dist) * delta);
 356         }
 357
 358         if (double temp = cool())
 359         {
 360             // Update positions
 361             BGL_FORALL_VERTICES_T(v, g, Graph)
 362             {
 363                 BOOST_USING_STD_MIN();
 364                 BOOST_USING_STD_MAX();
 365                 double disp_size = topology.norm(get(displacement, v));
 366                 put(position, v,
 367                     topology.adjust(get(position, v),
 368                         get(displacement, v)
 369                             * (min BOOST_PREVENT_MACRO_SUBSTITUTION(
 370                                    disp_size, temp)
 371                                 / disp_size)));
 372                 put(position, v, topology.bound(get(position, v)));
 373             }
 374         }
 375         else
 376         {
 377             break;
 378         }
 379     } while (true);
 380 }
 381
 382 namespace detail
 383 {
 384     template < typename DisplacementMap > struct fr_force_directed_layout
 385     {
 386         template < typename Topology, typename Graph, typename PositionMap,
 387             typename AttractiveForce, typename RepulsiveForce,
 388             typename ForcePairs, typename Cooling, typename Param, typename Tag,
 389             typename Rest >
 390         static void run(const Graph& g, PositionMap position,
 391             const Topology& topology, AttractiveForce attractive_force,
 392             RepulsiveForce repulsive_force, ForcePairs force_pairs,
 393             Cooling cool, DisplacementMap displacement,
 394             const bgl_named_params< Param, Tag, Rest >&)
 395         {
 396             fruchterman_reingold_force_directed_layout(g, position, topology,
 397                 attractive_force, repulsive_force, force_pairs, cool,
 398                 displacement);
 399         }
 400     };
 401
 402     template <> struct fr_force_directed_layout< param_not_found >
 403     {
 404         template < typename Topology, typename Graph, typename PositionMap,
 405             typename AttractiveForce, typename RepulsiveForce,
 406             typename ForcePairs, typename Cooling, typename Param, typename Tag,
 407             typename Rest >
 408         static void run(const Graph& g, PositionMap position,
 409             const Topology& topology, AttractiveForce attractive_force,
 410             RepulsiveForce repulsive_force, ForcePairs force_pairs,
 411             Cooling cool, param_not_found,
 412             const bgl_named_params< Param, Tag, Rest >& params)
 413         {
 414             typedef typename Topology::point_difference_type PointDiff;
 415             std::vector< PointDiff > displacements(num_vertices(g));
 416             fruchterman_reingold_force_directed_layout(g, position, topology,
 417                 attractive_force, repulsive_force, force_pairs, cool,
 418                 make_iterator_property_map(displacements.begin(),
 419                     choose_const_pmap(
 420                         get_param(params, vertex_index), g, vertex_index),
 421                     PointDiff()));
 422         }
 423     };
 424
 425 } // end namespace detail
 426
 427 template < typename Topology, typename Graph, typename PositionMap,
 428     typename Param, typename Tag, typename Rest >
 429 void fruchterman_reingold_force_directed_layout(const Graph& g,
 430     PositionMap position, const Topology& topology,
 431     const bgl_named_params< Param, Tag, Rest >& params)
 432 {
 433     typedef typename get_param_type< vertex_displacement_t,
 434         bgl_named_params< Param, Tag, Rest > >::type D;
 435
 436     detail::fr_force_directed_layout< D >::run(g, position, topology,
 437         choose_param(get_param(params, attractive_force_t()),
 438             square_distance_attractive_force()),
 439         choose_param(get_param(params, repulsive_force_t()),
 440             square_distance_repulsive_force()),
 441         choose_param(get_param(params, force_pairs_t()),
 442             make_grid_force_pairs(topology, position, g)),
 443         choose_param(
 444             get_param(params, cooling_t()), linear_cooling< double >(100)),
 445         get_param(params, vertex_displacement_t()), params);
 446 }
 447
 448 template < typename Topology, typename Graph, typename PositionMap >
 449 void fruchterman_reingold_force_directed_layout(
 450     const Graph& g, PositionMap position, const Topology& topology)
 451 {
 452     fruchterman_reingold_force_directed_layout(g, position, topology,
 453         attractive_force(square_distance_attractive_force()));
 454 }
 455
 456 } // end namespace boost
 457
 458 #include BOOST_GRAPH_MPI_INCLUDE(<boost/graph/distributed/fruchterman_reingold.hpp>)
 459
 460 #endif // BOOST_GRAPH_FRUCHTERMAN_REINGOLD_FORCE_DIRECTED_LAYOUT_HPP