Re: [ublas] Issue calculating matrix product with diagonal matrix

Subject: Re: [ublas] Issue calculating matrix product with diagonal matrix
From: sguazt (marco.guazzone_at_[hidden])
Date: 2011-03-15 16:58:30

On Tue, Mar 15, 2011 at 7:17 PM, Steven Butlin <steven_at_[hidden]> wrote:
> Dear all,
>
>
>
> I have encountered an issue calculating matrix products when the diagonal
> matrix is used.  I have written a simple test case to illustrate the problem
> (see code).  The calculation is fairly straight forward.   However, a static
> assertion fails in matrix_expression.hpp at line 4815.
>
>
>
> Has anyone else encountered this problem?  Can anyone suggest a solution?
>
>

See here:

http://www.crystalclearsoftware.com/cgi-bin/boost_wiki/wiki.pl?Effective_UBLAS

at the section titled "Controlling the complexity of nested products".

As suggested in that section you can solve your problem by writing,
for instance, the following program:
--- [code] ---
#include <boost/numeric/ublas/banded.hpp>
#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/io.hpp>
#include <iostream>

namespace ublas = ::boost::numeric::ublas;

int main()
{
    ublas::matrix<double> A(3,3);

    A(0,0) = 4; A(0,1) = 1; A(0,2) = 1;
    A(1,0) = 1; A(1,1) = 2; A(1,2) = 1;
    A(2,0) = 1; A(2,1) = 1; A(2,2) = 3;

    ublas::diagonal_matrix<double> D(3,3);

    D(0,0) = 1; D(1,1) = 2; D(2,2) = 3;

// Don't use this ...
//  ::std::cout << prod(prod(A,D),A) << ::std::endl;

// Instead use this...
    ublas::matrix<double> T;
    T = prod(A,D);
    T = prod(T,A);
    ::std::cout << "T=" << T << ::std::endl;
}
--- [/code] ---

Best,

-- Marco