Thread Subject: Convert array of class 'int8' to sparse

Hi,
Am new to sparse matrices and hoping someone can help me with this. I have a large (150MB) binary file of 'int8' type data which I am loading into an array A. The data has a large number of 0's, so I would like to hold it in a sparse matrix. However, if I type the command B=sparse(A), I get the error message: "??? Undefined function or method 'sparse' for input arguments of type 'int8'". I have read through the documentation on sparse matrices and can't see anywhere where it says one can only convert matrices of type 'double' to sparse (unless I missed it). In any case, is there any reason why it would not be possible to convert integer type arrays to sparse? Or am I doing something wrong?

Thanks in advance for any suggestions.

Bobane

"Bobane" <bobane73@hotmail.com> wrote in message <gfhh6c$dvr$1@fred.mathworks.com>...
> Hi,
> Am new to sparse matrices and hoping someone can help me with this. I have a large (150MB) binary file of 'int8' type data which I am loading into an array A. The data has a large number of 0's, so I would like to hold it in a sparse matrix. However, if I type the command B=sparse(A), I get the error message: "??? Undefined function or method 'sparse' for input arguments of type 'int8'". I have read through the documentation on sparse matrices and can't see anywhere where it says one can only convert matrices of type 'double' to sparse (unless I missed it). In any case, is there any reason why it would not be possible to convert integer type arrays to sparse? Or am I doing something wrong?
>
> Thanks in advance for any suggestions.
>
> Bobane

The sparse data type does indeed only support doubles, much to the chagrin of many users. So all you can do is use double() to pre-convert your data.

I too am looking for this feature to come to matlab's sparse matrix tools.

"Gus Lott" <lottg.nospam@janelia.hhmi.org> wrote in message <gq9bg5$3ra$1@fred.mathworks.com>...
> I too am looking for this feature to come to matlab's sparse matrix tools.

All,

You could always just modify the MATLAB mex example fulltosparse.c of filling a sparse matrix and change the input from double to int8. The Mathworks comments seem to imply that this is the very purpose of the example code. E.g., the code is below. To mex it just do this:

>> mex -setup
    (then select a C/C++ compiler, such as lcc)
>> mex int8tosparse.c

Now you have a function int8tosparse that takes a 2-dimensional int8 input and directly converts it into a sparse double matrix without taking the intermediate step of converting to a full double matrix first. If there is some interest in this, I could work this up into a complete mex function that works for *all* integer classes and post it to the FEX. Let me know.

James Tursa

Caveat, I only tested this for a single 100x100 example. It got the same answer as converting to a double matrix first, so I would assume the code is correct. However, be advised this is the Mathworks example code, not mine. I just modified it for int8. If you have a later version of MATLAB you might want to change the int j,k,etc to mwSize j,k,etc.

-----------------------------------------------------------------------------------------------------

/*
 * =============================================================
 * int8tosparse.c
 * This example demonstrates how to populate a sparse
 * matrix. For the purpose of this example, you must pass in a
 * non-sparse 2-dimensional argument of type int8.
 *
 * Mathworks fulltosparse.c modified by J. Tursa for int8.
 *
=============================================================
 */

#include <math.h> /* Needed for the ceil() prototype. */
#include "mex.h"

