ScaLAPACK Archives

[Scalapack] matrix inverse


Hello Desheng,

yes there seems to have a problem with PxGETRI with rectangular grids when 
one uses workspace query to get the optimal size of IWORK. The value 
returned in IWORK(1) after a workspace query (LWORK=-1) is too small. The 
documentation gives the same workspace size as the one returned, so the 
documentation of the routine seems to be incorrect as well.

Here is a temporary fix. It will need more testing to be validated but 
this works fine on one of our cluster. Let me know if this works as well 
for you.

Replace Line 221-222 of pdgetri.f 
-               LIWMIN = NQ + MAX( ICEIL( ICEIL( MP, DESCA( MB_ ) ),
-     $                            LCM / NPROW ), DESCA( NB_ ) )
with the following seven lines
+               LIWMIN = NUMROC( DESCA( M_ ) + DESCA( MB_ ) * NPROW
+     $                  + MOD ( IA - 1, MB_P ), DESCA ( NB_ ), MYCOL,
+     $                  DESCA( CSRC_ ), NPCOL ) + MAX ( DESCA( MB_ )
+     $                  * ICEIL ( ICEIL( NUMROC( DESCA( M_ )
+     $                  + DESCA( MB_ ) * NPROW, DESCA( MB_ ), MYROW,
+     $                  DESCA( RSRC_ ), NPROW ), DESCA( MB_ ) ),
+     $                  LCM / NPROW ), DESCA( NB_ ) )

I have tested matrices of size n=[550], with nb=[13 32] for all the grids 
possible with less than 16 processors ( 1-1, 1-2, 2-1, 1-3, 3-1, 1-4, 2-2, 
4-1, 1-5, .... , 16-1), the routine PDGETRI with this modification works 
fine, whereas the previous one was only working on square matrices if one 
was using the minimum LIWORK. You were right.

So, temporary fix, you compile pdgetri.f with your fortran compiler, link 
with scalapack but put pdgetri.o before the libscalapack.a library and 
this should work.

We'll take a bit of time here to have a second look and certainly release 
a more official patch soon.

Feel free for question, remarks, comments.

Julien



<Prev in Thread] Current Thread [Next in Thread>


For additional information you may use the LAPACK/ScaLAPACK Forum.
Or one of the mailing lists, or