d9/de4/cla__lin__berr_8f_source.html

00001       SUBROUTINE CLA_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       REAL               AYB( N, NRHS ), BERR( NRHS )
00018       COMPLEX            RES( N, NRHS )
00019 *     ..
00020 *
00021 *  Purpose
00022 *  =======
00023 *
00024 *     CLA_LIN_BERR computes componentwise 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 componentwise absolute value of the matrix
00028 *     or vector Z.
00029 *
00030 *     N       (input) INTEGER
00031 *     The number of linear equations, i.e., the order of the
00032 *     matrix A.  N >= 0.
00033 *
00034 *     NZ      (input) INTEGER
00035 *     We add (NZ+1)*SLAMCH( 'Safe minimum' ) to R(i) in the numerator to
00036 *     guard against spuriously zero residuals. Default value is N.
00037 *
00038 *     NRHS    (input) INTEGER
00039 *     The number of right hand sides, i.e., the number of columns
00040 *     of the matrices AYB, RES, and BERR.  NRHS >= 0.
00041 *
00042 *     RES    (input) DOUBLE PRECISION array, dimension (N,NRHS)
00043 *     The residual matrix, i.e., the matrix R in the relative backward
00044 *     error formula above.
00045 *
00046 *     AYB    (input) DOUBLE PRECISION array, dimension (N, NRHS)
00047 *     The denominator in the relative backward error formula above, i.e.,
00048 *     the matrix abs(op(A_s))*abs(Y) + abs(B_s). The matrices A, Y, and B
00049 *     are from iterative refinement (see cla_gerfsx_extended.f).
00050 *
00051 *     BERR   (output) COMPLEX array, dimension (NRHS)
00052 *     The componentwise relative backward error from the formula above.
00053 *
00054 *  =====================================================================
00055 *
00056 *     .. Local Scalars ..
00057       REAL               TMP
00058       INTEGER            I, J
00059       COMPLEX            CDUM
00060 *     ..
00061 *     .. Intrinsic Functions ..
00062       INTRINSIC          ABS, REAL, AIMAG, MAX
00063 *     ..
00064 *     .. External Functions ..
00065       EXTERNAL           SLAMCH
00066       REAL               SLAMCH
00067       REAL               SAFE1
00068 *     ..
00069 *     .. Statement Functions ..
00070       COMPLEX            CABS1
00071 *     ..
00072 *     .. Statement Function Definitions ..
00073       CABS1( CDUM ) = ABS( REAL( CDUM ) ) + ABS( AIMAG( CDUM ) )
00074 *     ..
00075 *     .. Executable Statements ..
00076 *
00077 *     Adding SAFE1 to the numerator guards against spuriously zero
00078 *     residuals.  A similar safeguard is in the CLA_yyAMV routine used
00079 *     to compute AYB.
00080 *
00081       SAFE1 = SLAMCH( 'Safe minimum' )
00082       SAFE1 = (NZ+1)*SAFE1
00083
00084       DO J = 1, NRHS
00085          BERR(J) = 0.0
00086          DO I = 1, N
00087             IF (AYB(I,J) .NE. 0.0) THEN
00088                TMP = (SAFE1 + CABS1(RES(I,J)))/AYB(I,J)
00089                BERR(J) = MAX( BERR(J), TMP )
00090             END IF
00091 *
00092 *     If AYB is exactly 0.0 (and if computed by CLA_yyAMV), then we know
00093 *     the true residual also must be exactly 0.0.
00094 *
00095          END DO
00096       END DO
00097       END