227 lines
7.3 KiB
C
227 lines
7.3 KiB
C
/* ========================================================================== */
|
|
/* === ccolamd mexFunction ================================================== */
|
|
/* ========================================================================== */
|
|
|
|
/* ----------------------------------------------------------------------------
|
|
* CCOLAMD, Copyright (C), Univ. of Florida. Authors: Timothy A. Davis,
|
|
* Sivasankaran Rajamanickam, and Stefan Larimore
|
|
* -------------------------------------------------------------------------- */
|
|
|
|
/*
|
|
* Usage:
|
|
* p = ccolamd (A) ;
|
|
* [p stats] = ccolamd (A, knobs, cmember) ;
|
|
*
|
|
* See ccolamd.m for a description.
|
|
*/
|
|
|
|
/* ========================================================================== */
|
|
/* === Include files ======================================================== */
|
|
/* ========================================================================== */
|
|
|
|
#include "ccolamd.h"
|
|
#include "mex.h"
|
|
#include "matrix.h"
|
|
#include <stdlib.h>
|
|
#include <string.h>
|
|
#define Long SuiteSparse_long
|
|
|
|
/* ========================================================================== */
|
|
/* === ccolamd mexFunction ================================================== */
|
|
/* ========================================================================== */
|
|
|
|
void mexFunction
|
|
(
|
|
/* === Parameters ======================================================= */
|
|
|
|
int nargout,
|
|
mxArray *pargout [ ],
|
|
int nargin,
|
|
const mxArray *pargin [ ]
|
|
)
|
|
{
|
|
/* === Local variables ================================================== */
|
|
|
|
Long *A ; /* ccolamd's copy of the matrix and workspace */
|
|
Long *cmember ; /* ccolamd's copy of the constraint set */
|
|
double *in_cmember ; /* input constraint set */
|
|
Long *p ; /* ccolamd's copy of the column pointers */
|
|
Long Alen ; /* size of A */
|
|
Long cslen ; /* size of CS */
|
|
Long n_col ; /* number of columns of A */
|
|
Long n_row ; /* number of rows of A */
|
|
Long nnz ; /* number of entries in A */
|
|
Long full ; /* TRUE if input matrix full, FALSE if sparse */
|
|
double knobs [CCOLAMD_KNOBS] ; /* ccolamd user-controllable parameters */
|
|
double *out_perm ; /* output permutation vector */
|
|
double *out_stats ; /* output stats vector */
|
|
double *in_knobs ; /* input knobs vector */
|
|
Long i ; /* loop counter */
|
|
mxArray *Ainput ; /* input matrix handle */
|
|
Long spumoni ; /* verbosity variable */
|
|
Long stats [CCOLAMD_STATS] ;/* stats for ccolamd */
|
|
|
|
/* === Check inputs ===================================================== */
|
|
|
|
if (nargin < 1 || nargin > 3 || nargout < 0 || nargout > 2)
|
|
{
|
|
mexErrMsgTxt ("Usage: [p stats] = ccolamd (A, knobs, cmember)") ;
|
|
}
|
|
|
|
/* === Get cmember ====================================================== */
|
|
|
|
cmember = NULL ;
|
|
cslen = 0 ;
|
|
if (nargin > 2)
|
|
{
|
|
in_cmember = mxGetPr (pargin [2]) ;
|
|
cslen = mxGetNumberOfElements (pargin [2]) ;
|
|
if (cslen != 0)
|
|
{
|
|
cmember = (Long *) mxCalloc (cslen, sizeof (Long)) ;
|
|
for (i = 0 ; i < cslen ; i++)
|
|
{
|
|
/* convert cmember from 1-based to 0-based */
|
|
cmember[i] = ((Long) in_cmember [i] - 1) ;
|
|
}
|
|
}
|
|
}
|
|
|
|
/* === Get knobs ======================================================== */
|
|
|
|
ccolamd_l_set_defaults (knobs) ;
|
|
spumoni = 0 ;
|
|
|
|
/* check for user-passed knobs */
|
|
if (nargin > 1)
|
|
{
|
|
in_knobs = mxGetPr (pargin [1]) ;
|
|
i = mxGetNumberOfElements (pargin [1]) ;
|
|
if (i > 0) knobs [CCOLAMD_LU] = (in_knobs [0] != 0) ;
|
|
if (i > 1) knobs [CCOLAMD_DENSE_ROW] = in_knobs [1] ;
|
|
if (i > 2) knobs [CCOLAMD_DENSE_COL] = in_knobs [2] ;
|
|
if (i > 3) knobs [CCOLAMD_AGGRESSIVE] = (in_knobs [3] != 0) ;
|
|
if (i > 4) spumoni = (in_knobs [4] != 0) ;
|
|
}
|
|
|
|
/* print knob settings if spumoni is set */
|
|
if (spumoni)
|
|
{
|
|
mexPrintf ("\nccolamd version %d.%d, %s:\nknobs(1): %g, order for %s\n",
|
|
CCOLAMD_MAIN_VERSION, CCOLAMD_SUB_VERSION, CCOLAMD_DATE,
|
|
in_knobs [0],
|
|
(knobs [CCOLAMD_LU] != 0) ? "lu(A)" : "chol(A'*A)") ;
|
|
if (knobs [CCOLAMD_DENSE_ROW] >= 0)
|
|
{
|
|
mexPrintf ("knobs(2): %g, rows with > max(16,%g*sqrt(size(A,2)))"
|
|
" entries removed\n", in_knobs [1], knobs [CCOLAMD_DENSE_ROW]) ;
|
|
}
|
|
else
|
|
{
|
|
mexPrintf ("knobs(2): %g, no dense rows removed\n", in_knobs [1]) ;
|
|
}
|
|
if (knobs [CCOLAMD_DENSE_COL] >= 0)
|
|
{
|
|
mexPrintf ("knobs(3): %g, cols with > max(16,%g*sqrt(min(size(A)))"
|
|
" entries removed\n", in_knobs [2], knobs [CCOLAMD_DENSE_COL]) ;
|
|
}
|
|
else
|
|
{
|
|
mexPrintf ("knobs(3): no dense columns removed\n", in_knobs [2]) ;
|
|
}
|
|
mexPrintf ("knobs(4): %g, aggressive absorption: %s\n",
|
|
in_knobs [3], (knobs [CCOLAMD_AGGRESSIVE] != 0) ? "yes" : "no") ;
|
|
mexPrintf ("knobs(5): %g, statistics and knobs printed\n",
|
|
in_knobs [4]) ;
|
|
}
|
|
|
|
/* === If A is full, convert to a sparse matrix ========================= */
|
|
|
|
Ainput = (mxArray *) pargin [0] ;
|
|
if (mxGetNumberOfDimensions (Ainput) != 2)
|
|
{
|
|
mexErrMsgTxt ("ccolamd: input matrix must be 2-dimensional") ;
|
|
}
|
|
full = !mxIsSparse (Ainput) ;
|
|
if (full)
|
|
{
|
|
mexCallMATLAB (1, &Ainput, 1, (mxArray **) pargin, "sparse") ;
|
|
}
|
|
|
|
/* === Allocate workspace for ccolamd =================================== */
|
|
|
|
/* get size of matrix */
|
|
n_row = mxGetM (Ainput) ;
|
|
n_col = mxGetN (Ainput) ;
|
|
|
|
/* get column pointer vector */
|
|
p = (Long *) mxCalloc (n_col+1, sizeof (Long)) ;
|
|
(void) memcpy (p, mxGetJc (Ainput), (n_col+1)*sizeof (Long)) ;
|
|
nnz = p [n_col] ;
|
|
Alen = (Long) ccolamd_l_recommended (nnz, n_row, n_col) ;
|
|
if (Alen == 0)
|
|
{
|
|
mexErrMsgTxt ("ccolamd: problem too large") ;
|
|
}
|
|
|
|
/* === Copy input matrix into workspace ================================= */
|
|
|
|
A = (Long *) mxCalloc (Alen, sizeof (Long)) ;
|
|
(void) memcpy (A, mxGetIr (Ainput), nnz*sizeof (Long)) ;
|
|
|
|
if (full)
|
|
{
|
|
mxDestroyArray (Ainput) ;
|
|
}
|
|
|
|
/* Check constraint set size */
|
|
if (cmember != NULL && cslen != n_col)
|
|
{
|
|
mexErrMsgTxt ("ccolamd: cmember must be of length equal to #cols of A");
|
|
}
|
|
|
|
/* === Order the columns (destroys A) =================================== */
|
|
|
|
if (!ccolamd_l (n_row, n_col, Alen, A, p, knobs, stats, cmember))
|
|
{
|
|
ccolamd_l_report (stats) ;
|
|
mexErrMsgTxt ("ccolamd error!") ;
|
|
}
|
|
mxFree (A) ;
|
|
mxFree (cmember) ;
|
|
|
|
/* === Return the permutation vector ==================================== */
|
|
|
|
pargout [0] = mxCreateDoubleMatrix (1, n_col, mxREAL) ;
|
|
out_perm = mxGetPr (pargout [0]) ;
|
|
for (i = 0 ; i < n_col ; i++)
|
|
{
|
|
/* ccolamd is 0-based, but MATLAB expects this to be 1-based */
|
|
out_perm [i] = p [i] + 1 ;
|
|
}
|
|
mxFree (p) ;
|
|
|
|
/* === Return the stats vector ========================================== */
|
|
|
|
/* print stats if spumoni is set */
|
|
if (spumoni)
|
|
{
|
|
ccolamd_l_report (stats) ;
|
|
}
|
|
|
|
if (nargout == 2)
|
|
{
|
|
pargout [1] = mxCreateDoubleMatrix (1, CCOLAMD_STATS, mxREAL) ;
|
|
out_stats = mxGetPr (pargout [1]) ;
|
|
for (i = 0 ; i < CCOLAMD_STATS ; i++)
|
|
{
|
|
out_stats [i] = stats [i] ;
|
|
}
|
|
|
|
/* fix stats (5) and (6), for 1-based information on jumbled matrix. */
|
|
/* note that this correction doesn't occur if symamd returns FALSE */
|
|
out_stats [CCOLAMD_INFO1] ++ ;
|
|
out_stats [CCOLAMD_INFO2] ++ ;
|
|
}
|
|
}
|