 LAPACK  3.10.1 LAPACK: Linear Algebra PACKage

## ◆ ztbrfs()

 subroutine ztbrfs ( 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, double precision, dimension( * ) FERR, double precision, dimension( * ) BERR, complex*16, dimension( * ) WORK, double precision, dimension( * ) RWORK, integer INFO )

ZTBRFS

Purpose:
``` ZTBRFS provides error bounds and backward error estimates for the
solution to a system of linear equations with a triangular band
coefficient matrix.

The solution matrix X must be computed by ZTBTRS or some other
means before entering this routine.  ZTBRFS does not do iterative
refinement because doing so cannot improve the backward error.```
Parameters
 [in] UPLO ``` UPLO is CHARACTER*1 = 'U': A is upper triangular; = 'L': A is lower triangular.``` [in] TRANS ``` TRANS is CHARACTER*1 Specifies the form of the system of equations: = 'N': A * X = B (No transpose) = 'T': A**T * X = B (Transpose) = 'C': A**H * X = B (Conjugate transpose)``` [in] DIAG ``` DIAG is CHARACTER*1 = 'N': A is non-unit triangular; = 'U': A is unit triangular.``` [in] N ``` N is INTEGER The order of the matrix A. N >= 0.``` [in] KD ``` KD is INTEGER The number of superdiagonals or subdiagonals of the triangular band matrix A. KD >= 0.``` [in] NRHS ``` NRHS is INTEGER The number of right hand sides, i.e., the number of columns of the matrices B and X. 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 matrix B.``` [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 solution matrix X.``` [in] LDX ``` LDX is INTEGER The leading dimension of the array X. LDX >= max(1,N).``` [out] FERR ``` FERR is DOUBLE PRECISION array, dimension (NRHS) The estimated forward error bound for each solution vector X(j) (the j-th column of the solution matrix X). If XTRUE is the true solution corresponding to X(j), FERR(j) is an estimated upper bound for the magnitude of the largest element in (X(j) - XTRUE) divided by the magnitude of the largest element in X(j). The estimate is as reliable as the estimate for RCOND, and is almost always a slight overestimate of the true error.``` [out] BERR ``` BERR is DOUBLE PRECISION array, dimension (NRHS) The componentwise relative backward error of each solution vector X(j) (i.e., the smallest relative change in any element of A or B that makes X(j) an exact solution).``` [out] WORK ` WORK is COMPLEX*16 array, dimension (2*N)` [out] RWORK ` RWORK is DOUBLE PRECISION array, dimension (N)` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value```

Definition at line 186 of file ztbrfs.f.

188 *
189 * -- LAPACK computational routine --
190 * -- LAPACK is a software package provided by Univ. of Tennessee, --
191 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
192 *
193 * .. Scalar Arguments ..
194  CHARACTER DIAG, TRANS, UPLO
195  INTEGER INFO, KD, LDAB, LDB, LDX, N, NRHS
196 * ..
197 * .. Array Arguments ..
198  DOUBLE PRECISION BERR( * ), FERR( * ), RWORK( * )
199  COMPLEX*16 AB( LDAB, * ), B( LDB, * ), WORK( * ),
200  \$ X( LDX, * )
201 * ..
202 *
203 * =====================================================================
204 *
205 * .. Parameters ..
206  DOUBLE PRECISION ZERO
207  parameter( zero = 0.0d+0 )
208  COMPLEX*16 ONE
209  parameter( one = ( 1.0d+0, 0.0d+0 ) )
210 * ..
211 * .. Local Scalars ..
212  LOGICAL NOTRAN, NOUNIT, UPPER
213  CHARACTER TRANSN, TRANST
214  INTEGER I, J, K, KASE, NZ
215  DOUBLE PRECISION EPS, LSTRES, S, SAFE1, SAFE2, SAFMIN, XK
216  COMPLEX*16 ZDUM
217 * ..
218 * .. Local Arrays ..
219  INTEGER ISAVE( 3 )
220 * ..
221 * .. External Subroutines ..
222  EXTERNAL xerbla, zaxpy, zcopy, zlacn2, ztbmv, ztbsv
223 * ..
224 * .. Intrinsic Functions ..
225  INTRINSIC abs, dble, dimag, max, min
226 * ..
227 * .. External Functions ..
228  LOGICAL LSAME
229  DOUBLE PRECISION DLAMCH
230  EXTERNAL lsame, dlamch
231 * ..
232 * .. Statement Functions ..
233  DOUBLE PRECISION CABS1
234 * ..
235 * .. Statement Function definitions ..
236  cabs1( zdum ) = abs( dble( zdum ) ) + abs( dimag( zdum ) )
237 * ..
238 * .. Executable Statements ..
239 *
240 * Test the input parameters.
241 *
242  info = 0
243  upper = lsame( uplo, 'U' )
244  notran = lsame( trans, 'N' )
245  nounit = lsame( diag, 'N' )
246 *
247  IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
248  info = -1
249  ELSE IF( .NOT.notran .AND. .NOT.lsame( trans, 'T' ) .AND. .NOT.
250  \$ lsame( trans, 'C' ) ) THEN
251  info = -2
252  ELSE IF( .NOT.nounit .AND. .NOT.lsame( diag, 'U' ) ) THEN
253  info = -3
254  ELSE IF( n.LT.0 ) THEN
255  info = -4
256  ELSE IF( kd.LT.0 ) THEN
257  info = -5
258  ELSE IF( nrhs.LT.0 ) THEN
259  info = -6
260  ELSE IF( ldab.LT.kd+1 ) THEN
261  info = -8
262  ELSE IF( ldb.LT.max( 1, n ) ) THEN
263  info = -10
264  ELSE IF( ldx.LT.max( 1, n ) ) THEN
265  info = -12
266  END IF
267  IF( info.NE.0 ) THEN
268  CALL xerbla( 'ZTBRFS', -info )
269  RETURN
270  END IF
271 *
272 * Quick return if possible
273 *
274  IF( n.EQ.0 .OR. nrhs.EQ.0 ) THEN
275  DO 10 j = 1, nrhs
276  ferr( j ) = zero
277  berr( j ) = zero
278  10 CONTINUE
279  RETURN
280  END IF
281 *
282  IF( notran ) THEN
283  transn = 'N'
284  transt = 'C'
285  ELSE
286  transn = 'C'
287  transt = 'N'
288  END IF
289 *
290 * NZ = maximum number of nonzero elements in each row of A, plus 1
291 *
292  nz = kd + 2
293  eps = dlamch( 'Epsilon' )
294  safmin = dlamch( 'Safe minimum' )
295  safe1 = nz*safmin
296  safe2 = safe1 / eps
297 *
298 * Do for each right hand side
299 *
300  DO 250 j = 1, nrhs
301 *
302 * Compute residual R = B - op(A) * X,
303 * where op(A) = A, A**T, or A**H, depending on TRANS.
304 *
305  CALL zcopy( n, x( 1, j ), 1, work, 1 )
306  CALL ztbmv( uplo, trans, diag, n, kd, ab, ldab, work, 1 )
307  CALL zaxpy( n, -one, b( 1, j ), 1, work, 1 )
308 *
309 * Compute componentwise relative backward error from formula
310 *
311 * max(i) ( abs(R(i)) / ( abs(op(A))*abs(X) + abs(B) )(i) )
312 *
313 * where abs(Z) is the componentwise absolute value of the matrix
314 * or vector Z. If the i-th component of the denominator is less
315 * than SAFE2, then SAFE1 is added to the i-th components of the
316 * numerator and denominator before dividing.
317 *
318  DO 20 i = 1, n
319  rwork( i ) = cabs1( b( i, j ) )
320  20 CONTINUE
321 *
322  IF( notran ) THEN
323 *
324 * Compute abs(A)*abs(X) + abs(B).
325 *
326  IF( upper ) THEN
327  IF( nounit ) THEN
328  DO 40 k = 1, n
329  xk = cabs1( x( k, j ) )
330  DO 30 i = max( 1, k-kd ), k
331  rwork( i ) = rwork( i ) +
332  \$ cabs1( ab( kd+1+i-k, k ) )*xk
333  30 CONTINUE
334  40 CONTINUE
335  ELSE
336  DO 60 k = 1, n
337  xk = cabs1( x( k, j ) )
338  DO 50 i = max( 1, k-kd ), k - 1
339  rwork( i ) = rwork( i ) +
340  \$ cabs1( ab( kd+1+i-k, k ) )*xk
341  50 CONTINUE
342  rwork( k ) = rwork( k ) + xk
343  60 CONTINUE
344  END IF
345  ELSE
346  IF( nounit ) THEN
347  DO 80 k = 1, n
348  xk = cabs1( x( k, j ) )
349  DO 70 i = k, min( n, k+kd )
350  rwork( i ) = rwork( i ) +
351  \$ cabs1( ab( 1+i-k, k ) )*xk
352  70 CONTINUE
353  80 CONTINUE
354  ELSE
355  DO 100 k = 1, n
356  xk = cabs1( x( k, j ) )
357  DO 90 i = k + 1, min( n, k+kd )
358  rwork( i ) = rwork( i ) +
359  \$ cabs1( ab( 1+i-k, k ) )*xk
360  90 CONTINUE
361  rwork( k ) = rwork( k ) + xk
362  100 CONTINUE
363  END IF
364  END IF
365  ELSE
366 *
367 * Compute abs(A**H)*abs(X) + abs(B).
368 *
369  IF( upper ) THEN
370  IF( nounit ) THEN
371  DO 120 k = 1, n
372  s = zero
373  DO 110 i = max( 1, k-kd ), k
374  s = s + cabs1( ab( kd+1+i-k, k ) )*
375  \$ cabs1( x( i, j ) )
376  110 CONTINUE
377  rwork( k ) = rwork( k ) + s
378  120 CONTINUE
379  ELSE
380  DO 140 k = 1, n
381  s = cabs1( x( k, j ) )
382  DO 130 i = max( 1, k-kd ), k - 1
383  s = s + cabs1( ab( kd+1+i-k, k ) )*
384  \$ cabs1( x( i, j ) )
385  130 CONTINUE
386  rwork( k ) = rwork( k ) + s
387  140 CONTINUE
388  END IF
389  ELSE
390  IF( nounit ) THEN
391  DO 160 k = 1, n
392  s = zero
393  DO 150 i = k, min( n, k+kd )
394  s = s + cabs1( ab( 1+i-k, k ) )*
395  \$ cabs1( x( i, j ) )
396  150 CONTINUE
397  rwork( k ) = rwork( k ) + s
398  160 CONTINUE
399  ELSE
400  DO 180 k = 1, n
401  s = cabs1( x( k, j ) )
402  DO 170 i = k + 1, min( n, k+kd )
403  s = s + cabs1( ab( 1+i-k, k ) )*
404  \$ cabs1( x( i, j ) )
405  170 CONTINUE
406  rwork( k ) = rwork( k ) + s
407  180 CONTINUE
408  END IF
409  END IF
410  END IF
411  s = zero
412  DO 190 i = 1, n
413  IF( rwork( i ).GT.safe2 ) THEN
414  s = max( s, cabs1( work( i ) ) / rwork( i ) )
415  ELSE
416  s = max( s, ( cabs1( work( i ) )+safe1 ) /
417  \$ ( rwork( i )+safe1 ) )
418  END IF
419  190 CONTINUE
420  berr( j ) = s
421 *
422 * Bound error from formula
423 *
424 * norm(X - XTRUE) / norm(X) .le. FERR =
425 * norm( abs(inv(op(A)))*
426 * ( abs(R) + NZ*EPS*( abs(op(A))*abs(X)+abs(B) ))) / norm(X)
427 *
428 * where
429 * norm(Z) is the magnitude of the largest component of Z
430 * inv(op(A)) is the inverse of op(A)
431 * abs(Z) is the componentwise absolute value of the matrix or
432 * vector Z
433 * NZ is the maximum number of nonzeros in any row of A, plus 1
434 * EPS is machine epsilon
435 *
436 * The i-th component of abs(R)+NZ*EPS*(abs(op(A))*abs(X)+abs(B))
437 * is incremented by SAFE1 if the i-th component of
438 * abs(op(A))*abs(X) + abs(B) is less than SAFE2.
439 *
440 * Use ZLACN2 to estimate the infinity-norm of the matrix
441 * inv(op(A)) * diag(W),
442 * where W = abs(R) + NZ*EPS*( abs(op(A))*abs(X)+abs(B) )))
443 *
444  DO 200 i = 1, n
445  IF( rwork( i ).GT.safe2 ) THEN
446  rwork( i ) = cabs1( work( i ) ) + nz*eps*rwork( i )
447  ELSE
448  rwork( i ) = cabs1( work( i ) ) + nz*eps*rwork( i ) +
449  \$ safe1
450  END IF
451  200 CONTINUE
452 *
453  kase = 0
454  210 CONTINUE
455  CALL zlacn2( n, work( n+1 ), work, ferr( j ), kase, isave )
456  IF( kase.NE.0 ) THEN
457  IF( kase.EQ.1 ) THEN
458 *
459 * Multiply by diag(W)*inv(op(A)**H).
460 *
461  CALL ztbsv( uplo, transt, diag, n, kd, ab, ldab, work,
462  \$ 1 )
463  DO 220 i = 1, n
464  work( i ) = rwork( i )*work( i )
465  220 CONTINUE
466  ELSE
467 *
468 * Multiply by inv(op(A))*diag(W).
469 *
470  DO 230 i = 1, n
471  work( i ) = rwork( i )*work( i )
472  230 CONTINUE
473  CALL ztbsv( uplo, transn, diag, n, kd, ab, ldab, work,
474  \$ 1 )
475  END IF
476  GO TO 210
477  END IF
478 *
479 * Normalize error.
480 *
481  lstres = zero
482  DO 240 i = 1, n
483  lstres = max( lstres, cabs1( x( i, j ) ) )
484  240 CONTINUE
485  IF( lstres.NE.zero )
486  \$ ferr( j ) = ferr( j ) / lstres
487 *
488  250 CONTINUE
489 *
490  RETURN
491 *
492 * End of ZTBRFS
493 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine zaxpy(N, ZA, ZX, INCX, ZY, INCY)
ZAXPY
Definition: zaxpy.f:88
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
Definition: zcopy.f:81
subroutine ztbsv(UPLO, TRANS, DIAG, N, K, A, LDA, X, INCX)
ZTBSV
Definition: ztbsv.f:189
subroutine ztbmv(UPLO, TRANS, DIAG, N, K, A, LDA, X, INCX)
ZTBMV
Definition: ztbmv.f:186
subroutine zlacn2(N, V, X, EST, KASE, ISAVE)
ZLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vec...
Definition: zlacn2.f:133
Here is the call graph for this function:
Here is the caller graph for this function: