4
$\begingroup$

I have a large $d\times d$ SparseArray $P$, and I want to multiply a vector (also a SparseArray) $v$ of length $d$ elementwise with each row of $P$. The value of $d$ is currently $2^{16}=65536$. The matrix $P$ has $d$ entries with the non-default value (which is 1, default is 0). The vector $v$ has $d/2$ non-default entries (again 1's, default 0).

I tried two approaches:

Transpose[v * Transpose[P]]

and

# * v & /@ P

The first version seems to explode in memory, as it immediately brings my computer to a crawl and memory usage goes beyond 90% (i7, 16GB RAM). In the end, the computation officially runs out of memory.

The second version works, but is very slow and uses no more than 4% of my memory. I should note that if I use numerical values in the matrix, the second version seems to parallelize automatically as CPU usage goes to 400%, but it is still very slow. ParallelMap didn't work apparently, because it is a SparseArray.

I wonder if there is not a middle way to speed up the computation without blowing my memory?

$\endgroup$
  • 7
    $\begingroup$ How about the more formal P.SparseArray[Band[{1, 1}] -> v]? $\endgroup$ – J. M. will be back soon Feb 26 '16 at 22:24
  • $\begingroup$ @J.M. Amazing! This is basically instantaneous! So the element wise multiplication should generally be avoided and replaced with matrix multiplication in this way? $\endgroup$ – tortortor Feb 26 '16 at 22:39
  • 3
    $\begingroup$ Well, it's more natural to use dot products rather than Hadamard products when dealing with sparse arrays. In this case, you are fortunate that the Hadamard product you want can be expressed as an appropriate product with a diagonal matrix. $\endgroup$ – J. M. will be back soon Feb 26 '16 at 22:56
7
$\begingroup$

I'll use a smaller example for this answer, since my computer really isn't up to such large sizes:

BlockRandom[SeedRandom[42];
            mat = SparseArray[{j_, k_} :> RandomInteger[], {512, 512}];
            vec = SparseArray[RandomInteger[1, 512]];]

One way to (roughly) look at the memory consumption of an evaluation is to use the function MaxMemoryUsed[], which returns a result in bytes. Let's try it out on the two Hadamard-based variants in the OP:

MaxMemoryUsed[Transpose[vec Transpose[mat]]]
   2642360

MaxMemoryUsed[(# vec) & /@ mat]
   2462792

So, a bit more than 2 MB for these two. How about the dot product version?

MaxMemoryUsed[mat.SparseArray[Band[{1, 1}] -> vec]]
   543656

Quite a drop in memory consumption! Now, I am not at all privy to how SparseArray[] is actually implemented under the hood, but I will speculate that dot products are more natural to do for the internal storage format used than Hadamard products. (At the very least, the sparse matrix formats I am familiar with are certainly optimized for taking dot products with.) To do a Hadamard product, probably some extra conversion/processing has to be done, which increases the memory consumption and effort.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.