void mexFunction(
        int nlhs, mxArray *plhs[],
        int nrhs, const mxArray *prhs[]
        )
{
  /* Declare variables. */
  int j,k,m,n,nzmax,*irs,*jcs,cmplx,isfull;
  double *si,*sr;
  double percent_sparse;
  signed char *pr, *pi;

  /* Check for proper number of input and output arguments. */
  if (nrhs != 1) {
    mexErrMsgTxt("One input argument required.");
  }
  if (nlhs > 1) {
    mexErrMsgTxt("Too many output arguments.");
  }

  /* Check data type of input argument. */
  if (!(mxIsInt8(prhs[0]))) {
    mexErrMsgTxt("Input argument must be of type int8.");
  }

  if (mxGetNumberOfDimensions(prhs[0]) != 2) {
    mexErrMsgTxt("Input argument must be two dimensional\n");
  }

  /* Get the size and pointers to input data. */
  m = mxGetM(prhs[0]);
  n = mxGetN(prhs[0]);
  pr = mxGetPr(prhs[0]);
  pi = mxGetPi(prhs[0]);
  cmplx = (pi == NULL ? 0 : 1);

  /* Allocate space for sparse matrix.
   * NOTE: Assume at most 20% of the data is sparse. Use ceil
   * to cause it to round up.
   */

  percent_sparse = 0.2;
  nzmax = (int)ceil((double)m*(double)n*percent_sparse);

  plhs[0] = mxCreateSparse(m,n,nzmax,cmplx);
  sr = mxGetPr(plhs[0]);
  si = mxGetPi(plhs[0]);
  irs = mxGetIr(plhs[0]);
  jcs = mxGetJc(plhs[0]);
    
  /* Copy nonzeros. */
  k = 0;
  isfull = 0;
  for (j = 0; (j < n); j++) {
    int i;
    jcs[j] = k;
    for (i = 0; (i < m); i++) {
      if (pr[i] || (cmplx && pi[i])) {

        /* Check to see if non-zero element will fit in
         * allocated output array. If not, increase
         * percent_sparse by 10%, recalculate nzmax, and augment
         * the sparse array.
         */
        if (k >= nzmax) {
          int oldnzmax = nzmax;
          percent_sparse += 0.1;
          nzmax = (int)ceil((double)m*(double)n*percent_sparse);

          /* Make sure nzmax increases atleast by 1. */
          if (oldnzmax == nzmax)
            nzmax++;

          mxSetNzmax(plhs[0], nzmax);
          mxSetPr(plhs[0], mxRealloc(sr, nzmax*sizeof(double)));
          if (si != NULL)
          mxSetPi(plhs[0], mxRealloc(si, nzmax*sizeof(double)));
          mxSetIr(plhs[0], mxRealloc(irs, nzmax*sizeof(int)));

          sr = mxGetPr(plhs[0]);
          si = mxGetPi(plhs[0]);
          irs = mxGetIr(plhs[0]);
        }
        sr[k] = pr[i];
        if (cmplx) {
          si[k] = pi[i];
        }
        irs[k] = i;
        k++;
      }
    }
    pr += m;
    pi += m;
  }
  jcs[n] = k;
}

"James Tursa" <aclassyguywithaknotac@hotmail.com> wrote in message <gqa8th$7j0$1@fred.mathworks.com>...
> "Gus Lott" <lottg.nospam@janelia.hhmi.org> wrote in message <gq9bg5$3ra$1@fred.mathworks.com>...
> > I too am looking for this feature to come to matlab's sparse matrix tools.
>
> All,
>
> You could always just modify the MATLAB mex example fulltosparse.c of filling a sparse matrix and change the input from double to int8. The Mathworks comments seem to imply that this is the very purpose of the example code. E.g., the code is below. To mex it just do this:
>
> >> mex -setup
> (then select a C/C++ compiler, such as lcc)
> >> mex int8tosparse.c
>
> Now you have a function int8tosparse that takes a 2-dimensional int8 input and directly converts it into a sparse double matrix without taking the intermediate step of converting to a full double matrix first. If there is some interest in this, I could work this up into a complete mex function that works for *all* integer classes and post it to the FEX. Let me know.
>
> James Tursa
>
> Caveat, I only tested this for a single 100x100 example. It got the same answer as converting to a double matrix first, so I would assume the code is correct. However, be advised this is the Mathworks example code, not mine. I just modified it for int8. If you have a later version of MATLAB you might want to change the int j,k,etc to mwSize j,k,etc.
>
> -----------------------------------------------------------------------------------------------------
>
> /*
> * =============================================================
> * int8tosparse.c
> * This example demonstrates how to populate a sparse
> * matrix. For the purpose of this example, you must pass in a
> * non-sparse 2-dimensional argument of type int8.
> *
> * Mathworks fulltosparse.c modified by J. Tursa for int8.
> *
> =============================================================
> */
>
> #include <math.h> /* Needed for the ceil() prototype. */
> #include "mex.h"
>
> void mexFunction(
> int nlhs, mxArray *plhs[],
> int nrhs, const mxArray *prhs[]
> )
> {
> /* Declare variables. */
> int j,k,m,n,nzmax,*irs,*jcs,cmplx,isfull;
> double *si,*sr;
> double percent_sparse;
> signed char *pr, *pi;
>
> /* Check for proper number of input and output arguments. */
> if (nrhs != 1) {
> mexErrMsgTxt("One input argument required.");
> }
> if (nlhs > 1) {
> mexErrMsgTxt("Too many output arguments.");
> }
>
> /* Check data type of input argument. */
> if (!(mxIsInt8(prhs[0]))) {
> mexErrMsgTxt("Input argument must be of type int8.");
> }
>
> if (mxGetNumberOfDimensions(prhs[0]) != 2) {
> mexErrMsgTxt("Input argument must be two dimensional\n");
> }
>
> /* Get the size and pointers to input data. */
> m = mxGetM(prhs[0]);
> n = mxGetN(prhs[0]);
> pr = mxGetPr(prhs[0]);
> pi = mxGetPi(prhs[0]);
> cmplx = (pi == NULL ? 0 : 1);
>
> /* Allocate space for sparse matrix.
> * NOTE: Assume at most 20% of the data is sparse. Use ceil
> * to cause it to round up.
> */
>
> percent_sparse = 0.2;
> nzmax = (int)ceil((double)m*(double)n*percent_sparse);
>
> plhs[0] = mxCreateSparse(m,n,nzmax,cmplx);
> sr = mxGetPr(plhs[0]);
> si = mxGetPi(plhs[0]);
> irs = mxGetIr(plhs[0]);
> jcs = mxGetJc(plhs[0]);
>
> /* Copy nonzeros. */
> k = 0;
> isfull = 0;
> for (j = 0; (j < n); j++) {
> int i;
> jcs[j] = k;
> for (i = 0; (i < m); i++) {
> if (pr[i] || (cmplx && pi[i])) {
>
> /* Check to see if non-zero element will fit in
> * allocated output array. If not, increase
> * percent_sparse by 10%, recalculate nzmax, and augment
> * the sparse array.
> */
> if (k >= nzmax) {
> int oldnzmax = nzmax;
> percent_sparse += 0.1;
> nzmax = (int)ceil((double)m*(double)n*percent_sparse);
>
> /* Make sure nzmax increases atleast by 1. */
> if (oldnzmax == nzmax)
> nzmax++;
>
> mxSetNzmax(plhs[0], nzmax);
> mxSetPr(plhs[0], mxRealloc(sr, nzmax*sizeof(double)));
> if (si != NULL)
> mxSetPi(plhs[0], mxRealloc(si, nzmax*sizeof(double)));
> mxSetIr(plhs[0], mxRealloc(irs, nzmax*sizeof(int)));
>
> sr = mxGetPr(plhs[0]);
> si = mxGetPi(plhs[0]);
> irs = mxGetIr(plhs[0]);
> }
> sr[k] = pr[i];
> if (cmplx) {
> si[k] = pi[i];
> }
> irs[k] = i;
> k++;
> }
> }
> pr += m;
> pi += m;
> }
> jcs[n] = k;
> }



