source: vis_dev/glu-2.3/src/sparse/matrix.c @ 84

Last change on this file since 84 was 13, checked in by cecile, 13 years ago

library glu 2.3

File size: 11.5 KB
Line 
1/*
2 * Revision Control Information
3 *
4 * $Id: matrix.c,v 1.4 2002/09/10 00:01:02 fabio Exp $
5 *
6 */
7#include <stdio.h>
8#include "sparse_int.h"
9
10/*
11 *  free-lists are only used if 'FAST_AND_LOOSE' is set; this is because
12 *  we lose the debugging capability of libmm_t which trashes objects when
13 *  they are free'd.  However, FAST_AND_LOOSE is much faster if matrices
14 *  are created and freed frequently.
15 */
16
17#ifdef FAST_AND_LOOSE
18sm_element *sm_element_freelist;
19sm_row *sm_row_freelist;
20sm_col *sm_col_freelist;
21#endif
22
23sm_matrix *
24sm_allocate(void)
25{
26    register sm_matrix *A;
27
28    A = ALLOC(sm_matrix, 1);
29    A->rows = NIL(sm_row *);
30    A->cols = NIL(sm_col *);
31    A->nrows = A->ncols = 0;
32    A->rows_size = A->cols_size = 0;
33    A->first_row = A->last_row = NIL(sm_row);
34    A->first_col = A->last_col = NIL(sm_col);
35    A->user_word = NIL(char);           /* for our user ... */
36    return A;
37}
38
39
40sm_matrix *
41sm_alloc_size(int row, int col)
42{
43    register sm_matrix *A;
44
45    A = sm_alloc();
46    sm_resize(A, row, col);
47    return A;
48}
49
50
51void
52sm_free_space(sm_matrix *A)
53{
54#ifdef FAST_AND_LOOSE
55    register sm_row *prow;
56
57    if (A->first_row != 0) {
58        for(prow = A->first_row; prow != 0; prow = prow->next_row) {
59            /* add the elements to the free list of elements */
60            prow->last_col->next_col = sm_element_freelist;
61            sm_element_freelist = prow->first_col;
62        }
63
64        /* Add the linked list of rows to the row-free-list */
65        A->last_row->next_row = sm_row_freelist;
66        sm_row_freelist = A->first_row;
67
68        /* Add the linked list of cols to the col-free-list */
69        A->last_col->next_col = sm_col_freelist;
70        sm_col_freelist = A->first_col;
71    }
72#else
73    register sm_row *prow, *pnext_row;
74    register sm_col *pcol, *pnext_col;
75
76    for(prow = A->first_row; prow != 0; prow = pnext_row) {
77        pnext_row = prow->next_row;
78        sm_row_free(prow);
79    }
80    for(pcol = A->first_col; pcol != 0; pcol = pnext_col) {
81        pnext_col = pcol->next_col;
82        pcol->first_row = pcol->last_row = NIL(sm_element);
83        sm_col_free(pcol);
84    }
85#endif
86
87    /* Free the arrays to map row/col numbers into pointers */
88    FREE(A->rows);
89    FREE(A->cols);
90    FREE(A);
91}
92
93
94sm_matrix *
95sm_dup(sm_matrix *A)
96{
97    register sm_row *prow;
98    register sm_element *p;
99    register sm_matrix *B;
100
101    B = sm_alloc();
102    if (A->last_row != 0) {
103        sm_resize(B, A->last_row->row_num, A->last_col->col_num);
104        for(prow = A->first_row; prow != 0; prow = prow->next_row) {
105            for(p = prow->first_col; p != 0; p = p->next_col) {
106                (void) sm_insert(B, p->row_num, p->col_num);
107            }
108        }
109    }
110    return B;
111}
112
113
114void 
115sm_resize(sm_matrix *A, int row, int col)
116{
117    register int i, new_size;
118
119    if (row >= A->rows_size) {
120        new_size = MAX(A->rows_size*2, row+1);
121        A->rows = REALLOC(sm_row *, A->rows, new_size);
122        for(i = A->rows_size; i < new_size; i++) {
123            A->rows[i] = NIL(sm_row);
124        }
125        A->rows_size = new_size;
126    }
127
128    if (col >= A->cols_size) {
129        new_size = MAX(A->cols_size*2, col+1);
130        A->cols = REALLOC(sm_col *, A->cols, new_size);
131        for(i = A->cols_size; i < new_size; i++) {
132            A->cols[i] = NIL(sm_col);
133        }
134        A->cols_size = new_size;
135    }
136}
137
138
139/* 
140 *  insert -- insert a value into the matrix
141 */
142sm_element *
143sm_insert(sm_matrix *A, int row, int col)
144{
145    register sm_row *prow;
146    register sm_col *pcol;
147    register sm_element *element;
148    sm_element *save_element;
149
150    if (row >= A->rows_size || col >= A->cols_size) {
151        sm_resize(A, row, col);
152    }
153
154    prow = A->rows[row];
155    if (prow == NIL(sm_row)) {
156        prow = A->rows[row] = sm_row_alloc();
157        prow->row_num = row;
158        sorted_insert(sm_row, A->first_row, A->last_row, A->nrows, 
159                        next_row, prev_row, row_num, row, prow);
160    }
161
162    pcol = A->cols[col];
163    if (pcol == NIL(sm_col)) {
164        pcol = A->cols[col] = sm_col_alloc();
165        pcol->col_num = col;
166        sorted_insert(sm_col, A->first_col, A->last_col, A->ncols, 
167                        next_col, prev_col, col_num, col, pcol);
168    }
169
170    /* get a new item, save its address */
171    sm_element_alloc(element);
172    save_element = element;
173
174    /* insert it into the row list */
175    sorted_insert(sm_element, prow->first_col, prow->last_col, 
176                prow->length, next_col, prev_col, col_num, col, element);
177
178    /* if it was used, also insert it into the column list */
179    if (element == save_element) {
180        sorted_insert(sm_element, pcol->first_row, pcol->last_row, 
181                pcol->length, next_row, prev_row, row_num, row, element);
182    } else {
183        /* otherwise, it was already in matrix -- free element we allocated */
184        sm_element_free(save_element);
185    }
186    return element;
187}
188
189
190sm_element *
191sm_find(sm_matrix *A, int rownum, int colnum)
192{
193    sm_row *prow;
194    sm_col *pcol;
195
196    prow = sm_get_row(A, rownum);
197    if (prow == NIL(sm_row)) {
198        return NIL(sm_element);
199    } else {
200        pcol = sm_get_col(A, colnum);
201        if (pcol == NIL(sm_col)) {
202            return NIL(sm_element);
203        }
204        if (prow->length < pcol->length) {
205            return sm_row_find(prow, colnum);
206        } else {
207            return sm_col_find(pcol, rownum);
208        }
209    }
210}
211
212
213void
214sm_remove(sm_matrix *A, int rownum, int colnum)
215{
216    sm_remove_element(A, sm_find(A, rownum, colnum));
217}
218
219
220
221void
222sm_remove_element(sm_matrix *A, sm_element *p)
223{
224    register sm_row *prow;
225    register sm_col *pcol;
226
227    if (p == 0) return;
228
229    /* Unlink the element from its row */
230    prow = sm_get_row(A, p->row_num);
231    dll_unlink(p, prow->first_col, prow->last_col, 
232                        next_col, prev_col, prow->length);
233
234    /* if no more elements in the row, discard the row header */
235    if (prow->first_col == NIL(sm_element)) {
236        sm_delrow(A, p->row_num);
237    }
238
239    /* Unlink the element from its column */
240    pcol = sm_get_col(A, p->col_num);
241    dll_unlink(p, pcol->first_row, pcol->last_row, 
242                        next_row, prev_row, pcol->length);
243
244    /* if no more elements in the column, discard the column header */
245    if (pcol->first_row == NIL(sm_element)) {
246        sm_delcol(A, p->col_num);
247    }
248
249    sm_element_free(p);
250}
251
252void 
253sm_delrow(sm_matrix *A, int i)
254{
255    register sm_element *p, *pnext;
256    sm_col *pcol;
257    sm_row *prow;
258
259    prow = sm_get_row(A, i);
260    if (prow != NIL(sm_row)) {
261        /* walk across the row */
262        for(p = prow->first_col; p != 0; p = pnext) {
263            pnext = p->next_col;
264
265            /* unlink the item from the column (and delete it) */
266            pcol = sm_get_col(A, p->col_num);
267            sm_col_remove_element(pcol, p);
268
269            /* discard the column if it is now empty */
270            if (pcol->first_row == NIL(sm_element)) {
271                sm_delcol(A, pcol->col_num);
272            }
273        }
274
275        /* discard the row -- we already threw away the elements */ 
276        A->rows[i] = NIL(sm_row);
277        dll_unlink(prow, A->first_row, A->last_row, 
278                                next_row, prev_row, A->nrows);
279        prow->first_col = prow->last_col = NIL(sm_element);
280        sm_row_free(prow);
281    }
282}
283
284void 
285sm_delcol(sm_matrix *A, int i)
286{
287    register sm_element *p, *pnext;
288    sm_row *prow;
289    sm_col *pcol;
290
291    pcol = sm_get_col(A, i);
292    if (pcol != NIL(sm_col)) {
293        /* walk down the column */
294        for(p = pcol->first_row; p != 0; p = pnext) {
295            pnext = p->next_row;
296
297            /* unlink the element from the row (and delete it) */
298            prow = sm_get_row(A, p->row_num);
299            sm_row_remove_element(prow, p);
300
301            /* discard the row if it is now empty */
302            if (prow->first_col == NIL(sm_element)) {
303                sm_delrow(A, prow->row_num);
304            }
305        }
306
307        /* discard the column -- we already threw away the elements */ 
308        A->cols[i] = NIL(sm_col);
309        dll_unlink(pcol, A->first_col, A->last_col, 
310                            next_col, prev_col, A->ncols);
311        pcol->first_row = pcol->last_row = NIL(sm_element);
312        sm_col_free(pcol);
313    }
314}
315
316void
317sm_copy_row(sm_matrix *dest, int dest_row, sm_row *prow)
318{
319    register sm_element *p;
320
321    for(p = prow->first_col; p != 0; p = p->next_col) {
322        (void) sm_insert(dest, dest_row, p->col_num);
323    }
324}
325
326
327void
328sm_copy_col(sm_matrix *dest, int dest_col, sm_col *pcol)
329{
330    register sm_element *p;
331
332    for(p = pcol->first_row; p != 0; p = p->next_row) {
333        (void) sm_insert(dest, dest_col, p->row_num);
334    }
335}
336
337
338sm_row *
339sm_longest_row(sm_matrix *A)
340{
341    register sm_row *large_row, *prow;
342    register int max_length;
343
344    max_length = 0;
345    large_row = NIL(sm_row);
346    for(prow = A->first_row; prow != 0; prow = prow->next_row) {
347        if (prow->length > max_length) {
348            max_length = prow->length;
349            large_row = prow;
350        }
351    }
352    return large_row;
353}
354
355
356sm_col *
357sm_longest_col(sm_matrix *A)
358{
359    register sm_col *large_col, *pcol;
360    register int max_length;
361
362    max_length = 0;
363    large_col = NIL(sm_col);
364    for(pcol = A->first_col; pcol != 0; pcol = pcol->next_col) {
365        if (pcol->length > max_length) {
366            max_length = pcol->length;
367            large_col = pcol;
368        }
369    }
370    return large_col;
371}
372
373int
374sm_num_elements(sm_matrix *A)
375{
376    register sm_row *prow;
377    register int num;
378
379    num = 0;
380    sm_foreach_row(A, prow) {
381        num += prow->length;
382    }
383    return num;
384}
385
386int 
387sm_read(FILE *fp, sm_matrix **A)
388{
389    int i, j, err;
390
391    *A = sm_alloc();
392    while (! feof(fp)) {
393        err = fscanf(fp, "%d %d", &i, &j);
394        if (err == EOF) {
395            return 1;
396        } else if (err != 2) {
397            return 0;
398        }
399        (void) sm_insert(*A, i, j);
400    }
401    return 1;
402}
403
404
405int 
406sm_read_compressed(FILE *fp, sm_matrix **A)
407{
408    int i, j, k, nrows, ncols;
409    unsigned long x;
410
411    *A = sm_alloc();
412    if (fscanf(fp, "%d %d", &nrows, &ncols) != 2) {
413        return 0;
414    }
415    sm_resize(*A, nrows, ncols);
416
417    for(i = 0; i < nrows; i++) {
418        if (fscanf(fp, "%lx", &x) != 1) {
419            return 0;
420        }
421        for(j = 0; j < ncols; j += 32) {
422            if (fscanf(fp, "%lx", &x) != 1) {
423                return 0;
424            }
425            for(k = j; x != 0; x >>= 1, k++) {
426                if (x & 1) {
427                    (void) sm_insert(*A, i, k);
428                }
429            }
430        }
431    }
432    return 1;
433}
434
435
436void 
437sm_write(FILE *fp, sm_matrix *A)
438{
439    register sm_row *prow;
440    register sm_element *p;
441
442    for(prow = A->first_row; prow != 0; prow = prow->next_row) {
443        for(p = prow->first_col; p != 0; p = p->next_col) {
444            (void) fprintf(fp, "%d %d\n", p->row_num, p->col_num);
445        }
446    }
447}
448
449void 
450sm_print(FILE *fp, sm_matrix *A)
451{
452    register sm_row *prow;
453    register sm_col *pcol;
454    int c;
455
456    if (A->last_col->col_num >= 100) {
457        (void) fprintf(fp, "    ");
458        for(pcol = A->first_col; pcol != 0; pcol = pcol->next_col) {
459            (void) fprintf(fp, "%d", (pcol->col_num / 100)%10);
460        }
461        putc('\n', fp);
462    }
463
464    if (A->last_col->col_num >= 10) {
465        (void) fprintf(fp, "    ");
466        for(pcol = A->first_col; pcol != 0; pcol = pcol->next_col) {
467            (void) fprintf(fp, "%d", (pcol->col_num / 10)%10);
468        }
469        putc('\n', fp);
470    }
471
472    (void) fprintf(fp, "    ");
473    for(pcol = A->first_col; pcol != 0; pcol = pcol->next_col) {
474        (void) fprintf(fp, "%d", pcol->col_num % 10);
475    }
476    putc('\n', fp);
477
478    (void) fprintf(fp, "    ");
479    for(pcol = A->first_col; pcol != 0; pcol = pcol->next_col) {
480        (void) fprintf(fp, "-");
481    }
482    putc('\n', fp);
483
484    for(prow = A->first_row; prow != 0; prow = prow->next_row) {
485        (void) fprintf(fp, "%3d:", prow->row_num);
486
487        for(pcol = A->first_col; pcol != 0; pcol = pcol->next_col) {
488            c = sm_row_find(prow, pcol->col_num) ? '1' : '.';
489            putc(c, fp);
490        }
491        putc('\n', fp);
492    }
493}
494
495
496void 
497sm_dump(sm_matrix *A, char *s, int max)
498{
499    FILE *fp = stdout;
500
501    (void) fprintf(fp, "%s %d rows by %d cols\n", s, A->nrows, A->ncols);
502    if (A->nrows < max) {
503        sm_print(fp, A);
504    }
505}
506
507void
508sm_cleanup(void)
509{
510#ifdef FAST_AND_LOOSE
511    register sm_element *p, *pnext;
512    register sm_row *prow, *pnextrow;
513    register sm_col *pcol, *pnextcol;
514
515    for(p = sm_element_freelist; p != 0; p = pnext) {
516        pnext = p->next_col;
517        FREE(p);
518    }
519    sm_element_freelist = 0;
520
521    for(prow = sm_row_freelist; prow != 0; prow = pnextrow) {
522        pnextrow = prow->next_row;
523        FREE(prow);
524    }
525    sm_row_freelist = 0;
526
527    for(pcol = sm_col_freelist; pcol != 0; pcol = pnextcol) {
528        pnextcol = pcol->next_col;
529        FREE(pcol);
530    }
531    sm_col_freelist = 0;
532#endif
533}
Note: See TracBrowser for help on using the repository browser.