Community
cancel
Showing results for 
Search instead for 
Did you mean: 
charlesqhappy
Novice
449 Views

Why the function “LAPACKE_spbtrf” of Intel MKL library keeps returning 1?

Jump to solution

guys, I'm using the function "LAPACKE_spbtrs" of MKL library to solve symmetric-positive-definite-band matrix equations. The 1st step is factorizing the matrix by "LAPACKE_spbtrf" according to the office website: https://software.intel.com/content/www/us/en/develop/documentation/mkl-developer-reference-c/top/lap....

But when I run the "LAPACKE_spbtrf", it always return 1, which indicates that the input matrix is not positive-definite. But I'm pretty sure it is.

Here is the full code:

/* Parameters */

#define N 5

#define NRHS 1

#define LDA N

#define LDB NRHS

/* Main program */ int main() {

//convert a symmetric band matrix A to lower band storage mode LA
/*A =  1 1
       1 3
*/

float ab[4] = { 0,1,
                1,3 };
const int n = 2;
int kd = 1;
int lda = kd + 1;
int info = LAPACKE_spbtrf(LAPACK_ROW_MAJOR, 'L', n, kd, ab, lda);
std::cout << "info:" << info << std::endl;

}

What am I missing?

0 Kudos
1 Solution
mecej4
Black Belt
433 Views

If you wanted to provide the lower triangle of the symmetric matrix, following the conventions at Netlib , you should set the values into the matrix using

 

float ab[4] = { 1,3,
                1,0 };

 

That is, the main diagonal should be in the first row, followed by the sub-diagonal in the second row. The zero is just a placeholder, since the subdiagonal has only one element, and that zero corresponds to the '*' in the diagram at the link that I listed.

View solution in original post

5 Replies
charlesqhappy
Novice
441 Views

Intersting, when I set the lower mode  "L" to upper mode "U", the function will return 0, which indicate sucess. But I suppose I do store the lower part of the original matrix, right?

mecej4
Black Belt
434 Views

If you wanted to provide the lower triangle of the symmetric matrix, following the conventions at Netlib , you should set the values into the matrix using

 

float ab[4] = { 1,3,
                1,0 };

 

That is, the main diagonal should be in the first row, followed by the sub-diagonal in the second row. The zero is just a placeholder, since the subdiagonal has only one element, and that zero corresponds to the '*' in the diagram at the link that I listed.

View solution in original post

charlesqhappy
Novice
427 Views

Hey Mecej, thanks very very much, it works! This problem has annoyed me 3 days. But I still have a confusion. It is that why I should use the band storage mode of Netlib instead of the MKL's?  After all, I'm calling the function of MKL. And the MKl website gives a different storage mode of the band matrix, of which the diagonal elements are in a column :(I use row major layout as you can see in my funcion "LAPACKE_spbtrf")

微信截图_20200926130741.png

But in Netlib, they are in a row:

n.png

mecej4
Black Belt
414 Views

Let me answer in a broad way since, if I try to analyze what went wrong with your original attempt, it is easy to get confused by C versus Fortran, row versus column versus band storage, and the version of MKL that you use versus the version of the MKL documentation that I refer to, etc.

Lapack is open source, and is maintained and contributed to by a large number of people from around the world. That fact can tend to increase inconsistencies between the software and the documentation, but there is a counter-influence: people worldwide use Lapack, and discrepancies get reported and fixed quickly.

MKL is different. It is closed source, and belongs to Intel. Lapack is a small subset of what MKL covers, and the MKL documentation has to cover much more than Lapack. Secondly, C users will probably use the LapackE interface to Lapack, and LapackE is somewhat of a step-child.

In the specific case of SPBTRF, the MKL  documentation and the Netlib documentation differ regarding the argument AB. In such cases, my inclination is to follow the Netlib documentation and see if that works. You also made LDAB = KD+1, which makes the 1-D array that you passed for AB completely equivalent to the packed 2-D array that should have been used.

During the last few years, it did not help that at times the MKL documentation for C and Fortran usage was merged into a single document. The result was, at least for the BLAS and Lapack parts, that usually half the monolingual (i.e., C/Fortran) readers were left confused. This problem persists in some places in the current MKL documentation, even though the C and Fortran versions have been "separated". One could say, "they have taken the C out of the Fortran, but they failed to take the Fortran out of the C". There are parts of MKL where the opposite is true.

charlesqhappy
Novice
400 Views

Hi , Mecej, thanks! But I got the same problem with another matrix using the Netlib band matrix storage mode. Please help me with it. The link is :

https://community.intel.com/t5/Intel-oneAPI-Math-Kernel-Library/Why-the-quot-LAPACKE-spbtrf-quot-alw...

Reply