topic Re: false results with zgesvd in Intel® oneAPI Math Kernel Library

false results with zgesvd

mouna_hammami — Tue, 22 Nov 2005 16:55:56 GMT

I try to compute the SVD of complex matrices using zgesvd but when I compare the results with those of matlab, I found a big difference!

this is my code written in C language:

#include

#include

#include "include/mkl.h"

#define SIZE 3

void ZGESVD(char *,char *,int *,int *, MKL_Complex16 *,int *,double *, MKL_Complex16 *,int *, MKL_Complex16 *,int *, MKL_Complex16 *,int *,double *,int *);

int main (){

char JOBU;

char JOBVT;

int M = SIZE;

int N = SIZE;

int LDA = M;

int LDU = M;

int LDVT = N;

int LWORK;

double RWORK;

int INFO;

int i,j;

MKL_Complex16 *a ;

a = (MKL_Complex16 *)calloc( M*N, sizeof(MKL_Complex16) );

if ( a == NULL ) {

printf(" Can't allocate memory arrays");

return 0;

}

for( i = 0; i M; i++ ){

for (j=0 ; j

{

a[N*(i-1)+j].real = (double)i;

a[N*(i-1)+j].imag = (double)i * 2.0;

}

}

double s[SIZE];

MKL_Complex16 wk[201];

MKL_Complex16 uu[SIZE*SIZE];

MKL_Complex16 vt[SIZE*SIZE];

JOBU = 'A';

JOBVT = 'A';

LWORK = 201;

ZGESVD(&JOBU, &JOBVT, &M, &N, a, &LDA, s, uu, &LDU, vt, &LDVT,wk,&LWORK, &RWORK, &INFO );

printf(" INFO=%d", INFO );

for ( i= 0; i SIZE; i++ ) {

printf(" s[ %d ] = %f", i, s[ i ] ); }

free(a);

return 1;

}

what's the problem? any help will be welcome :) thanx!

Re: false results with zgesvd

Gregory_H_Intel — Thu, 01 Dec 2005 07:00:30 GMT

What does your MATLAB code look like? Is it something like this:

A =

1.0000 + 2.0000i 1.0000 + 2.0000i 1.0000 + 2.0000i
2.0000 + 4.0000i 2.0000 + 4.0000i 2.0000 + 4.0000i
3.0000 + 6.0000i 3.0000 + 6.0000i 3.0000 + 6.0000i

>> svd(A)

ans =

14.4914
0.0000
0.0000

Also, can you do an attachment of your MATLAB code (and/or a .mat)and an attachment of your C code? Some of the characters in your message seemed to have gotten lost, whereas perhaps with an attachment it might not get lost.

Thanks,

Greg Henry (Intel MKL)

Re: false results with zgesvd

meretikoesra — Wed, 29 Apr 2009 13:26:32 GMT

Quoting - Gregory Henry (Intel)

What does your MATLAB code look like? Is it something like this:

A =

1.0000 + 2.0000i 1.0000 + 2.0000i 1.0000 + 2.0000i
2.0000 + 4.0000i 2.0000 + 4.0000i 2.0000 + 4.0000i
3.0000 + 6.0000i 3.0000 + 6.0000i 3.0000 + 6.0000i

>> svd(A)

ans =

14.4914
0.0000
0.0000

Also, can you do an attachment of your MATLAB code (and/or a .mat)and an attachment of your C code? Some of the characters in your message seemed to have gotten lost, whereas perhaps with an attachment it might not get lost.

Thanks,

Greg Henry (Intel MKL)

hi i have a similar problem

here is my source code

int main (int argc, char *argv[]) {
MKL_Complex16 *a,*u,*vt,*work, *carp, *carp2;
MKL_INT nm = 5;
int lw = 201;
int rw = 201;

u = (MKL_Complex16*)calloc(sizeof(MKL_Complex16), nm*nm);
vt = (MKL_Complex16*)calloc(sizeof(MKL_Complex16), nm*nm);
a = (MKL_Complex16*)calloc(sizeof(MKL_Complex16), nm*nm);
work = (MKL_Complex16*)calloc(sizeof(MKL_Complex16), rw);

a[0].real = 0.4499; a[0].imag = -0.1267;
a[1].real = -0.2361; a[1].imag = 0.0775;
a[2].real = 0.5008; a[2].imag = -0.7799;
a[3].real = -0.7702; a[3].imag = -0.3143;
a[4].real = -0.5084; a[4].imag = 0.5621;

a[5].real = -0.8432; a[5].imag = -0.3443;
a[6].real = 0.0186; a[6].imag = -0.6332;
a[7].real = 0.3505; a[7].imag = 0.0166;
a[8].real = -0.8892; a[8].imag = 0.2657;
a[9].real = -0.5066; a[9].imag = 0.6005;

a[10].real = 0.5896; a[10].imag = 0.2601;
a[11].real = 0.8928; a[11].imag = 0.0103;
a[12].real = 0.4350; a[12].imag = -0.2666;
a[13].real = -0.7799; a[13].imag = -0.5512;
a[14].real = -0.1022; a[14].imag = -0.4500;

a[15].real = 0.3911; a[15].imag = 0.3234;
a[16].real = -0.5605; a[16].imag = 0.7638;
a[17].real = -0.2009; a[17].imag = 0.2707;
a[18].real = 0.0314; a[18].imag = -0.4336;
a[19].real = -0.1046; a[19].imag = 0.3267;

a[20].real = -0.1395; a[20].imag = -0.1561;
a[21].real = 0.3064; a[21].imag = 0.0263;
a[22].real = -0.0975; a[22].imag = 0.9019;
a[23].real = 0.2491; a[23].imag = 0.5787;
a[24].real = 0.4106; a[24].imag = -0.8102;

double *s;
s = (double*)calloc(sizeof(double), nm);
int info, ldu = nm, ldvt = nm, lwork = lw;
double *rwork = (double*)calloc(rw,sizeof(double));
zgesvd("A","A",&nm,&nm,a,&nm,s,u,&ldu,vt,&ldvt,work,&lwork,rwork,&info);

return 0;
}

results are:

S vector: 2.4 1.9 1.2 0.71 0.25

U matrix:

-0.11 0.25i -0.17 0.27i -0.25 0.38i 0.61 0.14i 0.42 -0.23i

0.23 0.37i -0.43 0.44i 0.0014 -0.33i 0.015 -0.36i 0.024 0.44i

-0.5 -0.33i -0.074 0.31i 0.026 0.54i -0.19 -0.24i -0.094 0.38i

-0.41 0.25i 0.16 -0.34i -0.49 -0.17i -0.11 -0.55i 0.17 -0.088i

-0.18 -0.33i -0.52 0.091i 0.25 -0.23i -0.25 -0.13i 0.36 -0.51i

VT matrix:

-0.53 0i 0.41 0i -0.064 0i -0.38 -0i 0.63 0i

-0.45 0.26i -0.24 0.28i 0.63 -0.33i 0.16 -0.21i -0.068 -0.12i

-0.35 -0.39i -0.038 -0.34i -0.4 -0.4i 0.29 -0.28i -0.14 -0.32i

0.14 -0.11i 0.51 -0.081i 0.3 0.26i -0.26 -0.51i -0.34 -0.32i

0.36 -0.065i -0.54 -0.14i 0.055 -0.099i -0.37 -0.4i 0.44 -0.22i

Matlab output:
Singular values (Diagonal of diagonal S matrix):

2.43287685397618
1.94346744841323
1.15707710132509
0.70578696643348
0.24840203289605

U:

0.2242-0.4836 -0.2151-0.3438 0.0537+0.0351 -0.3236+0.2017 0.3038+0.5564
-0.0470-0.5195 0.3622+0.0569 -0.7086-0.0724 0.0222-0.2599 0.0695-0.1153
0.5033-0.1529 -0.2686+0.2132 0.1156+0.5523 0.0944-0.3938 0.2180-0.2752
0.0390+0.1680 -0.3397-0.3904 -0.1128-0.3801 -0.4899-0.2989 0.1194-0.4491
-0.0943+0.3583 0.1731+0.5345 -0.1006+0.0530 -0.5220-0.1439 0.4013+0.2825

V:

0.2731+0.0000 -0.4362+0.0000 0.5986+0.0000 -0.4818+0.0000 -0.3803+0.0000
0.3175-0.0394 -0.1459+0.5963 -0.1079+0.3004 0.3115+0.2061 -0.1693-0.5007
0.4497-0.0720 0.2791-0.1763 -0.3194+0.4403 -0.3277+0.4089 -0.0848+0.3256
-0.1286+0.6084 0.2932-0.2010 0.2858-0.0985 0.2001+0.5267 -0.2325-0.1550
-0.3843+0.2852 -0.3888+0.2152 -0.1289+0.3687 0.1898-0.0146 -0.2733+0.5568

thanks in advance..
esra

Re: false results with zgesvd

meretikoesra — Thu, 30 Apr 2009 13:33:55 GMT

Quoting - meretikoesra

ok problem solved..
U transpose is returned from c zgesvd function..
a = u' * s * v' gives the original a matrix.