From: Dima Sorkin (dsorkin_at_[hidden])
Date: 2005-03-13 04:46:42
Quoting Gunter Winkler <guwi17_at_[hidden]>:
> Am Freitag, 11. März 2005 18:15 schrieb Dima Sorkin:
> > Hi.
> > I have a following code:
> > 1: sparse_matrix<double> m1(N,N);
> > 2: // fill m1 ...
> > 3: compressed_matrix<double> m2(m1);
>
> An Alternative would be to use a coordinate_matrix to fill the data. The
> conversion from coordinate matrix to a compressed_matrix is very fast.
> The third way is to precompute the structure of the compressed_matrix
> and fill it using push_back(). See
> http://www.bauv.unibw-muenchen.de/~winkler/ublas/sparse_comparison.html
> for more details.
>
Hi.
Thank you.
Please read this message till its end, since I believe I have some serious claim
here.
There is currently no way to detect the sparse structure or even number of
non-zeros of matrix_expression, since it is unified to all (sparse and dense
matrix types).This leads to the following results (because of use of
Expr.Templates):
(I must remark that I use boost release 1.31.1 ,but I have examined the 1.32
sources and it seems that the results will be the same)
matrix_sparse<double> mt(N,N);
compressed_matrix<double> ms1(N,N),ms2(N,N);//Later filled with N
non-zeros,each
compressed_matrix<double> ms3(-ms1); // O(N^2) time
ms1 = mt ; //O(N^2) time
mt = prod(ms1,ms2) // seems O(N^2) ,I didn't have time to check it individually
Because of these ( and similar cases ) examples,I think that use of expression
templates in sparse operations, as it done now, is serious performance bug of
ublas. I may, however , miss some point. Can someone explain my mistake ?
If I am right , then some kind of "sparse_matrix_expression" class should be
introduced in ublas. For now I use some workaround: hand written functions that
access individual elements of ublas sparse metrices , and give O(N),O(N*logN)
performance in all cases I have. But for sure, they are far from being
optimal.
Dima.
P.S.
Minor remarks:
By the way it was undocumented in ublas help that compressed_matrix has direct
constructor from coordinate_matrix. According to ublas help it can be
constructed from matrix_expression only.
Second minor remark : in ublas documentation the arguments of sparse matrices
constructors written in reverse order - nnzeros comes first.