PCGESVX(3)    ScaLAPACK routine of NEC Numeric Library Collection   PCGESVX(3)



NAME
       PCGESVX - use the LU factorization to compute the solution to a complex
       system   of   linear   equations    A(IA:IA+N-1,JA:JA+N-1)   *   X    =
       B(IB:IB+N-1,JB:JB+NRHS-1),

SYNOPSIS
       SUBROUTINE PCGESVX( FACT,  TRANS,  N,  NRHS, A, IA, JA, DESCA, AF, IAF,
                           JAF, DESCAF, IPIV, EQUED, R, C, B, IB,  JB,  DESCB,
                           X,  IX,  JX, DESCX, RCOND, FERR, BERR, WORK, LWORK,
                           RWORK, LRWORK, INFO )

           CHARACTER       EQUED, FACT, TRANS

           INTEGER         IA, IAF, IB, INFO, IX, JA,  JAF,  JB,  JX,  LRWORK,
                           LWORK, N, NRHS

           REAL            RCOND

           INTEGER         DESCA(  *  ),  DESCAF( * ), DESCB( * ), DESCX( * ),
                           IPIV( * )

           REAL            BERR( * ), C( * ), FERR( * ), R( * ), RWORK( * )

           COMPLEX         A( * ), AF( * ), B( * ), WORK( * ), X( * )

PURPOSE
       PCGESVX uses the LU factorization to compute the solution to a  complex
       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 pro-
       vided.


       Notes
       =====

       Each global data object is described by an associated description  vec-
       tor.  This vector stores the information required to establish the map-
       ping 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 col-
       umn.
       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
       In the  following  description,  A  denotes  A(IA:IA+N-1,JA:JA+N-1),  B
       denotes B(IB:IB+N-1,JB:JB+NRHS-1) and X denotes
       X(IX:IX+N-1,JX:JX+NRHS-1).

       The following steps are performed:

       1. If FACT = 'E', real scaling factors are computed to equilibrate
          the system:
             TRANS = 'N':  diag(R)*A*diag(C)     *inv(diag(C))*X = diag(R)*B
             TRANS = 'T': (diag(R)*A*diag(C))**T *inv(diag(R))*X = diag(C)*B
             TRANS = 'C': (diag(R)*A*diag(C))**H *inv(diag(R))*X = diag(C)*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(R)*A*diag(C) and B by diag(R)*B (if TRANS='N')
          or diag(C)*B (if TRANS = 'T' or 'C').

       2. If FACT = 'N' or 'E', the LU decomposition is used to factor the
          matrix A (after equilibration if FACT = 'E') as
             A = P * L * U,
          where P is a permutation matrix, L is a unit lower triangular
          matrix, and U is upper triangular.

       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 FACT = 'E' and equilibration was used, the matrix X is
          premultiplied by diag(C) (if TRANS = 'N') or diag(R) (if
          TRANS = 'T' or 'C') 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(IA:IA+N-1,JA:JA+N-1) is supplied on entry, and if not,
               whether the matrix  A(IA:IA+N-1,JA:JA+N-1)  should  be  equili-
               brated   before   it   is   factored.    =   'F':    On  entry,
               AF(IAF:IAF+N-1,JAF:JAF+N-1) and IPIV con-
               tain the factored form of A(IA:IA+N-1,JA:JA+N-1).  If EQUED  is
               not  'N',  the  matrix  A(IA:IA+N-1,JA:JA+N-1) has been equili-
               brated   with   scaling   factors   given   by   R    and    C.
               A(IA:IA+N-1,JA:JA+N-1),  AF(IAF:IAF+N-1,JAF:JAF+N-1),  and IPIV
               are not modified.  = 'N':   The  matrix  A(IA:IA+N-1,JA:JA+N-1)
               will be copied to
               AF(IAF:IAF+N-1,JAF:JAF+N-1) and factored.
               =  'E':   The  matrix  A(IA:IA+N-1,JA:JA+N-1)  will  be equili-
               brated if necessary, then copied to AF(IAF:IAF+N-1,JAF:JAF+N-1)
               and factored.

       TRANS   (global input) CHARACTER
               Specifies the form of the system of equations:
               = 'N':  A(IA:IA+N-1,JA:JA+N-1) * X(IX:IX+N-1,JX:JX+NRHS-1)
               = B(IB:IB+N-1,JB:JB+NRHS-1)     (No transpose)
               = 'T':  A(IA:IA+N-1,JA:JA+N-1)**T * X(IX:IX+N-1,JX:JX+NRHS-1)
               = B(IB:IB+N-1,JB:JB+NRHS-1)  (Transpose)
               = 'C':  A(IA:IA+N-1,JA:JA+N-1)**H * X(IX:IX+N-1,JX:JX+NRHS-1)
               = B(IB:IB+N-1,JB:JB+NRHS-1)  (Conjugate transpose)

       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(IB:IB+N-1,JB:JB+NRHS-1) and
               X(IX:IX+N-1,JX:JX+NRHS-1).  NRHS >= 0.

       A       (local input/local output) COMPLEX pointer into
               the   local   memory   to   an   array   of   local   dimension
               (LLD_A,LOCc(JA+N-1)).     On    entry,    the   N-by-N   matrix
               A(IA:IA+N-1,JA:JA+N-1).  If FACT = 'F' and EQUED is not 'N',
               then A(IA:IA+N-1,JA:JA+N-1) must have been equilibrated by
               the scaling factors in R and/or C.   A(IA:IA+N-1,JA:JA+N-1)  is
               not  modified if FACT = 'F' or  'N', or if FACT = 'E' and EQUED
               = 'N' on exit.

               On exit, if EQUED .ne. 'N', A(IA:IA+N-1,JA:JA+N-1) is scaled as
               follows:
               EQUED = 'R':  A(IA:IA+N-1,JA:JA+N-1) :=
               diag(R) * A(IA:IA+N-1,JA:JA+N-1)
               EQUED = 'C':  A(IA:IA+N-1,JA:JA+N-1) :=
               A(IA:IA+N-1,JA:JA+N-1) * diag(C)
               EQUED = 'B':  A(IA:IA+N-1,JA:JA+N-1) :=
               diag(R) * A(IA:IA+N-1,JA:JA+N-1) * diag(C).

       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) COMPLEX pointer
               into   the   local  memory  to  an  array  of  local  dimension
               (LLD_AF,LOCc(JA+N-1)).      If     FACT     =     'F',     then
               AF(IAF:IAF+N-1,JAF:JAF+N-1)  is  an input argument and on entry
               contains  the  factors  L  and   U   from   the   factorization
               A(IA:IA+N-1,JA:JA+N-1)  =  P*L*U  as  computed  by PCGETRF.  If
               EQUED .ne. 'N', then AF is the factored  form  of  the  equili-
               brated matrix A(IA:IA+N-1,JA:JA+N-1).

               If  FACT  =  'N', then AF(IAF:IAF+N-1,JAF:JAF+N-1) is an output
               argument and on exit returns the factors L and U from the  fac-
               torization A(IA:IA+N-1,JA:JA+N-1) = P*L*U of the original
               matrix A(IA:IA+N-1,JA:JA+N-1).

               If  FACT  =  'E', then AF(IAF:IAF+N-1,JAF:JAF+N-1) is an output
               argument and on exit returns the factors L and U from the  fac-
               torization A(IA:IA+N-1,JA:JA+N-1) = P*L*U of the equili-
               brated matrix A(IA:IA+N-1,JA:JA+N-1) (see the description of
               A(IA:IA+N-1,JA:JA+N-1)   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.

       IPIV    (local input or local output) INTEGER array, dimension
               LOCr(M_A)+MB_A. If FACT = 'F', then IPIV is an input argu- ment
               and  on  entry contains the pivot indices from the fac- toriza-
               tion A(IA:IA+N-1,JA:JA+N-1) = P*L*U  as  computed  by  PCGETRF;
               IPIV(i)  ->  The global row local row i was swapped with.  This
               array must be aligned with A( IA:IA+N-1, * ).

               If FACT = 'N', then IPIV is an output argument and on exit con-
               tains    the    pivot    indices    from    the   factorization
               A(IA:IA+N-1,JA:JA+N-1) = P*L*U of the original matrix
               A(IA:IA+N-1,JA:JA+N-1).

               If FACT = 'E', then IPIV is an output argument and on exit con-
               tains    the    pivot    indices    from    the   factorization
               A(IA:IA+N-1,JA:JA+N-1) = P*L*U of the equilibrated matrix
               A(IA:IA+N-1,JA:JA+N-1).

       EQUED   (global input or global output) CHARACTER
               Specifies the form of equilibration that was done.  = 'N':   No
               equilibration (always true if FACT = 'N').
               =  'R':   Row  equilibration,  i.e., A(IA:IA+N-1,JA:JA+N-1) has
               been premultiplied by diag(R).  = 'C':   Column  equilibration,
               i.e.,   A(IA:IA+N-1,JA:JA+N-1)   has   been  postmultiplied  by
               diag(C).  = 'B':  Both row and column equilibration, i.e.,
               A(IA:IA+N-1,JA:JA+N-1) has been replaced by
               diag(R) * A(IA:IA+N-1,JA:JA+N-1) * diag(C).  EQUED is an  input
               variable if FACT = 'F'; otherwise, it is an output variable.

       R       (local input or local output) REAL array,
               dimension    LOCr(M_A).     The    row    scale   factors   for
               A(IA:IA+N-1,JA:JA+N-1).
               If EQUED = 'R' or 'B', A(IA:IA+N-1,JA:JA+N-1) is multiplied  on
               the  left by diag(R); if EQUED='N' or 'C', R is not acces- sed.
               R is an input variable if FACT = 'F'; otherwise, R is an output
               variable.   If  FACT = 'F' and EQUED = 'R' or 'B', each element
               of R must be  positive.   R  is  replicated  in  every  process
               column, and is aligned with the distributed matrix A.

       C       (local input or local output) REAL array,
               dimension    LOCc(N_A).    The   column   scale   factors   for
               A(IA:IA+N-1,JA:JA+N-1).
               If EQUED = 'C' or 'B', A(IA:IA+N-1,JA:JA+N-1) is multiplied  on
               the right by diag(C); if EQUED = 'N' or 'R', C is not accessed.
               C is an input variable if FACT = 'F'; otherwise, C is an output
               variable.   If FACT = 'F' and EQUED = 'C' or C is replicated in
               every process row, and is aligned with the  distributed  matrix
               A.

       B       (local input/local output) COMPLEX 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(IB:IB+N-1,JB:JB+NRHS-1). On exit, if
               EQUED  =  'N',  B(IB:IB+N-1,JB:JB+NRHS-1)  is  not modified; if
               TRANS = 'N' and EQUED  =  'R'  or  'B',  B  is  overwritten  by
               diag(R)*B(IB:IB+N-1,JB:JB+NRHS-1); if TRANS = 'T' or 'C'
               and EQUED = 'C' or 'B', B(IB:IB+N-1,JB:JB+NRHS-1) is over-
               written by diag(C)*B(IB:IB+N-1,JB:JB+NRHS-1).

       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) COMPLEX 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(IX:IX+N-1,JX:JX+NRHS-1) to the original
               system of equations.  Note that A(IA:IA+N-1,JA:JA+N-1) and
               B(IB:IB+N-1,JB:JB+NRHS-1)  are  modified  on exit if EQUED .ne.
               'N',  and  the  solution  to   the   equilibrated   system   is
               inv(diag(C))*X(IX:IX+N-1,JX:JX+NRHS-1) if TRANS = 'N' and EQUED
               = 'C'  or  'B',  or  inv(diag(R))*X(IX:IX+N-1,JX:JX+NRHS-1)  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(IA:IA+N-1,JA:JA+N-1) 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  condi-
               tion is indicated by a return code of INFO > 0.

       FERR    (local output) REAL array, dimension LOCc(N_B)
               The  estimated  forward  error  bounds for each solution vector
               X(j)   (the   j-th    column    of    the    solution    matrix
               X(IX:IX+N-1,JX:JX+NRHS-1).  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 estimate is as reliable as the estimate for RCOND,  and  is
               almost always a slight overestimate of the true error.  FERR is
               replicated in every process row, and is aligned with the matri-
               ces B and X.

       BERR    (local output) REAL array, dimension LOCc(N_B).
               The componentwise relative backward error of each solution vec-
               tor X(j) (i.e., the smallest relative change in  any  entry  of
               A(IA:IA+N-1,JA:JA+N-1) or
               B(IB:IB+N-1,JB:JB+NRHS-1)  that  makes X(j) an exact solution).
               BERR is replicated in every process row, and  is  aligned  with
               the matrices B and X.

       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 = MAX( PCGECON( LWORK ), PCGERFS( LWORK ) ) +
               LOCr( N_A ).

               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_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 INFO = -i, the i-th argument had an illegal value
               > 0:  if INFO = i, and i is
               <= N:  U(IA+I-1,IA+I-1) is exactly zero.  The factorization has
               been  completed,  but  the factor U is exactly singular, so 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.



ScaLAPACK routine               31 October 2017                     PCGESVX(3)