Routine PSPOCON() Documented By Universal Report

Routine: PSPOCON()

File: SRC\pspocon.f

# lines:	405
# code:	405
# comment:	0
# blank:	0
# Variables:	73
# Callers:	1
# Callings:	4
# Words:	223
# Keywords:	140

..
     .. Array Arguments ..
     ..
  Purpose
  =======
  PSPOCON estimates the reciprocal of the condition number (in the
  1-norm) of a real symmetric positive definite distributed matrix
  using the Cholesky factorization A = U**T*U or A = L*L**T computed by
  PSPOTRF.
  An estimate is obtained for norm(inv(A(IA:IA+N-1,JA:JA+N-1))), and
  the reciprocal of the condition number is computed as
             RCOND = 1 / ( norm( A(IA:IA+N-1,JA:JA+N-1)      ) *
                           norm( inv(A(IA:IA+N-1,JA:JA+N-1)) ) ).
  Notes
  =====
  Each global data object is described by an associated description
  vector.  This vector stores the information required to establish
  the mapping between an object element and its corresponding process
  and memory location.
  Let A be a generic term for any 2D block cyclicly distributed array.
  Such a global array has an associated description vector DESCA.
  In the following comments, the character _ should be read as
  "of the global array".
  NOTATION        STORED IN      EXPLANATION
  --------------- -------------- --------------------------------------
  DTYPE_A(global) DESCA( DTYPE_ )The descriptor type.  In this case,
                                 DTYPE_A = 1.
  CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
                                 the BLACS process grid A is distribu-
                                 ted over. The context itself is glo-
                                 bal, but the handle (the integer
                                 value) may vary.
  M_A    (global) DESCA( M_ )    The number of rows in the global
                                 array A.
  N_A    (global) DESCA( N_ )    The number of columns in the global
                                 array A.
  MB_A   (global) DESCA( MB_ )   The blocking factor used to distribute
                                 the rows of the array.
  NB_A   (global) DESCA( NB_ )   The blocking factor used to distribute
                                 the columns of the array.
  RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
                                 row of the array A is distributed.
  CSRC_A (global) DESCA( CSRC_ ) The process column over which the
                                 first column of the array A is
                                 distributed.
  LLD_A  (local)  DESCA( LLD_ )  The leading dimension of the local
                                 array.  LLD_A >= MAX(1,LOCr(M_A)).
  Let K be the number of rows or columns of a distributed matrix,
  and assume that its process grid has dimension p x q.
  LOCr( K ) denotes the number of elements of K that a process
  would receive if K were distributed over the p processes of its
  process column.
  Similarly, LOCc( K ) denotes the number of elements of K that a
  process would receive if K were distributed over the q processes of
  its process row.
  The values of LOCr() and LOCc() may be determined via a call to the
  ScaLAPACK tool function, NUMROC:
          LOCr( M ) = NUMROC( M, MB_A, MYROW, RSRC_A, NPROW ),
          LOCc( N ) = NUMROC( N, NB_A, MYCOL, CSRC_A, NPCOL ).
  An upper bound for these quantities may be computed by:
          LOCr( M ) <= ceil( ceil(M/MB_A)/NPROW )*MB_A
          LOCc( N ) <= ceil( ceil(N/NB_A)/NPCOL )*NB_A
  Arguments
  =========
  UPLO    (global input) CHARACTER
          Specifies whether the factor stored in
          A(IA:IA+N-1,JA:JA+N-1) is upper or lower triangular.
          = 'U':  Upper triangular
          = 'L':  Lower triangular
  N       (global input) INTEGER
          The order of the distributed matrix A(IA:IA+N-1,JA:JA+N-1).
          N >= 0.
  A       (local input) REAL pointer into the local memory to
          an array of dimension ( LLD_A, LOCc(JA+N-1) ). On entry, this
          array contains the local pieces of the factors L or U from
          the Cholesky factorization A(IA:IA+N-1,JA:JA+N-1) = U'*U or
          L*L', as computed by PSPOTRF.
  IA      (global input) INTEGER
          The row index in the global array A indicating the first
          row of sub( A ).
  JA      (global input) INTEGER
          The column index in the global array A indicating the
          first column of sub( A ).
  DESCA   (global and local input) INTEGER array of dimension DLEN_.
          The array descriptor for the distributed matrix A.
  ANORM   (global input) REAL
          The 1-norm (or infinity-norm) of the symmetric distributed
          matrix A(IA:IA+N-1,JA:JA+N-1).
  RCOND   (global output) REAL
          The reciprocal of the condition number of the distributed
          matrix A(IA:IA+N-1,JA:JA+N-1), computed as
             RCOND = 1 / ( norm( A(IA:IA+N-1,JA:JA+N-1)      ) *
                           norm( inv(A(IA:IA+N-1,JA:JA+N-1)) ) ).
  WORK    (local workspace/local output) REAL array,
                                                   dimension (LWORK)
          On exit, WORK(1) returns the minimal and optimal LWORK.
  LWORK   (local or global input) INTEGER
          The dimension of the array WORK.
          LWORK is local input and must be at least
          LWORK >= 2*LOCr(N+MOD(IA-1,MB_A)) + 2*LOCc(N+MOD(JA-1,NB_A))+
          MAX( 2, MAX(NB_A*CEIL(NPROW-1,NPCOL),LOCc(N+MOD(JA-1,NB_A)) +
          NB_A*CEIL(NPCOL-1,NPROW)) ).
          If LWORK = -1, then LWORK is global input and a workspace
          query is assumed; the routine only calculates the minimum
          and optimal size for all work arrays. Each of these
          values is returned in the first entry of the corresponding
          work array, and no error message is issued by PXERBLA.
  IWORK   (local workspace/local output) INTEGER array,
                                                 dimension (LIWORK)
          On exit, IWORK(1) returns the minimal and optimal LIWORK.
  LIWORK  (local or global input) INTEGER
          The dimension of the array IWORK.
          LIWORK is local input and must be at least
          LIWORK >= LOCr(N+MOD(IA-1,MB_A)).
          If LIWORK = -1, then LIWORK is global input and a workspace
          query is assumed; the routine only calculates the minimum
          and optimal size for all work arrays. Each of these
          values is returned in the first entry of the corresponding
          work array, and no error message is issued by PXERBLA.
  INFO    (global output) INTEGER
          = 0:  successful exit
          < 0:  If the i-th argument is an array and the j-entry had
                an illegal value, then INFO = -(i*100+j), if the i-th
                argument is a scalar and had an illegal value, then
                INFO = -i.
  =====================================================================
     .. Parameters ..

Display dynamic version

FontSize: - +

001        SUBROUTINE PSPOCON( UPLO , N , A , IA , JA , DESCA , ANORM , RCOND , WORK ,
002       $LWORK , IWORK , LIWORK , INFO )
003  
004  *     -- ScaLAPACK routine(version 1.7) --
005  *     University of Tennessee , Knoxville , Oak Ridge National Laboratory ,
006  *     and University of California , Berkeley.
007  *     May 25 , 2001
008  
009  *     .. Scalar Arguments ..
010        CHARACTER UPLO
011        INTEGER IA , INFO , JA , LIWORK , LWORK , N
012        REAL ANORM , RCOND
013        INTEGER BLOCK_CYCLIC_2D , CSRC_ , CTXT_ , DLEN_ , DTYPE_ ,
014       $LLD_ , MB_ , M_ , NB_ , N_ , RSRC_
015        PARAMETER( BLOCK_CYCLIC_2D = 1 , DLEN_ = 9 , DTYPE_ = 1 ,
016       $CTXT_ = 2 , M_ = 3 , N_ = 4 , MB_ = 5 , NB_ = 6 ,
017       $RSRC_ = 7 , CSRC_ = 8 , LLD_ = 9 )
018        REAL ONE , ZERO
019        PARAMETER( ONE = 1.0E + 0 , ZERO = 0.0E + 0 )
020  *     ..
021  *     .. Local Scalars ..
022        LOGICAL LQUERY , UPPER
023        CHARACTER CBTOP , COLCTOP , NORMIN , ROWCTOP
024        INTEGER IACOL , IAROW , ICOFF , ICTXT , IIA , IPNL , IPNU ,
025       $IPV , IPW , IPX , IROFF , IV , IX , IXX , JJA , JV ,
026       $JX , KASE , LIWMIN , LWMIN , MYCOL , MYROW , NP ,
027       $NPCOL , NPROW , NPMOD , NQ , NQMOD
028        REAL AINVNM , SCALE , SL , SU , SMLNUM
029        REAL WMAX
030  *     ..
031  *     .. Local Arrays ..
032        INTEGER DESCV( DLEN_ ) , DESCX( DLEN_ ) , IDUM1( 3 ) ,
033       $IDUM2( 3 )
034  *     ..
035  *     .. External Subroutines ..
036        EXTERNAL BLACS_GRIDINFO , CHK1MAT , DESCSET , INFOG2L ,
037       $PCHK1MAT , PSAMAX , PSLATRS , PSLACON ,
038       $PSRSCL , PB_TOPGET , PB_TOPSET , PXERBLA , SGEBR2D ,
039       $SGEBS2D
040  *     ..
041  *     .. External Functions ..
042        LOGICAL LSAME
043        INTEGER ICEIL , INDXG2P , NUMROC
044        REAL PSLAMCH
045        EXTERNAL ICEIL , INDXG2P , LSAME , NUMROC , PSLAMCH
046  *     ..
047  *     .. Intrinsic Functions ..
048        INTRINSIC ABS , ICHAR , MAX , MOD , REAL
049  *     ..
050  *     .. Executable Statements ..
051  
052  *     Get grid parameters
053  
054        ICTXT = DESCA( CTXT_ )
055        CALL BLACS_GRIDINFO( ICTXT , NPROW , NPCOL , MYROW , MYCOL )
056  
057  *     Test the input parameters
058  
059        INFO = 0
060        IF( NPROW.EQ. - 1 ) THEN
060               
061            INFO = - (600 + CTXT_)
062        ELSE
062               
063            CALL CHK1MAT( N , 2 , N , 2 , IA , JA , DESCA , 6 , INFO )
064            IF( INFO.EQ.0 ) THEN
064                   
065                UPPER = LSAME( UPLO , 'U' )
066                IAROW = INDXG2P( IA , DESCA( MB_ ) , MYROW , DESCA( RSRC_ ) ,
067       $        NPROW )
068                IACOL = INDXG2P( JA , DESCA( NB_ ) , MYCOL , DESCA( CSRC_ ) ,
069       $        NPCOL )
070                NPMOD = NUMROC( N + MOD( IA - 1 , DESCA( MB_ ) ) , DESCA( MB_ ) ,
071       $        MYROW , IAROW , NPROW )
072                NQMOD = NUMROC( N + MOD( JA - 1 , DESCA( NB_ ) ) , DESCA( NB_ ) ,
073       $        MYCOL , IACOL , NPCOL )
074                LWMIN = 2*NPMOD + 2*NQMOD +
075       $        MAX( 2 , MAX( DESCA( NB_ )*
076       $        MAX( 1 , ICEIL( NPROW - 1 , NPCOL ) ) , NQMOD +
077       $        DESCA( NB_ )*
078       $        MAX( 1 , ICEIL( NPCOL - 1 , NPROW ) ) ) )
079                WORK( 1 ) = REAL( LWMIN )
080                LIWMIN = NPMOD
081                IWORK( 1 ) = LIWMIN
082                LQUERY =( LWORK.EQ. - 1 .OR. LIWORK.EQ. - 1 )
083  
084                IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO , 'L' ) ) THEN
084                       
085                    INFO = - 1
086                ELSE IF( ANORM.LT.ZERO ) THEN
086                       
087                    INFO = - 7
088                ELSE IF( LWORK.LT.LWMIN .AND. .NOT.LQUERY ) THEN
088                       
089                    INFO = - 10
090                ELSE IF( LIWORK.LT.LIWMIN .AND. .NOT.LQUERY ) THEN
090                       
091                    IWORK( 1 ) = LIWMIN
092                    INFO = - 12
093                END IF
094            END IF
095  
096            IF( UPPER ) THEN
096                   
097                IDUM1( 1 ) = ICHAR( 'U' )
098            ELSE
098                   
099                IDUM1( 1 ) = ICHAR( 'L' )
100            END IF
101            IDUM2( 1 ) = 1
102            IF( LWORK.EQ. - 1 ) THEN
102                   
103                IDUM1( 2 ) = - 1
104            ELSE
104                   
105                IDUM1( 2 ) = 1
106            END IF
107            IDUM2( 2 ) = 10
108            IF( LIWORK.EQ. - 1 ) THEN
108                   
109                IDUM1( 3 ) = - 1
110            ELSE
110                   
111                IDUM1( 3 ) = 1
112            END IF
113            IDUM2( 3 ) = 12
114            CALL PCHK1MAT( N , 2 , N , 2 , IA , JA , DESCA , 6 , 3 , IDUM1 , IDUM2 ,
115       $    INFO )
116        END IF
117  
118        IF( INFO.NE.0 ) THEN
118               
119            CALL PXERBLA( ICTXT , 'PSPOCON' , - INFO )
120            RETURN
121        ELSE IF( LQUERY ) THEN
121               
122            RETURN
123        END IF
124  
125  *     Quick return if possible
126  
127        RCOND = ZERO
128        IF( N.EQ.0 ) THEN
128               
129            RCOND = ONE
130            RETURN
131        ELSE IF( ANORM.EQ.ZERO ) THEN
131               
132            RETURN
133        ELSE IF( N.EQ.1 ) THEN
133               
134            RCOND = ONE
135            RETURN
136        END IF
137  
138        CALL PB_TOPGET( ICTXT , 'Combine' , 'Columnwise' , COLCTOP )
139        CALL PB_TOPGET( ICTXT , 'Combine' , 'Rowwise' , ROWCTOP )
140        CALL PB_TOPSET( ICTXT , 'Combine' , 'Columnwise' , '1 - tree' )
141        CALL PB_TOPSET( ICTXT , 'Combine' , 'Rowwise' , '1 - tree' )
142  
143        SMLNUM = PSLAMCH( ICTXT , 'Safe minimum' )
144        IROFF = MOD( IA - 1 , DESCA( MB_ ) )
145        ICOFF = MOD( JA - 1 , DESCA( NB_ ) )
146        CALL INFOG2L( IA , JA , DESCA , NPROW , NPCOL , MYROW , MYCOL , IIA , JJA ,
147       $IAROW , IACOL )
148        NP = NUMROC( N + IROFF , DESCA( MB_ ) , MYROW , IAROW , NPROW )
149        NQ = NUMROC( N + ICOFF , DESCA( NB_ ) , MYCOL , IACOL , NPCOL )
150        IV = IROFF + 1
151        IX = IV
152        JV = ICOFF + 1
153        JX = JV
154  
155        IPX = 1
156        IPV = IPX + NP
157        IPNL = IPV + NP
158        IPNU = IPNL + NQ
159        IPW = IPNU + NQ
160  
161        CALL DESCSET( DESCV , N + IROFF , 1 , DESCA( MB_ ) , 1 , IAROW , MYCOL ,
162       $ICTXT , MAX( 1 , NP ) )
163        CALL DESCSET( DESCX , N + IROFF , 1 , DESCA( MB_ ) , 1 , IAROW , MYCOL ,
164       $ICTXT , MAX( 1 , NP ) )
165  
166  *     Estimate the 1 - norm(or I - norm) of inv(A).
167  
168        AINVNM = ZERO
169        KASE = 0
170        NORMIN = 'N'
171  
172     10 CONTINUE
173        CALL PSLACON ( N , WORK( IPV ) , IV , JV , DESCV , WORK( IPX ) , IX , JX ,
174       $DESCX , IWORK , AINVNM , KASE )
175        IF( KASE.NE.0 ) THEN
175               
176            IF( UPPER ) THEN
177  
178  *             Multiply by inv(U').
179  
179                   
180                DESCX( CSRC_ ) = IACOL
181                CALL PSLATRS ( 'Upper' , 'Transpose' , 'Non - unit' , NORMIN ,
182       $        N , A , IA , JA , DESCA , WORK( IPX ) , IX , JX ,
183       $        DESCX , SL , WORK( IPNL ) , WORK( IPW ) )
184                DESCX( CSRC_ ) = MYCOL
185                NORMIN = 'Y'
186  
187  *             Multiply by inv(U).
188  
189                DESCX( CSRC_ ) = IACOL
190                CALL PSLATRS ( 'Upper' , 'No transpose' , 'Non - unit' , NORMIN ,
191       $        N , A , IA , JA , DESCA , WORK( IPX ) , IX , JX ,
192       $        DESCX , SU , WORK( IPNU ) , WORK( IPW ) )
193                DESCX( CSRC_ ) = MYCOL
194            ELSE
195  
196  *             Multiply by inv(L).
197  
197                   
198                DESCX( CSRC_ ) = IACOL
199                CALL PSLATRS ( 'Lower' , 'No transpose' , 'Non - unit' , NORMIN ,
200       $        N , A , IA , JA , DESCA , WORK( IPX ) , IX , JX ,
201       $        DESCX , SL , WORK( IPNL ) , WORK( IPW ) )
202                DESCX( CSRC_ ) = MYCOL
203                NORMIN = 'Y'
204  
205  *             Multiply by inv(L').
206  
207                DESCX( CSRC_ ) = IACOL
208                CALL PSLATRS ( 'Lower' , 'Transpose' , 'Non - unit' , NORMIN ,
209       $        N , A , IA , JA , DESCA , WORK( IPX ) , IX , JX ,
210       $        DESCX , SU , WORK( IPNU ) , WORK( IPW ) )
211                DESCX( CSRC_ ) = MYCOL
212            END IF
213  
214  *         Multiply by 1 / SCALE if doing so will not cause overflow.
215  
216            SCALE = SL*SU
217            IF( SCALE.NE.ONE ) THEN
217                   
218                CALL PSAMAX( N , WMAX , IXX , WORK( IPX ) , IX , JX , DESCX , 1 )
219                IF( DESCX( M_ ).EQ.1 .AND. N.EQ.1 ) THEN
219                       
220                    CALL PB_TOPGET( ICTXT , 'Broadcast' , 'Columnwise' , CBTOP )
221                    IF( MYROW.EQ.IAROW ) THEN
221                           
222                        CALL SGEBS2D( ICTXT , 'Column' , CBTOP , 1 , 1 , WMAX , 1 )
223                    ELSE
223                           
224                        CALL SGEBR2D( ICTXT , 'Column' , CBTOP , 1 , 1 , WMAX , 1 ,
225       $                IAROW , MYCOL )
226                    END IF
227                END IF
228                IF( SCALE.LT.ABS( WMAX )*SMLNUM .OR. SCALE.EQ.ZERO )
228                       
229       $            GO TO 20
230                    CALL PSRSCL ( N , SCALE , WORK( IPX ) , IX , JX , DESCX , 1 )
231                END IF
232                GO TO 10
233            END IF
234  
235  *         Compute the estimate of the reciprocal condition number.
236  
237            IF( AINVNM.NE.ZERO )
237                   
238       $        RCOND =( ONE / AINVNM ) / ANORM
239  
240     20 CONTINUE
241  
242        CALL PB_TOPSET( ICTXT , 'Combine' , 'Columnwise' , COLCTOP )
243        CALL PB_TOPSET( ICTXT , 'Combine' , 'Rowwise' , ROWCTOP )
244  
245        RETURN
246  
247  *     End of PSPOCON
248  
249        END2835

Variables in Routine PSPOCON()

Summary Report
Data Type	Quantity	Size(byte)
CHARACTER	5	5
INTEGER	53	232
LOGICAL	3	3
REAL	12	48
TOTAL	73	288

List of Variables

CHARACTER

CBTOP COLCTOP NORMIN ROWCTOP UPLO

INTEGER

BLOCK_CYCLIC_2D CSRC_ CTXT_ DESCV( DLEN_ ) DESCX( DLEN_ )
DLEN_ DTYPE_ IA IACOL IAROW
ICEIL ICOFF ICTXT IDUM1( 3 ) IDUM2( 3 )
IIA INDXG2P INFO IPNL IPNU
IPV IPW IPX IROFF IV
IWORK IX IXX JA JJA
JV JX KASE LIWMIN LIWORK
LLD_ LWMIN LWORK M_ MB_
MYCOL MYROW N N_ NB_
NP NPCOL NPMOD NPROW NQ
NQMOD NUMROC RSRC_

LOGICAL

LQUERY LSAME UPPER

REAL

AINVNM ANORM ONE PSLAMCH RCOND
SCALE SL SMLNUM SU WMAX
WORK ZERO

Variables Dependence Graph
Put the mouse over a right hand side variable to display the corresponding line of the dependence
-
- -
AINVNM <--- ZEROAINVNM = ZERO

DESCX <--- IACOLDESCX( CSRC_ ) = IACOL{2DESCX( CSRC_ ) = IACOL, 3DESCX( CSRC_ ) = IACOL, 4DESCX( CSRC_ ) = IACOL}, MYCOLDESCX( CSRC_ ) = MYCOL{2DESCX( CSRC_ ) = MYCOL, 3DESCX( CSRC_ ) = MYCOL, 4DESCX( CSRC_ ) = MYCOL}

IACOL <--- INDXG2P, JA, MYCOL, NB_, NPCOL, CSRC_

IAROW <--- IA, INDXG2P, MB_, MYROW, NPROW, RSRC_

ICOFF <--- JAICOFF = MOD( JA-1, DESCA( NB_ ) ), NB_ICOFF = MOD( JA-1, DESCA( NB_ ) )

ICTXT <--- CTXT_ICTXT = DESCA( CTXT_ )

INFO <--- CTXT_INFO = -(600+CTXT_)

IPNL <--- IPVIPNL = IPV + NP, NPIPNL = IPV + NP

IPNU <--- IPNLIPNU = IPNL + NQ, NQIPNU = IPNL + NQ

IPV <--- IPXIPV = IPX + NP, NPIPV = IPX + NP

IPW <--- IPNUIPW = IPNU + NQ, NQIPW = IPNU + NQ

IROFF <--- IAIROFF = MOD( IA-1, DESCA( MB_ ) ), MB_IROFF = MOD( IA-1, DESCA( MB_ ) )

IV <--- IROFFIV = IROFF + 1

IWORK <--- LIWMINIWORK( 1 ) = LIWMIN{2IWORK( 1 ) = LIWMIN}

IX <--- IVIX = IV

JV <--- ICOFFJV = ICOFF + 1

JX <--- JVJX = JV

LIWMIN <--- NPMODLIWMIN = NPMOD

LWMIN <--- ICEILLWMIN = 2*NPMOD + 2*NQMOD +, NB_LWMIN = 2*NPMOD + 2*NQMOD +, NPCOLLWMIN = 2*NPMOD + 2*NQMOD +, NPMODLWMIN = 2*NPMOD + 2*NQMOD +, NPROWLWMIN = 2*NPMOD + 2*NQMOD +, NQMODLWMIN = 2*NPMOD + 2*NQMOD +

NORMIN <--- NNORMIN = 'N'

NP <--- IAROW, IROFF, MB_, MYROW, N, NPROW, NUMROC

NPMOD <--- IA, IAROW, MB_, MYROW, N, NPROW, NUMROC

NQ <--- IACOL, ICOFF, MYCOL, N, NB_, NPCOL, NUMROC

NQMOD <--- IACOL, JA, MYCOL, N, NB_, NPCOL, NUMROC

RCOND <--- ONERCOND = ONE{2RCOND = ONE}, ZERORCOND = ZERO

SCALE <--- SLSCALE = SL*SU, SUSCALE = SL*SU

SMLNUM <--- ICTXT, PSLAMCH

UPPER <--- LSAMEUPPER = LSAME( UPLO, 'U' ), UPLOUPPER = LSAME( UPLO, 'U' )

WORK <--- LWMINWORK( 1 ) = REAL( LWMIN )