1 /* -*- mode: C -*- */
2 /* vim:set ts=4 sw=4 sts=4 et: */
3 /*
4 IGraph R package.
5 Copyright (C) 2014 Gabor Csardi <csardi.gabor@gmail.com>
6 334 Harvard street, Cambridge, MA 02139 USA
7
8 This program is free software; you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation; either version 2 of the License, or
11 (at your option) any later version.
12
13 This program is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with this program; if not, write to the Free Software
20 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
21 02110-1301 USA
22
23 */
24
25 #include "igraph_layout.h"
26 #include "igraph_random.h"
27 #include "igraph_interface.h"
28 #include "igraph_components.h"
29 #include "igraph_types_internal.h"
30
igraph_layout_i_fr(const igraph_t * graph,igraph_matrix_t * res,igraph_bool_t use_seed,igraph_integer_t niter,igraph_real_t start_temp,const igraph_vector_t * weight,const igraph_vector_t * minx,const igraph_vector_t * maxx,const igraph_vector_t * miny,const igraph_vector_t * maxy)31 static int igraph_layout_i_fr(const igraph_t *graph,
32 igraph_matrix_t *res,
33 igraph_bool_t use_seed,
34 igraph_integer_t niter,
35 igraph_real_t start_temp,
36 const igraph_vector_t *weight,
37 const igraph_vector_t *minx,
38 const igraph_vector_t *maxx,
39 const igraph_vector_t *miny,
40 const igraph_vector_t *maxy) {
41
42 igraph_integer_t no_nodes = igraph_vcount(graph);
43 igraph_integer_t no_edges = igraph_ecount(graph);
44 igraph_integer_t i;
45 igraph_vector_float_t dispx, dispy;
46 igraph_real_t temp = start_temp;
47 igraph_real_t difftemp = start_temp / niter;
48 float width = sqrtf(no_nodes), height = width;
49 igraph_bool_t conn = 1;
50 float C;
51
52 igraph_is_connected(graph, &conn, IGRAPH_WEAK);
53 if (!conn) {
54 C = no_nodes * sqrtf(no_nodes);
55 }
56
57 RNG_BEGIN();
58
59 if (!use_seed) {
60 IGRAPH_CHECK(igraph_matrix_resize(res, no_nodes, 2));
61 for (i = 0; i < no_nodes; i++) {
62 igraph_real_t x1 = minx ? VECTOR(*minx)[i] : -width / 2;
63 igraph_real_t x2 = maxx ? VECTOR(*maxx)[i] : width / 2;
64 igraph_real_t y1 = miny ? VECTOR(*miny)[i] : -height / 2;
65 igraph_real_t y2 = maxy ? VECTOR(*maxy)[i] : height / 2;
66 if (!igraph_finite(x1)) {
67 x1 = -sqrt(no_nodes) / 2;
68 }
69 if (!igraph_finite(x2)) {
70 x2 = sqrt(no_nodes) / 2;
71 }
72 if (!igraph_finite(y1)) {
73 y1 = -sqrt(no_nodes) / 2;
74 }
75 if (!igraph_finite(y2)) {
76 y2 = sqrt(no_nodes) / 2;
77 }
78 MATRIX(*res, i, 0) = RNG_UNIF(x1, x2);
79 MATRIX(*res, i, 1) = RNG_UNIF(y1, y2);
80 }
81 }
82
83 IGRAPH_CHECK(igraph_vector_float_init(&dispx, no_nodes));
84 IGRAPH_FINALLY(igraph_vector_float_destroy, &dispx);
85 IGRAPH_CHECK(igraph_vector_float_init(&dispy, no_nodes));
86 IGRAPH_FINALLY(igraph_vector_float_destroy, &dispy);
87
88 for (i = 0; i < niter; i++) {
89 igraph_integer_t v, u, e;
90
91 /* calculate repulsive forces, we have a special version
92 for unconnected graphs */
93 igraph_vector_float_null(&dispx);
94 igraph_vector_float_null(&dispy);
95 if (conn) {
96 for (v = 0; v < no_nodes; v++) {
97 for (u = v + 1; u < no_nodes; u++) {
98 float dx = MATRIX(*res, v, 0) - MATRIX(*res, u, 0);
99 float dy = MATRIX(*res, v, 1) - MATRIX(*res, u, 1);
100 float dlen = dx * dx + dy * dy;
101
102 if (dlen == 0) {
103 dx = RNG_UNIF01() * 1e-9;
104 dy = RNG_UNIF01() * 1e-9;
105 dlen = dx * dx + dy * dy;
106 }
107
108 VECTOR(dispx)[v] += dx / dlen;
109 VECTOR(dispy)[v] += dy / dlen;
110 VECTOR(dispx)[u] -= dx / dlen;
111 VECTOR(dispy)[u] -= dy / dlen;
112 }
113 }
114 } else {
115 for (v = 0; v < no_nodes; v++) {
116 for (u = v + 1; u < no_nodes; u++) {
117 float dx = MATRIX(*res, v, 0) - MATRIX(*res, u, 0);
118 float dy = MATRIX(*res, v, 1) - MATRIX(*res, u, 1);
119 float dlen, rdlen;
120
121 dlen = dx * dx + dy * dy;
122 if (dlen == 0) {
123 dx = RNG_UNIF(0, 1e-6);
124 dy = RNG_UNIF(0, 1e-6);
125 dlen = dx * dx + dy * dy;
126 }
127
128 rdlen = sqrt(dlen);
129
130 VECTOR(dispx)[v] += dx * (C - dlen * rdlen) / (dlen * C);
131 VECTOR(dispy)[v] += dy * (C - dlen * rdlen) / (dlen * C);
132 VECTOR(dispx)[u] -= dx * (C - dlen * rdlen) / (dlen * C);
133 VECTOR(dispy)[u] -= dy * (C - dlen * rdlen) / (dlen * C);
134 }
135 }
136 }
137
138 /* calculate attractive forces */
139 for (e = 0; e < no_edges; e++) {
140 /* each edges is an ordered pair of vertices v and u */
141 igraph_integer_t v = IGRAPH_FROM(graph, e);
142 igraph_integer_t u = IGRAPH_TO(graph, e);
143 igraph_real_t dx = MATRIX(*res, v, 0) - MATRIX(*res, u, 0);
144 igraph_real_t dy = MATRIX(*res, v, 1) - MATRIX(*res, u, 1);
145 igraph_real_t w = weight ? VECTOR(*weight)[e] : 1.0;
146 igraph_real_t dlen = sqrt(dx * dx + dy * dy) * w;
147 VECTOR(dispx)[v] -= (dx * dlen);
148 VECTOR(dispy)[v] -= (dy * dlen);
149 VECTOR(dispx)[u] += (dx * dlen);
150 VECTOR(dispy)[u] += (dy * dlen);
151 }
152
153 /* limit max displacement to temperature t and prevent from
154 displacement outside frame */
155 for (v = 0; v < no_nodes; v++) {
156 igraph_real_t dx = VECTOR(dispx)[v] + RNG_UNIF01() * 1e-9;
157 igraph_real_t dy = VECTOR(dispy)[v] + RNG_UNIF01() * 1e-9;
158 igraph_real_t displen = sqrt(dx * dx + dy * dy);
159 igraph_real_t mx = fabs(dx) < temp ? dx : temp;
160 igraph_real_t my = fabs(dy) < temp ? dy : temp;
161 if (displen > 0) {
162 MATRIX(*res, v, 0) += (dx / displen) * mx;
163 MATRIX(*res, v, 1) += (dy / displen) * my;
164 }
165 if (minx && MATRIX(*res, v, 0) < VECTOR(*minx)[v]) {
166 MATRIX(*res, v, 0) = VECTOR(*minx)[v];
167 }
168 if (maxx && MATRIX(*res, v, 0) > VECTOR(*maxx)[v]) {
169 MATRIX(*res, v, 0) = VECTOR(*maxx)[v];
170 }
171 if (miny && MATRIX(*res, v, 1) < VECTOR(*miny)[v]) {
172 MATRIX(*res, v, 1) = VECTOR(*miny)[v];
173 }
174 if (maxy && MATRIX(*res, v, 1) > VECTOR(*maxy)[v]) {
175 MATRIX(*res, v, 1) = VECTOR(*maxy)[v];
176 }
177 }
178
179 temp -= difftemp;
180 }
181
182 RNG_END();
183
184 igraph_vector_float_destroy(&dispx);
185 igraph_vector_float_destroy(&dispy);
186 IGRAPH_FINALLY_CLEAN(2);
187
188 return 0;
189 }
190
igraph_layout_i_grid_fr(const igraph_t * graph,igraph_matrix_t * res,igraph_bool_t use_seed,igraph_integer_t niter,igraph_real_t start_temp,const igraph_vector_t * weight,const igraph_vector_t * minx,const igraph_vector_t * maxx,const igraph_vector_t * miny,const igraph_vector_t * maxy)191 static int igraph_layout_i_grid_fr(
192 const igraph_t *graph,
193 igraph_matrix_t *res, igraph_bool_t use_seed,
194 igraph_integer_t niter, igraph_real_t start_temp,
195 const igraph_vector_t *weight, const igraph_vector_t *minx,
196 const igraph_vector_t *maxx, const igraph_vector_t *miny,
197 const igraph_vector_t *maxy) {
198
199 igraph_integer_t no_nodes = igraph_vcount(graph);
200 igraph_integer_t no_edges = igraph_ecount(graph);
201 float width = sqrtf(no_nodes), height = width;
202 igraph_2dgrid_t grid;
203 igraph_vector_float_t dispx, dispy;
204 igraph_real_t temp = start_temp;
205 igraph_real_t difftemp = start_temp / niter;
206 igraph_2dgrid_iterator_t vidit;
207 igraph_integer_t i;
208 const float cellsize = 2.0;
209
210 RNG_BEGIN();
211
212 if (!use_seed) {
213 IGRAPH_CHECK(igraph_matrix_resize(res, no_nodes, 2));
214 for (i = 0; i < no_nodes; i++) {
215 igraph_real_t x1 = minx ? VECTOR(*minx)[i] : -width / 2;
216 igraph_real_t x2 = maxx ? VECTOR(*maxx)[i] : width / 2;
217 igraph_real_t y1 = miny ? VECTOR(*miny)[i] : -height / 2;
218 igraph_real_t y2 = maxy ? VECTOR(*maxy)[i] : height / 2;
219 if (!igraph_finite(x1)) {
220 x1 = -sqrt(no_nodes) / 2;
221 }
222 if (!igraph_finite(x2)) {
223 x2 = sqrt(no_nodes) / 2;
224 }
225 if (!igraph_finite(y1)) {
226 y1 = -sqrt(no_nodes) / 2;
227 }
228 if (!igraph_finite(y2)) {
229 y2 = sqrt(no_nodes) / 2;
230 }
231 MATRIX(*res, i, 0) = RNG_UNIF(x1, x2);
232 MATRIX(*res, i, 1) = RNG_UNIF(y1, y2);
233 }
234 }
235
236 /* make grid */
237 IGRAPH_CHECK(igraph_2dgrid_init(&grid, res, -width / 2, width / 2, cellsize,
238 -height / 2, height / 2, cellsize));
239 IGRAPH_FINALLY(igraph_2dgrid_destroy, &grid);
240
241 /* place vertices on grid */
242 for (i = 0; i < no_nodes; i++) {
243 igraph_2dgrid_add2(&grid, i);
244 }
245
246 IGRAPH_CHECK(igraph_vector_float_init(&dispx, no_nodes));
247 IGRAPH_FINALLY(igraph_vector_float_destroy, &dispx);
248 IGRAPH_CHECK(igraph_vector_float_init(&dispy, no_nodes));
249 IGRAPH_FINALLY(igraph_vector_float_destroy, &dispy);
250
251 for (i = 0; i < niter; i++) {
252 igraph_integer_t v, u, e;
253
254 igraph_vector_float_null(&dispx);
255 igraph_vector_float_null(&dispy);
256
257 /* repulsion */
258 igraph_2dgrid_reset(&grid, &vidit);
259 while ( (v = igraph_2dgrid_next(&grid, &vidit) - 1) != -1) {
260 while ( (u = igraph_2dgrid_next_nei(&grid, &vidit) - 1) != -1) {
261 float dx = MATRIX(*res, v, 0) - MATRIX(*res, u, 0);
262 float dy = MATRIX(*res, v, 1) - MATRIX(*res, u, 1);
263 float dlen = dx * dx + dy * dy;
264 if (dlen < cellsize * cellsize) {
265 VECTOR(dispx)[v] += dx / dlen;
266 VECTOR(dispy)[v] += dy / dlen;
267 VECTOR(dispx)[u] -= dx / dlen;
268 VECTOR(dispy)[u] -= dy / dlen;
269 }
270 }
271 }
272
273 /* attraction */
274 for (e = 0; e < no_edges; e++) {
275 igraph_integer_t v = IGRAPH_FROM(graph, e);
276 igraph_integer_t u = IGRAPH_TO(graph, e);
277 igraph_real_t dx = MATRIX(*res, v, 0) - MATRIX(*res, u, 0);
278 igraph_real_t dy = MATRIX(*res, v, 1) - MATRIX(*res, u, 1);
279 igraph_real_t w = weight ? VECTOR(*weight)[e] : 1.0;
280 igraph_real_t dlen = sqrt(dx * dx + dy * dy) * w;
281 VECTOR(dispx)[v] -= (dx * dlen);
282 VECTOR(dispy)[v] -= (dy * dlen);
283 VECTOR(dispx)[u] += (dx * dlen);
284 VECTOR(dispy)[u] += (dy * dlen);
285 }
286
287 /* update */
288 for (v = 0; v < no_nodes; v++) {
289 igraph_real_t dx = VECTOR(dispx)[v] + RNG_UNIF01() * 1e-9;
290 igraph_real_t dy = VECTOR(dispy)[v] + RNG_UNIF01() * 1e-9;
291 igraph_real_t displen = sqrt(dx * dx + dy * dy);
292 igraph_real_t mx = fabs(dx) < temp ? dx : temp;
293 igraph_real_t my = fabs(dy) < temp ? dy : temp;
294 if (displen > 0) {
295 MATRIX(*res, v, 0) += (dx / displen) * mx;
296 MATRIX(*res, v, 1) += (dy / displen) * my;
297 }
298 if (minx && MATRIX(*res, v, 0) < VECTOR(*minx)[v]) {
299 MATRIX(*res, v, 0) = VECTOR(*minx)[v];
300 }
301 if (maxx && MATRIX(*res, v, 0) > VECTOR(*maxx)[v]) {
302 MATRIX(*res, v, 0) = VECTOR(*maxx)[v];
303 }
304 if (miny && MATRIX(*res, v, 1) < VECTOR(*miny)[v]) {
305 MATRIX(*res, v, 1) = VECTOR(*miny)[v];
306 }
307 if (maxy && MATRIX(*res, v, 1) > VECTOR(*maxy)[v]) {
308 MATRIX(*res, v, 1) = VECTOR(*maxy)[v];
309 }
310 }
311
312 temp -= difftemp;
313 }
314
315 igraph_vector_float_destroy(&dispx);
316 igraph_vector_float_destroy(&dispy);
317 igraph_2dgrid_destroy(&grid);
318 IGRAPH_FINALLY_CLEAN(3);
319 return 0;
320 }
321
322 /**
323 * \ingroup layout
324 * \function igraph_layout_fruchterman_reingold
325 * \brief Places the vertices on a plane according to the Fruchterman-Reingold algorithm.
326 *
327 * </para><para>
328 * This is a force-directed layout, see Fruchterman, T.M.J. and
329 * Reingold, E.M.: Graph Drawing by Force-directed Placement.
330 * Software -- Practice and Experience, 21/11, 1129--1164,
331 * 1991.
332 * \param graph Pointer to an initialized graph object.
333 * \param res Pointer to an initialized matrix object. This will
334 * contain the result and will be resized as needed.
335 * \param use_seed Logical, if true the supplied values in the
336 * \p res argument are used as an initial layout, if
337 * false a random initial layout is used.
338 * \param niter The number of iterations to do. A reasonable
339 * default value is 500.
340 * \param start_temp Start temperature. This is the maximum amount
341 * of movement alloved along one axis, within one step, for a
342 * vertex. Currently it is decreased linearly to zero during
343 * the iteration.
344 * \param grid Whether to use the (fast but less accurate) grid based
345 * version of the algorithm. Possible values: \c
346 * IGRAPH_LAYOUT_GRID, \c IGRAPH_LAYOUT_NOGRID, \c
347 * IGRAPH_LAYOUT_AUTOGRID. The last one uses the grid based
348 * version only for large graphs, currently the ones with
349 * more than 1000 vertices.
350 * \param weight Pointer to a vector containing edge weights,
351 * the attraction along the edges will be multiplied by these.
352 * It will be ignored if it is a null-pointer.
353 * \param minx Pointer to a vector, or a \c NULL pointer. If not a
354 * \c NULL pointer then the vector gives the minimum
355 * \quote x \endquote coordinate for every vertex.
356 * \param maxx Same as \p minx, but the maximum \quote x \endquote
357 * coordinates.
358 * \param miny Pointer to a vector, or a \c NULL pointer. If not a
359 * \c NULL pointer then the vector gives the minimum
360 * \quote y \endquote coordinate for every vertex.
361 * \param maxy Same as \p miny, but the maximum \quote y \endquote
362 * coordinates.
363 * \return Error code.
364 *
365 * Time complexity: O(|V|^2) in each
366 * iteration, |V| is the number of
367 * vertices in the graph.
368 */
369
igraph_layout_fruchterman_reingold(const igraph_t * graph,igraph_matrix_t * res,igraph_bool_t use_seed,igraph_integer_t niter,igraph_real_t start_temp,igraph_layout_grid_t grid,const igraph_vector_t * weight,const igraph_vector_t * minx,const igraph_vector_t * maxx,const igraph_vector_t * miny,const igraph_vector_t * maxy)370 int igraph_layout_fruchterman_reingold(const igraph_t *graph,
371 igraph_matrix_t *res,
372 igraph_bool_t use_seed,
373 igraph_integer_t niter,
374 igraph_real_t start_temp,
375 igraph_layout_grid_t grid,
376 const igraph_vector_t *weight,
377 const igraph_vector_t *minx,
378 const igraph_vector_t *maxx,
379 const igraph_vector_t *miny,
380 const igraph_vector_t *maxy) {
381
382 igraph_integer_t no_nodes = igraph_vcount(graph);
383
384 if (niter < 0) {
385 IGRAPH_ERROR("Number of iterations must be non-negative in "
386 "Fruchterman-Reingold layout", IGRAPH_EINVAL);
387 }
388
389 if (use_seed && (igraph_matrix_nrow(res) != no_nodes ||
390 igraph_matrix_ncol(res) != 2)) {
391 IGRAPH_ERROR("Invalid start position matrix size in "
392 "Fruchterman-Reingold layout", IGRAPH_EINVAL);
393 }
394
395 if (weight && igraph_vector_size(weight) != igraph_ecount(graph)) {
396 IGRAPH_ERROR("Invalid weight vector length", IGRAPH_EINVAL);
397 }
398
399 if (minx && igraph_vector_size(minx) != no_nodes) {
400 IGRAPH_ERROR("Invalid minx vector length", IGRAPH_EINVAL);
401 }
402 if (maxx && igraph_vector_size(maxx) != no_nodes) {
403 IGRAPH_ERROR("Invalid maxx vector length", IGRAPH_EINVAL);
404 }
405 if (minx && maxx && !igraph_vector_all_le(minx, maxx)) {
406 IGRAPH_ERROR("minx must not be greater than maxx", IGRAPH_EINVAL);
407 }
408 if (miny && igraph_vector_size(miny) != no_nodes) {
409 IGRAPH_ERROR("Invalid miny vector length", IGRAPH_EINVAL);
410 }
411 if (maxy && igraph_vector_size(maxy) != no_nodes) {
412 IGRAPH_ERROR("Invalid maxy vector length", IGRAPH_EINVAL);
413 }
414 if (miny && maxy && !igraph_vector_all_le(miny, maxy)) {
415 IGRAPH_ERROR("miny must not be greater than maxy", IGRAPH_EINVAL);
416 }
417
418 if (grid == IGRAPH_LAYOUT_AUTOGRID) {
419 if (no_nodes > 1000) {
420 grid = IGRAPH_LAYOUT_GRID;
421 } else {
422 grid = IGRAPH_LAYOUT_NOGRID;
423 }
424 }
425
426 if (grid == IGRAPH_LAYOUT_GRID) {
427 return igraph_layout_i_grid_fr(graph, res, use_seed, niter, start_temp,
428 weight, minx, maxx, miny, maxy);
429 } else {
430 return igraph_layout_i_fr(graph, res, use_seed, niter, start_temp,
431 weight, minx, maxx, miny, maxy);
432 }
433 }
434
435 /**
436 * \function igraph_layout_fruchterman_reingold_3d
437 * \brief 3D Fruchterman-Reingold algorithm.
438 *
439 * This is the 3D version of the force based
440 * Fruchterman-Reingold layout (see \ref
441 * igraph_layout_fruchterman_reingold for the 2D version
442 *
443 * \param graph Pointer to an initialized graph object.
444 * \param res Pointer to an initialized matrix object. This will
445 * contain the result and will be resized as needed.
446 * \param use_seed Logical, if true the supplied values in the
447 * \p res argument are used as an initial layout, if
448 * false a random initial layout is used.
449 * \param niter The number of iterations to do. A reasonable
450 * default value is 500.
451 * \param start_temp Start temperature. This is the maximum amount
452 * of movement alloved along one axis, within one step, for a
453 * vertex. Currently it is decreased linearly to zero during
454 * the iteration.
455 * \param weight Pointer to a vector containing edge weights,
456 * the attraction along the edges will be multiplied by these.
457 * It will be ignored if it is a null-pointer.
458 * \param minx Pointer to a vector, or a \c NULL pointer. If not a
459 * \c NULL pointer then the vector gives the minimum
460 * \quote x \endquote coordinate for every vertex.
461 * \param maxx Same as \p minx, but the maximum \quote x \endquote
462 * coordinates.
463 * \param miny Pointer to a vector, or a \c NULL pointer. If not a
464 * \c NULL pointer then the vector gives the minimum
465 * \quote y \endquote coordinate for every vertex.
466 * \param maxy Same as \p miny, but the maximum \quote y \endquote
467 * coordinates.
468 * \param minz Pointer to a vector, or a \c NULL pointer. If not a
469 * \c NULL pointer then the vector gives the minimum
470 * \quote z \endquote coordinate for every vertex.
471 * \param maxz Same as \p minz, but the maximum \quote z \endquote
472 * coordinates.
473 * \return Error code.
474 *
475 * Added in version 0.2.</para><para>
476 *
477 * Time complexity: O(|V|^2) in each
478 * iteration, |V| is the number of
479 * vertices in the graph.
480 *
481 */
482
igraph_layout_fruchterman_reingold_3d(const igraph_t * graph,igraph_matrix_t * res,igraph_bool_t use_seed,igraph_integer_t niter,igraph_real_t start_temp,const igraph_vector_t * weight,const igraph_vector_t * minx,const igraph_vector_t * maxx,const igraph_vector_t * miny,const igraph_vector_t * maxy,const igraph_vector_t * minz,const igraph_vector_t * maxz)483 int igraph_layout_fruchterman_reingold_3d(const igraph_t *graph,
484 igraph_matrix_t *res,
485 igraph_bool_t use_seed,
486 igraph_integer_t niter,
487 igraph_real_t start_temp,
488 const igraph_vector_t *weight,
489 const igraph_vector_t *minx,
490 const igraph_vector_t *maxx,
491 const igraph_vector_t *miny,
492 const igraph_vector_t *maxy,
493 const igraph_vector_t *minz,
494 const igraph_vector_t *maxz) {
495
496 igraph_integer_t no_nodes = igraph_vcount(graph);
497 igraph_integer_t no_edges = igraph_ecount(graph);
498 igraph_integer_t i;
499 igraph_vector_float_t dispx, dispy, dispz;
500 igraph_real_t temp = start_temp;
501 igraph_real_t difftemp = start_temp / niter;
502 float width = sqrtf(no_nodes), height = width, depth = width;
503 igraph_bool_t conn = 1;
504 float C;
505
506 if (niter < 0) {
507 IGRAPH_ERROR("Number of iterations must be non-negative in "
508 "Fruchterman-Reingold layout", IGRAPH_EINVAL);
509 }
510
511 if (use_seed && (igraph_matrix_nrow(res) != no_nodes ||
512 igraph_matrix_ncol(res) != 3)) {
513 IGRAPH_ERROR("Invalid start position matrix size in "
514 "Fruchterman-Reingold layout", IGRAPH_EINVAL);
515 }
516
517 if (weight && igraph_vector_size(weight) != igraph_ecount(graph)) {
518 IGRAPH_ERROR("Invalid weight vector length", IGRAPH_EINVAL);
519 }
520
521 if (minx && igraph_vector_size(minx) != no_nodes) {
522 IGRAPH_ERROR("Invalid minx vector length", IGRAPH_EINVAL);
523 }
524 if (maxx && igraph_vector_size(maxx) != no_nodes) {
525 IGRAPH_ERROR("Invalid maxx vector length", IGRAPH_EINVAL);
526 }
527 if (minx && maxx && !igraph_vector_all_le(minx, maxx)) {
528 IGRAPH_ERROR("minx must not be greater than maxx", IGRAPH_EINVAL);
529 }
530 if (miny && igraph_vector_size(miny) != no_nodes) {
531 IGRAPH_ERROR("Invalid miny vector length", IGRAPH_EINVAL);
532 }
533 if (maxy && igraph_vector_size(maxy) != no_nodes) {
534 IGRAPH_ERROR("Invalid maxy vector length", IGRAPH_EINVAL);
535 }
536 if (miny && maxy && !igraph_vector_all_le(miny, maxy)) {
537 IGRAPH_ERROR("miny must not be greater than maxy", IGRAPH_EINVAL);
538 }
539 if (minz && igraph_vector_size(minz) != no_nodes) {
540 IGRAPH_ERROR("Invalid minz vector length", IGRAPH_EINVAL);
541 }
542 if (maxz && igraph_vector_size(maxz) != no_nodes) {
543 IGRAPH_ERROR("Invalid maxz vector length", IGRAPH_EINVAL);
544 }
545 if (minz && maxz && !igraph_vector_all_le(minz, maxz)) {
546 IGRAPH_ERROR("minz must not be greater than maxz", IGRAPH_EINVAL);
547 }
548
549 igraph_is_connected(graph, &conn, IGRAPH_WEAK);
550 if (!conn) {
551 C = no_nodes * sqrtf(no_nodes);
552 }
553
554 RNG_BEGIN();
555
556 if (!use_seed) {
557 IGRAPH_CHECK(igraph_matrix_resize(res, no_nodes, 3));
558 for (i = 0; i < no_nodes; i++) {
559 igraph_real_t x1 = minx ? VECTOR(*minx)[i] : -width / 2;
560 igraph_real_t x2 = maxx ? VECTOR(*maxx)[i] : width / 2;
561 igraph_real_t y1 = miny ? VECTOR(*miny)[i] : -height / 2;
562 igraph_real_t y2 = maxy ? VECTOR(*maxy)[i] : height / 2;
563 igraph_real_t z1 = minz ? VECTOR(*minz)[i] : -depth / 2;
564 igraph_real_t z2 = maxz ? VECTOR(*maxz)[i] : depth / 2;
565 MATRIX(*res, i, 0) = RNG_UNIF(x1, x2);
566 MATRIX(*res, i, 1) = RNG_UNIF(y1, y2);
567 MATRIX(*res, i, 2) = RNG_UNIF(z1, z2);
568 }
569 }
570
571 IGRAPH_CHECK(igraph_vector_float_init(&dispx, no_nodes));
572 IGRAPH_FINALLY(igraph_vector_float_destroy, &dispx);
573 IGRAPH_CHECK(igraph_vector_float_init(&dispy, no_nodes));
574 IGRAPH_FINALLY(igraph_vector_float_destroy, &dispy);
575 IGRAPH_CHECK(igraph_vector_float_init(&dispz, no_nodes));
576 IGRAPH_FINALLY(igraph_vector_float_destroy, &dispz);
577
578 for (i = 0; i < niter; i++) {
579 igraph_integer_t v, u, e;
580
581 /* calculate repulsive forces, we have a special version
582 for unconnected graphs */
583 igraph_vector_float_null(&dispx);
584 igraph_vector_float_null(&dispy);
585 igraph_vector_float_null(&dispz);
586 if (conn) {
587 for (v = 0; v < no_nodes; v++) {
588 for (u = v + 1; u < no_nodes; u++) {
589 float dx = MATRIX(*res, v, 0) - MATRIX(*res, u, 0);
590 float dy = MATRIX(*res, v, 1) - MATRIX(*res, u, 1);
591 float dz = MATRIX(*res, v, 2) - MATRIX(*res, u, 2);
592 float dlen = dx * dx + dy * dy + dz * dz;
593
594 if (dlen == 0) {
595 dx = RNG_UNIF01() * 1e-9;
596 dy = RNG_UNIF01() * 1e-9;
597 dz = RNG_UNIF01() * 1e-9;
598 dlen = dx * dx + dy * dy + dz * dz;
599 }
600
601 VECTOR(dispx)[v] += dx / dlen;
602 VECTOR(dispy)[v] += dy / dlen;
603 VECTOR(dispz)[v] += dz / dlen;
604 VECTOR(dispx)[u] -= dx / dlen;
605 VECTOR(dispy)[u] -= dy / dlen;
606 VECTOR(dispz)[u] -= dz / dlen;
607 }
608 }
609 } else {
610 for (v = 0; v < no_nodes; v++) {
611 for (u = v + 1; u < no_nodes; u++) {
612 float dx = MATRIX(*res, v, 0) - MATRIX(*res, u, 0);
613 float dy = MATRIX(*res, v, 1) - MATRIX(*res, u, 1);
614 float dz = MATRIX(*res, v, 2) - MATRIX(*res, u, 2);
615 float dlen, rdlen;
616
617 dlen = dx * dx + dy * dy + dz * dz;
618 if (dlen == 0) {
619 dx = RNG_UNIF01() * 1e-9;
620 dy = RNG_UNIF01() * 1e-9;
621 dz = RNG_UNIF01() * 1e-9;
622 dlen = dx * dx + dy * dy + dz * dz;
623 }
624
625 rdlen = sqrt(dlen);
626
627 VECTOR(dispx)[v] += dx * (C - dlen * rdlen) / (dlen * C);
628 VECTOR(dispy)[v] += dy * (C - dlen * rdlen) / (dlen * C);
629 VECTOR(dispy)[v] += dz * (C - dlen * rdlen) / (dlen * C);
630 VECTOR(dispx)[u] -= dx * (C - dlen * rdlen) / (dlen * C);
631 VECTOR(dispy)[u] -= dy * (C - dlen * rdlen) / (dlen * C);
632 VECTOR(dispz)[u] -= dz * (C - dlen * rdlen) / (dlen * C);
633 }
634 }
635 }
636
637 /* calculate attractive forces */
638 for (e = 0; e < no_edges; e++) {
639 /* each edges is an ordered pair of vertices v and u */
640 igraph_integer_t v = IGRAPH_FROM(graph, e);
641 igraph_integer_t u = IGRAPH_TO(graph, e);
642 igraph_real_t dx = MATRIX(*res, v, 0) - MATRIX(*res, u, 0);
643 igraph_real_t dy = MATRIX(*res, v, 1) - MATRIX(*res, u, 1);
644 igraph_real_t dz = MATRIX(*res, v, 2) - MATRIX(*res, u, 2);
645 igraph_real_t w = weight ? VECTOR(*weight)[e] : 1.0;
646 igraph_real_t dlen = sqrt(dx * dx + dy * dy + dz * dz) * w;
647 VECTOR(dispx)[v] -= (dx * dlen);
648 VECTOR(dispy)[v] -= (dy * dlen);
649 VECTOR(dispz)[v] -= (dz * dlen);
650 VECTOR(dispx)[u] += (dx * dlen);
651 VECTOR(dispy)[u] += (dy * dlen);
652 VECTOR(dispz)[u] += (dz * dlen);
653 }
654
655 /* limit max displacement to temperature t and prevent from
656 displacement outside frame */
657 for (v = 0; v < no_nodes; v++) {
658 igraph_real_t dx = VECTOR(dispx)[v] + RNG_UNIF01() * 1e-9;
659 igraph_real_t dy = VECTOR(dispy)[v] + RNG_UNIF01() * 1e-9;
660 igraph_real_t dz = VECTOR(dispz)[v] + RNG_UNIF01() * 1e-9;
661 igraph_real_t displen = sqrt(dx * dx + dy * dy + dz * dz);
662 igraph_real_t mx = fabs(dx) < temp ? dx : temp;
663 igraph_real_t my = fabs(dy) < temp ? dy : temp;
664 igraph_real_t mz = fabs(dz) < temp ? dz : temp;
665 if (displen > 0) {
666 MATRIX(*res, v, 0) += (dx / displen) * mx;
667 MATRIX(*res, v, 1) += (dy / displen) * my;
668 MATRIX(*res, v, 2) += (dz / displen) * mz;
669 }
670 if (minx && MATRIX(*res, v, 0) < VECTOR(*minx)[v]) {
671 MATRIX(*res, v, 0) = VECTOR(*minx)[v];
672 }
673 if (maxx && MATRIX(*res, v, 0) > VECTOR(*maxx)[v]) {
674 MATRIX(*res, v, 0) = VECTOR(*maxx)[v];
675 }
676 if (miny && MATRIX(*res, v, 1) < VECTOR(*miny)[v]) {
677 MATRIX(*res, v, 1) = VECTOR(*miny)[v];
678 }
679 if (maxy && MATRIX(*res, v, 1) > VECTOR(*maxy)[v]) {
680 MATRIX(*res, v, 1) = VECTOR(*maxy)[v];
681 }
682 if (minz && MATRIX(*res, v, 2) < VECTOR(*minz)[v]) {
683 MATRIX(*res, v, 2) = VECTOR(*minz)[v];
684 }
685 if (maxz && MATRIX(*res, v, 2) > VECTOR(*maxz)[v]) {
686 MATRIX(*res, v, 2) = VECTOR(*maxz)[v];
687 }
688 }
689
690 temp -= difftemp;
691 }
692
693 RNG_END();
694
695 igraph_vector_float_destroy(&dispx);
696 igraph_vector_float_destroy(&dispy);
697 igraph_vector_float_destroy(&dispz);
698 IGRAPH_FINALLY_CLEAN(3);
699
700 return 0;
701 }
702