4
Copyright (C) 2003, 2004, 2005, 2006 Gabor Csardi <csardi@rmki.kfki.hu>
5
MTA RMKI, Konkoly-Thege Miklos st. 29-33, Budapest 1121, Hungary
7
This program is free software; you can redistribute it and/or modify
8
it under the terms of the GNU General Public License as published by
9
the Free Software Foundation; either version 2 of the License, or
10
(at your option) any later version.
12
This program is distributed in the hope that it will be useful,
13
but WITHOUT ANY WARRANTY; without even the implied warranty of
14
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15
GNU General Public License for more details.
17
You should have received a copy of the GNU General Public License
18
along with this program; if not, write to the Free Software
19
Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
29
#include "igraph_math.h"
31
int igraph_i_layout_sphere_2d(igraph_matrix_t *coords, igraph_real_t *x, igraph_real_t *y,
33
int igraph_i_layout_sphere_3d(igraph_matrix_t *coords, igraph_real_t *x, igraph_real_t *y,
34
igraph_real_t *z, igraph_real_t *r);
37
* \section about_layouts
39
* <para>Layout generator functions (or at least most of them) try to place the
40
* vertices and edges of a graph on a 2D plane or in 3D space in a way
41
* which visually pleases the human eye.</para>
43
* <para>They take a graph object and a number of parameters as arguments
44
* and return an \type igraph_matrix_t, in which each row gives the
45
* coordinates of a vertex.</para>
50
* \function igraph_layout_random
51
* \brief Places the vertices uniform randomly on a plane.
53
* \param graph Pointer to an initialized graph object.
54
* \param res Pointer to an initialized graph object. This will
55
* contain the result and will be resized in needed.
56
* \return Error code. The current implementation always returns with
59
* Time complexity: O(|V|), the
63
int igraph_layout_random(const igraph_t *graph, igraph_matrix_t *res) {
65
long int no_of_nodes=igraph_vcount(graph);
68
IGRAPH_CHECK(igraph_matrix_resize(res, no_of_nodes, 2));
72
for (i=0; i<no_of_nodes; i++) {
73
MATRIX(*res, i, 0)=RNG_UNIF(-1, 1);
74
MATRIX(*res, i, 1)=RNG_UNIF(-1, 1);
83
* \function igraph_layout_random_3d
84
* \brief Random layout in 3D
86
* \param graph The graph to place.
87
* \param res Pointer to an initialized matrix object. It will be
88
* resized to hold the result.
89
* \return Error code. The current implementation always returns with
92
* Added in version 0.2.</para><para>
94
* Time complexity: O(|V|), the number of vertices.
97
int igraph_layout_random_3d(const igraph_t *graph, igraph_matrix_t *res) {
99
long int no_of_nodes=igraph_vcount(graph);
102
IGRAPH_CHECK(igraph_matrix_resize(res, no_of_nodes, 3));
106
for (i=0; i<no_of_nodes; i++) {
107
MATRIX(*res, i, 0)=RNG_UNIF(-1, 1);
108
MATRIX(*res, i, 1)=RNG_UNIF(-1, 1);
109
MATRIX(*res, i, 2)=RNG_UNIF(-1, 1);
119
* \function igraph_layout_circle
120
* \brief Places the vertices uniformly on a circle, in the order of vertex ids.
122
* \param graph Pointer to an initialized graph object.
123
* \param res Pointer to an initialized graph object. This will
124
* contain the result and will be resized in needed.
125
* \return Error code.
127
* Time complexity: O(|V|), the
128
* number of vertices.
131
int igraph_layout_circle(const igraph_t *graph, igraph_matrix_t *res) {
133
long int no_of_nodes=igraph_vcount(graph);
136
IGRAPH_CHECK(igraph_matrix_resize(res, no_of_nodes, 2));
138
for (i=0; i<no_of_nodes; i++) {
139
igraph_real_t phi=2*M_PI/no_of_nodes*i;
140
MATRIX(*res, i, 0)=cos(phi);
141
MATRIX(*res, i, 1)=sin(phi);
148
* \function igraph_layout_sphere
149
* \brief Places vertices (more or less) uniformly on a sphere.
152
* The algorithm was described in the following paper:
153
* Distributing many points on a sphere by E.B. Saff and
154
* A.B.J. Kuijlaars, \emb Mathematical Intelligencer \eme 19.1 (1997)
157
* \param graph Pointer to an initialized graph object.
158
* \param res Pointer to an initialized matrix object, the will be
159
* stored here. It will be resized.
160
* \return Error code. The current implementation always returns with
163
* Added in version 0.2.</para><para>
165
* Time complexity: O(|V|), the number of vertices in the graph.
168
int igraph_layout_sphere(const igraph_t *graph, igraph_matrix_t *res) {
170
long int no_of_nodes=igraph_vcount(graph);
174
IGRAPH_CHECK(igraph_matrix_resize(res, no_of_nodes, 3));
176
if (no_of_nodes != 0) {
177
MATRIX(*res, 0, 0)=M_PI;
178
MATRIX(*res, 0, 1)=0;
180
for (i=1; i<no_of_nodes-1; i++) {
181
h = -1 + 2*i/(double)(no_of_nodes-1);
182
MATRIX(*res, i, 0) = acos(h);
183
MATRIX(*res, i, 1) = fmod((MATRIX(*res, i-1, 1) +
184
3.6/sqrt(no_of_nodes*(1-h*h))), 2*M_PI);
185
IGRAPH_ALLOW_INTERRUPTION();
187
if (no_of_nodes >=2) {
188
MATRIX(*res, no_of_nodes-1, 0)=0;
189
MATRIX(*res, no_of_nodes-1, 1)=0;
192
for (i=0; i<no_of_nodes; i++) {
193
igraph_real_t x=cos(MATRIX(*res, i, 1))*sin(MATRIX(*res, i, 0));
194
igraph_real_t y=sin(MATRIX(*res, i, 1))*sin(MATRIX(*res, i, 0));
195
igraph_real_t z=cos(MATRIX(*res, i, 0));
196
MATRIX(*res, i, 0)=x;
197
MATRIX(*res, i, 1)=y;
198
MATRIX(*res, i, 2)=z;
199
IGRAPH_ALLOW_INTERRUPTION();
207
* \function igraph_layout_fruchterman_reingold
208
* \brief Places the vertices on a plane according to the Fruchterman-Reingold algorithm.
211
* This is a force-directed layout, see Fruchterman, T.M.J. and
212
* Reingold, E.M.: Graph Drawing by Force-directed Placement.
213
* Software -- Practice and Experience, 21/11, 1129--1164,
215
* This function was ported from the SNA R package.
216
* \param graph Pointer to an initialized graph object.
217
* \param res Pointer to an initialized matrix object. This will
218
* contain the result and will be resized in needed.
219
* \param niter The number of iterations to do.
220
* \param maxdelta The maximum distance to move a vertex in an
222
* \param area The area parameter of the algorithm.
223
* \param coolexp The cooling exponent of the simulated annealing.
224
* \param repulserad Determines the radius at which
225
* vertex-vertex repulsion cancels out attraction of
227
* \param use_seed Logical, if true the supplied values in the
228
* \p res argument are used as an initial layout, if
229
* false a random initial layout is used.
230
* \param weight Pointer to a vector containing edge weights,
231
* the attraction along the edges will be multiplied by these.
232
* It will be ignored if it is a null-pointer.
233
* \return Error code.
235
* Time complexity: O(|V|^2) in each
236
* iteration, |V| is the number of
237
* vertices in the graph.
240
int igraph_layout_fruchterman_reingold(const igraph_t *graph, igraph_matrix_t *res,
241
igraph_integer_t niter, igraph_real_t maxdelta,
242
igraph_real_t area, igraph_real_t coolexp,
243
igraph_real_t repulserad, igraph_bool_t use_seed,
244
const igraph_vector_t *weight) {
245
igraph_real_t frk,t,ded,xd,yd;
249
long int no_of_nodes=igraph_vcount(graph);
250
igraph_matrix_t dxdy=IGRAPH_MATRIX_NULL;
252
igraph_integer_t from, to;
254
if (weight && igraph_vector_size(weight) != igraph_ecount(graph)) {
255
IGRAPH_ERROR("Invalid weight vector length", IGRAPH_EINVAL);
258
IGRAPH_CHECK(igraph_matrix_resize(res, no_of_nodes, 2));
260
IGRAPH_CHECK(igraph_layout_random(graph, res));
262
IGRAPH_MATRIX_INIT_FINALLY(&dxdy, no_of_nodes, 2);
264
IGRAPH_CHECK(igraph_eit_create(graph, igraph_ess_all(0), &edgeit));
265
IGRAPH_FINALLY(igraph_eit_destroy, &edgeit);
267
frk=sqrt(area/no_of_nodes);
269
for(i=niter;i>0;i--) {
270
/* Report progress in approx. every 100th step */
272
IGRAPH_PROGRESS("Fruchterman-Reingold layout: ",
273
100.0-100.0*i/niter, NULL);
275
/* Set the temperature (maximum move/iteration) */
276
t=maxdelta*pow(i/(double)niter,coolexp);
277
/* Clear the deltas */
278
igraph_matrix_null(&dxdy);
279
/* Increment deltas for each undirected pair */
280
for(j=0;j<no_of_nodes;j++) {
281
IGRAPH_ALLOW_INTERRUPTION();
282
for(k=j+1;k<no_of_nodes;k++){
283
/* Obtain difference vector */
284
xd=MATRIX(*res, j, 0)-MATRIX(*res, k, 0);
285
yd=MATRIX(*res, j, 1)-MATRIX(*res, k, 1);
286
ded=sqrt(xd*xd+yd*yd); /* Get dyadic euclidean distance */
287
xd/=ded; /* Rescale differences to length 1 */
289
/* Calculate repulsive "force" */
290
rf=frk*frk*(1.0/ded-ded*ded/repulserad);
291
MATRIX(dxdy, j, 0)+=xd*rf; /* Add to the position change vector */
292
MATRIX(dxdy, k, 0)-=xd*rf;
293
MATRIX(dxdy, j, 1)+=yd*rf;
294
MATRIX(dxdy, k, 1)-=yd*rf;
297
/* Calculate the attractive "force" */
298
IGRAPH_EIT_RESET(edgeit);
299
while (!IGRAPH_EIT_END(edgeit)) {
300
long int edge=IGRAPH_EIT_GET(edgeit);
301
igraph_real_t w= weight ? VECTOR(*weight)[ edge ] : 1.0;
302
igraph_edge(graph, edge, &from, &to);
305
xd=MATRIX(*res, j, 0)-MATRIX(*res, k, 0);
306
yd=MATRIX(*res, j, 1)-MATRIX(*res, k, 1);
307
ded=sqrt(xd*xd+yd*yd); /* Get dyadic euclidean distance */
309
xd/=ded; /* Rescale differences to length 1 */
313
MATRIX(dxdy, j, 0)-=xd*af; /* Add to the position change vector */
314
MATRIX(dxdy, k, 0)+=xd*af;
315
MATRIX(dxdy, j, 1)-=yd*af;
316
MATRIX(dxdy, k, 1)+=yd*af;
317
IGRAPH_EIT_NEXT(edgeit);
320
/* Dampen motion, if needed, and move the points */
321
for(j=0;j<no_of_nodes;j++){
322
ded=sqrt(MATRIX(dxdy, j, 0)*MATRIX(dxdy, j, 0)+
323
MATRIX(dxdy, j, 1)*MATRIX(dxdy, j, 1));
324
if(ded>t){ /* Dampen to t */
326
MATRIX(dxdy, j, 0)*=ded;
327
MATRIX(dxdy, j, 1)*=ded;
329
MATRIX(*res, j, 0)+=MATRIX(dxdy, j, 0); /* Update positions */
330
MATRIX(*res, j, 1)+=MATRIX(dxdy, j, 1);
334
IGRAPH_PROGRESS("Fruchterman-Reingold layout: ", 100.0, NULL);
336
igraph_eit_destroy(&edgeit);
337
igraph_matrix_destroy(&dxdy);
338
IGRAPH_FINALLY_CLEAN(2);
344
* \function igraph_layout_fruchterman_reingold_3d
345
* \brief 3D Fruchterman-Reingold algorithm.
347
* This is the 3D version of the force based
348
* Fruchterman-Reingold layout (see \ref
349
* igraph_layout_fruchterman_reingold for the 2D version
352
* This function was ported from the SNA R package.
353
* \param graph Pointer to an initialized graph object.
354
* \param res Pointer to an initialized matrix object. This will
355
* contain the result and will be resized in needed.
356
* \param niter The number of iterations to do.
357
* \param maxdelta The maximum distance to move a vertex in an
359
* \param volume The volume parameter of the algorithm.
360
* \param coolexp The cooling exponent of the simulated annealing.
361
* \param repulserad Determines the radius at which
362
* vertex-vertex repulsion cancels out attraction of
364
* \param use_seed Logical, if true the supplied values in the
365
* \p res argument are used as an initial layout, if
366
* false a random initial layout is used.
367
* \param weight Pointer to a vector containing edge weights,
368
* the attraction along the edges will be multiplied by these.
369
* It will be ignored if it is a null-pointer.
370
* \return Error code.
372
* Added in version 0.2.</para><para>
374
* Time complexity: O(|V|^2) in each
375
* iteration, |V| is the number of
376
* vertices in the graph.
380
int igraph_layout_fruchterman_reingold_3d(const igraph_t *graph,
381
igraph_matrix_t *res,
382
igraph_integer_t niter, igraph_real_t maxdelta,
383
igraph_real_t volume, igraph_real_t coolexp,
384
igraph_real_t repulserad,
385
igraph_bool_t use_seed,
386
const igraph_vector_t *weight) {
388
igraph_real_t frk, t, ded, xd, yd, zd;
389
igraph_matrix_t dxdydz;
390
igraph_real_t rf, af;
393
long int no_of_nodes=igraph_vcount(graph);
395
igraph_integer_t from, to;
397
if (weight && igraph_vector_size(weight) != igraph_ecount(graph)) {
398
IGRAPH_ERROR("Invalid weight vector length", IGRAPH_EINVAL);
401
IGRAPH_CHECK(igraph_matrix_init(&dxdydz, no_of_nodes, 3));
402
IGRAPH_FINALLY(igraph_matrix_destroy, &dxdydz);
404
IGRAPH_CHECK(igraph_matrix_resize(res, no_of_nodes, 3));
406
IGRAPH_CHECK(igraph_layout_random_3d(graph, res));
409
IGRAPH_CHECK(igraph_eit_create(graph, igraph_ess_all(0), &edgeit));
410
IGRAPH_FINALLY(igraph_eit_destroy, &edgeit);
412
frk=pow(volume/(double)no_of_nodes,1.0/3.0); /*Define the F-R constant*/
414
/*Run the annealing loop*/
415
for(i=niter;i>=0;i--){
417
IGRAPH_PROGRESS("3D Fruchterman-Reingold layout: ",
418
100.0-100.0*i/niter, NULL);
420
/*Set the temperature (maximum move/iteration)*/
421
t=maxdelta*pow(i/(double)niter,coolexp);
423
igraph_matrix_null(&dxdydz);
424
/*Increment deltas for each undirected pair*/
425
for(j=0;j<no_of_nodes;j++) {
426
IGRAPH_ALLOW_INTERRUPTION();
427
for(k=j+1;k<no_of_nodes;k++){
428
/*Obtain difference vector*/
429
xd=MATRIX(*res, j, 0)-MATRIX(*res, k, 0);
430
yd=MATRIX(*res, j, 1)-MATRIX(*res, k, 1);
431
zd=MATRIX(*res, j, 2)-MATRIX(*res, k, 2);
432
ded=sqrt(xd*xd+yd*yd+zd*zd); /*Get dyadic euclidean distance*/
434
xd/=ded; /*Rescale differences to length 1*/
438
/*Calculate repulsive "force"*/
439
rf=frk*frk*(1.0/ded-ded*ded/repulserad);
440
MATRIX(dxdydz, j, 0)+=xd*rf; /*Add to the position change vector*/
441
MATRIX(dxdydz, k, 0)-=xd*rf;
442
MATRIX(dxdydz, j, 1)+=yd*rf;
443
MATRIX(dxdydz, k, 1)-=yd*rf;
444
MATRIX(dxdydz, j, 2)+=zd*rf;
445
MATRIX(dxdydz, k, 2)-=zd*rf;
449
/*Calculate the attractive "force"*/
450
IGRAPH_EIT_RESET(edgeit);
451
while (!IGRAPH_EIT_END(edgeit)) {
452
long int edge=IGRAPH_EIT_GET(edgeit);
453
igraph_real_t w= weight ? VECTOR(*weight)[edge] : 1.0;
454
igraph_edge(graph, edge, &from, &to);
457
xd=MATRIX(*res, j, 0)-MATRIX(*res, k, 0);
458
yd=MATRIX(*res, j, 1)-MATRIX(*res, k, 1);
459
zd=MATRIX(*res, j, 2)-MATRIX(*res, k, 2);
460
ded=sqrt(xd*xd+yd*yd+zd*zd); /*Get dyadic euclidean distance*/
462
xd/=ded; /*Rescale differences to length 1*/
467
MATRIX(dxdydz, j, 0)-=xd*af; /*Add to the position change vector*/
468
MATRIX(dxdydz, k, 0)+=xd*af;
469
MATRIX(dxdydz, j, 1)-=yd*af;
470
MATRIX(dxdydz, k, 1)+=yd*af;
471
MATRIX(dxdydz, j, 2)-=zd*af;
472
MATRIX(dxdydz, k, 2)+=zd*af;
473
IGRAPH_EIT_NEXT(edgeit);
476
/*Dampen motion, if needed, and move the points*/
477
for(j=0;j<no_of_nodes;j++){
478
ded=sqrt(MATRIX(dxdydz, j, 0)*MATRIX(dxdydz, j, 0)+
479
MATRIX(dxdydz, j, 1)*MATRIX(dxdydz, j, 1)+
480
MATRIX(dxdydz, j, 2)*MATRIX(dxdydz, j, 2));
481
if(ded>t){ /*Dampen to t*/
483
MATRIX(dxdydz, j, 0)*=ded;
484
MATRIX(dxdydz, j, 1)*=ded;
485
MATRIX(dxdydz, j, 2)*=ded;
487
MATRIX(*res, j, 0)+=MATRIX(dxdydz, j, 0); /*Update positions*/
488
MATRIX(*res, j, 1)+=MATRIX(dxdydz, j, 1);
489
MATRIX(*res, j, 2)+=MATRIX(dxdydz, j, 2);
493
IGRAPH_PROGRESS("3D Fruchterman-Reingold layout: ", 100.0, NULL);
495
igraph_matrix_destroy(&dxdydz);
496
igraph_eit_destroy(&edgeit);
497
IGRAPH_FINALLY_CLEAN(2);
503
* \function igraph_layout_kamada_kawai
504
* \brief Places the vertices on a plane according the Kamada-Kawai algorithm.
507
* This is a force directed layout, see Kamada, T. and Kawai, S.: An
508
* Algorithm for Drawing General Undirected Graphs. Information
509
* Processing Letters, 31/1, 7--15, 1989.
510
* This function was ported from the SNA R package.
511
* \param graph A graph object.
512
* \param res Pointer to an initialized matrix object. This will
513
* contain the result and will be resized if needed.
514
* \param niter The number of iterations to perform.
515
* \param sigma Sets the base standard deviation of position
517
* \param initemp Sets the initial temperature for the annealing.
518
* \param coolexp The cooling exponent of the annealing.
519
* \param kkconst The Kamada-Kawai vertex attraction constant.
520
* \param use_seed Boolean, whether to use the values cupplied in the \p res
521
* argument as the initial configuration. If zero then a random initial
522
* configuration is used.
523
* \return Error code.
525
* Time complexity: O(|V|^2) for each
526
* iteration, |V| is the number of
527
* vertices in the graph.
530
int igraph_layout_kamada_kawai(const igraph_t *graph, igraph_matrix_t *res,
531
igraph_integer_t niter, igraph_real_t sigma,
532
igraph_real_t initemp, igraph_real_t coolexp,
533
igraph_real_t kkconst, igraph_bool_t use_seed) {
535
igraph_real_t temp, candx, candy;
536
igraph_real_t dpot, odis, ndis, osqd, nsqd;
538
igraph_matrix_t elen;
540
/* Define various things */
541
n=igraph_vcount(graph);
543
/* Calculate elen, initial x & y */
547
IGRAPH_CHECK(igraph_matrix_resize(res, n, 2));
548
IGRAPH_MATRIX_INIT_FINALLY(&elen, n, n);
549
IGRAPH_CHECK(igraph_shortest_paths(graph, &elen, igraph_vss_all(),
552
for (i=0; i<n; i++) {
553
MATRIX(elen, i, i) = sqrt(n);
554
MATRIX(*res, i, 0) = RNG_NORMAL(0, n/4.0);
555
MATRIX(*res, i, 1) = RNG_NORMAL(0, n/4.0);
559
/*Perform the annealing loop*/
561
for(i=0;i<niter;i++){
562
/* Report progress in approx. every 100th step */
564
IGRAPH_PROGRESS("Kamada-Kawai layout: ",
565
100.0*i/niter, NULL);
566
/*Update each vertex*/
568
IGRAPH_ALLOW_INTERRUPTION();
569
/*Draw the candidate via a gaussian perturbation*/
570
candx=RNG_NORMAL(MATRIX(*res, j, 0),sigma*temp/initemp);
571
candy=RNG_NORMAL(MATRIX(*res, j, 1),sigma*temp/initemp);
572
/*Calculate the potential difference for the new position*/
574
for(k=0;k<n;k++) /*Potential differences for pairwise effects*/
576
odis=sqrt((MATRIX(*res, j, 0)-MATRIX(*res, k, 0))*
577
(MATRIX(*res, j, 0)-MATRIX(*res, k, 0))+
578
(MATRIX(*res, j, 1)-MATRIX(*res, k, 1))*
579
(MATRIX(*res, j, 1)-MATRIX(*res, k, 1)));
580
ndis=sqrt((candx-MATRIX(*res, k, 0))*(candx-MATRIX(*res, k, 0))+
581
(candy-MATRIX(*res, k, 1))*(candy-MATRIX(*res, k, 1)));
582
osqd=(odis-MATRIX(elen, j, k))*(odis-MATRIX(elen, j, k));
583
nsqd=(ndis-MATRIX(elen, j, k))*(ndis-MATRIX(elen, j, k));
584
dpot+=kkconst*(osqd-nsqd)/(MATRIX(elen, j, k)*MATRIX(elen, j, k));
586
/*Make a keep/reject decision*/
587
if(log(RNG_UNIF(0.0,1.0))<dpot/temp){
588
MATRIX(*res, j, 0)=candx;
589
MATRIX(*res, j, 1)=candy;
596
IGRAPH_PROGRESS("Kamada-Kawai layout: ", 100.0, NULL);
599
igraph_matrix_destroy(&elen);
600
IGRAPH_FINALLY_CLEAN(1);
606
* \function igraph_layout_kamada_kawai_3d
607
* \brief 3D version of the force based Kamada-Kawai layout.
609
* The pair of the \ref igraph_layout_kamada_kawai 2D layout generator
612
* This function was ported from the SNA R package.
613
* \param graph A graph object.
614
* \param res Pointer to an initialized matrix object. This will
615
* contain the result and will be resized if needed.
616
* \param niter The number of iterations to perform.
617
* \param sigma Sets the base standard deviation of position
619
* \param initemp Sets the initial temperature for the annealing.
620
* \param coolexp The cooling exponent of the annealing.
621
* \param kkconst The Kamada-Kawai vertex attraction constant.
622
* \param use_seed Boolean, whether to use the values cupplied in the \p res
623
* argument as the initial configuration. If zero then a random initial
624
* configuration is used.
625
* \return Error code.
627
* Added in version 0.2.</para><para>
629
* Time complexity: O(|V|^2) for each
630
* iteration, |V| is the number of
631
* vertices in the graph.
634
int igraph_layout_kamada_kawai_3d(const igraph_t *graph, igraph_matrix_t *res,
635
igraph_integer_t niter, igraph_real_t sigma,
636
igraph_real_t initemp, igraph_real_t coolexp,
637
igraph_real_t kkconst, igraph_bool_t use_seed) {
638
igraph_real_t temp, candx, candy, candz;
639
igraph_real_t dpot, odis, ndis, osqd, nsqd;
641
long int no_of_nodes=igraph_vcount(graph);
642
igraph_matrix_t elen;
646
IGRAPH_CHECK(igraph_matrix_resize(res, no_of_nodes, 3));
647
IGRAPH_MATRIX_INIT_FINALLY(&elen, no_of_nodes, no_of_nodes);
648
IGRAPH_CHECK(igraph_shortest_paths(graph, &elen, igraph_vss_all(),
652
for (i=0; i<no_of_nodes; i++) {
653
MATRIX(elen, i, i) = sqrt(no_of_nodes);
654
MATRIX(*res, i, 0) = RNG_NORMAL(0, no_of_nodes/4.0);
655
MATRIX(*res, i, 1) = RNG_NORMAL(0, no_of_nodes/4.0);
656
MATRIX(*res, i, 2) = RNG_NORMAL(0, no_of_nodes/4.0);
661
for(i=0;i<niter;i++){
663
IGRAPH_PROGRESS("3D Kamada-Kawai layout: ",
664
100.0*i/niter, NULL);
666
/*Update each vertex*/
667
for(j=0;j<no_of_nodes;j++){
668
IGRAPH_ALLOW_INTERRUPTION();
669
/*Draw the candidate via a gaussian perturbation*/
670
candx=RNG_NORMAL(MATRIX(*res, j, 0) ,sigma*temp/initemp);
671
candy=RNG_NORMAL(MATRIX(*res, j, 1) ,sigma*temp/initemp);
672
candz=RNG_NORMAL(MATRIX(*res, j, 2) ,sigma*temp/initemp);
673
/*Calculate the potential difference for the new position*/
675
for(k=0;k<no_of_nodes;k++) /*Potential differences for pairwise effects*/
677
odis=sqrt((MATRIX(*res, j, 0)-MATRIX(*res, k, 0))*
678
(MATRIX(*res, j, 0)-MATRIX(*res, k, 0))+
679
(MATRIX(*res, j, 1)-MATRIX(*res, k, 1))*
680
(MATRIX(*res, j, 1)-MATRIX(*res, k, 1))+
681
(MATRIX(*res, j, 2)-MATRIX(*res, k, 2))*
682
(MATRIX(*res, j, 2)-MATRIX(*res, k, 2)));
683
ndis=sqrt((candx-MATRIX(*res, k, 0))*(candx-MATRIX(*res, k, 0))+
684
(candy-MATRIX(*res, k, 1))*(candy-MATRIX(*res, k, 1))+
685
(candz-MATRIX(*res, k, 2))*(candz-MATRIX(*res, k, 2)));
686
osqd=(odis-MATRIX(elen, j, k))*(odis-MATRIX(elen, j, k));
687
nsqd=(ndis-MATRIX(elen, j, k))*(ndis-MATRIX(elen, j, k));
688
dpot+=kkconst*(osqd-nsqd)/(MATRIX(elen, j, k)*MATRIX(elen, j, k));
690
/*Make a keep/reject decision*/
691
if(log(RNG_UNIF(0.0,1.0))<dpot/temp){
692
MATRIX(*res, j, 0)=candx;
693
MATRIX(*res, j, 1)=candy;
694
MATRIX(*res, j, 2)=candz;
701
IGRAPH_PROGRESS("3D Kamada-Kawai layout: ", 100.0, NULL);
704
igraph_matrix_destroy(&elen);
705
IGRAPH_FINALLY_CLEAN(1);
710
int igraph_layout_springs(const igraph_t *graph, igraph_matrix_t *res,
711
igraph_real_t mass, igraph_real_t equil, igraph_real_t k,
712
igraph_real_t repeqdis, igraph_real_t kfr, igraph_bool_t repulse) {
717
void igraph_i_norm2d(igraph_real_t *x, igraph_real_t *y) {
718
igraph_real_t len=sqrt((*x)*(*x) + (*y)*(*y));
726
* \function igraph_layout_lgl
727
* \brief Force based layout algorithm for large graphs.
730
* This is a layout generator similar to the Large Graph Layout
731
* algorithm and program
732
* (http://bioinformatics.icmb.utexas.edu/lgl/). But unlike LGL, this
733
* version uses a Fruchterman-Reingold style simulated annealing
734
* algorithm for placing the vertices. The speedup is achived by
735
* placing the vertices on a grid and calculating the repulsion only
736
* for vertices which are closer to each other than a limit.
738
* \param graph The (initialized) graph object to place.
739
* \param res Pointer to an initialized matrix object to hold the
740
* result. It will be resized if needed.
741
* \param maxit The maximum number of cooling iterations to perform
742
* for each layout step.
743
* \param maxdelta The maximum length of the move allowed for a vertex
744
* in a single iteration.
745
* \param area This parameter gives the area of the square on which
746
* the vertices will be placed.
747
* \param coolexp The cooling exponent.
748
* \param repulserad Determines the radius at which vertex-vertex
749
* repulsion cancels out attraction of adjacenct vertices.
750
* \param cellsize The size of the grid cells, one side of the
752
* \param proot The root vertex, this is placed first, its neighbors
753
* in the first iteration, second neighbors in the second, etc. If
754
* negative then a random vertex is chosen.
755
* \return Error code.
757
* Added in version 0.2.</para><para>
759
* Time complexity: ideally O(dia*maxit*(|V|+|E|)), |V| is the number
761
* dia is the diameter of the graph, worst case complexity is still
762
* O(dia*maxit*(|V|^2+|E|)), this is the case when all vertices happen to be
763
* in the same grid cell.
766
int igraph_layout_lgl(const igraph_t *graph, igraph_matrix_t *res,
767
igraph_integer_t maxit, igraph_real_t maxdelta,
768
igraph_real_t area, igraph_real_t coolexp,
769
igraph_real_t repulserad, igraph_real_t cellsize,
770
igraph_integer_t proot) {
773
long int no_of_nodes=igraph_vcount(graph);
774
long int no_of_edges=igraph_ecount(graph);
777
long int no_of_layers, actlayer=0;
778
igraph_vector_t vids;
779
igraph_vector_t layers;
780
igraph_vector_t parents;
781
igraph_vector_t edges;
782
igraph_2dgrid_t grid;
783
igraph_vector_t eids;
784
igraph_vector_t forcex;
785
igraph_vector_t forcey;
787
igraph_real_t frk=sqrt(area/no_of_nodes);
790
IGRAPH_CHECK(igraph_minimum_spanning_tree_unweighted(graph, &mst));
791
IGRAPH_FINALLY(igraph_destroy, &mst);
793
IGRAPH_CHECK(igraph_matrix_resize(res, no_of_nodes, 2));
795
/* Determine the root vertex, random pick right now */
797
root=RNG_INTEGER(0, no_of_nodes-1);
802
/* Assign the layers */
803
IGRAPH_VECTOR_INIT_FINALLY(&vids, 0);
804
IGRAPH_VECTOR_INIT_FINALLY(&layers, 0);
805
IGRAPH_VECTOR_INIT_FINALLY(&parents, 0);
806
IGRAPH_CHECK(igraph_bfs(&mst, root, IGRAPH_ALL, &vids, &layers, &parents));
807
no_of_layers=igraph_vector_size(&layers)-1;
809
/* We don't need the mst any more */
810
igraph_destroy(&mst);
811
igraph_empty(&mst, 0, IGRAPH_UNDIRECTED); /* to make finalization work */
813
IGRAPH_VECTOR_INIT_FINALLY(&edges, 0);
814
IGRAPH_CHECK(igraph_vector_reserve(&edges, no_of_edges));
815
IGRAPH_VECTOR_INIT_FINALLY(&eids, 0);
816
IGRAPH_VECTOR_INIT_FINALLY(&forcex, no_of_nodes);
817
IGRAPH_VECTOR_INIT_FINALLY(&forcey, no_of_nodes);
819
/* Place the vertices randomly */
820
IGRAPH_CHECK(igraph_layout_random(graph, res));
821
igraph_matrix_scale(res, 1e6);
823
/* This is the grid for calculating the vertices near to a given vertex */
824
IGRAPH_CHECK(igraph_2dgrid_init(&grid, res,
825
-sqrt(area/M_PI),sqrt(area/M_PI), cellsize,
826
-sqrt(area/M_PI),sqrt(area/M_PI), cellsize));
827
IGRAPH_FINALLY(igraph_2dgrid_destroy, &grid);
829
/* Place the root vertex */
830
igraph_2dgrid_add(&grid, root, 0, 0);
832
for (actlayer=1; actlayer<no_of_layers; actlayer++) {
836
for (actlayer=1; actlayer<no_of_layers; actlayer++) {
840
igraph_real_t massx, massy;
841
igraph_real_t px, py;
842
igraph_real_t sx, sy;
845
igraph_real_t epsilon=10e-6;
846
igraph_real_t maxchange=epsilon+1;
848
igraph_real_t sconst=sqrt(area/M_PI) / H_n;
849
igraph_2dgrid_iterator_t vidit;
851
/* printf("Layer %li:\n", actlayer); */
853
/*-----------------------------------------*/
854
/* Step 1: place the next layer on spheres */
855
/*-----------------------------------------*/
859
j=VECTOR(layers)[actlayer];
860
for (i=VECTOR(layers)[actlayer-1]; i<VECTOR(layers)[actlayer]; i++) {
862
long int vid=VECTOR(vids)[i];
863
long int par=VECTOR(parents)[vid];
864
IGRAPH_ALLOW_INTERRUPTION();
865
igraph_2dgrid_getcenter(&grid, &massx, &massy);
866
igraph_i_norm2d(&massx, &massy);
867
px=MATRIX(*res, vid, 0)-MATRIX(*res, par, 0);
868
py=MATRIX(*res, vid, 1)-MATRIX(*res, par, 1);
869
igraph_i_norm2d(&px, &py);
870
sx=c*(massx+px)+MATRIX(*res, vid, 0);
871
sy=c*(massy+py)+MATRIX(*res, vid, 1);
873
/* The neighbors of 'vid' */
874
while (j < VECTOR(layers)[actlayer+1] &&
875
VECTOR(parents)[(long int)VECTOR(vids)[j]]==vid) {
876
igraph_real_t rx, ry;
878
igraph_real_t phi=2*M_PI/(VECTOR(layers)[2]-1)*(j-1);
885
igraph_i_norm2d(&rx, &ry);
886
rx = rx / actlayer * sconst;
887
ry = ry / actlayer * sconst;
888
igraph_2dgrid_add(&grid, VECTOR(vids)[j], sx+rx, sy+ry);
895
/*-----------------------------------------*/
896
/* Step 2: add the edges of the next layer */
897
/*-----------------------------------------*/
899
for (j=VECTOR(layers)[actlayer]; j<VECTOR(layers)[actlayer+1]; j++) {
900
long int vid=VECTOR(vids)[j];
902
IGRAPH_ALLOW_INTERRUPTION();
903
IGRAPH_CHECK(igraph_adjacent(graph, &eids, vid, IGRAPH_ALL));
904
for (k=0;k<igraph_vector_size(&eids);k++) {
905
long int eid=VECTOR(eids)[k];
906
igraph_integer_t from, to;
907
igraph_edge(graph, eid, &from, &to);
908
if ((from != vid && igraph_2dgrid_in(&grid, from)) ||
909
(to != vid && igraph_2dgrid_in(&grid, to))) {
910
igraph_vector_push_back(&edges, eid);
915
/*-----------------------------------------*/
916
/* Step 3: let the springs spring */
917
/*-----------------------------------------*/
920
while (it < maxit && maxchange > epsilon) {
922
igraph_real_t t=maxdelta*pow((maxit-it)/(double)maxit, coolexp);
925
IGRAPH_PROGRESS("Large graph layout",
926
100.0*((actlayer-1.0)/(no_of_layers-1.0)+((float)it)/(maxit*(no_of_layers-1.0))),
930
igraph_vector_null(&forcex);
931
igraph_vector_null(&forcey);
934
/* attractive "forces" along the edges */
935
for (j=0; j<igraph_vector_size(&edges); j++) {
936
igraph_integer_t from, to;
937
igraph_real_t xd, yd, dist, force;
938
IGRAPH_ALLOW_INTERRUPTION();
939
igraph_edge(graph, VECTOR(edges)[j], &from, &to);
940
xd=MATRIX(*res, (long int)from, 0)-MATRIX(*res, (long int)to, 0);
941
yd=MATRIX(*res, (long int)from, 1)-MATRIX(*res, (long int)to, 1);
942
dist=sqrt(xd*xd+yd*yd);
943
if (dist != 0) { xd/=dist; yd/=dist; }
945
VECTOR(forcex)[(long int)from] -= xd*force;
946
VECTOR(forcex)[(long int)to] += xd*force;
947
VECTOR(forcey)[(long int)from] -= yd*force;
948
VECTOR(forcey)[(long int)to] += yd*force;
951
/* repulsive "forces" of the vertices nearby */
953
igraph_2dgrid_reset(&grid, &vidit);
954
while ( (vid=igraph_2dgrid_next(&grid, &vidit)-1) != -1) {
955
while ( (nei=igraph_2dgrid_next_nei(&grid, &vidit)-1) != -1) {
956
igraph_real_t xd=MATRIX(*res, (long int)vid, 0)-
957
MATRIX(*res, (long int)nei, 0);
958
igraph_real_t yd=MATRIX(*res, (long int)vid, 1)-
959
MATRIX(*res, (long int)nei, 1);
960
igraph_real_t dist=sqrt(xd*xd+yd*yd);
962
if (dist < cellsize) {
964
if (dist==0) { dist=epsilon; };
966
force=frk*frk*(1.0/dist-dist*dist/repulserad);
967
VECTOR(forcex)[(long int)vid] += xd*force;
968
VECTOR(forcex)[(long int)nei] -= xd*force;
969
VECTOR(forcey)[(long int)vid] += yd*force;
970
VECTOR(forcey)[(long int)nei] -= yd*force;
975
/* printf("verties: %li iterations: %li\n", */
976
/* (long int) VECTOR(layers)[actlayer+1], pairs); */
978
/* apply the changes */
979
for (j=0; j<VECTOR(layers)[actlayer+1]; j++) {
980
long int vid=VECTOR(vids)[j];
981
igraph_real_t fx=VECTOR(forcex)[vid];
982
igraph_real_t fy=VECTOR(forcey)[vid];
983
igraph_real_t ded=sqrt(fx*fx+fy*fy);
988
igraph_2dgrid_move(&grid, vid, fx, fy);
989
if (fx > maxchange) { maxchange=fx; }
990
if (fy > maxchange) { maxchange=fy; }
993
/* printf("%li iterations, maxchange: %f\n", it, (double)maxchange); */
997
IGRAPH_PROGRESS("Large graph layout", 100.0, 0);
998
igraph_destroy(&mst);
999
igraph_vector_destroy(&vids);
1000
igraph_vector_destroy(&layers);
1001
igraph_vector_destroy(&parents);
1002
igraph_vector_destroy(&edges);
1003
igraph_2dgrid_destroy(&grid);
1004
igraph_vector_destroy(&eids);
1005
igraph_vector_destroy(&forcex);
1006
igraph_vector_destroy(&forcey);
1007
IGRAPH_FINALLY_CLEAN(9);
1013
* \function igraph_layout_grid_fruchterman_reingold
1014
* \brief Force based layout generator for large graphs.
1017
* This algorithm is the same as the Fruchterman-Reingold layout
1018
* generator, but it partitions the 2d space to a grid and and vertex
1019
* repulsion is calculated only for vertices nearby.
1021
* \param graph The graph object.
1022
* \param res The result, the coordinates in a matrix. The parameter
1023
* should point to an initialized matrix object and will be resized.
1024
* \param niter Number of iterations to perform.
1025
* \param maxdelta Maximum distance to move a vertex in an iteration.
1026
* \param area The area of the square on which the vertices will be
1028
* \param coolexp The cooling exponent.
1029
* \param repulserad Determines the radius at which vertex-vertex
1030
* repulsion cancels out attraction of adjacenct vertices.
1031
* \param cellsize The size of the grid cells.
1032
* \param use_seed Logical, if true, the coordinates passed in \p res
1033
* (should have the appropriate size) will be used for the first
1035
* \return Error code.
1037
* Added in version 0.2.</para><para>
1039
* Time complexity: ideally (constant number of vertices in each cell)
1040
* O(niter*(|V|+|E|)), in the worst case O(niter*(|V|^2+|E|)).
1043
int igraph_layout_grid_fruchterman_reingold(const igraph_t *graph,
1044
igraph_matrix_t *res,
1045
igraph_integer_t niter, igraph_real_t maxdelta,
1046
igraph_real_t area, igraph_real_t coolexp,
1047
igraph_real_t repulserad,
1048
igraph_real_t cellsize,
1049
igraph_bool_t use_seed) {
1051
long int no_of_nodes=igraph_vcount(graph);
1052
long int no_of_edges=igraph_ecount(graph);
1053
igraph_2dgrid_t grid;
1054
igraph_vector_t forcex;
1055
igraph_vector_t forcey;
1057
igraph_2dgrid_iterator_t vidit;
1059
igraph_real_t frk=sqrt(area/no_of_nodes);
1061
IGRAPH_CHECK(igraph_matrix_resize(res, no_of_nodes, 2));
1062
IGRAPH_VECTOR_INIT_FINALLY(&forcex, no_of_nodes);
1063
IGRAPH_VECTOR_INIT_FINALLY(&forcey, no_of_nodes);
1065
/* initial layout */
1067
IGRAPH_CHECK(igraph_layout_random(graph, res));
1068
igraph_matrix_scale(res, sqrt(area/M_PI));
1072
IGRAPH_CHECK(igraph_2dgrid_init(&grid, res,
1073
-sqrt(area/M_PI),sqrt(area/M_PI), cellsize,
1074
-sqrt(area/M_PI),sqrt(area/M_PI), cellsize));
1075
IGRAPH_FINALLY(igraph_2dgrid_destroy, &grid);
1077
/* place vertices on grid */
1078
for (i=0; i<no_of_nodes; i++) {
1079
igraph_2dgrid_add2(&grid, i);
1084
igraph_real_t t=maxdelta*pow((niter-it)/(double)niter, coolexp);
1087
/* Report progress */
1089
IGRAPH_PROGRESS("Grid based Fruchterman-Reingold layout: ",
1090
(100.0*it)/niter, NULL);
1093
igraph_vector_null(&forcex);
1094
igraph_vector_null(&forcey);
1097
for (j=0; j<no_of_edges; j++) {
1098
igraph_integer_t from, to;
1099
igraph_real_t xd, yd, dist, force;
1100
igraph_edge(graph, j, &from, &to);
1101
xd=MATRIX(*res, (long int)from, 0)-MATRIX(*res, (long int)to, 0);
1102
yd=MATRIX(*res, (long int)from, 1)-MATRIX(*res, (long int)to, 1);
1103
dist=sqrt(xd*xd+yd*yd);
1104
if (dist != 0) { xd/=dist; yd/=dist; }
1105
force=dist*dist/frk;
1106
VECTOR(forcex)[(long int)from] -= xd*force;
1107
VECTOR(forcex)[(long int)to] += xd*force;
1108
VECTOR(forcey)[(long int)from] -= yd*force;
1109
VECTOR(forcey)[(long int)to] += yd*force;
1113
igraph_2dgrid_reset(&grid, &vidit);
1114
while ( (vid=igraph_2dgrid_next(&grid, &vidit)-1) != -1) {
1115
IGRAPH_ALLOW_INTERRUPTION();
1116
while ( (nei=igraph_2dgrid_next_nei(&grid, &vidit)-1) != -1) {
1117
igraph_real_t xd=MATRIX(*res, (long int)vid, 0)-
1118
MATRIX(*res, (long int)nei, 0);
1119
igraph_real_t yd=MATRIX(*res, (long int)vid, 1)-
1120
MATRIX(*res, (long int)nei, 1);
1121
igraph_real_t dist=sqrt(xd*xd+yd*yd);
1122
igraph_real_t force;
1123
if (dist < cellsize) {
1124
if (dist==0) { dist=1e-6; };
1126
force=frk*frk*(1.0/dist-dist*dist/repulserad);
1127
VECTOR(forcex)[(long int)vid] += xd*force;
1128
VECTOR(forcex)[(long int)nei] -= xd*force;
1129
VECTOR(forcey)[(long int)vid] += yd*force;
1130
VECTOR(forcey)[(long int)nei] -= yd*force;
1136
for (j=0; j<no_of_nodes; j++) {
1138
igraph_real_t fx=VECTOR(forcex)[vid];
1139
igraph_real_t fy=VECTOR(forcey)[vid];
1140
igraph_real_t ded=sqrt(fx*fx+fy*fy);
1145
igraph_2dgrid_move(&grid, vid, fx, fy);
1151
IGRAPH_PROGRESS("Grid based Fruchterman-Reingold layout: ",
1154
igraph_vector_destroy(&forcex);
1155
igraph_vector_destroy(&forcey);
1156
igraph_2dgrid_destroy(&grid);
1157
IGRAPH_FINALLY_CLEAN(3);
1161
/* Internal structure for Reingold-Tilford layout */
1162
struct igraph_i_reingold_tilford_vertex {
1163
long int parent; /* Parent node index */
1164
long int level; /* Level of the node */
1165
igraph_real_t offset; /* X offset from parent node */
1166
long int left_contour; /* Next left node of the contour
1167
of the subtree rooted at this node */
1168
long int right_contour; /* Next right node of the contour
1169
of the subtree rooted at this node */
1170
igraph_real_t offset_follow_lc; /* X offset when following the left contour */
1171
igraph_real_t offset_follow_rc; /* X offset when following the right contour */
1174
int igraph_i_layout_reingold_tilford_postorder(struct igraph_i_reingold_tilford_vertex *vdata,
1175
long int node, long int vcount);
1176
int igraph_i_layout_reingold_tilford_calc_coords(struct igraph_i_reingold_tilford_vertex *vdata,
1177
igraph_matrix_t *res, long int node,
1178
long int vcount, igraph_real_t xpos);
1181
* \function igraph_layout_reingold_tilford
1182
* \brief Reingold-Tilford layout for tree graphs
1185
* Arranges the nodes in a tree where the given node is used as the root.
1186
* The tree is directed downwards and the parents are centered above its
1187
* children. For the exact algorithm, see:
1190
* Reingold, E and Tilford, J: Tidier drawing of trees.
1191
* IEEE Trans. Softw. Eng., SE-7(2):223--228, 1981
1194
* If the given graph is not a tree, a breadth-first search is executed
1195
* first to obtain a possible spanning tree.
1197
* \param graph The graph object.
1198
* \param res The result, the coordinates in a matrix. The parameter
1199
* should point to an initialized matrix object and will be resized.
1200
* \param root The index of the root vertex.
1201
* \return Error code.
1203
* Added in version 0.2.
1206
* TODO: decompose and merge for not fully connected graphs
1207
* TODO: possible speedup could be achieved if we use a table for storing
1208
* the children of each node in the tree. (Now the implementation uses a
1209
* single array containing the parent of each node and a node's children
1210
* are determined by looking for other nodes that have this node as parent)
1212
* \sa \ref igraph_layout_reingold_tilford_circular().
1214
int igraph_layout_reingold_tilford(const igraph_t *graph,
1215
igraph_matrix_t *res, long int root) {
1216
long int no_of_nodes=igraph_vcount(graph);
1218
igraph_dqueue_t q=IGRAPH_DQUEUE_NULL;
1219
igraph_adjlist_t allneis;
1220
igraph_vector_t *neis;
1221
struct igraph_i_reingold_tilford_vertex *vdata;
1223
if (root<0 || root>=no_of_nodes) {
1224
IGRAPH_ERROR("invalid vertex id", IGRAPH_EINVVID);
1227
IGRAPH_CHECK(igraph_matrix_resize(res, no_of_nodes, 2));
1228
IGRAPH_DQUEUE_INIT_FINALLY(&q, 100);
1230
IGRAPH_CHECK(igraph_adjlist_init(graph, &allneis, IGRAPH_ALL));
1231
IGRAPH_FINALLY(igraph_adjlist_destroy, &allneis);
1233
vdata=igraph_Calloc(no_of_nodes, struct igraph_i_reingold_tilford_vertex);
1235
IGRAPH_ERROR("igraph_layout_reingold_tilford failed", IGRAPH_ENOMEM);
1237
IGRAPH_FINALLY(igraph_free, vdata);
1239
for (i=0; i<no_of_nodes; i++) {
1242
vdata[i].offset=0.0;
1243
vdata[i].left_contour=-1;
1244
vdata[i].right_contour=-1;
1245
vdata[i].offset_follow_lc=0.0;
1246
vdata[i].offset_follow_rc=0.0;
1248
vdata[root].parent=root;
1249
vdata[root].level=0;
1250
MATRIX(*res, root, 1) = 0;
1252
/* Step 1: assign Y coordinates based on BFS and setup parents vector */
1253
IGRAPH_CHECK(igraph_dqueue_push(&q, root));
1254
IGRAPH_CHECK(igraph_dqueue_push(&q, 0));
1255
while (!igraph_dqueue_empty(&q)) {
1256
long int actnode=igraph_dqueue_pop(&q);
1257
long int actdist=igraph_dqueue_pop(&q);
1258
neis=igraph_adjlist_get(&allneis, actnode);
1259
n=igraph_vector_size(neis);
1260
for (j=0; j<n; j++) {
1261
long int neighbor=VECTOR(*neis)[j];
1262
if (vdata[neighbor].parent >= 0) { continue; }
1263
MATRIX(*res, neighbor, 1)=actdist+1;
1264
IGRAPH_CHECK(igraph_dqueue_push(&q, neighbor));
1265
IGRAPH_CHECK(igraph_dqueue_push(&q, actdist+1));
1266
vdata[neighbor].parent = actnode;
1267
vdata[neighbor].level = actdist+1;
1271
/* Step 2: postorder tree traversal, determines the appropriate X
1272
* offsets for every node */
1273
igraph_i_layout_reingold_tilford_postorder(vdata, root, no_of_nodes);
1275
/* Step 3: calculate real coordinates based on X offsets */
1276
igraph_i_layout_reingold_tilford_calc_coords(vdata, res, root, no_of_nodes, vdata[root].offset);
1278
igraph_dqueue_destroy(&q);
1279
igraph_adjlist_destroy(&allneis);
1281
IGRAPH_FINALLY_CLEAN(3);
1283
IGRAPH_PROGRESS("Reingold-Tilford tree layout", 100.0, NULL);
1288
int igraph_i_layout_reingold_tilford_calc_coords(struct igraph_i_reingold_tilford_vertex *vdata,
1289
igraph_matrix_t *res, long int node,
1290
long int vcount, igraph_real_t xpos) {
1292
MATRIX(*res, node, 0) = xpos;
1293
for (i=0, n=0; i<vcount; i++) {
1294
if (i == node) continue;
1295
if (vdata[i].parent == node) {
1296
igraph_i_layout_reingold_tilford_calc_coords(vdata, res, i, vcount,
1297
xpos+vdata[i].offset);
1303
int igraph_i_layout_reingold_tilford_postorder(struct igraph_i_reingold_tilford_vertex *vdata,
1304
long int node, long int vcount) {
1305
long int i, j, childcount, leftroot, leftrootidx, leftleftroot;
1308
/* printf("Starting visiting node %d\n", node); */
1310
/* Check whether this node is a leaf node */
1312
for (i=0; i<vcount; i++) {
1313
if (i == node) continue;
1314
if (vdata[i].parent == node) {
1315
/* Node i is a child, so visit it recursively */
1317
igraph_i_layout_reingold_tilford_postorder(vdata, i, vcount);
1321
if (childcount == 0) return 0;
1323
/* Here we can assume that all of the subtrees have been placed and their
1324
* left and right contours are calculated. Let's place them next to each
1325
* other as close as we can.
1326
* We will take each subtree in an arbitrary order. The root of the
1327
* first one will be placed at offset 0, the next ones will be placed
1328
* as close to each other as possible. leftroot stores the root of the
1329
* rightmost subtree of the already placed subtrees - its right contour
1330
* will be checked against the left contour of the next subtree */
1331
leftleftroot=leftroot=leftrootidx=-1;
1333
/*printf("Visited node %d and arranged its subtrees\n", node);*/
1334
for (i=0, j=0; i<vcount; i++) {
1335
if (i == node) continue;
1336
if (vdata[i].parent == node) {
1337
/*printf(" Placing child %d on level %d\n", i, vdata[i].level);*/
1338
if (leftroot >= 0) {
1339
/* Now we will follow the right contour of leftroot and the
1340
* left contour of the subtree rooted at i */
1342
igraph_real_t loffset, roffset, minsep, rootsep;
1343
lnode = leftroot; rnode = i;
1345
rootsep = vdata[leftroot].offset + minsep;
1346
loffset = 0; roffset = minsep;
1347
/*printf(" Contour: [%d, %d], offsets: [%lf, %lf], rootsep: %lf\n",
1348
lnode, rnode, loffset, roffset, rootsep);*/
1349
while ((lnode >= 0) && (rnode >= 0)) {
1350
/* Step to the next level on the right contour of the left subtree */
1351
if (vdata[lnode].right_contour >= 0) {
1352
loffset += vdata[lnode].offset_follow_rc;
1353
lnode = vdata[lnode].right_contour;
1355
/* Left subtree ended there. The right contour of the left subtree
1356
* will continue to the next step on the right subtree. */
1357
if (vdata[rnode].left_contour >= 0) {
1358
/*printf(" Left subtree ended, continuing left subtree's left and right contour on right subtree (node %ld)\n", vdata[rnode].left_contour);*/
1359
vdata[lnode].left_contour = vdata[rnode].left_contour;
1360
vdata[lnode].right_contour = vdata[rnode].left_contour;
1361
vdata[lnode].offset_follow_lc = vdata[lnode].offset_follow_rc =
1362
(roffset-loffset)+vdata[rnode].offset_follow_lc;
1363
/*printf(" vdata[lnode].offset_follow_* = %.4f\n", vdata[lnode].offset_follow_lc);*/
1367
/* Step to the next level on the left contour of the right subtree */
1368
if (vdata[rnode].left_contour >= 0) {
1369
roffset += vdata[rnode].offset_follow_lc;
1370
rnode = vdata[rnode].left_contour;
1372
/* Right subtree ended here. The left contour of the right
1373
* subtree will continue to the next step on the left subtree.
1374
* Note that lnode has already been advanced here */
1376
/*printf(" Right subtree ended, continuing right subtree's left and right contour on left subtree (node %ld)\n", lnode);*/
1377
vdata[rnode].left_contour = lnode;
1378
vdata[rnode].right_contour = lnode;
1379
vdata[rnode].offset_follow_lc = vdata[rnode].offset_follow_rc =
1380
(loffset-roffset); /* loffset has also been increased earlier */
1381
/*printf(" vdata[rnode].offset_follow_* = %.4f\n", vdata[rnode].offset_follow_lc);*/
1385
/*printf(" Contour: [%d, %d], offsets: [%lf, %lf], rootsep: %lf\n",
1386
lnode, rnode, loffset, roffset, rootsep);*/
1388
/* Push subtrees away if necessary */
1389
if ((lnode >= 0) && (rnode >= 0) && (roffset - loffset < minsep)) {
1390
/*printf(" Pushing right subtree away by %lf\n", minsep-roffset+loffset);*/
1391
rootsep += minsep-roffset+loffset;
1392
roffset = loffset+minsep;
1396
/*printf(" Offset of subtree with root node %d will be %lf\n", i, rootsep);*/
1397
vdata[i].offset = rootsep;
1398
vdata[node].right_contour = i;
1399
vdata[node].offset_follow_rc = rootsep;
1400
avg = (avg*j)/(j+1) + rootsep/(j+1);
1406
vdata[node].left_contour=i;
1407
vdata[node].right_contour=i;
1408
vdata[node].offset_follow_lc = 0.0;
1409
vdata[node].offset_follow_rc = 0.0;
1410
avg = vdata[i].offset;
1415
/*printf("Shifting node to be centered above children. Shift amount: %lf\n", avg);*/
1416
vdata[node].offset_follow_lc -= avg;
1417
vdata[node].offset_follow_rc -= avg;
1418
for (i=0, j=0; i<vcount; i++) {
1419
if (i == node) continue;
1420
if (vdata[i].parent == node) vdata[i].offset -= avg;
1427
* \function igraph_layout_reingold_tilford_circular
1428
* \brief Circular Reingold-Tilford layout for trees
1431
* This layout is almost the same as \ref igraph_layout_reingold_tilford(), but
1432
* the tree is drawn in a circular way, with the root vertex in the center.
1434
* \param graph The graph object.
1435
* \param res The result, the coordinates in a matrix. The parameter
1436
* should point to an initialized matrix object and will be resized.
1437
* \param root The index of the root vertex.
1438
* \return Error code.
1440
* \sa \ref igraph_layout_reingold_tilford().
1443
int igraph_layout_reingold_tilford_circular(const igraph_t *graph,
1444
igraph_matrix_t *res, long int root) {
1446
long int no_of_nodes=igraph_vcount(graph);
1448
igraph_real_t ratio=2*M_PI*(no_of_nodes-1.0)/no_of_nodes;
1449
igraph_real_t minx, maxx;
1451
IGRAPH_CHECK(igraph_layout_reingold_tilford(graph, res, root));
1453
if (no_of_nodes == 0) return 0;
1455
minx = maxx = MATRIX(*res, 0, 0);
1456
for (i=1; i<no_of_nodes; i++) {
1457
if (MATRIX(*res, i, 0) > maxx) maxx=MATRIX(*res, i, 0);
1458
if (MATRIX(*res, i, 0) < minx) minx=MATRIX(*res, i, 0);
1460
ratio /= (maxx-minx);
1461
for (i=0; i<no_of_nodes; i++) {
1462
igraph_real_t phi=(MATRIX(*res, i, 0)-minx)*ratio;
1463
igraph_real_t r=MATRIX(*res, i, 1);
1464
MATRIX(*res, i, 0) = r*cos(phi);
1465
MATRIX(*res, i, 1) = r*sin(phi);
1471
#define COULOMBS_CONSTANT 8987500000.0
1473
igraph_real_t igraph_i_distance_between(const igraph_matrix_t *c, long int a,
1475
igraph_real_t diffx=MATRIX(*c, a, 0)-MATRIX(*c, b, 0);
1476
igraph_real_t diffy=MATRIX(*c, a, 1)-MATRIX(*c, b, 1);
1477
return sqrt( diffx*diffx + diffy*diffy );
1480
int igraph_i_determine_electric_axal_forces(const igraph_matrix_t *pos,
1483
igraph_real_t directed_force,
1484
igraph_real_t distance,
1485
long int other_node,
1486
long int this_node) {
1488
// We know what the directed force is. We now need to translate it
1489
// into the appropriate x and y componenets.
1493
// directed_force / |
1498
// other_node.x > this_node.x
1499
// other_node.y > this_node.y
1500
// the force will be on this_node away from other_node
1502
// the proportion (distance/y_distance) is equal to the proportion
1503
// (directed_force/y_force), as the two triangles are similar.
1504
// therefore, the magnitude of y_force = (directed_force*y_distance)/distance
1505
// the sign of y_force is negative, away from other_node
1507
igraph_real_t x_distance, y_distance;
1508
y_distance = MATRIX(*pos, other_node, 1)-MATRIX(*pos, this_node, 1);
1509
if (y_distance < 0) { y_distance = -y_distance; }
1510
*y = -1 * ((directed_force * y_distance) / distance);
1512
// the x component works in exactly the same way.
1513
x_distance = MATRIX(*pos, other_node, 0)-MATRIX(*pos, this_node, 0);
1514
if (x_distance < 0) { x_distance = -x_distance; }
1515
*x = -1 * ((directed_force * x_distance) / distance);
1517
// Now we need to reverse the polarity of our answers based on the falsness
1518
// of our assumptions.
1519
if (MATRIX(*pos, other_node, 0) < MATRIX(*pos, this_node, 0)) {
1522
if (MATRIX(*pos, other_node, 1) < MATRIX(*pos, this_node, 1)) {
1529
int igraph_i_apply_electrical_force(const igraph_matrix_t *pos,
1530
igraph_vector_t *pending_forces_x,
1531
igraph_vector_t *pending_forces_y,
1532
long int other_node, long int this_node,
1533
igraph_real_t node_charge,
1534
igraph_real_t distance) {
1536
igraph_real_t directed_force = COULOMBS_CONSTANT *
1537
((node_charge * node_charge)/(distance * distance));
1539
igraph_real_t x_force, y_force;
1540
igraph_i_determine_electric_axal_forces(pos, &x_force, &y_force,
1541
directed_force, distance,
1542
other_node, this_node);
1544
VECTOR(*pending_forces_x)[this_node] += x_force;
1545
VECTOR(*pending_forces_y)[this_node] += y_force;
1546
VECTOR(*pending_forces_x)[other_node] -= x_force;
1547
VECTOR(*pending_forces_y)[other_node] -= y_force;
1552
int igraph_i_determine_spring_axal_forces(const igraph_matrix_t *pos,
1553
igraph_real_t *x, igraph_real_t *y,
1554
igraph_real_t directed_force,
1555
igraph_real_t distance,
1557
long int other_node, long int this_node) {
1559
// if the spring is just the right size, the forces will be 0, so we can
1560
// skip the computation.
1562
// if the spring is too long, our forces will be identical to those computed
1563
// by determine_electrical_axal_forces() (this_node will be pulled toward
1566
// if the spring is too short, our forces will be the opposite of those
1567
// computed by determine_electrical_axal_forces() (this_node will be pushed
1568
// away from other_node)
1570
// finally, since both nodes are movable, only one-half of the total force
1571
// should be applied to each node, so half the forces for our answer.
1573
if (distance == spring_length) {
1577
igraph_i_determine_electric_axal_forces(pos, x, y, directed_force, distance,
1578
other_node, this_node);
1579
if (distance < spring_length) {
1590
int igraph_i_apply_spring_force(const igraph_matrix_t *pos,
1591
igraph_vector_t *pending_forces_x,
1592
igraph_vector_t *pending_forces_y,
1593
long int other_node,
1594
long int this_node, int spring_length,
1595
igraph_real_t spring_constant) {
1597
// determined using Hooke's Law:
1600
// k = spring constant
1601
// x = displacement from ideal length in meters
1603
igraph_real_t distance, displacement, directed_force, x_force, y_force;
1604
distance = igraph_i_distance_between(pos, other_node, this_node);
1605
// let's protect ourselves from division by zero by ignoring two nodes that
1606
// happen to be in the same place. Since we separate all nodes before we
1607
// work on any of them, this will only happen in extremely rare circumstances,
1608
// and when it does, electrical force will probably push one or both of them
1609
// one way or another anyway.
1610
if (distance == 0.0) {
1614
displacement = distance - spring_length;
1615
if (displacement < 0) {
1616
displacement = -displacement;
1618
directed_force = -1 * spring_constant * displacement;
1619
// remember, this is force directed away from the spring;
1620
// a negative number is back towards the spring (or, in our case, back towards
1623
// get the force that should be applied to >this< node
1624
igraph_i_determine_spring_axal_forces(pos, &x_force, &y_force,
1625
directed_force, distance, spring_length,
1626
other_node, this_node);
1628
VECTOR(*pending_forces_x)[this_node] += x_force;
1629
VECTOR(*pending_forces_y)[this_node] += y_force;
1630
VECTOR(*pending_forces_x)[other_node] -= x_force;
1631
VECTOR(*pending_forces_y)[other_node] -= y_force;
1636
int igraph_i_move_nodes(igraph_matrix_t *pos,
1637
const igraph_vector_t *pending_forces_x,
1638
const igraph_vector_t *pending_forces_y,
1639
igraph_real_t node_mass,
1640
igraph_real_t max_sa_movement) {
1642
// Since each iteration is isolated, time is constant at 1.
1644
// Force effects acceleration.
1645
// acceleration (d(velocity)/time) = velocity
1646
// velocity (d(displacement)/time) = displacement
1647
// displacement = acceleration
1649
// determined using Newton's second law:
1652
// acceleration = force / mass
1653
// velocity = force / mass
1654
// displacement = force / mass
1656
long int this_node, no_of_nodes=igraph_vector_size(pending_forces_x);
1658
for (this_node=0; this_node < no_of_nodes; this_node++) {
1660
igraph_real_t x_movement, y_movement;
1662
x_movement = VECTOR(*pending_forces_x)[this_node] / node_mass;
1663
if (x_movement > max_sa_movement) {
1664
x_movement = max_sa_movement;
1665
} else if (x_movement < -max_sa_movement) {
1666
x_movement = -max_sa_movement;
1669
y_movement = VECTOR(*pending_forces_y)[this_node] / node_mass;
1670
if (y_movement > max_sa_movement) {
1671
y_movement = max_sa_movement;
1672
} else if (y_movement < -max_sa_movement) {
1673
y_movement = -max_sa_movement;
1676
MATRIX(*pos, this_node, 0) += x_movement;
1677
MATRIX(*pos, this_node, 1) += y_movement;
1684
* \function igraph_layout_graphopt
1685
* \brief Optimizes vertex layout via the graphopt algorithm.
1688
* This is a port of the graphopt layout algorithm by Michael Schmuhl.
1689
* graphopt version 0.4.1 was rewritten in C and the support for
1690
* layers was removed (might be added later) and a code was a bit
1691
* reorganized to avoid some unneccessary steps is the node charge (see below)
1695
* graphopt uses physical analogies for defining attracting and repelling
1696
* forces among the vertices and then the physical system is simulated
1697
* until it reaches an equilibrium. (There is no simulated annealing or
1698
* anything like that, so a stable fixed point is not guaranteed.)
1701
* See also http://www.schmuhl.org/graphopt/ for the original graphopt.
1702
* \param graph The input graph.
1703
* \param res Pointer to an initialized matrix, the result will be stored here
1704
* and its initial contents is used the starting point of the simulation
1705
* if the \p use_seed argument is true. Note that in this case the
1706
* matrix should have the proper size, otherwise a warning is issued and
1707
* the supplied values are ignored. If no starting positions are given
1708
* (or they are invalid) then a random staring position is used.
1709
* The matrix will be resized if needed.
1710
* \param niter Integer constant, the number of iterations to perform.
1711
* Should be a couple of hundred in general. If you have a large graph
1712
* then you might want to only do a few iterations and then check the
1713
* result. If it is not good enough you can feed it in again in
1714
* the \p res argument. The original graphopt default if 500.
1715
* \param node_charge The charge of the vertices, used to calculate electric
1716
* repulsion. The original graphopt default is 0.001.
1717
* \param node_mass The mass of the vertices, used for the spring forces.
1718
* The original graphopt defaults to 30.
1719
* \param spring_length The length of the springs, an integer number.
1720
* The original graphopt defaults to zero.
1721
* \param spring_constant The spring constant, the original graphopt defaults
1723
* \param max_sa_movement Real constant, it gives the maximum amount of movement
1724
* allowed in a single step along a single axis. The original graphopt
1726
* \param use_seed Logical scalar, whether to use the positions in \p res as
1727
* a starting configuration. See also \p res above.
1728
* \return Error code.
1730
* Time complexity: O(n (|V|^2+|E|) ), n is the number of iterations,
1731
* |V| is the number of vertices, |E| the number
1732
* of edges. If \p node_charge is zero then it is only O(n|E|).
1735
int igraph_layout_graphopt(const igraph_t *graph, igraph_matrix_t *res,
1736
igraph_integer_t niter,
1737
igraph_real_t node_charge, igraph_real_t node_mass,
1738
igraph_integer_t spring_length,
1739
igraph_real_t spring_constant,
1740
igraph_real_t max_sa_movement,
1741
igraph_bool_t use_seed) {
1743
long int no_of_nodes=igraph_vcount(graph);
1744
long int no_of_edges=igraph_ecount(graph);
1745
int my_spring_length=spring_length;
1746
igraph_vector_t pending_forces_x, pending_forces_y;
1747
/* Set a flag to calculate (or not) the electrical forces that the nodes */
1748
/* apply on each other based on if both node types' charges are zero. */
1749
igraph_bool_t apply_electric_charges= (node_charge!=0);
1751
long int this_node, other_node, edge;
1752
igraph_real_t distance;
1755
IGRAPH_VECTOR_INIT_FINALLY(&pending_forces_x, no_of_nodes);
1756
IGRAPH_VECTOR_INIT_FINALLY(&pending_forces_y, no_of_nodes);
1759
if (igraph_matrix_nrow(res) != no_of_nodes ||
1760
igraph_matrix_ncol(res) != 2) {
1761
IGRAPH_WARNING("Invalid size for initial matrix, starting from random layout");
1762
IGRAPH_CHECK(igraph_layout_random(graph, res));
1765
IGRAPH_CHECK(igraph_layout_random(graph, res));
1768
IGRAPH_PROGRESS("Graphopt layout", 0, NULL);
1769
for(i=niter;i>0;i--) {
1770
/* Report progress in approx. every 100th step */
1772
IGRAPH_PROGRESS("Graphopt layout", 100.0-100.0*i/niter, NULL);
1775
/* Clear pending forces on all nodes */
1776
igraph_vector_null(&pending_forces_x);
1777
igraph_vector_null(&pending_forces_y);
1779
// Apply electrical force applied by all other nodes
1780
if (apply_electric_charges) {
1781
// Iterate through all nodes
1782
for (this_node = 0; this_node < no_of_nodes; this_node++) {
1783
IGRAPH_ALLOW_INTERRUPTION();
1784
for (other_node = this_node + 1;
1785
other_node < no_of_nodes;
1787
distance = igraph_i_distance_between(res, this_node, other_node);
1788
// let's protect ourselves from division by zero by ignoring
1789
// two nodes that happen to be in the same place. Since we
1790
// separate all nodes before we work on any of them, this
1791
// will only happen in extremely rare circumstances, and when
1792
// it does, springs will probably pull them apart anyway.
1793
// also, if we are more than 50 away, the electric force
1794
// will be negligable.
1795
// ***** may not always be desirable ****
1796
if ((distance != 0.0) && (distance < 500.0)) {
1797
// if (distance != 0.0) {
1798
// Apply electrical force from node(counter2) on
1800
igraph_i_apply_electrical_force(res, &pending_forces_x,
1802
other_node, this_node,
1810
// Apply force from springs
1811
for (edge = 0; edge < no_of_edges; edge++) {
1812
long int this_node=IGRAPH_FROM(graph, edge);
1813
long int other_node=IGRAPH_TO(graph, edge);
1814
// Apply spring force on both nodes
1815
igraph_i_apply_spring_force(res, &pending_forces_x, &pending_forces_y,
1816
other_node, this_node, my_spring_length,
1820
// Effect the movement of the nodes based on all pending forces
1821
igraph_i_move_nodes(res, &pending_forces_x, &pending_forces_y, node_mass,
1824
IGRAPH_PROGRESS("Graphopt layout", 100, NULL);
1826
igraph_vector_destroy(&pending_forces_y);
1827
igraph_vector_destroy(&pending_forces_x);
1828
IGRAPH_FINALLY_CLEAN(2);
1833
int igraph_i_layout_merge_dla(igraph_i_layout_mergegrid_t *grid,
1834
long int actg, igraph_real_t *x, igraph_real_t *y, igraph_real_t r,
1835
igraph_real_t cx, igraph_real_t cy, igraph_real_t startr,
1836
igraph_real_t killr);
1839
* \function igraph_layout_merge_dla
1840
* \brief Merge multiple layouts by using a DLA algorithm
1843
* First each layout is covered by a circle. Then the layout of the
1844
* largest graph is placed at the origin. Then the other layouts are
1845
* placed by the DLA algorithm, larger ones first and smaller ones
1847
* \param thegraphs Pointer vector containing the graph object of
1848
* which the layouts will be merged.
1849
* \param coords Pointer vector containing matrix objects with the 2d
1850
* layouts of the graphs in \p thegraphs.
1851
* \param res Pointer to an initialized matrix object, the result will
1852
* be stored here. It will be resized if needed.
1853
* \return Error code.
1855
* Added in version 0.2. This function is experimental.
1858
* Time complexity: TODO.
1861
int igraph_layout_merge_dla(igraph_vector_ptr_t *thegraphs,
1862
igraph_vector_ptr_t *coords,
1863
igraph_matrix_t *res) {
1864
long int graphs=igraph_vector_ptr_size(coords);
1865
igraph_vector_t sizes;
1866
igraph_vector_t x, y, r;
1867
igraph_vector_t nx, ny, nr;
1868
long int allnodes=0;
1871
igraph_i_layout_mergegrid_t grid;
1873
igraph_real_t minx, maxx, miny, maxy;
1874
igraph_real_t area=0;
1875
igraph_real_t maxr=0;
1878
IGRAPH_VECTOR_INIT_FINALLY(&sizes, graphs);
1879
IGRAPH_VECTOR_INIT_FINALLY(&x, graphs);
1880
IGRAPH_VECTOR_INIT_FINALLY(&y, graphs);
1881
IGRAPH_VECTOR_INIT_FINALLY(&r, graphs);
1882
IGRAPH_VECTOR_INIT_FINALLY(&nx, graphs);
1883
IGRAPH_VECTOR_INIT_FINALLY(&ny, graphs);
1884
IGRAPH_VECTOR_INIT_FINALLY(&nr, graphs);
1886
for (i=0; i<igraph_vector_ptr_size(coords); i++) {
1887
igraph_matrix_t *mat=VECTOR(*coords)[i];
1888
long int size=igraph_matrix_nrow(mat);
1889
IGRAPH_ALLOW_INTERRUPTION();
1891
VECTOR(sizes)[i]=size;
1892
VECTOR(r)[i]=pow(size, .75);
1893
area+=VECTOR(r)[i] * VECTOR(r)[i];
1894
if (VECTOR(r)[i] > maxr) {
1898
igraph_i_layout_sphere_2d(mat,
1899
igraph_vector_e_ptr(&nx, i),
1900
igraph_vector_e_ptr(&ny, i),
1901
igraph_vector_e_ptr(&nr, i));
1904
igraph_vector_order2(&sizes); /* largest first */
1906
/* 0. create grid */
1907
minx=miny=-sqrt(5*area);
1908
maxx=maxy=sqrt(5*area);
1909
igraph_i_layout_mergegrid_init(&grid, minx, maxx, 200,
1911
IGRAPH_FINALLY(igraph_i_layout_mergegrid_destroy, &grid);
1913
/* fprintf(stderr, "Ok, starting DLA\n"); */
1915
/* 1. place the largest */
1916
actg=VECTOR(sizes)[jpos++];
1917
igraph_i_layout_merge_place_sphere(&grid, 0, 0, VECTOR(r)[actg], actg);
1919
IGRAPH_PROGRESS("Merging layouts via DLA", 0.0, NULL);
1920
while (jpos<graphs) {
1921
IGRAPH_ALLOW_INTERRUPTION();
1922
/* fprintf(stderr, "comp: %li", jpos); */
1923
IGRAPH_PROGRESS("Merging layouts via DLA", (100.0*jpos)/graphs, NULL);
1925
actg=VECTOR(sizes)[jpos++];
1926
/* 2. random walk, TODO: tune parameters */
1927
igraph_i_layout_merge_dla(&grid, actg,
1928
igraph_vector_e_ptr(&x, actg),
1929
igraph_vector_e_ptr(&y, actg),
1930
VECTOR(r)[actg], 0, 0,
1931
maxx-maxr, maxx-maxr+5);
1933
/* 3. place sphere */
1934
igraph_i_layout_merge_place_sphere(&grid, VECTOR(x)[actg], VECTOR(y)[actg],
1935
VECTOR(r)[actg], actg);
1937
IGRAPH_PROGRESS("Merging layouts via DLA", 100.0, NULL);
1939
/* Create the result */
1940
IGRAPH_CHECK(igraph_matrix_resize(res, allnodes, 2));
1942
for (i=0; i<graphs; i++) {
1943
long int size=igraph_matrix_nrow(VECTOR(*coords)[i]);
1944
igraph_real_t xx=VECTOR(x)[i];
1945
igraph_real_t yy=VECTOR(y)[i];
1946
igraph_real_t rr=VECTOR(r)[i]/VECTOR(nr)[i];
1947
igraph_matrix_t *mat=VECTOR(*coords)[i];
1948
IGRAPH_ALLOW_INTERRUPTION();
1949
if (VECTOR(nr)[i]==0) { rr=1; }
1950
for (j=0; j<size; j++) {
1951
MATRIX(*res, respos, 0)=rr*(MATRIX(*mat, j, 0)-VECTOR(nx)[i]);
1952
MATRIX(*res, respos, 1)=rr*(MATRIX(*mat, j, 1)-VECTOR(ny)[i]);
1953
MATRIX(*res, respos, 0)+=xx;
1954
MATRIX(*res, respos, 1)+=yy;
1959
igraph_i_layout_mergegrid_destroy(&grid);
1960
igraph_vector_destroy(&sizes);
1961
igraph_vector_destroy(&x);
1962
igraph_vector_destroy(&y);
1963
igraph_vector_destroy(&r);
1964
igraph_vector_destroy(&nx);
1965
igraph_vector_destroy(&ny);
1966
igraph_vector_destroy(&nr);
1967
IGRAPH_FINALLY_CLEAN(8);
1971
int igraph_i_layout_sphere_2d(igraph_matrix_t *coords, igraph_real_t *x, igraph_real_t *y,
1973
long int nodes=igraph_matrix_nrow(coords);
1975
igraph_real_t xmin, xmax, ymin, ymax;
1977
xmin=xmax=MATRIX(*coords,0,0);
1978
ymin=ymax=MATRIX(*coords,0,1);
1979
for (i=1; i<nodes; i++) {
1981
if (MATRIX(*coords,i,0) < xmin) {
1982
xmin=MATRIX(*coords,i,0);
1983
} else if (MATRIX(*coords,i,0)>xmax) {
1984
xmax=MATRIX(*coords,i,0);
1987
if (MATRIX(*coords,i,1) < ymin) {
1988
ymin=MATRIX(*coords,i,1);
1989
} else if (MATRIX(*coords,i,1)>ymax) {
1990
ymax=MATRIX(*coords,i,1);
1997
*r=sqrt( (xmax-xmin)*(xmax-xmin)+(ymax-ymin)*(ymax-ymin) ) / 2;
2002
int igraph_i_layout_sphere_3d(igraph_matrix_t *coords, igraph_real_t *x, igraph_real_t *y,
2003
igraph_real_t *z, igraph_real_t *r) {
2004
long int nodes=igraph_matrix_nrow(coords);
2006
igraph_real_t xmin, xmax, ymin, ymax, zmin, zmax;
2008
xmin=xmax=MATRIX(*coords,0,0);
2009
ymin=ymax=MATRIX(*coords,0,1);
2010
zmin=zmax=MATRIX(*coords,0,2);
2011
for (i=1; i<nodes; i++) {
2013
if (MATRIX(*coords,i,0) < xmin) {
2014
xmin=MATRIX(*coords,i,0);
2015
} else if (MATRIX(*coords,i,0)>xmax) {
2016
xmax=MATRIX(*coords,i,0);
2019
if (MATRIX(*coords,i,1) < ymin) {
2020
ymin=MATRIX(*coords,i,1);
2021
} else if (MATRIX(*coords,i,1)>ymax) {
2022
ymax=MATRIX(*coords,i,1);
2025
if (MATRIX(*coords,i,2) < zmin) {
2026
zmin=MATRIX(*coords,i,2);
2027
} else if (MATRIX(*coords,i,2)>zmax) {
2028
zmax=MATRIX(*coords,i,2);
2036
*r=sqrt( (xmax-xmin)*(xmax-xmin)+(ymax-ymin)*(ymax-ymin)+
2037
(zmax-zmin)*(zmax-zmin) ) / 2;
2042
#define DIST(x,y) (sqrt(pow((x)-cx,2)+pow((y)-cy,2)))
2044
int igraph_i_layout_merge_dla(igraph_i_layout_mergegrid_t *grid,
2045
long int actg, igraph_real_t *x, igraph_real_t *y, igraph_real_t r,
2046
igraph_real_t cx, igraph_real_t cy, igraph_real_t startr,
2047
igraph_real_t killr) {
2049
igraph_real_t angle, len;
2053
/* start particle */
2056
angle=RNG_UNIF(0,2*M_PI);
2057
len=RNG_UNIF(.5*startr, startr);
2058
*x=cx+len*cos(angle);
2059
*y=cy+len*sin(angle);
2060
sp=igraph_i_layout_mergegrid_get_sphere(grid, *x, *y, r);
2063
while (sp < 0 && DIST(*x,*y)<killr) {
2064
igraph_real_t nx, ny;
2066
angle=RNG_UNIF(0,2*M_PI);
2067
len=RNG_UNIF(0, startr/100);
2068
nx= *x + len * cos(angle);
2069
ny= *y + len * sin(angle);
2070
sp=igraph_i_layout_mergegrid_get_sphere(grid, nx, ny, r);
2077
/* fprintf(stderr, "%li ", steps); */
added
removed
src/layout.c
