~ubuntu-branches/ubuntu/trusty/igraph/trusty-proposed

« back to all changes in this revision

Viewing changes to src/decomposition.c

  • Committer: Package Import Robot
  • Author(s): Mathieu Malaterre
  • Date: 2013-05-27 14:01:54 UTC
  • mfrom: (4.1.1 experimental)
  • Revision ID: package-import@ubuntu.com-20130527140154-oxwwmr0gj3kdy4ol
Tags: 0.6.5-2
Upload to sid

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
/* -*- mode: C -*-  */
 
2
/* 
 
3
   IGraph library.
 
4
   Copyright (C) 2008-2012  Gabor Csardi <csardi.gabor@gmail.com>
 
5
   334 Harvard street, Cambridge, MA 02139 USA
 
6
   
 
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.
 
11
   
 
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.
 
16
   
 
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 
 
20
   02110-1301 USA
 
21
 
 
22
*/
 
23
 
 
24
#include "igraph_structural.h"
 
25
#include "igraph_error.h"
 
26
#include "igraph_adjlist.h"
 
27
#include "igraph_interface.h"
 
28
 
 
29
/**
 
30
 * \function igraph_maximum_cardinality_search
 
31
 * Maximum cardinality search
 
32
 * 
 
33
 * This function implements the maximum cardinality search algorithm
 
34
 * discussed in 
 
35
 * Robert E Tarjan and Mihalis Yannakakis: Simple linear-time
 
36
 * algorithms to test chordality of graphs, test acyclicity of
 
37
 * hypergraphs, and selectively reduce acyclic hypergraphs.
 
38
 * SIAM Journal of Computation 13, 566--579, 1984.
 
39
 * 
 
40
 * \param graph The input graph. Can be directed, but the direction
 
41
 *   of the edges is ignored.
 
42
 * \param alpha Pointer to an initialized vector, the result is stored here. 
 
43
 *   It will be resized, as needed. Upon return it contains
 
44
 *   the rank of the each vertex.
 
45
 * \param alpham1 Pointer to an initialized vector or a \c NULL
 
46
 *   pointer. If not \c NULL, then the inverse of \p alpha is stored
 
47
 *   here.
 
48
 * \return Error code.
 
49
 * 
 
50
 * Time complexity: O(|V|+|E|), linear in terms of the number of
 
51
 * vertices and edges.  
 
52
 * 
 
53
 * \sa \ref igraph_is_chordal().
 
54
 */
 
