Possible bug in MTL?

Hello.

We are using MTL to develop a MultiGrid method, and we have  =

encountered a problem that we cannot identify.

The problem can be shown by the following simple example:

#include <iostream>
#include <mtl/matrix.h>
#include <mtl/mtl.h>
#include <mtl/utils.h>
#include <mtl/linalg_vec.h>

using namespace mtl;
using namespace std;

typedef matrix< double, symmetric<upper>,
dense<>, row_major>::type Matrix;

typedef dense1D<double> Vec;

int main()
{
     const Matrix::size_type matrix_size=3D5;

     Matrix A(matrix_size,matrix_size);
     Vec y(matrix_size,0), x(matrix_size, 1), z(matrix_size);

     int r;
     for (r =3D 0;r < matrix_size;++r) {
         (A)(r,(r)) =3D 1;
         if (r + 1 < matrix_size)
             (A)(r,(r+1)) =3D 1;
     }
     mult(A,x,y,z);

     return 0;
}

By following with the debugger the execution of the mult() function,  =

we have found that the inner iterator in the version of mult() that  =

finally gets called goes out of bounds. It looks like it always  =

iterates to the entire size of the matirx, but the matrix is upper  =

triangular. As a result, we get writes outside the bounds of the  =

allocated vectors.

We experienced the same problem compiling under linux and under  =

Windows, so it does not seem to be a system-specific bug.

Could it be that we are choosing the wrong declaration for the matrix/ =

vector types?

				Thank you in advance

				  Antonio Frangioni



Antonio Frangioni    -    Associate Professor    -    Universit=E0 di Pisa
Dipartimento di Informatica - Largo B. Pontecorvo, 3 - 56127 Pisa, Italy
Ph: +39 050 2212789, Fax: +39 050 2212726, Web: www.di.unipi.it/~frangio




_______________________________________________
This list is archived at http://www.osl.iu.edu/MailArchives/mtl-devel/