[ublas] Troubles with triangular matrix - vector product

Subject: [ublas] Troubles with triangular matrix - vector product
From: Florent Teichteil (florent.teichteil_at_[hidden])
Date: 2012-10-05 18:57:45

Hi,

I'm trying to multiply a triangular matrix adapted from a banded
matrix, by a vector.
I'm using the latest blas bindings from boost sandbox.

First problem: it seems that the boost::numeric::bindings::trans
function does not work properly with triangular matrices.
For instance, the following piece of code correctly compiles, but
makes the backend BLAS library generate a runtime error (with ATLAS
backend, I get "On entry to DTRMV  parameter number  6 had an illegal
value"):

#include <boost/numeric/bindings/ublas.hpp>
#include <boost/numeric/bindings/blas.hpp>
#include <boost/numeric/bindings/trans.hpp>
using namespace boost::numeric::ublas;

int main () {
        banded_matrix<double, column_major> m(1000, 1000, 1, 1000);
        triangular_adaptor<banded_matrix<double, column_major>, upper> mT(m);
        vector<double> v(1000);
        boost::numeric::bindings::blas::trmv(boost::numeric::bindings::trans(mT), v);
        return 0;
}

If we do not transpose mT in the previous code, it executes correctly.
(By the way, I am not sure that boost::numeric::bindings::trans does
not change the transpose tag of mT for the next operations, since it
looks like it modified some internal tags of mT)

Second (independent) problem: it seems that we can not use trmv in
conjunction with column-major matrix proxies. The following code does
not compile:

#include <boost/numeric/bindings/ublas.hpp>
#include <boost/numeric/bindings/blas.hpp>
#include <boost/numeric/bindings/trans.hpp>
using namespace boost::numeric::ublas;

int main () {
        banded_matrix<double, column_major> m(1000, 1000, 1, 1000);
        matrix_range<banded_matrix<double, column_major> > mR = project(m,
range(100, 1000), range(100, 1000));
        triangular_adaptor<matrix_range<banded_matrix<double, column_major>
>, upper> mT(mR);
        vector<double> v(900);
        boost::numeric::bindings::blas::trmv(mT, v);
        return 0;
}

Yet, if we replace all occurrences of column_major by row_major in the
previous program, then it compiles correctly. It's really related to
using matrix proxies AND trmv, because the code compiles either if we
remove the the use of mR (see 1st program of the e-mail) or if we
remove the call to trmv.

Does it sound like bugs in the numeric bindings, or am I missing something?

Cheers,
Florent