#18 Equalling a matrix to a submatrix does not give a new matrix

[Anonymous]

Equalling a just declared matrix to a submatrix does not give a newly allocated matrix (modifying the new matrix modifies the old matrix). I mean:

// it just creates a random (gaussian) matrix
LaGenMatDouble A = StatUtil::RandnMatrix(2,3,0.0,1.0);

LaGenMatDouble B = A(LaIndex(0,1),LaIndex(1,2));

cout << "A is" << endl << A << "B is" << endl << B;
B(0,0) = 10;
cout << "A is" << endl << A << "B is" << endl << B;

A is
-0.531466 -1.82688 0.302208
0.849202 -0.983648 -0.00441152
B is
-1.82688 0.302208
-0.983648 -0.00441152
A is
-0.531466 10 0.302208
0.849202 -0.983648 -0.00441152
B is
10 0.302208
-0.983648 -0.00441152

I don't know if that's the behavior you were looking for, but it's not what I expected.

Regards,
Manuel A. Vázquez (mani001@terra.es)

[Christian Stimming]

Logged In: YES
user_id=531034
Originator: NO

Thanks for reporting this issue. As you can see in the example program below, the unwanted reference copy is created only when writing LaGenMatDouble B = A(LaIndex,LaIndex); i.e. only when assigning directly at declaration. It does not appear when assigning later, C=A(LaIndex,LaIndex), where a deep copy as expected is created. I'm actually a bit at a loss as to why this happens - the constructor and operator= of that class clearly will create a deep copy . Maybe I can think of a solution, but so far I don't have an idea.

For the record, providing full example programs for such bugs would make bugfixing much easier.

// Compile as follows: gcc bla.cpp -o bla -I/usr/include/lapackpp -llapackpp

include "gmd.h"

using namespace std;
int main(int argc, char* argv[]){
// Just an example matrix
LaGenMatDouble A(2,3);
A(0,0)=0; A(0,1)=1; A(0,2)=2;
A(1,0)=3; A(1,1)=4; A(1,2)=5;
A.debug(1);

LaGenMatDouble B = A(LaIndex(0,1),LaIndex(1,2)); // huh? Is a reference to A

LaGenMatDouble C;
C = A(LaIndex(0,1),LaIndex(1,2)); // correctly creates a copy of A

cout << "A is" << endl << A << "B is" << endl << B << "C is" << endl << C;
B(0,0) = 10;
cout << "A is" << endl << A << "B is" << endl << B << "C is" << endl << C;
return 0;
}

[Christian Stimming]

Logged In: YES
user_id=531034
Originator: NO

Thanks for reporting this issue. It turns out this is caused by an unexpected optimization of C++, and this cannot easily be fixed in the current lapackpp structure. A workaround is explained in the updated documentation here http://lapackpp.sourceforge.net/html/classLaGenMatDouble.html#32d24b40c52b38db9f4bccd0ac8df2f5 ; basically, you must write

LaGenMatDouble B;
B = A(LaIndex(0,1),LaIndex(1,2));

All other possibilities will not work with the current lapackpp matrix classes. Sorry for that. This bug will stay open until someone can come up with a permanent solution.