Routine: PSPOSVX()  File: SRC\psposvx.f

 
 
# lines: 668
  # code: 668
  # comment: 0
  # blank:0
# Variables:89
# Callers:0
# Callings:9
# Words:319
# Keywords:209
 

 

..
     .. Array Arguments ..
     ..
  Purpose
  =======
  PSPOSVX uses the Cholesky factorization A = U**T*U or A = L*L**T to
  compute the solution to a real system of linear equations
        A(IA:IA+N-1,JA:JA+N-1) * X = B(IB:IB+N-1,JB:JB+NRHS-1),
  where A(IA:IA+N-1,JA:JA+N-1) is an N-by-N matrix and X and
  B(IB:IB+N-1,JB:JB+NRHS-1) are N-by-NRHS matrices.
  Error bounds on the solution and a condition estimate are also
  provided.  In the following comments Y denotes Y(IY:IY+M-1,JY:JY+K-1)
  a M-by-K matrix where Y can be A, AF, B and X.
  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
  Description
  ===========
  The following steps are performed:
  1. If FACT = 'E', real scaling factors are computed to equilibrate
     the system:
        diag(SR) * A * diag(SC) * inv(diag(SC)) * X = diag(SR) * B
     Whether or not the system will be equilibrated depends on the
     scaling of the matrix A, but if equilibration is used, A is
     overwritten by diag(SR)*A*diag(SC) and B by diag(SR)*B.
  2. If FACT = 'N' or 'E', the Cholesky decomposition is used to
     factor the matrix A (after equilibration if FACT = 'E') as
        A = U**T* U,  if UPLO = 'U', or
        A = L * L**T,  if UPLO = 'L',
     where U is an upper triangular matrix and L is a lower triangular
     matrix.
  3. The factored form of A is used to estimate the condition number
     of the matrix A.  If the reciprocal of the condition number is
     less than machine precision, steps 4-6 are skipped.
  4. The system of equations is solved for X using the factored form
     of A.
  5. Iterative refinement is applied to improve the computed solution
     matrix and calculate error bounds and backward error estimates
     for it.
  6. If equilibration was used, the matrix X is premultiplied by
     diag(SR) so that it solves the original system before
     equilibration.
  Arguments
  =========
  FACT    (global input) CHARACTER
          Specifies whether or not the factored form of the matrix A is
          supplied on entry, and if not, whether the matrix A should be
          equilibrated before it is factored.
          = 'F':  On entry, AF contains the factored form of A.
                  If EQUED = 'Y', the matrix A has been equilibrated
                  with scaling factors given by S.  A and AF will not
                  be modified.
          = 'N':  The matrix A will be copied to AF and factored.
          = 'E':  The matrix A will be equilibrated if necessary, then
                  copied to AF and factored.
  UPLO    (global input) CHARACTER
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
  N       (global input) INTEGER
          The number of rows and columns to be operated on, i.e. the
          order of the distributed submatrix A(IA:IA+N-1,JA:JA+N-1).
          N >= 0.
  NRHS    (global input) INTEGER
          The number of right hand sides, i.e., the number of columns
          of the distributed submatrices B and X.  NRHS >= 0.
  A       (local input/local output) REAL pointer into
          the local memory to an array of local dimension
          ( LLD_A, LOCc(JA+N-1) ).
          On entry, the symmetric matrix A, except if FACT = 'F' and
          EQUED = 'Y', then A must contain the equilibrated matrix
          diag(SR)*A*diag(SC).  If UPLO = 'U', the leading
          N-by-N upper triangular part of A contains the upper
          triangular part of the matrix A, and the strictly lower
          triangular part of A is not referenced.  If UPLO = 'L', the
          leading N-by-N lower triangular part of A contains the lower
          triangular part of the matrix A, and the strictly upper
          triangular part of A is not referenced.  A is not modified if
          FACT = 'F' or 'N', or if FACT = 'E' and EQUED = 'N' on exit.
          On exit, if FACT = 'E' and EQUED = 'Y', A is overwritten by
          diag(SR)*A*diag(SC).
  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.
  AF      (local input or local output) REAL pointer
          into the local memory to an array of local dimension
          ( LLD_AF, LOCc(JA+N-1)).
          If FACT = 'F', then AF is an input argument and on entry
          contains the triangular factor U or L from the Cholesky
          factorization A = U**T*U or A = L*L**T, in the same storage
          format as A.  If EQUED .ne. 'N', then AF is the factored form
          of the equilibrated matrix diag(SR)*A*diag(SC).
          If FACT = 'N', then AF is an output argument and on exit
          returns the triangular factor U or L from the Cholesky
          factorization A = U**T*U or A = L*L**T of the original
          matrix A.
          If FACT = 'E', then AF is an output argument and on exit
          returns the triangular factor U or L from the Cholesky
          factorization A = U**T*U or A = L*L**T of the equilibrated
          matrix A (see the description of A for the form of the
          equilibrated matrix).
  IAF     (global input) INTEGER
          The row index in the global array AF indicating the first
          row of sub( AF ).
  JAF     (global input) INTEGER
          The column index in the global array AF indicating the
          first column of sub( AF ).
  DESCAF  (global and local input) INTEGER array of dimension DLEN_.
          The array descriptor for the distributed matrix AF.
  EQUED   (global input/global output) CHARACTER
          Specifies the form of equilibration that was done.
          = 'N':  No equilibration (always true if FACT = 'N').
          = 'Y':  Equilibration was done, i.e., A has been replaced by
                  diag(SR) * A * diag(SC).
          EQUED is an input variable if FACT = 'F'; otherwise, it is an
          output variable.
  SR      (local input/local output) REAL array,
                                                    dimension (LLD_A)
          The scale factors for A distributed across process rows;
          not accessed if EQUED = 'N'.  SR is an input variable if
          FACT = 'F'; otherwise, SR is an output variable.
          If FACT = 'F' and EQUED = 'Y', each element of SR must be
          positive.
  SC      (local input/local output) REAL array,
                                              dimension (LOC(N_A))
          The scale factors for A distributed across
          process columns; not accessed if EQUED = 'N'. SC is an input
          variable if FACT = 'F'; otherwise, SC is an output variable.
          If FACT = 'F' and EQUED = 'Y', each element of SC must be
          positive.
  B       (local input/local output) REAL pointer into
          the local memory to an array of local dimension
          ( LLD_B, LOCc(JB+NRHS-1) ).
          On entry, the N-by-NRHS right-hand side matrix B.
          On exit, if EQUED = 'N', B is not modified; if TRANS = 'N'
          and EQUED = 'R' or 'B', B is overwritten by diag(R)*B; if
          TRANS = 'T' or 'C' and EQUED = 'C' or 'B', B is overwritten
          by diag(C)*B.
  IB      (global input) INTEGER
          The row index in the global array B indicating the first
          row of sub( B ).
  JB      (global input) INTEGER
          The column index in the global array B indicating the
          first column of sub( B ).
  DESCB   (global and local input) INTEGER array of dimension DLEN_.
          The array descriptor for the distributed matrix B.
  X       (local input/local output) REAL pointer into
          the local memory to an array of local dimension
          ( LLD_X, LOCc(JX+NRHS-1) ).
          If INFO = 0, the N-by-NRHS solution matrix X to the original
          system of equations.  Note that A and B are modified on exit
          if EQUED .ne. 'N', and the solution to the equilibrated
          system is inv(diag(SC))*X if TRANS = 'N' and EQUED = 'C' or
          'B', or inv(diag(SR))*X if TRANS = 'T' or 'C' and EQUED = 'R'
          or 'B'.
  IX      (global input) INTEGER
          The row index in the global array X indicating the first
          row of sub( X ).
  JX      (global input) INTEGER
          The column index in the global array X indicating the
          first column of sub( X ).
  DESCX   (global and local input) INTEGER array of dimension DLEN_.
          The array descriptor for the distributed matrix X.
  RCOND   (global output) REAL
          The estimate of the reciprocal condition number of the matrix
          A after equilibration (if done).  If RCOND is less than the
          machine precision (in particular, if RCOND = 0), the matrix
          is singular to working precision.  This condition is
          indicated by a return code of INFO > 0, and the solution and
          error bounds are not computed.
  FERR    (local output) REAL array, dimension (LOC(N_B))
          The estimated forward error bounds for each solution vector
          X(j) (the j-th column of the solution matrix X).
          If XTRUE is the true solution, FERR(j) bounds the magnitude
          of the largest entry in (X(j) - XTRUE) divided by
          the magnitude of the largest entry in X(j).  The quality of
          the error bound depends on the quality of the estimate of
          norm(inv(A)) computed in the code; if the estimate of
          norm(inv(A)) is accurate, the error bound is guaranteed.
  BERR    (local output) REAL array, dimension (LOC(N_B))
          The componentwise relative backward error of each solution
          vector X(j) (i.e., the smallest relative change in
          any entry of A or B that makes X(j) an exact solution).
  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 = MAX( PSPOCON( LWORK ), PSPORFS( LWORK ) )
                  + LOCr( N_A ).
          LWORK = 3*DESCA( LLD_ )
          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 = DESCA( LLD_ )
          LIWORK = LOCr(N_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 INFO = -i, the i-th argument had an illegal value
          > 0: if INFO = i, and i is
               <= N: if INFO = i, the leading minor of order i of A
                     is not positive definite, so the factorization
                     could not be completed, and the solution and error
                     bounds could not be computed.
               = N+1: RCOND is less than machine precision.  The
                     factorization has been completed, but the matrix
                     is singular to working precision, and the solution
                     and error bounds have not been computed.
  =====================================================================
     .. Parameters ..

 
Display dynamic version Find AutoScroll Reload FontSize: - + Hide Comments Hide Blanks Frame FullScreen MailPrint

 
001        SUBROUTINE PSPOSVX( FACT , UPLO , N , NRHS , A , IA , JA , DESCA , AF ,
002       $IAF , JAF , DESCAF , EQUED , SR , SC , B , IB , JB ,
003       $DESCB , X , IX , JX , DESCX , RCOND , FERR , BERR ,
004       $WORK , LWORK , IWORK , LIWORK , INFO )
005  
006  *     -- ScaLAPACK routine(version 1.7) --
007  *     University of Tennessee , Knoxville , Oak Ridge National Laboratory ,
008  *     and University of California , Berkeley.
009  *     December 31 , 1998
010  
011  *     .. Scalar Arguments ..
012        CHARACTER EQUED , FACT , UPLO
013        INTEGER IA , IAF , IB , INFO , IX , JA , JAF , JB , JX , LIWORK ,
014       $LWORK , N , NRHS
015        REAL RCOND
016        INTEGER BLOCK_CYCLIC_2D , CSRC_ , CTXT_ , DLEN_ , DTYPE_ ,
017       $LLD_ , MB_ , M_ , NB_ , N_ , RSRC_
018        PARAMETER( BLOCK_CYCLIC_2D = 1 , DLEN_ = 9 , DTYPE_ = 1 ,
019       $CTXT_ = 2 , M_ = 3 , N_ = 4 , MB_ = 5 , NB_ = 6 ,
020       $RSRC_ = 7 , CSRC_ = 8 , LLD_ = 9 )
021        REAL ONE , ZERO
022        PARAMETER( ONE = 1.0E + 0 , ZERO = 0.0E + 0 )
023  *     ..
024  *     .. Local Scalars ..
025        LOGICAL EQUIL , LQUERY , NOFACT , RCEQU
026        INTEGER I , IACOL , IAROW , IAFROW , IBROW , IBCOL , ICOFF ,
027       $ICOFFA , ICTXT , IDUMM , IIA , IIB , IIX , INFEQU ,
028       $IROFF , IROFFA , IROFFAF , IROFFB , IROFFX , IXCOL ,
029       $IXROW , J , JJA , JJB , JJX , LDB , LDX , LIWMIN ,
030       $LWMIN , MYCOL , MYROW , NP , NPCOL , NPROW , NRHSQ ,
031       $NQ
032        REAL AMAX , ANORM , BIGNUM , SCOND , SMAX , SMIN , SMLNUM
033  *     ..
034  *     .. Local Arrays ..
035        INTEGER IDUM1( 5 ) , IDUM2( 5 )
036  *     ..
037  *     .. External Subroutines ..
038        EXTERNAL BLACS_GRIDINFO , CHK1MAT , PCHK2MAT , INFOG2L ,
039       $PSPOCON , PSPOEQU , PSPORFS ,
040       $PSPOTRF , PSPOTRS ,
041       $PSLACPY , PSLAQSY , PXERBLA ,
042       $SGAMN2D , SGAMX2D
043  *     ..
044  *     .. External Functions ..
045        LOGICAL LSAME
046        INTEGER INDXG2P , NUMROC
047        REAL PSLANSY , PSLAMCH
048        EXTERNAL INDXG2P , LSAME , NUMROC , PSLANSY , PSLAMCH
049  *     ..
050  *     .. Intrinsic Functions ..
051        INTRINSIC ICHAR , MAX , MIN , MOD
052  *     ..
053  *     .. Executable Statements ..
054  
055  *     Get grid parameters
056  
057        ICTXT = DESCA( CTXT_ )
058        CALL BLACS_GRIDINFO( ICTXT , NPROW , NPCOL , MYROW , MYCOL )
059  
060  *     Test the input parameters
061  
062        INFO = 0
063        IF( NPROW.EQ. - 1 ) THEN
064            INFO = - (800 + CTXT_)
065        ELSE
066            CALL CHK1MAT( N , 3 , N , 3 , IA , JA , DESCA , 8 , INFO )
067            IF( LSAME( FACT , 'F' ) )
068       $        CALL CHK1MAT( N , 3 , N , 3 , IAF , JAF , DESCAF , 12 , INFO )
069                CALL CHK1MAT( N , 3 , NRHS , 4 , IB , JB , DESCB , 20 , INFO )
070                IF( INFO.EQ.0 ) THEN
071                    IAROW = INDXG2P( IA , DESCA( MB_ ) , MYROW , DESCA( RSRC_ ) ,
072       $            NPROW )
073                    IAFROW = INDXG2P( IAF , DESCAF( MB_ ) , MYROW ,
074       $            DESCAF( RSRC_ ) , NPROW )
075                    IBROW = INDXG2P( IB , DESCB( MB_ ) , MYROW , DESCB( RSRC_ ) ,
076       $            NPROW )
077                    IXROW = INDXG2P( IX , DESCX( MB_ ) , MYROW , DESCX( RSRC_ ) ,
078       $            NPROW )
079                    IROFFA = MOD( IA - 1 , DESCA( MB_ ) )
080                    IROFFAF = MOD( IAF - 1 , DESCAF( MB_ ) )
081                    ICOFFA = MOD( JA - 1 , DESCA( NB_ ) )
082                    IROFFB = MOD( IB - 1 , DESCB( MB_ ) )
083                    IROFFX = MOD( IX - 1 , DESCX( MB_ ) )
084                    CALL INFOG2L( IA , JA , DESCA , NPROW , NPCOL , MYROW ,
085       $            MYCOL , IIA , JJA , IAROW , IACOL )
086                    NP = NUMROC( N + IROFFA , DESCA( MB_ ) , MYROW , IAROW , NPROW )
087                    IF( MYROW.EQ.IAROW )
088       $                NP = NP - IROFFA
089                        NQ = NUMROC( N + ICOFFA , DESCA( NB_ ) , MYCOL , IACOL , NPCOL )
090                        IF( MYCOL.EQ.IACOL )
091       $                    NQ = NQ - ICOFFA
092                            LWMIN = 3*DESCA( LLD_ )
093                            LIWMIN = NP
094                            NOFACT = LSAME( FACT , 'N' )
095                            EQUIL = LSAME( FACT , 'E' )
096                            IF( NOFACT .OR. EQUIL ) THEN
097                                EQUED = 'N'
098                                RCEQU = .FALSE.
099                            ELSE
100                                RCEQU = LSAME( EQUED , 'Y' )
101                                SMLNUM = PSLAMCH( ICTXT , 'Safe minimum' )
102                                BIGNUM = ONE / SMLNUM
103                            END IF
104                            IF( .NOT.NOFACT .AND. .NOT.EQUIL .AND.
105       $                        .NOT.LSAME( FACT , 'F' ) ) THEN
106                                INFO = - 1
107                            ELSE IF( .NOT.LSAME( UPLO , 'U' ) .AND.
108       $                        .NOT.LSAME( UPLO , 'L' ) ) THEN
109                                INFO = - 2
110                            ELSE IF( IROFFA.NE.0 ) THEN
111                                INFO = - 6
112                            ELSE IF( ICOFFA.NE.0 .OR. IROFFA.NE.ICOFFA ) THEN
113                                INFO = - 7
114                            ELSE IF( DESCA( MB_ ).NE.DESCA( NB_ ) ) THEN
115                                INFO = - (800 + NB_)
116                            ELSE IF( IAFROW.NE.IAROW ) THEN
117                                INFO = - 10
118                            ELSE IF( IROFFAF.NE.0 ) THEN
119                                INFO = - 10
120                            ELSE IF( ICTXT.NE.DESCAF( CTXT_ ) ) THEN
121                                INFO = - (1200 + CTXT_)
122                            ELSE IF( LSAME( FACT , 'F' ) .AND. .NOT.
123       $( RCEQU .OR. LSAME( EQUED , 'N' ) ) ) THEN
124                                INFO = - 13
125                            ELSE
126                                IF( RCEQU ) THEN
127  
128                                    SMIN = BIGNUM
129                                    SMAX = ZERO
130                                    DO 10 J = IIA , IIA + NP - 1
131                                        SMIN = MIN( SMIN , SR( J ) )
132                                        SMAX = MAX( SMAX , SR( J ) )
133     10                             CONTINUE
134                                    CALL SGAMN2D( ICTXT , 'Columnwise' , ' ' , 1 , 1 , SMIN ,
135       $                            1 , IDUMM , IDUMM , - 1 , - 1 , MYCOL )
136                                    CALL SGAMX2D( ICTXT , 'Columnwise' , ' ' , 1 , 1 , SMAX ,
137       $                            1 , IDUMM , IDUMM , - 1 , - 1 , MYCOL )
138                                    IF( SMIN.LE.ZERO ) THEN
139                                        INFO = - 14
140                                    ELSE IF( N.GT.0 ) THEN
141                                        SCOND = MAX( SMIN , SMLNUM ) / MIN( SMAX , BIGNUM )
142                                    ELSE
143                                        SCOND = ONE
144                                    END IF
145                                END IF
146                            END IF
147                        END IF
148  
149                        WORK( 1 ) = REAL( LWMIN )
150                        IWORK( 1 ) = LIWMIN
151                        LQUERY =( LWORK.EQ. - 1 .OR. LIWORK.EQ. - 1 )
152                        IF( INFO.EQ.0 ) THEN
153                            IF( IBROW.NE.IAROW ) THEN
154                                INFO = - 18
155                            ELSE IF( IXROW.NE.IBROW ) THEN
156                                INFO = - 22
157                            ELSE IF( DESCB( MB_ ).NE.DESCA( NB_ ) ) THEN
158                                INFO = - (2000 + NB_)
159                            ELSE IF( ICTXT.NE.DESCB( CTXT_ ) ) THEN
160                                INFO = - (2000 + CTXT_)
161                            ELSE IF( ICTXT.NE.DESCX( CTXT_ ) ) THEN
162                                INFO = - (2400 + CTXT_)
163                            ELSE IF( LWORK.LT.LWMIN .AND. .NOT.LQUERY ) THEN
164                                INFO = - 28
165                            ELSE IF( LIWORK.LT.LIWMIN .AND. .NOT.LQUERY ) THEN
166                                INFO = - 30
167                            END IF
168                            IDUM1( 1 ) = ICHAR( FACT )
169                            IDUM2( 1 ) = 1
170                            IDUM1( 2 ) = ICHAR( UPLO )
171                            IDUM2( 2 ) = 2
172                            IF( LSAME( FACT , 'F' ) ) THEN
173                                IDUM1( 3 ) = ICHAR( EQUED )
174                                IDUM2( 3 ) = 13
175                                IF( LWORK.EQ. - 1 ) THEN
176                                    IDUM1( 4 ) = - 1
177                                ELSE
178                                    IDUM1( 4 ) = 1
179                                END IF
180                                IDUM2( 4 ) = 28
181                                IF( LIWORK.EQ. - 1 ) THEN
182                                    IDUM1( 5 ) = - 1
183                                ELSE
184                                    IDUM1( 5 ) = 1
185                                END IF
186                                IDUM2( 5 ) = 30
187                                CALL PCHK2MAT( N , 3 , N , 3 , IA , JA , DESCA , 8 , N , 3 , NRHS ,
188       $                        4 , IB , JB , DESCB , 19 , 5 , IDUM1 , IDUM2 ,
189       $                        INFO )
190                            ELSE
191                                IF( LWORK.EQ. - 1 ) THEN
192                                    IDUM1( 3 ) = - 1
193                                ELSE
194                                    IDUM1( 3 ) = 1
195                                END IF
196                                IDUM2( 3 ) = 28
197                                IF( LIWORK.EQ. - 1 ) THEN
198                                    IDUM1( 4 ) = - 1
199                                ELSE
200                                    IDUM1( 4 ) = 1
201                                END IF
202                                IDUM2( 4 ) = 30
203                                CALL PCHK2MAT( N , 3 , N , 3 , IA , JA , DESCA , 8 , N , 3 , NRHS ,
204       $                        4 , IB , JB , DESCB , 19 , 4 , IDUM1 , IDUM2 ,
205       $                        INFO )
206                            END IF
207                        END IF
208                    END IF
209  
210                    IF( INFO.NE.0 ) THEN
211                        CALL PXERBLA( ICTXT , 'PSPOSVX' , - INFO )
212                        RETURN
213                    ELSE IF( LQUERY ) THEN
214                        RETURN
215                    END IF
216  
217                    IF( EQUIL ) THEN
218  
219  *                     Compute row and column scalings to equilibrate the matrix A.
220  
221                        CALL PSPOEQU ( N , A , IA , JA , DESCA , SR , SC , SCOND , AMAX ,
222       $                INFEQU )
223  
224                        IF( INFEQU.EQ.0 ) THEN
225  
226  *                         Equilibrate the matrix.
227  
228                            CALL PSLAQSY ( UPLO , N , A , IA , JA , DESCA , SR , SC , SCOND ,
229       $                    AMAX , EQUED )
230                            RCEQU = LSAME( EQUED , 'Y' )
231                        END IF
232                    END IF
233  
234  *                 Scale the right - hand side.
235  
236                    CALL INFOG2L( IB , JB , DESCB , NPROW , NPCOL , MYROW , MYCOL , IIB ,
237       $            JJB , IBROW , IBCOL )
238                    LDB = DESCB( LLD_ )
239                    IROFF = MOD( IB - 1 , DESCB( MB_ ) )
240                    ICOFF = MOD( JB - 1 , DESCB( NB_ ) )
241                    NP = NUMROC( N + IROFF , DESCB( MB_ ) , MYROW , IBROW , NPROW )
242                    NRHSQ = NUMROC( NRHS + ICOFF , DESCB( NB_ ) , MYCOL , IBCOL , NPCOL )
243                    IF( MYROW.EQ.IBROW ) NP = NP - IROFF
244                    IF( MYCOL.EQ.IBCOL ) NRHSQ = NRHSQ - ICOFF
245  
246                    IF( RCEQU ) THEN
247                        DO 30 J = JJB , JJB + NRHSQ - 1
248                            DO 20 I = IIB , IIB + NP - 1
249                                B( I + ( J - 1 )*LDB ) = SR( I )*B( I + ( J - 1 )*LDB )
250     20                     CONTINUE
251     30                 CONTINUE
252                    END IF
253  
254                    IF( NOFACT .OR. EQUIL ) THEN
255  
256  *                     Compute the Cholesky factorization A = U'*U or A = L*L'.
257  
258                        CALL PSLACPY ( 'Full' , N , N , A , IA , JA , DESCA , AF , IAF , JAF ,
259       $                DESCAF )
260                        CALL PSPOTRF ( UPLO , N , AF , IAF , JAF , DESCAF , INFO )
261  
262  *                     Return if INFO is non - zero.
263  
264                        IF( INFO.NE.0 ) THEN
265                            IF( INFO.GT.0 )
266       $                        RCOND = ZERO
267                                RETURN
268                            END IF
269                        END IF
270  
271  *                     Compute the norm of the matrix A.
272  
273                        ANORM = PSLANSY( '1' , UPLO , N , A , IA , JA , DESCA , WORK )
274  
275  *                     Compute the reciprocal of the condition number of A.
276  
277                        CALL PSPOCON ( UPLO , N , AF , IAF , JAF , DESCAF , ANORM , RCOND , WORK ,
278       $                LWORK , IWORK , LIWORK , INFO )
279  
280  *                     Return if the matrix is singular to working precision.
281  
282                        IF( RCOND.LT.PSLAMCH( ICTXT , 'Epsilon' ) ) THEN
283                            INFO = IA + N
284                            RETURN
285                        END IF
286  
287  *                     Compute the solution matrix X.
288  
289                        CALL PSLACPY ( 'Full' , N , NRHS , B , IB , JB , DESCB , X , IX , JX ,
290       $                DESCX )
291                        CALL PSPOTRS ( UPLO , N , NRHS , AF , IAF , JAF , DESCAF , X , IX , JX ,
292       $                DESCX , INFO )
293  
294  *                     Use iterative refinement to improve the computed solution and
295  *                     compute error bounds and backward error estimates for it.
296  
297                        CALL PSPORFS ( UPLO , N , NRHS , A , IA , JA , DESCA , AF , IAF , JAF ,
298       $                DESCAF , B , IB , JB , DESCB , X , IX , JX , DESCX , FERR ,
299       $                BERR , WORK , LWORK , IWORK , LIWORK , INFO )
300  
301  *                     Transform the solution matrix X to a solution of the original
302  *                     system.
303  
304                        CALL INFOG2L( IX , JX , DESCX , NPROW , NPCOL , MYROW , MYCOL , IIX ,
305       $                JJX , IXROW , IXCOL )
306                        LDX = DESCX( LLD_ )
307                        IROFF = MOD( IX - 1 , DESCX( MB_ ) )
308                        ICOFF = MOD( JX - 1 , DESCX( NB_ ) )
309                        NP = NUMROC( N + IROFF , DESCX( MB_ ) , MYROW , IXROW , NPROW )
310                        NRHSQ = NUMROC( NRHS + ICOFF , DESCX( NB_ ) , MYCOL , IXCOL , NPCOL )
311                        IF( MYROW.EQ.IBROW ) NP = NP - IROFF
312                        IF( MYCOL.EQ.IBCOL ) NRHSQ = NRHSQ - ICOFF
313  
314                        IF( RCEQU ) THEN
315                            DO 50 J = JJX , JJX + NRHSQ - 1
316                                DO 40 I = IIX , IIX + NP - 1
317                                    X( I + ( J - 1 )*LDX ) = SR( I )*X( I + ( J - 1 )*LDX )
318     40                         CONTINUE
319     50                     CONTINUE
320                            DO 60 J = JJX , JJX + NRHSQ - 1
321                                FERR( J ) = FERR( J ) / SCOND
322     60                     CONTINUE
323                        END IF
324  
325                        WORK( 1 ) = REAL( LWMIN )
326                        IWORK( 1 ) = LIWMIN
327                        RETURN
328  
329  *                     End of PSPOSVX
330  
331                    END