Go to main content
Oracle Developer Studio 12.5 Man Pages

Exit Print View

Updated: June 2017
 
 

cgbsvx (3p)

Name

cgbsvx - use the LU factorization to compute the solution to a complex system of linear equations A*X=B, A**T*X=B, or A**H*X=B, where A is a band matrix

Synopsis

SUBROUTINE CGBSVX(FACT, TRANSA, N, KL, KU, NRHS, A, LDA, AF,
LDAF, IPIVOT, EQUED, R, C, B, LDB, X, LDX, RCOND, FERR,
BERR, WORK, WORK2, INFO)

CHARACTER*1 FACT, TRANSA, EQUED
COMPLEX A(LDA,*), AF(LDAF,*), B(LDB,*), X(LDX,*), WORK(*)
INTEGER N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
INTEGER IPIVOT(*)
REAL RCOND
REAL R(*), C(*), FERR(*), BERR(*), WORK2(*)

SUBROUTINE CGBSVX_64(FACT, TRANSA, N, KL, KU, NRHS, A, LDA, AF,
LDAF, IPIVOT, EQUED, R, C, B, LDB, X, LDX, RCOND, FERR,
BERR, WORK, WORK2, INFO)

CHARACTER*1 FACT, TRANSA, EQUED
COMPLEX A(LDA,*), AF(LDAF,*), B(LDB,*), X(LDX,*), WORK(*)
INTEGER*8 N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
INTEGER*8 IPIVOT(*)
REAL RCOND
REAL R(*), C(*), FERR(*), BERR(*), WORK2(*)




F95 INTERFACE
SUBROUTINE GBSVX(FACT, TRANSA, N, KL, KU, NRHS, A, LDA,
AF, LDAF, IPIVOT, EQUED, R, C, B, LDB, X, LDX,
RCOND, FERR, BERR, WORK, WORK2, INFO)

CHARACTER(LEN=1) :: FACT, TRANSA, EQUED
COMPLEX, DIMENSION(:) :: WORK
COMPLEX, DIMENSION(:,:) :: A, AF, B, X
INTEGER :: N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
INTEGER, DIMENSION(:) :: IPIVOT
REAL :: RCOND
REAL, DIMENSION(:) :: R, C, FERR, BERR, WORK2

SUBROUTINE GBSVX_64(FACT, TRANSA, N, KL, KU, NRHS, A,
LDA, AF, LDAF, IPIVOT, EQUED, R, C, B, LDB, X, LDX,
RCOND, FERR, BERR, WORK, WORK2, INFO)

CHARACTER(LEN=1) :: FACT, TRANSA, EQUED
COMPLEX, DIMENSION(:) :: WORK
COMPLEX, DIMENSION(:,:) :: A, AF, B, X
INTEGER(8) :: N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
INTEGER(8), DIMENSION(:) :: IPIVOT
REAL :: RCOND
REAL, DIMENSION(:) :: R, C, FERR, BERR, WORK2




C INTERFACE
#include <sunperf.h>

void  cgbsvx(char  fact,  char transa, int n, int kl, int ku, int nrhs,
complex *a, int lda, complex *af, int ldaf, int *ipivot, char
*equed,  float *r, float *c, complex *b, int ldb, complex *x,
int ldx,  float  *rcond,  float  *ferr,  float  *berr,  float
*work2, int *info);

void  cgbsvx_64(char  fact, char transa, long n, long kl, long ku, long
nrhs, complex *a, long lda,  complex  *af,  long  ldaf,  long
*ipivot,  char  *equed,  float *r, float *c, complex *b, long
ldb, complex *x, long ldx, float *rcond, float  *ferr,  float
*berr, float *work2, long *info);

Description

Oracle Solaris Studio Performance Library                           cgbsvx(3P)



NAME
       cgbsvx  - use the LU factorization to compute the solution to a complex
       system of linear equations A*X=B, A**T*X=B, or A**H*X=B, where A  is  a
       band matrix


