1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/src/amd/amd_preprocess.c	Mon Dec 06 13:09:21 2010 +0100
     1.3 @@ -0,0 +1,119 @@
     1.4 +/* ========================================================================= */
     1.5 +/* === AMD_preprocess ====================================================== */
     1.6 +/* ========================================================================= */
     1.7 +
     1.8 +/* ------------------------------------------------------------------------- */
     1.9 +/* AMD, Copyright (c) Timothy A. Davis,                                      */
    1.10 +/* Patrick R. Amestoy, and Iain S. Duff.  See ../README.txt for License.     */
    1.11 +/* email: davis at cise.ufl.edu    CISE Department, Univ. of Florida.        */
    1.12 +/* web: http://www.cise.ufl.edu/research/sparse/amd                          */
    1.13 +/* ------------------------------------------------------------------------- */
    1.14 +
    1.15 +/* Sorts, removes duplicate entries, and transposes from the nonzero pattern of
    1.16 + * a column-form matrix A, to obtain the matrix R.  The input matrix can have
    1.17 + * duplicate entries and/or unsorted columns (AMD_valid (n,Ap,Ai) must not be
    1.18 + * AMD_INVALID).
    1.19 + *
    1.20 + * This input condition is NOT checked.  This routine is not user-callable.
    1.21 + */
    1.22 +
    1.23 +#include "amd_internal.h"
    1.24 +
    1.25 +/* ========================================================================= */
    1.26 +/* === AMD_preprocess ====================================================== */
    1.27 +/* ========================================================================= */
    1.28 +
    1.29 +/* AMD_preprocess does not check its input for errors or allocate workspace.
    1.30 + * On input, the condition (AMD_valid (n,n,Ap,Ai) != AMD_INVALID) must hold.
    1.31 + */
    1.32 +
    1.33 +GLOBAL void AMD_preprocess
    1.34 +(
    1.35 +    Int n,              /* input matrix: A is n-by-n */
    1.36 +    const Int Ap [ ],   /* size n+1 */
    1.37 +    const Int Ai [ ],   /* size nz = Ap [n] */
    1.38 +
    1.39 +    /* output matrix R: */
    1.40 +    Int Rp [ ],         /* size n+1 */
    1.41 +    Int Ri [ ],         /* size nz (or less, if duplicates present) */
    1.42 +
    1.43 +    Int W [ ],          /* workspace of size n */
    1.44 +    Int Flag [ ]        /* workspace of size n */
    1.45 +)
    1.46 +{
    1.47 +
    1.48 +    /* --------------------------------------------------------------------- */
    1.49 +    /* local variables */
    1.50 +    /* --------------------------------------------------------------------- */
    1.51 +
    1.52 +    Int i, j, p, p2 ;
    1.53 +
    1.54 +    ASSERT (AMD_valid (n, n, Ap, Ai) != AMD_INVALID) ;
    1.55 +
    1.56 +    /* --------------------------------------------------------------------- */
    1.57 +    /* count the entries in each row of A (excluding duplicates) */
    1.58 +    /* --------------------------------------------------------------------- */
    1.59 +
    1.60 +    for (i = 0 ; i < n ; i++)
    1.61 +    {
    1.62 +        W [i] = 0 ;             /* # of nonzeros in row i (excl duplicates) */
    1.63 +        Flag [i] = EMPTY ;      /* Flag [i] = j if i appears in column j */
    1.64 +    }
    1.65 +    for (j = 0 ; j < n ; j++)
    1.66 +    {
    1.67 +        p2 = Ap [j+1] ;
    1.68 +        for (p = Ap [j] ; p < p2 ; p++)
    1.69 +        {
    1.70 +            i = Ai [p] ;
    1.71 +            if (Flag [i] != j)
    1.72 +            {
    1.73 +                /* row index i has not yet appeared in column j */
    1.74 +                W [i]++ ;           /* one more entry in row i */
    1.75 +                Flag [i] = j ;      /* flag row index i as appearing in col j*/
    1.76 +            }
    1.77 +        }
    1.78 +    }
    1.79 +
    1.80 +    /* --------------------------------------------------------------------- */
    1.81 +    /* compute the row pointers for R */
    1.82 +    /* --------------------------------------------------------------------- */
    1.83 +
    1.84 +    Rp [0] = 0 ;
    1.85 +    for (i = 0 ; i < n ; i++)
    1.86 +    {
    1.87 +        Rp [i+1] = Rp [i] + W [i] ;
    1.88 +    }
    1.89 +    for (i = 0 ; i < n ; i++)
    1.90 +    {
    1.91 +        W [i] = Rp [i] ;
    1.92 +        Flag [i] = EMPTY ;
    1.93 +    }
    1.94 +
    1.95 +    /* --------------------------------------------------------------------- */
    1.96 +    /* construct the row form matrix R */
    1.97 +    /* --------------------------------------------------------------------- */
    1.98 +
    1.99 +    /* R = row form of pattern of A */
   1.100 +    for (j = 0 ; j < n ; j++)
   1.101 +    {
   1.102 +        p2 = Ap [j+1] ;
   1.103 +        for (p = Ap [j] ; p < p2 ; p++)
   1.104 +        {
   1.105 +            i = Ai [p] ;
   1.106 +            if (Flag [i] != j)
   1.107 +            {
   1.108 +                /* row index i has not yet appeared in column j */
   1.109 +                Ri [W [i]++] = j ;  /* put col j in row i */
   1.110 +                Flag [i] = j ;      /* flag row index i as appearing in col j*/
   1.111 +            }
   1.112 +        }
   1.113 +    }
   1.114 +
   1.115 +#ifndef NDEBUG
   1.116 +    ASSERT (AMD_valid (n, n, Rp, Ri) == AMD_OK) ;
   1.117 +    for (j = 0 ; j < n ; j++)
   1.118 +    {
   1.119 +        ASSERT (W [j] == Rp [j+1]) ;
   1.120 +    }
   1.121 +#endif
   1.122 +}