sparse matrix operations to gain efficiency

7 次查看（过去 30 天）
AM
AM 2019-7-15
编辑： Cedric 2019-7-16
[EDITED]
Hello,
I've never really used sparse matrices before but I was wondering if it would be useful in my case to imporve efficiency. My code is as follows
reac=25;spec=20;
kf=rand(1,reac); kb=rand(1,reac);
C=rand(1,spec);
fwd=zeros(reac,spec);bcwd=zeros(reac,spec);
for i=1:reac
a=randi([1 3],1,5);b=randi([1 spec],1,5);
fwd(i,b)=a;
c=randi([1 3],1,5);d=randi([1 spec],1,5);
bcwd(i,d)=c;
end
i=1:reac;q=kf(i)'.*prod(C.^fwd(i,:),2)-kb(i)'.*prod(C.^bcwd(i,:),2);
sto=fwd-bcwd;
i=1:spec;rate(i)=sum(sto(:,i).*q,1);
I've used an example of reac=25 and spec=20 but in reality I could have reac up to 1500 and spec up to 500. I use the values of rate in differential equations which I solve using ode15s. This is the part of my code that takes the longest, where fwd and bcwd are sparse matrices. Would truning them into sparse matrices using sparse help me? if so, how can I rewrite my code?
An help is appreciated,
Thank you!
  3 个评论
显示 1更早的评论
AM
AM 2019-7-16
I didn't know how to force the elements of b and d to be unique in my example but since it's maximum 5 non-zero values per row it is not a problem. I also realize now I can simply use the sparse command and the rest of my code does not change.
Cedric
Cedric 2019-7-16
编辑：Cedric 2019-7-16
It is likely more complex than this.
The gain in performance obtained from using sparse matrices generally involves building vectors I, J, and V of row and column indices and values of non-zero elements, and building a sparse matrix in one shot with a call like
S = sparse( I, J, V, nRows, nCols ) ;
You can easily build I using REPELEM and V using RANDI. But building J is not direct and you must enforce unique values (per row) otherwise the accumlating behavior of SPARSE will produce a side effect that you may not want.
To illustrate this behavior, look at the following:
>> S = sparse( [1 2 1], [1 2 1], [10 20 30], 2, 2 )
S =
(1,1) 40
(2,2) 20
Here you can see that as there were two elements (10 and 30) with same indices (1,1), the values were accumulated (summed) and S(1,1)=10+30=40.
PS: if you wanted to enforce uniqueness, you could use RANDPERM, but it won't allow you to generate J in one shot as far as I can tell.

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

John D'Errico
John D'Errico 2019-7-16
Are those matrices mostly zero? No. It looks like they are 50% zeros. Sparse will not gain anything, because thre is also overhead in working with computations on sparse matrices.
Next, sparse matrices gain when you are doing matrix multiplications, and solving linear systems of eequations, but mainly only then if they are really sparse. 50% zeros does not cut it.
  3 个评论
显示 1更早的评论
John D'Errico
John D'Errico 2019-7-16
Sparse matrices work just like other matrices. So you use them just like you use ANY other normal matrix.
Are you really creating them randomly? If so, then just use the sparse command to make them into sparse matrices. Calling sparse on a full matrix will return a new matrix that is sparse version of that matrix. (The function full does the opposite, as you would expect.)
You could create them directly using the sparse commnd also, but that takes more skill and understanding of how to build a sparse matrix to do it well, and far more than I will teach you.
But you should also recognize that CREATING the sparse matrix as such willl then add some time. You could also find that you won't gain much, if at all. So if some of the time, your matrices are not that sparse, expect it all to take more time when using sparse matrices. Only you can test it to see. But converting to sparse is simple as I said.
AM
AM 2019-7-16
I realize now there is no difficulty in using sparse matrices here, I will test out my code and see which way is more efficient, thank you for your time!
(I had tried to use sparse matrices but had made a mistake in my code so I thought they were handled differently, sorry)

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Sparse Matrices 的更多信息

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by