1*> \brief \b DTBT02
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*  Definition:
9*  ===========
10*
11*       SUBROUTINE DTBT02( UPLO, TRANS, DIAG, N, KD, NRHS, AB, LDAB, X,
12*                          LDX, B, LDB, WORK, RESID )
13*
14*       .. Scalar Arguments ..
15*       CHARACTER          DIAG, TRANS, UPLO
16*       INTEGER            KD, LDAB, LDB, LDX, N, NRHS
17*       DOUBLE PRECISION   RESID
18*       ..
19*       .. Array Arguments ..
20*       DOUBLE PRECISION   AB( LDAB, * ), B( LDB, * ), WORK( * ),
21*      $                   X( LDX, * )
22*       ..
23*
24*
25*> \par Purpose:
26*  =============
27*>
28*> \verbatim
29*>
30*> DTBT02 computes the residual for the computed solution to a
31*> triangular system of linear equations  A*x = b  or  A' *x = b when
32*> A is a triangular band matrix.  Here A' is the transpose of A and
33*> x and b are N by NRHS matrices.  The test ratio is the maximum over
34*> the number of right hand sides of
35*>    norm(b - op(A)*x) / ( norm(op(A)) * norm(x) * EPS ),
36*> where op(A) denotes A or A' and EPS is the machine epsilon.
37*> \endverbatim
38*
39*  Arguments:
40*  ==========
41*
42*> \param[in] UPLO
43*> \verbatim
44*>          UPLO is CHARACTER*1
45*>          Specifies whether the matrix A is upper or lower triangular.
46*>          = 'U':  Upper triangular
47*>          = 'L':  Lower triangular
48*> \endverbatim
49*>
50*> \param[in] TRANS
51*> \verbatim
52*>          TRANS is CHARACTER*1
53*>          Specifies the operation applied to A.
54*>          = 'N':  A *x = b  (No transpose)
55*>          = 'T':  A'*x = b  (Transpose)
56*>          = 'C':  A'*x = b  (Conjugate transpose = Transpose)
57*> \endverbatim
58*>
59*> \param[in] DIAG
60*> \verbatim
61*>          DIAG is CHARACTER*1
62*>          Specifies whether or not the matrix A is unit triangular.
63*>          = 'N':  Non-unit triangular
64*>          = 'U':  Unit triangular
65*> \endverbatim
66*>
67*> \param[in] N
68*> \verbatim
69*>          N is INTEGER
70*>          The order of the matrix A.  N >= 0.
71*> \endverbatim
72*>
73*> \param[in] KD
74*> \verbatim
75*>          KD is INTEGER
76*>          The number of superdiagonals or subdiagonals of the
77*>          triangular band matrix A.  KD >= 0.
78*> \endverbatim
79*>
80*> \param[in] NRHS
81*> \verbatim
82*>          NRHS is INTEGER
83*>          The number of right hand sides, i.e., the number of columns
84*>          of the matrices X and B.  NRHS >= 0.
85*> \endverbatim
86*>
87*> \param[in] AB
88*> \verbatim
89*>          AB is DOUBLE PRECISION array, dimension (LDAB,N)
90*>          The upper or lower triangular band matrix A, stored in the
91*>          first kd+1 rows of the array. The j-th column of A is stored
92*>          in the j-th column of the array AB as follows:
93*>          if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
94*>          if UPLO = 'L', AB(1+i-j,j)    = A(i,j) for j<=i<=min(n,j+kd).
95*> \endverbatim
96*>
97*> \param[in] LDAB
98*> \verbatim
99*>          LDAB is INTEGER
100*>          The leading dimension of the array AB.  LDAB >= KD+1.
101*> \endverbatim
102*>
103*> \param[in] X
104*> \verbatim
105*>          X is DOUBLE PRECISION array, dimension (LDX,NRHS)
106*>          The computed solution vectors for the system of linear
107*>          equations.
108*> \endverbatim
109*>
110*> \param[in] LDX
111*> \verbatim
112*>          LDX is INTEGER
113*>          The leading dimension of the array X.  LDX >= max(1,N).
114*> \endverbatim
115*>
116*> \param[in] B
117*> \verbatim
118*>          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
119*>          The right hand side vectors for the system of linear
120*>          equations.
121*> \endverbatim
122*>
123*> \param[in] LDB
124*> \verbatim
125*>          LDB is INTEGER
126*>          The leading dimension of the array B.  LDB >= max(1,N).
127*> \endverbatim
128*>
129*> \param[out] WORK
130*> \verbatim
131*>          WORK is DOUBLE PRECISION array, dimension (N)
132*> \endverbatim
133*>
134*> \param[out] RESID
135*> \verbatim
136*>          RESID is DOUBLE PRECISION
137*>          The maximum over the number of right hand sides of
138*>          norm(op(A)*x - b) / ( norm(op(A)) * norm(x) * EPS ).
139*> \endverbatim
140*
141*  Authors:
142*  ========
143*
144*> \author Univ. of Tennessee
145*> \author Univ. of California Berkeley
146*> \author Univ. of Colorado Denver
147*> \author NAG Ltd.
148*
149*> \date November 2011
150*
151*> \ingroup double_lin
152*
153*  =====================================================================
154      SUBROUTINE DTBT02( UPLO, TRANS, DIAG, N, KD, NRHS, AB, LDAB, X,
155     $                   LDX, B, LDB, WORK, RESID )
156*
157*  -- LAPACK test routine (version 3.4.0) --
158*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
159*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
160*     November 2011
161*
162*     .. Scalar Arguments ..
163      CHARACTER          DIAG, TRANS, UPLO
164      INTEGER            KD, LDAB, LDB, LDX, N, NRHS
165      DOUBLE PRECISION   RESID
166*     ..
167*     .. Array Arguments ..
168      DOUBLE PRECISION   AB( LDAB, * ), B( LDB, * ), WORK( * ),
169     $                   X( LDX, * )
170*     ..
171*
172*  =====================================================================
173*
174*     .. Parameters ..
175      DOUBLE PRECISION   ZERO, ONE
176      PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0 )
177*     ..
178*     .. Local Scalars ..
179      INTEGER            J
180      DOUBLE PRECISION   ANORM, BNORM, EPS, XNORM
181*     ..
182*     .. External Functions ..
183      LOGICAL            LSAME
184      DOUBLE PRECISION   DASUM, DLAMCH, DLANTB
185      EXTERNAL           LSAME, DASUM, DLAMCH, DLANTB
186*     ..
187*     .. External Subroutines ..
188      EXTERNAL           DAXPY, DCOPY, DTBMV
189*     ..
190*     .. Intrinsic Functions ..
191      INTRINSIC          MAX
192*     ..
193*     .. Executable Statements ..
194*
195*     Quick exit if N = 0 or NRHS = 0
196*
197      IF( N.LE.0 .OR. NRHS.LE.0 ) THEN
198         RESID = ZERO
199         RETURN
200      END IF
201*
202*     Compute the 1-norm of A or A'.
203*
204      IF( LSAME( TRANS, 'N' ) ) THEN
205         ANORM = DLANTB( '1', UPLO, DIAG, N, KD, AB, LDAB, WORK )
206      ELSE
207         ANORM = DLANTB( 'I', UPLO, DIAG, N, KD, AB, LDAB, WORK )
208      END IF
209*
210*     Exit with RESID = 1/EPS if ANORM = 0.
211*
212      EPS = DLAMCH( 'Epsilon' )
213      IF( ANORM.LE.ZERO ) THEN
214         RESID = ONE / EPS
215         RETURN
216      END IF
217*
218*     Compute the maximum over the number of right hand sides of
219*        norm(op(A)*x - b) / ( norm(op(A)) * norm(x) * EPS ).
220*
221      RESID = ZERO
222      DO 10 J = 1, NRHS
223         CALL DCOPY( N, X( 1, J ), 1, WORK, 1 )
224         CALL DTBMV( UPLO, TRANS, DIAG, N, KD, AB, LDAB, WORK, 1 )
225         CALL DAXPY( N, -ONE, B( 1, J ), 1, WORK, 1 )
226         BNORM = DASUM( N, WORK, 1 )
227         XNORM = DASUM( N, X( 1, J ), 1 )
228         IF( XNORM.LE.ZERO ) THEN
229            RESID = ONE / EPS
230         ELSE
231            RESID = MAX( RESID, ( ( BNORM / ANORM ) / XNORM ) / EPS )
232         END IF
233   10 CONTINUE
234*
235      RETURN
236*
237*     End of DTBT02
238*
239      END
240