Атомарная разреженная матрица, умножение и сравнение, чтобы избежать ошибки нехватки памяти

У меня есть две разреженные матрицы A (логические, 80274 x 80274) и B (неотрицательное целое, 21018 x 80274) и вектор c (положительное целое, 21018 x 1).

Я хотел бы найти разрешение результата (логическое, 21018 x 80274)

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

Проблема в том, что B * A имеет достаточно ненулевых значений (я думаю, 60824321, что не кажется большим, но каким-то образом вычисление spones (mat) использует более гигабайта памяти до сбоя октавы), чтобы исчерпать всю память моей машины хотя большинство из них не превышает с-1.

Есть ли способ сделать это без вычисления промежуточной матрицы mat = B * A ?

ПОЯСНЕНИЕ: Это, вероятно, не имеет значения, но B и c на самом деле являются двойными матрицами, которые содержат только целые значения (и B является разреженным).


person dspyz    schedule 14.10.2013    source источник
comment
Что такое spdiag? Я могу найти только spdiags в документе, и все равно не имеет особого смысла применять его к вектору?   -  person Luis Mendo    schedule 15.10.2013
comment
Извините, spdiag — устаревшая функция. Это эквивалентно sparse(diag(x)), но я предпочитаю его, потому что он не пытается создать полную матрицу, когда x является логическим вектором. Это не имеет отношения к этой проблеме, поэтому я изменю его.   -  person dspyz    schedule 15.10.2013


Ответы (2)


Разве вы не можете просто работать с ненулевыми значениями mat? (для нулевых значений, которые вы знаете, результат будет 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
person Luis Mendo    schedule 14.10.2013
comment
Могу, но это не решит проблему. Проблема в том, что у mat слишком много ненулевых значений, чтобы поместиться в памяти. Я хочу обойти часть, где я вычисляю мат, или, по крайней мере, делать это частями. - person dspyz; 15.10.2013
comment
Тогда вам, вероятно, понадобится цикл (скажем, по n), чтобы вычислить каждую строку mat как B(n,1)*A, а затем сравнить это с c(n)-1. Или вы можете работать в пакетах с таким количеством строк, с которым может справиться ваша память: bsxfun(@gt, B(1:1000,1)*A, c(1:1000)-1), затем bsxfun(@gt, B(1001:2000,1)*A, c(1001:2000)-1) и т. д. - person Luis Mendo; 15.10.2013

Вы можете попробовать обновиться до Octave 3.8.1, bsxfun был обновлен, чтобы быть в курсе, это значительно повышает производительность!

person gaborous    schedule 20.06.2014