SYNOPSIS
       SUBROUTINE CGBSVX(FACT, TRANSA, N, KL, KU, NRHS, A, LDA, AF,
             LDAF, IPIVOT, EQUED, R, C, B, LDB, X, LDX, RCOND, FERR,
             BERR, WORK, WORK2, INFO)

       CHARACTER*1 FACT, TRANSA, EQUED
       COMPLEX A(LDA,*), AF(LDAF,*), B(LDB,*), X(LDX,*), WORK(*)
       INTEGER N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
       INTEGER IPIVOT(*)
       REAL RCOND
       REAL R(*), C(*), FERR(*), BERR(*), WORK2(*)

       SUBROUTINE CGBSVX_64(FACT, TRANSA, N, KL, KU, NRHS, A, LDA, AF,
             LDAF, IPIVOT, EQUED, R, C, B, LDB, X, LDX, RCOND, FERR,
             BERR, WORK, WORK2, INFO)

       CHARACTER*1 FACT, TRANSA, EQUED
       COMPLEX A(LDA,*), AF(LDAF,*), B(LDB,*), X(LDX,*), WORK(*)
       INTEGER*8 N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
       INTEGER*8 IPIVOT(*)
       REAL RCOND
       REAL R(*), C(*), FERR(*), BERR(*), WORK2(*)




   F95 INTERFACE
       SUBROUTINE GBSVX(FACT, TRANSA, N, KL, KU, NRHS, A, LDA,
              AF, LDAF, IPIVOT, EQUED, R, C, B, LDB, X, LDX,
              RCOND, FERR, BERR, WORK, WORK2, INFO)

       CHARACTER(LEN=1) :: FACT, TRANSA, EQUED
       COMPLEX, DIMENSION(:) :: WORK
       COMPLEX, DIMENSION(:,:) :: A, AF, B, X
       INTEGER :: N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
       INTEGER, DIMENSION(:) :: IPIVOT
       REAL :: RCOND
       REAL, DIMENSION(:) :: R, C, FERR, BERR, WORK2

       SUBROUTINE GBSVX_64(FACT, TRANSA, N, KL, KU, NRHS, A,
              LDA, AF, LDAF, IPIVOT, EQUED, R, C, B, LDB, X, LDX,
              RCOND, FERR, BERR, WORK, WORK2, INFO)

       CHARACTER(LEN=1) :: FACT, TRANSA, EQUED
       COMPLEX, DIMENSION(:) :: WORK
       COMPLEX, DIMENSION(:,:) :: A, AF, B, X
       INTEGER(8) :: N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
       INTEGER(8), DIMENSION(:) :: IPIVOT
       REAL :: RCOND
       REAL, DIMENSION(:) :: R, C, FERR, BERR, WORK2




   C INTERFACE
       #include <sunperf.h>

       void  cgbsvx(char  fact,  char transa, int n, int kl, int ku, int nrhs,
                 complex *a, int lda, complex *af, int ldaf, int *ipivot, char
                 *equed,  float *r, float *c, complex *b, int ldb, complex *x,
                 int ldx,  float  *rcond,  float  *ferr,  float  *berr,  float
                 *work2, int *info);

       void  cgbsvx_64(char  fact, char transa, long n, long kl, long ku, long
                 nrhs, complex *a, long lda,  complex  *af,  long  ldaf,  long
                 *ipivot,  char  *equed,  float *r, float *c, complex *b, long
                 ldb, complex *x, long ldx, float *rcond, float  *ferr,  float
                 *berr, float *work2, long *info);



PURPOSE
       cgbsvx  uses  the LU factorization to compute the solution to a complex
       system of linear equations A*X=B, A**T *X=B, or A**H *X=B, where A is a
       band  matrix of order N with KL subdiagonals and KU superdiagonals, and
       X and B are N-by-NRHS matrices.

       Error bounds on the solution and a condition  estimate  are  also  pro-
       vided.

       The following steps are performed by this subroutine:

       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 = L * U,
       where L is a product of permutation and unit lower triangular  matrices
       with  KL subdiagonals, and U is upper triangular with KL+KU superdiago-
       nals.

       3. If some U(i,i)=0, so that U is exactly singular,  then  the  routine
       returns  with  INFO  =  i. Otherwise, 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, INFO = N+1 is
       returned as a warning, but the routine still goes on to solve for X and
       compute error bounds as described below.

       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(C)
       (if  TRANS  = 'N') or diag(R) (if TRANS = 'T' or 'C') so that it solves
       the original system before equilibration.


