~ubuntu-branches/ubuntu/natty/suitesparse/natty

« back to all changes in this revision

Viewing changes to CXSparse_newfiles/MATLAB/CSparse/cs_thumb_mex.c

  • Committer: Bazaar Package Importer
  • Author(s): Christophe Prud'homme
  • Date: 2007-05-29 09:36:29 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20070529093629-zowquo0b7slkk6nc
Tags: 3.0.0-2
* suitesparse builds properly twice in a row
* Bug fix: "suitesparse - FTBFS: Broken build depens: libgfortran1-dev",
  thanks to Bastian Blank (Closes: #426349).
* Bug fix: "suitesparse_3.0.0-1: FTBFS: build-depends on
  libgfortran1-dev", thanks to Steve Langasek (Closes: #426354).

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
#include "cs_mex.h"
 
2
/* cs_thumb: convert a sparse matrix to a dense 2D thumbnail matrix of size
 
3
 * at most k-by-k.  k defaults to 256.  A helper mexFunction for cspy. */
 
4
 
 
5
#define INDEX(i,j,lda) ((i)+(j)*(lda))
 
6
#define ISNAN(x) ((x) != (x))
 
7
#ifdef DBL_MAX
 
8
#define BIG_VALUE DBL_MAX
 
9
#else
 
10
#define BIG_VALUE 1.7e308
 
11
#endif
 
12
 
 
13
void mexFunction
 
14
(
 
15
    int nargout,
 
16
    mxArray *pargout [ ],
 
17
    int nargin,
 
18
    const mxArray *pargin [ ]
 
19
)
 
20
{
 
21
    CS_INT m, n, mn, m2, n2, k, s, j, ij, sj, si, p, *Ap, *Ai ;
 
22
    double aij, ax, az, *S, *Ax, *Az ;
 
23
    if (nargout > 1 || nargin < 1 || nargin > 2)
 
24
    {
 
25
        mexErrMsgTxt ("Usage: S = cs_thumb(A,k)") ;
 
26
    }
 
27
    cs_mex_check (0, -1, -1, 0, 1, 1, pargin [0]) ;
 
28
    m = mxGetM (pargin [0]) ;
 
29
    n = mxGetN (pargin [0]) ;
 
30
    mn = CS_MAX (m,n) ;
 
31
    k = (nargin == 1) ? 256 : mxGetScalar (pargin [1]) ;    /* get k */
 
32
    /* s = size of each submatrix; A(1:s,1:s) maps to S(1,1) */
 
33
    s = (mn < k) ? 1 : (CS_INT) ceil ((double) mn / (double) k) ;
 
34
    m2 = (CS_INT) ceil ((double) m / (double) s) ;
 
35
    n2 = (CS_INT) ceil ((double) n / (double) s) ;
 
36
    /* create S */
 
37
    pargout [0] = mxCreateDoubleMatrix (m2, n2, mxREAL) ;
 
38
    S = mxGetPr (pargout [0]) ;
 
39
    Ap = (CS_INT *) mxGetJc (pargin [0]) ;
 
40
    Ai = (CS_INT *) mxGetIr (pargin [0]) ;
 
41
    Ax = mxGetPr (pargin [0]) ;
 
42
    Az = (mxIsComplex (pargin [0])) ? mxGetPi (pargin [0]) : NULL ;
 
43
    for (j = 0 ; j < n ; j++)
 
44
    {
 
45
        sj = j/s ;
 
46
        for (p = Ap [j] ; p < Ap [j+1] ; p++)
 
47
        {
 
48
            si = Ai [p] / s ;
 
49
            ij = INDEX (si,sj,m2) ;
 
50
            ax = Ax [p] ;
 
51
            az = Az ? Az [p] : 0 ;
 
52
            if (az == 0)
 
53
            {
 
54
                aij = fabs (ax) ;
 
55
            }
 
56
            else
 
57
            {
 
58
                aij = sqrt (ax*ax + az*az) ;
 
59
            }
 
60
            if (ISNAN (aij)) aij = BIG_VALUE ;
 
61
            aij = CS_MIN (BIG_VALUE, aij) ;
 
62
            S [ij] = CS_MAX (S [ij], aij) ;
 
63
        }
 
64
    }
 
65
}