LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ ctbt05()

subroutine ctbt05 ( character  UPLO,
character  TRANS,
character  DIAG,
integer  N,
integer  KD,
integer  NRHS,
complex, dimension( ldab, * )  AB,
integer  LDAB,
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 
)

CTBT05

Purpose:
 CTBT05 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 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 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 187 of file ctbt05.f.

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