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

◆ dla_gbrfsx_extended()

subroutine dla_gbrfsx_extended ( integer prec_type,
integer trans_type,
integer n,
integer kl,
integer ku,
integer nrhs,
double precision, dimension( ldab, * ) ab,
integer ldab,
double precision, dimension( ldafb, * ) afb,
integer ldafb,
integer, dimension( * ) ipiv,
logical colequ,
double precision, dimension( * ) c,
double precision, dimension( ldb, * ) b,
integer ldb,
double precision, dimension( ldy, * ) y,
integer ldy,
double precision, dimension(*) berr_out,
integer n_norms,
double precision, dimension( nrhs, * ) err_bnds_norm,
double precision, dimension( nrhs, * ) err_bnds_comp,
double precision, dimension(*) res,
double precision, dimension(*) ayb,
double precision, dimension(*) dy,
double precision, dimension(*) y_tail,
double precision rcond,
integer ithresh,
double precision rthresh,
double precision dz_ub,
logical ignore_cwise,
integer info )

DLA_GBRFSX_EXTENDED improves the computed solution to a system of linear equations for general banded matrices by performing extra-precise iterative refinement and provides error bounds and backward error estimates for the solution.

Download DLA_GBRFSX_EXTENDED + dependencies [TGZ] [ZIP] [TXT]

