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

◆ dgbt05()

subroutine dgbt05 ( character trans,
integer n,
integer kl,
integer ku,
integer nrhs,
double precision, dimension( ldab, * ) ab,
integer ldab,
double precision, dimension( ldb, * ) b,
integer ldb,
double precision, dimension( ldx, * ) x,
integer ldx,
double precision, dimension( ldxact, * ) xact,
integer ldxact,
double precision, dimension( * ) ferr,
double precision, dimension( * ) berr,
double precision, dimension( * ) reslts )

DGBT05

Purpose:
!>
!> DGBT05 tests the error bounds from iterative refinement for the
!> computed solution to a system of equations op(A)*X = B, where A is a
!> general band matrix of order n with kl subdiagonals and ku
!> superdiagonals and op(A) = A or A**T, depending on TRANS.
!>
!> 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(op(A))*abs(X) +abs(b))_i )
!>             and NZ = max. number of nonzeros in any row of A, plus 1
!>
Parameters
[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 = Transpose)
!>
[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]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 columns of the matrices X, B, and XACT.
!>          NRHS >= 0.
!>
[in]AB
!>          AB is DOUBLE PRECISION array, dimension (LDAB,N)
!>          The original band matrix A, stored in rows 1 to KL+KU+1.
!>          The j-th column of A is stored in the j-th column of the
!>          array AB as follows:
!>          AB(ku+1+i-j,j) = A(i,j) for max(1,j-ku)<=i<=min(n,j+kl).
!>
[in]LDAB
!>          LDAB is INTEGER
!>          The leading dimension of the array AB.  LDAB >= KL+KU+1.
!>
[in]B
!>          B is DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 + (*) )
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 174 of file dgbt05.f.

176*
177* -- LAPACK test routine --
178* -- LAPACK is a software package provided by Univ. of Tennessee, --
179* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
180*
181* .. Scalar Arguments ..
182 CHARACTER TRANS
183 INTEGER KL, KU, LDAB, LDB, LDX, LDXACT, N, NRHS
184* ..
185* .. Array Arguments ..
186 DOUBLE PRECISION AB( LDAB, * ), B( LDB, * ), BERR( * ),
187 $ FERR( * ), RESLTS( * ), X( LDX, * ),
188 $ XACT( LDXACT, * )
189* ..
190*
191* =====================================================================
192*
193* .. Parameters ..
194 DOUBLE PRECISION ZERO, ONE
195 parameter( zero = 0.0d+0, one = 1.0d+0 )
196* ..
197* .. Local Scalars ..
198 LOGICAL NOTRAN
199 INTEGER I, IMAX, J, K, NZ
200 DOUBLE PRECISION AXBI, DIFF, EPS, ERRBND, OVFL, TMP, UNFL, XNORM
201* ..
202* .. External Functions ..
203 LOGICAL LSAME
204 INTEGER IDAMAX
205 DOUBLE PRECISION DLAMCH
206 EXTERNAL lsame, idamax, dlamch
207* ..
208* .. Intrinsic Functions ..
209 INTRINSIC abs, max, min
210* ..
211* .. Executable Statements ..
212*
213* Quick exit if N = 0 or NRHS = 0.
214*
215 IF( n.LE.0 .OR. nrhs.LE.0 ) THEN
216 reslts( 1 ) = zero
217 reslts( 2 ) = zero
218 RETURN
219 END IF
220*
221 eps = dlamch( 'Epsilon' )
222 unfl = dlamch( 'Safe minimum' )
223 ovfl = one / unfl
224 notran = lsame( trans, 'N' )
225 nz = min( kl+ku+2, n+1 )
226*
227* Test 1: Compute the maximum of
228* norm(X - XACT) / ( norm(X) * FERR )
229* over all the vectors X and XACT using the infinity-norm.
230*
231 errbnd = zero
232 DO 30 j = 1, nrhs
233 imax = idamax( n, x( 1, j ), 1 )
234 xnorm = max( abs( x( imax, j ) ), unfl )
235 diff = zero
236 DO 10 i = 1, n
237 diff = max( diff, abs( x( i, j )-xact( i, j ) ) )
238 10 CONTINUE
239*
240 IF( xnorm.GT.one ) THEN
241 GO TO 20
242 ELSE IF( diff.LE.ovfl*xnorm ) THEN
243 GO TO 20
244 ELSE
245 errbnd = one / eps
246 GO TO 30
247 END IF
248*
249 20 CONTINUE
250 IF( diff / xnorm.LE.ferr( j ) ) THEN
251 errbnd = max( errbnd, ( diff / xnorm ) / ferr( j ) )
252 ELSE
253 errbnd = one / eps
254 END IF
255 30 CONTINUE
256 reslts( 1 ) = errbnd
257*
258* Test 2: Compute the maximum of BERR / ( NZ*EPS + (*) ), where
259* (*) = NZ*UNFL / (min_i (abs(op(A))*abs(X) +abs(b))_i )
260*
261 DO 70 k = 1, nrhs
262 DO 60 i = 1, n
263 tmp = abs( b( i, k ) )
264 IF( notran ) THEN
265 DO 40 j = max( i-kl, 1 ), min( i+ku, n )
266 tmp = tmp + abs( ab( ku+1+i-j, j ) )*abs( x( j, k ) )
267 40 CONTINUE
268 ELSE
269 DO 50 j = max( i-ku, 1 ), min( i+kl, n )
270 tmp = tmp + abs( ab( ku+1+j-i, i ) )*abs( x( j, k ) )
271 50 CONTINUE
272 END IF
273 IF( i.EQ.1 ) THEN
274 axbi = tmp
275 ELSE
276 axbi = min( axbi, tmp )
277 END IF
278 60 CONTINUE
279 tmp = berr( k ) / ( nz*eps+nz*unfl / max( axbi, nz*unfl ) )
280 IF( k.EQ.1 ) THEN
281 reslts( 2 ) = tmp
282 ELSE
283 reslts( 2 ) = max( reslts( 2 ), tmp )
284 END IF
285 70 CONTINUE
286*
287 RETURN
288*
289* End of DGBT05
290*
idamax
integer function idamax(n, dx, incx)
IDAMAX
Definition idamax.f:71
dlamch
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
Here is the caller graph for this function: