LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE DLA_LIN_BERR ( N, NZ, NRHS, RES, AYB, BERR ) 00002 * 00003 * -- LAPACK routine (version 3.2.2) -- 00004 * -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and -- 00005 * -- Jason Riedy of Univ. of California Berkeley. -- 00006 * -- June 2010 -- 00007 * 00008 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00009 * -- Univ. of California Berkeley and NAG Ltd. -- 00010 * 00011 IMPLICIT NONE 00012 * .. 00013 * .. Scalar Arguments .. 00014 INTEGER N, NZ, NRHS 00015 * .. 00016 * .. Array Arguments .. 00017 DOUBLE PRECISION AYB( N, NRHS ), BERR( NRHS ) 00018 DOUBLE PRECISION RES( N, NRHS ) 00019 * .. 00020 * 00021 * Purpose 00022 * ======= 00023 * 00024 * DLA_LIN_BERR computes component-wise relative backward error from 00025 * the formula 00026 * max(i) ( abs(R(i)) / ( abs(op(A_s))*abs(Y) + abs(B_s) )(i) ) 00027 * where abs(Z) is the component-wise absolute value of the matrix 00028 * or vector Z. 00029 * 00030 * Arguments 00031 * ========== 00032 * 00033 * N (input) INTEGER 00034 * The number of linear equations, i.e., the order of the 00035 * matrix A. N >= 0. 00036 * 00037 * NZ (input) INTEGER 00038 * We add (NZ+1)*SLAMCH( 'Safe minimum' ) to R(i) in the numerator to 00039 * guard against spuriously zero residuals. Default value is N. 00040 * 00041 * NRHS (input) INTEGER 00042 * The number of right hand sides, i.e., the number of columns 00043 * of the matrices AYB, RES, and BERR. NRHS >= 0. 00044 * 00045 * RES (input) DOUBLE PRECISION array, dimension (N,NRHS) 00046 * The residual matrix, i.e., the matrix R in the relative backward 00047 * error formula above. 00048 * 00049 * AYB (input) DOUBLE PRECISION array, dimension (N, NRHS) 00050 * The denominator in the relative backward error formula above, i.e., 00051 * the matrix abs(op(A_s))*abs(Y) + abs(B_s). The matrices A, Y, and B 00052 * are from iterative refinement (see dla_gerfsx_extended.f). 00053 * 00054 * BERR (output) DOUBLE PRECISION array, dimension (NRHS) 00055 * The component-wise relative backward error from the formula above. 00056 * 00057 * ===================================================================== 00058 * 00059 * .. Local Scalars .. 00060 DOUBLE PRECISION TMP 00061 INTEGER I, J 00062 * .. 00063 * .. Intrinsic Functions .. 00064 INTRINSIC ABS, MAX 00065 * .. 00066 * .. External Functions .. 00067 EXTERNAL DLAMCH 00068 DOUBLE PRECISION DLAMCH 00069 DOUBLE PRECISION SAFE1 00070 * .. 00071 * .. Executable Statements .. 00072 * 00073 * Adding SAFE1 to the numerator guards against spuriously zero 00074 * residuals. A similar safeguard is in the SLA_yyAMV routine used 00075 * to compute AYB. 00076 * 00077 SAFE1 = DLAMCH( 'Safe minimum' ) 00078 SAFE1 = (NZ+1)*SAFE1 00079 00080 DO J = 1, NRHS 00081 BERR(J) = 0.0D+0 00082 DO I = 1, N 00083 IF (AYB(I,J) .NE. 0.0D+0) THEN 00084 TMP = (SAFE1+ABS(RES(I,J)))/AYB(I,J) 00085 BERR(J) = MAX( BERR(J), TMP ) 00086 END IF 00087 * 00088 * If AYB is exactly 0.0 (and if computed by SLA_yyAMV), then we know 00089 * the true residual also must be exactly 0.0. 00090 * 00091 END DO 00092 END DO 00093 END