Purpose:
!>
!>
!> DLA_GBRFSX_EXTENDED improves the computed solution to a system of
!> linear equations by performing extra-precise iterative refinement
!> and provides error bounds and backward error estimates for the solution.
!> This subroutine is called by DGBRFSX to perform iterative refinement.
!> In addition to normwise error bound, the code provides maximum
!> componentwise error bound if possible. See comments for ERR_BNDS_NORM
!> and ERR_BNDS_COMP for details of the error bounds. Note that this
!> subroutine is only responsible for setting the second fields of
!> ERR_BNDS_NORM and ERR_BNDS_COMP.
!> 
Parameters
[in]PREC_TYPE
!>          PREC_TYPE is INTEGER
!>     Specifies the intermediate precision to be used in refinement.
!>     The value is defined by ILAPREC(P) where P is a CHARACTER and P
!>          = 'S':  Single
!>          = 'D':  Double
!>          = 'I':  Indigenous
!>          = 'X' or 'E':  Extra
!> 
[in]TRANS_TYPE
!>          TRANS_TYPE is INTEGER
!>     Specifies the transposition operation on A.
!>     The value is defined by ILATRANS(T) where T is a CHARACTER and T
!>          = 'N':  No transpose
!>          = 'T':  Transpose
!>          = 'C':  Conjugate transpose
!> 
[in]N
!>          N is INTEGER
!>     The number of linear equations, i.e., the order of the
!>     matrix A.  N >= 0.
!> 
[in]KL
!>          KL is INTEGER
!>     The number of subdiagonals within the band of A.  KL >= 0.
!> 
[in]KU
!>          KU is INTEGER
!>     The number of superdiagonals within the band of A.  KU >= 0
!> 
[in]NRHS
!>          NRHS is INTEGER
!>     The number of right-hand-sides, i.e., the number of columns of the
!>     matrix B.
!> 
[in]AB
!>          AB is DOUBLE PRECISION array, dimension (LDAB,N)
!>          On entry, the N-by-N matrix AB.
!> 
[in]LDAB
!>          LDAB is INTEGER
!>          The leading dimension of the array AB.  LDBA >= max(1,N).
!> 
[in]AFB
!>          AFB is DOUBLE PRECISION array, dimension (LDAFB,N)
!>     The factors L and U from the factorization
!>     A = P*L*U as computed by DGBTRF.
!> 
[in]LDAFB
!>          LDAFB is INTEGER
!>     The leading dimension of the array AF.  LDAFB >= max(1,N).
!> 
[in]IPIV
!>          IPIV is INTEGER array, dimension (N)
!>     The pivot indices from the factorization A = P*L*U
!>     as computed by DGBTRF; row i of the matrix was interchanged
!>     with row IPIV(i).
!> 
[in]COLEQU
!>          COLEQU is LOGICAL
!>     If .TRUE. then column equilibration was done to A before calling
!>     this routine. This is needed to compute the solution and error
!>     bounds correctly.
!> 
[in]C
!>          C is DOUBLE PRECISION array, dimension (N)
!>     The column scale factors for A. If COLEQU = .FALSE., C
!>     is not accessed. If C is input, each element of C should be a power
!>     of the radix to ensure a reliable solution and error estimates.
!>     Scaling by powers of the radix does not cause rounding errors unless
!>     the result underflows or overflows. Rounding errors during scaling
!>     lead to refining with a matrix that is not equivalent to the
!>     input matrix, producing error estimates that may not be
!>     reliable.
!> 
[in]B
!>          B is DOUBLE PRECISION 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,out]Y
!>          Y is DOUBLE PRECISION array, dimension (LDY,NRHS)
!>     On entry, the solution matrix X, as computed by DGBTRS.
!>     On exit, the improved solution matrix Y.
!> 
[in]LDY
!>          LDY is INTEGER
!>     The leading dimension of the array Y.  LDY >= max(1,N).
!> 
[out]BERR_OUT
!>          BERR_OUT is DOUBLE PRECISION array, dimension (NRHS)
!>     On exit, BERR_OUT(j) contains the componentwise relative backward
!>     error for right-hand-side j from the formula
!>         max(i) ( abs(RES(i)) / ( abs(op(A_s))*abs(Y) + abs(B_s) )(i) )
!>     where abs(Z) is the componentwise absolute value of the matrix
!>     or vector Z. This is computed by DLA_LIN_BERR.
!> 
[in]N_NORMS
!>          N_NORMS is INTEGER
!>     Determines which error bounds to return (see ERR_BNDS_NORM
!>     and ERR_BNDS_COMP).
!>     If N_NORMS >= 1 return normwise error bounds.
!>     If N_NORMS >= 2 return componentwise error bounds.
!> 
[in,out]ERR_BNDS_NORM
!>          ERR_BNDS_NORM is DOUBLE PRECISION array, dimension (NRHS, N_ERR_BNDS)
!>     For each right-hand side, this array contains information about
!>     various error bounds and condition numbers corresponding to the
!>     normwise relative error, which is defined as follows:
!>
!>     Normwise relative error in the ith solution vector:
!>             max_j (abs(XTRUE(j,i) - X(j,i)))
!>            ------------------------------
!>                  max_j abs(X(j,i))
!>
!>     The array is indexed by the type of error information as described
!>     below. There currently are up to three pieces of information
!>     returned.
!>
!>     The first index in ERR_BNDS_NORM(i,:) corresponds to the ith
!>     right-hand side.
!>
!>     The second index in ERR_BNDS_NORM(:,err) contains the following
!>     three fields:
!>     err = 1  boolean. Trust the answer if the
!>              reciprocal condition number is less than the threshold
!>              sqrt(n) * slamch('Epsilon').
!>
!>     err = 2  error bound: The estimated forward error,
!>              almost certainly within a factor of 10 of the true error
!>              so long as the next entry is greater than the threshold
!>              sqrt(n) * slamch('Epsilon'). This error bound should only
!>              be trusted if the previous boolean is true.
!>
!>     err = 3  Reciprocal condition number: Estimated normwise
!>              reciprocal condition number.  Compared with the threshold
!>              sqrt(n) * slamch('Epsilon') to determine if the error
!>              estimate is . These reciprocal condition
!>              numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some
!>              appropriately scaled matrix Z.
!>              Let Z = S*A, where S scales each row by a power of the
!>              radix so all absolute row sums of Z are approximately 1.
!>
!>     This subroutine is only responsible for setting the second field
!>     above.
!>     See Lapack Working Note 165 for further details and extra
!>     cautions.
!> 
[in,out]ERR_BNDS_COMP
!>          ERR_BNDS_COMP is DOUBLE PRECISION array, dimension (NRHS, N_ERR_BNDS)
!>     For each right-hand side, this array contains information about
!>     various error bounds and condition numbers corresponding to the
!>     componentwise relative error, which is defined as follows:
!>
!>     Componentwise relative error in the ith solution vector:
!>                    abs(XTRUE(j,i) - X(j,i))
!>             max_j ----------------------
!>                         abs(X(j,i))
!>
!>     The array is indexed by the right-hand side i (on which the
!>     componentwise relative error depends), and the type of error
!>     information as described below. There currently are up to three
!>     pieces of information returned for each right-hand side. If
!>     componentwise accuracy is not requested (PARAMS(3) = 0.0), then
!>     ERR_BNDS_COMP is not accessed.  If N_ERR_BNDS < 3, then at most
!>     the first (:,N_ERR_BNDS) entries are returned.
!>
!>     The first index in ERR_BNDS_COMP(i,:) corresponds to the ith
!>     right-hand side.
!>
!>     The second index in ERR_BNDS_COMP(:,err) contains the following
!>     three fields:
!>     err = 1  boolean. Trust the answer if the
!>              reciprocal condition number is less than the threshold
!>              sqrt(n) * slamch('Epsilon').
!>
!>     err = 2  error bound: The estimated forward error,
!>              almost certainly within a factor of 10 of the true error
!>              so long as the next entry is greater than the threshold
!>              sqrt(n) * slamch('Epsilon'). This error bound should only
!>              be trusted if the previous boolean is true.
!>
!>     err = 3  Reciprocal condition number: Estimated componentwise
!>              reciprocal condition number.  Compared with the threshold
!>              sqrt(n) * slamch('Epsilon') to determine if the error
!>              estimate is . These reciprocal condition
!>              numbers are 1 / (norm(Z^{-1},inf) * norm(Z,inf)) for some
!>              appropriately scaled matrix Z.
!>              Let Z = S*(A*diag(x)), where x is the solution for the
!>              current right-hand side and S scales each row of
!>              A*diag(x) by a power of the radix so all absolute row
!>              sums of Z are approximately 1.
!>
!>     This subroutine is only responsible for setting the second field
!>     above.
!>     See Lapack Working Note 165 for further details and extra
!>     cautions.
!> 
[in]RES
!>          RES is DOUBLE PRECISION array, dimension (N)
!>     Workspace to hold the intermediate residual.
!> 
[in]AYB
!>          AYB is DOUBLE PRECISION array, dimension (N)
!>     Workspace. This can be the same workspace passed for Y_TAIL.
!> 
[in]DY
!>          DY is DOUBLE PRECISION array, dimension (N)
!>     Workspace to hold the intermediate solution.
!> 
[in]Y_TAIL
!>          Y_TAIL is DOUBLE PRECISION array, dimension (N)
!>     Workspace to hold the trailing bits of the intermediate solution.
!> 
[in]RCOND
!>          RCOND is DOUBLE PRECISION
!>     Reciprocal scaled condition number.  This is an estimate of the
!>     reciprocal Skeel condition number of the matrix A after
!>     equilibration (if done).  If this is less than the machine
!>     precision (in particular, if it is zero), the matrix is singular
!>     to working precision.  Note that the error may still be small even
!>     if this number is very small and the matrix appears ill-
!>     conditioned.
!> 
[in]ITHRESH
!>          ITHRESH is INTEGER
!>     The maximum number of residual computations allowed for
!>     refinement. The default is 10. For 'aggressive' set to 100 to
!>     permit convergence using approximate factorizations or
!>     factorizations other than LU. If the factorization uses a
!>     technique other than Gaussian elimination, the guarantees in
!>     ERR_BNDS_NORM and ERR_BNDS_COMP may no longer be trustworthy.
!> 
[in]RTHRESH
!>          RTHRESH is DOUBLE PRECISION
!>     Determines when to stop refinement if the error estimate stops
!>     decreasing. Refinement will stop when the next solution no longer
!>     satisfies norm(dx_{i+1}) < RTHRESH * norm(dx_i) where norm(Z) is
!>     the infinity norm of Z. RTHRESH satisfies 0 < RTHRESH <= 1. The
!>     default value is 0.5. For 'aggressive' set to 0.9 to permit
!>     convergence on extremely ill-conditioned matrices. See LAWN 165
!>     for more details.
!> 
[in]DZ_UB
!>          DZ_UB is DOUBLE PRECISION
!>     Determines when to start considering componentwise convergence.
!>     Componentwise convergence is only considered after each component
!>     of the solution Y is stable, which we define as the relative
!>     change in each component being less than DZ_UB. The default value
!>     is 0.25, requiring the first bit to be stable. See LAWN 165 for
!>     more details.
!> 
[in]IGNORE_CWISE
!>          IGNORE_CWISE is LOGICAL
!>     If .TRUE. then ignore componentwise convergence. Default value
!>     is .FALSE..
!> 
[out]INFO
!>          INFO is INTEGER
!>       = 0:  Successful exit.
!>       < 0:  if INFO = -i, the ith argument to DGBTRS had an illegal
!>             value
!> 
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 402 of file dla_gbrfsx_extended.f.

410*
411* -- LAPACK computational routine --
412* -- LAPACK is a software package provided by Univ. of Tennessee, --
413* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
414*
415* .. Scalar Arguments ..
416 INTEGER INFO, LDAB, LDAFB, LDB, LDY, N, KL, KU, NRHS,
417 $ PREC_TYPE, TRANS_TYPE, N_NORMS, ITHRESH
418 LOGICAL COLEQU, IGNORE_CWISE
419 DOUBLE PRECISION RTHRESH, DZ_UB
420* ..
421* .. Array Arguments ..
422 INTEGER IPIV( * )
423 DOUBLE PRECISION AB( LDAB, * ), AFB( LDAFB, * ), B( LDB, * ),
424 $ Y( LDY, * ), RES(*), DY(*), Y_TAIL(*)
425 DOUBLE PRECISION C( * ), AYB(*), RCOND, BERR_OUT(*),
426 $ ERR_BNDS_NORM( NRHS, * ),
427 $ ERR_BNDS_COMP( NRHS, * )
428* ..
429*
430* =====================================================================
431*
432* .. Local Scalars ..
433 CHARACTER TRANS
434 INTEGER CNT, I, J, M, X_STATE, Z_STATE, Y_PREC_STATE
435 DOUBLE PRECISION YK, DYK, YMIN, NORMY, NORMX, NORMDX, DXRAT,
436 $ DZRAT, PREVNORMDX, PREV_DZ_Z, DXRATMAX,
437 $ DZRATMAX, DX_X, DZ_Z, FINAL_DX_X, FINAL_DZ_Z,
438 $ EPS, HUGEVAL, INCR_THRESH
439 LOGICAL INCR_PREC
440* ..
441* .. Parameters ..
442 INTEGER UNSTABLE_STATE, WORKING_STATE, CONV_STATE,
443 $ NOPROG_STATE, BASE_RESIDUAL, EXTRA_RESIDUAL,
444 $ EXTRA_Y
445 parameter( unstable_state = 0, working_state = 1,
446 $ conv_state = 2, noprog_state = 3 )
447 parameter( base_residual = 0, extra_residual = 1,
448 $ extra_y = 2 )
449 INTEGER FINAL_NRM_ERR_I, FINAL_CMP_ERR_I, BERR_I
450 INTEGER RCOND_I, NRM_RCOND_I, NRM_ERR_I, CMP_RCOND_I
451 INTEGER CMP_ERR_I, PIV_GROWTH_I
452 parameter( final_nrm_err_i = 1, final_cmp_err_i = 2,
453 $ berr_i = 3 )
454 parameter( rcond_i = 4, nrm_rcond_i = 5, nrm_err_i = 6 )
455 parameter( cmp_rcond_i = 7, cmp_err_i = 8,
456 $ piv_growth_i = 9 )
457 INTEGER LA_LINRX_ITREF_I, LA_LINRX_ITHRESH_I,
458 $ LA_LINRX_CWISE_I
459 parameter( la_linrx_itref_i = 1,
460 $ la_linrx_ithresh_i = 2 )
461 parameter( la_linrx_cwise_i = 3 )
462 INTEGER LA_LINRX_TRUST_I, LA_LINRX_ERR_I,
463 $ LA_LINRX_RCOND_I
464 parameter( la_linrx_trust_i = 1, la_linrx_err_i = 2 )
465 parameter( la_linrx_rcond_i = 3 )
466* ..
467* .. External Subroutines ..
468 EXTERNAL daxpy, dcopy, dgbtrs, dgbmv,
469 $ blas_dgbmv_x,
470 $ blas_dgbmv2_x, dla_gbamv, dla_wwaddw, dlamch,
472 DOUBLE PRECISION DLAMCH
473 CHARACTER CHLA_TRANSTYPE
474* ..
475* .. Intrinsic Functions ..
476 INTRINSIC abs, max, min
477* ..
478* .. Executable Statements ..
479*
480 IF (info.NE.0) RETURN
481 trans = chla_transtype(trans_type)
482 eps = dlamch( 'Epsilon' )
483 hugeval = dlamch( 'Overflow' )
484* Force HUGEVAL to Inf
485 hugeval = hugeval * hugeval
486* Using HUGEVAL may lead to spurious underflows.
487 incr_thresh = dble( n ) * eps
488 m = kl+ku+1
489
490 DO j = 1, nrhs
491 y_prec_state = extra_residual
492 IF ( y_prec_state .EQ. extra_y ) THEN
493 DO i = 1, n
494 y_tail( i ) = 0.0d+0
495 END DO
496 END IF
497
498 dxrat = 0.0d+0
499 dxratmax = 0.0d+0
500 dzrat = 0.0d+0
501 dzratmax = 0.0d+0
502 final_dx_x = hugeval
503 final_dz_z = hugeval
504 prevnormdx = hugeval
505 prev_dz_z = hugeval
506 dz_z = hugeval
507 dx_x = hugeval
508
509 x_state = working_state
510 z_state = unstable_state
511 incr_prec = .false.
512
513 DO cnt = 1, ithresh
514*
515* Compute residual RES = B_s - op(A_s) * Y,
516* op(A) = A, A**T, or A**H depending on TRANS (and type).
517*
518 CALL dcopy( n, b( 1, j ), 1, res, 1 )
519 IF ( y_prec_state .EQ. base_residual ) THEN
520 CALL dgbmv( trans, m, n, kl, ku, -1.0d+0, ab, ldab,
521 $ y( 1, j ), 1, 1.0d+0, res, 1 )
522 ELSE IF ( y_prec_state .EQ. extra_residual ) THEN
523 CALL blas_dgbmv_x( trans_type, n, n, kl, ku,
524 $ -1.0d+0, ab, ldab, y( 1, j ), 1, 1.0d+0, res, 1,
525 $ prec_type )
526 ELSE
527 CALL blas_dgbmv2_x( trans_type, n, n, kl, ku, -1.0d+0,
528 $ ab, ldab, y( 1, j ), y_tail, 1, 1.0d+0, res, 1,
529 $ prec_type )
530 END IF
531
532! XXX: RES is no longer needed.
533 CALL dcopy( n, res, 1, dy, 1 )
534 CALL dgbtrs( trans, n, kl, ku, 1, afb, ldafb, ipiv, dy,
535 $ n,
536 $ info )
537*
538* Calculate relative changes DX_X, DZ_Z and ratios DXRAT, DZRAT.
539*
540 normx = 0.0d+0
541 normy = 0.0d+0
542 normdx = 0.0d+0
543 dz_z = 0.0d+0
544 ymin = hugeval
545
546 DO i = 1, n
547 yk = abs( y( i, j ) )
548 dyk = abs( dy( i ) )
549
550 IF ( yk .NE. 0.0d+0 ) THEN
551 dz_z = max( dz_z, dyk / yk )
552 ELSE IF ( dyk .NE. 0.0d+0 ) THEN
553 dz_z = hugeval
554 END IF
555
556 ymin = min( ymin, yk )
557
558 normy = max( normy, yk )
559
560 IF ( colequ ) THEN
561 normx = max( normx, yk * c( i ) )
562 normdx = max( normdx, dyk * c( i ) )
563 ELSE
564 normx = normy
565 normdx = max( normdx, dyk )
566 END IF
567 END DO
568
569 IF ( normx .NE. 0.0d+0 ) THEN
570 dx_x = normdx / normx
571 ELSE IF ( normdx .EQ. 0.0d+0 ) THEN
572 dx_x = 0.0d+0
573 ELSE
574 dx_x = hugeval
575 END IF
576
577 dxrat = normdx / prevnormdx
578 dzrat = dz_z / prev_dz_z
579*
580* Check termination criteria.
581*
582 IF ( .NOT.ignore_cwise
583 $ .AND. ymin*rcond .LT. incr_thresh*normy
584 $ .AND. y_prec_state .LT. extra_y )
585 $ incr_prec = .true.
586
587 IF ( x_state .EQ. noprog_state .AND. dxrat .LE. rthresh )
588 $ x_state = working_state
589 IF ( x_state .EQ. working_state ) THEN
590 IF ( dx_x .LE. eps ) THEN
591 x_state = conv_state
592 ELSE IF ( dxrat .GT. rthresh ) THEN
593 IF ( y_prec_state .NE. extra_y ) THEN
594 incr_prec = .true.
595 ELSE
596 x_state = noprog_state
597 END IF
598 ELSE
599 IF ( dxrat .GT. dxratmax ) dxratmax = dxrat
600 END IF
601 IF ( x_state .GT. working_state ) final_dx_x = dx_x
602 END IF
603
604 IF ( z_state .EQ. unstable_state .AND. dz_z .LE. dz_ub )
605 $ z_state = working_state
606 IF ( z_state .EQ. noprog_state .AND. dzrat .LE. rthresh )
607 $ z_state = working_state
608 IF ( z_state .EQ. working_state ) THEN
609 IF ( dz_z .LE. eps ) THEN
610 z_state = conv_state
611 ELSE IF ( dz_z .GT. dz_ub ) THEN
612 z_state = unstable_state
613 dzratmax = 0.0d+0
614 final_dz_z = hugeval
615 ELSE IF ( dzrat .GT. rthresh ) THEN
616 IF ( y_prec_state .NE. extra_y ) THEN
617 incr_prec = .true.
618 ELSE
619 z_state = noprog_state
620 END IF
621 ELSE
622 IF ( dzrat .GT. dzratmax ) dzratmax = dzrat
623 END IF
624 IF ( z_state .GT. working_state ) final_dz_z = dz_z
625 END IF
626*
627* Exit if both normwise and componentwise stopped working,
628* but if componentwise is unstable, let it go at least two
629* iterations.
630*
631 IF ( x_state.NE.working_state ) THEN
632 IF ( ignore_cwise ) GOTO 666
633 IF ( z_state.EQ.noprog_state .OR. z_state.EQ.conv_state )
634 $ GOTO 666
635 IF ( z_state.EQ.unstable_state .AND. cnt.GT.1 ) GOTO 666
636 END IF
637
638 IF ( incr_prec ) THEN
639 incr_prec = .false.
640 y_prec_state = y_prec_state + 1
641 DO i = 1, n
642 y_tail( i ) = 0.0d+0
643 END DO
644 END IF
645
646 prevnormdx = normdx
647 prev_dz_z = dz_z
648*
649* Update solution.
650*
651 IF (y_prec_state .LT. extra_y) THEN
652 CALL daxpy( n, 1.0d+0, dy, 1, y(1,j), 1 )
653 ELSE
654 CALL dla_wwaddw( n, y(1,j), y_tail, dy )
655 END IF
656
657 END DO
658* Target of "IF (Z_STOP .AND. X_STOP)". Sun's f77 won't EXIT.
659 666 CONTINUE
660*
661* Set final_* when cnt hits ithresh.
662*
663 IF ( x_state .EQ. working_state ) final_dx_x = dx_x
664 IF ( z_state .EQ. working_state ) final_dz_z = dz_z
665*
666* Compute error bounds.
667*
668 IF ( n_norms .GE. 1 ) THEN
669 err_bnds_norm( j, la_linrx_err_i ) =
670 $ final_dx_x / (1 - dxratmax)
671 END IF
672 IF (n_norms .GE. 2) THEN
673 err_bnds_comp( j, la_linrx_err_i ) =
674 $ final_dz_z / (1 - dzratmax)
675 END IF
676*
677* Compute componentwise relative backward error from formula
678* max(i) ( abs(R(i)) / ( abs(op(A_s))*abs(Y) + abs(B_s) )(i) )
679* where abs(Z) is the componentwise absolute value of the matrix
680* or vector Z.
681*
682* Compute residual RES = B_s - op(A_s) * Y,
683* op(A) = A, A**T, or A**H depending on TRANS (and type).
684*
685 CALL dcopy( n, b( 1, j ), 1, res, 1 )
686 CALL dgbmv(trans, n, n, kl, ku, -1.0d+0, ab, ldab, y(1,j),
687 $ 1, 1.0d+0, res, 1 )
688
689 DO i = 1, n
690 ayb( i ) = abs( b( i, j ) )
691 END DO
692*
693* Compute abs(op(A_s))*abs(Y) + abs(B_s).
694*
695 CALL dla_gbamv( trans_type, n, n, kl, ku, 1.0d+0,
696 $ ab, ldab, y(1, j), 1, 1.0d+0, ayb, 1 )
697
698 CALL dla_lin_berr( n, n, 1, res, ayb, berr_out( j ) )
699*
700* End of loop for each RHS
701*
702 END DO
703*
704 RETURN
705*
706* End of DLA_GBRFSX_EXTENDED
707*
subroutine daxpy(n, da, dx, incx, dy, incy)
DAXPY
Definition daxpy.f:89
subroutine dcopy(n, dx, incx, dy, incy)
DCOPY
Definition dcopy.f:82
subroutine dgbmv(trans, m, n, kl, ku, alpha, a, lda, x, incx, beta, y, incy)
DGBMV
Definition dgbmv.f:188
subroutine dgbtrs(trans, n, kl, ku, nrhs, ab, ldab, ipiv, b, ldb, info)
DGBTRS
Definition dgbtrs.f:137
subroutine dla_gbamv(trans, m, n, kl, ku, alpha, ab, ldab, x, incx, beta, y, incy)
DLA_GBAMV performs a matrix-vector operation to calculate error bounds.
Definition dla_gbamv.f:185
subroutine dla_lin_berr(n, nz, nrhs, res, ayb, berr)
DLA_LIN_BERR computes a component-wise relative backward error.
character *1 function chla_transtype(trans)
CHLA_TRANSTYPE
subroutine dla_wwaddw(n, x, y, w)
DLA_WWADDW adds a vector into a doubled-single vector.
Definition dla_wwaddw.f:79
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
Here is the call graph for this function:
Here is the caller graph for this function: