Auto-suggest helps you quickly narrow down your search results by suggesting possible matches as you type.

Showing results for

- Intel Community
- Software
- Software Development SDKs and Libraries
- Intel® oneAPI Math Kernel Library
- LAPACKE_zgesvd has a bug

- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Mute
- Printer Friendly Page

Chunyu_W_

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-06-2014
01:59 PM

76 Views

LAPACKE_zgesvd has a bug

I think LAPACKE_zgesvd has a bug.

When I use this function to calculate the SVD of a 2*1 matrix [1;3], the result seems to be wrong.

Thanks.

1 Solution

mecej4

Black Belt

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-08-2014
01:36 AM

76 Views

Link Copied

9 Replies

mecej4

Black Belt

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-06-2014
02:41 PM

76 Views

Please show us the result, and tell us why it "seems to be wrong".

I tried the same 2 X 1 matrix example with MKL 11.1 and 11.2, and the results agreed with those from Matlab.

Given that your matrix is real, why do you call the complex SVD routine?

Chunyu_W_

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-07-2014
07:25 AM

76 Views

Thank you for your reply!

Do you use LAPACKE_zgesvd to calculate?

So you mean with different data type, we must use corresponding function. I can not call complex SVD routine when I use real matrix.

Oh, I think that's the problem. I always treat the real number as a special case of complex number and I just let the imaginary part to be 0 to call the complex routine.

mecej4

Black Belt

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-07-2014
11:29 AM

76 Views

Chunyu W. wrote:

Do you use LAPACKE_zgesvd to calculate?

Yes. I used a modified version of the file lapacke_zgesvd_row.c that is distributed in the examples subdirectory of MKL.

So you mean with different data type, we must use corresponding function. I can not call complex SVD routine when I use real matrix.

You should know the properties of the three matrices that the SVD produces. If the original matrix is real, the U, S and V^{T} matrices are all real. When the routines to handle this case (sgesvd and dgesvd) are provided, why would you resort to using the more general routine for complex matrices? All the imaginary parts in the results will be zero, and the calculation will be much slower, without any purpose.

Chunyu_W_

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-07-2014
04:19 PM

76 Views

I still think that complex SVD routine can calculate real matrix and it should be compatible with real matrix.

I use LAPACKE_zgesvd because sometimes I get complex matrix and sometimes I get real matrix, so I call this function overall, you know it is troublesome to judge what matrix I get.

I think the problem might have something to do with the superb parameter.

Below is an explanation about superb from the documentation(https://software.intel.com/sites/products/documentation/hpc/mkl/mklman/)

On exit, `superb(0:min( m,n)-2)` contains the unconverged superdiagonal elements of an upper bidiagonal matrix

So my question is that when min(m,n) = 1, what should be the size of superb? This is just the case when I caculate [1;3].

mecej4

Black Belt

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-07-2014
06:25 PM

76 Views

Chunyu W. wrote:

I think the problem might have something to do with the superb parameter.

Since you have not yet stated what exactly the problem is, I cannot comment on that. The matrix B is of the same size as S. However, S is diagonal whereas B is upper bi-diagonal.

In your example, S is 1 X 1, so the size of the first upper diagonal is zero. However, C does not allow you to declare an array of size zero, so you could declare "double superb[max(1,min(M,N)-1)]. Alternatively, pass NULL as the last argument to LAPACKE_zgesvd.

Chunyu_W_

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-07-2014
07:40 PM

76 Views

Thank you for your patience!

I might paste my code here as well.

/**********************code*******************************/

void svd(int m, int n, FemasDComplex * a, FemasDComplex * u, FemasDComplex * vt, double * s, double * superb)

{

(void)LAPACKE_zgesvd(LAPACK_ROW_MAJOR, 'A', 'A', m, n, (MKL_Complex16 *)a, m, s, (MKL_Complex16 *)u, m, (MKL_Complex16 *)vt, n, superb);

}

int _tmain(int argc, _TCHAR* argv[])

{

FemasDComplex a[2], u[4], vt[1];

double s[1];

a[0].real = 1;

a[0].imag = 0;

a[1].real = 3;

a[1].imag = 0;

for(int i=0; i<2; i++)

cout<<a*.real<<'+'<<a .imag<<"i\n";*

cout<<endl;

svd(2, 1, a, u, vt, s, NULL);

//print the result

cout<<s[0]<<endl;

for(int i=0; i<4; i++)

cout<<u*.real<<'+'<<u .imag<<"i\n";*

cout<<endl;

for(int i=0; i<1; i++)

cout<<vt.real<<'+'<<vt.imag<<"i\n";

cout<<endl;

while(1);

return 0;

}

/****************************code*******************************************/

above is c++ result

above is matlab result.

Chunyu_W_

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-07-2014
07:45 PM

76 Views

FemasDComplex is a structure the same as MKL_Complex16.

typedef struct FemasDComplex {

double real;

double imag;

} FemasDComplex;

I test some other matrix, they give the same result as matlab, just this one show a different result. [1 ; 3]

mecej4

Black Belt

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-08-2014
01:36 AM

77 Views

Chunyu_W_

Beginner

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Email to a Friend
- Report Inappropriate Content

09-08-2014
04:17 PM

76 Views

Thank you. That's exactly the reason.

I also notice that svd is not unique actually, so it is sometimes possible that the result may not be the same as matlab.

- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Printer Friendly Page

For more complete information about compiler optimizations, see our Optimization Notice.