Boost logo

Ublas :

Subject: [ublas] Slow (seeming) sparse matrix multiplication?
From: James Davidson (james_at_[hidden])
Date: 2010-01-27 17:56:25


Howdy,

I am attempting to multiply two 1000x1000 compressed_matrix<double> matrices
and the computation time seems excessive. When I compile the code attached
at the end of this message I get the results:

> g++ -DNDEBUG -I"/home/james/local/include/boost-1_37" -O3
test_sparse_mult.cpp -o test_sparse_mult
> time ./test_sparse_mult

real 0m10.609s
user 0m10.381s
sys 0m0.012s

I also compiled the code without the NDEBUG flag and obtained:

> g++ -I"/home/james/local/include/boost-1_37" -O3 test_sparse_mult.cpp -o
test_sparse_mult
> time ./test_sparse_mult

real 0m10.004s
user 0m9.373s
sys 0m0.568s

As a reference, I ran substantially similar code in Matlab and achieved many
orders of magnitude improvement in running time--reducing the execution time
to 0.000262 seconds.
Does this seem correct? Is there something I am missing? Any help on this
matter will be greatly appreciated. Thanks.

--James

//
-------------------------------------------------------------------------------
#include <boost/numeric/ublas/matrix_sparse.hpp>

#include <boost/random.hpp>
#include <boost/random/variate_generator.hpp>
#include <boost/generator_iterator.hpp>
#include <boost/random/uniform_01.hpp>
#include <ctime>
//define uniform sampling function
typedef boost::lagged_fibonacci607 base_generator_type;
extern base_generator_type range;
extern boost::uniform_01 < base_generator_type > uni_01;
base_generator_type range(static_cast<int>(std::time(0)));
boost::uniform_01 < base_generator_type > uni_01(range);

using boost::numeric::ublas::compressed_matrix;

int main(int argc, char* argv[])
{
    unsigned int size = 1000;
    compressed_matrix<double> A(size,size);
    for(unsigned int i = 0; i < size; ++i)
    {
        unsigned int x = floor(uni_01()*size);
        unsigned int y = floor(uni_01()*size);
        A(x,y) = uni_01();
    }
    compressed_matrix<double> B(size,size);
    for(unsigned int i = 0; i < size; ++i)
    {
        unsigned int x = floor(uni_01()*size);
        unsigned int y = floor(uni_01()*size);
        B(x,y) = uni_01();
    }

    compressed_matrix<double> C = prod(B,A);

}
//-------------------------------------------------------------------------------------