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*> \date November 2011
146*
147*> \ingroup complex_lin
148*
149*  =====================================================================
150      SUBROUTINE CPTT05( N, NRHS, D, E, B, LDB, X, LDX, XACT, LDXACT,
151     $                   FERR, BERR, RESLTS )
152*
153*  -- LAPACK test routine (version 3.4.0) --
154*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
155*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
156*     November 2011
157*
158*     .. Scalar Arguments ..
159      INTEGER            LDB, LDX, LDXACT, N, NRHS
160*     ..
161*     .. Array Arguments ..
162      REAL               BERR( * ), D( * ), FERR( * ), RESLTS( * )
163      COMPLEX            B( LDB, * ), E( * ), X( LDX, * ),
164     $                   XACT( LDXACT, * )
165*     ..
166*
167*  =====================================================================
168*
169*     .. Parameters ..
170      REAL               ZERO, ONE
171      PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0 )
172*     ..
173*     .. Local Scalars ..
174      INTEGER            I, IMAX, J, K, NZ
175      REAL               AXBI, DIFF, EPS, ERRBND, OVFL, TMP, UNFL, XNORM
176      COMPLEX            ZDUM
177*     ..
178*     .. External Functions ..
179      INTEGER            ICAMAX
180      REAL               SLAMCH
181      EXTERNAL           ICAMAX, SLAMCH
182*     ..
183*     .. Intrinsic Functions ..
184      INTRINSIC          ABS, AIMAG, MAX, MIN, REAL
185*     ..
186*     .. Statement Functions ..
187      REAL               CABS1
188*     ..
189*     .. Statement Function definitions ..
190      CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
191*     ..
192*     .. Executable Statements ..
193*
194*     Quick exit if N = 0 or NRHS = 0.
195*
196      IF( N.LE.0 .OR. NRHS.LE.0 ) THEN
197         RESLTS( 1 ) = ZERO
198         RESLTS( 2 ) = ZERO
199         RETURN
200      END IF
201*
202      EPS = SLAMCH( 'Epsilon' )
203      UNFL = SLAMCH( 'Safe minimum' )
204      OVFL = ONE / UNFL
205      NZ = 4
206*
207*     Test 1:  Compute the maximum of
208*        norm(X - XACT) / ( norm(X) * FERR )
209*     over all the vectors X and XACT using the infinity-norm.
210*
211      ERRBND = ZERO
212      DO 30 J = 1, NRHS
213         IMAX = ICAMAX( N, X( 1, J ), 1 )
214         XNORM = MAX( CABS1( X( IMAX, J ) ), UNFL )
215         DIFF = ZERO
216         DO 10 I = 1, N
217            DIFF = MAX( DIFF, CABS1( X( I, J )-XACT( I, J ) ) )
218   10    CONTINUE
219*
220         IF( XNORM.GT.ONE ) THEN
221            GO TO 20
222         ELSE IF( DIFF.LE.OVFL*XNORM ) THEN
223            GO TO 20
224         ELSE
225            ERRBND = ONE / EPS
226            GO TO 30
227         END IF
228*
229   20    CONTINUE
230         IF( DIFF / XNORM.LE.FERR( J ) ) THEN
231            ERRBND = MAX( ERRBND, ( DIFF / XNORM ) / FERR( J ) )
232         ELSE
233            ERRBND = ONE / EPS
234         END IF
235   30 CONTINUE
236      RESLTS( 1 ) = ERRBND
237*
238*     Test 2:  Compute the maximum of BERR / ( NZ*EPS + (*) ), where
239*     (*) = NZ*UNFL / (min_i (abs(A)*abs(X) +abs(b))_i )
240*
241      DO 50 K = 1, NRHS
242         IF( N.EQ.1 ) THEN
243            AXBI = CABS1( B( 1, K ) ) + CABS1( D( 1 )*X( 1, K ) )
244         ELSE
245            AXBI = CABS1( B( 1, K ) ) + CABS1( D( 1 )*X( 1, K ) ) +
246     $             CABS1( E( 1 ) )*CABS1( X( 2, K ) )
247            DO 40 I = 2, N - 1
248               TMP = CABS1( B( I, K ) ) + CABS1( E( I-1 ) )*
249     $               CABS1( X( I-1, K ) ) + CABS1( D( I )*X( I, K ) ) +
250     $               CABS1( E( I ) )*CABS1( X( I+1, K ) )
251               AXBI = MIN( AXBI, TMP )
252   40       CONTINUE
253            TMP = CABS1( B( N, K ) ) + CABS1( E( N-1 ) )*
254     $            CABS1( X( N-1, K ) ) + CABS1( D( N )*X( N, K ) )
255            AXBI = MIN( AXBI, TMP )
256         END IF
257         TMP = BERR( K ) / ( NZ*EPS+NZ*UNFL / MAX( AXBI, NZ*UNFL ) )
258         IF( K.EQ.1 ) THEN
259            RESLTS( 2 ) = TMP
260         ELSE
261            RESLTS( 2 ) = MAX( RESLTS( 2 ), TMP )
262         END IF
263   50 CONTINUE
264*
265      RETURN
266*
267*     End of CPTT05
268*
269      END
270