ARGUMENTS
       FACT (input)
                 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 and IPIVOT contain the factored form  of
                 A.   If  EQUED is not 'N', the matrix A has been equilibrated
                 with scaling factors given by R and C.  A, AF, and IPIVOT are
                 not 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.


       TRANSA (input)
                 Specifies the form of the system of equations.  = 'N':  A * X
                 = B     (No transpose)
                 = 'T':  A**T * X = B  (Transpose)
                 = 'C':  A**H * X = B  (Conjugate transpose)


       N (input) The number of linear equations, i.e., the order of the matrix
                 A.  N >= 0.


       KL (input)
                 The number of subdiagonals within the band of A.  KL >= 0.


       KU (input)
                 The number of superdiagonals within the band of A.  KU >=  0.


       NRHS (input)
                 The  number  of right hand sides, i.e., the number of columns
                 of the matrices B and X.  NRHS >= 0.


       A (input/output)
                 COMPLEX array, dimension (LDA,N) On entry, the  matrix  A  in
                 band  storage, in rows 1 to KL+KU+1.  The j-th column of A is
                 stored in the j-th column of the array A as follows:
                 A(KU+1+i-j,j)=A(i,j) for max(1,j-KU)<=i<=min(N,j+kl)
                 If FACT = 'F' and EQUED is not 'N', then  A  must  have  been
                 equilibrated  by  the scaling factors in R and/or C. A is not
                 modified if FACT = 'F' or 'N', or if FACT = 'E' and  EQUED  =
                 'N' on exit.
                 On exit, if EQUED .ne. 'N', A is scaled as follows:
                 EQUED = 'R':  A := diag(R) * A
                 EQUED = 'C':  A := A * diag(C)
                 EQUED = 'B':  A := diag(R) * A * diag(C).


       LDA (input)
                 The leading dimension of the array A.
                 LDA >= KL+KU+1.


       AF (input or output)
                 COMPLEX array, dimension (LDAF,N)
                 If FACT = 'F', then AF is an input argument and on entry con-
                 tains details of the LU factorization of the band  matrix  A,
                 as  computed  by  CGBTRF.  U is stored as an upper triangular
                 band matrix with KL+KU superdiagonals in rows 1  to  KL+KU+1,
                 and  the multipliers used during the factorization are stored
                 in rows KL+KU+2 to 2*KL+KU+1. If EQUED .ne. 'N', then  AF  is
                 the factored form of the equilibrated matrix A.
                 If  FACT  =  'N',  then  AF is an output argument and on exit
                 returns details of the LU factorization of A.
                 If FACT = 'E', then AF is an  output  argument  and  on  exit
                 returns  details  of the LU factorization of the equilibrated
                 matrix A (see the description of A for the form of the  equi-
                 librated matrix).


       LDAF (input)
                 The leading dimension of the array AF.
                 LDAF >= 2*KL+KU+1.


       IPIVOT (input or output)
                 INTEGER array, dimension (N) If FACT = 'F', then IPIVOT is an
                 input argument and on entry contains the pivot  indices  from
                 the  factorization  A=L*U as computed by CGBTRF; row i of the
                 matrix was interchanged with row IPIVOT(i).
                 If FACT = 'N', then IPIVOT is an output argument and on  exit
                 contains  the pivot indices from the factorization A = L*U of
                 the original matrix A.
                 If FACT = 'E', then IPIVOT is an output argument and on  exit
                 contains  the pivot indices from the factorization A = L*U of
                 the equilibrated matrix A.


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


       R (input or output)
                 REAL array, dimension (N) The row scale  factors  for  A.  If
                 EQUED  =  'R' or 'B', A is multiplied on the left by diag(R);
                 if EQUED = 'N' or 'C', R is not accessed. R is an input argu-
                 ment  if  FACT  = 'F'; otherwise, R is an output argument. If
                 FACT = 'F' and EQUED = 'R' or 'B', each element of R must  be
                 positive.


       C (input or output)
                 REAL  array, dimension (N) The column scale factors for A. If
                 EQUED = 'C' or 'B', A is multiplied on the right by  diag(C);
                 if EQUED = 'N' or 'R', C is not accessed. C is an input argu-
                 ment if FACT = 'F'; otherwise, C is an  output  argument.  If
                 FACT  = 'F' and EQUED = 'C' or 'B', each element of C must be
                 positive.


       B (input/output)
                 COMPLEX array, dimension (LDB,NRHS) On entry, the right  hand
                 side matrix B.
                 On exit,
                 if EQUED = 'N', B is not modified;
                 if  TRANSA  = 'N' and EQUED = 'R' or 'B', B is overwritten by
                 diag(R)*B;
                 if TRANSA = 'T' or 'C' and EQUED = 'C' or 'B', B is overwrit-
                 ten by diag(C)*B.


       LDB (input)
                 The leading dimension of the array B.
                 LDB >= max(1,N).


       X (output)
                 COMPLEX  array,  dimension  (LDX,NRHS)  If INFO = 0 or INFO =
                 N+1, 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(C))*X  if  TRANSA  =  'N' and EQUED = 'C' or 'B', or
                 inv(diag(R))*X if TRANSA = 'T' or 'C' and EQUED = 'R' or 'B'.


       LDX (input)
                 The leading dimension of the array X.
                 LDX >= max(1,N).


       RCOND (output)
                 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.


       FERR (output)
                 The  estimated  forward  error bound for each solution vector
                 X(j) (the j-th column of the solution matrix X).  If XTRUE is
                 the  true solution corresponding to X(j), FERR(j) is an esti-
                 mated upper bound for the magnitude of the largest element in
                 (X(j)  -  XTRUE) divided by the magnitude of the largest ele-
                 ment in X(j). The estimate is as reliable as the estimate for
                 RCOND, and is almost always a slight overestimate of the true
                 error.


       BERR (output)
                 The componentwise relative backward error  of  each  solution
                 vector  X(j)  (i.e., the smallest relative change in any ele-
                 ment of A or B that makes X(j) an exact solution).


       WORK (workspace)
                 WORK is COMPLEX array, dimension(2*N)


       WORK2 (output)
                 WORK2 is REAL array, dimension(N) On exit, WORK2(1)  contains
                 the  reciprocal pivot growth factor norm(A)/norm(U). The "max
                 absolute element" norm is used. If WORK2(1) is much less than
                 1, then the stability of the LU factorization of the (equili-
                 brated) matrix A could be poor.  This  also  means  that  the
                 solution  X,  condition  estimator  RCOND,  and forward error
                 bound FERR could be unreliable. If factorization  fails  with
                 0<INFO<=N, then WORK2(1) contains the reciprocal pivot growth
                 factor for the leading INFO columns of A.


       INFO (output)
                 = 0:  successful exit;
                 < 0:  if INFO = -i, the i-th argument had an illegal value;
                 > 0:  if INFO = i, and i is
                 <= N:  U(i,i) is exactly zero.  The  factorization  has  been
                 completed, but the factor U is exactly singular, so the solu-
                 tion and error bounds could not be computed.  RCOND  =  0  is
                 returned.
                 =  N+1: U is nonsingular, but RCOND is less than machine pre-
                 cision, meaning that the matrix is singular to working preci-
                 sion.  Nevertheless,  the  solution and error bounds are com-
                 puted because there are a number of situations where the com-
                 puted  solution  can be more accurate than the value of RCOND
                 would suggest.




                                  7 Nov 2015                        cgbsvx(3P)