TESTING/EIG/zlatm4.f

   1 *> \brief \b ZLATM4
   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 ZLATM4( ITYPE, N, NZ1, NZ2, RSIGN, AMAGN, RCOND,
  12 *                          TRIANG, IDIST, ISEED, A, LDA )
  13 *
  14 *       .. Scalar Arguments ..
  15 *       LOGICAL            RSIGN
  16 *       INTEGER            IDIST, ITYPE, LDA, N, NZ1, NZ2
  17 *       DOUBLE PRECISION   AMAGN, RCOND, TRIANG
  18 *       ..
  19 *       .. Array Arguments ..
  20 *       INTEGER            ISEED( 4 )
  21 *       COMPLEX*16         A( LDA, * )
  22 *       ..
  23 *
  24 *
  25 *> \par Purpose:
  26 *  =============
  27 *>
  28 *> \verbatim
  29 *>
  30 *> ZLATM4 generates basic square matrices, which may later be
  31 *> multiplied by others in order to produce test matrices.  It is
  32 *> intended mainly to be used to test the generalized eigenvalue
  33 *> routines.
  34 *>
  35 *> It first generates the diagonal and (possibly) subdiagonal,
  36 *> according to the value of ITYPE, NZ1, NZ2, RSIGN, AMAGN, and RCOND.
  37 *> It then fills in the upper triangle with random numbers, if TRIANG is
  38 *> non-zero.
  39 *> \endverbatim
  40 *
  41 *  Arguments:
  42 *  ==========
  43 *
  44 *> \param[in] ITYPE
  45 *> \verbatim
  46 *>          ITYPE is INTEGER
  47 *>          The "type" of matrix on the diagonal and sub-diagonal.
  48 *>          If ITYPE < 0, then type abs(ITYPE) is generated and then
  49 *>             swapped end for end (A(I,J) := A'(N-J,N-I).)  See also
  50 *>             the description of AMAGN and RSIGN.
  51 *>
  52 *>          Special types:
  53 *>          = 0:  the zero matrix.
  54 *>          = 1:  the identity.
  55 *>          = 2:  a transposed Jordan block.
  56 *>          = 3:  If N is odd, then a k+1 x k+1 transposed Jordan block
  57 *>                followed by a k x k identity block, where k=(N-1)/2.
  58 *>                If N is even, then k=(N-2)/2, and a zero diagonal entry
  59 *>                is tacked onto the end.
  60 *>
  61 *>          Diagonal types.  The diagonal consists of NZ1 zeros, then
  62 *>             k=N-NZ1-NZ2 nonzeros.  The subdiagonal is zero.  ITYPE
  63 *>             specifies the nonzero diagonal entries as follows:
  64 *>          = 4:  1, ..., k
  65 *>          = 5:  1, RCOND, ..., RCOND
  66 *>          = 6:  1, ..., 1, RCOND
  67 *>          = 7:  1, a, a^2, ..., a^(k-1)=RCOND
  68 *>          = 8:  1, 1-d, 1-2*d, ..., 1-(k-1)*d=RCOND
  69 *>          = 9:  random numbers chosen from (RCOND,1)
  70 *>          = 10: random numbers with distribution IDIST (see ZLARND.)
  71 *> \endverbatim
  72 *>
  73 *> \param[in] N
  74 *> \verbatim
  75 *>          N is INTEGER
  76 *>          The order of the matrix.
  77 *> \endverbatim
  78 *>
  79 *> \param[in] NZ1
  80 *> \verbatim
  81 *>          NZ1 is INTEGER
  82 *>          If abs(ITYPE) > 3, then the first NZ1 diagonal entries will
  83 *>          be zero.
  84 *> \endverbatim
  85 *>
  86 *> \param[in] NZ2
  87 *> \verbatim
  88 *>          NZ2 is INTEGER
  89 *>          If abs(ITYPE) > 3, then the last NZ2 diagonal entries will
  90 *>          be zero.
  91 *> \endverbatim
  92 *>
  93 *> \param[in] RSIGN
  94 *> \verbatim
  95 *>          RSIGN is LOGICAL
  96 *>          = .TRUE.:  The diagonal and subdiagonal entries will be
  97 *>                     multiplied by random numbers of magnitude 1.
  98 *>          = .FALSE.: The diagonal and subdiagonal entries will be
  99 *>                     left as they are (usually non-negative real.)
 100 *> \endverbatim
 101 *>
 102 *> \param[in] AMAGN
 103 *> \verbatim
 104 *>          AMAGN is DOUBLE PRECISION
 105 *>          The diagonal and subdiagonal entries will be multiplied by
 106 *>          AMAGN.
 107 *> \endverbatim
 108 *>
 109 *> \param[in] RCOND
 110 *> \verbatim
 111 *>          RCOND is DOUBLE PRECISION
 112 *>          If abs(ITYPE) > 4, then the smallest diagonal entry will be
 113 *>          RCOND.  RCOND must be between 0 and 1.
 114 *> \endverbatim
 115 *>
 116 *> \param[in] TRIANG
 117 *> \verbatim
 118 *>          TRIANG is DOUBLE PRECISION
 119 *>          The entries above the diagonal will be random numbers with
 120 *>          magnitude bounded by TRIANG (i.e., random numbers multiplied
 121 *>          by TRIANG.)
 122 *> \endverbatim
 123 *>
 124 *> \param[in] IDIST
 125 *> \verbatim
 126 *>          IDIST is INTEGER
 127 *>          On entry, DIST specifies the type of distribution to be used
 128 *>          to generate a random matrix .
 129 *>          = 1: real and imaginary parts each UNIFORM( 0, 1 )
 130 *>          = 2: real and imaginary parts each UNIFORM( -1, 1 )
 131 *>          = 3: real and imaginary parts each NORMAL( 0, 1 )
 132 *>          = 4: complex number uniform in DISK( 0, 1 )
 133 *> \endverbatim
 134 *>
 135 *> \param[in,out] ISEED
 136 *> \verbatim
 137 *>          ISEED is INTEGER array, dimension (4)
 138 *>          On entry ISEED specifies the seed of the random number
 139 *>          generator.  The values of ISEED are changed on exit, and can
 140 *>          be used in the next call to ZLATM4 to continue the same
 141 *>          random number sequence.
 142 *>          Note: ISEED(4) should be odd, for the random number generator
 143 *>          used at present.
 144 *> \endverbatim
 145 *>
 146 *> \param[out] A
 147 *> \verbatim
 148 *>          A is COMPLEX*16 array, dimension (LDA, N)
 149 *>          Array to be computed.
 150 *> \endverbatim
 151 *>
 152 *> \param[in] LDA
 153 *> \verbatim
 154 *>          LDA is INTEGER
 155 *>          Leading dimension of A.  Must be at least 1 and at least N.
 156 *> \endverbatim
 157 *
 158 *  Authors:
 159 *  ========
 160 *
 161 *> \author Univ. of Tennessee
 162 *> \author Univ. of California Berkeley
 163 *> \author Univ. of Colorado Denver
 164 *> \author NAG Ltd.
 165 *
 166 *> \date November 2011
 167 *
 168 *> \ingroup complex16_eig
 169 *
 170 *  =====================================================================
 171       SUBROUTINE ZLATM4( ITYPE, N, NZ1, NZ2, RSIGN, AMAGN, RCOND,
 172      $                   TRIANG, IDIST, ISEED, A, LDA )
 173 *
 174 *  -- LAPACK test routine (version 3.4.0) --
 175 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 176 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 177 *     November 2011
 178 *
 179 *     .. Scalar Arguments ..
 180       LOGICAL            RSIGN
 181       INTEGER            IDIST, ITYPE, LDA, N, NZ1, NZ2
 182       DOUBLE PRECISION   AMAGN, RCOND, TRIANG
 183 *     ..
 184 *     .. Array Arguments ..
 185       INTEGER            ISEED( 4 )
 186       COMPLEX*16         A( LDA, * )
 187 *     ..
 188 *
 189 *  =====================================================================
 190 *
 191 *     .. Parameters ..
 192       DOUBLE PRECISION   ZERO, ONE
 193       PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0 )
 194       COMPLEX*16         CZERO, CONE
 195       PARAMETER          ( CZERO = ( 0.0D+0, 0.0D+0 ),
 196      $                   CONE = ( 1.0D+0, 0.0D+0 ) )
 197 *     ..
 198 *     .. Local Scalars ..
 199       INTEGER            I, ISDB, ISDE, JC, JD, JR, K, KBEG, KEND, KLEN
 200       DOUBLE PRECISION   ALPHA
 201       COMPLEX*16         CTEMP
 202 *     ..
 203 *     .. External Functions ..
 204       DOUBLE PRECISION   DLARAN
 205       COMPLEX*16         ZLARND
 206       EXTERNAL           DLARAN, ZLARND
 207 *     ..
 208 *     .. External Subroutines ..
 209       EXTERNAL           ZLASET
 210 *     ..
 211 *     .. Intrinsic Functions ..
 212       INTRINSIC          ABS, DBLE, DCMPLX, EXP, LOG, MAX, MIN, MOD
 213 *     ..
 214 *     .. Executable Statements ..
 215 *
 216       IF( N.LE.0 )
 217      $   RETURN
 218       CALL ZLASET( 'Full', N, N, CZERO, CZERO, A, LDA )
 219 *
 220 *     Insure a correct ISEED
 221 *
 222       IF( MOD( ISEED( 4 ), 2 ).NE.1 )
 223      $   ISEED( 4 ) = ISEED( 4 ) + 1
 224 *
 225 *     Compute diagonal and subdiagonal according to ITYPE, NZ1, NZ2,
 226 *     and RCOND
 227 *
 228       IF( ITYPE.NE.0 ) THEN
 229          IF( ABS( ITYPE ).GE.4 ) THEN
 230             KBEG = MAX( 1, MIN( N, NZ1+1 ) )
 231             KEND = MAX( KBEG, MIN( N, N-NZ2 ) )
 232             KLEN = KEND + 1 - KBEG
 233          ELSE
 234             KBEG = 1
 235             KEND = N
 236             KLEN = N
 237          END IF
 238          ISDB = 1
 239          ISDE = 0
 240          GO TO ( 10, 30, 50, 80, 100, 120, 140, 160,
 241      $           180, 200 )ABS( ITYPE )
 242 *
 243 *        abs(ITYPE) = 1: Identity
 244 *
 245    10    CONTINUE
 246          DO 20 JD = 1, N
 247             A( JD, JD ) = CONE
 248    20    CONTINUE
 249          GO TO 220
 250 *
 251 *        abs(ITYPE) = 2: Transposed Jordan block
 252 *
 253    30    CONTINUE
 254          DO 40 JD = 1, N - 1
 255             A( JD+1, JD ) = CONE
 256    40    CONTINUE
 257          ISDB = 1
 258          ISDE = N - 1
 259          GO TO 220
 260 *
 261 *        abs(ITYPE) = 3: Transposed Jordan block, followed by the
 262 *                        identity.
 263 *
 264    50    CONTINUE
 265          K = ( N-1 ) / 2
 266          DO 60 JD = 1, K
 267             A( JD+1, JD ) = CONE
 268    60    CONTINUE
 269          ISDB = 1
 270          ISDE = K
 271          DO 70 JD = K + 2, 2*K + 1
 272             A( JD, JD ) = CONE
 273    70    CONTINUE
 274          GO TO 220
 275 *
 276 *        abs(ITYPE) = 4: 1,...,k
 277 *
 278    80    CONTINUE
 279          DO 90 JD = KBEG, KEND
 280             A( JD, JD ) = DCMPLX( JD-NZ1 )
 281    90    CONTINUE
 282          GO TO 220
 283 *
 284 *        abs(ITYPE) = 5: One large D value:
 285 *
 286   100    CONTINUE
 287          DO 110 JD = KBEG + 1, KEND
 288             A( JD, JD ) = DCMPLX( RCOND )
 289   110    CONTINUE
 290          A( KBEG, KBEG ) = CONE
 291          GO TO 220
 292 *
 293 *        abs(ITYPE) = 6: One small D value:
 294 *
 295   120    CONTINUE
 296          DO 130 JD = KBEG, KEND - 1
 297             A( JD, JD ) = CONE
 298   130    CONTINUE
 299          A( KEND, KEND ) = DCMPLX( RCOND )
 300          GO TO 220
 301 *
 302 *        abs(ITYPE) = 7: Exponentially distributed D values:
 303 *
 304   140    CONTINUE
 305          A( KBEG, KBEG ) = CONE
 306          IF( KLEN.GT.1 ) THEN
 307             ALPHA = RCOND**( ONE / DBLE( KLEN-1 ) )
 308             DO 150 I = 2, KLEN
 309                A( NZ1+I, NZ1+I ) = DCMPLX( ALPHA**DBLE( I-1 ) )
 310   150       CONTINUE
 311          END IF
 312          GO TO 220
 313 *
 314 *        abs(ITYPE) = 8: Arithmetically distributed D values:
 315 *
 316   160    CONTINUE
 317          A( KBEG, KBEG ) = CONE
 318          IF( KLEN.GT.1 ) THEN
 319             ALPHA = ( ONE-RCOND ) / DBLE( KLEN-1 )
 320             DO 170 I = 2, KLEN
 321                A( NZ1+I, NZ1+I ) = DCMPLX( DBLE( KLEN-I )*ALPHA+RCOND )
 322   170       CONTINUE
 323          END IF
 324          GO TO 220
 325 *
 326 *        abs(ITYPE) = 9: Randomly distributed D values on ( RCOND, 1):
 327 *
 328   180    CONTINUE
 329          ALPHA = LOG( RCOND )
 330          DO 190 JD = KBEG, KEND
 331             A( JD, JD ) = EXP( ALPHA*DLARAN( ISEED ) )
 332   190    CONTINUE
 333          GO TO 220
 334 *
 335 *        abs(ITYPE) = 10: Randomly distributed D values from DIST
 336 *
 337   200    CONTINUE
 338          DO 210 JD = KBEG, KEND
 339             A( JD, JD ) = ZLARND( IDIST, ISEED )
 340   210    CONTINUE
 341 *
 342   220    CONTINUE
 343 *
 344 *        Scale by AMAGN
 345 *
 346          DO 230 JD = KBEG, KEND
 347             A( JD, JD ) = AMAGN*DBLE( A( JD, JD ) )
 348   230    CONTINUE
 349          DO 240 JD = ISDB, ISDE
 350             A( JD+1, JD ) = AMAGN*DBLE( A( JD+1, JD ) )
 351   240    CONTINUE
 352 *
 353 *        If RSIGN = .TRUE., assign random signs to diagonal and
 354 *        subdiagonal
 355 *
 356          IF( RSIGN ) THEN
 357             DO 250 JD = KBEG, KEND
 358                IF( DBLE( A( JD, JD ) ).NE.ZERO ) THEN
 359                   CTEMP = ZLARND( 3, ISEED )
 360                   CTEMP = CTEMP / ABS( CTEMP )
 361                   A( JD, JD ) = CTEMP*DBLE( A( JD, JD ) )
 362                END IF
 363   250       CONTINUE
 364             DO 260 JD = ISDB, ISDE
 365                IF( DBLE( A( JD+1, JD ) ).NE.ZERO ) THEN
 366                   CTEMP = ZLARND( 3, ISEED )
 367                   CTEMP = CTEMP / ABS( CTEMP )
 368                   A( JD+1, JD ) = CTEMP*DBLE( A( JD+1, JD ) )
 369                END IF
 370   260       CONTINUE
 371          END IF
 372 *
 373 *        Reverse if ITYPE < 0
 374 *
 375          IF( ITYPE.LT.0 ) THEN
 376             DO 270 JD = KBEG, ( KBEG+KEND-1 ) / 2
 377                CTEMP = A( JD, JD )
 378                A( JD, JD ) = A( KBEG+KEND-JD, KBEG+KEND-JD )
 379                A( KBEG+KEND-JD, KBEG+KEND-JD ) = CTEMP
 380   270       CONTINUE
 381             DO 280 JD = 1, ( N-1 ) / 2
 382                CTEMP = A( JD+1, JD )
 383                A( JD+1, JD ) = A( N+1-JD, N-JD )
 384                A( N+1-JD, N-JD ) = CTEMP
 385   280       CONTINUE
 386          END IF
 387 *
 388       END IF
 389 *
 390 *     Fill in upper triangle
 391 *
 392       IF( TRIANG.NE.ZERO ) THEN
 393          DO 300 JC = 2, N
 394             DO 290 JR = 1, JC - 1
 395                A( JR, JC ) = TRIANG*ZLARND( IDIST, ISEED )
 396   290       CONTINUE
 397   300    CONTINUE
 398       END IF
 399 *
 400       RETURN
 401 *
 402 *     End of ZLATM4
 403 *
 404       END