Creating a sparse vector inside a for loop.

Suppose I have a vector u that does not change. I want to create another vector v = M*u, where the matrix M depends on the index i . For example,

v = sparse(N,1);
for i=1:N
  M = makeSomeMatrix(i);
  v(i) = (M*u)'*u;
end

Is there a faster way to fill in v? This is the slowest part of my code and it gives me the warning " This sparse indexing expression is likely to be slow."

[Amit]

How big is N?

[Alex]

N is between 2000 and 4000.

Answer 1

[Amit]

The sparse matrix is slow but effective for huge matrixes with a lot of zeros.

but MATLAB can handle big matrixes in normal way quite effectively. In you case, it seems like probably you can work without a sparse matrix. I said that simple because you're looping through N and the matrix v is initiated as Nx1.

[Amit]

if N is between 2000 and 4000, dont use sparse. simply initiate v as

v = zeros(N,1);

[Amit]

With the new code, I would suggest that sparse is unnecessary. It will be much faster to use normal matrix for such a small matrix size.

You can understand the difference between sparse matrix and a normal matrix. In normal matrix, the whole matrix is stored at a location in memory. So replacing a value of a given index, does not produces significant overhead. Now take it for a sparse matrix. Sparse matrix only stores non-zero value. So, adding a new non-zero matrix over and over again makes MATLAB store the new sparse matrix in a new location which is significant overhead in a loop.

Answer 2

[Matt J]

Based on what you've shown, M is a row vector (that's the only way the scalar location v(i) can accept the result). You should really construct/concatenate all rows of the matrix in advance and do a single matrix-vector multiplication.

 v=M_all*u

It may help if you show us makeSomeMatrix, to see if/why each row needs to be built separately and if/why it requires a for-loop.

[Alex]

Hi Matt, you're right. That was a mistake. I've changed the example code so that everything is dimensionally correct. The matrix M changes every time and the number (M*u)'*u is stored in an element of the vector v .

[Matt J]

Alex,

Your change doesn't affect my remarks in any substantive way. You are still generating one row/column (M*u) of a bigger matrix in every iteration of the loop, and the computation of v could be vectorized if you saved and post concatenated them, e.g.,

 c=cell(1,N);

 for i=1:N
  c{i} = makeSomeMatrix(i)*u;
 end

 v=u.'*[c{:}];

The question of whether you could avoid N calls to makesomeMatrix is still open as well, and will be until we understand what's going on inside that function.