TESTING/LIN/cqrt17.f

   1 *> \brief \b CQRT17
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *  Definition:
   9 *  ===========
  10 *
  11 *       REAL             FUNCTION CQRT17( TRANS, IRESID, M, N, NRHS, A,
  12 *                        LDA, X, LDX, B, LDB, C, WORK, LWORK )
  13 *
  14 *       .. Scalar Arguments ..
  15 *       CHARACTER          TRANS
  16 *       INTEGER            IRESID, LDA, LDB, LDX, LWORK, M, N, NRHS
  17 *       ..
  18 *       .. Array Arguments ..
  19 *       COMPLEX            A( LDA, * ), B( LDB, * ), C( LDB, * ),
  20 *      $                   WORK( LWORK ), X( LDX, * )
  21 *       ..
  22 *
  23 *
  24 *> \par Purpose:
  25 *  =============
  26 *>
  27 *> \verbatim
  28 *>
  29 *> CQRT17 computes the ratio
  30 *>
  31 *>    || R'*op(A) ||/(||A||*alpha*max(M,N,NRHS)*eps)
  32 *>
  33 *> where R = op(A)*X - B, op(A) is A or A', and
  34 *>
  35 *>    alpha = ||B|| if IRESID = 1 (zero-residual problem)
  36 *>    alpha = ||R|| if IRESID = 2 (otherwise).
  37 *> \endverbatim
  38 *
  39 *  Arguments:
  40 *  ==========
  41 *
  42 *> \param[in] TRANS
  43 *> \verbatim
  44 *>          TRANS is CHARACTER*1
  45 *>          Specifies whether or not the transpose of A is used.
  46 *>          = 'N':  No transpose, op(A) = A.
  47 *>          = 'C':  Conjugate transpose, op(A) = A'.
  48 *> \endverbatim
  49 *>
  50 *> \param[in] IRESID
  51 *> \verbatim
  52 *>          IRESID is INTEGER
  53 *>          IRESID = 1 indicates zero-residual problem.
  54 *>          IRESID = 2 indicates non-zero residual.
  55 *> \endverbatim
  56 *>
  57 *> \param[in] M
  58 *> \verbatim
  59 *>          M is INTEGER
  60 *>          The number of rows of the matrix A.
  61 *>          If TRANS = 'N', the number of rows of the matrix B.
  62 *>          If TRANS = 'C', the number of rows of the matrix X.
  63 *> \endverbatim
  64 *>
  65 *> \param[in] N
  66 *> \verbatim
  67 *>          N is INTEGER
  68 *>          The number of columns of the matrix  A.
  69 *>          If TRANS = 'N', the number of rows of the matrix X.
  70 *>          If TRANS = 'C', the number of rows of the matrix B.
  71 *> \endverbatim
  72 *>
  73 *> \param[in] NRHS
  74 *> \verbatim
  75 *>          NRHS is INTEGER
  76 *>          The number of columns of the matrices X and B.
  77 *> \endverbatim
  78 *>
  79 *> \param[in] A
  80 *> \verbatim
  81 *>          A is COMPLEX array, dimension (LDA,N)
  82 *>          The m-by-n matrix A.
  83 *> \endverbatim
  84 *>
  85 *> \param[in] LDA
  86 *> \verbatim
  87 *>          LDA is INTEGER
  88 *>          The leading dimension of the array A. LDA >= M.
  89 *> \endverbatim
  90 *>
  91 *> \param[in] X
  92 *> \verbatim
  93 *>          X is COMPLEX array, dimension (LDX,NRHS)
  94 *>          If TRANS = 'N', the n-by-nrhs matrix X.
  95 *>          If TRANS = 'C', the m-by-nrhs matrix X.
  96 *> \endverbatim
  97 *>
  98 *> \param[in] LDX
  99 *> \verbatim
 100 *>          LDX is INTEGER
 101 *>          The leading dimension of the array X.
 102 *>          If TRANS = 'N', LDX >= N.
 103 *>          If TRANS = 'C', LDX >= M.
 104 *> \endverbatim
 105 *>
 106 *> \param[in] B
 107 *> \verbatim
 108 *>          B is COMPLEX array, dimension (LDB,NRHS)
 109 *>          If TRANS = 'N', the m-by-nrhs matrix B.
 110 *>          If TRANS = 'C', the n-by-nrhs matrix B.
 111 *> \endverbatim
 112 *>
 113 *> \param[in] LDB
 114 *> \verbatim
 115 *>          LDB is INTEGER
 116 *>          The leading dimension of the array B.
 117 *>          If TRANS = 'N', LDB >= M.
 118 *>          If TRANS = 'C', LDB >= N.
 119 *> \endverbatim
 120 *>
 121 *> \param[out] C
 122 *> \verbatim
 123 *>          C is COMPLEX array, dimension (LDB,NRHS)
 124 *> \endverbatim
 125 *>
 126 *> \param[out] WORK
 127 *> \verbatim
 128 *>          WORK is COMPLEX array, dimension (LWORK)
 129 *> \endverbatim
 130 *>
 131 *> \param[in] LWORK
 132 *> \verbatim
 133 *>          LWORK is INTEGER
 134 *>          The length of the array WORK.  LWORK >= NRHS*(M+N).
 135 *> \endverbatim
 136 *
 137 *  Authors:
 138 *  ========
 139 *
 140 *> \author Univ. of Tennessee
 141 *> \author Univ. of California Berkeley
 142 *> \author Univ. of Colorado Denver
 143 *> \author NAG Ltd.
 144 *
 145 *> \date November 2015
 146 *
 147 *> \ingroup complex_lin
 148 *
 149 *  =====================================================================
 150       REAL             FUNCTION CQRT17( TRANS, IRESID, M, N, NRHS, A,
 151      $                 LDA, X, LDX, B, LDB, C, WORK, LWORK )
 152 *
 153 *  -- LAPACK test routine (version 3.6.0) --
 154 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 155 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 156 *     November 2015
 157 *
 158 *     .. Scalar Arguments ..
 159       CHARACTER          TRANS
 160       INTEGER            IRESID, LDA, LDB, LDX, LWORK, M, N, NRHS
 161 *     ..
 162 *     .. Array Arguments ..
 163       COMPLEX            A( LDA, * ), B( LDB, * ), C( LDB, * ),
 164      $                   WORK( LWORK ), X( LDX, * )
 165 *     ..
 166 *
 167 *  =====================================================================
 168 *
 169 *     .. Parameters ..
 170       REAL               ZERO, ONE
 171       PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0 )
 172 *     ..
 173 *     .. Local Scalars ..
 174       INTEGER            INFO, ISCL, NCOLS, NROWS
 175       REAL               BIGNUM, ERR, NORMA, NORMB, NORMRS, NORMX,
 176      $                   SMLNUM
 177 *     ..
 178 *     .. Local Arrays ..
 179       REAL               RWORK( 1 )
 180 *     ..
 181 *     .. External Functions ..
 182       LOGICAL            LSAME
 183       REAL               CLANGE, SLAMCH
 184       EXTERNAL           LSAME, CLANGE, SLAMCH
 185 *     ..
 186 *     .. External Subroutines ..
 187       EXTERNAL           CGEMM, CLACPY, CLASCL, XERBLA
 188 *     ..
 189 *     .. Intrinsic Functions ..
 190       INTRINSIC          CMPLX, MAX, REAL
 191 *     ..
 192 *     .. Executable Statements ..
 193 *
 194       CQRT17 = ZERO
 195 *
 196       IF( LSAME( TRANS, 'N' ) ) THEN
 197          NROWS = M
 198          NCOLS = N
 199       ELSE IF( LSAME( TRANS, 'C' ) ) THEN
 200          NROWS = N
 201          NCOLS = M
 202       ELSE
 203          CALL XERBLA( 'CQRT17', 1 )
 204          RETURN
 205       END IF
 206 *
 207       IF( LWORK.LT.NCOLS*NRHS ) THEN
 208          CALL XERBLA( 'CQRT17', 13 )
 209          RETURN
 210       END IF
 211 *
 212       IF( M.LE.0 .OR. N.LE.0 .OR. NRHS.LE.0 )
 213      $   RETURN
 214 *
 215       NORMA = CLANGE( 'One-norm', M, N, A, LDA, RWORK )
 216       SMLNUM = SLAMCH( 'Safe minimum' ) / SLAMCH( 'Precision' )
 217       BIGNUM = ONE / SMLNUM
 218       ISCL = 0
 219 *
 220 *     compute residual and scale it
 221 *
 222       CALL CLACPY( 'All', NROWS, NRHS, B, LDB, C, LDB )
 223       CALL CGEMM( TRANS, 'No transpose', NROWS, NRHS, NCOLS,
 224      $            CMPLX( -ONE ), A, LDA, X, LDX, CMPLX( ONE ), C, LDB )
 225       NORMRS = CLANGE( 'Max', NROWS, NRHS, C, LDB, RWORK )
 226       IF( NORMRS.GT.SMLNUM ) THEN
 227          ISCL = 1
 228          CALL CLASCL( 'General', 0, 0, NORMRS, ONE, NROWS, NRHS, C, LDB,
 229      $                INFO )
 230       END IF
 231 *
 232 *     compute R'*A
 233 *
 234       CALL CGEMM( 'Conjugate transpose', TRANS, NRHS, NCOLS, NROWS,
 235      $            CMPLX( ONE ), C, LDB, A, LDA, CMPLX( ZERO ), WORK,
 236      $            NRHS )
 237 *
 238 *     compute and properly scale error
 239 *
 240       ERR = CLANGE( 'One-norm', NRHS, NCOLS, WORK, NRHS, RWORK )
 241       IF( NORMA.NE.ZERO )
 242      $   ERR = ERR / NORMA
 243 *
 244       IF( ISCL.EQ.1 )
 245      $   ERR = ERR*NORMRS
 246 *
 247       IF( IRESID.EQ.1 ) THEN
 248          NORMB = CLANGE( 'One-norm', NROWS, NRHS, B, LDB, RWORK )
 249          IF( NORMB.NE.ZERO )
 250      $      ERR = ERR / NORMB
 251       ELSE
 252          IF( NORMRS.NE.ZERO )
 253      $      ERR = ERR / NORMRS
 254       END IF
 255 *
 256       CQRT17 = ERR / ( SLAMCH( 'Epsilon' )*REAL( MAX( M, N, NRHS ) ) )
 257       RETURN
 258 *
 259 *     End of CQRT17
 260 *
 261       END