Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.

Pivoting procedure in LU Factorization

wugaxp
Beginner
353 Views
These days I'm working on the partial LU factorization of some matrix. Then I findsome interesting things.
When we use the usual Gaussian elimination to factorize the matrix, we will compare the diagonal elements and the off-diagonal elements below them to determine the permutation matrix. if the matrix is a complex matrix, we must compare the absolute value of the elements. I find that there are two kinds of methods to calculate the "absolute value" in routine zgetrf. The first one is the normal definition of the absolute value of complex number, sqrt(a*a+b*b), if the complex number is z=a+b*i. However, in MKL and the netlib routine,another definition of theabsolute value of complex number is adopted, namely abs(a)+abs(b). The two kinds of normal absolute valuesare calculated with routines CDABS and DCABS1, respectively. So is there any reason besides of the performance consideration here? It can be easily proven that sqrt(a*a+b*b)I find some researchers have compared the Pivoting Strategies in LU Factorization, http://archives.math.utk.edu/ICTCM/i/08/C006.html, and they found that none of the five strategies is absolutely more accurate than the others. Then is it the same here for different definition of absolute value of complex number?
0 Kudos
2 Replies
Michael_C_Intel4
Employee
353 Views

MKL and NetlibLU both use the normal definition of the absolute value of thecomplex number, namely abs(z)=sqrt(a*a+b*b). Please refer to cgetf2/zgetf2 functions at www.netlib.org/lapack.

Michael.
0 Kudos
wugaxp
Beginner
353 Views

MKL and NetlibLU both use the normal definition of the absolute value of thecomplex number, namely abs(z)=sqrt(a*a+b*b). Please refer to cgetf2/zgetf2 functions at www.netlib.org/lapack.

Michael.

Are you sure?
Please check:
http://www.netlib.org/lapack/explore-html/zgetf2.f.html

In the routine, function 'IZAMAX' is used to find the max complex number (See Line 110).
And in http://www.netlib.org/lapack/explore-html/izamax.f.html, you will find that only DCABS1 is used.

In fact, I also write a small routine to verify my guess. In some cases, the DCABS1 and CDABS gives different Pivoting vector. My numerical result indicates the DCABS1 routine is used.
0 Kudos
Reply