Amesos Package Browser (Single Doxygen Collection) Development
amesos_cholmod_l_camd.c
Go to the documentation of this file.
00001 /* ========================================================================== */
00002 /* === Partition/cholmod_camd =============================================== */
00003 /* ========================================================================== */
00004 
00005 /* -----------------------------------------------------------------------------
00006  * CHOLMOD/Partition Module.  Copyright (C) 2005-2006, Timothy A. Davis
00007  * The CHOLMOD/Partition Module is licensed under Version 2.1 of the GNU
00008  * Lesser General Public License.  See lesser.txt for a text of the license.
00009  * CHOLMOD is also available under other licenses; contact authors for details.
00010  * http://www.cise.ufl.edu/research/sparse
00011  * -------------------------------------------------------------------------- */
00012 
00013 /* CHOLMOD interface to the CAMD ordering routine.  Orders A if the matrix is
00014  * symmetric.  On output, Perm [k] = i if row/column i of A is the kth
00015  * row/column of P*A*P'.  This corresponds to A(p,p) in MATLAB notation.
00016  *
00017  * If A is unsymmetric, cholmod_camd orders A*A'.  On output, Perm [k] = i if
00018  * row/column i of A*A' is the kth row/column of P*A*A'*P'.  This corresponds to
00019  * A(p,:)*A(p,:)' in MATLAB notation.  If f is present, A(p,f)*A(p,f)' is
00020  * ordered.
00021  *
00022  * Computes the flop count for a subsequent LL' factorization, the number
00023  * of nonzeros in L, and the number of nonzeros in the matrix ordered (A,
00024  * A*A' or A(:,f)*A(:,f)').
00025  *
00026  * workspace: Iwork (4*nrow). Head (nrow).
00027  *
00028  * Allocates a temporary copy of A+A' or A*A' (with
00029  * both upper and lower triangular parts) as input to CAMD.
00030  * Also allocates 3*(n+1) additional integer workspace (not in Common).
00031  *
00032  * Supports any xtype (pattern, real, complex, or zomplex)
00033  */
00034 
00035 #ifndef NPARTITION
00036 
00037 /* This file should make the long int version of CHOLMOD */
00038 #define DLONG 1
00039 
00040 #include "amesos_cholmod_internal.h"
00041 #include "amesos_camd.h"
00042 #include "amesos_cholmod_partition.h"
00043 
00044 #if (CAMD_VERSION < CAMD_VERSION_CODE (2,0))
00045 #error "CAMD v2.0 or later is required"
00046 #endif
00047 
00048 /* ========================================================================== */
00049 /* === cholmod_camd ========================================================= */
00050 /* ========================================================================== */
00051 
00052 int CHOLMOD(camd)
00053 (
00054     /* ---- input ---- */
00055     cholmod_sparse *A,  /* matrix to order */
00056     Int *fset,    /* subset of 0:(A->ncol)-1 */
00057     size_t fsize, /* size of fset */
00058     Int *Cmember, /* size nrow.  see cholmod_ccolamd.c for description.*/
00059     /* ---- output ---- */
00060     Int *Perm,    /* size A->nrow, output permutation */
00061     /* --------------- */
00062     cholmod_common *Common
00063 )
00064 {
00065     double Info [CAMD_INFO], Control2 [CAMD_CONTROL], *Control ;
00066     Int *Cp, *Len, *Nv, *Head, *Elen, *Degree, *Wi, *Next, *BucketSet,
00067   *Work3n, *p ;
00068     cholmod_sparse *C ;
00069     Int j, n, cnz ;
00070     size_t s ;
00071     int ok = TRUE ;
00072 
00073     /* ---------------------------------------------------------------------- */
00074     /* get inputs */
00075     /* ---------------------------------------------------------------------- */
00076 
00077     RETURN_IF_NULL_COMMON (FALSE) ;
00078     RETURN_IF_NULL (A, FALSE) ;
00079     n = A->nrow ;
00080 
00081     /* s = 4*n */
00082     s = CHOLMOD(mult_size_t) (n, 4, &ok) ;
00083     if (!ok)
00084     {
00085   ERROR (CHOLMOD_TOO_LARGE, "problem too large") ;
00086   return (FALSE) ;
00087     }
00088 
00089     RETURN_IF_NULL (Perm, FALSE) ;
00090     RETURN_IF_XTYPE_INVALID (A, CHOLMOD_PATTERN, CHOLMOD_ZOMPLEX, FALSE) ;
00091     Common->status = CHOLMOD_OK ;
00092     if (n == 0)
00093     {
00094   /* nothing to do */
00095   Common->fl = 0 ;
00096   Common->lnz = 0 ;
00097   Common->anz = 0 ;
00098   return (TRUE) ;
00099     }
00100 
00101     /* ---------------------------------------------------------------------- */
00102     /* get workspace */
00103     /* ---------------------------------------------------------------------- */
00104 
00105     /* cholmod_analyze has allocated Cmember at Common->Iwork + 5*n+uncol, and
00106      * CParent at Common->Iwork + 4*n+uncol, where uncol is 0 if A is symmetric
00107      * or A->ncol otherwise.  Thus, only the first 4n integers in Common->Iwork
00108      * can be used here. */
00109 
00110     CHOLMOD(allocate_work) (n, s, 0, Common) ;
00111     if (Common->status < CHOLMOD_OK)
00112     {
00113   return (FALSE) ;
00114     }
00115 
00116     p = Common->Iwork ;
00117     Degree = p ; p += n ; /* size n */
00118     Elen   = p ; p += n ; /* size n */
00119     Len    = p ; p += n ; /* size n */
00120     Nv     = p ; p += n ; /* size n */
00121 
00122     Work3n = CHOLMOD(malloc) (n+1, 3*sizeof (Int), Common) ;
00123     if (Common->status < CHOLMOD_OK)
00124     {
00125   return (FALSE) ;
00126     }
00127     p = Work3n ;
00128     Next = p ; p += n ;   /* size n */
00129     Wi   = p ; p += (n+1) ; /* size n+1 */
00130     BucketSet = p ;   /* size n */
00131 
00132     Head = Common->Head ; /* size n+1 */
00133 
00134     /* ---------------------------------------------------------------------- */
00135     /* construct the input matrix for CAMD */
00136     /* ---------------------------------------------------------------------- */
00137 
00138     if (A->stype == 0)
00139     {
00140   /* C = A*A' or A(:,f)*A(:,f)', add extra space of nnz(C)/2+n to C */
00141   C = CHOLMOD(aat) (A, fset, fsize, -2, Common) ;
00142     }
00143     else
00144     {
00145   /* C = A+A', but use only the upper triangular part of A if A->stype = 1
00146    * and only the lower part of A if A->stype = -1.  Add extra space of
00147    * nnz(C)/2+n to C. */
00148   C = CHOLMOD(copy) (A, 0, -2, Common) ;
00149     }
00150 
00151     if (Common->status < CHOLMOD_OK)
00152     {
00153   /* out of memory, fset invalid, or other error */
00154   CHOLMOD(free) (n+1, 3*sizeof (Int), Work3n, Common) ;
00155   return (FALSE) ;
00156     }
00157 
00158     Cp = C->p ;
00159     for (j = 0 ; j < n ; j++)
00160     {
00161   Len [j] = Cp [j+1] - Cp [j] ;
00162     }
00163 
00164     /* C does not include the diagonal, and both upper and lower parts.
00165      * Common->anz includes the diagonal, and just the lower part of C */
00166     cnz = Cp [n] ;
00167     Common->anz = cnz / 2 + n ;
00168 
00169     /* ---------------------------------------------------------------------- */
00170     /* order C using CAMD */
00171     /* ---------------------------------------------------------------------- */
00172 
00173     /* get parameters */
00174     if (Common->current < 0 || Common->current >= CHOLMOD_MAXMETHODS)
00175     {
00176   /* use CAMD defaults */
00177   Control = NULL ;
00178     }
00179     else
00180     {
00181   Control = Control2 ;
00182   Control [CAMD_DENSE] = Common->method [Common->current].prune_dense ;
00183   Control [CAMD_AGGRESSIVE] = Common->method [Common->current].aggressive;
00184     }
00185 
00186     /* CAMD_2 does not use camd_malloc and camd_free, but set these pointers
00187      * just be safe. */
00188     amesos_camd_malloc = Common->malloc_memory ;
00189     amesos_camd_free = Common->free_memory ;
00190     amesos_camd_calloc = Common->calloc_memory ;
00191     amesos_camd_realloc = Common->realloc_memory ;
00192 
00193     /* CAMD_2 doesn't print anything either, but future versions might,
00194      * so set the camd_printf pointer too. */
00195     amesos_camd_printf = Common->print_function ;
00196 
00197 #ifdef LONG
00198     /* DEBUG (camd_l_debug_init ("cholmod_l_camd")) ; */
00199     amesos_camd_l2 (n, C->p,  C->i, Len, C->nzmax, cnz, Nv, Next, Perm, Head, Elen,
00200       Degree, Wi, Control, Info, Cmember, BucketSet) ;
00201 #else
00202     /* DEBUG (camd_debug_init ("cholmod_camd")) ; */
00203     amesos_camd_2 (n, C->p,  C->i, Len, C->nzmax, cnz, Nv, Next, Perm, Head, Elen,
00204       Degree, Wi, Control, Info, Cmember, BucketSet) ;
00205 #endif
00206 
00207     /* LL' flop count.  Need to subtract n for LL' flop count.  Note that this
00208      * is a slight upper bound which is often exact (see CAMD/Source/camd_2.c
00209      * for details).  cholmod_analyze computes an exact flop count and
00210      * fill-in. */
00211     Common->fl = Info [CAMD_NDIV] + 2 * Info [CAMD_NMULTSUBS_LDL] + n ;
00212 
00213     /* Info [CAMD_LNZ] excludes the diagonal */
00214     Common->lnz = n + Info [CAMD_LNZ] ;
00215 
00216     /* ---------------------------------------------------------------------- */
00217     /* free the CAMD workspace and clear the persistent workspace in Common */
00218     /* ---------------------------------------------------------------------- */
00219 
00220     ASSERT (IMPLIES (Common->status == CHOLMOD_OK,
00221     CHOLMOD(dump_perm) (Perm, n, n, "CAMD2 perm", Common))) ;
00222     CHOLMOD(free_sparse) (&C, Common) ;
00223     for (j = 0 ; j <= n ; j++)
00224     {
00225   Head [j] = EMPTY ;
00226     }
00227     CHOLMOD(free) (n+1, 3*sizeof (Int), Work3n, Common) ;
00228     return (TRUE) ;
00229 }
00230 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines