Re: [ublas] ublas Digest, Vol 63, Issue 1

Subject: Re: [ublas] ublas Digest, Vol 63, Issue 1
From: Jörn Ungermann (j.ungermann_at_[hidden])
Date: 2010-02-01 17:57:45

> Message: 1
> Date: Mon, 1 Feb 2010 09:54:36 +1100
> From: "Paul C. Leopardi" <paul.leopardi_at_[hidden]>
> To: ublas_at_[hidden]
> Subject: Re: [ublas] Slow (seeming) sparse matrix multiplication?
>       (James  Davidson)
> Message-ID: <201002010954.36293.paul.leopardi_at_[hidden]>
> Content-Type: Text/Plain;  charset="utf-8"
>
> On Saturday 30 January 2010 23:30:13 J?rn Ungermann wrote:
> > Hello James,
> >
> > we had the same problem. You might want to have a look at the
> > "axpy_prod" and "sparse_prod" functions, which offer specialized
> > implementations for sparse matrices. However, we found all sparse matrix
> > products severely lacking and finally implemented the sparse matrix
> > multiplications we needed ourselves (which are mostly several orders of
> > magnitude faster than the best ublas routine). Sadly, we could not spend
> > the effort to do a implementation general enough and worthy to be
> > included into ublas.
> >
> > But if you need something along the lines of compressed_matrix *
> > compressed_matrix -> matrix, send me an email and I can send you our
> > code. We have not implemented all products with respect to combinations
> > of transposedness/row-column-majority, but the missing ones could be
> > easily added given the existing implementations. As the product of two
> > sparse matrices is seldom sparse itself (and certainly not in our
> > use-case) we have not thought about a compressed_matrix as target, which
> > certainly requires a special order in which the resulting matrix is to
> > filled.
>
> Hi J?rn ,
> GluCat ( http://glucat.sf.net ) also uses uBLAS compressed matrices. It has
> special products for monomial matrices (one non-zero per row and column) and
> uses ublas::sparse_prod otherwise. See glucat/matrix_imp.h. I'd also be
> interested in seeing your code.
> Best, Paul

Hi Paul,

the sparse_prod optimizes IMHO only according to the target matrix
majority (i.e. column or row major) and not according to the "source"
matrices.
By using the same techniques that were used for axpy_prod for
compressed_matrix-vector products it is very simple, if tedious, to
write specializations for all different combinations of majority and
transposedness of the source matrices, which we therefore did only for
those combinations needed by us and also only for ublas::matrix target
matrices (as the product of sparse matrices is seldom sparse itself and
certainly not in our use case). However, one can probably also optimize
for different target matrices (i.e. random-access or some kind of
compressed_matrix) at the cost of *even more* specializations and
probably run time.

I'll send you our code separately, but as I said, it is not very pretty
(else I'd have tried to get it into uBLAS).

Kind regards,
Jörn Ungermann
--
Jörn Ungermann
Dipl.-Mathematiker
Institut für Chemie und Dynamik der Geosphäre: Stratosphäre
Forschungszentrum Jülich
Telefon: [phone]
------------------------------------------------------------------------------------------------
------------------------------------------------------------------------------------------------
Forschungszentrum Juelich GmbH
52425 Juelich
Sitz der Gesellschaft: Juelich
Eingetragen im Handelsregister des Amtsgerichts Dueren Nr. HR B 3498
Vorsitzende des Aufsichtsrats: MinDir'in Baerbel Brumme-Bothe
Geschaeftsfuehrung: Prof. Dr. Achim Bachem (Vorsitzender),
Dr. Ulrich Krafft (stellv. Vorsitzender), Prof. Dr.-Ing. Harald Bolt,
Prof. Dr. Sebastian M. Schmidt
------------------------------------------------------------------------------------------------
------------------------------------------------------------------------------------------------