Dear John,

With
----------------------------------------------------
>> ver
MATLAB Version 7.6.0.324 (R2008a)
Operating System: Microsoft Windows Vista Version 6.0
 (Build 6001: Service Pack 1)
Java VM Version: Java 1.6.0 with Sun Microsystems Inc.
Java HotSpot(TM) 64-Bit Server VM mixed mode
-----------------------------------------------------

I ran

-----------------------------------------------------
>> mex -largeArrayDims int8tosparse.c
Microsoft (R) C/C++ Optimizing Compiler Version
15.00.21022.08 for x64 Copyright (C) Microsoft
Corporation. All rights reserved.
-----------------------------------------------------

which generated the file int8tosparse.mexw64

I then ran

--------------------------------------------------------
>> n = 1000; Sp16 = sprand(n,n,0.001); D16 = full(Sp16);...
               D8 = int8(D16); Sp8 = int8tosparse(D8);...
               whos

  Name Size Bytes Class Attributes


  Sp16 1000x1000 24,008 double sparse
  D16 1000x1000 8,000,000 double
  D8 1000x1000 1,000,000 int8
  Sp8 1000x1000 3,208,008 double sparse
  n 1x1 8 double

commas added so I can read the sizes
-----------------------------------------------------------

I know nothing about C/C++ programming and compiling. Any suggestions?

A complete mex function that works for *all* integer classes would be very useful, especially when working with large sparse graphs.

Regards,

Derek O'Connor

"Derek O'Connor" <derekroconnor@eircom.net> wrote in message <gqae67$vg$1@fred.mathworks.com>...
>
> I then ran
>
> --------------------------------------------------------
> >> n = 1000; Sp16 = sprand(n,n,0.001); D16 = full(Sp16);...
> D8 = int8(D16); Sp8 = int8tosparse(D8);...
> whos
>
> Name Size Bytes Class Attributes
>
>
> Sp16 1000x1000 24,008 double sparse
> D16 1000x1000 8,000,000 double
> D8 1000x1000 1,000,000 int8
> Sp8 1000x1000 3,208,008 double sparse
> n 1x1 8 double
>
> commas added so I can read the sizes
> -----------------------------------------------------------
>
> I know nothing about C/C++ programming and compiling. Any suggestions?
>
> A complete mex function that works for *all* integer classes would be very useful, especially when working with large sparse graphs.
>
> Regards,
>
> Derek O'Connor

