[Lapackpp-devel] bugs in operator= for LaIndex/LaGenMatComplex/LaGenMatDouble?

[Brian W. <bw...@cs...>]

Hello,

I think I've uncovered a few bugs.

In LaIndex::oeprator=(const LaIndex &i), the following definition
ignores the increment of the LaIndex it is copying:

        inline LaIndex& operator=(const LaIndex& i){
            start_=i.start_; inc_=1; end_=i.end_;
            return *this;}

I think inc_=1 should be inc_=i.inc_:

        inline LaIndex& operator=(const LaIndex& i){
            start_=i.start_; inc_=i.inc_; end_=i.end_;
            return *this;}

Also, LaGenMatComplex::operator=(COMPLEX s) and
LaGenMatDouble::operator=(double s) seem to be problematic since both
assign to the vectorDouble v without accounting for the active view
of the LaGenMatDouble object.  I believe this means that if we assign
(via operator=) through a submatrix view, the entire matrix will be
effected.

So, instead of doing

    v = s;

in both LaGenMatComplex::operator=(COMPLEX s) and
LaGenMatDouble::operator=(double s), we should instead have:

    int i,j;

    for (j=0; j<size(1); j++)
    {
        for (i=0; i<size(0); i++)
        {
            (*this)(i,j) = s;
        }
    }

and something similar for the complex case.  Note that this
mirrors operator+= (obviously, with += replaced by =).

Below is some example code (taken essentially from the Lapack++ 1.1
manual when it was still being developed by NIST/Tennessee/Oak Ridge).
I've included what I believe to be the correct output and also
the output I get before making the above changes.

By the way, this differs from the example in the manual only
in so far as assignments are replaced by references.
i.e., B = A(II,II) is replaced below by B.ref(A(II,II)).
If I'm reading that manual correctly, B = A(II,II) should
establish B as a reference to elements in A.  Instead it
appears to perform a (deep) copy.

Thanks,
Brian

int main(int argc, char **argv)
{

  LaGenMatDouble A(10,10), B, C;
  LaIndex II(1,9,2);
  LaIndex JJ(1,1,1);

  B.ref(A(II,II));   // B references the following elements of A:
                     //    (1, 1), (1, 3), (1, 5), (1, 7), (1, 9)
                     //    (3, 1), (3, 3), (3, 5), (3, 7), (3, 9)
                     //    (5, 1), (5, 3), (5, 5), (5, 7), (5, 9)
                     //    (7, 1), (7, 3), (7, 5), (7, 7), (7, 9)
                     //    (9, 1), (9, 3), (9, 5), (9, 7), (9, 9)
  B(2,3) = 3.1;      // Set A(5,7) = 3.1

  C.ref(B(LaIndex(2,2,4), JJ));
                     // C references the following elements of B:
                     //    (2, 1)
                     // C references the following elements of A:
                     //    (5, 3)

  C = 1.1;           // Set B(2,1) = A(5,3) = 1.1

  cout << A << endl;
  cout << B << endl;
  cout << C << endl;

  return 0;
}

Correct output:

0  0  0  0  0  0  0  0  0  0
7.94687e-311  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0
0  0  0  1.1  0  0  0  3.1  0  0
0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0
0  0  0  0  0  0  0  0  0  0

// Note the following matrix is 5x5
0  0  0  0  0
0  0  0  0  0
0  1.1  0  3.1  0
0  0  0  0  0
0  0  0  0  0

1.1


Incorrect output:

1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1

// The following matrix is 9x9 because the
// view (II=(1,9,2),II) is incorrectly copied to
// (I=(1,9,1),I) when the increment is set to 1
// irrespective of the LaIndex being copied.

// Also, we have assigned to the entire matrix,
// rather than simply the view specified in
//   C.ref(B(LaIndex(2,2,4), JJ));
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1
1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1  1.1

1.1