Routine PDPOCON() Documented By Universal Report

Routine: PDPOCON()

File: SRC\pdpocon.f

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

..
     .. Array Arguments ..
     ..
  Purpose
  =======
  PDPOCON 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
  PDPOTRF.
  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) DOUBLE PRECISION 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 PDPOTRF.
  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) DOUBLE PRECISION
          The 1-norm (or infinity-norm) of the symmetric distributed
          matrix A(IA:IA+N-1,JA:JA+N-1).
  RCOND   (global output) DOUBLE PRECISION
          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) DOUBLE PRECISION 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 PDPOCON( 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        DOUBLE PRECISION 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        DOUBLE PRECISION ONE , ZERO
019        PARAMETER( ONE = 1.0D + 0 , ZERO = 0.0D + 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        DOUBLE PRECISION AINVNM , SCALE , SL , SU , SMLNUM
029        DOUBLE PRECISION 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 , DGEBR2D ,
037       $DGEBS2D , INFOG2L , PCHK1MAT , PDAMAX ,
038       $PDLATRS , PDLACON , PDRSCL , PB_TOPGET ,
039       $PB_TOPSET , PXERBLA
040  *     ..
041  *     .. External Functions ..
042        LOGICAL LSAME
043        INTEGER ICEIL , INDXG2P , NUMROC
044        DOUBLE PRECISION PDLAMCH
045        EXTERNAL ICEIL , INDXG2P , LSAME , NUMROC , PDLAMCH
046  *     ..
047  *     .. Intrinsic Functions ..
048        INTRINSIC ABS , DBLE , ICHAR , MAX , MOD
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 ) = DBLE( 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 , 'PDPOCON' , - 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 = PDLAMCH( 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 PDLACON ( 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 PDLATRS ( '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 PDLATRS ( '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 PDLATRS ( '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 PDLATRS ( '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 PDAMAX( 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 DGEBS2D( ICTXT , 'Column' , CBTOP , 1 , 1 , WMAX , 1 )
223                    ELSE
223                           
224                        CALL DGEBR2D( 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 PDRSCL ( 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 PDPOCON
248  
249        END2835

Variables in Routine PDPOCON()

Summary Report
Data Type	Quantity	Size(byte)
CHARACTER	5	5
DOUBLE PRECISION	11	44
INTEGER	53	232
LOGICAL	3	3
REAL	1	4
TOTAL	73	288

List of Variables

CHARACTER

CBTOP COLCTOP NORMIN ROWCTOP UPLO

DOUBLE PRECISION

AINVNM ANORM ONE PDLAMCH RCOND
SCALE SL SMLNUM SU WMAX
ZERO

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

WORK

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, PDLAMCH

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

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