lu for complex<double> bug?
Mateusz Koreń <[email protected]> Mon, 08 May 2006 11:38:55 +0200
| Newsgroups | gmane.comp.lib.mtl.devel |
|---|---|
| Organization | Nowa Poczta Wirtualnej Polski S.A. http://www.wp.pl/ |
| Message-ID | <[email protected]> |
I'm trying to use mtl as a inverter for my tvmet matrices and i found that something gets wrong when matrix elements have imaginary parts.
I checked how it works with lu_factorization.cc from your website and it seems to give something like this (sorry for an inelegant code):
#include <iostream>
#include <complex>
using namespace std;
#include <mtl/lu.h>
#include <mtl/matrix.h>
#include <mtl/dense1D.h>
#include <mtl/utils.h>
#include <mtl/blais.h>
int main() {
using namespace mtl;
typedef complex<double> comp;
typedef matrix<comp,
rectangle<>,
dense<external>,
column_major>::type Matrix;
const int N = 3;
comp da [] = {3,comp(1,1),0,0,3,0,comp(2,0),0,3};
comp da2 [] = {3,comp(1,1),0,0,3,0,comp(2,0),0,3};
comp db [] = {3,comp(1,1),0,0,3,0,comp(2,2),0,3};
comp db2 [] = {3,comp(1,1),0,0,3,0,comp(2,2),0,3};
comp resulta[9];
comp resultb[9];
comp onea[9];
comp oneb[9];
Matrix A2(da2,N,N);
Matrix A(da,N,N);
Matrix B2(db2,N,N);
Matrix B(db,N,N);
Matrix A_inv(resulta,N,N);
Matrix B_inv(resultb,N,N);
Matrix identa(onea,N,N);
Matrix identb(oneb,N,N);
dense1D<int> pivotsa(N,0);
dense1D<int> pivotsb(N,0);
lu_factor(A, pivotsa);
lu_inverse(A,pivotsa,A_inv);
lu_factor(B, pivotsb);
lu_inverse(B,pivotsb,B_inv);
mult(A2,A_inv,identa);
mult(B2,B_inv,identb);
cout<<"A*A^-1 (should be identity):"<<endl;
print_all_matrix(identa);
cout<<"B*B^-1 (should be identity):"<<endl;
print_all_matrix(identb);
return 0;
}
resulting in:
A*A^-1 (should be identity):
3x3
[
[(1,0),(0,0),(0,0)],
[(0,0),(1,0),(0,0)],
[(0,0),(0,0),(1,0)]
]
B*B^-1 (should be identity):
3x3
[
[(1,0),(0,0),(0,0)],
[(0,0),(1,0),(0.444444,-0.444444)],
[(0,0),(0,0),(1,0)]
]
and so on with more "imaginary" matrices. Seems like a bug to me. Or maybe I'm doing something wrong? I've got Debian/GNU Linux 3.1 with gcc 4.0
Best regards,
Mateusz
PS I also checked the internal storage way - with the same result.
----------------------------------------------------
Koncert zespołu TOOL!
24 czerwca w katowickim Spodku!
http://klik.wp.pl/?adr=http%3A%2F%2Fadv.reklama.wp.pl%2Fas%2Ftool.html&sid=742
_______________________________________________
This list is archived at http://www.osl.iu.edu/MailArchives/mtl-devel/