001:       SUBROUTINE CGBSVXX( FACT, TRANS, N, KL, KU, NRHS, AB, LDAB, AFB,
002:      $                    LDAFB, IPIV, EQUED, R, C, B, LDB, X, LDX,
003:      $                    RCOND, RPVGRW, BERR, N_ERR_BNDS,
004:      $                    ERR_BNDS_NORM, ERR_BNDS_COMP, NPARAMS, PARAMS,
005:      $                    WORK, RWORK, INFO )
006: *
007: *     -- LAPACK driver routine (version 3.2)                          --
008: *     -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
009: *     -- Jason Riedy of Univ. of California Berkeley.                 --
010: *     -- November 2008                                                --
011: *
012: *     -- LAPACK is a software package provided by Univ. of Tennessee, --
013: *     -- Univ. of California Berkeley and NAG Ltd.                    --
014: *
015:       IMPLICIT NONE
016: *     ..
017: *     .. Scalar Arguments ..
018:       CHARACTER          EQUED, FACT, TRANS
019:       INTEGER            INFO, LDAB, LDAFB, LDB, LDX, N, NRHS, NPARAMS,
020:      $                   N_ERR_BNDS
021:       REAL               RCOND, RPVGRW
022: *     ..
023: *     .. Array Arguments ..
024:       INTEGER            IPIV( * )
025:       COMPLEX            AB( LDAB, * ), AFB( LDAFB, * ), B( LDB, * ),
026:      $                   X( LDX , * ),WORK( * )
027:       REAL               R( * ), C( * ), PARAMS( * ), BERR( * ),
028:      $                   ERR_BNDS_NORM( NRHS, * ),
029:      $                   ERR_BNDS_COMP( NRHS, * ), RWORK( * )
030: *     ..
031: *
032: *     Purpose
033: *     =======
034: *
035: *     CGBSVXX uses the LU factorization to compute the solution to a
036: *     complex system of linear equations  A * X = B,  where A is an
037: *     N-by-N matrix and X and B are N-by-NRHS matrices.
038: *
039: *     If requested, both normwise and maximum componentwise error bounds
040: *     are returned. CGBSVXX will return a solution with a tiny
041: *     guaranteed error (O(eps) where eps is the working machine
042: *     precision) unless the matrix is very ill-conditioned, in which
043: *     case a warning is returned. Relevant condition numbers also are
044: *     calculated and returned.
045: *
046: *     CGBSVXX accepts user-provided factorizations and equilibration
047: *     factors; see the definitions of the FACT and EQUED options.
048: *     Solving with refinement and using a factorization from a previous
049: *     CGBSVXX call will also produce a solution with either O(eps)
050: *     errors or warnings, but we cannot make that claim for general
051: *     user-provided factorizations and equilibration factors if they
052: *     differ from what CGBSVXX would itself produce.
053: *
054: *     Description
055: *     ===========
056: *
057: *     The following steps are performed:
058: *
059: *     1. If FACT = 'E', real scaling factors are computed to equilibrate
060: *     the system:
061: *
062: *       TRANS = 'N':  diag(R)*A*diag(C)     *inv(diag(C))*X = diag(R)*B
063: *       TRANS = 'T': (diag(R)*A*diag(C))**T *inv(diag(R))*X = diag(C)*B
064: *       TRANS = 'C': (diag(R)*A*diag(C))**H *inv(diag(R))*X = diag(C)*B
065: *
066: *     Whether or not the system will be equilibrated depends on the
067: *     scaling of the matrix A, but if equilibration is used, A is
068: *     overwritten by diag(R)*A*diag(C) and B by diag(R)*B (if TRANS='N')
069: *     or diag(C)*B (if TRANS = 'T' or 'C').
070: *
071: *     2. If FACT = 'N' or 'E', the LU decomposition is used to factor
072: *     the matrix A (after equilibration if FACT = 'E') as
073: *
074: *       A = P * L * U,
075: *
076: *     where P is a permutation matrix, L is a unit lower triangular
077: *     matrix, and U is upper triangular.
078: *
079: *     3. If some U(i,i)=0, so that U is exactly singular, then the
080: *     routine returns with INFO = i. Otherwise, the factored form of A
081: *     is used to estimate the condition number of the matrix A (see
082: *     argument RCOND). If the reciprocal of the condition number is less
083: *     than machine precision, the routine still goes on to solve for X
084: *     and compute error bounds as described below.
085: *
086: *     4. The system of equations is solved for X using the factored form
087: *     of A.
088: *
089: *     5. By default (unless PARAMS(LA_LINRX_ITREF_I) is set to zero),
090: *     the routine will use iterative refinement to try to get a small
091: *     error and error bounds.  Refinement calculates the residual to at
092: *     least twice the working precision.
093: *
094: *     6. If equilibration was used, the matrix X is premultiplied by
095: *     diag(C) (if TRANS = 'N') or diag(R) (if TRANS = 'T' or 'C') so
096: *     that it solves the original system before equilibration.
097: *
098: *     Arguments
099: *     =========
100: *
101: *     Some optional parameters are bundled in the PARAMS array.  These
102: *     settings determine how refinement is performed, but often the
103: *     defaults are acceptable.  If the defaults are acceptable, users
104: *     can pass NPARAMS = 0 which prevents the source code from accessing
105: *     the PARAMS argument.
106: *
107: *     FACT    (input) CHARACTER*1
108: *     Specifies whether or not the factored form of the matrix A is
109: *     supplied on entry, and if not, whether the matrix A should be
110: *     equilibrated before it is factored.
111: *       = 'F':  On entry, AF and IPIV contain the factored form of A.
112: *               If EQUED is not 'N', the matrix A has been
113: *               equilibrated with scaling factors given by R and C.
114: *               A, AF, and IPIV are not modified.
115: *       = 'N':  The matrix A will be copied to AF and factored.
116: *       = 'E':  The matrix A will be equilibrated if necessary, then
117: *               copied to AF and factored.
118: *
119: *     TRANS   (input) CHARACTER*1
120: *     Specifies the form of the system of equations:
121: *       = 'N':  A * X = B     (No transpose)
122: *       = 'T':  A**T * X = B  (Transpose)
123: *       = 'C':  A**H * X = B  (Conjugate Transpose = Transpose)
124: *
125: *     N       (input) INTEGER
126: *     The number of linear equations, i.e., the order of the
127: *     matrix A.  N >= 0.
128: *
129: *     KL      (input) INTEGER
130: *     The number of subdiagonals within the band of A.  KL >= 0.
131: *
132: *     KU      (input) INTEGER
133: *     The number of superdiagonals within the band of A.  KU >= 0.
134: *
135: *     NRHS    (input) INTEGER
136: *     The number of right hand sides, i.e., the number of columns
137: *     of the matrices B and X.  NRHS >= 0.
138: *
139: *     AB      (input/output) REAL array, dimension (LDAB,N)
140: *     On entry, the matrix A in band storage, in rows 1 to KL+KU+1.
141: *     The j-th column of A is stored in the j-th column of the
142: *     array AB as follows:
143: *     AB(KU+1+i-j,j) = A(i,j) for max(1,j-KU)<=i<=min(N,j+kl)
144: *
145: *     If FACT = 'F' and EQUED is not 'N', then AB must have been
146: *     equilibrated by the scaling factors in R and/or C.  AB is not
147: *     modified if FACT = 'F' or 'N', or if FACT = 'E' and
148: *     EQUED = 'N' on exit.
149: *
150: *     On exit, if EQUED .ne. 'N', A is scaled as follows:
151: *     EQUED = 'R':  A := diag(R) * A
152: *     EQUED = 'C':  A := A * diag(C)
153: *     EQUED = 'B':  A := diag(R) * A * diag(C).
154: *
155: *     LDAB    (input) INTEGER
156: *     The leading dimension of the array AB.  LDAB >= KL+KU+1.
157: *
158: *     AFB     (input or output) REAL array, dimension (LDAFB,N)
159: *     If FACT = 'F', then AFB is an input argument and on entry
160: *     contains details of the LU factorization of the band matrix
161: *     A, as computed by CGBTRF.  U is stored as an upper triangular
162: *     band matrix with KL+KU superdiagonals in rows 1 to KL+KU+1,
163: *     and the multipliers used during the factorization are stored
164: *     in rows KL+KU+2 to 2*KL+KU+1.  If EQUED .ne. 'N', then AFB is
165: *     the factored form of the equilibrated matrix A.
166: *
167: *     If FACT = 'N', then AF is an output argument and on exit
168: *     returns the factors L and U from the factorization A = P*L*U
169: *     of the original matrix A.
170: *
171: *     If FACT = 'E', then AF is an output argument and on exit
172: *     returns the factors L and U from the factorization A = P*L*U
173: *     of the equilibrated matrix A (see the description of A for
174: *     the form of the equilibrated matrix).
175: *
176: *     LDAFB   (input) INTEGER
177: *     The leading dimension of the array AFB.  LDAFB >= 2*KL+KU+1.
178: *
179: *     IPIV    (input or output) INTEGER array, dimension (N)
180: *     If FACT = 'F', then IPIV is an input argument and on entry
181: *     contains the pivot indices from the factorization A = P*L*U
182: *     as computed by SGETRF; row i of the matrix was interchanged
183: *     with row IPIV(i).
184: *
185: *     If FACT = 'N', then IPIV is an output argument and on exit
186: *     contains the pivot indices from the factorization A = P*L*U
187: *     of the original matrix A.
188: *
189: *     If FACT = 'E', then IPIV is an output argument and on exit
190: *     contains the pivot indices from the factorization A = P*L*U
191: *     of the equilibrated matrix A.
192: *
193: *     EQUED   (input or output) CHARACTER*1
194: *     Specifies the form of equilibration that was done.
195: *       = 'N':  No equilibration (always true if FACT = 'N').
196: *       = 'R':  Row equilibration, i.e., A has been premultiplied by
197: *               diag(R).
198: *       = 'C':  Column equilibration, i.e., A has been postmultiplied
199: *               by diag(C).
200: *       = 'B':  Both row and column equilibration, i.e., A has been
201: *               replaced by diag(R) * A * diag(C).
202: *     EQUED is an input argument if FACT = 'F'; otherwise, it is an
203: *     output argument.
204: *
205: *     R       (input or output) REAL array, dimension (N)
206: *     The row scale factors for A.  If EQUED = 'R' or 'B', A is
207: *     multiplied on the left by diag(R); if EQUED = 'N' or 'C', R
208: *     is not accessed.  R is an input argument if FACT = 'F';
209: *     otherwise, R is an output argument.  If FACT = 'F' and
210: *     EQUED = 'R' or 'B', each element of R must be positive.
211: *     If R is output, each element of R is a power of the radix.
212: *     If R is input, each element of R should be a power of the radix
213: *     to ensure a reliable solution and error estimates. Scaling by
214: *     powers of the radix does not cause rounding errors unless the
215: *     result underflows or overflows. Rounding errors during scaling
216: *     lead to refining with a matrix that is not equivalent to the
217: *     input matrix, producing error estimates that may not be
218: *     reliable.
219: *
220: *     C       (input or output) REAL array, dimension (N)
221: *     The column scale factors for A.  If EQUED = 'C' or 'B', A is
222: *     multiplied on the right by diag(C); if EQUED = 'N' or 'R', C
223: *     is not accessed.  C is an input argument if FACT = 'F';
224: *     otherwise, C is an output argument.  If FACT = 'F' and
225: *     EQUED = 'C' or 'B', each element of C must be positive.
226: *     If C is output, each element of C is a power of the radix.
227: *     If C is input, each element of C should be a power of the radix
228: *     to ensure a reliable solution and error estimates. Scaling by
229: *     powers of the radix does not cause rounding errors unless the
230: *     result underflows or overflows. Rounding errors during scaling
231: *     lead to refining with a matrix that is not equivalent to the
232: *     input matrix, producing error estimates that may not be
233: *     reliable.
234: *
235: *     B       (input/output) REAL array, dimension (LDB,NRHS)
236: *     On entry, the N-by-NRHS right hand side matrix B.
237: *     On exit,
238: *     if EQUED = 'N', B is not modified;
239: *     if TRANS = 'N' and EQUED = 'R' or 'B', B is overwritten by
240: *        diag(R)*B;
241: *     if TRANS = 'T' or 'C' and EQUED = 'C' or 'B', B is
242: *        overwritten by diag(C)*B.
243: *
244: *     LDB     (input) INTEGER
245: *     The leading dimension of the array B.  LDB >= max(1,N).
246: *
247: *     X       (output) REAL array, dimension (LDX,NRHS)
248: *     If INFO = 0, the N-by-NRHS solution matrix X to the original
249: *     system of equations.  Note that A and B are modified on exit
250: *     if EQUED .ne. 'N', and the solution to the equilibrated system is
251: *     inv(diag(C))*X if TRANS = 'N' and EQUED = 'C' or 'B', or
252: *     inv(diag(R))*X if TRANS = 'T' or 'C' and EQUED = 'R' or 'B'.
253: *
254: *     LDX     (input) INTEGER
255: *     The leading dimension of the array X.  LDX >= max(1,N).
256: *
257: *     RCOND   (output) REAL
258: *     Reciprocal scaled condition number.  This is an estimate of the
259: *     reciprocal Skeel condition number of the matrix A after
260: *     equilibration (if done).  If this is less than the machine
261: *     precision (in particular, if it is zero), the matrix is singular
262: *     to working precision.  Note that the error may still be small even
263: *     if this number is very small and the matrix appears ill-
264: *     conditioned.
265: *
266: *     RPVGRW  (output) REAL
267: *     Reciprocal pivot growth.  On exit, this contains the reciprocal
268: *     pivot growth factor norm(A)/norm(U). The "max absolute element"
269: *     norm is used.  If this is much less than 1, then the stability of
270: *     the LU factorization of the (equilibrated) matrix A could be poor.
271: *     This also means that the solution X, estimated condition numbers,
272: *     and error bounds could be unreliable. If factorization fails with
273: *     0<INFO<=N, then this contains the reciprocal pivot growth factor
274: *     for the leading INFO columns of A.  In SGESVX, this quantity is
275: *     returned in WORK(1).
276: *
277: *     BERR    (output) REAL array, dimension (NRHS)
278: *     Componentwise relative backward error.  This is the
279: *     componentwise relative backward error of each solution vector X(j)
280: *     (i.e., the smallest relative change in any element of A or B that
281: *     makes X(j) an exact solution).
282: *
283: *     N_ERR_BNDS (input) INTEGER
284: *     Number of error bounds to return for each right hand side
285: *     and each type (normwise or componentwise).  See ERR_BNDS_NORM and
286: *     ERR_BNDS_COMP below.
287: *
288: *     ERR_BNDS_NORM  (output) REAL array, dimension (NRHS, N_ERR_BNDS)
289: *     For each right-hand side, this array contains information about
290: *     various error bounds and condition numbers corresponding to the
291: *     normwise relative error, which is defined as follows:
292: *
293: *     Normwise relative error in the ith solution vector:
294: *             max_j (abs(XTRUE(j,i) - X(j,i)))
295: *            ------------------------------
296: *                  max_j abs(X(j,i))
297: *
298: *     The array is indexed by the type of error information as described
299: *     below. There currently are up to three pieces of information
300: *     returned.
301: *
302: *     The first index in ERR_BNDS_NORM(i,:) corresponds to the ith
303: *     right-hand side.
304: *
305: *     The second index in ERR_BNDS_NORM(:,err) contains the following
306: *     three fields:
307: *     err = 1 "Trust/don't trust" boolean. Trust the answer if the
308: *              reciprocal condition number is less than the threshold
309: *              sqrt(n) * slamch('Epsilon').
310: *
311: *     err = 2 "Guaranteed" error bound: The estimated forward error,
312: *              almost certainly within a factor of 10 of the true error
313: *              so long as the next entry is greater than the threshold
314: *              sqrt(n) * slamch('Epsilon'). This error bound should only
315: *              be trusted if the previous boolean is true.
316: *
317: *     err = 3  Reciprocal condition number: Estimated normwise
318: *              reciprocal condition number.  Compared with the threshold
319: *              sqrt(n) * slamch('Epsilon') to determine if the error
320: *              estimate is "guaranteed". These reciprocal condition
321: *              numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some
322: *              appropriately scaled matrix Z.
323: *              Let Z = S*A, where S scales each row by a power of the
324: *              radix so all absolute row sums of Z are approximately 1.
325: *
326: *     See Lapack Working Note 165 for further details and extra
327: *     cautions.
328: *
329: *     ERR_BNDS_COMP  (output) REAL array, dimension (NRHS, N_ERR_BNDS)
330: *     For each right-hand side, this array contains information about
331: *     various error bounds and condition numbers corresponding to the
332: *     componentwise relative error, which is defined as follows:
333: *
334: *     Componentwise relative error in the ith solution vector:
335: *                    abs(XTRUE(j,i) - X(j,i))
336: *             max_j ----------------------
337: *                         abs(X(j,i))
338: *
339: *     The array is indexed by the right-hand side i (on which the
340: *     componentwise relative error depends), and the type of error
341: *     information as described below. There currently are up to three
342: *     pieces of information returned for each right-hand side. If
343: *     componentwise accuracy is not requested (PARAMS(3) = 0.0), then
344: *     ERR_BNDS_COMP is not accessed.  If N_ERR_BNDS .LT. 3, then at most
345: *     the first (:,N_ERR_BNDS) entries are returned.
346: *
347: *     The first index in ERR_BNDS_COMP(i,:) corresponds to the ith
348: *     right-hand side.
349: *
350: *     The second index in ERR_BNDS_COMP(:,err) contains the following
351: *     three fields:
352: *     err = 1 "Trust/don't trust" boolean. Trust the answer if the
353: *              reciprocal condition number is less than the threshold
354: *              sqrt(n) * slamch('Epsilon').
355: *
356: *     err = 2 "Guaranteed" error bound: The estimated forward error,
357: *              almost certainly within a factor of 10 of the true error
358: *              so long as the next entry is greater than the threshold
359: *              sqrt(n) * slamch('Epsilon'). This error bound should only
360: *              be trusted if the previous boolean is true.
361: *
362: *     err = 3  Reciprocal condition number: Estimated componentwise
363: *              reciprocal condition number.  Compared with the threshold
364: *              sqrt(n) * slamch('Epsilon') to determine if the error
365: *              estimate is "guaranteed". These reciprocal condition
366: *              numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some
367: *              appropriately scaled matrix Z.
368: *              Let Z = S*(A*diag(x)), where x is the solution for the
369: *              current right-hand side and S scales each row of
370: *              A*diag(x) by a power of the radix so all absolute row
371: *              sums of Z are approximately 1.
372: *
373: *     See Lapack Working Note 165 for further details and extra
374: *     cautions.
375: *
376: *     NPARAMS (input) INTEGER
377: *     Specifies the number of parameters set in PARAMS.  If .LE. 0, the
378: *     PARAMS array is never referenced and default values are used.
379: *
380: *     PARAMS  (input / output) REAL array, dimension NPARAMS
381: *     Specifies algorithm parameters.  If an entry is .LT. 0.0, then
382: *     that entry will be filled with default value used for that
383: *     parameter.  Only positions up to NPARAMS are accessed; defaults
384: *     are used for higher-numbered parameters.
385: *
386: *       PARAMS(LA_LINRX_ITREF_I = 1) : Whether to perform iterative
387: *            refinement or not.
388: *         Default: 1.0
389: *            = 0.0 : No refinement is performed, and no error bounds are
390: *                    computed.
391: *            = 1.0 : Use the double-precision refinement algorithm,
392: *                    possibly with doubled-single computations if the
393: *                    compilation environment does not support DOUBLE
394: *                    PRECISION.
395: *              (other values are reserved for future use)
396: *
397: *       PARAMS(LA_LINRX_ITHRESH_I = 2) : Maximum number of residual
398: *            computations allowed for refinement.
399: *         Default: 10
400: *         Aggressive: Set to 100 to permit convergence using approximate
401: *                     factorizations or factorizations other than LU. If
402: *                     the factorization uses a technique other than
403: *                     Gaussian elimination, the guarantees in
404: *                     err_bnds_norm and err_bnds_comp may no longer be
405: *                     trustworthy.
406: *
407: *       PARAMS(LA_LINRX_CWISE_I = 3) : Flag determining if the code
408: *            will attempt to find a solution with small componentwise
409: *            relative error in the double-precision algorithm.  Positive
410: *            is true, 0.0 is false.
411: *         Default: 1.0 (attempt componentwise convergence)
412: *
413: *     WORK    (workspace) COMPLEX array, dimension (2*N)
414: *
415: *     RWORK   (workspace) REAL array, dimension (2*N)
416: *
417: *     INFO    (output) INTEGER
418: *       = 0:  Successful exit. The solution to every right-hand side is
419: *         guaranteed.
420: *       < 0:  If INFO = -i, the i-th argument had an illegal value
421: *       > 0 and <= N:  U(INFO,INFO) is exactly zero.  The factorization
422: *         has been completed, but the factor U is exactly singular, so
423: *         the solution and error bounds could not be computed. RCOND = 0
424: *         is returned.
425: *       = N+J: The solution corresponding to the Jth right-hand side is
426: *         not guaranteed. The solutions corresponding to other right-
427: *         hand sides K with K > J may not be guaranteed as well, but
428: *         only the first such right-hand side is reported. If a small
429: *         componentwise error is not requested (PARAMS(3) = 0.0) then
430: *         the Jth right-hand side is the first with a normwise error
431: *         bound that is not guaranteed (the smallest J such
432: *         that ERR_BNDS_NORM(J,1) = 0.0). By default (PARAMS(3) = 1.0)
433: *         the Jth right-hand side is the first with either a normwise or
434: *         componentwise error bound that is not guaranteed (the smallest
435: *         J such that either ERR_BNDS_NORM(J,1) = 0.0 or
436: *         ERR_BNDS_COMP(J,1) = 0.0). See the definition of
437: *         ERR_BNDS_NORM(:,1) and ERR_BNDS_COMP(:,1). To get information
438: *         about all of the right-hand sides check ERR_BNDS_NORM or
439: *         ERR_BNDS_COMP.
440: *
441: *     ==================================================================
442: *
443: *     .. Parameters ..
444:       REAL               ZERO, ONE
445:       PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0 )
446:       INTEGER            FINAL_NRM_ERR_I, FINAL_CMP_ERR_I, BERR_I
447:       INTEGER            RCOND_I, NRM_RCOND_I, NRM_ERR_I, CMP_RCOND_I
448:       INTEGER            CMP_ERR_I, PIV_GROWTH_I
449:       PARAMETER          ( FINAL_NRM_ERR_I = 1, FINAL_CMP_ERR_I = 2,
450:      $                   BERR_I = 3 )
451:       PARAMETER          ( RCOND_I = 4, NRM_RCOND_I = 5, NRM_ERR_I = 6 )
452:       PARAMETER          ( CMP_RCOND_I = 7, CMP_ERR_I = 8,
453:      $                   PIV_GROWTH_I = 9 )
454: *     ..
455: *     .. Local Scalars ..
456:       LOGICAL            COLEQU, EQUIL, NOFACT, NOTRAN, ROWEQU
457:       INTEGER            INFEQU, I, J, KL, KU
458:       REAL               AMAX, BIGNUM, COLCND, RCMAX, RCMIN,
459:      $                   ROWCND, SMLNUM
460: *     ..
461: *     .. External Functions ..
462:       EXTERNAL           LSAME, SLAMCH, CLA_GBRPVGRW
463:       LOGICAL            LSAME
464:       REAL               SLAMCH, CLA_GBRPVGRW
465: *     ..
466: *     .. External Subroutines ..
467:       EXTERNAL           CGBEQUB, CGBTRF, CGBTRS, CLACPY, CLAQGB,
468:      $                   XERBLA, CLASCL2, CGBRFSX
469: *     ..
470: *     .. Intrinsic Functions ..
471:       INTRINSIC          MAX, MIN
472: *     ..
473: *     .. Executable Statements ..
474: *
475:       INFO = 0
476:       NOFACT = LSAME( FACT, 'N' )
477:       EQUIL = LSAME( FACT, 'E' )
478:       NOTRAN = LSAME( TRANS, 'N' )
479:       SMLNUM = SLAMCH( 'Safe minimum' )
480:       BIGNUM = ONE / SMLNUM
481:       IF( NOFACT .OR. EQUIL ) THEN
482:          EQUED = 'N'
483:          ROWEQU = .FALSE.
484:          COLEQU = .FALSE.
485:       ELSE
486:          ROWEQU = LSAME( EQUED, 'R' ) .OR. LSAME( EQUED, 'B' )
487:          COLEQU = LSAME( EQUED, 'C' ) .OR. LSAME( EQUED, 'B' )
488:       END IF
489: *
490: *     Default is failure.  If an input parameter is wrong or
491: *     factorization fails, make everything look horrible.  Only the
492: *     pivot growth is set here, the rest is initialized in CGBRFSX.
493: *
494:       RPVGRW = ZERO
495: *
496: *     Test the input parameters.  PARAMS is not tested until SGERFSX.
497: *
498:       IF( .NOT.NOFACT .AND. .NOT.EQUIL .AND. .NOT.
499:      $     LSAME( FACT, 'F' ) ) THEN
500:          INFO = -1
501:       ELSE IF( .NOT.NOTRAN .AND. .NOT.LSAME( TRANS, 'T' ) .AND. .NOT.
502:      $        LSAME( TRANS, 'C' ) ) THEN
503:          INFO = -2
504:       ELSE IF( N.LT.0 ) THEN
505:          INFO = -3
506:       ELSE IF( KL.LT.0 ) THEN
507:          INFO = -4
508:       ELSE IF( KU.LT.0 ) THEN
509:          INFO = -5
510:       ELSE IF( NRHS.LT.0 ) THEN
511:          INFO = -6
512:       ELSE IF( LDAB.LT.KL+KU+1 ) THEN
513:          INFO = -8
514:       ELSE IF( LDAFB.LT.2*KL+KU+1 ) THEN
515:          INFO = -10
516:       ELSE IF( LSAME( FACT, 'F' ) .AND. .NOT.
517:      $        ( ROWEQU .OR. COLEQU .OR. LSAME( EQUED, 'N' ) ) ) THEN
518:          INFO = -12
519:       ELSE
520:          IF( ROWEQU ) THEN
521:             RCMIN = BIGNUM
522:             RCMAX = ZERO
523:             DO 10 J = 1, N
524:                RCMIN = MIN( RCMIN, R( J ) )
525:                RCMAX = MAX( RCMAX, R( J ) )
526:  10         CONTINUE
527:             IF( RCMIN.LE.ZERO ) THEN
528:                INFO = -13
529:             ELSE IF( N.GT.0 ) THEN
530:                ROWCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
531:             ELSE
532:                ROWCND = ONE
533:             END IF
534:          END IF
535:          IF( COLEQU .AND. INFO.EQ.0 ) THEN
536:             RCMIN = BIGNUM
537:             RCMAX = ZERO
538:             DO 20 J = 1, N
539:                RCMIN = MIN( RCMIN, C( J ) )
540:                RCMAX = MAX( RCMAX, C( J ) )
541:  20         CONTINUE
542:             IF( RCMIN.LE.ZERO ) THEN
543:                INFO = -14
544:             ELSE IF( N.GT.0 ) THEN
545:                COLCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
546:             ELSE
547:                COLCND = ONE
548:             END IF
549:          END IF
550:          IF( INFO.EQ.0 ) THEN
551:             IF( LDB.LT.MAX( 1, N ) ) THEN
552:                INFO = -15
553:             ELSE IF( LDX.LT.MAX( 1, N ) ) THEN
554:                INFO = -16
555:             END IF
556:          END IF
557:       END IF
558: *
559:       IF( INFO.NE.0 ) THEN
560:          CALL XERBLA( 'CGBSVXX', -INFO )
561:          RETURN
562:       END IF
563: *
564:       IF( EQUIL ) THEN
565: *
566: *     Compute row and column scalings to equilibrate the matrix A.
567: *
568:          CALL CGBEQUB( N, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND,
569:      $        AMAX, INFEQU )
570:          IF( INFEQU.EQ.0 ) THEN
571: *
572: *     Equilibrate the matrix.
573: *
574:             CALL CLAQGB( N, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND,
575:      $           AMAX, EQUED )
576:             ROWEQU = LSAME( EQUED, 'R' ) .OR. LSAME( EQUED, 'B' )
577:             COLEQU = LSAME( EQUED, 'C' ) .OR. LSAME( EQUED, 'B' )
578:          END IF
579: *
580: *     If the scaling factors are not applied, set them to 1.0.
581: *
582:          IF ( .NOT.ROWEQU ) THEN
583:             DO J = 1, N
584:                R( J ) = 1.0
585:             END DO
586:          END IF
587:          IF ( .NOT.COLEQU ) THEN
588:             DO J = 1, N
589:                C( J ) = 1.0
590:             END DO
591:          END IF
592:       END IF
593: *
594: *     Scale the right-hand side.
595: *
596:       IF( NOTRAN ) THEN
597:          IF( ROWEQU ) CALL CLASCL2( N, NRHS, R, B, LDB )
598:       ELSE
599:          IF( COLEQU ) CALL CLASCL2( N, NRHS, C, B, LDB )
600:       END IF
601: *
602:       IF( NOFACT .OR. EQUIL ) THEN
603: *
604: *        Compute the LU factorization of A.
605: *
606:          DO 40, J = 1, N
607:             DO 30, I = KL+1, 2*KL+KU+1
608:                AFB( I, J ) = AB( I-KL, J )
609:  30         CONTINUE
610:  40      CONTINUE
611:          CALL CGBTRF( N, N, KL, KU, AFB, LDAFB, IPIV, INFO )
612: *
613: *        Return if INFO is non-zero.
614: *
615:          IF( INFO.GT.0 ) THEN
616: *
617: *           Pivot in column INFO is exactly 0
618: *           Compute the reciprocal pivot growth factor of the
619: *           leading rank-deficient INFO columns of A.
620: *
621:             RPVGRW = CLA_GBRPVGRW( N, KL, KU, INFO, AB, LDAB, AFB,
622:      $           LDAFB )
623:             RETURN
624:          END IF
625:       END IF
626: *
627: *     Compute the reciprocal pivot growth factor RPVGRW.
628: *
629:       RPVGRW = CLA_GBRPVGRW( N, KL, KU, N, AB, LDAB, AFB, LDAFB )
630: *
631: *     Compute the solution matrix X.
632: *
633:       CALL CLACPY( 'Full', N, NRHS, B, LDB, X, LDX )
634:       CALL CGBTRS( TRANS, N, KL, KU, NRHS, AFB, LDAFB, IPIV, X, LDX,
635:      $     INFO )
636: *
637: *     Use iterative refinement to improve the computed solution and
638: *     compute error bounds and backward error estimates for it.
639: *
640:       CALL CGBRFSX( TRANS, EQUED, N, KL, KU, NRHS, AB, LDAB, AFB, LDAFB,
641:      $     IPIV, R, C, B, LDB, X, LDX, RCOND, BERR,
642:      $     N_ERR_BNDS, ERR_BNDS_NORM, ERR_BNDS_COMP, NPARAMS, PARAMS,
643:      $     WORK, RWORK, INFO )
644: 
645: *
646: *     Scale solutions.
647: *
648:       IF ( COLEQU .AND. NOTRAN ) THEN
649:          CALL CLASCL2( N, NRHS, C, X, LDX )
650:       ELSE IF ( ROWEQU .AND. .NOT.NOTRAN ) THEN
651:          CALL CLASCL2( N, NRHS, R, X, LDX )
652:       END IF
653: *
654:       RETURN
655: *
656: *     End of CGBSVXX
657: *
658:       END
659: