Why does MATLAB calculates when using speye?
13 views (last 30 days)
Show older comments
Hello,
i just noticed a weird thing
>> tic;speye(10000)*speye(10000);toc
Elapsed time is 0.000753 seconds.
>> tic;speye(100000)*speye(100000);toc
Elapsed time is 0.007663 seconds.
>> tic;speye(1000000)*speye(1000000);toc
Elapsed time is 0.073734 seconds.
>> tic;speye(10000000)*speye(10000000);toc
Elapsed time is 0.737171 seconds.
Where is the sense in using an implemented sparse unit matrix for calculation? Same goes for a multiplication between a dense matrix and a sparse unit matrix.
>> tic;eye(1000)*speye(1000);toc
Elapsed time is 0.032003 seconds.
>> tic;eye(10000)*speye(10000);toc
Elapsed time is 1.331953 seconds.
and between zero matrices ...
>> tic;sparse(10000,10000)*sparse(10000,10000);toc
Elapsed time is 0.000416 seconds.
>> tic;sparse(100000,100000)*sparse(100000,100000);toc
Elapsed time is 0.001833 seconds.
>> tic;sparse(1000000,1000000)*sparse(1000000,1000000);toc
Elapsed time is 0.012763 seconds.
>> tic;sparse(10000000,10000000)*sparse(10000000,10000000);toc
Elapsed time is 0.111892 seconds.
>> tic;eye(1000,1000)*sparse(1000,1000);toc
Elapsed time is 0.013433 seconds.
>> tic;eye(10000,10000)*sparse(10000,10000);toc
Elapsed time is 0.753174 seconds.
I thought this datatype is also designed to avoid such unnecessary performance lacks.
This does really matter for me, because in my calculation some matrices will be unit matrices in special cases. So i created an independent calculation function which does not care about the matrices, it only calculates them. In the other part i check the conditions and set some matrices to unit sparse matrices, because i want to avoid the calculation time and send them to the calculation function. Unfortunately this does not work as expected and i get an unnecessary performance lack.
0 Comments
Answers (1)
Jan
on 11 Jul 2012
You do not measure the time for the multiplication, but the time for the creation of the matrix is included also. Please post the reuslts of:
tic;speye(10000)*speye(10000);toc
A = speye(10000);
tic; A * A; toc
tic;speye(10000000)*speye(10000000);toc
A = speye(10000000);
tic; A * A; toc
The sense of the sparse arrays is that e.g. the multiplication considers the 0 elements and does not perform a multiplication for them. In addition the storage of sparse matrices is more dense than a full matrix, such that speye(10000000) is possible even when eye(10000000) would cause an out-of-memory already. It is not the idea of sparse matrices to short-cut multiplications with the identity matrix and it is easy to avoid such actions programmatically.
4 Comments
Jan
on 12 Jul 2012
Edited: Jan
on 12 Jul 2012
Of course using sparse arrays is common. When I started programming in the early 80th, sparse data types have been even more important than today, because 64MB RAM have been very expensive.
I understand your point of premature optimization. Exactly the same argument matchs, when the handling of special values is included in the libraries of a basic data type. "Sparse" matrices are defined such, that their values are not stored sequentially, but only the non-zero values are stored together with a list of their indices. Incorporating specific operations for rarely used exceptional values would mean to embed high-level functionality in the low-level data type.
Exploiting the block-structure of an array is definitely a high-level function and it is usual to catch such exceptions manually. My code example would be very ugly, when implemented directly. But when encapsulated in a user-defined class as "SmartSparse", which has e.g. the fields "isEye", "isNull" and "data", and these flags are checked in dedicated inv(), times(), mtimes(), rdiv(), ldiv() etc., you can use your "km = ..." line without changes. In addtion you can compare the timings and results of your function just by defining B as standard sparse double array, or as your SmartSparse class.
I'd call the including this feature in the basic datatype "premature", because it matters your 300*15*8 mutliplications and divisions, while the billions of multiplications of users with standard cases would waste time for the additional tiny checks.
Adjusted linear algebra solvers, which consider the specific block-sparse structure of the system, are common in scientific software for e.g. the simulation of multi-body-systems. E.g. "natural coordinates" (Garcia de Jalon) lead to large sparse mass matrices, with a 1 in the half of the diagonal elements. Then the "isEye" and "isNull" flags would be too trivial, but only hand-coded LA solvers can be efficient. Therefore I repeat: User-specific exceptions should be managed in user-specific code.
Btw, how should speye()-eps be handled? Or 2*speye()-speye()? The number of exceptions will grow exponentially also, when you really catch all possible optimizations.
See Also
Categories
Find more on Loops and Conditional Statements in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!