ceph/src/boost/boost/graph/kamada_kawai_spring_layout.hpp

   1 // Copyright 2004 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_KAMADA_KAWAI_SPRING_LAYOUT_HPP
  10 #define BOOST_GRAPH_KAMADA_KAWAI_SPRING_LAYOUT_HPP
  11
  12 #include <boost/graph/graph_traits.hpp>
  13 #include <boost/graph/topology.hpp>
  14 #include <boost/graph/iteration_macros.hpp>
  15 #include <boost/graph/johnson_all_pairs_shortest.hpp>
  16 #include <boost/type_traits/is_convertible.hpp>
  17 #include <utility>
  18 #include <iterator>
  19 #include <vector>
  20 #include <iostream>
  21 #include <boost/limits.hpp>
  22 #include <boost/config/no_tr1/cmath.hpp>
  23
  24 namespace boost
  25 {
  26 namespace detail
  27 {
  28     namespace graph
  29     {
  30         /**
  31          * Denotes an edge or display area side length used to scale a
  32          * Kamada-Kawai drawing.
  33          */
  34         template < bool Edge, typename T > struct edge_or_side
  35         {
  36             explicit edge_or_side(T value) : value(value) {}
  37
  38             T value;
  39         };
  40
  41         /**
  42          * Compute the edge length from an edge length. This is trivial.
  43          */
  44         template < typename Graph, typename DistanceMap, typename IndexMap,
  45             typename T >
  46         T compute_edge_length(
  47             const Graph&, DistanceMap, IndexMap, edge_or_side< true, T > length)
  48         {
  49             return length.value;
  50         }
  51
  52         /**
  53          * Compute the edge length based on the display area side
  54            length. We do this by dividing the side length by the largest
  55            shortest distance between any two vertices in the graph.
  56          */
  57         template < typename Graph, typename DistanceMap, typename IndexMap,
  58             typename T >
  59         T compute_edge_length(const Graph& g, DistanceMap distance,
  60             IndexMap index, edge_or_side< false, T > length)
  61         {
  62             T result(0);
  63
  64             typedef
  65                 typename graph_traits< Graph >::vertex_iterator vertex_iterator;
  66
  67             for (vertex_iterator ui = vertices(g).first,
  68                                  end = vertices(g).second;
  69                  ui != end; ++ui)
  70             {
  71                 vertex_iterator vi = ui;
  72                 for (++vi; vi != end; ++vi)
  73                 {
  74                     T dij = distance[get(index, *ui)][get(index, *vi)];
  75                     if (dij > result)
  76                         result = dij;
  77                 }
  78             }
  79             return length.value / result;
  80         }
  81
  82         /**
  83          * Dense linear solver for fixed-size matrices.
  84          */
  85         template < std::size_t Size > struct linear_solver
  86         {
  87             // Indices in mat are (row, column)
  88             // template <typename Vec>
  89             // static Vec solve(double mat[Size][Size], Vec rhs);
  90         };
  91
  92         template <> struct linear_solver< 1 >
  93         {
  94             template < typename Vec >
  95             static Vec solve(double mat[1][1], Vec rhs)
  96             {
  97                 return rhs / mat[0][0];
  98             }
  99         };
 100
 101         // These are from http://en.wikipedia.org/wiki/Cramer%27s_rule
 102         template <> struct linear_solver< 2 >
 103         {
 104             template < typename Vec >
 105             static Vec solve(double mat[2][2], Vec rhs)
 106             {
 107                 double denom = mat[0][0] * mat[1][1] - mat[1][0] * mat[0][1];
 108                 double x_num = rhs[0] * mat[1][1] - rhs[1] * mat[0][1];
 109                 double y_num = mat[0][0] * rhs[1] - mat[1][0] * rhs[0];
 110                 Vec result;
 111                 result[0] = x_num / denom;
 112                 result[1] = y_num / denom;
 113                 return result;
 114             }
 115         };
 116
 117         template <> struct linear_solver< 3 >
 118         {
 119             template < typename Vec >
 120             static Vec solve(double mat[3][3], Vec rhs)
 121             {
 122                 double denom = mat[0][0]
 123                         * (mat[1][1] * mat[2][2] - mat[2][1] * mat[1][2])
 124                     - mat[1][0]
 125                         * (mat[0][1] * mat[2][2] - mat[2][1] * mat[0][2])
 126                     + mat[2][0]
 127                         * (mat[0][1] * mat[1][2] - mat[1][1] * mat[0][2]);
 128                 double x_num
 129                     = rhs[0] * (mat[1][1] * mat[2][2] - mat[2][1] * mat[1][2])
 130                     - rhs[1] * (mat[0][1] * mat[2][2] - mat[2][1] * mat[0][2])
 131                     + rhs[2] * (mat[0][1] * mat[1][2] - mat[1][1] * mat[0][2]);
 132                 double y_num
 133                     = mat[0][0] * (rhs[1] * mat[2][2] - rhs[2] * mat[1][2])
 134                     - mat[1][0] * (rhs[0] * mat[2][2] - rhs[2] * mat[0][2])
 135                     + mat[2][0] * (rhs[0] * mat[1][2] - rhs[1] * mat[0][2]);
 136                 double z_num
 137                     = mat[0][0] * (mat[1][1] * rhs[2] - mat[2][1] * rhs[1])
 138                     - mat[1][0] * (mat[0][1] * rhs[2] - mat[2][1] * rhs[0])
 139                     + mat[2][0] * (mat[0][1] * rhs[1] - mat[1][1] * rhs[0]);
 140                 Vec result;
 141                 result[0] = x_num / denom;
 142                 result[1] = y_num / denom;
 143                 result[2] = z_num / denom;
 144                 return result;
 145             }
 146         };
 147
 148         /**
 149          * Implementation of the Kamada-Kawai spring layout algorithm.
 150          */
 151         template < typename Topology, typename Graph, typename PositionMap,
 152             typename WeightMap, typename EdgeOrSideLength, typename Done,
 153             typename VertexIndexMap, typename DistanceMatrix,
 154             typename SpringStrengthMatrix, typename PartialDerivativeMap >
 155         struct kamada_kawai_spring_layout_impl
 156         {
 157             typedef
 158                 typename property_traits< WeightMap >::value_type weight_type;
 159             typedef typename Topology::point_type Point;
 160             typedef
 161                 typename Topology::point_difference_type point_difference_type;
 162             typedef point_difference_type deriv_type;
 163             typedef
 164                 typename graph_traits< Graph >::vertex_iterator vertex_iterator;
 165             typedef typename graph_traits< Graph >::vertex_descriptor
 166                 vertex_descriptor;
 167
 168             kamada_kawai_spring_layout_impl(const Topology& topology,
 169                 const Graph& g, PositionMap position, WeightMap weight,
 170                 EdgeOrSideLength edge_or_side_length, Done done,
 171                 weight_type spring_constant, VertexIndexMap index,
 172                 DistanceMatrix distance, SpringStrengthMatrix spring_strength,
 173                 PartialDerivativeMap partial_derivatives)
 174             : topology(topology)
 175             , g(g)
 176             , position(position)
 177             , weight(weight)
 178             , edge_or_side_length(edge_or_side_length)
 179             , done(done)
 180             , spring_constant(spring_constant)
 181             , index(index)
 182             , distance(distance)
 183             , spring_strength(spring_strength)
 184             , partial_derivatives(partial_derivatives)
 185             {
 186             }
 187
 188             // Compute contribution of vertex i to the first partial
 189             // derivatives (dE/dx_m, dE/dy_m) (for vertex m)
 190             deriv_type compute_partial_derivative(
 191                 vertex_descriptor m, vertex_descriptor i)
 192             {
 193 #ifndef BOOST_NO_STDC_NAMESPACE
 194                 using std::sqrt;
 195 #endif // BOOST_NO_STDC_NAMESPACE
 196
 197                 deriv_type result;
 198                 if (i != m)
 199                 {
 200                     point_difference_type diff
 201                         = topology.difference(position[m], position[i]);
 202                     weight_type dist = topology.norm(diff);
 203                     result = spring_strength[get(index, m)][get(index, i)]
 204                         * (diff
 205                             - distance[get(index, m)][get(index, i)] / dist
 206                                 * diff);
 207                 }
 208
 209                 return result;
 210             }
 211
 212             // Compute partial derivatives dE/dx_m and dE/dy_m
 213             deriv_type compute_partial_derivatives(vertex_descriptor m)
 214             {
 215 #ifndef BOOST_NO_STDC_NAMESPACE
 216                 using std::sqrt;
 217 #endif // BOOST_NO_STDC_NAMESPACE
 218
 219                 deriv_type result;
 220
 221                 // TBD: looks like an accumulate to me
 222                 BGL_FORALL_VERTICES_T(i, g, Graph)
 223                 {
 224                     deriv_type deriv = compute_partial_derivative(m, i);
 225                     result += deriv;
 226                 }
 227
 228                 return result;
 229             }
 230
 231             // The actual Kamada-Kawai spring layout algorithm implementation
 232             bool run()
 233             {
 234 #ifndef BOOST_NO_STDC_NAMESPACE
 235                 using std::sqrt;
 236 #endif // BOOST_NO_STDC_NAMESPACE
 237
 238                 // Compute d_{ij} and place it in the distance matrix
 239                 if (!johnson_all_pairs_shortest_paths(
 240                         g, distance, index, weight, weight_type(0)))
 241                     return false;
 242
 243                 // Compute L based on side length (if needed), or retrieve L
 244                 weight_type edge_length = detail::graph::compute_edge_length(
 245                     g, distance, index, edge_or_side_length);
 246
 247                 // std::cerr << "edge_length = " << edge_length << std::endl;
 248
 249                 // Compute l_{ij} and k_{ij}
 250                 const weight_type K = spring_constant;
 251                 vertex_iterator ui, end;
 252                 for (ui = vertices(g).first, end = vertices(g).second;
 253                      ui != end; ++ui)
 254                 {
 255                     vertex_iterator vi = ui;
 256                     for (++vi; vi != end; ++vi)
 257                     {
 258                         weight_type dij
 259                             = distance[get(index, *ui)][get(index, *vi)];
 260                         if (dij == (std::numeric_limits< weight_type >::max)())
 261                             return false;
 262                         distance[get(index, *ui)][get(index, *vi)]
 263                             = edge_length * dij;
 264                         distance[get(index, *vi)][get(index, *ui)]
 265                             = edge_length * dij;
 266                         spring_strength[get(index, *ui)][get(index, *vi)]
 267                             = K / (dij * dij);
 268                         spring_strength[get(index, *vi)][get(index, *ui)]
 269                             = K / (dij * dij);
 270                     }
 271                 }
 272
 273                 // Compute Delta_i and find max
 274                 vertex_descriptor p = *vertices(g).first;
 275                 weight_type delta_p(0);
 276
 277                 for (ui = vertices(g).first, end = vertices(g).second;
 278                      ui != end; ++ui)
 279                 {
 280                     deriv_type deriv = compute_partial_derivatives(*ui);
 281                     put(partial_derivatives, *ui, deriv);
 282
 283                     weight_type delta = topology.norm(deriv);
 284
 285                     if (delta > delta_p)
 286                     {
 287                         p = *ui;
 288                         delta_p = delta;
 289                     }
 290                 }
 291
 292                 while (!done(delta_p, p, g, true))
 293                 {
 294                     // The contribution p makes to the partial derivatives of
 295                     // each vertex. Computing this (at O(n) cost) allows us to
 296                     // update the delta_i values in O(n) time instead of O(n^2)
 297                     // time.
 298                     std::vector< deriv_type > p_partials(num_vertices(g));
 299                     for (ui = vertices(g).first, end = vertices(g).second;
 300                          ui != end; ++ui)
 301                     {
 302                         vertex_descriptor i = *ui;
 303                         p_partials[get(index, i)]
 304                             = compute_partial_derivative(i, p);
 305                     }
 306
 307                     do
 308                     {
 309                         // For debugging, compute the energy value E
 310                         double E = 0.;
 311                         for (ui = vertices(g).first, end = vertices(g).second;
 312                              ui != end; ++ui)
 313                         {
 314                             vertex_iterator vi = ui;
 315                             for (++vi; vi != end; ++vi)
 316                             {
 317                                 double dist = topology.distance(
 318                                     position[*ui], position[*vi]);
 319                                 weight_type k_ij = spring_strength[get(
 320                                     index, *ui)][get(index, *vi)];
 321                                 weight_type l_ij = distance[get(index, *ui)]
 322                                                            [get(index, *vi)];
 323                                 E += .5 * k_ij * (dist - l_ij) * (dist - l_ij);
 324                             }
 325                         }
 326                         // std::cerr << "E = " << E << std::endl;
 327
 328                         // Compute the elements of the Jacobian
 329                         // From
 330                         // http://www.cs.panam.edu/~rfowler/papers/1994_kumar_fowler_A_Spring_UTPACSTR.pdf
 331                         // with the bugs fixed in the off-diagonal case
 332                         weight_type dE_d_d[Point::dimensions]
 333                                           [Point::dimensions];
 334                         for (std::size_t i = 0; i < Point::dimensions; ++i)
 335                             for (std::size_t j = 0; j < Point::dimensions; ++j)
 336                                 dE_d_d[i][j] = 0.;
 337                         for (ui = vertices(g).first, end = vertices(g).second;
 338                              ui != end; ++ui)
 339                         {
 340                             vertex_descriptor i = *ui;
 341                             if (i != p)
 342                             {
 343                                 point_difference_type diff
 344                                     = topology.difference(
 345                                         position[p], position[i]);
 346                                 weight_type dist = topology.norm(diff);
 347                                 weight_type dist_squared = dist * dist;
 348                                 weight_type inv_dist_cubed
 349                                     = 1. / (dist_squared * dist);
 350                                 weight_type k_mi = spring_strength[get(
 351                                     index, p)][get(index, i)];
 352                                 weight_type l_mi
 353                                     = distance[get(index, p)][get(index, i)];
 354                                 for (std::size_t i = 0; i < Point::dimensions;
 355                                      ++i)
 356                                 {
 357                                     for (std::size_t j = 0;
 358                                          j < Point::dimensions; ++j)
 359                                     {
 360                                         if (i == j)
 361                                         {
 362                                             dE_d_d[i][i] += k_mi
 363                                                 * (1
 364                                                     + (l_mi
 365                                                         * (diff[i] * diff[i]
 366                                                             - dist_squared)
 367                                                         * inv_dist_cubed));
 368                                         }
 369                                         else
 370                                         {
 371                                             dE_d_d[i][j] += k_mi * l_mi
 372                                                 * diff[i] * diff[j]
 373                                                 * inv_dist_cubed;
 374                                             // dE_d_d[i][j] += k_mi * l_mi *
 375                                             // sqrt(hypot(diff[i], diff[j])) *
 376                                             // inv_dist_cubed;
 377                                         }
 378                                     }
 379                                 }
 380                             }
 381                         }
 382
 383                         deriv_type dE_d = get(partial_derivatives, p);
 384
 385                         // Solve dE_d_d * delta = -dE_d to get delta
 386                         point_difference_type delta
 387                             = -linear_solver< Point::dimensions >::solve(
 388                                 dE_d_d, dE_d);
 389
 390                         // Move p by delta
 391                         position[p] = topology.adjust(position[p], delta);
 392
 393                         // Recompute partial derivatives and delta_p
 394                         deriv_type deriv = compute_partial_derivatives(p);
 395                         put(partial_derivatives, p, deriv);
 396
 397                         delta_p = topology.norm(deriv);
 398                     } while (!done(delta_p, p, g, false));
 399
 400                     // Select new p by updating each partial derivative and
 401                     // delta
 402                     vertex_descriptor old_p = p;
 403                     for (ui = vertices(g).first, end = vertices(g).second;
 404                          ui != end; ++ui)
 405                     {
 406                         deriv_type old_deriv_p = p_partials[get(index, *ui)];
 407                         deriv_type old_p_partial
 408                             = compute_partial_derivative(*ui, old_p);
 409                         deriv_type deriv = get(partial_derivatives, *ui);
 410
 411                         deriv += old_p_partial - old_deriv_p;
 412
 413                         put(partial_derivatives, *ui, deriv);
 414                         weight_type delta = topology.norm(deriv);
 415
 416                         if (delta > delta_p)
 417                         {
 418                             p = *ui;
 419                             delta_p = delta;
 420                         }
 421                     }
 422                 }
 423
 424                 return true;
 425             }
 426
 427             const Topology& topology;
 428             const Graph& g;
 429             PositionMap position;
 430             WeightMap weight;
 431             EdgeOrSideLength edge_or_side_length;
 432             Done done;
 433             weight_type spring_constant;
 434             VertexIndexMap index;
 435             DistanceMatrix distance;
 436             SpringStrengthMatrix spring_strength;
 437             PartialDerivativeMap partial_derivatives;
 438         };
 439     }
 440 } // end namespace detail::graph
 441
 442 /// States that the given quantity is an edge length.
 443 template < typename T >
 444 inline detail::graph::edge_or_side< true, T > edge_length(T x)
 445 {
 446     return detail::graph::edge_or_side< true, T >(x);
 447 }
 448
 449 /// States that the given quantity is a display area side length.
 450 template < typename T >
 451 inline detail::graph::edge_or_side< false, T > side_length(T x)
 452 {
 453     return detail::graph::edge_or_side< false, T >(x);
 454 }
 455
 456 /**
 457  * \brief Determines when to terminate layout of a particular graph based
 458  * on a given relative tolerance.
 459  */
 460 template < typename T = double > struct layout_tolerance
 461 {
 462     layout_tolerance(const T& tolerance = T(0.001))
 463     : tolerance(tolerance)
 464     , last_energy((std::numeric_limits< T >::max)())
 465     , last_local_energy((std::numeric_limits< T >::max)())
 466     {
 467     }
 468
 469     template < typename Graph >
 470     bool operator()(T delta_p,
 471         typename boost::graph_traits< Graph >::vertex_descriptor p,
 472         const Graph& g, bool global)
 473     {
 474         if (global)
 475         {
 476             if (last_energy == (std::numeric_limits< T >::max)())
 477             {
 478                 last_energy = delta_p;
 479                 return false;
 480             }
 481
 482             T diff = last_energy - delta_p;
 483             if (diff < T(0))
 484                 diff = -diff;
 485             bool done = (delta_p == T(0) || diff / last_energy < tolerance);
 486             last_energy = delta_p;
 487             return done;
 488         }
 489         else
 490         {
 491             if (last_local_energy == (std::numeric_limits< T >::max)())
 492             {
 493                 last_local_energy = delta_p;
 494                 return delta_p == T(0);
 495             }
 496
 497             T diff = last_local_energy - delta_p;
 498             bool done
 499                 = (delta_p == T(0) || (diff / last_local_energy) < tolerance);
 500             last_local_energy = delta_p;
 501             return done;
 502         }
 503     }
 504
 505 private:
 506     T tolerance;
 507     T last_energy;
 508     T last_local_energy;
 509 };
 510
 511 /** \brief Kamada-Kawai spring layout for undirected graphs.
 512  *
 513  * This algorithm performs graph layout (in two dimensions) for
 514  * connected, undirected graphs. It operates by relating the layout
 515  * of graphs to a dynamic spring system and minimizing the energy
 516  * within that system. The strength of a spring between two vertices
 517  * is inversely proportional to the square of the shortest distance
 518  * (in graph terms) between those two vertices. Essentially,
 519  * vertices that are closer in the graph-theoretic sense (i.e., by
 520  * following edges) will have stronger springs and will therefore be
 521  * placed closer together.
 522  *
 523  * Prior to invoking this algorithm, it is recommended that the
 524  * vertices be placed along the vertices of a regular n-sided
 525  * polygon.
 526  *
 527  * \param g (IN) must be a model of Vertex List Graph, Edge List
 528  * Graph, and Incidence Graph and must be undirected.
 529  *
 530  * \param position (OUT) must be a model of Lvalue Property Map,
 531  * where the value type is a class containing fields @c x and @c y
 532  * that will be set to the @c x and @c y coordinates of each vertex.
 533  *
 534  * \param weight (IN) must be a model of Readable Property Map,
 535  * which provides the weight of each edge in the graph @p g.
 536  *
 537  * \param topology (IN) must be a topology object (see topology.hpp),
 538  * which provides operations on points and differences between them.
 539  *
 540  * \param edge_or_side_length (IN) provides either the unit length
 541  * @c e of an edge in the layout or the length of a side @c s of the
 542  * display area, and must be either @c boost::edge_length(e) or @c
 543  * boost::side_length(s), respectively.
 544  *
 545  * \param done (IN) is a 4-argument function object that is passed
 546  * the current value of delta_p (i.e., the energy of vertex @p p),
 547  * the vertex @p p, the graph @p g, and a boolean flag indicating
 548  * whether @p delta_p is the maximum energy in the system (when @c
 549  * true) or the energy of the vertex being moved. Defaults to @c
 550  * layout_tolerance instantiated over the value type of the weight
 551  * map.
 552  *
 553  * \param spring_constant (IN) is the constant multiplied by each
 554  * spring's strength. Larger values create systems with more energy
 555  * that can take longer to stabilize; smaller values create systems
 556  * with less energy that stabilize quickly but do not necessarily
 557  * result in pleasing layouts. The default value is 1.
 558  *
 559  * \param index (IN) is a mapping from vertices to index values
 560  * between 0 and @c num_vertices(g). The default is @c
 561  * get(vertex_index,g).
 562  *
 563  * \param distance (UTIL/OUT) will be used to store the distance
 564  * from every vertex to every other vertex, which is computed in the
 565  * first stages of the algorithm. This value's type must be a model
 566  * of BasicMatrix with value type equal to the value type of the
 567  * weight map. The default is a vector of vectors.
 568  *
 569  * \param spring_strength (UTIL/OUT) will be used to store the
 570  * strength of the spring between every pair of vertices. This
 571  * value's type must be a model of BasicMatrix with value type equal
 572  * to the value type of the weight map. The default is a vector of
 573  * vectors.
 574  *
 575  * \param partial_derivatives (UTIL) will be used to store the
 576  * partial derivates of each vertex with respect to the @c x and @c
 577  * y coordinates. This must be a Read/Write Property Map whose value
 578  * type is a pair with both types equivalent to the value type of
 579  * the weight map. The default is an iterator property map.
 580  *
 581  * \returns @c true if layout was successful or @c false if a
 582  * negative weight cycle was detected.
 583  */
 584 template < typename Topology, typename Graph, typename PositionMap,
 585     typename WeightMap, typename T, bool EdgeOrSideLength, typename Done,
 586     typename VertexIndexMap, typename DistanceMatrix,
 587     typename SpringStrengthMatrix, typename PartialDerivativeMap >
 588 bool kamada_kawai_spring_layout(const Graph& g, PositionMap position,
 589     WeightMap weight, const Topology& topology,
 590     detail::graph::edge_or_side< EdgeOrSideLength, T > edge_or_side_length,
 591     Done done,
 592     typename property_traits< WeightMap >::value_type spring_constant,
 593     VertexIndexMap index, DistanceMatrix distance,
 594     SpringStrengthMatrix spring_strength,
 595     PartialDerivativeMap partial_derivatives)
 596 {
 597     BOOST_STATIC_ASSERT(
 598         (is_convertible< typename graph_traits< Graph >::directed_category*,
 599             undirected_tag* >::value));
 600
 601     detail::graph::kamada_kawai_spring_layout_impl< Topology, Graph,
 602         PositionMap, WeightMap,
 603         detail::graph::edge_or_side< EdgeOrSideLength, T >, Done,
 604         VertexIndexMap, DistanceMatrix, SpringStrengthMatrix,
 605         PartialDerivativeMap >
 606         alg(topology, g, position, weight, edge_or_side_length, done,
 607             spring_constant, index, distance, spring_strength,
 608             partial_derivatives);
 609     return alg.run();
 610 }
 611
 612 /**
 613  * \overload
 614  */
 615 template < typename Topology, typename Graph, typename PositionMap,
 616     typename WeightMap, typename T, bool EdgeOrSideLength, typename Done,
 617     typename VertexIndexMap >
 618 bool kamada_kawai_spring_layout(const Graph& g, PositionMap position,
 619     WeightMap weight, const Topology& topology,
 620     detail::graph::edge_or_side< EdgeOrSideLength, T > edge_or_side_length,
 621     Done done,
 622     typename property_traits< WeightMap >::value_type spring_constant,
 623     VertexIndexMap index)
 624 {
 625     typedef typename property_traits< WeightMap >::value_type weight_type;
 626
 627     typename graph_traits< Graph >::vertices_size_type n = num_vertices(g);
 628     typedef std::vector< weight_type > weight_vec;
 629
 630     std::vector< weight_vec > distance(n, weight_vec(n));
 631     std::vector< weight_vec > spring_strength(n, weight_vec(n));
 632     std::vector< typename Topology::point_difference_type > partial_derivatives(
 633         n);
 634
 635     return kamada_kawai_spring_layout(g, position, weight, topology,
 636         edge_or_side_length, done, spring_constant, index, distance.begin(),
 637         spring_strength.begin(),
 638         make_iterator_property_map(partial_derivatives.begin(), index,
 639             typename Topology::point_difference_type()));
 640 }
 641
 642 /**
 643  * \overload
 644  */
 645 template < typename Topology, typename Graph, typename PositionMap,
 646     typename WeightMap, typename T, bool EdgeOrSideLength, typename Done >
 647 bool kamada_kawai_spring_layout(const Graph& g, PositionMap position,
 648     WeightMap weight, const Topology& topology,
 649     detail::graph::edge_or_side< EdgeOrSideLength, T > edge_or_side_length,
 650     Done done,
 651     typename property_traits< WeightMap >::value_type spring_constant)
 652 {
 653     return kamada_kawai_spring_layout(g, position, weight, topology,
 654         edge_or_side_length, done, spring_constant, get(vertex_index, g));
 655 }
 656
 657 /**
 658  * \overload
 659  */
 660 template < typename Topology, typename Graph, typename PositionMap,
 661     typename WeightMap, typename T, bool EdgeOrSideLength, typename Done >
 662 bool kamada_kawai_spring_layout(const Graph& g, PositionMap position,
 663     WeightMap weight, const Topology& topology,
 664     detail::graph::edge_or_side< EdgeOrSideLength, T > edge_or_side_length,
 665     Done done)
 666 {
 667     typedef typename property_traits< WeightMap >::value_type weight_type;
 668     return kamada_kawai_spring_layout(g, position, weight, topology,
 669         edge_or_side_length, done, weight_type(1));
 670 }
 671
 672 /**
 673  * \overload
 674  */
 675 template < typename Topology, typename Graph, typename PositionMap,
 676     typename WeightMap, typename T, bool EdgeOrSideLength >
 677 bool kamada_kawai_spring_layout(const Graph& g, PositionMap position,
 678     WeightMap weight, const Topology& topology,
 679     detail::graph::edge_or_side< EdgeOrSideLength, T > edge_or_side_length)
 680 {
 681     typedef typename property_traits< WeightMap >::value_type weight_type;
 682     return kamada_kawai_spring_layout(g, position, weight, topology,
 683         edge_or_side_length, layout_tolerance< weight_type >(),
 684         weight_type(1.0), get(vertex_index, g));
 685 }
 686 } // end namespace boost
 687
 688 #endif // BOOST_GRAPH_KAMADA_KAWAI_SPRING_LAYOUT_HPP