First, you are comparing apples to oranges. You need to make sure the original sparse matrix can be converted to int8 and still get the same values before any comparison will be meaningful. i.e., make sure the original values are all integral between -128 and +127. Second, you cannot just compare storage size to see if the conversion happened correctly, you need to use the isequal function to determine this. Third, the code guesses that the storage needed is 20% and then incrementally bumps up the memory requirement by 10% if needed. This may be too large of a guess for your case, and would need to be made smarter for general cases. Finally, for your later version of MATLAB I would change int to mwSize as I mentioned earlier. So change the c code lines:

int j,k,m,n,nzmax,*irs,*jcs,cmplx,isfull;

to

mwSize j,k,m,n,nzmax,*irs,*jcs,cmplx,isfull;

and

  percent_sparse = 0.2;

to

  percent_sparse = 0.001;

and

          percent_sparse += 0.1;

to

          percent_sparse += 0.001;

And change your example to this:

n = 1000;
Sp16 = sprand(n,n,0.10);
Sp16(Sp16>0) = 1./Sp16(Sp16>0);
Sp16(Sp16>127) = 127;
Sp16 = floor(Sp16);
D16 = full(Sp16);
D8 = int8(D16);
Sp8 = int8tosparse(D8);
whos
isequal(D16,D8)
isequal(Sp16,Sp8)


The percent_sparse changes are *only* for your particular example and are just a quick fix to get it working for this particular case. I will work on some more generic code later to handle arbitrary size matrices.

James Tursa

"James Tursa" <aclassyguywithaknotac@hotmail.com> wrote in message <gqaqc5$l2k$1@fred.mathworks.com>...
>
> And change your example to this:
>
> n = 1000;
> Sp16 = sprand(n,n,0.10);
> Sp16(Sp16>0) = 1./Sp16(Sp16>0);
> Sp16(Sp16>127) = 127;
> Sp16 = floor(Sp16);
> D16 = full(Sp16);
> D8 = int8(D16);
> Sp8 = int8tosparse(D8);
> whos
> isequal(D16,D8)
> isequal(Sp16,Sp8)
>

Oops .. typo holdover from my small test. Try this instead:

n = 1000;
Sp16 = sprand(n,n,0.001); % <-- changed 0.10 to 0.001
Sp16(Sp16>0) = 1./Sp16(Sp16>0);
Sp16(Sp16>127) = 127;
Sp16 = floor(Sp16);
D16 = full(Sp16);
D8 = int8(D16);
Sp8 = int8tosparse(D8);
whos
isequal(D16,D8)
isequal(Sp16,Sp8)

James Tursa

"James Tursa" <aclassyguywithaknotac@hotmail.com> wrote in message <gqarna$pv1$1@fred.mathworks.com>...
> "James Tursa" <aclassyguywithaknotac@hotmail.com> wrote in message <gqaqc5$l2k$1@fred.mathworks.com>...
> >
> > And change your example to this:
> >
> > n = 1000;
> > Sp16 = sprand(n,n,0.10);
> > Sp16(Sp16>0) = 1./Sp16(Sp16>0);
> > Sp16(Sp16>127) = 127;
> > Sp16 = floor(Sp16);
> > D16 = full(Sp16);
> > D8 = int8(D16);
> > Sp8 = int8tosparse(D8);
> > whos
> > isequal(D16,D8)
> > isequal(Sp16,Sp8)
> >
>
> Oops .. typo holdover from my small test. Try this instead:
>
> n = 1000;
> Sp16 = sprand(n,n,0.001); % <-- changed 0.10 to 0.001
> Sp16(Sp16>0) = 1./Sp16(Sp16>0);
> Sp16(Sp16>127) = 127;
> Sp16 = floor(Sp16);
> D16 = full(Sp16);
> D8 = int8(D16);
> Sp8 = int8tosparse(D8);
> whos
> isequal(D16,D8)
> isequal(Sp16,Sp8)
>
> James Tursa


Dear John,

Thanks for your very prompt reply. I made the changes you suggested and got this:

---------------------------------------
>> n = 1000;
Sp16 = sprand(n,n,0.001); % <-- changed 0.10 to 0.001
Sp16(Sp16>0) = 1./Sp16(Sp16>0);
Sp16(Sp16>127) = 127;
Sp16 = floor(Sp16);
D16 = full(Sp16);
D8 = int8(D16);
Sp8 = int8tosparse(D8);
>> whos
  Name Size Bytes Class Attributes

  D16 1000x1000 8000000 double
  D8 1000x1000 1000000 int8
  Sp16 1000x1000 24008 double sparse
  Sp8 1000x1000 24008 double sparse

