Go to main content
Oracle Developer Studio 12.5 Man Pages

Exit Print View

Updated: June 2017
 
 

dgbsvx (3p)

Name

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

Synopsis

SUBROUTINE DGBSVX(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
INTEGER N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
INTEGER IPIVOT(*), WORK2(*)
DOUBLE PRECISION RCOND
DOUBLE  PRECISION A(LDA,*), AF(LDAF,*), R(*), C(*), B(LDB,*), X(LDX,*),
FERR(*), BERR(*), WORK(*)

SUBROUTINE DGBSVX_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
INTEGER*8 N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
INTEGER*8 IPIVOT(*), WORK2(*)
DOUBLE PRECISION RCOND
DOUBLE PRECISION A(LDA,*), AF(LDAF,*), R(*), C(*), B(LDB,*),  X(LDX,*),
FERR(*), BERR(*), WORK(*)




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
INTEGER :: N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
INTEGER, DIMENSION(:) :: IPIVOT, WORK2
REAL(8) :: RCOND
REAL(8), DIMENSION(:) :: R, C, FERR, BERR, WORK
REAL(8), DIMENSION(:,:) :: A, AF, B, X

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
INTEGER(8) :: N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
INTEGER(8), DIMENSION(:) :: IPIVOT, WORK2
REAL(8) :: RCOND
REAL(8), DIMENSION(:) :: R, C, FERR, BERR, WORK
REAL(8), DIMENSION(:,:) :: A, AF, B, X




C INTERFACE
#include <sunperf.h>

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

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

Description

Oracle Solaris Studio Performance Library                           dgbsvx(3P)



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


SYNOPSIS
       SUBROUTINE DGBSVX(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
       INTEGER N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
       INTEGER IPIVOT(*), WORK2(*)
       DOUBLE PRECISION RCOND
       DOUBLE  PRECISION A(LDA,*), AF(LDAF,*), R(*), C(*), B(LDB,*), X(LDX,*),
       FERR(*), BERR(*), WORK(*)

       SUBROUTINE DGBSVX_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
       INTEGER*8 N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
       INTEGER*8 IPIVOT(*), WORK2(*)
       DOUBLE PRECISION RCOND
       DOUBLE PRECISION A(LDA,*), AF(LDAF,*), R(*), C(*), B(LDB,*),  X(LDX,*),
       FERR(*), BERR(*), WORK(*)




   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
       INTEGER :: N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
       INTEGER, DIMENSION(:) :: IPIVOT, WORK2
       REAL(8) :: RCOND
       REAL(8), DIMENSION(:) :: R, C, FERR, BERR, WORK
       REAL(8), DIMENSION(:,:) :: A, AF, B, X

       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
       INTEGER(8) :: N, KL, KU, NRHS, LDA, LDAF, LDB, LDX, INFO
       INTEGER(8), DIMENSION(:) :: IPIVOT, WORK2
       REAL(8) :: RCOND
       REAL(8), DIMENSION(:) :: R, C, FERR, BERR, WORK
       REAL(8), DIMENSION(:,:) :: A, AF, B, X




   C INTERFACE
       #include <sunperf.h>

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

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



PURPOSE
       dgbsvx uses the LU factorization to compute the solution to a real sys-
       tem 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 superdiago-
       nals, 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  (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)
                 DOUBLE PRECISION  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)
                 DOUBLE PRECISION 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  DGBTRF.  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 DGBTRF; 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 argument if
                 FACT = 'F'; otherwise, it is an output argument.


       R (input or output)
                 DOUBLE  PRECISION  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 argument if FACT = 'F'; otherwise, R is an output argu-
                 ment. If FACT = 'F' and EQUED = 'R' or 'B', each element of R
                 must be positive.


       C (input or output)
                 DOUBLE PRECISION array, dimension (N) The column  scale  fac-
                 tors  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 argument 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)
                 DOUBLE  PRECISION  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)
                 DOUBLE PRECISION 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 sys-
                 tem 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  element
                 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 (output)
                 WORK  is  DOUBLE  PRECISION  array,  dimension (3*N) On exit,
                 WORK(1)  contains  the   reciprocal   pivot   growth   factor
                 norm(A)/norm(U).  The "max absolute element" norm is used. If
                 WORK(1) is much less than 1, then the  stability  of  the  LU
                 factorization  of  the (equilibrated) 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 WORK(1) contains the
                 reciprocal  pivot  growth factor for the leading INFO columns
                 of A.


       WORK2 (workspace)
                 WORK2  is INTEGER array, dimension(N)


       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                        dgbsvx(3P)