LAPACK 3.3.0
|
00001 SUBROUTINE ZLATM4( ITYPE, N, NZ1, NZ2, RSIGN, AMAGN, RCOND, 00002 $ TRIANG, IDIST, ISEED, A, LDA ) 00003 * 00004 * -- LAPACK auxiliary test routine (version 3.1) -- 00005 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 LOGICAL RSIGN 00010 INTEGER IDIST, ITYPE, LDA, N, NZ1, NZ2 00011 DOUBLE PRECISION AMAGN, RCOND, TRIANG 00012 * .. 00013 * .. Array Arguments .. 00014 INTEGER ISEED( 4 ) 00015 COMPLEX*16 A( LDA, * ) 00016 * .. 00017 * 00018 * Purpose 00019 * ======= 00020 * 00021 * ZLATM4 generates basic square matrices, which may later be 00022 * multiplied by others in order to produce test matrices. It is 00023 * intended mainly to be used to test the generalized eigenvalue 00024 * routines. 00025 * 00026 * It first generates the diagonal and (possibly) subdiagonal, 00027 * according to the value of ITYPE, NZ1, NZ2, RSIGN, AMAGN, and RCOND. 00028 * It then fills in the upper triangle with random numbers, if TRIANG is 00029 * non-zero. 00030 * 00031 * Arguments 00032 * ========= 00033 * 00034 * ITYPE (input) INTEGER 00035 * The "type" of matrix on the diagonal and sub-diagonal. 00036 * If ITYPE < 0, then type abs(ITYPE) is generated and then 00037 * swapped end for end (A(I,J) := A'(N-J,N-I).) See also 00038 * the description of AMAGN and RSIGN. 00039 * 00040 * Special types: 00041 * = 0: the zero matrix. 00042 * = 1: the identity. 00043 * = 2: a transposed Jordan block. 00044 * = 3: If N is odd, then a k+1 x k+1 transposed Jordan block 00045 * followed by a k x k identity block, where k=(N-1)/2. 00046 * If N is even, then k=(N-2)/2, and a zero diagonal entry 00047 * is tacked onto the end. 00048 * 00049 * Diagonal types. The diagonal consists of NZ1 zeros, then 00050 * k=N-NZ1-NZ2 nonzeros. The subdiagonal is zero. ITYPE 00051 * specifies the nonzero diagonal entries as follows: 00052 * = 4: 1, ..., k 00053 * = 5: 1, RCOND, ..., RCOND 00054 * = 6: 1, ..., 1, RCOND 00055 * = 7: 1, a, a^2, ..., a^(k-1)=RCOND 00056 * = 8: 1, 1-d, 1-2*d, ..., 1-(k-1)*d=RCOND 00057 * = 9: random numbers chosen from (RCOND,1) 00058 * = 10: random numbers with distribution IDIST (see ZLARND.) 00059 * 00060 * N (input) INTEGER 00061 * The order of the matrix. 00062 * 00063 * NZ1 (input) INTEGER 00064 * If abs(ITYPE) > 3, then the first NZ1 diagonal entries will 00065 * be zero. 00066 * 00067 * NZ2 (input) INTEGER 00068 * If abs(ITYPE) > 3, then the last NZ2 diagonal entries will 00069 * be zero. 00070 * 00071 * RSIGN (input) LOGICAL 00072 * = .TRUE.: The diagonal and subdiagonal entries will be 00073 * multiplied by random numbers of magnitude 1. 00074 * = .FALSE.: The diagonal and subdiagonal entries will be 00075 * left as they are (usually non-negative real.) 00076 * 00077 * AMAGN (input) DOUBLE PRECISION 00078 * The diagonal and subdiagonal entries will be multiplied by 00079 * AMAGN. 00080 * 00081 * RCOND (input) DOUBLE PRECISION 00082 * If abs(ITYPE) > 4, then the smallest diagonal entry will be 00083 * RCOND. RCOND must be between 0 and 1. 00084 * 00085 * TRIANG (input) DOUBLE PRECISION 00086 * The entries above the diagonal will be random numbers with 00087 * magnitude bounded by TRIANG (i.e., random numbers multiplied 00088 * by TRIANG.) 00089 * 00090 * IDIST (input) INTEGER 00091 * On entry, DIST specifies the type of distribution to be used 00092 * to generate a random matrix . 00093 * = 1: real and imaginary parts each UNIFORM( 0, 1 ) 00094 * = 2: real and imaginary parts each UNIFORM( -1, 1 ) 00095 * = 3: real and imaginary parts each NORMAL( 0, 1 ) 00096 * = 4: complex number uniform in DISK( 0, 1 ) 00097 * 00098 * ISEED (input/output) INTEGER array, dimension (4) 00099 * On entry ISEED specifies the seed of the random number 00100 * generator. The values of ISEED are changed on exit, and can 00101 * be used in the next call to ZLATM4 to continue the same 00102 * random number sequence. 00103 * Note: ISEED(4) should be odd, for the random number generator 00104 * used at present. 00105 * 00106 * A (output) COMPLEX*16 array, dimension (LDA, N) 00107 * Array to be computed. 00108 * 00109 * LDA (input) INTEGER 00110 * Leading dimension of A. Must be at least 1 and at least N. 00111 * 00112 * ===================================================================== 00113 * 00114 * .. Parameters .. 00115 DOUBLE PRECISION ZERO, ONE 00116 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) 00117 COMPLEX*16 CZERO, CONE 00118 PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ), 00119 $ CONE = ( 1.0D+0, 0.0D+0 ) ) 00120 * .. 00121 * .. Local Scalars .. 00122 INTEGER I, ISDB, ISDE, JC, JD, JR, K, KBEG, KEND, KLEN 00123 DOUBLE PRECISION ALPHA 00124 COMPLEX*16 CTEMP 00125 * .. 00126 * .. External Functions .. 00127 DOUBLE PRECISION DLARAN 00128 COMPLEX*16 ZLARND 00129 EXTERNAL DLARAN, ZLARND 00130 * .. 00131 * .. External Subroutines .. 00132 EXTERNAL ZLASET 00133 * .. 00134 * .. Intrinsic Functions .. 00135 INTRINSIC ABS, DBLE, DCMPLX, EXP, LOG, MAX, MIN, MOD 00136 * .. 00137 * .. Executable Statements .. 00138 * 00139 IF( N.LE.0 ) 00140 $ RETURN 00141 CALL ZLASET( 'Full', N, N, CZERO, CZERO, A, LDA ) 00142 * 00143 * Insure a correct ISEED 00144 * 00145 IF( MOD( ISEED( 4 ), 2 ).NE.1 ) 00146 $ ISEED( 4 ) = ISEED( 4 ) + 1 00147 * 00148 * Compute diagonal and subdiagonal according to ITYPE, NZ1, NZ2, 00149 * and RCOND 00150 * 00151 IF( ITYPE.NE.0 ) THEN 00152 IF( ABS( ITYPE ).GE.4 ) THEN 00153 KBEG = MAX( 1, MIN( N, NZ1+1 ) ) 00154 KEND = MAX( KBEG, MIN( N, N-NZ2 ) ) 00155 KLEN = KEND + 1 - KBEG 00156 ELSE 00157 KBEG = 1 00158 KEND = N 00159 KLEN = N 00160 END IF 00161 ISDB = 1 00162 ISDE = 0 00163 GO TO ( 10, 30, 50, 80, 100, 120, 140, 160, 00164 $ 180, 200 )ABS( ITYPE ) 00165 * 00166 * abs(ITYPE) = 1: Identity 00167 * 00168 10 CONTINUE 00169 DO 20 JD = 1, N 00170 A( JD, JD ) = CONE 00171 20 CONTINUE 00172 GO TO 220 00173 * 00174 * abs(ITYPE) = 2: Transposed Jordan block 00175 * 00176 30 CONTINUE 00177 DO 40 JD = 1, N - 1 00178 A( JD+1, JD ) = CONE 00179 40 CONTINUE 00180 ISDB = 1 00181 ISDE = N - 1 00182 GO TO 220 00183 * 00184 * abs(ITYPE) = 3: Transposed Jordan block, followed by the 00185 * identity. 00186 * 00187 50 CONTINUE 00188 K = ( N-1 ) / 2 00189 DO 60 JD = 1, K 00190 A( JD+1, JD ) = CONE 00191 60 CONTINUE 00192 ISDB = 1 00193 ISDE = K 00194 DO 70 JD = K + 2, 2*K + 1 00195 A( JD, JD ) = CONE 00196 70 CONTINUE 00197 GO TO 220 00198 * 00199 * abs(ITYPE) = 4: 1,...,k 00200 * 00201 80 CONTINUE 00202 DO 90 JD = KBEG, KEND 00203 A( JD, JD ) = DCMPLX( JD-NZ1 ) 00204 90 CONTINUE 00205 GO TO 220 00206 * 00207 * abs(ITYPE) = 5: One large D value: 00208 * 00209 100 CONTINUE 00210 DO 110 JD = KBEG + 1, KEND 00211 A( JD, JD ) = DCMPLX( RCOND ) 00212 110 CONTINUE 00213 A( KBEG, KBEG ) = CONE 00214 GO TO 220 00215 * 00216 * abs(ITYPE) = 6: One small D value: 00217 * 00218 120 CONTINUE 00219 DO 130 JD = KBEG, KEND - 1 00220 A( JD, JD ) = CONE 00221 130 CONTINUE 00222 A( KEND, KEND ) = DCMPLX( RCOND ) 00223 GO TO 220 00224 * 00225 * abs(ITYPE) = 7: Exponentially distributed D values: 00226 * 00227 140 CONTINUE 00228 A( KBEG, KBEG ) = CONE 00229 IF( KLEN.GT.1 ) THEN 00230 ALPHA = RCOND**( ONE / DBLE( KLEN-1 ) ) 00231 DO 150 I = 2, KLEN 00232 A( NZ1+I, NZ1+I ) = DCMPLX( ALPHA**DBLE( I-1 ) ) 00233 150 CONTINUE 00234 END IF 00235 GO TO 220 00236 * 00237 * abs(ITYPE) = 8: Arithmetically distributed D values: 00238 * 00239 160 CONTINUE 00240 A( KBEG, KBEG ) = CONE 00241 IF( KLEN.GT.1 ) THEN 00242 ALPHA = ( ONE-RCOND ) / DBLE( KLEN-1 ) 00243 DO 170 I = 2, KLEN 00244 A( NZ1+I, NZ1+I ) = DCMPLX( DBLE( KLEN-I )*ALPHA+RCOND ) 00245 170 CONTINUE 00246 END IF 00247 GO TO 220 00248 * 00249 * abs(ITYPE) = 9: Randomly distributed D values on ( RCOND, 1): 00250 * 00251 180 CONTINUE 00252 ALPHA = LOG( RCOND ) 00253 DO 190 JD = KBEG, KEND 00254 A( JD, JD ) = EXP( ALPHA*DLARAN( ISEED ) ) 00255 190 CONTINUE 00256 GO TO 220 00257 * 00258 * abs(ITYPE) = 10: Randomly distributed D values from DIST 00259 * 00260 200 CONTINUE 00261 DO 210 JD = KBEG, KEND 00262 A( JD, JD ) = ZLARND( IDIST, ISEED ) 00263 210 CONTINUE 00264 * 00265 220 CONTINUE 00266 * 00267 * Scale by AMAGN 00268 * 00269 DO 230 JD = KBEG, KEND 00270 A( JD, JD ) = AMAGN*DBLE( A( JD, JD ) ) 00271 230 CONTINUE 00272 DO 240 JD = ISDB, ISDE 00273 A( JD+1, JD ) = AMAGN*DBLE( A( JD+1, JD ) ) 00274 240 CONTINUE 00275 * 00276 * If RSIGN = .TRUE., assign random signs to diagonal and 00277 * subdiagonal 00278 * 00279 IF( RSIGN ) THEN 00280 DO 250 JD = KBEG, KEND 00281 IF( DBLE( A( JD, JD ) ).NE.ZERO ) THEN 00282 CTEMP = ZLARND( 3, ISEED ) 00283 CTEMP = CTEMP / ABS( CTEMP ) 00284 A( JD, JD ) = CTEMP*DBLE( A( JD, JD ) ) 00285 END IF 00286 250 CONTINUE 00287 DO 260 JD = ISDB, ISDE 00288 IF( DBLE( A( JD+1, JD ) ).NE.ZERO ) THEN 00289 CTEMP = ZLARND( 3, ISEED ) 00290 CTEMP = CTEMP / ABS( CTEMP ) 00291 A( JD+1, JD ) = CTEMP*DBLE( A( JD+1, JD ) ) 00292 END IF 00293 260 CONTINUE 00294 END IF 00295 * 00296 * Reverse if ITYPE < 0 00297 * 00298 IF( ITYPE.LT.0 ) THEN 00299 DO 270 JD = KBEG, ( KBEG+KEND-1 ) / 2 00300 CTEMP = A( JD, JD ) 00301 A( JD, JD ) = A( KBEG+KEND-JD, KBEG+KEND-JD ) 00302 A( KBEG+KEND-JD, KBEG+KEND-JD ) = CTEMP 00303 270 CONTINUE 00304 DO 280 JD = 1, ( N-1 ) / 2 00305 CTEMP = A( JD+1, JD ) 00306 A( JD+1, JD ) = A( N+1-JD, N-JD ) 00307 A( N+1-JD, N-JD ) = CTEMP 00308 280 CONTINUE 00309 END IF 00310 * 00311 END IF 00312 * 00313 * Fill in upper triangle 00314 * 00315 IF( TRIANG.NE.ZERO ) THEN 00316 DO 300 JC = 2, N 00317 DO 290 JR = 1, JC - 1 00318 A( JR, JC ) = TRIANG*ZLARND( IDIST, ISEED ) 00319 290 CONTINUE 00320 300 CONTINUE 00321 END IF 00322 * 00323 RETURN 00324 * 00325 * End of ZLATM4 00326 * 00327 END