LAPACK 3.11.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ cptt05()

subroutine cptt05 ( integer  N,
integer  NRHS,
real, dimension( * )  D,
complex, dimension( * )  E,
complex, dimension( ldb, * )  B,
integer  LDB,
complex, dimension( ldx, * )  X,
integer  LDX,
complex, dimension( ldxact, * )  XACT,
integer  LDXACT,
real, dimension( * )  FERR,
real, dimension( * )  BERR,
real, dimension( * )  RESLTS 
)

CPTT05

Purpose:
 CPTT05 tests the error bounds from iterative refinement for the
 computed solution to a system of equations A*X = B, where A is a
 Hermitian tridiagonal matrix of order n.

 RESLTS(1) = test of the error bound
           = norm(X - XACT) / ( norm(X) * FERR )

 A large value is returned if this ratio is not less than one.

 RESLTS(2) = residual from the iterative refinement routine
           = the maximum of BERR / ( NZ*EPS + (*) ), where
             (*) = NZ*UNFL / (min_i (abs(A)*abs(X) +abs(b))_i )
             and NZ = max. number of nonzeros in any row of A, plus 1
Parameters
[in]N
          N is INTEGER
          The number of rows of the matrices X, B, and XACT, and the
          order of the matrix A.  N >= 0.
[in]NRHS
          NRHS is INTEGER
          The number of columns of the matrices X, B, and XACT.
          NRHS >= 0.
[in]D
          D is REAL array, dimension (N)
          The n diagonal elements of the tridiagonal matrix A.
[in]E
          E is COMPLEX array, dimension (N-1)
          The (n-1) subdiagonal elements of the tridiagonal matrix A.
[in]B
          B is COMPLEX array, dimension (LDB,NRHS)
          The right hand side vectors for the system of linear
          equations.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
[in]X
          X is COMPLEX array, dimension (LDX,NRHS)
          The computed solution vectors.  Each vector is stored as a
          column of the matrix X.
[in]LDX
          LDX is INTEGER
          The leading dimension of the array X.  LDX >= max(1,N).
[in]XACT
          XACT is COMPLEX array, dimension (LDX,NRHS)
          The exact solution vectors.  Each vector is stored as a
          column of the matrix XACT.
[in]LDXACT
          LDXACT is INTEGER
          The leading dimension of the array XACT.  LDXACT >= max(1,N).
[in]FERR
          FERR is REAL array, dimension (NRHS)
          The estimated forward error bounds for each solution vector
          X.  If XTRUE is the true solution, FERR bounds the magnitude
          of the largest entry in (X - XTRUE) divided by the magnitude
          of the largest entry in X.
[in]BERR
          BERR is REAL array, dimension (NRHS)
          The componentwise relative backward error of each solution
          vector (i.e., the smallest relative change in any entry of A
          or B that makes X an exact solution).
[out]RESLTS
          RESLTS is REAL array, dimension (2)
          The maximum over the NRHS solution vectors of the ratios:
          RESLTS(1) = norm(X - XACT) / ( norm(X) * FERR )
          RESLTS(2) = BERR / ( NZ*EPS + (*) )
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 148 of file cptt05.f.

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*
integer function icamax(N, CX, INCX)
ICAMAX
Definition: icamax.f:71
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:68
Here is the caller graph for this function: