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

◆ ctrt05()

subroutine ctrt05 ( character  uplo,
character  trans,
character  diag,
integer  n,
integer  nrhs,
complex, dimension( lda, * )  a,
integer  lda,
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 
)

CTRT05

Purpose:
 CTRT05 tests the error bounds from iterative refinement for the
 computed solution to a system of equations A*X = B, where A is a
 triangular n by n 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 / ( (n+1)*EPS + (*) ), where
             (*) = (n+1)*UNFL / (min_i (abs(A)*abs(X) +abs(b))_i )
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          Specifies whether the matrix A is upper or lower triangular.
          = 'U':  Upper triangular
          = 'L':  Lower triangular
[in]TRANS
          TRANS is CHARACTER*1
          Specifies the form of the system of equations.
          = 'N':  A * X = B  (No transpose)
          = 'T':  A'* X = B  (Transpose)
          = 'C':  A'* X = B  (Conjugate transpose = Transpose)
[in]DIAG
          DIAG is CHARACTER*1
          Specifies whether or not the matrix A is unit triangular.
          = 'N':  Non-unit triangular
          = 'U':  Unit 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]NRHS
          NRHS is INTEGER
          The number of columns of the matrices X, B, and XACT.
          NRHS >= 0.
[in]A
          A is COMPLEX array, dimension (LDA,N)
          The triangular matrix A.  If UPLO = 'U', the leading n by n
          upper triangular part of the array A contains the upper
          triangular matrix, and the strictly lower triangular part of
          A is not referenced.  If UPLO = 'L', the leading n by n lower
          triangular part of the array A contains the lower triangular
          matrix, and the strictly upper triangular part of A is not
          referenced.  If DIAG = 'U', the diagonal elements of A are
          also not referenced and are assumed to be 1.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[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 / ( (n+1)*EPS + (*) )
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 180 of file ctrt05.f.

182*
183* -- LAPACK test routine --
184* -- LAPACK is a software package provided by Univ. of Tennessee, --
185* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
186*
187* .. Scalar Arguments ..
188 CHARACTER DIAG, TRANS, UPLO
189 INTEGER LDA, LDB, LDX, LDXACT, N, NRHS
190* ..
191* .. Array Arguments ..
192 REAL BERR( * ), FERR( * ), RESLTS( * )
193 COMPLEX A( LDA, * ), B( LDB, * ), X( LDX, * ),
194 $ XACT( LDXACT, * )
195* ..
196*
197* =====================================================================
198*
199* .. Parameters ..
200 REAL ZERO, ONE
201 parameter( zero = 0.0e+0, one = 1.0e+0 )
202* ..
203* .. Local Scalars ..
204 LOGICAL NOTRAN, UNIT, UPPER
205 INTEGER I, IFU, IMAX, J, K
206 REAL AXBI, DIFF, EPS, ERRBND, OVFL, TMP, UNFL, XNORM
207 COMPLEX ZDUM
208* ..
209* .. External Functions ..
210 LOGICAL LSAME
211 INTEGER ICAMAX
212 REAL SLAMCH
213 EXTERNAL lsame, icamax, slamch
214* ..
215* .. Intrinsic Functions ..
216 INTRINSIC abs, aimag, max, min, real
217* ..
218* .. Statement Functions ..
219 REAL CABS1
220* ..
221* .. Statement Function definitions ..
222 cabs1( zdum ) = abs( real( zdum ) ) + abs( aimag( zdum ) )
223* ..
224* .. Executable Statements ..
225*
226* Quick exit if N = 0 or NRHS = 0.
227*
228 IF( n.LE.0 .OR. nrhs.LE.0 ) THEN
229 reslts( 1 ) = zero
230 reslts( 2 ) = zero
231 RETURN
232 END IF
233*
234 eps = slamch( 'Epsilon' )
235 unfl = slamch( 'Safe minimum' )
236 ovfl = one / unfl
237 upper = lsame( uplo, 'U' )
238 notran = lsame( trans, 'N' )
239 unit = lsame( diag, 'U' )
240*
241* Test 1: Compute the maximum of
242* norm(X - XACT) / ( norm(X) * FERR )
243* over all the vectors X and XACT using the infinity-norm.
244*
245 errbnd = zero
246 DO 30 j = 1, nrhs
247 imax = icamax( n, x( 1, j ), 1 )
248 xnorm = max( cabs1( x( imax, j ) ), unfl )
249 diff = zero
250 DO 10 i = 1, n
251 diff = max( diff, cabs1( x( i, j )-xact( i, j ) ) )
252 10 CONTINUE
253*
254 IF( xnorm.GT.one ) THEN
255 GO TO 20
256 ELSE IF( diff.LE.ovfl*xnorm ) THEN
257 GO TO 20
258 ELSE
259 errbnd = one / eps
260 GO TO 30
261 END IF
262*
263 20 CONTINUE
264 IF( diff / xnorm.LE.ferr( j ) ) THEN
265 errbnd = max( errbnd, ( diff / xnorm ) / ferr( j ) )
266 ELSE
267 errbnd = one / eps
268 END IF
269 30 CONTINUE
270 reslts( 1 ) = errbnd
271*
272* Test 2: Compute the maximum of BERR / ( (n+1)*EPS + (*) ), where
273* (*) = (n+1)*UNFL / (min_i (abs(A)*abs(X) +abs(b))_i )
274*
275 ifu = 0
276 IF( unit )
277 $ ifu = 1
278 DO 90 k = 1, nrhs
279 DO 80 i = 1, n
280 tmp = cabs1( b( i, k ) )
281 IF( upper ) THEN
282 IF( .NOT.notran ) THEN
283 DO 40 j = 1, i - ifu
284 tmp = tmp + cabs1( a( j, i ) )*cabs1( x( j, k ) )
285 40 CONTINUE
286 IF( unit )
287 $ tmp = tmp + cabs1( x( i, k ) )
288 ELSE
289 IF( unit )
290 $ tmp = tmp + cabs1( x( i, k ) )
291 DO 50 j = i + ifu, n
292 tmp = tmp + cabs1( a( i, j ) )*cabs1( x( j, k ) )
293 50 CONTINUE
294 END IF
295 ELSE
296 IF( notran ) THEN
297 DO 60 j = 1, i - ifu
298 tmp = tmp + cabs1( a( i, j ) )*cabs1( x( j, k ) )
299 60 CONTINUE
300 IF( unit )
301 $ tmp = tmp + cabs1( x( i, k ) )
302 ELSE
303 IF( unit )
304 $ tmp = tmp + cabs1( x( i, k ) )
305 DO 70 j = i + ifu, n
306 tmp = tmp + cabs1( a( j, i ) )*cabs1( x( j, k ) )
307 70 CONTINUE
308 END IF
309 END IF
310 IF( i.EQ.1 ) THEN
311 axbi = tmp
312 ELSE
313 axbi = min( axbi, tmp )
314 END IF
315 80 CONTINUE
316 tmp = berr( k ) / ( ( n+1 )*eps+( n+1 )*unfl /
317 $ max( axbi, ( n+1 )*unfl ) )
318 IF( k.EQ.1 ) THEN
319 reslts( 2 ) = tmp
320 ELSE
321 reslts( 2 ) = max( reslts( 2 ), tmp )
322 END IF
323 90 CONTINUE
324*
325 RETURN
326*
327* End of CTRT05
328*
icamax
integer function icamax(n, cx, incx)
ICAMAX
Definition icamax.f:71
slamch
real function slamch(cmach)
SLAMCH
Definition slamch.f:68
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
Here is the caller graph for this function: