TESTING/EIG/cget35.f

   1 *> \brief \b CGET35
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *  Definition:
   9 *  ===========
  10 *
  11 *       SUBROUTINE CGET35( RMAX, LMAX, NINFO, KNT, NIN )
  12 *
  13 *       .. Scalar Arguments ..
  14 *       INTEGER            KNT, LMAX, NIN, NINFO
  15 *       REAL               RMAX
  16 *       ..
  17 *
  18 *
  19 *> \par Purpose:
  20 *  =============
  21 *>
  22 *> \verbatim
  23 *>
  24 *> CGET35 tests CTRSYL, a routine for solving the Sylvester matrix
  25 *> equation
  26 *>
  27 *>    op(A)*X + ISGN*X*op(B) = scale*C,
  28 *>
  29 *> A and B are assumed to be in Schur canonical form, op() represents an
  30 *> optional transpose, and ISGN can be -1 or +1.  Scale is an output
  31 *> less than or equal to 1, chosen to avoid overflow in X.
  32 *>
  33 *> The test code verifies that the following residual is order 1:
  34 *>
  35 *>    norm(op(A)*X + ISGN*X*op(B) - scale*C) /
  36 *>        (EPS*max(norm(A),norm(B))*norm(X))
  37 *> \endverbatim
  38 *
  39 *  Arguments:
  40 *  ==========
  41 *
  42 *> \param[out] RMAX
  43 *> \verbatim
  44 *>          RMAX is REAL
  45 *>          Value of the largest test ratio.
  46 *> \endverbatim
  47 *>
  48 *> \param[out] LMAX
  49 *> \verbatim
  50 *>          LMAX is INTEGER
  51 *>          Example number where largest test ratio achieved.
  52 *> \endverbatim
  53 *>
  54 *> \param[out] NINFO
  55 *> \verbatim
  56 *>          NINFO is INTEGER
  57 *>          Number of examples where INFO is nonzero.
  58 *> \endverbatim
  59 *>
  60 *> \param[out] KNT
  61 *> \verbatim
  62 *>          KNT is INTEGER
  63 *>          Total number of examples tested.
  64 *> \endverbatim
  65 *>
  66 *> \param[in] NIN
  67 *> \verbatim
  68 *>          NIN is INTEGER
  69 *>          Input logical unit number.
  70 *> \endverbatim
  71 *
  72 *  Authors:
  73 *  ========
  74 *
  75 *> \author Univ. of Tennessee
  76 *> \author Univ. of California Berkeley
  77 *> \author Univ. of Colorado Denver
  78 *> \author NAG Ltd.
  79 *
  80 *> \date November 2011
  81 *
  82 *> \ingroup complex_eig
  83 *
  84 *  =====================================================================
  85       SUBROUTINE CGET35( RMAX, LMAX, NINFO, KNT, NIN )
  86 *
  87 *  -- LAPACK test routine (version 3.4.0) --
  88 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  89 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  90 *     November 2011
  91 *
  92 *     .. Scalar Arguments ..
  93       INTEGER            KNT, LMAX, NIN, NINFO
  94       REAL               RMAX
  95 *     ..
  96 *
  97 *  =====================================================================
  98 *
  99 *     .. Parameters ..
 100       INTEGER            LDT
 101       PARAMETER          ( LDT = 10 )
 102       REAL               ZERO, ONE, TWO
 103       PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0 )
 104       REAL               LARGE
 105       PARAMETER          ( LARGE = 1.0E6 )
 106       COMPLEX            CONE
 107       PARAMETER          ( CONE = 1.0E0 )
 108 *     ..
 109 *     .. Local Scalars ..
 110       CHARACTER          TRANA, TRANB
 111       INTEGER            I, IMLA, IMLAD, IMLB, IMLC, INFO, ISGN, ITRANA,
 112      $                   ITRANB, J, M, N
 113       REAL               BIGNUM, EPS, RES, RES1, SCALE, SMLNUM, TNRM,
 114      $                   XNRM
 115       COMPLEX            RMUL
 116 *     ..
 117 *     .. Local Arrays ..
 118       REAL               DUM( 1 ), VM1( 3 ), VM2( 3 )
 119       COMPLEX            A( LDT, LDT ), ATMP( LDT, LDT ), B( LDT, LDT ),
 120      $                   BTMP( LDT, LDT ), C( LDT, LDT ),
 121      $                   CSAV( LDT, LDT ), CTMP( LDT, LDT )
 122 *     ..
 123 *     .. External Functions ..
 124       REAL               CLANGE, SLAMCH
 125       EXTERNAL           CLANGE, SLAMCH
 126 *     ..
 127 *     .. External Subroutines ..
 128       EXTERNAL           CGEMM, CTRSYL
 129 *     ..
 130 *     .. Intrinsic Functions ..
 131       INTRINSIC          ABS, MAX, REAL, SQRT
 132 *     ..
 133 *     .. Executable Statements ..
 134 *
 135 *     Get machine parameters
 136 *
 137       EPS = SLAMCH( 'P' )
 138       SMLNUM = SLAMCH( 'S' ) / EPS
 139       BIGNUM = ONE / SMLNUM
 140       CALL SLABAD( SMLNUM, BIGNUM )
 141 *
 142 *     Set up test case parameters
 143 *
 144       VM1( 1 ) = SQRT( SMLNUM )
 145       VM1( 2 ) = ONE
 146       VM1( 3 ) = LARGE
 147       VM2( 1 ) = ONE
 148       VM2( 2 ) = ONE + TWO*EPS
 149       VM2( 3 ) = TWO
 150 *
 151       KNT = 0
 152       NINFO = 0
 153       LMAX = 0
 154       RMAX = ZERO
 155 *
 156 *     Begin test loop
 157 *
 158    10 CONTINUE
 159       READ( NIN, FMT = * )M, N
 160       IF( N.EQ.0 )
 161      $   RETURN
 162       DO 20 I = 1, M
 163          READ( NIN, FMT = * )( ATMP( I, J ), J = 1, M )
 164    20 CONTINUE
 165       DO 30 I = 1, N
 166          READ( NIN, FMT = * )( BTMP( I, J ), J = 1, N )
 167    30 CONTINUE
 168       DO 40 I = 1, M
 169          READ( NIN, FMT = * )( CTMP( I, J ), J = 1, N )
 170    40 CONTINUE
 171       DO 170 IMLA = 1, 3
 172          DO 160 IMLAD = 1, 3
 173             DO 150 IMLB = 1, 3
 174                DO 140 IMLC = 1, 3
 175                   DO 130 ITRANA = 1, 2
 176                      DO 120 ITRANB = 1, 2
 177                         DO 110 ISGN = -1, 1, 2
 178                            IF( ITRANA.EQ.1 )
 179      $                        TRANA = 'N'
 180                            IF( ITRANA.EQ.2 )
 181      $                        TRANA = 'C'
 182                            IF( ITRANB.EQ.1 )
 183      $                        TRANB = 'N'
 184                            IF( ITRANB.EQ.2 )
 185      $                        TRANB = 'C'
 186                            TNRM = ZERO
 187                            DO 60 I = 1, M
 188                               DO 50 J = 1, M
 189                                  A( I, J ) = ATMP( I, J )*VM1( IMLA )
 190                                  TNRM = MAX( TNRM, ABS( A( I, J ) ) )
 191    50                         CONTINUE
 192                               A( I, I ) = A( I, I )*VM2( IMLAD )
 193                               TNRM = MAX( TNRM, ABS( A( I, I ) ) )
 194    60                      CONTINUE
 195                            DO 80 I = 1, N
 196                               DO 70 J = 1, N
 197                                  B( I, J ) = BTMP( I, J )*VM1( IMLB )
 198                                  TNRM = MAX( TNRM, ABS( B( I, J ) ) )
 199    70                         CONTINUE
 200    80                      CONTINUE
 201                            IF( TNRM.EQ.ZERO )
 202      $                        TNRM = ONE
 203                            DO 100 I = 1, M
 204                               DO 90 J = 1, N
 205                                  C( I, J ) = CTMP( I, J )*VM1( IMLC )
 206                                  CSAV( I, J ) = C( I, J )
 207    90                         CONTINUE
 208   100                      CONTINUE
 209                            KNT = KNT + 1
 210                            CALL CTRSYL( TRANA, TRANB, ISGN, M, N, A,
 211      $                                  LDT, B, LDT, C, LDT, SCALE,
 212      $                                  INFO )
 213                            IF( INFO.NE.0 )
 214      $                        NINFO = NINFO + 1
 215                            XNRM = CLANGE( 'M', M, N, C, LDT, DUM )
 216                            RMUL = CONE
 217                            IF( XNRM.GT.ONE .AND. TNRM.GT.ONE ) THEN
 218                               IF( XNRM.GT.BIGNUM / TNRM ) THEN
 219                                  RMUL = MAX( XNRM, TNRM )
 220                                  RMUL = CONE / RMUL
 221                               END IF
 222                            END IF
 223                            CALL CGEMM( TRANA, 'N', M, N, M, RMUL, A,
 224      $                                 LDT, C, LDT, -SCALE*RMUL, CSAV,
 225      $                                 LDT )
 226                            CALL CGEMM( 'N', TRANB, M, N, N,
 227      $                                 REAL( ISGN )*RMUL, C, LDT, B,
 228      $                                 LDT, CONE, CSAV, LDT )
 229                            RES1 = CLANGE( 'M', M, N, CSAV, LDT, DUM )
 230                            RES = RES1 / MAX( SMLNUM, SMLNUM*XNRM,
 231      $                           ( ( ABS( RMUL )*TNRM )*EPS )*XNRM )
 232                            IF( RES.GT.RMAX ) THEN
 233                               LMAX = KNT
 234                               RMAX = RES
 235                            END IF
 236   110                   CONTINUE
 237   120                CONTINUE
 238   130             CONTINUE
 239   140          CONTINUE
 240   150       CONTINUE
 241   160    CONTINUE
 242   170 CONTINUE
 243       GO TO 10
 244 *
 245 *     End of CGET35
 246 *
 247       END