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