Boost logo

Ublas :

Subject: Re: [ublas] Problems with nested product compilation?
From: Jesse Perla (jesseperla_at_[hidden])
Date: 2009-04-23 18:58:20


On Thu, Apr 23, 2009 at 4:11 PM, Gunter Winkler <guwi17_at_[hidden]> wrote:

> maybe you have to use boost::numeric::ublas::prod in order to force the
> use of the correct namespace.

Thanks Gunter,
I messed with the namespaces to no avail. I did, however, read a bit of
matrix_expression.hpp to track down what was happening.

The following (given in the docs) doesn't work:
D = prod(A, prod<matrix<double>>(B,C))

However the following does when you explicitly list the 3 matrix types to
get that overload of "prod".
D = prod(A, prod<matrix<double>,matrix<double>,matrix<double>>(B,C))

I think I will stick with the following which works fine and is less clumsy:
D = prod(A, matrix<double>(prod(B,C)));

The following also works fine, which is also relatively pretty:
matrix<double> temp(2,2);
D = prod(A, prod(B, C, temp));

Unless I am missing something in the source code, it appears that the
"prod(A, prod(B,C)) " N^3 version of nested products with no temporaries has
been removed and there is a static assertion under all circumstances to
prevent this style. Am I reading this correctly?

I got confused on all of this stuff because of the current ublas docs:
Look at the bottom Q&A question:
http://www.boost.org/doc/libs/1_38_0/libs/numeric/ublas/doc/index.htm
And the last section of:
http://www.boost.org/doc/libs/1_38_0/libs/numeric/ublas/doc/operations_overview.htm

These appear to be out of date. If they are, then I suggest getting rid of
the Q&A question at the bottom of the index since it is no longer relavent,
and then changing the text at the bottom of operations overview to be
something like:

Nested matrix products

Since they have a higher algorithmic complexity, nested matrix products
without the use of a temporary are not supported in uBLAS. You will not be
able to use the following code

 R = prod(A, prod(B,C)); //Does not compile

Instead, uBLAS provides 3 alternative syntaxes for this purpose given a user
chosen temp_type:

 temp_type T = prod(B,C); R = prod(A,T); //Manual use of a temporary

 prod(A, temp_type(prod(B,C)); //Constructs a temp_type with the
prod(B,C) expression

 prod(A, prod(B,C,T)); //Uses the T value as a temporary for
calculation in preallocated T.

Users will want to manage their own temporary matrix types for nested matrix
products to ensure maximum efficiency. For example:

matrix<double> A(2,2), B(2,2), C(2,2); //... fill in

matrix<double> T(2,2); //A natural, preallocated temporary

prod(A, prod(B,C,T));

or to have dynamic allocation:

prod(A, matrix<double>(prod(B,C))); //Dynamically allocates temporary

or use a more efficient matrix type given the structure of A,B,C:

bounded_matrix<double, 2, 2> A, B, C; //... fill in

prod(A, bounded_matrix<double, 2, 2>(prod(B, C));