Mercurial > octave-nkf
changeset 5438:49ff3dd744ee
[project @ 2005-09-04 12:25:21 by dbateman]
| author | dbateman |
|---|---|
| date | Sun, 04 Sep 2005 12:25:21 +0000 |
| parents | 009606303874 |
| children | ca7ebff13a16 |
| files | liboctave/COLAMD.README liboctave/COLAMD.files liboctave/COLAMD/ChangeLog liboctave/COLAMD/Makefile liboctave/COLAMD/colamd.c liboctave/COLAMD/colamd.h liboctave/COLAMD/colamd.m liboctave/COLAMD/colamd_demo.m liboctave/COLAMD/colamd_example.c liboctave/COLAMD/colamd_example.out liboctave/COLAMD/colamd_test.m liboctave/COLAMD/colamdmex.c liboctave/COLAMD/colamdtestmex.c liboctave/COLAMD/symamd.m liboctave/COLAMD/symamdmex.c liboctave/COLAMD/symamdtestmex.c liboctave/ChangeLog |
| diffstat | 17 files changed, 1169 insertions(+), 1080 deletions(-) [+] |
line wrap: on
line diff
--- a/liboctave/COLAMD.README Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD.README Sun Sep 04 12:25:21 2005 +0000 @@ -1,4 +1,4 @@ -This directory contains an unmodified copy of COLAMD version 2.3 in +This directory contains an unmodified copy of COLAMD version 2.4 in the subdirectory COLAMD. COLAMD was written by Stefan I. Larimore and Timothy A. Davis (davis@cise.ufl.edu), University of Florida. The algorithm was developed in collaboration with John Gilbert, Xerox
--- a/liboctave/COLAMD.files Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD.files Sun Sep 04 12:25:21 2005 +0000 @@ -1,7 +1,9 @@ -COLAMD_SRC := colamd.c +COLAMD_SRC := colamd_global.c colamd.c -COLAMD_OBJ := $(COLAMD_SRC:.c=.o) +COLAMD_VIRTUAL_SRC := colamd_l.c + +COLAMD_OBJ := $(COLAMD_SRC:.c=.o) $(COLAMD_VIRTUAL_SRC:.c=.o) COLAMD_DEP := $(COLAMD_SRC:.c=.d) @@ -9,3 +11,9 @@ COLAMD_EXTRAS := COLAMD.files COLAMD.README +# Special rules for long version of colamd +colamd_l.o : colamd.c + $(CC) -c $(CPPFLAGS) $(ALL_CFLAGS) -DDLONG $< -o $@ + +pic/colamd_l.o : colamd.c + $(CC) -c $(CPPFLAGS) $(CPICFLAGS) $(ALL_CFLAGS) -DDLONG $< -o $@
--- a/liboctave/COLAMD/ChangeLog Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/ChangeLog Sun Sep 04 12:25:21 2005 +0000 @@ -1,3 +1,31 @@ +Changes from Version 2.3 to 2.4 (Aug 30, 2005) + + * Makefile now relies on ../UFconfig/UFconfig.mk + + * changed the dense row/col detection. The meaning of the knobs + has thus changed. + + * added an option to turn off aggressive absorption. It was + always on in versions 2.3 and earlier. + + * added a #define'd version number + + * added a function pointer (colamd_printf) for COLAMD's printing. + + * added a -DNPRINT option, to turn off printing at compile-time. + + * added a check for integer overflow in colamd_recommended + + * minor changes to allow for more simpler 100% test coverage + + * bug fix. If symamd v2.3 fails to allocate its copy of the input + matrix, then it erroneously frees a calloc'd workspace twice. + This bug has no effect on the MATLAB symamd mexFunction, since + mxCalloc terminates the mexFunction if it fails to allocate + memory. Similarly, UMFPACK is not affected because it does not + use symamd. The bug has no affect on the colamd ordering + routine in v2.3. + Changes from Version 2.2 to 2.3 (Sept. 8, 2003) * removed the call to the MATLAB spparms ('spumoni') function. @@ -38,7 +66,7 @@ No bugs were found in version 1.1. These changes merely add new functionality. - * added the COLAMD_RECOMMENDED (nnz, n_row, n_col) macro. + * added the COLAMD_RECOMMENDED (nnz, n_row, n_col) macro. * moved the output statistics, from A, to a separate output argument. The arguments changed for the C-callable routines.
--- a/liboctave/COLAMD/Makefile Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/Makefile Sun Sep 04 12:25:21 2005 +0000 @@ -1,19 +1,46 @@ -colamd_example: colamd_example.c colamd.c colamd.h - cc -O -o colamd_example colamd_example.c colamd.c - - ./colamd_example + +default: libcolamd.a colamd_example colamd_l_example + +include ../UFconfig/UFconfig.mk + +colamd_example: colamd_example.c libcolamd.a + $(CC) $(CFLAGS) -o colamd_example colamd_example.c libcolamd.a -lm + - ./colamd_example > my_colamd_example.out + - diff colamd_example.out my_colamd_example.out -clean: - - rm *.o colamd_example - - rm colamdmex.mex* symamdmex.mex* - - rm colamdtestmex.mex* symamdtestmex.mex* +colamd_l_example: colamd_l_example.c libcolamd.a + $(CC) $(CFLAGS) -o colamd_l_example colamd_l_example.c libcolamd.a -lm + - ./colamd_l_example > my_colamd_l_example.out + - diff colamd_example.out my_colamd_example.out + +purge: distclean + +distclean: clean2 + - $(RM) libcolamd.a + +clean2: clean + - $(RM) *.o *.dll colamd_example colamd_l_example + - $(RM) colamdmex.mex* symamdmex.mex* + - $(RM) colamdtestmex.mex* symamdtestmex.mex* + - $(RM) my_colamd_example.out my_colamd_l_example.out # Compiles the MATLAB-callable routines -matlab: colamdmex.c colamd.c colamd.h - mex -O colamdmex.c colamd.c - mex -O symamdmex.c colamd.c +matlab: colamdmex.c symamdmex.c libcolamd.a + $(MEX) colamdmex.c libcolamd.a + $(MEX) symamdmex.c libcolamd.a # Compiles the extensive test code -test: matlab colamdmex.c colamd.c colamd.h - mex -O colamdtestmex.c colamd.c - mex -O symamdtestmex.c colamd.c +test: matlab colamdtestmex.c symamdtestmex.c libcolamd.a + $(MEX) colamdtestmex.c libcolamd.a + $(MEX) symamdtestmex.c libcolamd.a +# creates libcolamd.a, a C-callable COLAMD library +libcolamd.a: colamd.c colamd_global.c colamd.h + $(CC) $(CFLAGS) -c colamd_global.c + $(CC) $(CFLAGS) -c colamd.c + $(CC) $(CFLAGS) -c colamd.c -DDLONG -o colamd_l.o + $(AR) libcolamd.a colamd.o colamd_l.o colamd_global.o + +ccode: libcolamd.a + +library: libcolamd.a
--- a/liboctave/COLAMD/colamd.c Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/colamd.c Sun Sep 04 12:25:21 2005 +0000 @@ -2,7 +2,8 @@ /* === colamd/symamd - a sparse matrix column ordering algorithm ============ */ /* ========================================================================== */ -/* +/* COLAMD Version 2.4. + colamd: an approximate minimum degree column ordering algorithm, for LU factorization of symmetric or unsymmetric matrices, QR factorization, least squares, interior point methods for @@ -34,14 +35,10 @@ Authors: The authors of the code itself are Stefan I. Larimore and Timothy A. - Davis (davis@cise.ufl.edu), University of Florida. The algorithm was + Davis (davis at cise.ufl.edu), University of Florida. The algorithm was developed in collaboration with John Gilbert, Xerox PARC, and Esmond Ng, Oak Ridge National Laboratory. - Date: - - September 8, 2003. Version 2.3. - Acknowledgements: This work was supported by the National Science Foundation, under @@ -49,17 +46,34 @@ Copyright and License: - Copyright (c) 1998-2003 by the University of Florida. - All Rights Reserved. - - THIS MATERIAL IS PROVIDED AS IS, WITH ABSOLUTELY NO WARRANTY - EXPRESSED OR IMPLIED. ANY USE IS AT YOUR OWN RISK. - - Permission is hereby granted to use, copy, modify, and/or distribute - this program, provided that the Copyright, this License, and the - Availability of the original version is retained on all copies and made - accessible to the end-user of any code or package that includes COLAMD - or any modified version of COLAMD. + Copyright (c) 1998-2005, Timothy A. Davis, All Rights Reserved. + COLAMD is also available under alternate licenses, contact T. Davis + for details. + + This library is free software; you can redistribute it and/or + modify it under the terms of the GNU Lesser General Public + License as published by the Free Software Foundation; either + version 2.1 of the License, or (at your option) any later version. + + This library is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + Lesser General Public License for more details. + + You should have received a copy of the GNU Lesser General Public + License along with this library; if not, write to the Free Software + Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 + USA + + Permission is hereby granted to use or copy this program under the + terms of the GNU LGPL, provided that the Copyright, this License, + and the Availability of the original version is retained on all copies. + User documentation of any code that uses this code or any modified + version of this code must cite the Copyright, this License, the + Availability note, and "Used by permission." Permission to modify + the code and to distribute modified code is granted, provided the + Copyright, this License, and the Availability note are retained, + and a notice that the code was modified is included. Availability: @@ -70,16 +84,31 @@ This is the http://www.cise.ufl.edu/research/sparse/colamd/colamd.c file. It requires the colamd.h file. It is required by the colamdmex.c and symamdmex.c files, for the MATLAB interface to colamd and symamd. + Appears as ACM Algorithm 836. See the ChangeLog file for changes since Version 1.0. + References: + + T. A. Davis, J. R. Gilbert, S. Larimore, E. Ng, An approximate column + minimum degree ordering algorithm, ACM Transactions on Mathematical + Software, vol. 30, no. 3., pp. 353-376, 2004. + + T. A. Davis, J. R. Gilbert, S. Larimore, E. Ng, Algorithm 836: COLAMD, + an approximate column minimum degree ordering algorithm, ACM + Transactions on Mathematical Software, vol. 30, no. 3., pp. 377-380, + 2004. + */ /* ========================================================================== */ /* === Description of user-callable routines ================================ */ /* ========================================================================== */ -/* +/* COLAMD includes both int and long versions of all its routines. The + * description below os fir the int version. For long, all int arguments + * become long. + ---------------------------------------------------------------------------- colamd_recommended: ---------------------------------------------------------------------------- @@ -88,11 +117,7 @@ #include "colamd.h" int colamd_recommended (int nnz, int n_row, int n_col) ; - - or as a C macro - - #include "colamd.h" - Alen = COLAMD_RECOMMENDED (int nnz, int n_row, int n_col) ; + long colamd_l_recommended (long nnz, long n_row, long n_col) ; Purpose: @@ -122,6 +147,7 @@ #include "colamd.h" colamd_set_defaults (double knobs [COLAMD_KNOBS]) ; + colamd_l_set_defaults (double knobs [COLAMD_KNOBS]) ; Purpose: @@ -131,19 +157,26 @@ double knobs [COLAMD_KNOBS] ; Output only. - Colamd: rows with more than (knobs [COLAMD_DENSE_ROW] * n_col) + NOTE: the meaning of the dense row/col knobs has changed in v2.4 + + knobs [0] and knobs [1] control dense row and col detection: + + Colamd: rows with more than + max (16, knobs [COLAMD_DENSE_ROW] * sqrt (n_col)) entries are removed prior to ordering. Columns with more than - (knobs [COLAMD_DENSE_COL] * n_row) entries are removed prior to + max (16, knobs [COLAMD_DENSE_COL] * sqrt (MIN (n_row,n_col))) + entries are removed prior to ordering, and placed last in the output column ordering. Symamd: uses only knobs [COLAMD_DENSE_ROW], which is knobs [0]. - Rows and columns with more than (knobs [COLAMD_DENSE_ROW] * n) + Rows and columns with more than + max (16, knobs [COLAMD_DENSE_ROW] * sqrt (n)) entries are removed prior to ordering, and placed last in the output ordering. COLAMD_DENSE_ROW and COLAMD_DENSE_COL are defined as 0 and 1, respectively, in colamd.h. Default values of these two knobs - are both 0.5. Currently, only knobs [0] and knobs [1] are + are both 10. Currently, only knobs [0] and knobs [1] are used, but future versions may use more knobs. If so, they will be properly set to their defaults by the future version of colamd_set_defaults, so that the code that calls colamd will @@ -151,6 +184,12 @@ colamd_set_defaults, or pass a (double *) NULL pointer as the knobs array to colamd or symamd. + knobs [2]: aggressive absorption + + knobs [COLAMD_AGGRESSIVE] controls whether or not to do + aggressive absorption during the ordering. Default is TRUE. + + ---------------------------------------------------------------------------- colamd: ---------------------------------------------------------------------------- @@ -160,6 +199,8 @@ #include "colamd.h" int colamd (int n_row, int n_col, int Alen, int *A, int *p, double knobs [COLAMD_KNOBS], int stats [COLAMD_STATS]) ; + long colamd_l (long n_row, long n_col, long Alen, long *A, long *p, + double knobs [COLAMD_KNOBS], long stats [COLAMD_STATS]) ; Purpose: @@ -202,7 +243,9 @@ Alen >= COLAMD_RECOMMENDED (nnz, n_row, n_col) - will be sufficient. + will be sufficient. Note: the macro version does not check + for integer overflow, and thus is not recommended. Use + the colamd_recommended routine instead. int A [Alen] ; Input argument, undefined on output. @@ -352,8 +395,8 @@ with default knobs and no output statistics, do the following: #include "colamd.h" - #define ALEN COLAMD_RECOMMENDED (11, 5, 4) - int A [ALEN] = {1, 2, 5, 3, 5, 1, 2, 3, 4, 2, 4} ; + #define ALEN 100 + int A [ALEN] = {0, 1, 4, 2, 4, 0, 1, 2, 3, 1, 3} ; int p [ ] = {0, 3, 5, 9, 11} ; int stats [COLAMD_STATS] ; colamd (5, 4, ALEN, A, p, (double *) NULL, stats) ; @@ -370,6 +413,9 @@ int symamd (int n, int *A, int *p, int *perm, double knobs [COLAMD_KNOBS], int stats [COLAMD_STATS], void (*allocate) (size_t, size_t), void (*release) (void *)) ; + long symamd_l (long n, long *A, long *p, long *perm, + double knobs [COLAMD_KNOBS], long stats [COLAMD_STATS], + void (*allocate) (size_t, size_t), void (*release) (void *)) ; Purpose: @@ -510,13 +556,6 @@ workspace for M or count arrays using the "allocate" routine passed into symamd). - -999 internal error. colamd failed to order the - matrix M, when it should have succeeded. This - indicates a bug. If this (and *only* this) - error code occurs, please contact the authors. - Don't contact the authors if you get any other - error code. - Future versions may return more statistics in the stats array. void * (*allocate) (size_t, size_t) @@ -543,6 +582,7 @@ #include "colamd.h" colamd_report (int stats [COLAMD_STATS]) ; + colamd_l_report (long stats [COLAMD_STATS]) ; Purpose: @@ -563,6 +603,7 @@ #include "colamd.h" symamd_report (int stats [COLAMD_STATS]) ; + symamd_l_report (long stats [COLAMD_STATS]) ; Purpose: @@ -584,7 +625,11 @@ /* Ensure that debugging is turned off: */ #ifndef NDEBUG #define NDEBUG -#endif /* NDEBUG */ +#endif + +/* turn on debugging by uncommenting the following line + #undef NDEBUG +*/ /* Our "scaffolding code" philosophy: In our opinion, well-written library @@ -603,9 +648,7 @@ (3) (gasp!) for actually finding bugs. This code has been heavily tested and "should" be fully functional and bug-free ... but you never know... - To enable debugging, comment out the "#define NDEBUG" above. For a MATLAB - mexFunction, you will also need to modify mexopts.sh to remove the -DNDEBUG - definition. The code will become outrageously slow when debugging is + The code will become outrageously slow when debugging is enabled. To control the level of debugging output, set an environment variable D to 0 (little), 1 (some), 2, 3, or 4 (lots). When debugging, you should see the following message on the standard output: @@ -623,14 +666,115 @@ #include "colamd.h" #include <limits.h> +#include <math.h> #ifdef MATLAB_MEX_FILE #include "mex.h" #include "matrix.h" -#else +#endif /* MATLAB_MEX_FILE */ + +#if !defined (NPRINT) || !defined (NDEBUG) #include <stdio.h> -#include <assert.h> -#endif /* MATLAB_MEX_FILE */ +#endif + +#ifndef NULL +#define NULL ((void *) 0) +#endif + +/* ========================================================================== */ +/* === int or long ========================================================== */ +/* ========================================================================== */ + +#ifdef DLONG + +#define Int long +#define ID "%ld" +#define Int_MAX LONG_MAX +#define COLAMD_recommended colamd_l_recommended +#define COLAMD_set_defaults colamd_l_set_defaults +#define COLAMD_MAIN colamd_l +#define SYMAMD_MAIN symamd_l +#define COLAMD_report colamd_l_report +#define SYMAMD_report symamd_l_report + +#else + +#define Int int +#define ID "%d" +#define Int_MAX INT_MAX +#define COLAMD_recommended colamd_recommended +#define COLAMD_set_defaults colamd_set_defaults +#define COLAMD_MAIN colamd +#define SYMAMD_MAIN symamd +#define COLAMD_report colamd_report +#define SYMAMD_report symamd_report + +#endif + +/* ========================================================================== */ +/* === Row and Column structures ============================================ */ +/* ========================================================================== */ + +/* User code that makes use of the colamd/symamd routines need not directly */ +/* reference these structures. They are used only for the COLAMD_RECOMMENDED */ +/* macro. */ + +typedef struct Colamd_Col_struct +{ + Int start ; /* index for A of first row in this column, or DEAD */ + /* if column is dead */ + Int length ; /* number of rows in this column */ + union + { + Int thickness ; /* number of original columns represented by this */ + /* col, if the column is alive */ + Int parent ; /* parent in parent tree super-column structure, if */ + /* the column is dead */ + } shared1 ; + union + { + Int score ; /* the score used to maintain heap, if col is alive */ + Int order ; /* pivot ordering of this column, if col is dead */ + } shared2 ; + union + { + Int headhash ; /* head of a hash bucket, if col is at the head of */ + /* a degree list */ + Int hash ; /* hash value, if col is not in a degree list */ + Int prev ; /* previous column in degree list, if col is in a */ + /* degree list (but not at the head of a degree list) */ + } shared3 ; + union + { + Int degree_next ; /* next column, if col is in a degree list */ + Int hash_next ; /* next column, if col is in a hash list */ + } shared4 ; + +} Colamd_Col ; + +typedef struct Colamd_Row_struct +{ + Int start ; /* index for A of first col in this row */ + Int length ; /* number of principal columns in this row */ + union + { + Int degree ; /* number of principal & non-principal columns in row */ + Int p ; /* used as a row pointer in init_rows_cols () */ + } shared1 ; + union + { + Int mark ; /* for computing set differences and marking dead rows*/ + Int first_column ;/* first column in row (used in garbage collection) */ + } shared2 ; + +} Colamd_Row ; + +/* size of the Col and Row structures */ +#define COLAMD_C(n_col) \ + ((Int) (((n_col) + 1) * sizeof (Colamd_Col) / sizeof (Int))) + +#define COLAMD_R(n_row) \ + ((Int) (((n_row) + 1) * sizeof (Colamd_Row) / sizeof (Int))) /* ========================================================================== */ /* === Definitions ========================================================== */ @@ -640,6 +784,9 @@ #define PUBLIC #define PRIVATE static +#define DENSE_DEGREE(alpha,n) \ + ((Int) MAX (16.0, (alpha) * sqrt ((double) (n)))) + #define MAX(a,b) (((a) > (b)) ? (a) : (b)) #define MIN(a,b) (((a) < (b)) ? (a) : (b)) @@ -684,111 +831,105 @@ /* === Colamd reporting mechanism =========================================== */ /* ========================================================================== */ -#ifdef MATLAB_MEX_FILE - -/* use mexPrintf in a MATLAB mexFunction, for debugging and statistics output */ -#define PRINTF mexPrintf - +#if defined (MATLAB_MEX_FILE) || defined (MATHWORKS) /* In MATLAB, matrices are 1-based to the user, but 0-based internally */ #define INDEX(i) ((i)+1) - #else - -/* Use printf in standard C environment, for debugging and statistics output. */ -/* Output is generated only if debugging is enabled at compile time, or if */ -/* the caller explicitly calls colamd_report or symamd_report. */ -#define PRINTF printf - /* In C, matrices are 0-based and indices are reported as such in *_report */ #define INDEX(i) (i) - -#endif /* MATLAB_MEX_FILE */ +#endif + +/* All output goes through the PRINTF macro. */ +#define PRINTF(params) { if (colamd_printf != NULL) (void) colamd_printf params ; } /* ========================================================================== */ /* === Prototypes of PRIVATE routines ======================================= */ /* ========================================================================== */ -PRIVATE int init_rows_cols +PRIVATE Int init_rows_cols ( - int n_row, - int n_col, + Int n_row, + Int n_col, Colamd_Row Row [], Colamd_Col Col [], - int A [], - int p [], - int stats [COLAMD_STATS] + Int A [], + Int p [], + Int stats [COLAMD_STATS] ) ; PRIVATE void init_scoring ( - int n_row, - int n_col, + Int n_row, + Int n_col, Colamd_Row Row [], Colamd_Col Col [], - int A [], - int head [], + Int A [], + Int head [], double knobs [COLAMD_KNOBS], - int *p_n_row2, - int *p_n_col2, - int *p_max_deg + Int *p_n_row2, + Int *p_n_col2, + Int *p_max_deg ) ; -PRIVATE int find_ordering +PRIVATE Int find_ordering ( - int n_row, - int n_col, - int Alen, + Int n_row, + Int n_col, + Int Alen, Colamd_Row Row [], Colamd_Col Col [], - int A [], - int head [], - int n_col2, - int max_deg, - int pfree + Int A [], + Int head [], + Int n_col2, + Int max_deg, + Int pfree, + Int aggressive ) ; PRIVATE void order_children ( - int n_col, + Int n_col, Colamd_Col Col [], - int p [] + Int p [] ) ; PRIVATE void detect_super_cols ( #ifndef NDEBUG - int n_col, + Int n_col, Colamd_Row Row [], #endif /* NDEBUG */ Colamd_Col Col [], - int A [], - int head [], - int row_start, - int row_length + Int A [], + Int head [], + Int row_start, + Int row_length ) ; -PRIVATE int garbage_collection +PRIVATE Int garbage_collection ( - int n_row, - int n_col, + Int n_row, + Int n_col, Colamd_Row Row [], Colamd_Col Col [], - int A [], - int *pfree + Int A [], + Int *pfree ) ; -PRIVATE int clear_mark +PRIVATE Int clear_mark ( - int n_row, + Int tag_mark, + Int max_mark, + Int n_row, Colamd_Row Row [] ) ; PRIVATE void print_report ( char *method, - int stats [COLAMD_STATS] + Int stats [COLAMD_STATS] ) ; /* ========================================================================== */ @@ -797,16 +938,18 @@ #ifndef NDEBUG +#include <assert.h> + /* colamd_debug is the *ONLY* global variable, and is only */ /* present when debugging */ -PRIVATE int colamd_debug ; /* debug print level */ - -#define DEBUG0(params) { (void) PRINTF params ; } -#define DEBUG1(params) { if (colamd_debug >= 1) (void) PRINTF params ; } -#define DEBUG2(params) { if (colamd_debug >= 2) (void) PRINTF params ; } -#define DEBUG3(params) { if (colamd_debug >= 3) (void) PRINTF params ; } -#define DEBUG4(params) { if (colamd_debug >= 4) (void) PRINTF params ; } +PRIVATE Int colamd_debug = 0 ; /* debug print level */ + +#define DEBUG0(params) { PRINTF (params) ; } +#define DEBUG1(params) { if (colamd_debug >= 1) PRINTF (params) ; } +#define DEBUG2(params) { if (colamd_debug >= 2) PRINTF (params) ; } +#define DEBUG3(params) { if (colamd_debug >= 3) PRINTF (params) ; } +#define DEBUG4(params) { if (colamd_debug >= 4) PRINTF (params) ; } #ifdef MATLAB_MEX_FILE #define ASSERT(expression) (mxAssert ((expression), "")) @@ -821,41 +964,41 @@ PRIVATE void debug_deg_lists ( - int n_row, - int n_col, + Int n_row, + Int n_col, Colamd_Row Row [], Colamd_Col Col [], - int head [], - int min_score, - int should, - int max_deg + Int head [], + Int min_score, + Int should, + Int max_deg ) ; PRIVATE void debug_mark ( - int n_row, + Int n_row, Colamd_Row Row [], - int tag_mark, - int max_mark + Int tag_mark, + Int max_mark ) ; PRIVATE void debug_matrix ( - int n_row, - int n_col, + Int n_row, + Int n_col, Colamd_Row Row [], Colamd_Col Col [], - int A [] + Int A [] ) ; PRIVATE void debug_structures ( - int n_row, - int n_col, + Int n_row, + Int n_col, Colamd_Row Row [], Colamd_Col Col [], - int A [], - int n_col2 + Int A [], + Int n_col2 ) ; #else /* NDEBUG */ @@ -868,19 +1011,14 @@ #define DEBUG3(params) ; #define DEBUG4(params) ; -#define ASSERT(expression) ((void) 0) +#define ASSERT(expression) #endif /* NDEBUG */ /* ========================================================================== */ - - - -/* ========================================================================== */ /* === USER-CALLABLE ROUTINES: ============================================== */ /* ========================================================================== */ - /* ========================================================================== */ /* === colamd_recommended =================================================== */ /* ========================================================================== */ @@ -892,17 +1030,54 @@ argument is negative, a -1 is returned as an error condition. This function is also available as a macro defined in colamd.h, so that you can use it for a statically-allocated array size. + + The recommended length Alen of the array A passed to colamd is given by + the COLAMD_RECOMMENDED (nnz, n_row, n_col) macro. It returns -1 if any + argument is negative. 2*nnz space is required for the row and column + indices of the matrix. COLAMD_C (n_col) + COLAMD_R (n_row) space is + required for the Col and Row arrays, respectively, which are internal to + colamd. An additional n_col space is the minimal amount of "elbow room", + and nnz/5 more space is recommended for run time efficiency. + + This macro is not needed when using symamd. + + Explicit typecast to Int added Sept. 23, 2002, COLAMD version 2.2, to avoid + gcc -pedantic warning messages. + + Change in version 2.4: + The COLAMD_RECOMMEND macro does not check for integer overflow, but the + routine colamd_recommended does. The macro is thus no longer available + to the user. */ -PUBLIC int colamd_recommended /* returns recommended value of Alen. */ +#define COLAMD_RECOMMENDED(nnz, n_row, n_col) \ +( \ +((nnz) < 0 || (n_row) < 0 || (n_col) < 0) \ +? \ + (-1) \ +: \ + (2 * (nnz) + COLAMD_C (n_col) + COLAMD_R (n_row) + (n_col) + ((nnz) / 5)) \ +) + +PUBLIC Int COLAMD_recommended /* returns recommended value of Alen. */ ( /* === Parameters ======================================================= */ - int nnz, /* number of nonzeros in A */ - int n_row, /* number of rows in A */ - int n_col /* number of columns in A */ + Int nnz, /* number of nonzeros in A */ + Int n_row, /* number of rows in A */ + Int n_col /* number of columns in A */ ) { + double xnz = nnz ; + double xncol = n_col ; + double xnrow = n_row ; + + /* change for version 2.4: check Alen for integer overflow */ + if (COLAMD_RECOMMENDED (xnz, xnrow, xncol) >= INT_MAX) + { + return (-1) ; + } + return (COLAMD_RECOMMENDED (nnz, n_row, n_col)) ; } @@ -928,7 +1103,7 @@ knobs [2..19] unused, but future versions might use this */ -PUBLIC void colamd_set_defaults +PUBLIC void COLAMD_set_defaults ( /* === Parameters ======================================================= */ @@ -937,7 +1112,7 @@ { /* === Local variables ================================================== */ - int i ; + Int i ; if (!knobs) { @@ -947,8 +1122,9 @@ { knobs [i] = 0 ; } - knobs [COLAMD_DENSE_ROW] = 0.5 ; /* ignore rows over 50% dense */ - knobs [COLAMD_DENSE_COL] = 0.5 ; /* ignore columns over 50% dense */ + knobs [COLAMD_DENSE_ROW] = 10 ; + knobs [COLAMD_DENSE_COL] = 10 ; + knobs [COLAMD_AGGRESSIVE] = TRUE ; /* default: do aggressive absorption*/ } @@ -956,16 +1132,16 @@ /* === symamd =============================================================== */ /* ========================================================================== */ -PUBLIC int symamd /* return TRUE if OK, FALSE otherwise */ +PUBLIC Int SYMAMD_MAIN /* return TRUE if OK, FALSE otherwise */ ( /* === Parameters ======================================================= */ - int n, /* number of rows and columns of A */ - int A [], /* row indices of A */ - int p [], /* column pointers of A */ - int perm [], /* output permutation, size n+1 */ + Int n, /* number of rows and columns of A */ + Int A [], /* row indices of A */ + Int p [], /* column pointers of A */ + Int perm [], /* output permutation, size n+1 */ double knobs [COLAMD_KNOBS], /* parameters (uses defaults if NULL) */ - int stats [COLAMD_STATS], /* output statistics and error codes */ + Int stats [COLAMD_STATS], /* output statistics and error codes */ void * (*allocate) (size_t, size_t), /* pointer to calloc (ANSI C) or */ /* mxCalloc (for MATLAB mexFunction) */ @@ -976,23 +1152,22 @@ { /* === Local variables ================================================== */ - int *count ; /* length of each column of M, and col pointer*/ - int *mark ; /* mark array for finding duplicate entries */ - int *M ; /* row indices of matrix M */ - int Mlen ; /* length of M */ - int n_row ; /* number of rows in M */ - int nnz ; /* number of entries in A */ - int i ; /* row index of A */ - int j ; /* column index of A */ - int k ; /* row index of M */ - int mnz ; /* number of nonzeros in M */ - int pp ; /* index into a column of A */ - int last_row ; /* last row seen in the current column */ - int length ; /* number of nonzeros in a column */ + Int *count ; /* length of each column of M, and col pointer*/ + Int *mark ; /* mark array for finding duplicate entries */ + Int *M ; /* row indices of matrix M */ + Int Mlen ; /* length of M */ + Int n_row ; /* number of rows in M */ + Int nnz ; /* number of entries in A */ + Int i ; /* row index of A */ + Int j ; /* column index of A */ + Int k ; /* row index of M */ + Int mnz ; /* number of nonzeros in M */ + Int pp ; /* index into a column of A */ + Int last_row ; /* last row seen in the current column */ + Int length ; /* number of nonzeros in a column */ double cknobs [COLAMD_KNOBS] ; /* knobs for colamd */ double default_knobs [COLAMD_KNOBS] ; /* default knobs for colamd */ - int cstats [COLAMD_STATS] ; /* colamd stats */ #ifndef NDEBUG colamd_get_debug ("symamd") ; @@ -1056,13 +1231,13 @@ if (!knobs) { - colamd_set_defaults (default_knobs) ; + COLAMD_set_defaults (default_knobs) ; knobs = default_knobs ; } /* === Allocate count and mark ========================================== */ - count = (int *) ((*allocate) (n+1, sizeof (int))) ; + count = (Int *) ((*allocate) (n+1, sizeof (Int))) ; if (!count) { stats [COLAMD_STATUS] = COLAMD_ERROR_out_of_memory ; @@ -1070,7 +1245,7 @@ return (FALSE) ; } - mark = (int *) ((*allocate) (n+1, sizeof (int))) ; + mark = (Int *) ((*allocate) (n+1, sizeof (Int))) ; if (!mark) { stats [COLAMD_STATUS] = COLAMD_ERROR_out_of_memory ; @@ -1146,11 +1321,7 @@ } } - if (stats [COLAMD_STATUS] == COLAMD_OK) - { - /* if there are no duplicate entries, then mark is no longer needed */ - (*release) ((void *) mark) ; - } + /* v2.4: removed free(mark) */ /* === Compute column pointers of M ===================================== */ @@ -1169,8 +1340,8 @@ mnz = perm [n] ; n_row = mnz / 2 ; - Mlen = colamd_recommended (mnz, n_row, n) ; - M = (int *) ((*allocate) (Mlen, sizeof (int))) ; + Mlen = COLAMD_recommended (mnz, n_row, n) ; + M = (Int *) ((*allocate) (Mlen, sizeof (Int))) ; DEBUG0 (("symamd: M is %d-by-%d with %d entries, Mlen = %d\n", n_row, n, mnz, Mlen)) ; @@ -1230,11 +1401,12 @@ } } } - (*release) ((void *) mark) ; + /* v2.4: free(mark) moved below */ } /* count and mark no longer needed */ (*release) ((void *) count) ; + (*release) ((void *) mark) ; /* v2.4: free (mark) moved here */ ASSERT (k == n_row) ; /* === Adjust the knobs for M =========================================== */ @@ -1245,41 +1417,20 @@ } /* there are no dense rows in M */ - cknobs [COLAMD_DENSE_ROW] = 1.0 ; - - if (n_row != 0 && n < n_row) - { - /* On input, the knob is a fraction of 1..n, the number of rows of A. */ - /* Convert it to a fraction of 1..n_row, of the number of rows of M. */ - cknobs [COLAMD_DENSE_COL] = (knobs [COLAMD_DENSE_ROW] * n) / n_row ; - } - else - { - /* no dense columns in M */ - cknobs [COLAMD_DENSE_COL] = 1.0 ; - } - - DEBUG0 (("symamd: dense col knob for M: %g\n", cknobs [COLAMD_DENSE_COL])) ; + cknobs [COLAMD_DENSE_ROW] = -1 ; + cknobs [COLAMD_DENSE_COL] = knobs [COLAMD_DENSE_ROW] ; /* === Order the columns of M =========================================== */ - if (!colamd (n_row, n, Mlen, M, perm, cknobs, cstats)) - { - /* This "cannot" happen, unless there is a bug in the code. */ - stats [COLAMD_STATUS] = COLAMD_ERROR_internal_error ; - (*release) ((void *) M) ; - DEBUG0 (("symamd: internal error!\n")) ; - return (FALSE) ; - } + /* v2.4: colamd cannot fail here, so the error check is removed */ + (void) COLAMD_MAIN (n_row, n, Mlen, M, perm, cknobs, stats) ; /* Note that the output permutation is now in perm */ /* === get the statistics for symamd from colamd ======================== */ - /* note that a dense column in colamd means a dense row and col in symamd */ - stats [COLAMD_DENSE_ROW] = cstats [COLAMD_DENSE_COL] ; - stats [COLAMD_DENSE_COL] = cstats [COLAMD_DENSE_COL] ; - stats [COLAMD_DEFRAG_COUNT] = cstats [COLAMD_DEFRAG_COUNT] ; + /* a dense column in colamd means a dense row and col in symamd */ + stats [COLAMD_DENSE_ROW] = stats [COLAMD_DENSE_COL] ; /* === Free M =========================================================== */ @@ -1301,33 +1452,34 @@ (AQ)'(AQ) = LL' remains sparse. */ -PUBLIC int colamd /* returns TRUE if successful, FALSE otherwise*/ +PUBLIC Int COLAMD_MAIN /* returns TRUE if successful, FALSE otherwise*/ ( /* === Parameters ======================================================= */ - int n_row, /* number of rows in A */ - int n_col, /* number of columns in A */ - int Alen, /* length of A */ - int A [], /* row indices of A */ - int p [], /* pointers to columns in A */ + Int n_row, /* number of rows in A */ + Int n_col, /* number of columns in A */ + Int Alen, /* length of A */ + Int A [], /* row indices of A */ + Int p [], /* pointers to columns in A */ double knobs [COLAMD_KNOBS],/* parameters (uses defaults if NULL) */ - int stats [COLAMD_STATS] /* output statistics and error codes */ + Int stats [COLAMD_STATS] /* output statistics and error codes */ ) { /* === Local variables ================================================== */ - int i ; /* loop index */ - int nnz ; /* nonzeros in A */ - int Row_size ; /* size of Row [], in integers */ - int Col_size ; /* size of Col [], in integers */ - int need ; /* minimum required length of A */ + Int i ; /* loop index */ + Int nnz ; /* nonzeros in A */ + Int Row_size ; /* size of Row [], in integers */ + Int Col_size ; /* size of Col [], in integers */ + Int need ; /* minimum required length of A */ Colamd_Row *Row ; /* pointer into A of Row [0..n_row] array */ Colamd_Col *Col ; /* pointer into A of Col [0..n_col] array */ - int n_col2 ; /* number of non-dense, non-empty columns */ - int n_row2 ; /* number of non-dense, non-empty rows */ - int ngarbage ; /* number of garbage collections performed */ - int max_deg ; /* maximum row degree */ + Int n_col2 ; /* number of non-dense, non-empty columns */ + Int n_row2 ; /* number of non-dense, non-empty rows */ + Int ngarbage ; /* number of garbage collections performed */ + Int max_deg ; /* maximum row degree */ double default_knobs [COLAMD_KNOBS] ; /* default knobs array */ + Int aggressive ; /* do aggressive absorption */ #ifndef NDEBUG colamd_get_debug ("colamd") ; @@ -1399,10 +1551,12 @@ if (!knobs) { - colamd_set_defaults (default_knobs) ; + COLAMD_set_defaults (default_knobs) ; knobs = default_knobs ; } + aggressive = (knobs [COLAMD_AGGRESSIVE] != FALSE) ; + /* === Allocate the Row and Col arrays from array A ===================== */ Col_size = COLAMD_C (n_col) ; @@ -1440,7 +1594,7 @@ /* === Order the supercolumns =========================================== */ ngarbage = find_ordering (n_row, n_col, Alen, Row, Col, A, p, - n_col2, max_deg, 2*nnz) ; + n_col2, max_deg, 2*nnz, aggressive) ; /* === Order the non-principal columns ================================== */ @@ -1460,9 +1614,9 @@ /* === colamd_report ======================================================== */ /* ========================================================================== */ -PUBLIC void colamd_report +PUBLIC void COLAMD_report ( - int stats [COLAMD_STATS] + Int stats [COLAMD_STATS] ) { print_report ("colamd", stats) ; @@ -1473,9 +1627,9 @@ /* === symamd_report ======================================================== */ /* ========================================================================== */ -PUBLIC void symamd_report +PUBLIC void SYMAMD_report ( - int stats [COLAMD_STATS] + Int stats [COLAMD_STATS] ) { print_report ("symamd", stats) ; @@ -1503,28 +1657,28 @@ TRUE otherwise. Not user-callable. */ -PRIVATE int init_rows_cols /* returns TRUE if OK, or FALSE otherwise */ +PRIVATE Int init_rows_cols /* returns TRUE if OK, or FALSE otherwise */ ( /* === Parameters ======================================================= */ - int n_row, /* number of rows of A */ - int n_col, /* number of columns of A */ + Int n_row, /* number of rows of A */ + Int n_col, /* number of columns of A */ Colamd_Row Row [], /* of size n_row+1 */ Colamd_Col Col [], /* of size n_col+1 */ - int A [], /* row indices of A, of size Alen */ - int p [], /* pointers to columns in A, of size n_col+1 */ - int stats [COLAMD_STATS] /* colamd statistics */ + Int A [], /* row indices of A, of size Alen */ + Int p [], /* pointers to columns in A, of size n_col+1 */ + Int stats [COLAMD_STATS] /* colamd statistics */ ) { /* === Local variables ================================================== */ - int col ; /* a column index */ - int row ; /* a row index */ - int *cp ; /* a column pointer */ - int *cp_end ; /* a pointer to the end of a column */ - int *rp ; /* a row pointer */ - int *rp_end ; /* a pointer to the end of a row */ - int last_row ; /* previous row */ + Int col ; /* a column index */ + Int row ; /* a row index */ + Int *cp ; /* a column pointer */ + Int *cp_end ; /* a pointer to the end of a column */ + Int *rp ; /* a row pointer */ + Int *rp_end ; /* a pointer to the end of a row */ + Int last_row ; /* previous row */ /* === Initialize columns, and check column pointers ==================== */ @@ -1744,44 +1898,63 @@ ( /* === Parameters ======================================================= */ - int n_row, /* number of rows of A */ - int n_col, /* number of columns of A */ + Int n_row, /* number of rows of A */ + Int n_col, /* number of columns of A */ Colamd_Row Row [], /* of size n_row+1 */ Colamd_Col Col [], /* of size n_col+1 */ - int A [], /* column form and row form of A */ - int head [], /* of size n_col+1 */ + Int A [], /* column form and row form of A */ + Int head [], /* of size n_col+1 */ double knobs [COLAMD_KNOBS],/* parameters */ - int *p_n_row2, /* number of non-dense, non-empty rows */ - int *p_n_col2, /* number of non-dense, non-empty columns */ - int *p_max_deg /* maximum row degree */ + Int *p_n_row2, /* number of non-dense, non-empty rows */ + Int *p_n_col2, /* number of non-dense, non-empty columns */ + Int *p_max_deg /* maximum row degree */ ) { /* === Local variables ================================================== */ - int c ; /* a column index */ - int r, row ; /* a row index */ - int *cp ; /* a column pointer */ - int deg ; /* degree of a row or column */ - int *cp_end ; /* a pointer to the end of a column */ - int *new_cp ; /* new column pointer */ - int col_length ; /* length of pruned column */ - int score ; /* current column score */ - int n_col2 ; /* number of non-dense, non-empty columns */ - int n_row2 ; /* number of non-dense, non-empty rows */ - int dense_row_count ; /* remove rows with more entries than this */ - int dense_col_count ; /* remove cols with more entries than this */ - int min_score ; /* smallest column score */ - int max_deg ; /* maximum row degree */ - int next_col ; /* Used to add to degree list.*/ + Int c ; /* a column index */ + Int r, row ; /* a row index */ + Int *cp ; /* a column pointer */ + Int deg ; /* degree of a row or column */ + Int *cp_end ; /* a pointer to the end of a column */ + Int *new_cp ; /* new column pointer */ + Int col_length ; /* length of pruned column */ + Int score ; /* current column score */ + Int n_col2 ; /* number of non-dense, non-empty columns */ + Int n_row2 ; /* number of non-dense, non-empty rows */ + Int dense_row_count ; /* remove rows with more entries than this */ + Int dense_col_count ; /* remove cols with more entries than this */ + Int min_score ; /* smallest column score */ + Int max_deg ; /* maximum row degree */ + Int next_col ; /* Used to add to degree list.*/ #ifndef NDEBUG - int debug_count ; /* debug only. */ + Int debug_count ; /* debug only. */ #endif /* NDEBUG */ /* === Extract knobs ==================================================== */ - dense_row_count = MAX (0, MIN (knobs [COLAMD_DENSE_ROW] * n_col, n_col)) ; - dense_col_count = MAX (0, MIN (knobs [COLAMD_DENSE_COL] * n_row, n_row)) ; + /* Note: if knobs contains a NaN, this is undefined: */ + if (knobs [COLAMD_DENSE_ROW] < 0) + { + /* only remove completely dense rows */ + dense_row_count = n_col-1 ; + } + else + { + dense_row_count = DENSE_DEGREE (knobs [COLAMD_DENSE_ROW], n_col) ; + } + if (knobs [COLAMD_DENSE_COL] < 0) + { + /* only remove completely dense columns */ + dense_col_count = n_row-1 ; + } + else + { + dense_col_count = + DENSE_DEGREE (knobs [COLAMD_DENSE_COL], MIN (n_row, n_col)) ; + } + DEBUG1 (("colamd: densecount: %d %d\n", dense_row_count, dense_col_count)) ; max_deg = 0 ; n_col2 = n_col ; @@ -1886,7 +2059,7 @@ score = MIN (score, n_col) ; } /* determine pruned column length */ - col_length = (int) (new_cp - &A [Col [c].start]) ; + col_length = (Int) (new_cp - &A [Col [c].start]) ; if (col_length == 0) { /* a newly-made null column (all rows in this col are "dense" */ @@ -1996,65 +2169,66 @@ degree ordering method. Not user-callable. */ -PRIVATE int find_ordering /* return the number of garbage collections */ +PRIVATE Int find_ordering /* return the number of garbage collections */ ( /* === Parameters ======================================================= */ - int n_row, /* number of rows of A */ - int n_col, /* number of columns of A */ - int Alen, /* size of A, 2*nnz + n_col or larger */ + Int n_row, /* number of rows of A */ + Int n_col, /* number of columns of A */ + Int Alen, /* size of A, 2*nnz + n_col or larger */ Colamd_Row Row [], /* of size n_row+1 */ Colamd_Col Col [], /* of size n_col+1 */ - int A [], /* column form and row form of A */ - int head [], /* of size n_col+1 */ - int n_col2, /* Remaining columns to order */ - int max_deg, /* Maximum row degree */ - int pfree /* index of first free slot (2*nnz on entry) */ + Int A [], /* column form and row form of A */ + Int head [], /* of size n_col+1 */ + Int n_col2, /* Remaining columns to order */ + Int max_deg, /* Maximum row degree */ + Int pfree, /* index of first free slot (2*nnz on entry) */ + Int aggressive ) { /* === Local variables ================================================== */ - int k ; /* current pivot ordering step */ - int pivot_col ; /* current pivot column */ - int *cp ; /* a column pointer */ - int *rp ; /* a row pointer */ - int pivot_row ; /* current pivot row */ - int *new_cp ; /* modified column pointer */ - int *new_rp ; /* modified row pointer */ - int pivot_row_start ; /* pointer to start of pivot row */ - int pivot_row_degree ; /* number of columns in pivot row */ - int pivot_row_length ; /* number of supercolumns in pivot row */ - int pivot_col_score ; /* score of pivot column */ - int needed_memory ; /* free space needed for pivot row */ - int *cp_end ; /* pointer to the end of a column */ - int *rp_end ; /* pointer to the end of a row */ - int row ; /* a row index */ - int col ; /* a column index */ - int max_score ; /* maximum possible score */ - int cur_score ; /* score of current column */ - unsigned int hash ; /* hash value for supernode detection */ - int head_column ; /* head of hash bucket */ - int first_col ; /* first column in hash bucket */ - int tag_mark ; /* marker value for mark array */ - int row_mark ; /* Row [row].shared2.mark */ - int set_difference ; /* set difference size of row with pivot row */ - int min_score ; /* smallest column score */ - int col_thickness ; /* "thickness" (no. of columns in a supercol) */ - int max_mark ; /* maximum value of tag_mark */ - int pivot_col_thickness ; /* number of columns represented by pivot col */ - int prev_col ; /* Used by Dlist operations. */ - int next_col ; /* Used by Dlist operations. */ - int ngarbage ; /* number of garbage collections performed */ + Int k ; /* current pivot ordering step */ + Int pivot_col ; /* current pivot column */ + Int *cp ; /* a column pointer */ + Int *rp ; /* a row pointer */ + Int pivot_row ; /* current pivot row */ + Int *new_cp ; /* modified column pointer */ + Int *new_rp ; /* modified row pointer */ + Int pivot_row_start ; /* pointer to start of pivot row */ + Int pivot_row_degree ; /* number of columns in pivot row */ + Int pivot_row_length ; /* number of supercolumns in pivot row */ + Int pivot_col_score ; /* score of pivot column */ + Int needed_memory ; /* free space needed for pivot row */ + Int *cp_end ; /* pointer to the end of a column */ + Int *rp_end ; /* pointer to the end of a row */ + Int row ; /* a row index */ + Int col ; /* a column index */ + Int max_score ; /* maximum possible score */ + Int cur_score ; /* score of current column */ + unsigned Int hash ; /* hash value for supernode detection */ + Int head_column ; /* head of hash bucket */ + Int first_col ; /* first column in hash bucket */ + Int tag_mark ; /* marker value for mark array */ + Int row_mark ; /* Row [row].shared2.mark */ + Int set_difference ; /* set difference size of row with pivot row */ + Int min_score ; /* smallest column score */ + Int col_thickness ; /* "thickness" (no. of columns in a supercol) */ + Int max_mark ; /* maximum value of tag_mark */ + Int pivot_col_thickness ; /* number of columns represented by pivot col */ + Int prev_col ; /* Used by Dlist operations. */ + Int next_col ; /* Used by Dlist operations. */ + Int ngarbage ; /* number of garbage collections performed */ #ifndef NDEBUG - int debug_d ; /* debug loop counter */ - int debug_step = 0 ; /* debug loop counter */ + Int debug_d ; /* debug loop counter */ + Int debug_step = 0 ; /* debug loop counter */ #endif /* NDEBUG */ /* === Initialization and clear mark ==================================== */ max_mark = INT_MAX - n_col ; /* INT_MAX defined in <limits.h> */ - tag_mark = clear_mark (n_row, Row) ; + tag_mark = clear_mark (0, max_mark, n_row, Row) ; min_score = 0 ; ngarbage = 0 ; DEBUG1 (("colamd: Ordering, n_col2=%d\n", n_col2)) ; @@ -2108,7 +2282,6 @@ } ASSERT (COL_IS_ALIVE (pivot_col)) ; - DEBUG3 (("Pivot col: %d\n", pivot_col)) ; /* remember score for defrag check */ pivot_col_score = Col [pivot_col].shared2.score ; @@ -2120,6 +2293,7 @@ pivot_col_thickness = Col [pivot_col].shared1.thickness ; k += pivot_col_thickness ; ASSERT (pivot_col_thickness > 0) ; + DEBUG3 (("Pivot col: %d thick %d\n", pivot_col, pivot_col_thickness)) ; /* === Garbage_collection, if necessary ============================= */ @@ -2131,7 +2305,7 @@ /* after garbage collection we will have enough */ ASSERT (pfree + needed_memory < Alen) ; /* garbage collection has wiped out the Row[].shared2.mark array */ - tag_mark = clear_mark (n_row, Row) ; + tag_mark = clear_mark (0, max_mark, n_row, Row) ; #ifndef NDEBUG debug_matrix (n_row, n_col, Row, Col, A) ; @@ -2159,26 +2333,25 @@ row = *cp++ ; DEBUG4 (("Pivot col pattern %d %d\n", ROW_IS_ALIVE (row), row)) ; /* skip if row is dead */ - if (ROW_IS_DEAD (row)) - { - continue ; - } - rp = &A [Row [row].start] ; - rp_end = rp + Row [row].length ; - while (rp < rp_end) + if (ROW_IS_ALIVE (row)) { - /* get a column */ - col = *rp++ ; - /* add the column, if alive and untagged */ - col_thickness = Col [col].shared1.thickness ; - if (col_thickness > 0 && COL_IS_ALIVE (col)) + rp = &A [Row [row].start] ; + rp_end = rp + Row [row].length ; + while (rp < rp_end) { - /* tag column in pivot row */ - Col [col].shared1.thickness = -col_thickness ; - ASSERT (pfree < Alen) ; - /* place column in pivot row */ - A [pfree++] = col ; - pivot_row_degree += col_thickness ; + /* get a column */ + col = *rp++ ; + /* add the column, if alive and untagged */ + col_thickness = Col [col].shared1.thickness ; + if (col_thickness > 0 && COL_IS_ALIVE (col)) + { + /* tag column in pivot row */ + Col [col].shared1.thickness = -col_thickness ; + ASSERT (pfree < Alen) ; + /* place column in pivot row */ + A [pfree++] = col ; + pivot_row_degree += col_thickness ; + } } } } @@ -2309,7 +2482,7 @@ set_difference -= col_thickness ; ASSERT (set_difference >= 0) ; /* absorb this row if the set difference becomes zero */ - if (set_difference == 0) + if (set_difference == 0 && aggressive) { DEBUG3 (("aggressive absorption. Row: %d\n", row)) ; KILL_ROW (row) ; @@ -2357,9 +2530,11 @@ /* skip if dead */ if (ROW_IS_MARKED_DEAD (row_mark)) { + DEBUG4 ((" Row %d, dead\n", row)) ; continue ; } - ASSERT (row_mark > tag_mark) ; + DEBUG4 ((" Row %d, set diff %d\n", row, row_mark-tag_mark)); + ASSERT (row_mark >= tag_mark) ; /* compact the column */ *new_cp++ = row ; /* compute hash function */ @@ -2371,7 +2546,7 @@ } /* recompute the column's length */ - Col [col].length = (int) (new_cp - &A [Col [col].start]) ; + Col [col].length = (Int) (new_cp - &A [Col [col].start]) ; /* === Further mass elimination ================================= */ @@ -2400,7 +2575,7 @@ hash %= n_col + 1 ; DEBUG4 ((" Hash = %d, n_col = %d.\n", hash, n_col)) ; - ASSERT (hash <= n_col) ; + ASSERT (((Int) hash) <= n_col) ; head_column = head [hash] ; if (head_column > EMPTY) @@ -2419,7 +2594,7 @@ Col [col].shared4.hash_next = first_col ; /* save hash function in Col [col].shared3.hash */ - Col [col].shared3.hash = (int) hash ; + Col [col].shared3.hash = (Int) hash ; ASSERT (COL_IS_ALIVE (col)) ; } } @@ -2444,12 +2619,7 @@ /* === Clear mark =================================================== */ - tag_mark += (max_deg + 1) ; - if (tag_mark >= max_mark) - { - DEBUG2 (("clearing tag_mark\n")) ; - tag_mark = clear_mark (n_row, Row) ; - } + tag_mark = clear_mark (tag_mark+max_deg+1, max_mark, n_row, Row) ; #ifndef NDEBUG DEBUG3 (("check3\n")) ; @@ -2530,10 +2700,14 @@ /* update pivot row length to reflect any cols that were killed */ /* during super-col detection and mass elimination */ Row [pivot_row].start = pivot_row_start ; - Row [pivot_row].length = (int) (new_rp - &A[pivot_row_start]) ; + Row [pivot_row].length = (Int) (new_rp - &A[pivot_row_start]) ; + ASSERT (Row [pivot_row].length > 0) ; Row [pivot_row].shared1.degree = pivot_row_degree ; Row [pivot_row].shared2.mark = 0 ; /* pivot row is no longer dead */ + + DEBUG1 (("Resurrect Pivot_row %d deg: %d\n", + pivot_row, pivot_row_degree)) ; } } @@ -2564,17 +2738,17 @@ ( /* === Parameters ======================================================= */ - int n_col, /* number of columns of A */ + Int n_col, /* number of columns of A */ Colamd_Col Col [], /* of size n_col+1 */ - int p [] /* p [0 ... n_col-1] is the column permutation*/ + Int p [] /* p [0 ... n_col-1] is the column permutation*/ ) { /* === Local variables ================================================== */ - int i ; /* loop counter for all columns */ - int c ; /* column index */ - int parent ; /* index of column's parent */ - int order ; /* column's order */ + Int i ; /* loop counter for all columns */ + Int c ; /* column index */ + Int parent ; /* index of column's parent */ + Int order ; /* column's order */ /* === Order each non-principal column ================================== */ @@ -2667,32 +2841,32 @@ #ifndef NDEBUG /* these two parameters are only needed when debugging is enabled: */ - int n_col, /* number of columns of A */ + Int n_col, /* number of columns of A */ Colamd_Row Row [], /* of size n_row+1 */ #endif /* NDEBUG */ Colamd_Col Col [], /* of size n_col+1 */ - int A [], /* row indices of A */ - int head [], /* head of degree lists and hash buckets */ - int row_start, /* pointer to set of columns to check */ - int row_length /* number of columns to check */ + Int A [], /* row indices of A */ + Int head [], /* head of degree lists and hash buckets */ + Int row_start, /* pointer to set of columns to check */ + Int row_length /* number of columns to check */ ) { /* === Local variables ================================================== */ - int hash ; /* hash value for a column */ - int *rp ; /* pointer to a row */ - int c ; /* a column index */ - int super_c ; /* column index of the column to absorb into */ - int *cp1 ; /* column pointer for column super_c */ - int *cp2 ; /* column pointer for column c */ - int length ; /* length of column super_c */ - int prev_c ; /* column preceding c in hash bucket */ - int i ; /* loop counter */ - int *rp_end ; /* pointer to the end of the row */ - int col ; /* a column index in the row to check */ - int head_column ; /* first column in hash bucket or degree list */ - int first_col ; /* first column in hash bucket */ + Int hash ; /* hash value for a column */ + Int *rp ; /* pointer to a row */ + Int c ; /* a column index */ + Int super_c ; /* column index of the column to absorb into */ + Int *cp1 ; /* column pointer for column super_c */ + Int *cp2 ; /* column pointer for column c */ + Int length ; /* length of column super_c */ + Int prev_c ; /* column preceding c in hash bucket */ + Int i ; /* loop counter */ + Int *rp_end ; /* pointer to the end of the row */ + Int col ; /* a column index in the row to check */ + Int head_column ; /* first column in hash bucket or degree list */ + Int first_col ; /* first column in hash bucket */ /* === Consider each column in the row ================================== */ @@ -2818,29 +2992,29 @@ Not user-callable. */ -PRIVATE int garbage_collection /* returns the new value of pfree */ +PRIVATE Int garbage_collection /* returns the new value of pfree */ ( /* === Parameters ======================================================= */ - int n_row, /* number of rows */ - int n_col, /* number of columns */ + Int n_row, /* number of rows */ + Int n_col, /* number of columns */ Colamd_Row Row [], /* row info */ Colamd_Col Col [], /* column info */ - int A [], /* A [0 ... Alen-1] holds the matrix */ - int *pfree /* &A [0] ... pfree is in use */ + Int A [], /* A [0 ... Alen-1] holds the matrix */ + Int *pfree /* &A [0] ... pfree is in use */ ) { /* === Local variables ================================================== */ - int *psrc ; /* source pointer */ - int *pdest ; /* destination pointer */ - int j ; /* counter */ - int r ; /* a row index */ - int c ; /* a column index */ - int length ; /* length of a row or column */ + Int *psrc ; /* source pointer */ + Int *pdest ; /* destination pointer */ + Int j ; /* counter */ + Int r ; /* a row index */ + Int c ; /* a column index */ + Int length ; /* length of a row or column */ #ifndef NDEBUG - int debug_rows ; + Int debug_rows ; DEBUG2 (("Defrag..\n")) ; for (psrc = &A[0] ; psrc < pfree ; psrc++) ASSERT (*psrc >= 0) ; debug_rows = 0 ; @@ -2857,7 +3031,7 @@ /* move and compact the column */ ASSERT (pdest <= psrc) ; - Col [c].start = (int) (pdest - &A [0]) ; + Col [c].start = (Int) (pdest - &A [0]) ; length = Col [c].length ; for (j = 0 ; j < length ; j++) { @@ -2867,7 +3041,7 @@ *pdest++ = r ; } } - Col [c].length = (int) (pdest - &A [Col [c].start]) ; + Col [c].length = (Int) (pdest - &A [Col [c].start]) ; } } @@ -2875,28 +3049,25 @@ for (r = 0 ; r < n_row ; r++) { - if (ROW_IS_ALIVE (r)) + if (ROW_IS_DEAD (r) || (Row [r].length == 0)) + { + /* This row is already dead, or is of zero length. Cannot compact + * a row of zero length, so kill it. NOTE: in the current version, + * there are no zero-length live rows. Kill the row (for the first + * time, or again) just to be safe. */ + KILL_ROW (r) ; + } + else { - if (Row [r].length == 0) - { - /* this row is of zero length. cannot compact it, so kill it */ - DEBUG3 (("Defrag row kill\n")) ; - KILL_ROW (r) ; - } - else - { - /* save first column index in Row [r].shared2.first_column */ - psrc = &A [Row [r].start] ; - Row [r].shared2.first_column = *psrc ; - ASSERT (ROW_IS_ALIVE (r)) ; - /* flag the start of the row with the one's complement of row */ - *psrc = ONES_COMPLEMENT (r) ; - + /* save first column index in Row [r].shared2.first_column */ + psrc = &A [Row [r].start] ; + Row [r].shared2.first_column = *psrc ; + ASSERT (ROW_IS_ALIVE (r)) ; + /* flag the start of the row with the one's complement of row */ + *psrc = ONES_COMPLEMENT (r) ; #ifndef NDEBUG - debug_rows++ ; + debug_rows++ ; #endif /* NDEBUG */ - - } } } @@ -2915,10 +3086,10 @@ /* restore first column index */ *psrc = Row [r].shared2.first_column ; ASSERT (ROW_IS_ALIVE (r)) ; - + ASSERT (Row [r].length > 0) ; /* move and compact the row */ ASSERT (pdest <= psrc) ; - Row [r].start = (int) (pdest - &A [0]) ; + Row [r].start = (Int) (pdest - &A [0]) ; length = Row [r].length ; for (j = 0 ; j < length ; j++) { @@ -2928,12 +3099,11 @@ *pdest++ = c ; } } - Row [r].length = (int) (pdest - &A [Row [r].start]) ; - + Row [r].length = (Int) (pdest - &A [Row [r].start]) ; + ASSERT (Row [r].length > 0) ; #ifndef NDEBUG debug_rows-- ; #endif /* NDEBUG */ - } } /* ensure we found all the rows */ @@ -2941,7 +3111,7 @@ /* === Return the new value of pfree ==================================== */ - return ((int) (pdest - &A [0])) ; + return ((Int) (pdest - &A [0])) ; } @@ -2954,26 +3124,34 @@ Return value is the new tag_mark. Not user-callable. */ -PRIVATE int clear_mark /* return the new value for tag_mark */ +PRIVATE Int clear_mark /* return the new value for tag_mark */ ( /* === Parameters ======================================================= */ - int n_row, /* number of rows in A */ + Int tag_mark, /* new value of tag_mark */ + Int max_mark, /* max allowed value of tag_mark */ + + Int n_row, /* number of rows in A */ Colamd_Row Row [] /* Row [0 ... n_row-1].shared2.mark is set to zero */ ) { /* === Local variables ================================================== */ - int r ; - - for (r = 0 ; r < n_row ; r++) + Int r ; + + if (tag_mark <= 0 || tag_mark >= max_mark) { - if (ROW_IS_ALIVE (r)) + for (r = 0 ; r < n_row ; r++) { - Row [r].shared2.mark = 0 ; + if (ROW_IS_ALIVE (r)) + { + Row [r].shared2.mark = 0 ; + } } + tag_mark = 1 ; } - return (1) ; + + return (tag_mark) ; } @@ -2984,15 +3162,18 @@ PRIVATE void print_report ( char *method, - int stats [COLAMD_STATS] + Int stats [COLAMD_STATS] ) { - int i1, i2, i3 ; + Int i1, i2, i3 ; + + PRINTF (("\n%s version %d.%d, %s: ", method, + COLAMD_MAIN_VERSION, COLAMD_SUB_VERSION, COLAMD_DATE)) ; if (!stats) { - PRINTF ("%s: No statistics available.\n", method) ; + PRINTF (("No statistics available.\n")) ; return ; } @@ -3002,11 +3183,11 @@ if (stats [COLAMD_STATUS] >= 0) { - PRINTF ("%s: OK. ", method) ; + PRINTF (("OK. ")) ; } else { - PRINTF ("%s: ERROR. ", method) ; + PRINTF (("ERROR. ")) ; } switch (stats [COLAMD_STATUS]) @@ -3014,95 +3195,90 @@ case COLAMD_OK_BUT_JUMBLED: - PRINTF ("Matrix has unsorted or duplicate row indices.\n") ; - - PRINTF ("%s: number of duplicate or out-of-order row indices: %d\n", - method, i3) ; - - PRINTF ("%s: last seen duplicate or out-of-order row index: %d\n", - method, INDEX (i2)) ; - - PRINTF ("%s: last seen in column: %d", - method, INDEX (i1)) ; + PRINTF(("Matrix has unsorted or duplicate row indices.\n")) ; + + PRINTF(("%s: number of duplicate or out-of-order row indices: %d\n", + method, i3)) ; + + PRINTF(("%s: last seen duplicate or out-of-order row index: %d\n", + method, INDEX (i2))) ; + + PRINTF(("%s: last seen in column: %d", + method, INDEX (i1))) ; /* no break - fall through to next case instead */ case COLAMD_OK: - PRINTF ("\n") ; - - PRINTF ("%s: number of dense or empty rows ignored: %d\n", - method, stats [COLAMD_DENSE_ROW]) ; - - PRINTF ("%s: number of dense or empty columns ignored: %d\n", - method, stats [COLAMD_DENSE_COL]) ; - - PRINTF ("%s: number of garbage collections performed: %d\n", - method, stats [COLAMD_DEFRAG_COUNT]) ; + PRINTF(("\n")) ; + + PRINTF(("%s: number of dense or empty rows ignored: %d\n", + method, stats [COLAMD_DENSE_ROW])) ; + + PRINTF(("%s: number of dense or empty columns ignored: %d\n", + method, stats [COLAMD_DENSE_COL])) ; + + PRINTF(("%s: number of garbage collections performed: %d\n", + method, stats [COLAMD_DEFRAG_COUNT])) ; break ; case COLAMD_ERROR_A_not_present: - PRINTF ("Array A (row indices of matrix) not present.\n") ; + PRINTF(("Array A (row indices of matrix) not present.\n")) ; break ; case COLAMD_ERROR_p_not_present: - PRINTF ("Array p (column pointers for matrix) not present.\n") ; + PRINTF(("Array p (column pointers for matrix) not present.\n")) ; break ; case COLAMD_ERROR_nrow_negative: - PRINTF ("Invalid number of rows (%d).\n", i1) ; + PRINTF(("Invalid number of rows (%d).\n", i1)) ; break ; case COLAMD_ERROR_ncol_negative: - PRINTF ("Invalid number of columns (%d).\n", i1) ; + PRINTF(("Invalid number of columns (%d).\n", i1)) ; break ; case COLAMD_ERROR_nnz_negative: - PRINTF ("Invalid number of nonzero entries (%d).\n", i1) ; + PRINTF(("Invalid number of nonzero entries (%d).\n", i1)) ; break ; case COLAMD_ERROR_p0_nonzero: - PRINTF ("Invalid column pointer, p [0] = %d, must be zero.\n", i1) ; + PRINTF(("Invalid column pointer, p [0] = %d, must be zero.\n", i1)); break ; case COLAMD_ERROR_A_too_small: - PRINTF ("Array A too small.\n") ; - PRINTF (" Need Alen >= %d, but given only Alen = %d.\n", - i1, i2) ; + PRINTF(("Array A too small.\n")) ; + PRINTF((" Need Alen >= %d, but given only Alen = %d.\n", + i1, i2)) ; break ; case COLAMD_ERROR_col_length_negative: PRINTF - ("Column %d has a negative number of nonzero entries (%d).\n", - INDEX (i1), i2) ; + (("Column %d has a negative number of nonzero entries (%d).\n", + INDEX (i1), i2)) ; break ; case COLAMD_ERROR_row_index_out_of_bounds: PRINTF - ("Row index (row %d) out of bounds (%d to %d) in column %d.\n", - INDEX (i2), INDEX (0), INDEX (i3-1), INDEX (i1)) ; + (("Row index (row %d) out of bounds (%d to %d) in column %d.\n", + INDEX (i2), INDEX (0), INDEX (i3-1), INDEX (i1))) ; break ; case COLAMD_ERROR_out_of_memory: - PRINTF ("Out of memory.\n") ; + PRINTF(("Out of memory.\n")) ; break ; - case COLAMD_ERROR_internal_error: - - /* if this happens, there is a bug in the code */ - PRINTF - ("Internal error! Please contact authors (davis@cise.ufl.edu).\n") ; - break ; + /* v2.4: internal-error case deleted */ } } @@ -3133,26 +3309,26 @@ ( /* === Parameters ======================================================= */ - int n_row, - int n_col, + Int n_row, + Int n_col, Colamd_Row Row [], Colamd_Col Col [], - int A [], - int n_col2 + Int A [], + Int n_col2 ) { /* === Local variables ================================================== */ - int i ; - int c ; - int *cp ; - int *cp_end ; - int len ; - int score ; - int r ; - int *rp ; - int *rp_end ; - int deg ; + Int i ; + Int c ; + Int *cp ; + Int *cp_end ; + Int len ; + Int score ; + Int r ; + Int *rp ; + Int *rp_end ; + Int deg ; /* === Check A, Row, and Col ============================================ */ @@ -3220,22 +3396,22 @@ ( /* === Parameters ======================================================= */ - int n_row, - int n_col, + Int n_row, + Int n_col, Colamd_Row Row [], Colamd_Col Col [], - int head [], - int min_score, - int should, - int max_deg + Int head [], + Int min_score, + Int should, + Int max_deg ) { /* === Local variables ================================================== */ - int deg ; - int col ; - int have ; - int row ; + Int deg ; + Int col ; + Int have ; + Int row ; /* === Check the degree lists =========================================== */ @@ -3294,15 +3470,15 @@ ( /* === Parameters ======================================================= */ - int n_row, + Int n_row, Colamd_Row Row [], - int tag_mark, - int max_mark + Int tag_mark, + Int max_mark ) { /* === Local variables ================================================== */ - int r ; + Int r ; /* === Check the Row marks ============================================== */ @@ -3330,21 +3506,21 @@ ( /* === Parameters ======================================================= */ - int n_row, - int n_col, + Int n_row, + Int n_col, Colamd_Row Row [], Colamd_Col Col [], - int A [] + Int A [] ) { /* === Local variables ================================================== */ - int r ; - int c ; - int *rp ; - int *rp_end ; - int *cp ; - int *cp_end ; + Int r ; + Int c ; + Int *rp ; + Int *rp_end ; + Int *cp ; + Int *cp_end ; /* === Dump the rows and columns of the matrix ========================== */ @@ -3396,17 +3572,20 @@ char *method ) { + FILE *f ; colamd_debug = 0 ; /* no debug printing */ - - /* get "D" environment variable, which gives the debug printing level */ - if (getenv ("D")) + f = fopen ("debug", "r") ; + if (f == (FILE *) NULL) { - colamd_debug = atoi (getenv ("D")) ; + colamd_debug = 0 ; } - + else + { + fscanf (f, "%d", &colamd_debug) ; + fclose (f) ; + } DEBUG0 (("%s: debug version, D = %d (THIS WILL BE SLOW!)\n", method, colamd_debug)) ; } #endif /* NDEBUG */ -
--- a/liboctave/COLAMD/colamd.h Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/colamd.h Sun Sep 04 12:25:21 2005 +0000 @@ -2,21 +2,18 @@ /* === colamd/symamd prototypes and definitions ============================= */ /* ========================================================================== */ -/* +/* COLAMD Version 2.4 + You must include this file (colamd.h) in any routine that uses colamd, symamd, or the related macros and definitions. Authors: The authors of the code itself are Stefan I. Larimore and Timothy A. - Davis (davis@cise.ufl.edu), University of Florida. The algorithm was + Davis (davis at cise.ufl.edu), University of Florida. The algorithm was developed in collaboration with John Gilbert, Xerox PARC, and Esmond Ng, Oak Ridge National Laboratory. - Date: - - September 8, 2003. Version 2.3. - Acknowledgements: This work was supported by the National Science Foundation, under @@ -24,8 +21,7 @@ Notice: - Copyright (c) 1998-2003 by the University of Florida. - All Rights Reserved. + Copyright (c) 1998-2005, Timothy A. Davis, All Rights Reserved. THIS MATERIAL IS PROVIDED AS IS, WITH ABSOLUTELY NO WARRANTY EXPRESSED OR IMPLIED. ANY USE IS AT YOUR OWN RISK. @@ -52,6 +48,11 @@ #ifndef COLAMD_H #define COLAMD_H +/* make it easy for C++ programs to include COLAMD */ +#ifdef __cplusplus +extern "C" { +#endif + /* ========================================================================== */ /* === Include files ======================================================== */ /* ========================================================================== */ @@ -59,6 +60,35 @@ #include <stdlib.h> /* ========================================================================== */ +/* === COLAMD version ======================================================= */ +/* ========================================================================== */ + +/* COLAMD Version 2.4 and later will include the following definitions. + * As an example, to test if the version you are using is 2.4 or later: + * + * #ifdef COLAMD_VERSION + * if (COLAMD_VERSION >= COLAMD_VERSION_CODE (2,4)) ... + * #endif + * + * This also works during compile-time: + * + * #if defined(COLAMD_VERSION) && (COLAMD_VERSION >= COLAMD_VERSION_CODE (2,4)) + * printf ("This is version 2.4 or later\n") ; + * #else + * printf ("This is an early version\n") ; + * #endif + * + * Versions 2.3 and earlier of COLAMD do not include a #define'd version number. + */ + +#define COLAMD_DATE "Aug. 30, 2005" +#define COLAMD_VERSION_CODE(main,sub) ((main) * 1000 + (sub)) +#define COLAMD_MAIN_VERSION 2 +#define COLAMD_SUB_VERSION 4 +#define COLAMD_VERSION \ + COLAMD_VERSION_CODE(COLAMD_MAIN_VERSION,COLAMD_SUB_VERSION) + +/* ========================================================================== */ /* === Knob and statistics definitions ====================================== */ /* ========================================================================== */ @@ -74,6 +104,9 @@ /* knobs [1] and stats [1]: dense column knob and output statistic. */ #define COLAMD_DENSE_COL 1 +/* knobs [2]: aggressive absorption */ +#define COLAMD_AGGRESSIVE 2 + /* stats [2]: memory defragmentation count output statistic */ #define COLAMD_DEFRAG_COUNT 2 @@ -100,94 +133,6 @@ #define COLAMD_ERROR_out_of_memory (-10) #define COLAMD_ERROR_internal_error (-999) -/* ========================================================================== */ -/* === Row and Column structures ============================================ */ -/* ========================================================================== */ - -/* User code that makes use of the colamd/symamd routines need not directly */ -/* reference these structures. They are used only for the COLAMD_RECOMMENDED */ -/* macro. */ - -typedef struct Colamd_Col_struct -{ - int start ; /* index for A of first row in this column, or DEAD */ - /* if column is dead */ - int length ; /* number of rows in this column */ - union - { - int thickness ; /* number of original columns represented by this */ - /* col, if the column is alive */ - int parent ; /* parent in parent tree super-column structure, if */ - /* the column is dead */ - } shared1 ; - union - { - int score ; /* the score used to maintain heap, if col is alive */ - int order ; /* pivot ordering of this column, if col is dead */ - } shared2 ; - union - { - int headhash ; /* head of a hash bucket, if col is at the head of */ - /* a degree list */ - int hash ; /* hash value, if col is not in a degree list */ - int prev ; /* previous column in degree list, if col is in a */ - /* degree list (but not at the head of a degree list) */ - } shared3 ; - union - { - int degree_next ; /* next column, if col is in a degree list */ - int hash_next ; /* next column, if col is in a hash list */ - } shared4 ; - -} Colamd_Col ; - -typedef struct Colamd_Row_struct -{ - int start ; /* index for A of first col in this row */ - int length ; /* number of principal columns in this row */ - union - { - int degree ; /* number of principal & non-principal columns in row */ - int p ; /* used as a row pointer in init_rows_cols () */ - } shared1 ; - union - { - int mark ; /* for computing set differences and marking dead rows*/ - int first_column ;/* first column in row (used in garbage collection) */ - } shared2 ; - -} Colamd_Row ; - -/* ========================================================================== */ -/* === Colamd recommended memory size ======================================= */ -/* ========================================================================== */ - -/* - The recommended length Alen of the array A passed to colamd is given by - the COLAMD_RECOMMENDED (nnz, n_row, n_col) macro. It returns -1 if any - argument is negative. 2*nnz space is required for the row and column - indices of the matrix. COLAMD_C (n_col) + COLAMD_R (n_row) space is - required for the Col and Row arrays, respectively, which are internal to - colamd. An additional n_col space is the minimal amount of "elbow room", - and nnz/5 more space is recommended for run time efficiency. - - This macro is not needed when using symamd. - - Explicit typecast to int added Sept. 23, 2002, COLAMD version 2.2, to avoid - gcc -pedantic warning messages. -*/ - -#define COLAMD_C(n_col) ((int) (((n_col) + 1) * sizeof (Colamd_Col) / sizeof (int))) -#define COLAMD_R(n_row) ((int) (((n_row) + 1) * sizeof (Colamd_Row) / sizeof (int))) - -#define COLAMD_RECOMMENDED(nnz, n_row, n_col) \ -( \ -((nnz) < 0 || (n_row) < 0 || (n_col) < 0) \ -? \ - (-1) \ -: \ - (2 * (nnz) + COLAMD_C (n_col) + COLAMD_R (n_row) + (n_col) + ((nnz) / 5)) \ -) /* ========================================================================== */ /* === Prototypes of user-callable routines ================================= */ @@ -201,11 +146,24 @@ int n_col /* number of columns in A */ ) ; +long colamd_l_recommended /* returns recommended value of Alen, */ + /* or (-1) if input arguments are erroneous */ +( + long nnz, /* nonzeros in A */ + long n_row, /* number of rows in A */ + long n_col /* number of columns in A */ +) ; + void colamd_set_defaults /* sets default parameters */ ( /* knobs argument is modified on output */ double knobs [COLAMD_KNOBS] /* parameter settings for colamd */ ) ; +void colamd_l_set_defaults /* sets default parameters */ +( /* knobs argument is modified on output */ + double knobs [COLAMD_KNOBS] /* parameter settings for colamd */ +) ; + int colamd /* returns (1) if successful, (0) otherwise*/ ( /* A and p arguments are modified on output */ int n_row, /* number of rows in A */ @@ -217,6 +175,17 @@ int stats [COLAMD_STATS] /* colamd output statistics and error codes */ ) ; +long colamd_l /* returns (1) if successful, (0) otherwise*/ +( /* A and p arguments are modified on output */ + long n_row, /* number of rows in A */ + long n_col, /* number of columns in A */ + long Alen, /* size of the array A */ + long A [], /* row indices of A, of size Alen */ + long p [], /* column pointers of A, of size n_col+1 */ + double knobs [COLAMD_KNOBS],/* parameter settings for colamd */ + long stats [COLAMD_STATS] /* colamd output statistics and error codes */ +) ; + int symamd /* return (1) if OK, (0) otherwise */ ( int n, /* number of rows and columns of A */ @@ -233,14 +202,46 @@ /* mxFree (for MATLAB mexFunction) */ ) ; +long symamd_l /* return (1) if OK, (0) otherwise */ +( + long n, /* number of rows and columns of A */ + long A [], /* row indices of A */ + long p [], /* column pointers of A */ + long perm [], /* output permutation, size n_col+1 */ + double knobs [COLAMD_KNOBS], /* parameters (uses defaults if NULL) */ + long stats [COLAMD_STATS], /* output statistics and error codes */ + void * (*allocate) (size_t, size_t), + /* pointer to calloc (ANSI C) or */ + /* mxCalloc (for MATLAB mexFunction) */ + void (*release) (void *) + /* pointer to free (ANSI C) or */ + /* mxFree (for MATLAB mexFunction) */ +) ; + void colamd_report ( int stats [COLAMD_STATS] ) ; +void colamd_l_report +( + long stats [COLAMD_STATS] +) ; + void symamd_report ( int stats [COLAMD_STATS] ) ; +void symamd_l_report +( + long stats [COLAMD_STATS] +) ; + +extern int (*colamd_printf) (const char *, ...) ; + +#ifdef __cplusplus +} +#endif + #endif /* COLAMD_H */
--- a/liboctave/COLAMD/colamd.m Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/colamd.m Sun Sep 04 12:25:21 2005 +0000 @@ -1,25 +1,70 @@ function [p,stats] = colamd (S, knobs) %COLAMD Column approximate minimum degree permutation. -% P = COLAMD (S) returns the column approximate minimum degree permutation -% vector for the sparse matrix S. For a non-symmetric matrix S, S (:,P) +% P = COLAMD(S) returns the column approximate minimum degree permutation +% vector for the sparse matrix S. For a non-symmetric matrix S, S(:,P) % tends to have sparser LU factors than S. The Cholesky factorization of -% S (:,P)' * S (:,P) also tends to be sparser than that of S'*S. COLAMD -% tends to be faster than COLMMD and tends to return a better ordering. +% S(:,P)'*S(:,P) also tends to be sparser than that of S'*S. The ordering +% is followed by a column elimination tree post-ordering. % -% See also COLMMD, COLPERM, SPPARMS, SYMAMD, SYMMMD, SYMRCM. +% See also AMD, CCOLAMD, CSYMAMD, SYMAMD, COLPERM, SYMRCM. % % Usage: P = colamd (S) -% P = colamd (S, knobs) -% [P, stats] = colamd (S) % [P, stats] = colamd (S, knobs) % -% knobs is an optional two-element input vector. If S is m-by-n, then -% rows with more than (knobs (1))*n entries are ignored. Columns with -% more than (knobs (2))*m entries are removed prior to ordering, and -% ordered last in the output permutation P. If the knobs parameter is not -% present, then 0.5 is used instead, for both knobs (1) and knobs (2). -% knobs (3) controls the printing of statistics and error messages. +% knobs is an optional one- to three-element input vector. If S is m-by-n, +% then rows with more than max(16,knobs(1)*sqrt(n)) entries are ignored. +% Columns with more than max(16,knobs(2)*sqrt(min(m,n))) entries are +% removed prior to ordering, and ordered last in the output permutation P. +% Only completely dense rows or columns are removed if knobs(1) and knobs(2) +% are < 0, respectively. If knobs(3) is nonzero, stats and knobs are +% printed. The default is knobs = [10 10 0]. Note that knobs differs from +% earlier versions of colamd. +% +% Type the command "type colamd" for a description of the optional stats +% output and for the copyright information. +% +% Authors: S. Larimore and T. Davis, University of Florida. Developed in +% collaboration with J. Gilbert and E. Ng. Version 2.4. +% +% Acknowledgements: This work was supported by the National Science +% Foundation, under grants DMS-9504974 and DMS-9803599. + +% Notice: +% +% Copyright (c) 1998-2005, Timothy A. Davis, All Rights Reserved. +% +% See http://www.cise.ufl.edu/research/sparse/colamd (the colamd.c +% file) for the License. +% +% Availability: % +% The colamd, symamd, amd, ccolamd, and csymamd are available at +% http://www.cise.ufl.edu/research/sparse + +%------------------------------------------------------------------------------- +% Perform the colamd ordering: +%------------------------------------------------------------------------------- + +if (nargout <= 1 && nargin == 1) + p = colamdmex (S) ; +elseif (nargout <= 1 && nargin == 2) + p = colamdmex (S, knobs) ; +elseif (nargout == 2 && nargin == 1) + [p, stats] = colamdmex (S) ; +elseif (nargout == 2 && nargin == 2) + [p, stats] = colamdmex (S, knobs) ; +else + error ('MATLAB:colamd:WrongInputOrOutputNumber',... + 'colamd: incorrect number of input and/or output arguments') ; +end + +%------------------------------------------------------------------------------- +% column elimination tree post-ordering: +%------------------------------------------------------------------------------- + +[ignore, q] = sparsfun ('coletree', S (:,p)) ; +p = p (q) ; + % stats is an optional 20-element output vector that provides data about the % ordering and the validity of the input matrix S. Ordering statistics are % in stats (1:3). stats (1) and stats (2) are the number of dense or empty @@ -50,60 +95,4 @@ % % stats (8:20) is always zero in the current version of COLAMD (reserved % for future use). -% -% The ordering is followed by a column elimination tree post-ordering. -% -% Authors: -% -% The authors of the code itself are Stefan I. Larimore and Timothy A. -% Davis (davis@cise.ufl.edu), University of Florida. The algorithm was -% developed in collaboration with John Gilbert, Xerox PARC, and Esmond -% Ng, Oak Ridge National Laboratory. -% -% Date: -% -% September 8, 2003. Version 2.3. -% -% Acknowledgements: -% -% This work was supported by the National Science Foundation, under -% grants DMS-9504974 and DMS-9803599. -% Notice: -% -% Copyright (c) 1998-2003 by the University of Florida. -% All Rights Reserved. -% -% See http://www.cise.ufl.edu/research/sparse/colamd (the colamd.c -% file) for the License. -% -% Availability: -% -% The colamd/symamd library is available at -% -% http://www.cise.ufl.edu/research/sparse/colamd/ -% - -%------------------------------------------------------------------------------- -% Perform the colamd ordering: -%------------------------------------------------------------------------------- - -if (nargout <= 1 & nargin == 1) - p = colamdmex (S) ; -elseif (nargout <= 1 & nargin == 2) - p = colamdmex (S, knobs) ; -elseif (nargout == 2 & nargin == 1) - [p, stats] = colamdmex (S) ; -elseif (nargout == 2 & nargin == 2) - [p, stats] = colamdmex (S, knobs) ; -else - error ('colamd: incorrect number of input and/or output arguments') ; -end - -%------------------------------------------------------------------------------- -% column elimination tree post-ordering: -%------------------------------------------------------------------------------- - -[ignore, q] = sparsfun ('coletree', S (:,p)) ; -p = p (q) ; -
--- a/liboctave/COLAMD/colamd_demo.m Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/colamd_demo.m Sun Sep 04 12:25:21 2005 +0000 @@ -14,16 +14,17 @@ % % For a description of the methods used, see the colamd.c file. % -% September 8, 2003. Version 2.3. +% COLAMD Version 2.4. % http://www.cise.ufl.edu/research/sparse/colamd/ % +% Minor changes: in MATLAB 7, symmmd and colmmd are flagged as "obsolete". +% This demo checks if they exist, so it should still work when they are removed. + %------------------------------------------------------------------------------- % Print the introduction, the help info, and compile the mexFunctions %------------------------------------------------------------------------------- -more on - fprintf (1, '\n-----------------------------------------------------------\n') ; fprintf (1, 'Colamd/symamd demo.') ; fprintf (1, '\n-----------------------------------------------------------\n') ; @@ -39,15 +40,6 @@ fprintf (1, '\n-----------------------------------------------------------\n') ; help symamd ; -fprintf (1, '\n-----------------------------------------------------------\n') ; -fprintf (1, 'colamd and symamd mexFunctions will now be compiled.') ; -fprintf (1, '\n-----------------------------------------------------------\n') ; -s = input ('\n\nHit any key (or type ''n'' to skip compilation): ', 's') ; -if (~strcmp (s, 'n')) - disp ('mex -O colamdmex.c colamd.c') ; mex -O colamdmex.c colamd.c - disp ('mex -O symamdmex.c colamd.c') ; mex -O symamdmex.c colamd.c -end - %------------------------------------------------------------------------------- % Solving Ax=b %------------------------------------------------------------------------------- @@ -77,6 +69,7 @@ fprintf (1, '\nFlop count for [L,U,P] = lu (A*Q): %d\n', fl) ; fprintf (1, 'residual: %e\n', norm (A*x-b)); +try fprintf (1, '\n\nSolving via lu (PAQ = LU), where Q is from colmmd:\n') ; q = colmmd (A) ; I = speye (n) ; @@ -86,6 +79,9 @@ x = Q * (U \ (L \ (P * b))) ; fprintf (1, '\nFlop count for [L,U,P] = lu (A*Q): %d\n', fl) ; fprintf (1, 'residual: %e\n', norm (A*x-b)); +catch +fprintf (1, 'colmmd is obsolete\n') ; +end fprintf (1, '\n\nSolving via lu (PA = LU), without regard for sparsity:\n') ; [L,U,P] = lu (A) ; @@ -123,6 +119,7 @@ fprintf (1, 'nz in Cholesky factors of A(:,p)''A(:,p): %d\n', sum (lnz)) ; fprintf (1, 'flop count for Cholesky of A(:,p)''A(:,p): %d\n', sum (lnz.^2)) ; +try tic ; p = colmmd (A) ; t = toc ; @@ -131,6 +128,9 @@ fprintf (1, 'colmmd ordering quality: \n') ; fprintf (1, 'nz in Cholesky factors of A(:,p)''A(:,p): %d\n', sum (lnz)) ; fprintf (1, 'flop count for Cholesky of A(:,p)''A(:,p): %d\n', sum (lnz.^2)) ; +catch +fprintf (1, 'colmmd is obsolete\n') ; +end %------------------------------------------------------------------------------- % Large demo for symamd @@ -157,6 +157,7 @@ fprintf (1, 'nz in Cholesky factors of A(p,p): %d\n', sum (lnz)) ; fprintf (1, 'flop count for Cholesky of A(p,p): %d\n', sum (lnz.^2)) ; +try tic ; p = symmmd (A) ; t = toc ; @@ -165,6 +166,6 @@ fprintf (1, 'symmmd ordering quality: \n') ; fprintf (1, 'nz in Cholesky factors of A(p,p): %d\n', sum (lnz)) ; fprintf (1, 'flop count for Cholesky of A(p,p): %d\n', sum (lnz.^2)) ; - -more off - +catch +fprintf (1, 'symmmd is obsolete\n') ; +end
--- a/liboctave/COLAMD/colamd_example.c Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/colamd_example.c Sun Sep 04 12:25:21 2005 +0000 @@ -2,7 +2,8 @@ /* === colamd and symamd example ============================================ */ /* ========================================================================== */ -/* +/* COLAMD Version 2.4. + colamd example of use, to order the columns of a 5-by-4 matrix with 11 nonzero entries in the following nonzero pattern, with default knobs. @@ -24,8 +25,6 @@ (where x denotes a nonzero value). - September 8, 2003. Version 2.3. - See http://www.cise.ufl.edu/research/sparse/colamd/ (the colamd.c file) for the routines this program calls, and for the License. */ @@ -38,12 +37,12 @@ #define A_NNZ 11 #define A_NROW 5 #define A_NCOL 4 -#define ALEN (COLAMD_RECOMMENDED (A_NNZ, A_NCOL, A_NROW)) +#define ALEN 150 #define B_NNZ 4 #define B_N 5 -int main (int argc, char **argv) +int main (void) { /* ====================================================================== */
--- a/liboctave/COLAMD/colamd_example.out Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/colamd_example.out Sun Sep 04 12:25:21 2005 +0000 @@ -14,15 +14,16 @@ Column 3, with 2 entries: row 1 row 3 -colamd: OK. -colamd: number of dense or empty rows ignored: 1 -colamd: number of dense or empty columns ignored: 2 + +colamd version 2.4, Aug. 30, 2005: OK. +colamd: number of dense or empty rows ignored: 0 +colamd: number of dense or empty columns ignored: 0 colamd: number of garbage collections performed: 0 colamd column ordering: 1st column: 1 -2nd column: 3 -3rd column: 0 -4th column: 2 +2nd column: 0 +3rd column: 2 +4th column: 3 symamd 5-by-5 input matrix: @@ -36,7 +37,8 @@ Column 3, with 1 entries: row 4 Column 4, with 0 entries: -symamd: OK. + +symamd version 2.4, Aug. 30, 2005: OK. symamd: number of dense or empty rows ignored: 0 symamd: number of dense or empty columns ignored: 0 symamd: number of garbage collections performed: 0
--- a/liboctave/COLAMD/colamd_test.m Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/colamd_test.m Sun Sep 04 12:25:21 2005 +0000 @@ -1,23 +1,17 @@ -function colamd_test (spumoni) +function colamd_test % -% colamd_test (spumoni) +% colamd_test % % COLAMD and SYMAMD testing function. Here we try to give colamd and symamd % every possible type of matrix and erroneous input that they may encounter. % We want either a valid permutation returned or we want them to fail -% gracefully. If the optional spumoni argument is present, additional -% diagnostic messages are printed during the test: -% -% spumoni: -% 0: very little, just a progress meter (default) -% 1: lots -% 2: far too much +% gracefully. % % You are prompted as to whether or not the colamd and symand routines and % the test mexFunctions are to be compiled. % % Tim Davis -% September 8, 2003. Version 2.3. +% COLAMD Version 2.4. % http://www.cise.ufl.edu/research/sparse/colamd/ help colamd_test @@ -55,31 +49,34 @@ rand ('state', 0) ; randn ('state', 0) ; -if (nargin < 1) - spumoni = 0 ; -end -old = spparms ('spumoni') ; -spparms ('spumoni', spumoni) ; + +A = sprandn (500,500,0.4) ; + +p = colamd (A, [10 10 1]) ; check_perm (p, A) ; +p = colamd (A, [2 7 1]) ; check_perm (p, A) ; +p = symamd (A, [10 1]) ; check_perm (p, A) ; +p = symamd (A, [7 1]) ; check_perm (p, A) ; +p = symamd (A, [4 1]) ; check_perm (p, A) ; fprintf ('Null matrices') ; A = zeros (0,0) ; A = sparse (A) ; -[p, stats] = colamd (A, [.5 .5 spumoni]) ; +[p, stats] = colamd (A, [10 10 0]) ; check_perm (p, A) ; -[p, stats] = symamd (A, [.5 spumoni]) ; +[p, stats] = symamd (A, [10 0]) ; check_perm (p, A) ; A = zeros (0, 100) ; A = sparse (A) ; -[p, stats] = colamd (A, [.5 .5 spumoni]) ; +[p, stats] = colamd (A, [10 10 0]) ; check_perm (p, A) ; A = zeros (100, 0) ; A = sparse (A) ; -[p, stats] = colamd (A, [.5 .5 spumoni]) ; +[p, stats] = colamd (A, [10 10 0]) ; check_perm (p, A) ; fprintf (' OK\n') ; @@ -92,20 +89,19 @@ A = rand_matrix (1000, 1000, 1, 10, 20) ; [m n] = size (A) ; + for tol = [0:.1:2 3:20 1e6] - for tol = [0 0.001 0.005 0.01:0.01:0.10 0.5:.1:1] - - [p, stats] = colamd (A, [tol tol spumoni]) ; + [p, stats] = colamd (A, [tol tol 0]) ; check_perm (p, A) ; B = A + A' ; - [p, stats] = symamd (B, [tol spumoni]) ; + [p, stats] = symamd (B, [tol 0]) ; check_perm (p, A) ; - [p, stats] = colamd (A, [tol 1 spumoni]) ; + [p, stats] = colamd (A, [tol 1 0]) ; check_perm (p, A) ; - [p, stats] = colamd (A, [1 tol spumoni]) ; + [p, stats] = colamd (A, [1 tol 0]) ; check_perm (p, A) ; fprintf ('.') ; @@ -142,31 +138,31 @@ for err = 1:13 - p = Tcolamd (A, [1 1 spumoni 0 err]) ; + p = Tcolamd (A, [n n 0 0 err]) ; if p ~= -1 check_perm (p, A) ; end if (err == 1) % check different (valid) input args to colamd - p = Acolamd (A, spumoni) ; + p = Acolamd (A) ; - p2 = Acolamd (A, spumoni, [.5 .5 spumoni 0 0]) ; + p2 = Acolamd (A, [10 10 0 0 0]) ; if (any (p ~= p2)) error ('colamd: mismatch 1!') ; end - [p2 stats] = Acolamd (A, spumoni) ; + [p2 stats] = Acolamd (A) ; if (any (p ~= p2)) error ('colamd: mismatch 2!') ; end - [p2 stats] = Acolamd (A, spumoni, [.5 .5 spumoni 0 0]) ; + [p2 stats] = Acolamd (A, [10 10 0 0 0]) ; if (any (p ~= p2)) error ('colamd: mismatch 3!') ; end end B = A'*A ; - p = Tsymamd (B, [1 spumoni err]) ; + p = Tsymamd (B, [n 0 err]) ; if p ~= -1 check_perm (p, A) ; end @@ -174,17 +170,17 @@ if (err == 1) % check different (valid) input args to symamd - p = Asymamd (B, spumoni) ; + p = Asymamd (B) ; check_perm (p, A) ; - p2 = Asymamd (B, spumoni, [.5 spumoni 0]) ; + p2 = Asymamd (B, [10 0 0]) ; if (any (p ~= p2)) error ('symamd: mismatch 1!') ; end - [p2 stats] = Asymamd (B, spumoni) ; + [p2 stats] = Asymamd (B) ; if (any (p ~= p2)) error ('symamd: mismatch 2!') ; end - [p2 stats] = Asymamd (B, spumoni, [.5 spumoni 0]) ; + [p2 stats] = Asymamd (B, [10 0 0]) ; if (any (p ~= p2)) error ('symamd: mismatch 3!') ; end @@ -213,12 +209,17 @@ A (:, null_col) = 0 ; % Order the matrix and make sure that the null columns are ordered last. - [p, stats] = colamd (A, [1 1 spumoni]) ; + [p, stats] = colamd (A, [1e6 1e6 0]) ; check_perm (p, A) ; - if (stats (2) ~= 5) - error ('colamd: wrong number of null columns') ; - end +% if (stats (2) ~= 5) +% stats (2) +% error ('colamd: wrong number of null columns') ; +% end + + % find all null columns in A + null_col = find (sum (spones (A), 1) == 0) ; + nnull = length (null_col) ; if (any (null_col ~= p ((n-4):n))) error ('colamd: Null cols are not ordered last in natural order') ; end @@ -246,7 +247,7 @@ A (null_col, :) = 0 ; % Order the matrix and make sure that the null rows/cols are ordered last. - [p,stats] = symamd (A, [1 spumoni]) ; + [p,stats] = symamd (A, [10 0]) ; check_perm (p, A) ; % find actual number of null rows and columns @@ -282,7 +283,7 @@ null_row = sort (null_row (1:5)) ; A (null_row, :) = 0 ; - p = colamd (A, [.5 .5 spumoni]) ; + p = colamd (A, [10 10 0]) ; check_perm (p, A) ; if (stats (1) ~= 5) error ('colamd: wrong number of null rows') ; @@ -293,25 +294,20 @@ fprintf ('\ncolamd and symamd: all tests passed\n\n') ; -spparms ('spumoni', old) ; %------------------------------------------------------------------------------- % Acolamd: compare colamd and Tcolamd results %------------------------------------------------------------------------------- -function [p,stats] = Acolamd (S, spumoni, knobs) - -if (nargin < 2) - spumoni = 0 ; -end +function [p,stats] = Acolamd (S, knobs) if (nargin < 3) if (nargout == 1) [p] = colamd (S) ; - [p1] = Tcolamd (S, [.5 .5 spumoni 0 0]) ; + [p1] = Tcolamd (S, [10 10 0 0 0]) ; else [p, stats] = colamd (S) ; - [p1, stats1] = Tcolamd (S, [.5 .5 spumoni 0 0]) ; + [p1, stats1] = Tcolamd (S, [10 10 0 0 0]) ; end else if (nargout == 1) @@ -335,19 +331,15 @@ % Asymamd: compare symamd and Tsymamd results %------------------------------------------------------------------------------- -function [p,stats] = Asymamd (S, spumoni, knobs) - -if (nargin < 2) - spumoni = 0 ; -end +function [p,stats] = Asymamd (S, knobs) if (nargin < 3) if (nargout == 1) [p] = symamd (S) ; - [p1] = Tsymamd (S, [.5 spumoni 0]) ; + [p1] = Tsymamd (S, [10 0 0]) ; else [p, stats] = symamd (S) ; - [p1, stats1] = Tsymamd (S, [.5 spumoni 0]) ; + [p1, stats1] = Tsymamd (S, [10 0 0]) ; end else if (nargout == 1) @@ -518,4 +510,3 @@ p = p (q) ; check_perm (p, S) ; end -
--- a/liboctave/COLAMD/colamdmex.c Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/colamdmex.c Sun Sep 04 12:25:21 2005 +0000 @@ -2,76 +2,22 @@ /* === colamd mexFunction =================================================== */ /* ========================================================================== */ -/* +/* COLAMD Version 2.4 + Usage: P = colamd (A) ; - - P = colamd (A, knobs) ; - - [ P, stats ] = colamd (A) ; - [ P, stats ] = colamd (A, knobs) ; - Returns a permutation vector P such that the LU factorization of A (:,P) - tends to be sparser than that of A. The Cholesky factorization of - (A (:,P))'*(A (:,P)) will also tend to be sparser than that of A'*A. - This routine provides the same functionality as COLMMD, but is much faster - and returns a better permutation vector. Note that the COLMMD m-file in - MATLAB 5.2 also performs a column elimination tree post-ordering. This - mexFunction does not do this post-ordering. This mexFunction is a - replacement for the p = sparsfun ('colmmd', A) operation. - - The knobs and stats vectors are optional: - - knobs (1) rows with more than (knobs (1))*n_col entries - are removed prior to ordering. If knobs is not present, - then the default is used (0.5). - - knobs (2) columns with more than (knobs (2))*n_row entries - are removed prior to ordering, and placed last in the - column permutation. If knobs is not present, - then the default is used (0.5). - - knobs (3) print level, similar to spparms ('spumoni') - - stats (1) the number of dense (or empty) rows ignored - - stats (2) the number of dense (or empty) columms. These - are ordered last, in their natural order. - - stats (3) the number of garbage collections performed. - - stats (4) return status: - - 0: matrix is a valid MATLAB matrix. - - 1: matrix has duplicate entries or unsorted columns. - This should not occur in a valid MATLAB matrix, - but the ordering proceeded anyway by sorting the - row indices in each column and by ignoring the - duplicate row indices in each column. See - stats (5:7) for more information. - - stats (5) highest numbered column that is unsorted or has - duplicate entries (zero if none) - - stats (6) last seen duplicate or unsorted row index - (zero if none) - - stats (7) number of duplicate or unsorted row indices + see colamd.m for a description. Authors: The authors of the code itself are Stefan I. Larimore and Timothy A. - Davis (davis@cise.ufl.edu), University of Florida. The algorithm was + Davis (davis at cise.ufl.edu), University of Florida. The algorithm was developed in collaboration with John Gilbert, Xerox PARC, and Esmond Ng, Oak Ridge National Laboratory. - Date: - - September 8, 2003. Version 2.3. - Acknowledgements: This work was supported by the National Science Foundation, under @@ -79,8 +25,7 @@ Notice: - Copyright (c) 1998-2003 by the University of Florida. - All Rights Reserved. + Copyright (c) 1998-2005, Timothy A. Davis, All Rights Reserved. See http://www.cise.ufl.edu/research/sparse/colamd (the colamd.c file) for the License. @@ -138,6 +83,8 @@ int spumoni ; /* verbosity variable */ int stats [COLAMD_STATS] ; /* stats for colamd */ + colamd_printf = mexPrintf ; /* COLAMD printf routine */ + /* === Check inputs ===================================================== */ if (nrhs < 1 || nrhs > 2 || nlhs < 0 || nlhs > 2) @@ -156,18 +103,38 @@ { in_knobs = mxGetPr (prhs [1]) ; i = mxGetNumberOfElements (prhs [1]) ; - if (i > 0) knobs [COLAMD_DENSE_ROW] = in_knobs [COLAMD_DENSE_ROW] ; - if (i > 1) knobs [COLAMD_DENSE_COL] = in_knobs [COLAMD_DENSE_COL] ; - if (i > 2) spumoni = (int) in_knobs [2] ; + if (i > 0) knobs [COLAMD_DENSE_ROW] = in_knobs [0] ; + if (i > 1) knobs [COLAMD_DENSE_COL] = in_knobs [1] ; + if (i > 2) spumoni = (int) (in_knobs [2] != 0) ; } /* print knob settings if spumoni is set */ - if (spumoni > 0) + if (spumoni) { - mexPrintf ("colamd: dense row fraction: %f\n", - knobs [COLAMD_DENSE_ROW]) ; - mexPrintf ("colamd: dense col fraction: %f\n", - knobs [COLAMD_DENSE_COL]) ; + mexPrintf ("\ncolamd version %d.%d, %s:\n", + COLAMD_MAIN_VERSION, COLAMD_SUB_VERSION, COLAMD_DATE) ; + if (knobs [COLAMD_DENSE_ROW] >= 0) + { + mexPrintf ("knobs(1): %g, rows with > max(16,%g*sqrt(size(A,2)))" + " entries removed\n", in_knobs [0], knobs [COLAMD_DENSE_ROW]) ; + } + else + { + mexPrintf ("knobs(1): %g, only completely dense rows removed\n", + in_knobs [0]) ; + } + if (knobs [COLAMD_DENSE_COL] >= 0) + { + mexPrintf ("knobs(2): %g, cols with > max(16,%g*sqrt(min(size(A)))" + " entries removed\n", in_knobs [1], knobs [COLAMD_DENSE_COL]) ; + } + else + { + mexPrintf ("knobs(2): %g, only completely dense columns removed\n", + in_knobs [1]) ; + } + mexPrintf ("knobs(3): %g, statistics and knobs printed\n", + in_knobs [2]) ; } /* === If A is full, convert to a sparse matrix ========================= */ @@ -227,8 +194,8 @@ /* === Return the stats vector ========================================== */ - /* print stats if spumoni > 0 */ - if (spumoni > 0) + /* print stats if spumoni is set */ + if (spumoni) { colamd_report (stats) ; }
--- a/liboctave/COLAMD/colamdtestmex.c Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/colamdtestmex.c Sun Sep 04 12:25:21 2005 +0000 @@ -2,7 +2,8 @@ /* === colamdtest mexFunction =============================================== */ /* ========================================================================== */ -/* +/* COLAMD Version 2.4. + This MATLAB mexFunction is for testing only. It is not meant for production use. See colamdmex.c instead. @@ -10,28 +11,11 @@ [ P, stats ] = colamdtest (A, knobs) ; - Returns a permutation vector P such that the LU factorization of A (:,P) - tends to be sparser than that of A. The Cholesky factorization of - (A (:,P))'*(A (:,P)) will also tend to be sparser than that of A'*A. - This routine provides the same functionality as COLMMD, but is much faster - and returns a better permutation vector. Note that the COLMMD m-file in - MATLAB 5.2 also performs a column elimination tree post-ordering. This - mexFunction does not do this post-ordering. This mexFunction is a - replacement for the p = sparsfun ('colmmd', A) operation. - - The knobs and stats vectors are optional: + See colamd.m for a description. knobs is required. - knobs (1) rows with more than (knobs (1))*n_col entries - are removed prior to ordering. If knobs is not present, - then the default is used (0.5). - - knobs (2) columns with more than (knobs (2))*n_row entries - are removed prior to ordering, and placed last in the - column permutation. If knobs is not present, - then the default is used (0.5). - - knobs (3) print level, similar to spparms ('spumoni') - + knobs (1) dense row control + knobs (2) dense column control + knobs (3) spumoni knobs (4) for testing only. Controls the workspace used by colamd. @@ -39,43 +23,13 @@ jumbled prior to calling colamd, to test its error handling capability. - stats (1) the number of dense (or empty) rows ignored - - stats (2) the number of dense (or empty) columms. These - are ordered last, in their natural order. - - stats (3) the number of garbage collections performed. - - stats (4) return status: - - 0: matrix is a valid MATLAB matrix. - - 1: matrix has duplicate entries or unsorted columns. - This should not occur in a valid MATLAB matrix, - but the ordering proceeded anyway by sorting the - row indices in each column and by ignoring the - duplicate row indices in each column. See - stats (5:7) for more information. - - stats (5) highest numbered column that is unsorted or has - duplicate entries (zero if none) - - stats (6) last seen duplicate or unsorted row index - (zero if none) - - stats (7) number of duplicate or unsorted row indices - Authors: The authors of the code itself are Stefan I. Larimore and Timothy A. - Davis (davis@cise.ufl.edu), University of Florida. The algorithm was + Davis (davis at cise.ufl.edu), University of Florida. The algorithm was developed in collaboration with John Gilbert, Xerox PARC, and Esmond Ng, Oak Ridge National Laboratory. - Date: - - September 8, 2003. Version 2.3. - Acknowledgements: This work was supported by the National Science Foundation, under @@ -83,8 +37,7 @@ Notice: - Copyright (c) 1998-2003 by the University of Florida. - All Rights Reserved. + Copyright (c) 1998-2005, Timothy A. Davis, All Rights Reserved. See http://www.cise.ufl.edu/research/sparse/colamd (the colamd.c file) for the License. @@ -157,6 +110,8 @@ int *stats ; stats = stats2 ; + colamd_printf = mexPrintf ; /* COLAMD printf routine */ + /* === Check inputs ===================================================== */ if (nrhs < 1 || nrhs > 2 || nlhs < 0 || nlhs > 2) @@ -170,7 +125,7 @@ mexErrMsgTxt ("colamdtest: knobs are required") ; } /* for testing we require all 5 knobs */ - if (mxGetNumberOfElements (prhs [1]) < 5) + if (mxGetNumberOfElements (prhs [1]) != 5) { mexErrMsgTxt ("colamd: must have all 5 knobs for testing") ; } @@ -185,18 +140,38 @@ { in_knobs = mxGetPr (prhs [1]) ; i = mxGetNumberOfElements (prhs [1]) ; - if (i > 0) knobs [COLAMD_DENSE_ROW] = in_knobs [COLAMD_DENSE_ROW] ; - if (i > 1) knobs [COLAMD_DENSE_COL] = in_knobs [COLAMD_DENSE_COL] ; + if (i > 0) knobs [COLAMD_DENSE_ROW] = in_knobs [0] ; + if (i > 1) knobs [COLAMD_DENSE_COL] = in_knobs [1] ; if (i > 2) spumoni = (int) in_knobs [2] ; } /* print knob settings if spumoni is set */ - if (spumoni > 0) + if (spumoni) { - mexPrintf ("colamd: dense row fraction: %f\n", - knobs [COLAMD_DENSE_ROW]) ; - mexPrintf ("colamd: dense col fraction: %f\n", - knobs [COLAMD_DENSE_COL]) ; + mexPrintf ("\ncolamd version %d.%d, %s:\n", + COLAMD_MAIN_VERSION, COLAMD_SUB_VERSION, COLAMD_DATE) ; + if (knobs [COLAMD_DENSE_ROW] >= 0) + { + mexPrintf ("knobs(1): %g, rows with > max(16,%g*sqrt(size(A,2)))" + " entries removed\n", in_knobs [0], knobs [COLAMD_DENSE_ROW]) ; + } + else + { + mexPrintf ("knobs(1): %g, only completely dense rows removed\n", + in_knobs [0]) ; + } + if (knobs [COLAMD_DENSE_COL] >= 0) + { + mexPrintf ("knobs(2): %g, cols with > max(16,%g*sqrt(min(size(A)))" + " entries removed\n", in_knobs [1], knobs [COLAMD_DENSE_COL]) ; + } + else + { + mexPrintf ("knobs(2): %g, only completely dense columns removed\n", + in_knobs [1]) ; + } + mexPrintf ("knobs(3): %g, statistics and knobs printed\n", + in_knobs [2]) ; } /* === If A is full, convert to a sparse matrix ========================= */
--- a/liboctave/COLAMD/symamd.m Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/symamd.m Sun Sep 04 12:25:21 2005 +0000 @@ -1,25 +1,71 @@ function [p, stats] = symamd (S, knobs) %SYMAMD Symmetric approximate minimum degree permutation. -% P = SYMAMD (S) for a symmetric positive definite matrix S, returns the +% P = SYMAMD(S) for a symmetric positive definite matrix S, returns the % permutation vector p such that S(p,p) tends to have a sparser Cholesky % factor than S. Sometimes SYMAMD works well for symmetric indefinite -% matrices too. SYMAMD tends to be faster than SYMMMD and tends to return -% a better ordering. The matrix S is assumed to be symmetric; only the +% matrices too. The matrix S is assumed to be symmetric; only the % strictly lower triangular part is referenced. S must be square. +% Note that p = amd(S) is much faster and generates comparable orderings. +% The ordering is followed by an elimination tree post-ordering. % -% See also COLMMD, COLPERM, SPPARMS, COLAMD, SYMMMD, SYMRCM. +% See also AMD, CCOLAMD, CSYMAMD, COLAMD, COLPERM, SYMRCM. % % Usage: P = symamd (S) -% P = symamd (S, knobs) -% [P, stats] = symamd (S) % [P, stats] = symamd (S, knobs) % -% knobs is an optional input argument. If S is n-by-n, then rows and -% columns with more than knobs(1)*n entries are removed prior to ordering, -% and ordered last in the output permutation P. If the knobs parameter is not -% present, then the default of 0.5 is used instead. knobs (2) controls the -% printing of statistics and error messages. +% knobs is an optional one- to two-element input vector. If S is n-by-n, +% then rows and columns with more than max(16,knobs(1)*sqrt(n)) entries are +% removed prior to ordering, and ordered last in the output permutation P. +% No rows/columns are removed if knobs(1)<0. If knobs(2) is nonzero, stats +% and knobs are printed. The default is knobs = [10 0]. Note that knobs +% differs from earlier versions of symamd. +% +% Type the command "type symamd" for a description of the optional stats +% output and for the copyright information. +% +% Authors: S. Larimore and T. Davis, University of Florida. Developed in +% collaboration with J. Gilbert and E. Ng. Version 2.4. +% +% Acknowledgements: This work was supported by the National Science +% Foundation, under grants DMS-9504974 and DMS-9803599. + +% Notice: +% +% Copyright (c) 1998-2005, Timothy A. Davis, All Rights Reserved. +% +% See http://www.cise.ufl.edu/research/sparse/colamd (the colamd.c +% file) for the License. +% +% Availability: % +% The colamd, symamd, amd, ccolamd, and csymamd are available at +% http://www.cise.ufl.edu/research/sparse + +%------------------------------------------------------------------------------- +% perform the symamd ordering: +%------------------------------------------------------------------------------- + +if (nargout <= 1 && nargin == 1) + p = symamdmex (S) ; +elseif (nargout <= 1 && nargin == 2) + p = symamdmex (S, knobs) ; +elseif (nargout == 2 && nargin == 1) + [p, stats] = symamdmex (S) ; +elseif (nargout == 2 && nargin == 2) + [p, stats] = symamdmex (S, knobs) ; +else + error('MATLAB:symamd:WrongInputOrOutputNumber',... + 'symamd: incorrect number of input and/or output arguments.') ; +end + +%------------------------------------------------------------------------------- +% symmetric elimination tree post-ordering: +%------------------------------------------------------------------------------- + +[ignore, q] = sparsfun ('symetree', S (p,p)) ; +p = p (q) ; + + % stats is an optional 20-element output vector that provides data about the % ordering and the validity of the input matrix S. Ordering statistics are % in stats (1:3). stats (1) = stats (2) is the number of dense or empty @@ -50,59 +96,3 @@ % % stats (8:20) is always zero in the current version of SYMAMD (reserved % for future use). -% -% The ordering is followed by a column elimination tree post-ordering. -% -% Authors: -% -% The authors of the code itself are Stefan I. Larimore and Timothy A. -% Davis (davis@cise.ufl.edu), University of Florida. The algorithm was -% developed in collaboration with John Gilbert, Xerox PARC, and Esmond -% Ng, Oak Ridge National Laboratory. -% -% Date: -% -% September 8, 2003. Version 2.3. -% -% Acknowledgements: -% -% This work was supported by the National Science Foundation, under -% grants DMS-9504974 and DMS-9803599. - -% Notice: -% Copyright (c) 1998-2003 by the University of Florida. -% All Rights Reserved. -% -% See http://www.cise.ufl.edu/research/sparse/colamd (the colamd.c -% file) for the License. -% -% Availability: -% -% The colamd/symamd library is available at -% -% http://www.cise.ufl.edu/research/sparse/colamd/ -% - -%------------------------------------------------------------------------------- -% perform the symamd ordering: -%------------------------------------------------------------------------------- - -if (nargout <= 1 & nargin == 1) - p = symamdmex (S) ; -elseif (nargout <= 1 & nargin == 2) - p = symamdmex (S, knobs) ; -elseif (nargout == 2 & nargin == 1) - [p, stats] = symamdmex (S) ; -elseif (nargout == 2 & nargin == 2) - [p, stats] = symamdmex (S, knobs) ; -else - error ('symamd: incorrect number of input and/or output arguments') ; -end - -%------------------------------------------------------------------------------- -% symmetric elimination tree post-ordering: -%------------------------------------------------------------------------------- - -[ignore, q] = sparsfun ('symetree', S (p,p)) ; -p = p (q) ; -
--- a/liboctave/COLAMD/symamdmex.c Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/symamdmex.c Sun Sep 04 12:25:21 2005 +0000 @@ -2,76 +2,22 @@ /* === symamd mexFunction =================================================== */ /* ========================================================================== */ -/* +/* COLAMD Version 2.4. + Usage: P = symamd (A) ; - - P = symamd (A, knobs) ; - - [ P, stats ] = symamd (A) ; - [ P, stats ] = symamd (A, knobs) ; - Returns a permutation vector P such that the Cholesky factorization of - A (P,P) tends to be sparser than that of A. This routine provides the same - functionality as SYMMMD, but tends to be much faster and tends to return a - better permutation vector. Note that the SYMMMD m-file in - MATLAB 5.2 also performs a symmetric elimination tree post-ordering. This - mexFunction does not do this post-ordering. This mexFunction is a - replacement for the p = sparsfun ('symmmd', A) operation. - - A must be square, and is assummed to have a symmetric nonzero pattern. - Only the nonzero pattern of the lower triangular portion of A is accessed. - This routine constructs a matrix M such that the nonzero pattern of M'M is - equal to A (assuming A has symmetric pattern), and then performs a column - ordering of M using colamd. - - The knobs and stats vectors are optional: - - knobs (1) rows and columns with more than (knobs(1))*n entries - are removed prior to ordering, and placed last in - the output ordering. If knobs is not present, then the - default of 0.5 is used. - - knobs (2) print level, similar to spparms ('spumoni') - - stats (1) the number of dense (or empty) rows and columms. These - are ordered last, in their natural order. - - stats (2) (same as stats (1)) - - stats (3) the number of garbage collections performed. - - stats (4) return status: - - 0: matrix is a valid MATLAB matrix. - - 1: matrix has duplicate entries or unsorted columns. - This should not occur in a valid MATLAB matrix, - but the ordering proceeded anyway by ignoring the - duplicate row indices in each column. See - stats (5:7) for more information. - - stats (5) highest numbered column that is unsorted or has - duplicate entries (zero if none) - - stats (6) last seen duplicate or unsorted row index - (zero if none) - - stats (7) number of duplicate or unsorted row indices + See symamd.m for a description. Authors: The authors of the code itself are Stefan I. Larimore and Timothy A. - Davis (davis@cise.ufl.edu), University of Florida. The algorithm was + Davis (davis at cise.ufl.edu), University of Florida. The algorithm was developed in collaboration with John Gilbert, Xerox PARC, and Esmond Ng, Oak Ridge National Laboratory. - Date: - - September 8, 2003. Version 2.3. - Acknowledgements: This work was supported by the National Science Foundation, under @@ -79,8 +25,7 @@ Notice: - Copyright (c) 1998-2003 by the University of Florida. - All Rights Reserved. + Copyright (c) 1998-2005, Timothy A. Davis. All Rights Reserved. See http://www.cise.ufl.edu/research/sparse/colamd (the colamd.c file) for the License. @@ -136,6 +81,8 @@ int spumoni ; /* verbosity variable */ int stats [COLAMD_STATS] ; /* stats for symamd */ + colamd_printf = mexPrintf ; /* COLAMD printf routine */ + /* === Check inputs ===================================================== */ if (nrhs < 1 || nrhs > 2 || nlhs < 0 || nlhs > 2) @@ -154,15 +101,27 @@ { in_knobs = mxGetPr (prhs [1]) ; i = mxGetNumberOfElements (prhs [1]) ; - if (i > 0) knobs [COLAMD_DENSE_ROW] = in_knobs [COLAMD_DENSE_ROW] ; - if (i > 1) spumoni = (int) in_knobs [1] ; + if (i > 0) knobs [COLAMD_DENSE_ROW] = in_knobs [0] ; + if (i > 1) spumoni = (int) (in_knobs [1] != 0) ; } - /* print knob settings if spumoni > 0 */ - if (spumoni > 0) + /* print knob settings if spumoni is set */ + if (spumoni) { - mexPrintf ("symamd: dense row/col fraction: %f\n", - knobs [COLAMD_DENSE_ROW]) ; + mexPrintf ("\nsymamd version %d.%d, %s:\n", + COLAMD_MAIN_VERSION, COLAMD_SUB_VERSION, COLAMD_DATE) ; + if (knobs [COLAMD_DENSE_ROW] >= 0) + { + mexPrintf ("knobs(1): %g, rows/cols with > " + "max(16,%g*sqrt(size(A,2))) entries removed\n", + in_knobs [0], knobs [COLAMD_DENSE_ROW]) ; + } + else + { + mexPrintf ("knobs(1): %g, no dense rows removed\n", in_knobs [0]) ; + } + mexPrintf ("knobs(2): %g, statistics and knobs printed\n", + in_knobs [1]) ; } /* === If A is full, convert to a sparse matrix ========================= */ @@ -218,8 +177,8 @@ /* === Return the stats vector ========================================== */ - /* print stats if spumoni > 0 */ - if (spumoni > 0) + /* print stats if spumoni is set */ + if (spumoni) { symamd_report (stats) ; }
--- a/liboctave/COLAMD/symamdtestmex.c Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/COLAMD/symamdtestmex.c Sun Sep 04 12:25:21 2005 +0000 @@ -2,7 +2,8 @@ /* === symamdtest mexFunction =============================================== */ /* ========================================================================== */ -/* +/* COLAMD Version 2.4. + This MATLAB mexFunction is for testing only. It is not meant for production use. See symamdmex.c instead. @@ -10,69 +11,21 @@ [ P, stats ] = symamdtest (A, knobs) ; - Returns a permutation vector P such that the Cholesky factorization of - A (P,P) tends to be sparser than that of A. This routine provides the same - functionality as SYMMMD, but tends to be much faster and tends to return a - better permutation vector. Note that the SYMMMD m-file in - MATLAB 5.2 also performs a symmetric elimination tree post-ordering. This - mexFunction does not do this post-ordering. This mexFunction is a - replacement for the p = sparsfun ('symmmd', A) operation. + See symamd.m for a description. knobs is required. - A must be square, and is assummed to have a symmetric nonzero pattern. - Only the nonzero pattern of the lower triangular portion of A is accessed. - This routine constructs a matrix M such that the nonzero pattern of M'M is - equal to A (assuming A has symmetric pattern), and then performs a column - ordering of M using colamd. - - The knobs and stats vectors are optional: - - knobs (1) rows and columns with more than (knobs(1))*n entries - are removed prior to ordering, and placed last in - the output ordering. If knobs is not present, then the - default of 0.5 is used. - - knobs (2) print level, similar to spparms ('spumoni') - + knobs (1) dense row control + knobs (2) spumoni knobs (3) for testing only. Controls how the input matrix is jumbled prior to calling symamd, to test its error handling capability. - stats (1) the number of dense (or empty) rows and columms. These - are ordered last, in their natural order. - - stats (2) (same as stats (1)) - - stats (3) the number of garbage collections performed. - - stats (4) return status: - - 0: matrix is a valid MATLAB matrix. - - 1: matrix has duplicate entries or unsorted columns. - This should not occur in a valid MATLAB matrix, - but the ordering proceeded anyway by ignoring the - duplicate row indices in each column. See - stats (5:7) for more information. - - stats (5) highest numbered column that is unsorted or has - duplicate entries (zero if none) - - stats (6) last seen duplicate or unsorted row index - (zero if none) - - stats (7) number of duplicate or unsorted row indices - Authors: The authors of the code itself are Stefan I. Larimore and Timothy A. - Davis (davis@cise.ufl.edu), University of Florida. The algorithm was + Davis (davis at cise.ufl.edu), University of Florida. The algorithm was developed in collaboration with John Gilbert, Xerox PARC, and Esmond Ng, Oak Ridge National Laboratory. - Date: - - September 8, 2003. Version 2.3. - Acknowledgements: This work was supported by the National Science Foundation, under @@ -80,8 +33,7 @@ Notice: - Copyright (c) 1998-2003 by the University of Florida. - All Rights Reserved. + Copyright (c) 1998-2005, Timothy A. Davis. All Rights Reserved. See http://www.cise.ufl.edu/research/sparse/colamd (the colamd.c file) for the License. @@ -153,6 +105,8 @@ int *stats ; stats = stats2 ; + colamd_printf = mexPrintf ; /* COLAMD printf routine */ + /* === Check inputs ===================================================== */ if (nrhs < 1 || nrhs > 2 || nlhs < 0 || nlhs > 2) @@ -166,7 +120,7 @@ mexErrMsgTxt ("symamdtest: knobs are required") ; } /* for testing we require all 3 knobs */ - if (mxGetNumberOfElements (prhs [1]) < 3) + if (mxGetNumberOfElements (prhs [1]) != 3) { mexErrMsgTxt ("symamdtest: must have all 3 knobs for testing") ; } @@ -181,15 +135,28 @@ { in_knobs = mxGetPr (prhs [1]) ; i = mxGetNumberOfElements (prhs [1]) ; - if (i > 0) knobs [COLAMD_DENSE_ROW] = in_knobs [COLAMD_DENSE_ROW] ; + if (i > 0) knobs [COLAMD_DENSE_ROW] = in_knobs [0] ; if (i > 1) spumoni = (int) in_knobs [1] ; } - /* print knob settings if spumoni > 0 */ - if (spumoni > 0) + /* print knob settings if spumoni is set */ + if (spumoni) { - mexPrintf ("symamd: dense row/col fraction: %f\n", - knobs [COLAMD_DENSE_ROW]) ; + mexPrintf ("\nsymamd version %d.%d, %s:\n", + COLAMD_MAIN_VERSION, COLAMD_SUB_VERSION, COLAMD_DATE) ; + if (knobs [COLAMD_DENSE_ROW] >= 0) + { + mexPrintf ("knobs(1): %g, rows/cols with > " + "max(16,%g*sqrt(size(A,2))) entries removed\n", + in_knobs [0], knobs [COLAMD_DENSE_ROW]) ; + } + else + { + mexPrintf ("knobs(1): %g, no dense rows removed\n", in_knobs [0]) ; + } + mexPrintf ("knobs(2): %g, statistics and knobs printed\n", + in_knobs [1]) ; + mexPrintf ("Testing %d\n", in_knobs [2]) ; } /* === If A is full, convert to a sparse matrix ========================= */
--- a/liboctave/ChangeLog Fri Sep 02 06:27:52 2005 +0000 +++ b/liboctave/ChangeLog Sun Sep 04 12:25:21 2005 +0000 @@ -1,3 +1,9 @@ +2005-09-04 David Bateman <dbateman@free.fr> + + * COLAMD: Update version of colamd to v2.4. + * COLAMD.files: Add colamd_global.c to COLAMD_SRC and second build of + colamd.c for long version. + 2005-08-25 David Bateman <dbateman@free.fr> * Sparse-op-defs.h (FULL_SPARSE_MUL, SPARSE_FULL_MUL): Macro for