sparse function is intrisically slow?

I am building the sparse matrix using the function sparse(I,J,V,N,N)

I, J and V are already a column vector. I checked previous questions, it seems already being the best way. But, it takes extreamely long time if the size goes to 1E5.

Any way out? Thanks!

[Bjorn Gustavsson]

Is N 1e5 or the number of non-zero elements 1e5?

[Jiangtao Lu]

noe-zero element V 1E5, with matrix size N also 1E5

[Matt J]

I suggest you attach I,J,V in a .mat file so we can test it.

[Bjorn Gustavsson]

I get nothing obviously "strange" when running this snippet:

Nall = 10.^[1:8];
for i1 = 1:numel(Nall),
  iRows = randi(Nall(i1),[Nall(i1),1]);
  iCols = randi(Nall(i1),[Nall(i1),1]);
  Vals = randn([Nall(i1),1]);
  tic,M = sparse(iRows,iCols,Vals,Nall(i1),Nall(i1));toc
end
Elapsed time is 0.014283 seconds.
Elapsed time is 0.000035 seconds.
Elapsed time is 0.000135 seconds.
Elapsed time is 0.001263 seconds.
Elapsed time is 0.014471 seconds.
Elapsed time is 0.167366 seconds.
Elapsed time is 1.779908 seconds.
Error using randn
Out of memory. Type HELP MEMORY for your options.

The time increases reasonably linearly with N until I run out of memory. Transposing the inputs doesn't change timings much either. (I know that the timing is only order-of-magnitude, I know that the indices generated are not necessarily unique).

1e5 is nothing in terms of nnz elements.

[Jiangtao Lu]

Hi Bjorn, thanks. I got similar time in my real code. Perhaps you have noticed, if we futher increase the size (I need to do so in my work), this function becomes really time consuming. Is there any faster way for sparse creation?

[Walter Roberson]

Is there any faster way for sparse creation?

NO

Question: are you assigning into the matrix after you build it with sparse? If you are, then be sure to use the sixth parameter to sparse() to tell it the maximum number of non-zero entries there will be.

[Jiangtao Lu]

Thanks Walter. For now, I newly create sparse matrix every time step. Could it help, if I initialize the sparse matrix at the first time step (store in memory) and update it every time step later? The number and position of non-zero elements are unchanged but the values are all updated. I still need to resort to "sparse"? Or is there any faster re-assign way?

[Walter Roberson]

I am not sure at the moment for large vectors whether it would be faster to sparse() a new matrix into existence or to update all of the entries using linear indexing.

[James Tursa]

For the newly created sparse matrix at each iteration, how are you generating the values? Are they in a full double vector that you use with the same indexing to do the rebuild? What version of MATLAB are you using? There may be a very fast unnoficial mex solution involving modifying the values in place if you are desperate for speed, but it is a bit tricky to do without clobbering another shared variable. How do you use this sparse matrix downstream in your code? Any chance there are shared data copies (e.g., from A = your_sparse_matrix) or reference copies (e.g., from C{1} = your_sparse_matrix) of it created? How many iterations are we talking about? Do all of the values of the sparse matrix change each iteration, or only some of them?

[Jiangtao Lu]

the values are in a double vector but not generated using the same index for sparse creation.

I am using 2018b.

I have several terms with the same unknown, and then the sparse matrix of each term is added to solve Ax=b.

I am not sure I understand your "shared data copy". For each time step, I calculate V and sparse it only once. All values are changed but positions are not.

[James Tursa]

"... not generated using the same index for sparse creation ..."

"... but positions are not. ..."

These two statements seem to say the opposite. Is the sparse pattern the same for each iteration or not? First sentence seems to say they are not the same, but second sentence seems to say they are the same. Which is it?

[Jiangtao Lu]

sorry, I misunderstood the first sentence. The sparse pattern does not change (I and J are pre-calculated and stored), only the values V are updated every time step.

[James Tursa]

OK. Are the values in the double vector stored in the same memory order as the sparse matrix columns? Can you give a small example of your indexing and values vector? I'm asking all of these questions to see if a mex routine could work for you.

[Jiangtao Lu]

I am not clear to that question. The calculation is like:

work.jacValueIndex = numC*(work.g(work.J)-1)+work.I;
jacValue = df(work.jacValueIndex)./dc(work.J);
jac = sparse(work.I,work.J,jacValue,numC,numC);

I attached an example of work and jacValue

Answer 1

[Robert]

I was able to speed my sparse matrix build for very large problems by using sparse2.m (mexFunction) from SuiteSparse.