Atomic sparse-matrix-multiply-and-compare to avoid out-of-memory error

Question

I have two sparse matrices A (logical, 80274 x 80274) and B (non-negative integer, 21018 x 80274) and a vector c (positive integer, 21018 x 1).

I'd like to find the result res (logical, 21018 x 80274) of

mat = B * A;
res = mat > sparse(diag(c - 1)) * spones(mat);

# Note that since c is positive, this is equivalent to
# res = bsxfun(@gt, mat, c-1)
# but octave's sparse bsxfun support is somewhat shoddy,
# so I'm doing this instead as a workaround

The problem is B * A has enough nonzero values (I think 60824321 which doesn't seem like a lot, but somehow, the computation of spones(mat) uses up over a gigabyte of memory before octave crashes) to exhaust all my machine's memory even though most of these do not exceed c-1.

Is there a way to do this without computing the intermediate matrix mat = B * A ?

CLARIFICATION: It probably doesn't matter, but B and c are actually double matrices that happen to only be holding integer values (and B is sparse).

What is spdiag? I can only find spdiags in the doc, and still it doesn't seem to make much sense to apply it to a vector? — Luis Mendo
Sorry, spdiag is a deprecated function. It's equivalent to sparse(diag(x)), but I prefer it because it doesn't try to create a full matrix when x is a logical vector. It's not relevant for this problem, so I'll change it. — dspyz

Luis Mendo Luis Mendo · Accepted Answer · 2013-10-14T21:19:04

Can't you just work on the nonzero values of mat? (for the zero values you know the result will be 0):

c = c(:); % make c a column
ind = find(mat>0); % linear index
[row, ~] = ind2sub(size(mat),ind); % row index within mat (for use in c)
res = mat(ind) > c(row)-1; % results for the nonzero values of mat

Atomic sparse-matrix-multiply-and-compare to avoid out-of-memory error

2 Answers