Re: [ublas] blocked products and relation UBLAS/GLAS

From: Gunter Winkler (guwi17_at_[hidden])
Date: 2006-02-15 06:53:13

On Tuesday 14 February 2006 19:37, Riccardo Rossi wrote:
> dear list,
>          i am more or less a Newbie at using ublas for sparse operations
> and i would love a little advice on some delicate questions.
>
> in particular i would like to work with a sparse matrix in which all of the
> single terms are vectors or matrices instead of double. All the containers
> are templated matrices which implies that i am free to do it
> nevertheless...

You can use the matrix and vector classes as a container for any type provided 
the type is default and zero constructable. That means you need a constructor 
like

my_type(const int zero) : ... 
{ assert( 0 == zero ); }

which is necessary to let the sparse types know what "zero element" means.

> how can i make the library to know how to perform the multiplication
> between all of the different terms? the case is for example of a
> SparseMatrix of Matrices to be multiplied by a vector of vectors.

the easiest way is to use a custom type as value type which behaves like a 
scalar. That means you have to overload the usual arithmetic operators for 
your type.

> Similarly i would like to create a matrix doing something like A =
> trans(B)*B  (i think it is called Symmetric Product)... is it available in
> ublas?

this is more complicated because ublas assumes that a scalar is symmetric and 
therefore never calls trans(x) for any scalar x. Maybe you can create a new 
specialization of 

matrix_unary2_traits (in matrix_expression.hpp, lines 1648ff)

that defines a new result type.

    // (trans m) [i] [j] = m [j] [i]
    template<class E>
    BOOST_UBLAS_INLINE
    typename matrix_unary2_traits<const E, 
          scalar_identity<typename E::value_type> >::result_type
    trans (const matrix_expression<E> &e) {
        typedef typename matrix_unary2_traits<const E,
                  scalar_identity<typename E::value_type> 
               >::expression_type expression_type;
        return expression_type (e ());
    }

(you have to replace scalar_identity by scalar_transpose and write you own 
scalar_transpose, see functional.hpp - ask me if you need help doing this.)

> finally i would like to know what is the relation of ublas to GLAS

GLAS is (will be) a general linear algebra library supporting different back 
ends, e.g. ublas.  Maybe Karl can better comment on that.

mfg
Gunter