LAPACK  3.10.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.

Definition at line 148 of file cptt05.f.

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