LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
cptt05.f
Go to the documentation of this file.
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
subroutine cptt05(n, nrhs, d, e, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
CPTT05
Definition cptt05.f:150