>> isequal(D16,D8)
isequal(Sp16,Sp8)
ans =
     1
ans =
     1
--------------------------------------

You can see that the size of Sp8 is the same as Sp16 so it
seems there is no advantage in converting.

Here is the type of problem I'm interested in:

A 10^8 x 10^8 distance (adjacency) matrix with 20 non-zeros per row,
and each non-zero is an integer in the range 1:127. Only addition
and subtraction performed on these, and the values are never changed.
This is typical of shortest path and network flow problems.

Storage for these values is 20 x 10^8 = 8 x 20 x 10^8 = 16 GB in double.
                                                                              = 2 GB in int8.

The storage for indexing information should be the same for both.

Is there some fundamental difference between sparse matrices whose
non-zeros are in the rang 1:127 and those in the range 10^(-308) to 10^(308) ?
Or am I missing something here?


Regards,

Derek O'Connor

"Derek O'Connor" <derekroconnor@eircom.net> wrote in message <gqb07m$aau$1@fred.mathworks.com>...
>
> Dear John,
>
> Thanks for your very prompt reply. I made the changes you suggested and got this:
>
> ---------------------------------------
> >> n = 1000;
> Sp16 = sprand(n,n,0.001); % <-- changed 0.10 to 0.001
> Sp16(Sp16>0) = 1./Sp16(Sp16>0);
> Sp16(Sp16>127) = 127;
> Sp16 = floor(Sp16);
> D16 = full(Sp16);
> D8 = int8(D16);
> Sp8 = int8tosparse(D8);
> >> whos
> Name Size Bytes Class Attributes
>
> D16 1000x1000 8000000 double
> D8 1000x1000 1000000 int8
> Sp16 1000x1000 24008 double sparse
> Sp8 1000x1000 24008 double sparse
>
> >> isequal(D16,D8)
> isequal(Sp16,Sp8)
> ans =
> 1
> ans =
> 1
> --------------------------------------
>
> You can see that the size of Sp8 is the same as Sp16 so it
> seems there is no advantage in converting.
>
> Here is the type of problem I'm interested in:
>
> A 10^8 x 10^8 distance (adjacency) matrix with 20 non-zeros per row,
> and each non-zero is an integer in the range 1:127. Only addition
> and subtraction performed on these, and the values are never changed.
> This is typical of shortest path and network flow problems.
>
> Storage for these values is 20 x 10^8 = 8 x 20 x 10^8 = 16 GB in double.
> = 2 GB in int8.
>
> The storage for indexing information should be the same for both.
>
> Is there some fundamental difference between sparse matrices whose
> non-zeros are in the rang 1:127 and those in the range 10^(-308) to 10^(308) ?
> Or am I missing something here?
>

Yes, it seems you have a basic misunderstanding of what I have been doing. MATLAB only supports sparse matrices of double ... it does *not* support sparse matrices of any other class. The code I provided to you simply converts the full int8 input matrix to a double sparse matrix without first converting to a full double matrix like you did originally. The C code is simply a way to avoid creating the full double matrix intermediate step ... but the final results are the same regardless of what class you start with ... you will always end up with a *double* sparse matrix. There will be *no* storage savings on the final result. It appears you were expecting to have the fundamental storage class behind the scenes of int8tosparse be int8, but that is not the case. The only way to realize this would be to create a custom class that contained the int8 raw data and indexing and then write all of
the sparse functions you need to support this custom class from scratch. That would be a huge undertaking.

James Tursa

"James Tursa" <aclassyguywithaknotac@hotmail.com> wrote in message <gqb751$a7b$1@fred.mathworks.com>...
>
> ...MATLAB only supports sparse matrices of double ... it does *not* support sparse matrices of any other class...

I meant for the operations you are interested in (addition, subtraction, etc.). MATLAB *does* of course support logical sparse matrices. Since the storage class behind the scenes for logical sparse matrices is 1-byte, I suppose this could be used as a starting point for a custom sparse int8 or uint8 class, and then write c-mex code to do the operations. But again, this would be quite a bit of work ...

James Tursa

"James Tursa" <aclassyguywithaknotac@hotmail.com> wrote in message <gqb751$a7b$1@fred.mathworks.com>...
>
> ...MATLAB only supports sparse matrices of double ... it does *not* support sparse matrices of any other class...

I meant for the operations you are interested in (addition, subtraction, etc.). MATLAB *does* of course support logical sparse matrices. Since the storage class behind the scenes for logical sparse matrices is 1-byte, I suppose this could be used as a starting point for a custom sparse int8 or uint8 class, and then write c-mex code to do the operations. But again, this would be quite a bit of work ...

James Tursa