LAPACK  3.6.1 LAPACK: Linear Algebra PACKage
 subroutine ztbt05 ( character UPLO, character TRANS, character DIAG, integer N, integer KD, integer NRHS, complex*16, dimension( ldab, * ) AB, integer LDAB, complex*16, dimension( ldb, * ) B, integer LDB, complex*16, dimension( ldx, * ) X, integer LDX, complex*16, dimension( ldxact, * ) XACT, integer LDXACT, double precision, dimension( * ) FERR, double precision, dimension( * ) BERR, double precision, dimension( * ) RESLTS )

ZTBT05

Purpose:
ZTBT05 tests the error bounds from iterative refinement for the
computed solution to a system of equations A*X = B, where A is a
triangular band matrix.

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] UPLO UPLO is CHARACTER*1 Specifies whether the matrix A is upper or lower triangular. = 'U': Upper triangular = 'L': Lower triangular [in] TRANS TRANS is CHARACTER*1 Specifies the form of the system of equations. = 'N': A * X = B (No transpose) = 'T': A'* X = B (Transpose) = 'C': A'* X = B (Conjugate transpose = Transpose) [in] DIAG DIAG is CHARACTER*1 Specifies whether or not the matrix A is unit triangular. = 'N': Non-unit triangular = 'U': Unit triangular [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] KD KD is INTEGER The number of super-diagonals of the matrix A if UPLO = 'U', or the number of sub-diagonals if UPLO = 'L'. KD >= 0. [in] NRHS NRHS is INTEGER The number of columns of the matrices X, B, and XACT. NRHS >= 0. [in] AB AB is COMPLEX*16 array, dimension (LDAB,N) The upper or lower triangular band matrix A, stored in the first kd+1 rows of the array. The j-th column of A is stored in the j-th column of the array AB as follows: if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j; if UPLO = 'L', AB(1+i-j,j) = A(i,j) for j<=i<=min(n,j+kd). If DIAG = 'U', the diagonal elements of A are not referenced and are assumed to be 1. [in] LDAB LDAB is INTEGER The leading dimension of the array AB. LDAB >= KD+1. [in] B B is COMPLEX*16 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*16 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*16 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 + (*) )
Date
November 2011

Definition at line 191 of file ztbt05.f.

191 *
192 * -- LAPACK test routine (version 3.4.0) --
193 * -- LAPACK is a software package provided by Univ. of Tennessee, --
194 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
195 * November 2011
196 *
197 * .. Scalar Arguments ..
198  CHARACTER diag, trans, uplo
199  INTEGER kd, ldab, ldb, ldx, ldxact, n, nrhs
200 * ..
201 * .. Array Arguments ..
202  DOUBLE PRECISION berr( * ), ferr( * ), reslts( * )
203  COMPLEX*16 ab( ldab, * ), b( ldb, * ), x( ldx, * ),
204  \$ xact( ldxact, * )
205 * ..
206 *
207 * =====================================================================
208 *
209 * .. Parameters ..
210  DOUBLE PRECISION zero, one
211  parameter ( zero = 0.0d+0, one = 1.0d+0 )
212 * ..
213 * .. Local Scalars ..
214  LOGICAL notran, unit, upper
215  INTEGER i, ifu, imax, j, k, nz
216  DOUBLE PRECISION axbi, diff, eps, errbnd, ovfl, tmp, unfl, xnorm
217  COMPLEX*16 zdum
218 * ..
219 * .. External Functions ..
220  LOGICAL lsame
221  INTEGER izamax
222  DOUBLE PRECISION dlamch
223  EXTERNAL lsame, izamax, dlamch
224 * ..
225 * .. Intrinsic Functions ..
226  INTRINSIC abs, dble, dimag, max, min
227 * ..
228 * .. Statement Functions ..
229  DOUBLE PRECISION cabs1
230 * ..
231 * .. Statement Function definitions ..
232  cabs1( zdum ) = abs( dble( zdum ) ) + abs( dimag( zdum ) )
233 * ..
234 * .. Executable Statements ..
235 *
236 * Quick exit if N = 0 or NRHS = 0.
237 *
238  IF( n.LE.0 .OR. nrhs.LE.0 ) THEN
239  reslts( 1 ) = zero
240  reslts( 2 ) = zero
241  RETURN
242  END IF
243 *
244  eps = dlamch( 'Epsilon' )
245  unfl = dlamch( 'Safe minimum' )
246  ovfl = one / unfl
247  upper = lsame( uplo, 'U' )
248  notran = lsame( trans, 'N' )
249  unit = lsame( diag, 'U' )
250  nz = min( kd, n-1 ) + 1
251 *
252 * Test 1: Compute the maximum of
253 * norm(X - XACT) / ( norm(X) * FERR )
254 * over all the vectors X and XACT using the infinity-norm.
255 *
256  errbnd = zero
257  DO 30 j = 1, nrhs
258  imax = izamax( n, x( 1, j ), 1 )
259  xnorm = max( cabs1( x( imax, j ) ), unfl )
260  diff = zero
261  DO 10 i = 1, n
262  diff = max( diff, cabs1( x( i, j )-xact( i, j ) ) )
263  10 CONTINUE
264 *
265  IF( xnorm.GT.one ) THEN
266  GO TO 20
267  ELSE IF( diff.LE.ovfl*xnorm ) THEN
268  GO TO 20
269  ELSE
270  errbnd = one / eps
271  GO TO 30
272  END IF
273 *
274  20 CONTINUE
275  IF( diff / xnorm.LE.ferr( j ) ) THEN
276  errbnd = max( errbnd, ( diff / xnorm ) / ferr( j ) )
277  ELSE
278  errbnd = one / eps
279  END IF
280  30 CONTINUE
281  reslts( 1 ) = errbnd
282 *
283 * Test 2: Compute the maximum of BERR / ( NZ*EPS + (*) ), where
284 * (*) = NZ*UNFL / (min_i (abs(A)*abs(X) +abs(b))_i )
285 *
286  ifu = 0
287  IF( unit )
288  \$ ifu = 1
289  DO 90 k = 1, nrhs
290  DO 80 i = 1, n
291  tmp = cabs1( b( i, k ) )
292  IF( upper ) THEN
293  IF( .NOT.notran ) THEN
294  DO 40 j = max( i-kd, 1 ), i - ifu
295  tmp = tmp + cabs1( ab( kd+1-i+j, i ) )*
296  \$ cabs1( x( j, k ) )
297  40 CONTINUE
298  IF( unit )
299  \$ tmp = tmp + cabs1( x( i, k ) )
300  ELSE
301  IF( unit )
302  \$ tmp = tmp + cabs1( x( i, k ) )
303  DO 50 j = i + ifu, min( i+kd, n )
304  tmp = tmp + cabs1( ab( kd+1+i-j, j ) )*
305  \$ cabs1( x( j, k ) )
306  50 CONTINUE
307  END IF
308  ELSE
309  IF( notran ) THEN
310  DO 60 j = max( i-kd, 1 ), i - ifu
311  tmp = tmp + cabs1( ab( 1+i-j, j ) )*
312  \$ cabs1( x( j, k ) )
313  60 CONTINUE
314  IF( unit )
315  \$ tmp = tmp + cabs1( x( i, k ) )
316  ELSE
317  IF( unit )
318  \$ tmp = tmp + cabs1( x( i, k ) )
319  DO 70 j = i + ifu, min( i+kd, n )
320  tmp = tmp + cabs1( ab( 1+j-i, i ) )*
321  \$ cabs1( x( j, k ) )
322  70 CONTINUE
323  END IF
324  END IF
325  IF( i.EQ.1 ) THEN
326  axbi = tmp
327  ELSE
328  axbi = min( axbi, tmp )
329  END IF
330  80 CONTINUE
331  tmp = berr( k ) / ( nz*eps+nz*unfl / max( axbi, nz*unfl ) )
332  IF( k.EQ.1 ) THEN
333  reslts( 2 ) = tmp
334  ELSE
335  reslts( 2 ) = max( reslts( 2 ), tmp )
336  END IF
337  90 CONTINUE
338 *
339  RETURN
340 *
341 * End of ZTBT05
342 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
integer function izamax(N, ZX, INCX)
IZAMAX
Definition: izamax.f:53
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the caller graph for this function: