Routine PCPOCON() Documented By Universal Report

Routine: PCPOCON()

File: SRC\pcpocon.f

# lines:	413
# code:	413
# comment:	0
# blank:	0
# Variables:	75
# Callers:	1
# Callings:	4
# Words:	229
# Keywords:	141

..
     .. Array Arguments ..
     ..
  Purpose
  =======
  PCPOCON estimates the reciprocal of the condition number (in the
  1-norm) of a complex Hermitian positive definite distributed matrix
  using the Cholesky factorization A = U**H*U or A = L*L**H computed by
  PCPOTRF.
  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) COMPLEX 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 PCPOTRF.
  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 hermitian 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) COMPLEX 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)) +
          MAX( 2, MAX(NB_A*MAX(1,CEIL(P-1,Q)),LOCc(N+MOD(JA-1,NB_A)) +
          NB_A*MAX(1,CEIL(Q-1,P))) ).
          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.
  RWORK   (local workspace/local output) REAL array,
                                                  dimension (LRWORK)
          On exit, RWORK(1) returns the minimal and optimal LRWORK.
  LRWORK  (local or global input) INTEGER
          The dimension of the array RWORK.
          LRWORK is local input and must be at least
          LRWORK >= 2*LOCc(N+MOD(JA-1,NB_A)).
          If LRWORK = -1, then LRWORK 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 PCPOCON( UPLO , N , A , IA , JA , DESCA , ANORM , RCOND , WORK ,
002       $LWORK , RWORK , LRWORK , 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 , LRWORK , 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 , LRWMIN , LWMIN , MYCOL , MYROW , NP ,
027       $NPCOL , NPROW , NPMOD , NQ , NQMOD
028        REAL AINVNM , SCALE , SL , SU , SMLNUM
029        COMPLEX WMAX , ZDUM
030  *     ..
031  *     .. Local Arrays ..
032        INTEGER DESCV( DLEN_ ) , DESCX( DLEN_ ) , IDUM1( 3 ) ,
033       $IDUM2( 3 )
034  *     ..
035  *     .. External Subroutines ..
036        EXTERNAL BLACS_GRIDINFO , CGEBR2D , CGEBS2D , CHK1MAT ,
037       $DESCSET , INFOG2L , PCAMAX , PCHK1MAT , PCLATRS ,
038       $PCLACON , PCSRSCL , PB_TOPGET , PB_TOPSET ,
039       $PXERBLA
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 , AIMAG , ICHAR , MAX , MOD , REAL
049  *     ..
050  *     .. Statement Functions ..
051        REAL CABS1
052  *     ..
053  *     .. Statement Function definitions ..
054        CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
055  *     ..
056  *     .. Executable Statements ..
057  
058  *     Get grid parameters
059  
060        ICTXT = DESCA( CTXT_ )
061        CALL BLACS_GRIDINFO( ICTXT , NPROW , NPCOL , MYROW , MYCOL )
062  
063  *     Test the input parameters
064  
065        INFO = 0
066        IF( NPROW.EQ. - 1 ) THEN
066               
067            INFO = - (600 + CTXT_)
068        ELSE
068               
069            CALL CHK1MAT( N , 2 , N , 2 , IA , JA , DESCA , 6 , INFO )
070            IF( INFO.EQ.0 ) THEN
070                   
071                UPPER = LSAME( UPLO , 'U' )
072                IAROW = INDXG2P( IA , DESCA( MB_ ) , MYROW , DESCA( RSRC_ ) ,
073       $        NPROW )
074                IACOL = INDXG2P( JA , DESCA( NB_ ) , MYCOL , DESCA( CSRC_ ) ,
075       $        NPCOL )
076                NPMOD = NUMROC( N + MOD( IA - 1 , DESCA( MB_ ) ) , DESCA( MB_ ) ,
077       $        MYROW , IAROW , NPROW )
078                NQMOD = NUMROC( N + MOD( JA - 1 , DESCA( NB_ ) ) , DESCA( NB_ ) ,
079       $        MYCOL , IACOL , NPCOL )
080                LWMIN = 2*NPMOD +
081       $        MAX( 2 , MAX( DESCA( NB_ )*
082       $        MAX( 1 , ICEIL( NPROW - 1 , NPCOL ) ) , NQMOD +
083       $        DESCA( NB_ )*
084       $        MAX( 1 , ICEIL( NPCOL - 1 , NPROW ) ) ) )
085                WORK( 1 ) = REAL( LWMIN )
086                LRWMIN = 2*NQMOD
087                RWORK( 1 ) = REAL( LRWMIN )
088                LQUERY =( LWORK.EQ. - 1 .OR. LRWORK.EQ. - 1 )
089  
090                IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO , 'L' ) ) THEN
090                       
091                    INFO = - 1
092                ELSE IF( ANORM.LT.ZERO ) THEN
092                       
093                    INFO = - 7
094                ELSE IF( LWORK.LT.LWMIN .AND. .NOT.LQUERY ) THEN
094                       
095                    INFO = - 10
096                ELSE IF( LRWORK.LT.LRWMIN .AND. .NOT.LQUERY ) THEN
096                       
097                    INFO = - 12
098                END IF
099            END IF
100  
101            IF( UPPER ) THEN
101                   
102                IDUM1( 1 ) = ICHAR( 'U' )
103            ELSE
103                   
104                IDUM1( 1 ) = ICHAR( 'L' )
105            END IF
106            IDUM2( 1 ) = 1
107            IF( LWORK.EQ. - 1 ) THEN
107                   
108                IDUM1( 2 ) = - 1
109            ELSE
109                   
110                IDUM1( 2 ) = 1
111            END IF
112            IDUM2( 2 ) = 10
113            IF( LRWORK.EQ. - 1 ) THEN
113                   
114                IDUM1( 3 ) = - 1
115            ELSE
115                   
116                IDUM1( 3 ) = 1
117            END IF
118            IDUM2( 3 ) = 12
119            CALL PCHK1MAT( N , 2 , N , 2 , IA , JA , DESCA , 6 , 3 , IDUM1 , IDUM2 ,
120       $    INFO )
121        END IF
122  
123        IF( INFO.NE.0 ) THEN
123               
124            CALL PXERBLA( ICTXT , 'PCPOCON' , - INFO )
125            RETURN
126        ELSE IF( LQUERY ) THEN
126               
127            RETURN
128        END IF
129  
130  *     Quick return if possible
131  
132        RCOND = ZERO
133        IF( N.EQ.0 ) THEN
133               
134            RCOND = ONE
135            RETURN
136        ELSE IF( ANORM.EQ.ZERO ) THEN
136               
137            RETURN
138        ELSE IF( N.EQ.1 ) THEN
138               
139            RCOND = ONE
140            RETURN
141        END IF
142  
143        CALL PB_TOPGET( ICTXT , 'Combine' , 'Columnwise' , COLCTOP )
144        CALL PB_TOPGET( ICTXT , 'Combine' , 'Rowwise' , ROWCTOP )
145        CALL PB_TOPSET( ICTXT , 'Combine' , 'Columnwise' , '1 - tree' )
146        CALL PB_TOPSET( ICTXT , 'Combine' , 'Rowwise' , '1 - tree' )
147  
148        SMLNUM = PSLAMCH( ICTXT , 'Safe minimum' )
149        IROFF = MOD( IA - 1 , DESCA( MB_ ) )
150        ICOFF = MOD( JA - 1 , DESCA( NB_ ) )
151        CALL INFOG2L( IA , JA , DESCA , NPROW , NPCOL , MYROW , MYCOL , IIA , JJA ,
152       $IAROW , IACOL )
153        NP = NUMROC( N + IROFF , DESCA( MB_ ) , MYROW , IAROW , NPROW )
154        NQ = NUMROC( N + ICOFF , DESCA( NB_ ) , MYCOL , IACOL , NPCOL )
155        IV = IROFF + 1
156        IX = IV
157        JV = ICOFF + 1
158        JX = JV
159  
160        IPX = 1
161        IPV = IPX + NP
162        IPW = IPV + NP
163        IPNL = 1
164        IPNU = IPNL + NQ
165  
166        CALL DESCSET( DESCV , N + IROFF , 1 , DESCA( MB_ ) , 1 , IAROW , MYCOL ,
167       $ICTXT , MAX( 1 , NP ) )
168        CALL DESCSET( DESCX , N + IROFF , 1 , DESCA( MB_ ) , 1 , IAROW , MYCOL ,
169       $ICTXT , MAX( 1 , NP ) )
170  
171  *     Estimate the 1 - norm(or I - norm) of inv(A).
172  
173        AINVNM = ZERO
174        KASE = 0
175        NORMIN = 'N'
176  
177     10 CONTINUE
178        CALL PCLACON ( N , WORK( IPV ) , IV , JV , DESCV , WORK( IPX ) , IX , JX ,
179       $DESCX , AINVNM , KASE )
180        IF( KASE.NE.0 ) THEN
180               
181            IF( UPPER ) THEN
182  
183  *             Multiply by inv(U').
184  
184                   
185                DESCX( CSRC_ ) = IACOL
186                CALL PCLATRS ( 'Upper' , 'Conjugate transpose' , 'Non - unit' ,
187       $        NORMIN , N , A , IA , JA , DESCA , WORK( IPX ) ,
188       $        IX , JX , DESCX , SL , RWORK( IPNL ) ,
189       $        WORK( IPW ) )
190                DESCX( CSRC_ ) = MYCOL
191                NORMIN = 'Y'
192  
193  *             Multiply by inv(U).
194  
195                DESCX( CSRC_ ) = IACOL
196                CALL PCLATRS ( 'Upper' , 'No transpose' , 'Non - unit' , NORMIN ,
197       $        N , A , IA , JA , DESCA , WORK( IPX ) , IX , JX ,
198       $        DESCX , SU , RWORK( IPNU ) , WORK( IPW ) )
199                DESCX( CSRC_ ) = MYCOL
200            ELSE
201  
202  *             Multiply by inv(L).
203  
203                   
204                DESCX( CSRC_ ) = IACOL
205                CALL PCLATRS ( 'Lower' , 'No transpose' , 'Non - unit' , NORMIN ,
206       $        N , A , IA , JA , DESCA , WORK( IPX ) , IX , JX ,
207       $        DESCX , SL , RWORK( IPNL ) , WORK( IPW ) )
208                DESCX( CSRC_ ) = MYCOL
209                NORMIN = 'Y'
210  
211  *             Multiply by inv(L').
212  
213                DESCX( CSRC_ ) = IACOL
214                CALL PCLATRS ( 'Lower' , 'Conjugate transpose' , 'Non - unit' ,
215       $        NORMIN , N , A , IA , JA , DESCA , WORK( IPX ) ,
216       $        IX , JX , DESCX , SU , RWORK( IPNU ) ,
217       $        WORK( IPW ) )
218                DESCX( CSRC_ ) = MYCOL
219            END IF
220  
221  *         Multiply by 1 / SCALE if doing so will not cause overflow.
222  
223            SCALE = SL*SU
224            IF( SCALE.NE.ONE ) THEN
224                   
225                CALL PCAMAX( N , WMAX , IXX , WORK( IPX ) , IX , JX , DESCX , 1 )
226                IF( DESCX( M_ ).EQ.1 .AND. N.EQ.1 ) THEN
226                       
227                    CALL PB_TOPGET( ICTXT , 'Broadcast' , 'Columnwise' , CBTOP )
228                    IF( MYROW.EQ.IAROW ) THEN
228                           
229                        CALL CGEBS2D( ICTXT , 'Column' , CBTOP , 1 , 1 , WMAX , 1 )
230                    ELSE
230                           
231                        CALL CGEBR2D( ICTXT , 'Column' , CBTOP , 1 , 1 , WMAX , 1 ,
232       $                IAROW , MYCOL )
233                    END IF
234                END IF
235                IF( SCALE.LT.CABS1( WMAX )*SMLNUM .OR. SCALE.EQ.ZERO )
235                       
236       $            GO TO 20
237                    CALL PCSRSCL ( N , SCALE , WORK( IPX ) , IX , JX , DESCX , 1 )
238                END IF
239                GO TO 10
240            END IF
241  
242  *         Compute the estimate of the reciprocal condition number.
243  
244            IF( AINVNM.NE.ZERO )
244                   
245       $        RCOND =( ONE / AINVNM ) / ANORM
246  
247     20 CONTINUE
248  
249        CALL PB_TOPSET( ICTXT , 'Combine' , 'Columnwise' , COLCTOP )
250        CALL PB_TOPSET( ICTXT , 'Combine' , 'Rowwise' , ROWCTOP )
251  
252        RETURN
253  
254  *     End of PCPOCON
255  
256        END3235

Variables in Routine PCPOCON()

Summary Report
Data Type	Quantity	Size(byte)
CHARACTER	5	5
COMPLEX	2	8
INTEGER	52	228
LOGICAL	3	3
REAL	13	52
TOTAL	75	296

List of Variables

CHARACTER

CBTOP COLCTOP NORMIN ROWCTOP UPLO

COMPLEX

WMAX ZDUM

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
IX IXX JA JJA JV
JX KASE LLD_ LRWMIN LRWORK
LWMIN LWORK M_ MB_ MYCOL
MYROW N N_ NB_ NP
NPCOL NPMOD NPROW NQ NQMOD
NUMROC RSRC_

LOGICAL

LQUERY LSAME UPPER

REAL

AINVNM ANORM CABS1 ONE PSLAMCH
RCOND RWORK SCALE SL SMLNUM
SU 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

CABS1 <--- ZDUM

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_)

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

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

IPW <--- IPVIPW = IPV + NP, NPIPW = IPV + NP

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

IV <--- IROFFIV = IROFF + 1

IX <--- IVIX = IV

JV <--- ICOFFJV = ICOFF + 1

JX <--- JVJX = JV

LRWMIN <--- NQMODLRWMIN = 2*NQMOD

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

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

RWORK <--- LRWMINRWORK( 1 ) = REAL( LRWMIN )

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

SMLNUM <--- ICTXT, PSLAMCH

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

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