I welcome any help for the following code optimization problem:
I have a collection of N sparse matrices of identical sizes ([s1 s2]
) stored in a cell array A
and a corresponding number of scalar weights stored in an vector w
. I want to compute the sum of all the matrices in A
weighted by the values stored in w
. Through the iterations of my program, only the values in w
change. I can therefore compute a priori the number of non-zero elements in my result and pre-allocate some memory for it using spalloc
.
For the moment I have something like:
result = spalloc(s1,s1,number_of_non_zero);
for i=1:N
result = result + w(i)*A{i};
end
I really need do optimize this part which for the moment takes most of the computing time in my program (checked with the profiling tools).
Some additional information:
-The above code runs millions of times so even minor improvements are welcome.
-The matrices in A
come from a finite element code (1D or 2D)
-I have no problem moving away from the cell structure if I can save some time (like using cell2mat(A)
)
Thank you for any hint on how to speed up this part of the code.
A.
I can't say for sure if this solution will be more computationally-efficient, but it's something else to try...
If you really have to represent your matrices as sparse (i.e. full matrices take up too much memory) and the built-in sparse representation in MATLAB isn't giving you the performance you desire, then you can try representing the sparse matrices in a different way. Specifically, you can represent them as N-by-3 matrices where the first two columns contain the row and column indices into the matrix for all the non-zero values and the third column contains the non-zero values. You can convert your data to this form using the function FIND like so:
for iMatrix = 1:numel(A)
[r,c,v] = find(A{iMatrix});
A{iMatrix} = [r c v];
end
Each time you need to compute the weighted sum of these matrices, you first need to multiply the values by the weights:
B = A; %# Store a temporary copy of A
for iMatrix = 1:numel(B)
B{iMatrix}(:,3) = w(iMatrix).*B{iMatrix}(:,3);
end
Then, you can compute the final sum using the function ACCUMARRAY:
B = vertcat(B{:}); %# Convert B from a cell array to an N-by-3 matrix
result = accumarray(B(:,1:2),B(:,3));
The variable result
in this case will be a full matrix. If you need result
to be a sparse matrix, you can add extra arguments in the call to ACCUMARRAY like so:
result = accumarray(B(:,1:2),B(:,3),[],[],[],true);