LAPACK  3.8.0
LAPACK: Linear Algebra PACKage

◆ spbt05()

subroutine spbt05 ( character  UPLO,
integer  N,
integer  KD,
integer  NRHS,
real, dimension( ldab, * )  AB,
integer  LDAB,
real, dimension( ldb, * )  B,
integer  LDB,
real, dimension( ldx, * )  X,
integer  LDX,
real, dimension( ldxact, * )  XACT,
integer  LDXACT,
real, dimension( * )  FERR,
real, dimension( * )  BERR,
real, dimension( * )  RESLTS 
)

SPBT05

Purpose:
 SPBT05 tests the error bounds from iterative refinement for the
 computed solution to a system of equations A*X = B, where A is a
 symmetric 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 upper or lower triangular part of the
          symmetric matrix A is stored.
          = 'U':  Upper triangular
          = 'L':  Lower 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 REAL array, dimension (LDAB,N)
          The upper or lower triangle of the symmetric 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).
[in]LDAB
          LDAB is INTEGER
          The leading dimension of the array AB.  LDAB >= KD+1.
[in]B
          B is REAL 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 REAL 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 REAL 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.
Date
December 2016

Definition at line 173 of file spbt05.f.

173 *
174 * -- LAPACK test routine (version 3.7.0) --
175 * -- LAPACK is a software package provided by Univ. of Tennessee, --
176 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
177 * December 2016
178 *
179 * .. Scalar Arguments ..
180  CHARACTER uplo
181  INTEGER kd, ldab, ldb, ldx, ldxact, n, nrhs
182 * ..
183 * .. Array Arguments ..
184  REAL ab( ldab, * ), b( ldb, * ), berr( * ),
185  $ ferr( * ), reslts( * ), x( ldx, * ),
186  $ xact( ldxact, * )
187 * ..
188 *
189 * =====================================================================
190 *
191 * .. Parameters ..
192  REAL zero, one
193  parameter( zero = 0.0e+0, one = 1.0e+0 )
194 * ..
195 * .. Local Scalars ..
196  LOGICAL upper
197  INTEGER i, imax, j, k, nz
198  REAL axbi, diff, eps, errbnd, ovfl, tmp, unfl, xnorm
199 * ..
200 * .. External Functions ..
201  LOGICAL lsame
202  INTEGER isamax
203  REAL slamch
204  EXTERNAL lsame, isamax, slamch
205 * ..
206 * .. Intrinsic Functions ..
207  INTRINSIC abs, max, min
208 * ..
209 * .. Executable Statements ..
210 *
211 * Quick exit if N = 0 or NRHS = 0.
212 *
213  IF( n.LE.0 .OR. nrhs.LE.0 ) THEN
214  reslts( 1 ) = zero
215  reslts( 2 ) = zero
216  RETURN
217  END IF
218 *
219  eps = slamch( 'Epsilon' )
220  unfl = slamch( 'Safe minimum' )
221  ovfl = one / unfl
222  upper = lsame( uplo, 'U' )
223  nz = 2*max( kd, n-1 ) + 1
224 *
225 * Test 1: Compute the maximum of
226 * norm(X - XACT) / ( norm(X) * FERR )
227 * over all the vectors X and XACT using the infinity-norm.
228 *
229  errbnd = zero
230  DO 30 j = 1, nrhs
231  imax = isamax( n, x( 1, j ), 1 )
232  xnorm = max( abs( x( imax, j ) ), unfl )
233  diff = zero
234  DO 10 i = 1, n
235  diff = max( diff, abs( x( i, j )-xact( i, j ) ) )
236  10 CONTINUE
237 *
238  IF( xnorm.GT.one ) THEN
239  GO TO 20
240  ELSE IF( diff.LE.ovfl*xnorm ) THEN
241  GO TO 20
242  ELSE
243  errbnd = one / eps
244  GO TO 30
245  END IF
246 *
247  20 CONTINUE
248  IF( diff / xnorm.LE.ferr( j ) ) THEN
249  errbnd = max( errbnd, ( diff / xnorm ) / ferr( j ) )
250  ELSE
251  errbnd = one / eps
252  END IF
253  30 CONTINUE
254  reslts( 1 ) = errbnd
255 *
256 * Test 2: Compute the maximum of BERR / ( NZ*EPS + (*) ), where
257 * (*) = NZ*UNFL / (min_i (abs(A)*abs(X) +abs(b))_i )
258 *
259  DO 90 k = 1, nrhs
260  DO 80 i = 1, n
261  tmp = abs( b( i, k ) )
262  IF( upper ) THEN
263  DO 40 j = max( i-kd, 1 ), i
264  tmp = tmp + abs( ab( kd+1-i+j, i ) )*abs( x( j, k ) )
265  40 CONTINUE
266  DO 50 j = i + 1, min( i+kd, n )
267  tmp = tmp + abs( ab( kd+1+i-j, j ) )*abs( x( j, k ) )
268  50 CONTINUE
269  ELSE
270  DO 60 j = max( i-kd, 1 ), i - 1
271  tmp = tmp + abs( ab( 1+i-j, j ) )*abs( x( j, k ) )
272  60 CONTINUE
273  DO 70 j = i, min( i+kd, n )
274  tmp = tmp + abs( ab( 1+j-i, i ) )*abs( x( j, k ) )
275  70 CONTINUE
276  END IF
277  IF( i.EQ.1 ) THEN
278  axbi = tmp
279  ELSE
280  axbi = min( axbi, tmp )
281  END IF
282  80 CONTINUE
283  tmp = berr( k ) / ( nz*eps+nz*unfl / max( axbi, nz*unfl ) )
284  IF( k.EQ.1 ) THEN
285  reslts( 2 ) = tmp
286  ELSE
287  reslts( 2 ) = max( reslts( 2 ), tmp )
288  END IF
289  90 CONTINUE
290 *
291  RETURN
292 *
293 * End of SPBT05
294 *
integer function isamax(N, SX, INCX)
ISAMAX
Definition: isamax.f:73
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
Here is the caller graph for this function: