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

permutation vector from getrf when matrix is rank deficient

may_ka
Beginner
191 Views

Hi,

 

When using "getrf" for factorizing squared symmetric matrices I have come across permutation vectors as output which contain repeated values and may not cover all columns/rows of the input matrix. This seems to happen only when the input matrix is rank deficient.

 

It appears as if a permutation matrix P usable in A=PLU cannot be obtained directly from such vector.

 

I could recover a usable permutation vector (i.e. a vector from which a permutation matrix can be build) such that A=PLU after studying the "asPerm" C code in "Matrix" package of "R", which makes use of the the system Lapack installation.

 

A C++ example rebuild from the would be:

#include <array>
int main(){
  //p is a potential output from getrf for a 10x10 matrix
  auto p=std::array<int,10>{2,9,7,9,7,1,5,5,6,9};
  //ip will eventually contain the permutation vector covering columns 0 to 9
  auto ip=std::array<int,10>{0,1,2,3,4,5,6,7,8,9};
  //r is the permutation vector sought
  auto r=std::array<int,10>{2,5,7,1,0,4,8,6,9,3};
  //build permutation vector
  for(int j=0;j<ip.size();++j){
    if(p[j]!=j) std::swap(ip[p[j]],ip[j]);
  }
  //check whether results againts vector sought
  for(int j=0;j<ip.size();++j){
    if(ip[j]!=r[j]) return(1);
  }
  return(0);
}

 

What remains unclear to me is the purpose the initial representation of the permutation using repeated numbers.

Google didn't give much results.

s there any documentation/explanation?

And is this representation limited to "getrf"?

 

Thank you.

 

0 Kudos
1 Reply
Gennady_F_Intel
Moderator
97 Views




Here are the lapack's developer comments:

"I think the user is erroneously expecting the IPIV array to be an explicit permutation of (1,…,n) where repetition of index would be an error. The rank-deficient comment is unrelated to the format LAPACK uses.


  • LAPACK uses the IPIV array to encode a permutation as a product of transpositions, such a representation naturally arises out of gaussian elimination with pivoting. Assuming Fortran style indexing, the permutation can be written in cycle notation (composition going from left to right) as
  • (1, ipiv[1]) (2, ipiv[2]) ... (n, ipiv[n])
  • Our documentation for getrf has a terse explanation for IPIV:
  • ipiv

Array, size at least max(1,min(m,n)). Contains the pivot indices; for 1 ≤ i ≤ min(m,n), row i was interchanged with row ipiv(i).

 

  • This representation is common in LAPACK, though we do not have any exhaustive list.

"



0 Kudos
Reply