55
 
 
56
int igraph_maximum_cardinality_search(const igraph_t *graph,
 
57
                                      igraph_vector_t *alpha,
 
58
                                      igraph_vector_t *alpham1) {
 
59
  
 
60
  long int no_of_nodes=igraph_vcount(graph);
 
61
  igraph_vector_long_t size;
 
62
  igraph_vector_long_t head, next, prev; /* doubly linked list with head */
 
63
  long int i;
 
64
  igraph_adjlist_t adjlist;
 
65
  
 
66
  /***************/
 
67
  /* local j, v; */
 
68
  /***************/
 
69
  
 
70
  long int j, v;
 
71
 
 
72
  IGRAPH_CHECK(igraph_vector_long_init(&size, no_of_nodes));
 
73
  IGRAPH_FINALLY(igraph_vector_long_destroy, &size);
 
74
  IGRAPH_CHECK(igraph_vector_long_init(&head, no_of_nodes));
 
75
  IGRAPH_FINALLY(igraph_vector_long_destroy, &head);
 
76
  IGRAPH_CHECK(igraph_vector_long_init(&next, no_of_nodes));
 
77
  IGRAPH_FINALLY(igraph_vector_long_destroy, &next);
 
78
  IGRAPH_CHECK(igraph_vector_long_init(&prev, no_of_nodes));
 
79
  IGRAPH_FINALLY(igraph_vector_long_destroy, &prev);
 
80
 
 
81
  IGRAPH_CHECK(igraph_adjlist_init(graph, &adjlist, IGRAPH_ALL));
 
82
  IGRAPH_FINALLY(igraph_adjlist_destroy, &adjlist);
 
83
  
 
84
  IGRAPH_CHECK(igraph_vector_resize(alpha, no_of_nodes));
 
85
  if (alpham1) {
 
86
    IGRAPH_CHECK(igraph_vector_resize(alpham1, no_of_nodes));
 
87
  }
 
88
  
 
89
  /***********************************************/
 
90
  /* for i in [0,n-1] -> set(i) := emptyset rof; */
 
91
  /***********************************************/
 
92
 
 
93
  /* nothing to do, 'head' contains all zeros */
 
94
 
 
95
  /*********************************************************/
 
96
  /* for v in vertices -> size(v):=0; add v to set(0) rof; */
 
97
  /*********************************************************/
 
98
 
 
99
  VECTOR(head)[0]=1;
 
100
  for (v=0; v<no_of_nodes; v++) {
 
101
    VECTOR(next)[v]=v+2;
 
102
    VECTOR(prev)[v]=v;
 
103
  }
 
104
  VECTOR(next)[no_of_nodes-1] = 0;
 
105
  /* size is already all zero */
 
106
 
 
107
  /***************/
 
108
  /* i:=n; j:=0; */
 
109
  /***************/
 
110
 
 
111
  i=no_of_nodes; j=0;
 
112
  
 
113
  /**************/
 
114
  /* do i>=1 -> */
 
115
  /**************/
 
116
 
 
117
  while (i>=1) {
 
118
    long int x, k, len;
 
119
    igraph_vector_t *neis;
 
120
 
 
121
    /********************************/
 
122
    /* v :=  delete any from set(j) */
 
123
    /********************************/
 
124
 
 
125
    v=VECTOR(head)[j]-1;
 
126
    x=VECTOR(next)[v];
 
127
    VECTOR(head)[j]=x;
 
128
    if (x != 0) {
 
129
      VECTOR(prev)[x-1]=0;
 
130
    }
 
131
    
 
132
    /*************************************************/
 
133
    /* alpha(v) := i; alpham1(i) := v; size(v) := -1 */
 
134
    /*************************************************/
 
135
 
 
136
    VECTOR(*alpha)[v]=i-1; 
 
137
    if (alpham1) { 
 
138
      VECTOR(*alpham1)[i-1]=v;
 
139
    }
 
140
    VECTOR(size)[v]=-1;
 
141
    
 
142
    /********************************************/
 
143
    /* for {v,w} in E such that size(w) >= 0 -> */
 
144
    /********************************************/
 
145
    
 
146
    neis=igraph_adjlist_get(&adjlist, v);
 
147
    len=igraph_vector_size(neis);
 
148
    for (k=0; k<len; k++) {
 
149
      long int w=VECTOR(*neis)[k];
 
150
      long int ws=VECTOR(size)[w];
 
151
      if (ws >= 0) {
 
152
 
 
153
        /******************************/
 
154
        /* delete w from set(size(w)) */
 
155
        /******************************/
 
156
        
 
157
        long int nw=VECTOR(next)[w];
 
158
        long int pw=VECTOR(prev)[w];
 
159
        if (nw != 0) {
 
160
          VECTOR(prev)[nw-1] = pw;
 
161
        }
 
162
        if (pw != 0) {
 
163
          VECTOR(next)[pw-1] = nw;
 
164
        } else {
 
165
          VECTOR(head)[ws]=nw;
 
166
        }
 
167
 
 
168
        /******************************/
 
169
        /* size(w) := size(w)+1       */
 
170
        /******************************/
 
171
 
 
172
        VECTOR(size)[w] += 1;
 
173
 
 
174
        /******************************/
 
175
        /* add w to set(size(w))      */
 
176
        /******************************/
 
177
 
 
178
        ws=VECTOR(size)[w];
 
179
        nw=VECTOR(head)[ws];
 
180
        VECTOR(next)[w]=nw;
 
181
        VECTOR(prev)[w]=0;
 
182
        if (nw != 0) {
 
183
          VECTOR(prev)[nw-1]=w+1;
 
184
        }
 
185
        VECTOR(head)[ws]=w+1;
 
186
 
 
187
      }
 
188
    }
 
189
    
 
190
    /***********************/
 
191
    /* i := i-1; j := j+1; */
 
192
    /***********************/
 
193
    
 
194
    i -= 1;
 
195
    j += 1;
 
196
    
 
197
    /*********************************************/
 
198
    /* do j>=0 and set(j)=emptyset -> j:=j-1; od */ 
 
199
    /*********************************************/
 
200
    
 
201
    while (j>=0 && VECTOR(head)[j]==0) j--;
 
202
    
 
203
  }
 
204
  
 
205
  igraph_adjlist_destroy(&adjlist);
 
206
  igraph_vector_long_destroy(&prev);
 
207
  igraph_vector_long_destroy(&next);
 
208
  igraph_vector_long_destroy(&head);
 
209
  igraph_vector_long_destroy(&size);
 
210
  IGRAPH_FINALLY_CLEAN(5);
 
211
  
 
212
  return 0;
 
213
}
 
214
 
 
215
/**
 
216
 * \function igraph_is_chordal
 
217
 * Decides whether a graph is chordal
 
218
 * 
 
219
 * A graph is chordal if each of its cycles of four or more nodes
 
220
 * has a chord, which is an edge joining two nodes that are not
 
221
 * adjacent in the cycle. An equivalent definition is that any
 
222
 * chordless cycles have at most three nodes.
 
223
 *
 
224
 * If either \p alpha or \p alpha1 is given, then the other is
 
225
 * calculated by taking simply the inverse. If neither are given,
 
226
 * then \ref igraph_maximum_cardinality_search() is called to calculate
 
227
 * them.
 
228
 * \param graph The input graph, it might be directed, but edge
 
229
 *    direction is ignored.
 
230
 * \param alpha Either an alpha vector coming from
 
231
 *    \ref igraph_maximum_cardinality_search() (on the same graph), or a
 
232
 *    null pointer. 
 
233
 * \param alpham1 Either an inverse alpha vector coming from \ref
 
234
 *    igraph_maximum_cardinality_search() (on the same graph) or a null
 
235
 *    pointer.
 
236
 * \param chordal Pointer to a boolean, the result is stored here.
 
237
 * \param fill_in Pointer to an initialized vector, or a null
 
238
 *    pointer. If not a null pointer, then the fill-in of the graph is
 
239
 *    stored here. The fill-in is the set of edges that are needed to
 
240
 *    make the graph chordal. The vector is resized as needed.
 
241
 * \param newgraph Pointer to an uninitialized graph, or a null
 
242
 *   pointer. If not a null pointer, then a new triangulated graph is
 
243
 *   created here. This essentially means adding the fill-in edges to 
 
244
 *   the original graph.
 
245
 * \return Error code.
 
246
 * 
 
247
 * Time complexity: O(n).
 
248
 * 
 
249
 * \sa \ref igraph_maximum_cardinality_search().
 
250
 */
 
251
 
 
252
int igraph_is_chordal(const igraph_t *graph,
 
253
                      const igraph_vector_t *alpha,
 
254
                      const igraph_vector_t *alpham1,
 
255
                      igraph_bool_t *chordal,
 
256
                      igraph_vector_t *fill_in,
 
257
                      igraph_t *newgraph) {
 
258
  
 
259
  long int no_of_nodes=igraph_vcount(graph);
 
260
  const igraph_vector_t *my_alpha=alpha, *my_alpham1=alpham1;
 
261
  igraph_vector_t v_alpha, v_alpham1;
 
262
  igraph_vector_long_t f, index;
 
263
  long int i;
 
264
  igraph_adjlist_t adjlist;
 
265
  igraph_vector_long_t mark;
 
266
  igraph_bool_t calc_edges= fill_in || newgraph;
 
267
  igraph_vector_t *my_fill_in=fill_in, v_fill_in;
 
268
 
 
269
  /*****************/
 
270
  /* local v, w, x */
 
271
  /*****************/
 
272
 
 
273
  long int v, w, x;
 
274
 
 
275
  if (!chordal && !calc_edges) {
 
276
    /* Nothing to calculate */
 
277
    return 0;
 
278
  }
 
279
  
 
280
  if (!alpha && !alpham1) {
 
281
    IGRAPH_VECTOR_INIT_FINALLY(&v_alpha, no_of_nodes);
 
282
    my_alpha=&v_alpha;
 
283
    IGRAPH_VECTOR_INIT_FINALLY(&v_alpham1, no_of_nodes);
 
284
    my_alpham1=&v_alpham1;
 
285
    IGRAPH_CHECK(igraph_maximum_cardinality_search(graph, 
 
286
                                            (igraph_vector_t*) my_alpha, 
 
287
                                            (igraph_vector_t*) my_alpham1));
 
288
  } else if (alpha && !alpham1) {
 
289
    long int v;
 
290
    IGRAPH_VECTOR_INIT_FINALLY(&v_alpham1, no_of_nodes);
 
291
    my_alpham1=&v_alpham1;
 
292
    for (v=0; v<no_of_nodes; v++) {
 
293
      long int i=VECTOR(*my_alpha)[v];
 
294
      VECTOR(*my_alpham1)[i]=v;
 
295
    }
 
296
  } else if (!alpha && alpham1) { 
 
297
    long int i;
 
298
    IGRAPH_VECTOR_INIT_FINALLY(&v_alpha, no_of_nodes);
 
299
    my_alpha=&v_alpha;
 
300
    for (i=0; i<no_of_nodes; i++) {
 
301
      long int v=VECTOR(*my_alpham1)[i];
 
302
      VECTOR(*my_alpha)[v]=i;
 
303
    }
 
304
  }
 
305
 
 
306
  if (!fill_in && newgraph) {
 
307
    IGRAPH_VECTOR_INIT_FINALLY(&v_fill_in, 0);
 
308
    my_fill_in=&v_fill_in;
 
309
  }
 
310
  
 
311
  IGRAPH_CHECK(igraph_vector_long_init(&f, no_of_nodes));
 
312
  IGRAPH_FINALLY(igraph_vector_long_destroy, &f);
 
313
  IGRAPH_CHECK(igraph_vector_long_init(&index, no_of_nodes));
 
314
  IGRAPH_FINALLY(igraph_vector_long_destroy, &index);
 
315
  IGRAPH_CHECK(igraph_adjlist_init(graph, &adjlist, IGRAPH_ALL));
 
316
  IGRAPH_FINALLY(igraph_adjlist_destroy, &adjlist);
 
317
  IGRAPH_CHECK(igraph_vector_long_init(&mark, no_of_nodes));
 
318
  IGRAPH_FINALLY(igraph_vector_long_destroy, &mark);
 
319
  if (my_fill_in) { igraph_vector_clear(my_fill_in); }
 
320
 
 
321
  if (chordal) { *chordal=1; }
 
322
 
 
323
  /*********************/
 
324
  /* for i in [1,n] -> */
 
325
  /*********************/
 
326
 
 
327
  for (i=0; i<no_of_nodes; i++) {
 
328
    igraph_vector_t *neis;
 
329
    long int j, len;
 
330
    
 
331
    /**********************************************/
 
332
    /* w := alpham1(i); f(w) := w; index(w) := i; */
 
333
    /**********************************************/
 
334
    
 
335
    w=VECTOR(*my_alpham1)[i]; 
 
336
    VECTOR(f)[w]=w;
 
337
    VECTOR(index)[w]=i;
 
338
   
 
339
    /******************************************/
 
340
    /* for {v,w} in E such that alpha(v)<i -> */
 
341
    /******************************************/
 
342
 
 
343
    neis=igraph_adjlist_get(&adjlist, w);
 
344
    len=igraph_vector_size(neis);
 
345
    for (j=0; j<len; j++) {
 
346
      v=VECTOR(*neis)[j];
 
347
      VECTOR(mark)[v] = w+1;
 
348
    }
 
349
    
 
350
    for (j=0; j<len; j++) {
 
351
      v=VECTOR(*neis)[j];
 
352
      if (VECTOR(*my_alpha)[v] >= i) { continue; }
 
353
      
 
354
      /**********/
 
355
      /* x := v */
 
356
      /**********/
 
357
 
 
358
      x=v;
 
359
      
 
360
      /********************/
 
361
      /* do index(x)<i -> */
 
362
      /********************/
 
363
      
 
364
      while (VECTOR(index)[x] < i) {
 
365
        
 
366
        /******************/
 
367
        /* index(x) := i; */
 
368
        /******************/
 
369
        
 
370
        VECTOR(index)[x] = i;
 
371
        
 
372
        /**********************************/
 
373
        /* add {x,w} to E union F(alpha); */
 
374
        /**********************************/
 
375
 
 
376
        if (VECTOR(mark)[x] != w+1) {
 
377
 
 
378
          if (chordal) {
 
379
            *chordal=0;
 
380
          }
 
381
          
 
382
          if (my_fill_in) {
 
383
            IGRAPH_CHECK(igraph_vector_push_back(my_fill_in, x));
 
384
            IGRAPH_CHECK(igraph_vector_push_back(my_fill_in, w));
 
385
          }
 
386
          
 
387
          if (!calc_edges) { 
 
388
            /* make sure that we exit from all loops */
 
389
            i=no_of_nodes;
 
390
            j=len; 
 
391
            break; 
 
392
          }
 
393
        }
 
394
        
 
395
        /*************/
 
396
        /* x := f(x) */
 
397
        /*************/
 
398
        
 
399
        x=VECTOR(f)[x];
 
400
        
 
401
      } /* while (VECTOR(index)[x] < i) */
 
402
      
 
403
      /*****************************/
 
404
      /* if (f(x)=x -> f(x):=w; fi */
 
405
      /*****************************/
 
406
      
 
407
      if (VECTOR(f)[x] == x) { 
 
408
        VECTOR(f)[x] = w;
 
409
      }
 
410
    }
 
411
  }
 
412
 
 
413
  igraph_vector_long_destroy(&mark);
 
414
  igraph_adjlist_destroy(&adjlist);
 
415
  igraph_vector_long_destroy(&index);
 
416
  igraph_vector_long_destroy(&f);
 
417
  IGRAPH_FINALLY_CLEAN(4);
 
418
 
 
419
  if (newgraph) {
 
420
    IGRAPH_CHECK(igraph_copy(newgraph, graph));
 
421
    IGRAPH_FINALLY(igraph_destroy, newgraph);
 
422
    IGRAPH_CHECK(igraph_add_edges(newgraph, my_fill_in, 0));
 
423
    IGRAPH_FINALLY_CLEAN(1);
 
424
  }
 
425
 
 
426
  if (!fill_in && newgraph) {
 
427
    igraph_vector_destroy(&v_fill_in);
 
428
    IGRAPH_FINALLY_CLEAN(1);
 
429
  }
 
430
 
 
431
  if (!alpha && !alpham1) {
 
432
    igraph_vector_destroy(&v_alpham1);
 
433
    igraph_vector_destroy(&v_alpha);
 
434
    IGRAPH_FINALLY_CLEAN(2);
 
435
  } else if (alpha && !alpham1) {
 
436
    igraph_vector_destroy(&v_alpham1);
 
437
    IGRAPH_FINALLY_CLEAN(1);
 
438
  } else if (!alpha && alpham1) { 
 
439
    igraph_vector_destroy(&v_alpha);
 
440
    IGRAPH_FINALLY_CLEAN(1);
 
441
  }
 
442
 
 
443
  return 0;
 
444
}