1*> \brief \b CPTT05
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 CPTT05( N, NRHS, D, E, B, LDB, X, LDX, XACT, LDXACT,
12*                          FERR, BERR, RESLTS )
13*
14*       .. Scalar Arguments ..
15*       INTEGER            LDB, LDX, LDXACT, N, NRHS
16*       ..
17*       .. Array Arguments ..
18*       REAL               BERR( * ), D( * ), FERR( * ), RESLTS( * )
19*       COMPLEX            B( LDB, * ), E( * ), X( LDX, * ),
20*      $                   XACT( LDXACT, * )
21*       ..
22*
23*
24*> \par Purpose:
25*  =============
26*>
27*> \verbatim
28*>
29*> CPTT05 tests the error bounds from iterative refinement for the
30*> computed solution to a system of equations A*X = B, where A is a
31*> Hermitian tridiagonal matrix of order n.
32*>
33*> RESLTS(1) = test of the error bound
34*>           = norm(X - XACT) / ( norm(X) * FERR )
35*>
36*> A large value is returned if this ratio is not less than one.
37*>
38*> RESLTS(2) = residual from the iterative refinement routine
39*>           = the maximum of BERR / ( NZ*EPS + (*) ), where
40*>             (*) = NZ*UNFL / (min_i (abs(A)*abs(X) +abs(b))_i )
41*>             and NZ = max. number of nonzeros in any row of A, plus 1
42*> \endverbatim
43*
44*  Arguments:
45*  ==========
46*
47*> \param[in] N
48*> \verbatim
49*>          N is INTEGER
50*>          The number of rows of the matrices X, B, and XACT, and the
51*>          order of the matrix A.  N >= 0.
52*> \endverbatim
53*>
54*> \param[in] NRHS
55*> \verbatim
56*>          NRHS is INTEGER
57*>          The number of columns of the matrices X, B, and XACT.
58*>          NRHS >= 0.
59*> \endverbatim
60*>
61*> \param[in] D
62*> \verbatim
63*>          D is REAL array, dimension (N)
64*>          The n diagonal elements of the tridiagonal matrix A.
65*> \endverbatim
66*>
67*> \param[in] E
68*> \verbatim
69*>          E is COMPLEX array, dimension (N-1)
70*>          The (n-1) subdiagonal elements of the tridiagonal matrix A.
71*> \endverbatim
72*>
73*> \param[in] B
74*> \verbatim
75*>          B is COMPLEX array, dimension (LDB,NRHS)
76*>          The right hand side vectors for the system of linear
77*>          equations.
78*> \endverbatim
79*>
80*> \param[in] LDB
81*> \verbatim
82*>          LDB is INTEGER
83*>          The leading dimension of the array B.  LDB >= max(1,N).
84*> \endverbatim
85*>
86*> \param[in] X
87*> \verbatim
88*>          X is COMPLEX array, dimension (LDX,NRHS)
89*>          The computed solution vectors.  Each vector is stored as a
90*>          column of the matrix X.
91*> \endverbatim
92*>
93*> \param[in] LDX
94*> \verbatim
95*>          LDX is INTEGER
96*>          The leading dimension of the array X.  LDX >= max(1,N).
97*> \endverbatim
98*>
99*> \param[in] XACT
100*> \verbatim
101*>          XACT is COMPLEX array, dimension (LDX,NRHS)
102*>          The exact solution vectors.  Each vector is stored as a
103*>          column of the matrix XACT.
104*> \endverbatim
105*>
106*> \param[in] LDXACT
107*> \verbatim
108*>          LDXACT is INTEGER
109*>          The leading dimension of the array XACT.  LDXACT >= max(1,N).
110*> \endverbatim
111*>
112*> \param[in] FERR
113*> \verbatim
114*>          FERR is REAL array, dimension (NRHS)
115*>          The estimated forward error bounds for each solution vector
116*>          X.  If XTRUE is the true solution, FERR bounds the magnitude
117*>          of the largest entry in (X - XTRUE) divided by the magnitude
118*>          of the largest entry in X.
119*> \endverbatim
120*>
121*> \param[in] BERR
122*> \verbatim
123*>          BERR is REAL array, dimension (NRHS)
124*>          The componentwise relative backward error of each solution
125*>          vector (i.e., the smallest relative change in any entry of A
126*>          or B that makes X an exact solution).
127*> \endverbatim
128*>
129*> \param[out] RESLTS
130*> \verbatim
131*>          RESLTS is REAL array, dimension (2)
132*>          The maximum over the NRHS solution vectors of the ratios:
133*>          RESLTS(1) = norm(X - XACT) / ( norm(X) * FERR )
134*>          RESLTS(2) = BERR / ( NZ*EPS + (*) )
135*> \endverbatim
136*
137*  Authors:
138*  ========
139*
140*> \author Univ. of Tennessee
141*> \author Univ. of California Berkeley
142*> \author Univ. of Colorado Denver
143*> \author NAG Ltd.
144*
145*> \ingroup complex_lin
146*
147*  =====================================================================
148      SUBROUTINE CPTT05( N, NRHS, D, E, B, LDB, X, LDX, XACT, LDXACT,
149     $                   FERR, BERR, RESLTS )
150*
151*  -- LAPACK test routine --
152*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
153*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
154*
155*     .. Scalar Arguments ..
156      INTEGER            LDB, LDX, LDXACT, N, NRHS
157*     ..
158*     .. Array Arguments ..
159      REAL               BERR( * ), D( * ), FERR( * ), RESLTS( * )
160      COMPLEX            B( LDB, * ), E( * ), X( LDX, * ),
161     $                   XACT( LDXACT, * )
162*     ..
163*
164*  =====================================================================
165*
166*     .. Parameters ..
167      REAL               ZERO, ONE
168      PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0 )
169*     ..
170*     .. Local Scalars ..
171      INTEGER            I, IMAX, J, K, NZ
172      REAL               AXBI, DIFF, EPS, ERRBND, OVFL, TMP, UNFL, XNORM
173      COMPLEX            ZDUM
174*     ..
175*     .. External Functions ..
176      INTEGER            ICAMAX
177      REAL               SLAMCH
178      EXTERNAL           ICAMAX, SLAMCH
179*     ..
180*     .. Intrinsic Functions ..
181      INTRINSIC          ABS, AIMAG, MAX, MIN, REAL
182*     ..
183*     .. Statement Functions ..
184      REAL               CABS1
185*     ..
186*     .. Statement Function definitions ..
187      CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
188*     ..
189*     .. Executable Statements ..
190*
191*     Quick exit if N = 0 or NRHS = 0.
192*
193      IF( N.LE.0 .OR. NRHS.LE.0 ) THEN
194         RESLTS( 1 ) = ZERO
195         RESLTS( 2 ) = ZERO
196         RETURN
197      END IF
198*
199      EPS = SLAMCH( 'Epsilon' )
200      UNFL = SLAMCH( 'Safe minimum' )
201      OVFL = ONE / UNFL
202      NZ = 4
203*
204*     Test 1:  Compute the maximum of
205*        norm(X - XACT) / ( norm(X) * FERR )
206*     over all the vectors X and XACT using the infinity-norm.
207*
208      ERRBND = ZERO
209      DO 30 J = 1, NRHS
210         IMAX = ICAMAX( N, X( 1, J ), 1 )
211         XNORM = MAX( CABS1( X( IMAX, J ) ), UNFL )
212         DIFF = ZERO
213         DO 10 I = 1, N
214            DIFF = MAX( DIFF, CABS1( X( I, J )-XACT( I, J ) ) )
215   10    CONTINUE
216*
217         IF( XNORM.GT.ONE ) THEN
218            GO TO 20
219         ELSE IF( DIFF.LE.OVFL*XNORM ) THEN
220            GO TO 20
221         ELSE
222            ERRBND = ONE / EPS
223            GO TO 30
224         END IF
225*
226   20    CONTINUE
227         IF( DIFF / XNORM.LE.FERR( J ) ) THEN
228            ERRBND = MAX( ERRBND, ( DIFF / XNORM ) / FERR( J ) )
229         ELSE
230            ERRBND = ONE / EPS
231         END IF
232   30 CONTINUE
233      RESLTS( 1 ) = ERRBND
234*
235*     Test 2:  Compute the maximum of BERR / ( NZ*EPS + (*) ), where
236*     (*) = NZ*UNFL / (min_i (abs(A)*abs(X) +abs(b))_i )
237*
238      DO 50 K = 1, NRHS
239         IF( N.EQ.1 ) THEN
240            AXBI = CABS1( B( 1, K ) ) + CABS1( D( 1 )*X( 1, K ) )
241         ELSE
242            AXBI = CABS1( B( 1, K ) ) + CABS1( D( 1 )*X( 1, K ) ) +
243     $             CABS1( E( 1 ) )*CABS1( X( 2, K ) )
244            DO 40 I = 2, N - 1
245               TMP = CABS1( B( I, K ) ) + CABS1( E( I-1 ) )*
246     $               CABS1( X( I-1, K ) ) + CABS1( D( I )*X( I, K ) ) +
247     $               CABS1( E( I ) )*CABS1( X( I+1, K ) )
248               AXBI = MIN( AXBI, TMP )
249   40       CONTINUE
250            TMP = CABS1( B( N, K ) ) + CABS1( E( N-1 ) )*
251     $            CABS1( X( N-1, K ) ) + CABS1( D( N )*X( N, K ) )
252            AXBI = MIN( AXBI, TMP )
253         END IF
254         TMP = BERR( K ) / ( NZ*EPS+NZ*UNFL / MAX( AXBI, NZ*UNFL ) )
255         IF( K.EQ.1 ) THEN
256            RESLTS( 2 ) = TMP
257         ELSE
258            RESLTS( 2 ) = MAX( RESLTS( 2 ), TMP )
259         END IF
260   50 CONTINUE
261*
262      RETURN
263*
264*     End of CPTT05
265*
266      END
267