On Thu, Apr 23, 2009 at 4:11 PM, Gunter Winkler <guwi17@gmx.de> 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));