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

## ◆ 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: