SRC/zla_lin_berr.f

   1 *> \brief \b ZLA_LIN_BERR computes a component-wise relative backward error.
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download ZLA_LIN_BERR + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zla_lin_berr.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zla_lin_berr.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zla_lin_berr.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE ZLA_LIN_BERR ( N, NZ, NRHS, RES, AYB, BERR )
  22 *
  23 *       .. Scalar Arguments ..
  24 *       INTEGER            N, NZ, NRHS
  25 *       ..
  26 *       .. Array Arguments ..
  27 *       DOUBLE PRECISION   AYB( N, NRHS ), BERR( NRHS )
  28 *       COMPLEX*16         RES( N, NRHS )
  29 *       ..
  30 *
  31 *
  32 *> \par Purpose:
  33 *  =============
  34 *>
  35 *> \verbatim
  36 *>
  37 *>    ZLA_LIN_BERR computes componentwise relative backward error from
  38 *>    the formula
  39 *>        max(i) ( abs(R(i)) / ( abs(op(A_s))*abs(Y) + abs(B_s) )(i) )
  40 *>    where abs(Z) is the componentwise absolute value of the matrix
  41 *>    or vector Z.
  42 *> \endverbatim
  43 *
  44 *  Arguments:
  45 *  ==========
  46 *
  47 *> \param[in] N
  48 *> \verbatim
  49 *>          N is INTEGER
  50 *>     The number of linear equations, i.e., the order of the
  51 *>     matrix A.  N >= 0.
  52 *> \endverbatim
  53 *>
  54 *> \param[in] NZ
  55 *> \verbatim
  56 *>          NZ is INTEGER
  57 *>     We add (NZ+1)*SLAMCH( 'Safe minimum' ) to R(i) in the numerator to
  58 *>     guard against spuriously zero residuals. Default value is N.
  59 *> \endverbatim
  60 *>
  61 *> \param[in] NRHS
  62 *> \verbatim
  63 *>          NRHS is INTEGER
  64 *>     The number of right hand sides, i.e., the number of columns
  65 *>     of the matrices AYB, RES, and BERR.  NRHS >= 0.
  66 *> \endverbatim
  67 *>
  68 *> \param[in] RES
  69 *> \verbatim
  70 *>          RES is COMPLEX*16 array, dimension (N,NRHS)
  71 *>     The residual matrix, i.e., the matrix R in the relative backward
  72 *>     error formula above.
  73 *> \endverbatim
  74 *>
  75 *> \param[in] AYB
  76 *> \verbatim
  77 *>          AYB is DOUBLE PRECISION array, dimension (N, NRHS)
  78 *>     The denominator in the relative backward error formula above, i.e.,
  79 *>     the matrix abs(op(A_s))*abs(Y) + abs(B_s). The matrices A, Y, and B
  80 *>     are from iterative refinement (see zla_gerfsx_extended.f).
  81 *> \endverbatim
  82 *>
  83 *> \param[out] BERR
  84 *> \verbatim
  85 *>          BERR is DOUBLE PRECISION array, dimension (NRHS)
  86 *>     The componentwise relative backward error from the formula above.
  87 *> \endverbatim
  88 *
  89 *  Authors:
  90 *  ========
  91 *
  92 *> \author Univ. of Tennessee
  93 *> \author Univ. of California Berkeley
  94 *> \author Univ. of Colorado Denver
  95 *> \author NAG Ltd.
  96 *
  97 *> \date June 2016
  98 *
  99 *> \ingroup complex16OTHERcomputational
 100 *
 101 *  =====================================================================
 102       SUBROUTINE ZLA_LIN_BERR ( N, NZ, NRHS, RES, AYB, BERR )
 103 *
 104 *  -- LAPACK computational routine (version 3.6.1) --
 105 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 106 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 107 *     June 2016
 108 *
 109 *     .. Scalar Arguments ..
 110       INTEGER            N, NZ, NRHS
 111 *     ..
 112 *     .. Array Arguments ..
 113       DOUBLE PRECISION   AYB( N, NRHS ), BERR( NRHS )
 114       COMPLEX*16         RES( N, NRHS )
 115 *     ..
 116 *
 117 *  =====================================================================
 118 *
 119 *     .. Local Scalars ..
 120       DOUBLE PRECISION   TMP
 121       INTEGER            I, J
 122       COMPLEX*16         CDUM
 123 *     ..
 124 *     .. Intrinsic Functions ..
 125       INTRINSIC          ABS, REAL, DIMAG, MAX
 126 *     ..
 127 *     .. External Functions ..
 128       EXTERNAL           DLAMCH
 129       DOUBLE PRECISION   DLAMCH
 130       DOUBLE PRECISION   SAFE1
 131 *     ..
 132 *     .. Statement Functions ..
 133       COMPLEX*16         CABS1
 134 *     ..
 135 *     .. Statement Function Definitions ..
 136       CABS1( CDUM ) = ABS( DBLE( CDUM ) ) + ABS( DIMAG( CDUM ) )
 137 *     ..
 138 *     .. Executable Statements ..
 139 *
 140 *     Adding SAFE1 to the numerator guards against spuriously zero
 141 *     residuals.  A similar safeguard is in the CLA_yyAMV routine used
 142 *     to compute AYB.
 143 *
 144       SAFE1 = DLAMCH( 'Safe minimum' )
 145       SAFE1 = (NZ+1)*SAFE1
 146
 147       DO J = 1, NRHS
 148          BERR(J) = 0.0D+0
 149          DO I = 1, N
 150             IF (AYB(I,J) .NE. 0.0D+0) THEN
 151                TMP = (SAFE1 + CABS1(RES(I,J)))/AYB(I,J)
 152                BERR(J) = MAX( BERR(J), TMP )
 153             END IF
 154 *
 155 *     If AYB is exactly 0.0 (and if computed by CLA_yyAMV), then we know
 156 *     the true residual also must be exactly 0.0.
 157 *
 158          END DO
 159       END DO
 160       END