LAPACK  3.8.0
LAPACK: Linear Algebra PACKage

◆ cptt05()

subroutine cptt05 ( integer  N,
integer  NRHS,
real, dimension( * )  D,
complex, dimension( * )  E,
complex, dimension( ldb, * )  B,
integer  LDB,
complex, dimension( ldx, * )  X,
integer  LDX,
complex, dimension( ldxact, * )  XACT,
integer  LDXACT,
real, dimension( * )  FERR,
real, dimension( * )  BERR,
real, dimension( * )  RESLTS 
)

CPTT05

Purpose:
 CPTT05 tests the error bounds from iterative refinement for the
 computed solution to a system of equations A*X = B, where A is a
 Hermitian tridiagonal matrix of order n.

 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]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]NRHS
          NRHS is INTEGER
          The number of columns of the matrices X, B, and XACT.
          NRHS >= 0.
[in]D
          D is REAL array, dimension (N)
          The n diagonal elements of the tridiagonal matrix A.
[in]E
          E is COMPLEX array, dimension (N-1)
          The (n-1) subdiagonal elements of the tridiagonal matrix A.
[in]B
          B is COMPLEX 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 COMPLEX 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 COMPLEX 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 152 of file cptt05.f.

152 *
153 * -- LAPACK test routine (version 3.7.0) --
154 * -- LAPACK is a software package provided by Univ. of Tennessee, --
155 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
156 * December 2016
157 *
158 * .. Scalar Arguments ..
159  INTEGER ldb, ldx, ldxact, n, nrhs
160 * ..
161 * .. Array Arguments ..
162  REAL berr( * ), d( * ), ferr( * ), reslts( * )
163  COMPLEX b( ldb, * ), e( * ), x( ldx, * ),
164  $ xact( ldxact, * )
165 * ..
166 *
167 * =====================================================================
168 *
169 * .. Parameters ..
170  REAL zero, one
171  parameter( zero = 0.0e+0, one = 1.0e+0 )
172 * ..
173 * .. Local Scalars ..
174  INTEGER i, imax, j, k, nz
175  REAL axbi, diff, eps, errbnd, ovfl, tmp, unfl, xnorm
176  COMPLEX zdum
177 * ..
178 * .. External Functions ..
179  INTEGER icamax
180  REAL slamch
181  EXTERNAL icamax, slamch
182 * ..
183 * .. Intrinsic Functions ..
184  INTRINSIC abs, aimag, max, min, real
185 * ..
186 * .. Statement Functions ..
187  REAL cabs1
188 * ..
189 * .. Statement Function definitions ..
190  cabs1( zdum ) = abs( REAL( ZDUM ) ) + abs( aimag( zdum ) )
191 * ..
192 * .. Executable Statements ..
193 *
194 * Quick exit if N = 0 or NRHS = 0.
195 *
196  IF( n.LE.0 .OR. nrhs.LE.0 ) THEN
197  reslts( 1 ) = zero
198  reslts( 2 ) = zero
199  RETURN
200  END IF
201 *
202  eps = slamch( 'Epsilon' )
203  unfl = slamch( 'Safe minimum' )
204  ovfl = one / unfl
205  nz = 4
206 *
207 * Test 1: Compute the maximum of
208 * norm(X - XACT) / ( norm(X) * FERR )
209 * over all the vectors X and XACT using the infinity-norm.
210 *
211  errbnd = zero
212  DO 30 j = 1, nrhs
213  imax = icamax( n, x( 1, j ), 1 )
214  xnorm = max( cabs1( x( imax, j ) ), unfl )
215  diff = zero
216  DO 10 i = 1, n
217  diff = max( diff, cabs1( x( i, j )-xact( i, j ) ) )
218  10 CONTINUE
219 *
220  IF( xnorm.GT.one ) THEN
221  GO TO 20
222  ELSE IF( diff.LE.ovfl*xnorm ) THEN
223  GO TO 20
224  ELSE
225  errbnd = one / eps
226  GO TO 30
227  END IF
228 *
229  20 CONTINUE
230  IF( diff / xnorm.LE.ferr( j ) ) THEN
231  errbnd = max( errbnd, ( diff / xnorm ) / ferr( j ) )
232  ELSE
233  errbnd = one / eps
234  END IF
235  30 CONTINUE
236  reslts( 1 ) = errbnd
237 *
238 * Test 2: Compute the maximum of BERR / ( NZ*EPS + (*) ), where
239 * (*) = NZ*UNFL / (min_i (abs(A)*abs(X) +abs(b))_i )
240 *
241  DO 50 k = 1, nrhs
242  IF( n.EQ.1 ) THEN
243  axbi = cabs1( b( 1, k ) ) + cabs1( d( 1 )*x( 1, k ) )
244  ELSE
245  axbi = cabs1( b( 1, k ) ) + cabs1( d( 1 )*x( 1, k ) ) +
246  $ cabs1( e( 1 ) )*cabs1( x( 2, k ) )
247  DO 40 i = 2, n - 1
248  tmp = cabs1( b( i, k ) ) + cabs1( e( i-1 ) )*
249  $ cabs1( x( i-1, k ) ) + cabs1( d( i )*x( i, k ) ) +
250  $ cabs1( e( i ) )*cabs1( x( i+1, k ) )
251  axbi = min( axbi, tmp )
252  40 CONTINUE
253  tmp = cabs1( b( n, k ) ) + cabs1( e( n-1 ) )*
254  $ cabs1( x( n-1, k ) ) + cabs1( d( n )*x( n, k ) )
255  axbi = min( axbi, tmp )
256  END IF
257  tmp = berr( k ) / ( nz*eps+nz*unfl / max( axbi, nz*unfl ) )
258  IF( k.EQ.1 ) THEN
259  reslts( 2 ) = tmp
260  ELSE
261  reslts( 2 ) = max( reslts( 2 ), tmp )
262  END IF
263  50 CONTINUE
264 *
265  RETURN
266 *
267 * End of CPTT05
268 *
integer function icamax(N, CX, INCX)
ICAMAX
Definition: icamax.f:73
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
Here is the caller graph for this function: