![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b ZLATMS 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 * Definition: 00009 * =========== 00010 * 00011 * SUBROUTINE ZLATMS( M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, 00012 * KL, KU, PACK, A, LDA, WORK, INFO ) 00013 * 00014 * .. Scalar Arguments .. 00015 * CHARACTER DIST, PACK, SYM 00016 * INTEGER INFO, KL, KU, LDA, M, MODE, N 00017 * DOUBLE PRECISION COND, DMAX 00018 * .. 00019 * .. Array Arguments .. 00020 * INTEGER ISEED( 4 ) 00021 * DOUBLE PRECISION D( * ) 00022 * COMPLEX*16 A( LDA, * ), WORK( * ) 00023 * .. 00024 * 00025 * 00026 *> \par Purpose: 00027 * ============= 00028 *> 00029 *> \verbatim 00030 *> 00031 *> ZLATMS generates random matrices with specified singular values 00032 *> (or hermitian with specified eigenvalues) 00033 *> for testing LAPACK programs. 00034 *> 00035 *> ZLATMS operates by applying the following sequence of 00036 *> operations: 00037 *> 00038 *> Set the diagonal to D, where D may be input or 00039 *> computed according to MODE, COND, DMAX, and SYM 00040 *> as described below. 00041 *> 00042 *> Generate a matrix with the appropriate band structure, by one 00043 *> of two methods: 00044 *> 00045 *> Method A: 00046 *> Generate a dense M x N matrix by multiplying D on the left 00047 *> and the right by random unitary matrices, then: 00048 *> 00049 *> Reduce the bandwidth according to KL and KU, using 00050 *> Householder transformations. 00051 *> 00052 *> Method B: 00053 *> Convert the bandwidth-0 (i.e., diagonal) matrix to a 00054 *> bandwidth-1 matrix using Givens rotations, "chasing" 00055 *> out-of-band elements back, much as in QR; then convert 00056 *> the bandwidth-1 to a bandwidth-2 matrix, etc. Note 00057 *> that for reasonably small bandwidths (relative to M and 00058 *> N) this requires less storage, as a dense matrix is not 00059 *> generated. Also, for hermitian or symmetric matrices, 00060 *> only one triangle is generated. 00061 *> 00062 *> Method A is chosen if the bandwidth is a large fraction of the 00063 *> order of the matrix, and LDA is at least M (so a dense 00064 *> matrix can be stored.) Method B is chosen if the bandwidth 00065 *> is small (< 1/2 N for hermitian or symmetric, < .3 N+M for 00066 *> non-symmetric), or LDA is less than M and not less than the 00067 *> bandwidth. 00068 *> 00069 *> Pack the matrix if desired. Options specified by PACK are: 00070 *> no packing 00071 *> zero out upper half (if hermitian) 00072 *> zero out lower half (if hermitian) 00073 *> store the upper half columnwise (if hermitian or upper 00074 *> triangular) 00075 *> store the lower half columnwise (if hermitian or lower 00076 *> triangular) 00077 *> store the lower triangle in banded format (if hermitian or 00078 *> lower triangular) 00079 *> store the upper triangle in banded format (if hermitian or 00080 *> upper triangular) 00081 *> store the entire matrix in banded format 00082 *> If Method B is chosen, and band format is specified, then the 00083 *> matrix will be generated in the band format, so no repacking 00084 *> will be necessary. 00085 *> \endverbatim 00086 * 00087 * Arguments: 00088 * ========== 00089 * 00090 *> \param[in] M 00091 *> \verbatim 00092 *> M is INTEGER 00093 *> The number of rows of A. Not modified. 00094 *> \endverbatim 00095 *> 00096 *> \param[in] N 00097 *> \verbatim 00098 *> N is INTEGER 00099 *> The number of columns of A. N must equal M if the matrix 00100 *> is symmetric or hermitian (i.e., if SYM is not 'N') 00101 *> Not modified. 00102 *> \endverbatim 00103 *> 00104 *> \param[in] DIST 00105 *> \verbatim 00106 *> DIST is CHARACTER*1 00107 *> On entry, DIST specifies the type of distribution to be used 00108 *> to generate the random eigen-/singular values. 00109 *> 'U' => UNIFORM( 0, 1 ) ( 'U' for uniform ) 00110 *> 'S' => UNIFORM( -1, 1 ) ( 'S' for symmetric ) 00111 *> 'N' => NORMAL( 0, 1 ) ( 'N' for normal ) 00112 *> Not modified. 00113 *> \endverbatim 00114 *> 00115 *> \param[in,out] ISEED 00116 *> \verbatim 00117 *> ISEED is INTEGER array, dimension ( 4 ) 00118 *> On entry ISEED specifies the seed of the random number 00119 *> generator. They should lie between 0 and 4095 inclusive, 00120 *> and ISEED(4) should be odd. The random number generator 00121 *> uses a linear congruential sequence limited to small 00122 *> integers, and so should produce machine independent 00123 *> random numbers. The values of ISEED are changed on 00124 *> exit, and can be used in the next call to ZLATMS 00125 *> to continue the same random number sequence. 00126 *> Changed on exit. 00127 *> \endverbatim 00128 *> 00129 *> \param[in] SYM 00130 *> \verbatim 00131 *> SYM is CHARACTER*1 00132 *> If SYM='H', the generated matrix is hermitian, with 00133 *> eigenvalues specified by D, COND, MODE, and DMAX; they 00134 *> may be positive, negative, or zero. 00135 *> If SYM='P', the generated matrix is hermitian, with 00136 *> eigenvalues (= singular values) specified by D, COND, 00137 *> MODE, and DMAX; they will not be negative. 00138 *> If SYM='N', the generated matrix is nonsymmetric, with 00139 *> singular values specified by D, COND, MODE, and DMAX; 00140 *> they will not be negative. 00141 *> If SYM='S', the generated matrix is (complex) symmetric, 00142 *> with singular values specified by D, COND, MODE, and 00143 *> DMAX; they will not be negative. 00144 *> Not modified. 00145 *> \endverbatim 00146 *> 00147 *> \param[in,out] D 00148 *> \verbatim 00149 *> D is DOUBLE PRECISION array, dimension ( MIN( M, N ) ) 00150 *> This array is used to specify the singular values or 00151 *> eigenvalues of A (see SYM, above.) If MODE=0, then D is 00152 *> assumed to contain the singular/eigenvalues, otherwise 00153 *> they will be computed according to MODE, COND, and DMAX, 00154 *> and placed in D. 00155 *> Modified if MODE is nonzero. 00156 *> \endverbatim 00157 *> 00158 *> \param[in] MODE 00159 *> \verbatim 00160 *> MODE is INTEGER 00161 *> On entry this describes how the singular/eigenvalues are to 00162 *> be specified: 00163 *> MODE = 0 means use D as input 00164 *> MODE = 1 sets D(1)=1 and D(2:N)=1.0/COND 00165 *> MODE = 2 sets D(1:N-1)=1 and D(N)=1.0/COND 00166 *> MODE = 3 sets D(I)=COND**(-(I-1)/(N-1)) 00167 *> MODE = 4 sets D(i)=1 - (i-1)/(N-1)*(1 - 1/COND) 00168 *> MODE = 5 sets D to random numbers in the range 00169 *> ( 1/COND , 1 ) such that their logarithms 00170 *> are uniformly distributed. 00171 *> MODE = 6 set D to random numbers from same distribution 00172 *> as the rest of the matrix. 00173 *> MODE < 0 has the same meaning as ABS(MODE), except that 00174 *> the order of the elements of D is reversed. 00175 *> Thus if MODE is positive, D has entries ranging from 00176 *> 1 to 1/COND, if negative, from 1/COND to 1, 00177 *> If SYM='H', and MODE is neither 0, 6, nor -6, then 00178 *> the elements of D will also be multiplied by a random 00179 *> sign (i.e., +1 or -1.) 00180 *> Not modified. 00181 *> \endverbatim 00182 *> 00183 *> \param[in] COND 00184 *> \verbatim 00185 *> COND is DOUBLE PRECISION 00186 *> On entry, this is used as described under MODE above. 00187 *> If used, it must be >= 1. Not modified. 00188 *> \endverbatim 00189 *> 00190 *> \param[in] DMAX 00191 *> \verbatim 00192 *> DMAX is DOUBLE PRECISION 00193 *> If MODE is neither -6, 0 nor 6, the contents of D, as 00194 *> computed according to MODE and COND, will be scaled by 00195 *> DMAX / max(abs(D(i))); thus, the maximum absolute eigen- or 00196 *> singular value (which is to say the norm) will be abs(DMAX). 00197 *> Note that DMAX need not be positive: if DMAX is negative 00198 *> (or zero), D will be scaled by a negative number (or zero). 00199 *> Not modified. 00200 *> \endverbatim 00201 *> 00202 *> \param[in] KL 00203 *> \verbatim 00204 *> KL is INTEGER 00205 *> This specifies the lower bandwidth of the matrix. For 00206 *> example, KL=0 implies upper triangular, KL=1 implies upper 00207 *> Hessenberg, and KL being at least M-1 means that the matrix 00208 *> has full lower bandwidth. KL must equal KU if the matrix 00209 *> is symmetric or hermitian. 00210 *> Not modified. 00211 *> \endverbatim 00212 *> 00213 *> \param[in] KU 00214 *> \verbatim 00215 *> KU is INTEGER 00216 *> This specifies the upper bandwidth of the matrix. For 00217 *> example, KU=0 implies lower triangular, KU=1 implies lower 00218 *> Hessenberg, and KU being at least N-1 means that the matrix 00219 *> has full upper bandwidth. KL must equal KU if the matrix 00220 *> is symmetric or hermitian. 00221 *> Not modified. 00222 *> \endverbatim 00223 *> 00224 *> \param[in] PACK 00225 *> \verbatim 00226 *> PACK is CHARACTER*1 00227 *> This specifies packing of matrix as follows: 00228 *> 'N' => no packing 00229 *> 'U' => zero out all subdiagonal entries (if symmetric 00230 *> or hermitian) 00231 *> 'L' => zero out all superdiagonal entries (if symmetric 00232 *> or hermitian) 00233 *> 'C' => store the upper triangle columnwise (only if the 00234 *> matrix is symmetric, hermitian, or upper triangular) 00235 *> 'R' => store the lower triangle columnwise (only if the 00236 *> matrix is symmetric, hermitian, or lower triangular) 00237 *> 'B' => store the lower triangle in band storage scheme 00238 *> (only if the matrix is symmetric, hermitian, or 00239 *> lower triangular) 00240 *> 'Q' => store the upper triangle in band storage scheme 00241 *> (only if the matrix is symmetric, hermitian, or 00242 *> upper triangular) 00243 *> 'Z' => store the entire matrix in band storage scheme 00244 *> (pivoting can be provided for by using this 00245 *> option to store A in the trailing rows of 00246 *> the allocated storage) 00247 *> 00248 *> Using these options, the various LAPACK packed and banded 00249 *> storage schemes can be obtained: 00250 *> GB - use 'Z' 00251 *> PB, SB, HB, or TB - use 'B' or 'Q' 00252 *> PP, SP, HB, or TP - use 'C' or 'R' 00253 *> 00254 *> If two calls to ZLATMS differ only in the PACK parameter, 00255 *> they will generate mathematically equivalent matrices. 00256 *> Not modified. 00257 *> \endverbatim 00258 *> 00259 *> \param[in,out] A 00260 *> \verbatim 00261 *> A is COMPLEX*16 array, dimension ( LDA, N ) 00262 *> On exit A is the desired test matrix. A is first generated 00263 *> in full (unpacked) form, and then packed, if so specified 00264 *> by PACK. Thus, the first M elements of the first N 00265 *> columns will always be modified. If PACK specifies a 00266 *> packed or banded storage scheme, all LDA elements of the 00267 *> first N columns will be modified; the elements of the 00268 *> array which do not correspond to elements of the generated 00269 *> matrix are set to zero. 00270 *> Modified. 00271 *> \endverbatim 00272 *> 00273 *> \param[in] LDA 00274 *> \verbatim 00275 *> LDA is INTEGER 00276 *> LDA specifies the first dimension of A as declared in the 00277 *> calling program. If PACK='N', 'U', 'L', 'C', or 'R', then 00278 *> LDA must be at least M. If PACK='B' or 'Q', then LDA must 00279 *> be at least MIN( KL, M-1) (which is equal to MIN(KU,N-1)). 00280 *> If PACK='Z', LDA must be large enough to hold the packed 00281 *> array: MIN( KU, N-1) + MIN( KL, M-1) + 1. 00282 *> Not modified. 00283 *> \endverbatim 00284 *> 00285 *> \param[out] WORK 00286 *> \verbatim 00287 *> WORK is COMPLEX*16 array, dimension ( 3*MAX( N, M ) ) 00288 *> Workspace. 00289 *> Modified. 00290 *> \endverbatim 00291 *> 00292 *> \param[out] INFO 00293 *> \verbatim 00294 *> INFO is INTEGER 00295 *> Error code. On exit, INFO will be set to one of the 00296 *> following values: 00297 *> 0 => normal return 00298 *> -1 => M negative or unequal to N and SYM='S', 'H', or 'P' 00299 *> -2 => N negative 00300 *> -3 => DIST illegal string 00301 *> -5 => SYM illegal string 00302 *> -7 => MODE not in range -6 to 6 00303 *> -8 => COND less than 1.0, and MODE neither -6, 0 nor 6 00304 *> -10 => KL negative 00305 *> -11 => KU negative, or SYM is not 'N' and KU is not equal to 00306 *> KL 00307 *> -12 => PACK illegal string, or PACK='U' or 'L', and SYM='N'; 00308 *> or PACK='C' or 'Q' and SYM='N' and KL is not zero; 00309 *> or PACK='R' or 'B' and SYM='N' and KU is not zero; 00310 *> or PACK='U', 'L', 'C', 'R', 'B', or 'Q', and M is not 00311 *> N. 00312 *> -14 => LDA is less than M, or PACK='Z' and LDA is less than 00313 *> MIN(KU,N-1) + MIN(KL,M-1) + 1. 00314 *> 1 => Error return from DLATM1 00315 *> 2 => Cannot scale to DMAX (max. sing. value is 0) 00316 *> 3 => Error return from ZLAGGE, CLAGHE or CLAGSY 00317 *> \endverbatim 00318 * 00319 * Authors: 00320 * ======== 00321 * 00322 *> \author Univ. of Tennessee 00323 *> \author Univ. of California Berkeley 00324 *> \author Univ. of Colorado Denver 00325 *> \author NAG Ltd. 00326 * 00327 *> \date November 2011 00328 * 00329 *> \ingroup complex16_matgen 00330 * 00331 * ===================================================================== 00332 SUBROUTINE ZLATMS( M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, 00333 $ KL, KU, PACK, A, LDA, WORK, INFO ) 00334 * 00335 * -- LAPACK computational routine (version 3.4.0) -- 00336 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00337 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00338 * November 2011 00339 * 00340 * .. Scalar Arguments .. 00341 CHARACTER DIST, PACK, SYM 00342 INTEGER INFO, KL, KU, LDA, M, MODE, N 00343 DOUBLE PRECISION COND, DMAX 00344 * .. 00345 * .. Array Arguments .. 00346 INTEGER ISEED( 4 ) 00347 DOUBLE PRECISION D( * ) 00348 COMPLEX*16 A( LDA, * ), WORK( * ) 00349 * .. 00350 * 00351 * ===================================================================== 00352 * 00353 * .. Parameters .. 00354 DOUBLE PRECISION ZERO 00355 PARAMETER ( ZERO = 0.0D+0 ) 00356 DOUBLE PRECISION ONE 00357 PARAMETER ( ONE = 1.0D+0 ) 00358 COMPLEX*16 CZERO 00359 PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ) ) 00360 DOUBLE PRECISION TWOPI 00361 PARAMETER ( TWOPI = 6.2831853071795864769252867663D+0 ) 00362 * .. 00363 * .. Local Scalars .. 00364 LOGICAL GIVENS, ILEXTR, ILTEMP, TOPDWN, ZSYM 00365 INTEGER I, IC, ICOL, IDIST, IENDCH, IINFO, IL, ILDA, 00366 $ IOFFG, IOFFST, IPACK, IPACKG, IR, IR1, IR2, 00367 $ IROW, IRSIGN, ISKEW, ISYM, ISYMPK, J, JC, JCH, 00368 $ JKL, JKU, JR, K, LLB, MINLDA, MNMIN, MR, NC, 00369 $ UUB 00370 DOUBLE PRECISION ALPHA, ANGLE, REALC, TEMP 00371 COMPLEX*16 C, CT, CTEMP, DUMMY, EXTRA, S, ST 00372 * .. 00373 * .. External Functions .. 00374 LOGICAL LSAME 00375 DOUBLE PRECISION DLARND 00376 COMPLEX*16 ZLARND 00377 EXTERNAL LSAME, DLARND, ZLARND 00378 * .. 00379 * .. External Subroutines .. 00380 EXTERNAL DLATM1, DSCAL, XERBLA, ZLAGGE, ZLAGHE, ZLAGSY, 00381 $ ZLAROT, ZLARTG, ZLASET 00382 * .. 00383 * .. Intrinsic Functions .. 00384 INTRINSIC ABS, COS, DBLE, DCMPLX, DCONJG, MAX, MIN, MOD, 00385 $ SIN 00386 * .. 00387 * .. Executable Statements .. 00388 * 00389 * 1) Decode and Test the input parameters. 00390 * Initialize flags & seed. 00391 * 00392 INFO = 0 00393 * 00394 * Quick return if possible 00395 * 00396 IF( M.EQ.0 .OR. N.EQ.0 ) 00397 $ RETURN 00398 * 00399 * Decode DIST 00400 * 00401 IF( LSAME( DIST, 'U' ) ) THEN 00402 IDIST = 1 00403 ELSE IF( LSAME( DIST, 'S' ) ) THEN 00404 IDIST = 2 00405 ELSE IF( LSAME( DIST, 'N' ) ) THEN 00406 IDIST = 3 00407 ELSE 00408 IDIST = -1 00409 END IF 00410 * 00411 * Decode SYM 00412 * 00413 IF( LSAME( SYM, 'N' ) ) THEN 00414 ISYM = 1 00415 IRSIGN = 0 00416 ZSYM = .FALSE. 00417 ELSE IF( LSAME( SYM, 'P' ) ) THEN 00418 ISYM = 2 00419 IRSIGN = 0 00420 ZSYM = .FALSE. 00421 ELSE IF( LSAME( SYM, 'S' ) ) THEN 00422 ISYM = 2 00423 IRSIGN = 0 00424 ZSYM = .TRUE. 00425 ELSE IF( LSAME( SYM, 'H' ) ) THEN 00426 ISYM = 2 00427 IRSIGN = 1 00428 ZSYM = .FALSE. 00429 ELSE 00430 ISYM = -1 00431 END IF 00432 * 00433 * Decode PACK 00434 * 00435 ISYMPK = 0 00436 IF( LSAME( PACK, 'N' ) ) THEN 00437 IPACK = 0 00438 ELSE IF( LSAME( PACK, 'U' ) ) THEN 00439 IPACK = 1 00440 ISYMPK = 1 00441 ELSE IF( LSAME( PACK, 'L' ) ) THEN 00442 IPACK = 2 00443 ISYMPK = 1 00444 ELSE IF( LSAME( PACK, 'C' ) ) THEN 00445 IPACK = 3 00446 ISYMPK = 2 00447 ELSE IF( LSAME( PACK, 'R' ) ) THEN 00448 IPACK = 4 00449 ISYMPK = 3 00450 ELSE IF( LSAME( PACK, 'B' ) ) THEN 00451 IPACK = 5 00452 ISYMPK = 3 00453 ELSE IF( LSAME( PACK, 'Q' ) ) THEN 00454 IPACK = 6 00455 ISYMPK = 2 00456 ELSE IF( LSAME( PACK, 'Z' ) ) THEN 00457 IPACK = 7 00458 ELSE 00459 IPACK = -1 00460 END IF 00461 * 00462 * Set certain internal parameters 00463 * 00464 MNMIN = MIN( M, N ) 00465 LLB = MIN( KL, M-1 ) 00466 UUB = MIN( KU, N-1 ) 00467 MR = MIN( M, N+LLB ) 00468 NC = MIN( N, M+UUB ) 00469 * 00470 IF( IPACK.EQ.5 .OR. IPACK.EQ.6 ) THEN 00471 MINLDA = UUB + 1 00472 ELSE IF( IPACK.EQ.7 ) THEN 00473 MINLDA = LLB + UUB + 1 00474 ELSE 00475 MINLDA = M 00476 END IF 00477 * 00478 * Use Givens rotation method if bandwidth small enough, 00479 * or if LDA is too small to store the matrix unpacked. 00480 * 00481 GIVENS = .FALSE. 00482 IF( ISYM.EQ.1 ) THEN 00483 IF( DBLE( LLB+UUB ).LT.0.3D0*DBLE( MAX( 1, MR+NC ) ) ) 00484 $ GIVENS = .TRUE. 00485 ELSE 00486 IF( 2*LLB.LT.M ) 00487 $ GIVENS = .TRUE. 00488 END IF 00489 IF( LDA.LT.M .AND. LDA.GE.MINLDA ) 00490 $ GIVENS = .TRUE. 00491 * 00492 * Set INFO if an error 00493 * 00494 IF( M.LT.0 ) THEN 00495 INFO = -1 00496 ELSE IF( M.NE.N .AND. ISYM.NE.1 ) THEN 00497 INFO = -1 00498 ELSE IF( N.LT.0 ) THEN 00499 INFO = -2 00500 ELSE IF( IDIST.EQ.-1 ) THEN 00501 INFO = -3 00502 ELSE IF( ISYM.EQ.-1 ) THEN 00503 INFO = -5 00504 ELSE IF( ABS( MODE ).GT.6 ) THEN 00505 INFO = -7 00506 ELSE IF( ( MODE.NE.0 .AND. ABS( MODE ).NE.6 ) .AND. COND.LT.ONE ) 00507 $ THEN 00508 INFO = -8 00509 ELSE IF( KL.LT.0 ) THEN 00510 INFO = -10 00511 ELSE IF( KU.LT.0 .OR. ( ISYM.NE.1 .AND. KL.NE.KU ) ) THEN 00512 INFO = -11 00513 ELSE IF( IPACK.EQ.-1 .OR. ( ISYMPK.EQ.1 .AND. ISYM.EQ.1 ) .OR. 00514 $ ( ISYMPK.EQ.2 .AND. ISYM.EQ.1 .AND. KL.GT.0 ) .OR. 00515 $ ( ISYMPK.EQ.3 .AND. ISYM.EQ.1 .AND. KU.GT.0 ) .OR. 00516 $ ( ISYMPK.NE.0 .AND. M.NE.N ) ) THEN 00517 INFO = -12 00518 ELSE IF( LDA.LT.MAX( 1, MINLDA ) ) THEN 00519 INFO = -14 00520 END IF 00521 * 00522 IF( INFO.NE.0 ) THEN 00523 CALL XERBLA( 'ZLATMS', -INFO ) 00524 RETURN 00525 END IF 00526 * 00527 * Initialize random number generator 00528 * 00529 DO 10 I = 1, 4 00530 ISEED( I ) = MOD( ABS( ISEED( I ) ), 4096 ) 00531 10 CONTINUE 00532 * 00533 IF( MOD( ISEED( 4 ), 2 ).NE.1 ) 00534 $ ISEED( 4 ) = ISEED( 4 ) + 1 00535 * 00536 * 2) Set up D if indicated. 00537 * 00538 * Compute D according to COND and MODE 00539 * 00540 CALL DLATM1( MODE, COND, IRSIGN, IDIST, ISEED, D, MNMIN, IINFO ) 00541 IF( IINFO.NE.0 ) THEN 00542 INFO = 1 00543 RETURN 00544 END IF 00545 * 00546 * Choose Top-Down if D is (apparently) increasing, 00547 * Bottom-Up if D is (apparently) decreasing. 00548 * 00549 IF( ABS( D( 1 ) ).LE.ABS( D( MNMIN ) ) ) THEN 00550 TOPDWN = .TRUE. 00551 ELSE 00552 TOPDWN = .FALSE. 00553 END IF 00554 * 00555 IF( MODE.NE.0 .AND. ABS( MODE ).NE.6 ) THEN 00556 * 00557 * Scale by DMAX 00558 * 00559 TEMP = ABS( D( 1 ) ) 00560 DO 20 I = 2, MNMIN 00561 TEMP = MAX( TEMP, ABS( D( I ) ) ) 00562 20 CONTINUE 00563 * 00564 IF( TEMP.GT.ZERO ) THEN 00565 ALPHA = DMAX / TEMP 00566 ELSE 00567 INFO = 2 00568 RETURN 00569 END IF 00570 * 00571 CALL DSCAL( MNMIN, ALPHA, D, 1 ) 00572 * 00573 END IF 00574 * 00575 CALL ZLASET( 'Full', LDA, N, CZERO, CZERO, A, LDA ) 00576 * 00577 * 3) Generate Banded Matrix using Givens rotations. 00578 * Also the special case of UUB=LLB=0 00579 * 00580 * Compute Addressing constants to cover all 00581 * storage formats. Whether GE, HE, SY, GB, HB, or SB, 00582 * upper or lower triangle or both, 00583 * the (i,j)-th element is in 00584 * A( i - ISKEW*j + IOFFST, j ) 00585 * 00586 IF( IPACK.GT.4 ) THEN 00587 ILDA = LDA - 1 00588 ISKEW = 1 00589 IF( IPACK.GT.5 ) THEN 00590 IOFFST = UUB + 1 00591 ELSE 00592 IOFFST = 1 00593 END IF 00594 ELSE 00595 ILDA = LDA 00596 ISKEW = 0 00597 IOFFST = 0 00598 END IF 00599 * 00600 * IPACKG is the format that the matrix is generated in. If this is 00601 * different from IPACK, then the matrix must be repacked at the 00602 * end. It also signals how to compute the norm, for scaling. 00603 * 00604 IPACKG = 0 00605 * 00606 * Diagonal Matrix -- We are done, unless it 00607 * is to be stored HP/SP/PP/TP (PACK='R' or 'C') 00608 * 00609 IF( LLB.EQ.0 .AND. UUB.EQ.0 ) THEN 00610 DO 30 J = 1, MNMIN 00611 A( ( 1-ISKEW )*J+IOFFST, J ) = DCMPLX( D( J ) ) 00612 30 CONTINUE 00613 * 00614 IF( IPACK.LE.2 .OR. IPACK.GE.5 ) 00615 $ IPACKG = IPACK 00616 * 00617 ELSE IF( GIVENS ) THEN 00618 * 00619 * Check whether to use Givens rotations, 00620 * Householder transformations, or nothing. 00621 * 00622 IF( ISYM.EQ.1 ) THEN 00623 * 00624 * Non-symmetric -- A = U D V 00625 * 00626 IF( IPACK.GT.4 ) THEN 00627 IPACKG = IPACK 00628 ELSE 00629 IPACKG = 0 00630 END IF 00631 * 00632 DO 40 J = 1, MNMIN 00633 A( ( 1-ISKEW )*J+IOFFST, J ) = DCMPLX( D( J ) ) 00634 40 CONTINUE 00635 * 00636 IF( TOPDWN ) THEN 00637 JKL = 0 00638 DO 70 JKU = 1, UUB 00639 * 00640 * Transform from bandwidth JKL, JKU-1 to JKL, JKU 00641 * 00642 * Last row actually rotated is M 00643 * Last column actually rotated is MIN( M+JKU, N ) 00644 * 00645 DO 60 JR = 1, MIN( M+JKU, N ) + JKL - 1 00646 EXTRA = CZERO 00647 ANGLE = TWOPI*DLARND( 1, ISEED ) 00648 C = COS( ANGLE )*ZLARND( 5, ISEED ) 00649 S = SIN( ANGLE )*ZLARND( 5, ISEED ) 00650 ICOL = MAX( 1, JR-JKL ) 00651 IF( JR.LT.M ) THEN 00652 IL = MIN( N, JR+JKU ) + 1 - ICOL 00653 CALL ZLAROT( .TRUE., JR.GT.JKL, .FALSE., IL, C, 00654 $ S, A( JR-ISKEW*ICOL+IOFFST, ICOL ), 00655 $ ILDA, EXTRA, DUMMY ) 00656 END IF 00657 * 00658 * Chase "EXTRA" back up 00659 * 00660 IR = JR 00661 IC = ICOL 00662 DO 50 JCH = JR - JKL, 1, -JKL - JKU 00663 IF( IR.LT.M ) THEN 00664 CALL ZLARTG( A( IR+1-ISKEW*( IC+1 )+IOFFST, 00665 $ IC+1 ), EXTRA, REALC, S, DUMMY ) 00666 DUMMY = ZLARND( 5, ISEED ) 00667 C = DCONJG( REALC*DUMMY ) 00668 S = DCONJG( -S*DUMMY ) 00669 END IF 00670 IROW = MAX( 1, JCH-JKU ) 00671 IL = IR + 2 - IROW 00672 CTEMP = CZERO 00673 ILTEMP = JCH.GT.JKU 00674 CALL ZLAROT( .FALSE., ILTEMP, .TRUE., IL, C, S, 00675 $ A( IROW-ISKEW*IC+IOFFST, IC ), 00676 $ ILDA, CTEMP, EXTRA ) 00677 IF( ILTEMP ) THEN 00678 CALL ZLARTG( A( IROW+1-ISKEW*( IC+1 )+IOFFST, 00679 $ IC+1 ), CTEMP, REALC, S, DUMMY ) 00680 DUMMY = ZLARND( 5, ISEED ) 00681 C = DCONJG( REALC*DUMMY ) 00682 S = DCONJG( -S*DUMMY ) 00683 * 00684 ICOL = MAX( 1, JCH-JKU-JKL ) 00685 IL = IC + 2 - ICOL 00686 EXTRA = CZERO 00687 CALL ZLAROT( .TRUE., JCH.GT.JKU+JKL, .TRUE., 00688 $ IL, C, S, A( IROW-ISKEW*ICOL+ 00689 $ IOFFST, ICOL ), ILDA, EXTRA, 00690 $ CTEMP ) 00691 IC = ICOL 00692 IR = IROW 00693 END IF 00694 50 CONTINUE 00695 60 CONTINUE 00696 70 CONTINUE 00697 * 00698 JKU = UUB 00699 DO 100 JKL = 1, LLB 00700 * 00701 * Transform from bandwidth JKL-1, JKU to JKL, JKU 00702 * 00703 DO 90 JC = 1, MIN( N+JKL, M ) + JKU - 1 00704 EXTRA = CZERO 00705 ANGLE = TWOPI*DLARND( 1, ISEED ) 00706 C = COS( ANGLE )*ZLARND( 5, ISEED ) 00707 S = SIN( ANGLE )*ZLARND( 5, ISEED ) 00708 IROW = MAX( 1, JC-JKU ) 00709 IF( JC.LT.N ) THEN 00710 IL = MIN( M, JC+JKL ) + 1 - IROW 00711 CALL ZLAROT( .FALSE., JC.GT.JKU, .FALSE., IL, C, 00712 $ S, A( IROW-ISKEW*JC+IOFFST, JC ), 00713 $ ILDA, EXTRA, DUMMY ) 00714 END IF 00715 * 00716 * Chase "EXTRA" back up 00717 * 00718 IC = JC 00719 IR = IROW 00720 DO 80 JCH = JC - JKU, 1, -JKL - JKU 00721 IF( IC.LT.N ) THEN 00722 CALL ZLARTG( A( IR+1-ISKEW*( IC+1 )+IOFFST, 00723 $ IC+1 ), EXTRA, REALC, S, DUMMY ) 00724 DUMMY = ZLARND( 5, ISEED ) 00725 C = DCONJG( REALC*DUMMY ) 00726 S = DCONJG( -S*DUMMY ) 00727 END IF 00728 ICOL = MAX( 1, JCH-JKL ) 00729 IL = IC + 2 - ICOL 00730 CTEMP = CZERO 00731 ILTEMP = JCH.GT.JKL 00732 CALL ZLAROT( .TRUE., ILTEMP, .TRUE., IL, C, S, 00733 $ A( IR-ISKEW*ICOL+IOFFST, ICOL ), 00734 $ ILDA, CTEMP, EXTRA ) 00735 IF( ILTEMP ) THEN 00736 CALL ZLARTG( A( IR+1-ISKEW*( ICOL+1 )+IOFFST, 00737 $ ICOL+1 ), CTEMP, REALC, S, 00738 $ DUMMY ) 00739 DUMMY = ZLARND( 5, ISEED ) 00740 C = DCONJG( REALC*DUMMY ) 00741 S = DCONJG( -S*DUMMY ) 00742 IROW = MAX( 1, JCH-JKL-JKU ) 00743 IL = IR + 2 - IROW 00744 EXTRA = CZERO 00745 CALL ZLAROT( .FALSE., JCH.GT.JKL+JKU, .TRUE., 00746 $ IL, C, S, A( IROW-ISKEW*ICOL+ 00747 $ IOFFST, ICOL ), ILDA, EXTRA, 00748 $ CTEMP ) 00749 IC = ICOL 00750 IR = IROW 00751 END IF 00752 80 CONTINUE 00753 90 CONTINUE 00754 100 CONTINUE 00755 * 00756 ELSE 00757 * 00758 * Bottom-Up -- Start at the bottom right. 00759 * 00760 JKL = 0 00761 DO 130 JKU = 1, UUB 00762 * 00763 * Transform from bandwidth JKL, JKU-1 to JKL, JKU 00764 * 00765 * First row actually rotated is M 00766 * First column actually rotated is MIN( M+JKU, N ) 00767 * 00768 IENDCH = MIN( M, N+JKL ) - 1 00769 DO 120 JC = MIN( M+JKU, N ) - 1, 1 - JKL, -1 00770 EXTRA = CZERO 00771 ANGLE = TWOPI*DLARND( 1, ISEED ) 00772 C = COS( ANGLE )*ZLARND( 5, ISEED ) 00773 S = SIN( ANGLE )*ZLARND( 5, ISEED ) 00774 IROW = MAX( 1, JC-JKU+1 ) 00775 IF( JC.GT.0 ) THEN 00776 IL = MIN( M, JC+JKL+1 ) + 1 - IROW 00777 CALL ZLAROT( .FALSE., .FALSE., JC+JKL.LT.M, IL, 00778 $ C, S, A( IROW-ISKEW*JC+IOFFST, 00779 $ JC ), ILDA, DUMMY, EXTRA ) 00780 END IF 00781 * 00782 * Chase "EXTRA" back down 00783 * 00784 IC = JC 00785 DO 110 JCH = JC + JKL, IENDCH, JKL + JKU 00786 ILEXTR = IC.GT.0 00787 IF( ILEXTR ) THEN 00788 CALL ZLARTG( A( JCH-ISKEW*IC+IOFFST, IC ), 00789 $ EXTRA, REALC, S, DUMMY ) 00790 DUMMY = ZLARND( 5, ISEED ) 00791 C = REALC*DUMMY 00792 S = S*DUMMY 00793 END IF 00794 IC = MAX( 1, IC ) 00795 ICOL = MIN( N-1, JCH+JKU ) 00796 ILTEMP = JCH + JKU.LT.N 00797 CTEMP = CZERO 00798 CALL ZLAROT( .TRUE., ILEXTR, ILTEMP, ICOL+2-IC, 00799 $ C, S, A( JCH-ISKEW*IC+IOFFST, IC ), 00800 $ ILDA, EXTRA, CTEMP ) 00801 IF( ILTEMP ) THEN 00802 CALL ZLARTG( A( JCH-ISKEW*ICOL+IOFFST, 00803 $ ICOL ), CTEMP, REALC, S, DUMMY ) 00804 DUMMY = ZLARND( 5, ISEED ) 00805 C = REALC*DUMMY 00806 S = S*DUMMY 00807 IL = MIN( IENDCH, JCH+JKL+JKU ) + 2 - JCH 00808 EXTRA = CZERO 00809 CALL ZLAROT( .FALSE., .TRUE., 00810 $ JCH+JKL+JKU.LE.IENDCH, IL, C, S, 00811 $ A( JCH-ISKEW*ICOL+IOFFST, 00812 $ ICOL ), ILDA, CTEMP, EXTRA ) 00813 IC = ICOL 00814 END IF 00815 110 CONTINUE 00816 120 CONTINUE 00817 130 CONTINUE 00818 * 00819 JKU = UUB 00820 DO 160 JKL = 1, LLB 00821 * 00822 * Transform from bandwidth JKL-1, JKU to JKL, JKU 00823 * 00824 * First row actually rotated is MIN( N+JKL, M ) 00825 * First column actually rotated is N 00826 * 00827 IENDCH = MIN( N, M+JKU ) - 1 00828 DO 150 JR = MIN( N+JKL, M ) - 1, 1 - JKU, -1 00829 EXTRA = CZERO 00830 ANGLE = TWOPI*DLARND( 1, ISEED ) 00831 C = COS( ANGLE )*ZLARND( 5, ISEED ) 00832 S = SIN( ANGLE )*ZLARND( 5, ISEED ) 00833 ICOL = MAX( 1, JR-JKL+1 ) 00834 IF( JR.GT.0 ) THEN 00835 IL = MIN( N, JR+JKU+1 ) + 1 - ICOL 00836 CALL ZLAROT( .TRUE., .FALSE., JR+JKU.LT.N, IL, 00837 $ C, S, A( JR-ISKEW*ICOL+IOFFST, 00838 $ ICOL ), ILDA, DUMMY, EXTRA ) 00839 END IF 00840 * 00841 * Chase "EXTRA" back down 00842 * 00843 IR = JR 00844 DO 140 JCH = JR + JKU, IENDCH, JKL + JKU 00845 ILEXTR = IR.GT.0 00846 IF( ILEXTR ) THEN 00847 CALL ZLARTG( A( IR-ISKEW*JCH+IOFFST, JCH ), 00848 $ EXTRA, REALC, S, DUMMY ) 00849 DUMMY = ZLARND( 5, ISEED ) 00850 C = REALC*DUMMY 00851 S = S*DUMMY 00852 END IF 00853 IR = MAX( 1, IR ) 00854 IROW = MIN( M-1, JCH+JKL ) 00855 ILTEMP = JCH + JKL.LT.M 00856 CTEMP = CZERO 00857 CALL ZLAROT( .FALSE., ILEXTR, ILTEMP, IROW+2-IR, 00858 $ C, S, A( IR-ISKEW*JCH+IOFFST, 00859 $ JCH ), ILDA, EXTRA, CTEMP ) 00860 IF( ILTEMP ) THEN 00861 CALL ZLARTG( A( IROW-ISKEW*JCH+IOFFST, JCH ), 00862 $ CTEMP, REALC, S, DUMMY ) 00863 DUMMY = ZLARND( 5, ISEED ) 00864 C = REALC*DUMMY 00865 S = S*DUMMY 00866 IL = MIN( IENDCH, JCH+JKL+JKU ) + 2 - JCH 00867 EXTRA = CZERO 00868 CALL ZLAROT( .TRUE., .TRUE., 00869 $ JCH+JKL+JKU.LE.IENDCH, IL, C, S, 00870 $ A( IROW-ISKEW*JCH+IOFFST, JCH ), 00871 $ ILDA, CTEMP, EXTRA ) 00872 IR = IROW 00873 END IF 00874 140 CONTINUE 00875 150 CONTINUE 00876 160 CONTINUE 00877 * 00878 END IF 00879 * 00880 ELSE 00881 * 00882 * Symmetric -- A = U D U' 00883 * Hermitian -- A = U D U* 00884 * 00885 IPACKG = IPACK 00886 IOFFG = IOFFST 00887 * 00888 IF( TOPDWN ) THEN 00889 * 00890 * Top-Down -- Generate Upper triangle only 00891 * 00892 IF( IPACK.GE.5 ) THEN 00893 IPACKG = 6 00894 IOFFG = UUB + 1 00895 ELSE 00896 IPACKG = 1 00897 END IF 00898 * 00899 DO 170 J = 1, MNMIN 00900 A( ( 1-ISKEW )*J+IOFFG, J ) = DCMPLX( D( J ) ) 00901 170 CONTINUE 00902 * 00903 DO 200 K = 1, UUB 00904 DO 190 JC = 1, N - 1 00905 IROW = MAX( 1, JC-K ) 00906 IL = MIN( JC+1, K+2 ) 00907 EXTRA = CZERO 00908 CTEMP = A( JC-ISKEW*( JC+1 )+IOFFG, JC+1 ) 00909 ANGLE = TWOPI*DLARND( 1, ISEED ) 00910 C = COS( ANGLE )*ZLARND( 5, ISEED ) 00911 S = SIN( ANGLE )*ZLARND( 5, ISEED ) 00912 IF( ZSYM ) THEN 00913 CT = C 00914 ST = S 00915 ELSE 00916 CTEMP = DCONJG( CTEMP ) 00917 CT = DCONJG( C ) 00918 ST = DCONJG( S ) 00919 END IF 00920 CALL ZLAROT( .FALSE., JC.GT.K, .TRUE., IL, C, S, 00921 $ A( IROW-ISKEW*JC+IOFFG, JC ), ILDA, 00922 $ EXTRA, CTEMP ) 00923 CALL ZLAROT( .TRUE., .TRUE., .FALSE., 00924 $ MIN( K, N-JC )+1, CT, ST, 00925 $ A( ( 1-ISKEW )*JC+IOFFG, JC ), ILDA, 00926 $ CTEMP, DUMMY ) 00927 * 00928 * Chase EXTRA back up the matrix 00929 * 00930 ICOL = JC 00931 DO 180 JCH = JC - K, 1, -K 00932 CALL ZLARTG( A( JCH+1-ISKEW*( ICOL+1 )+IOFFG, 00933 $ ICOL+1 ), EXTRA, REALC, S, DUMMY ) 00934 DUMMY = ZLARND( 5, ISEED ) 00935 C = DCONJG( REALC*DUMMY ) 00936 S = DCONJG( -S*DUMMY ) 00937 CTEMP = A( JCH-ISKEW*( JCH+1 )+IOFFG, JCH+1 ) 00938 IF( ZSYM ) THEN 00939 CT = C 00940 ST = S 00941 ELSE 00942 CTEMP = DCONJG( CTEMP ) 00943 CT = DCONJG( C ) 00944 ST = DCONJG( S ) 00945 END IF 00946 CALL ZLAROT( .TRUE., .TRUE., .TRUE., K+2, C, S, 00947 $ A( ( 1-ISKEW )*JCH+IOFFG, JCH ), 00948 $ ILDA, CTEMP, EXTRA ) 00949 IROW = MAX( 1, JCH-K ) 00950 IL = MIN( JCH+1, K+2 ) 00951 EXTRA = CZERO 00952 CALL ZLAROT( .FALSE., JCH.GT.K, .TRUE., IL, CT, 00953 $ ST, A( IROW-ISKEW*JCH+IOFFG, JCH ), 00954 $ ILDA, EXTRA, CTEMP ) 00955 ICOL = JCH 00956 180 CONTINUE 00957 190 CONTINUE 00958 200 CONTINUE 00959 * 00960 * If we need lower triangle, copy from upper. Note that 00961 * the order of copying is chosen to work for 'q' -> 'b' 00962 * 00963 IF( IPACK.NE.IPACKG .AND. IPACK.NE.3 ) THEN 00964 DO 230 JC = 1, N 00965 IROW = IOFFST - ISKEW*JC 00966 IF( ZSYM ) THEN 00967 DO 210 JR = JC, MIN( N, JC+UUB ) 00968 A( JR+IROW, JC ) = A( JC-ISKEW*JR+IOFFG, JR ) 00969 210 CONTINUE 00970 ELSE 00971 DO 220 JR = JC, MIN( N, JC+UUB ) 00972 A( JR+IROW, JC ) = DCONJG( A( JC-ISKEW*JR+ 00973 $ IOFFG, JR ) ) 00974 220 CONTINUE 00975 END IF 00976 230 CONTINUE 00977 IF( IPACK.EQ.5 ) THEN 00978 DO 250 JC = N - UUB + 1, N 00979 DO 240 JR = N + 2 - JC, UUB + 1 00980 A( JR, JC ) = CZERO 00981 240 CONTINUE 00982 250 CONTINUE 00983 END IF 00984 IF( IPACKG.EQ.6 ) THEN 00985 IPACKG = IPACK 00986 ELSE 00987 IPACKG = 0 00988 END IF 00989 END IF 00990 ELSE 00991 * 00992 * Bottom-Up -- Generate Lower triangle only 00993 * 00994 IF( IPACK.GE.5 ) THEN 00995 IPACKG = 5 00996 IF( IPACK.EQ.6 ) 00997 $ IOFFG = 1 00998 ELSE 00999 IPACKG = 2 01000 END IF 01001 * 01002 DO 260 J = 1, MNMIN 01003 A( ( 1-ISKEW )*J+IOFFG, J ) = DCMPLX( D( J ) ) 01004 260 CONTINUE 01005 * 01006 DO 290 K = 1, UUB 01007 DO 280 JC = N - 1, 1, -1 01008 IL = MIN( N+1-JC, K+2 ) 01009 EXTRA = CZERO 01010 CTEMP = A( 1+( 1-ISKEW )*JC+IOFFG, JC ) 01011 ANGLE = TWOPI*DLARND( 1, ISEED ) 01012 C = COS( ANGLE )*ZLARND( 5, ISEED ) 01013 S = SIN( ANGLE )*ZLARND( 5, ISEED ) 01014 IF( ZSYM ) THEN 01015 CT = C 01016 ST = S 01017 ELSE 01018 CTEMP = DCONJG( CTEMP ) 01019 CT = DCONJG( C ) 01020 ST = DCONJG( S ) 01021 END IF 01022 CALL ZLAROT( .FALSE., .TRUE., N-JC.GT.K, IL, C, S, 01023 $ A( ( 1-ISKEW )*JC+IOFFG, JC ), ILDA, 01024 $ CTEMP, EXTRA ) 01025 ICOL = MAX( 1, JC-K+1 ) 01026 CALL ZLAROT( .TRUE., .FALSE., .TRUE., JC+2-ICOL, 01027 $ CT, ST, A( JC-ISKEW*ICOL+IOFFG, 01028 $ ICOL ), ILDA, DUMMY, CTEMP ) 01029 * 01030 * Chase EXTRA back down the matrix 01031 * 01032 ICOL = JC 01033 DO 270 JCH = JC + K, N - 1, K 01034 CALL ZLARTG( A( JCH-ISKEW*ICOL+IOFFG, ICOL ), 01035 $ EXTRA, REALC, S, DUMMY ) 01036 DUMMY = ZLARND( 5, ISEED ) 01037 C = REALC*DUMMY 01038 S = S*DUMMY 01039 CTEMP = A( 1+( 1-ISKEW )*JCH+IOFFG, JCH ) 01040 IF( ZSYM ) THEN 01041 CT = C 01042 ST = S 01043 ELSE 01044 CTEMP = DCONJG( CTEMP ) 01045 CT = DCONJG( C ) 01046 ST = DCONJG( S ) 01047 END IF 01048 CALL ZLAROT( .TRUE., .TRUE., .TRUE., K+2, C, S, 01049 $ A( JCH-ISKEW*ICOL+IOFFG, ICOL ), 01050 $ ILDA, EXTRA, CTEMP ) 01051 IL = MIN( N+1-JCH, K+2 ) 01052 EXTRA = CZERO 01053 CALL ZLAROT( .FALSE., .TRUE., N-JCH.GT.K, IL, 01054 $ CT, ST, A( ( 1-ISKEW )*JCH+IOFFG, 01055 $ JCH ), ILDA, CTEMP, EXTRA ) 01056 ICOL = JCH 01057 270 CONTINUE 01058 280 CONTINUE 01059 290 CONTINUE 01060 * 01061 * If we need upper triangle, copy from lower. Note that 01062 * the order of copying is chosen to work for 'b' -> 'q' 01063 * 01064 IF( IPACK.NE.IPACKG .AND. IPACK.NE.4 ) THEN 01065 DO 320 JC = N, 1, -1 01066 IROW = IOFFST - ISKEW*JC 01067 IF( ZSYM ) THEN 01068 DO 300 JR = JC, MAX( 1, JC-UUB ), -1 01069 A( JR+IROW, JC ) = A( JC-ISKEW*JR+IOFFG, JR ) 01070 300 CONTINUE 01071 ELSE 01072 DO 310 JR = JC, MAX( 1, JC-UUB ), -1 01073 A( JR+IROW, JC ) = DCONJG( A( JC-ISKEW*JR+ 01074 $ IOFFG, JR ) ) 01075 310 CONTINUE 01076 END IF 01077 320 CONTINUE 01078 IF( IPACK.EQ.6 ) THEN 01079 DO 340 JC = 1, UUB 01080 DO 330 JR = 1, UUB + 1 - JC 01081 A( JR, JC ) = CZERO 01082 330 CONTINUE 01083 340 CONTINUE 01084 END IF 01085 IF( IPACKG.EQ.5 ) THEN 01086 IPACKG = IPACK 01087 ELSE 01088 IPACKG = 0 01089 END IF 01090 END IF 01091 END IF 01092 * 01093 * Ensure that the diagonal is real if Hermitian 01094 * 01095 IF( .NOT.ZSYM ) THEN 01096 DO 350 JC = 1, N 01097 IROW = IOFFST + ( 1-ISKEW )*JC 01098 A( IROW, JC ) = DCMPLX( DBLE( A( IROW, JC ) ) ) 01099 350 CONTINUE 01100 END IF 01101 * 01102 END IF 01103 * 01104 ELSE 01105 * 01106 * 4) Generate Banded Matrix by first 01107 * Rotating by random Unitary matrices, 01108 * then reducing the bandwidth using Householder 01109 * transformations. 01110 * 01111 * Note: we should get here only if LDA .ge. N 01112 * 01113 IF( ISYM.EQ.1 ) THEN 01114 * 01115 * Non-symmetric -- A = U D V 01116 * 01117 CALL ZLAGGE( MR, NC, LLB, UUB, D, A, LDA, ISEED, WORK, 01118 $ IINFO ) 01119 ELSE 01120 * 01121 * Symmetric -- A = U D U' or 01122 * Hermitian -- A = U D U* 01123 * 01124 IF( ZSYM ) THEN 01125 CALL ZLAGSY( M, LLB, D, A, LDA, ISEED, WORK, IINFO ) 01126 ELSE 01127 CALL ZLAGHE( M, LLB, D, A, LDA, ISEED, WORK, IINFO ) 01128 END IF 01129 END IF 01130 * 01131 IF( IINFO.NE.0 ) THEN 01132 INFO = 3 01133 RETURN 01134 END IF 01135 END IF 01136 * 01137 * 5) Pack the matrix 01138 * 01139 IF( IPACK.NE.IPACKG ) THEN 01140 IF( IPACK.EQ.1 ) THEN 01141 * 01142 * 'U' -- Upper triangular, not packed 01143 * 01144 DO 370 J = 1, M 01145 DO 360 I = J + 1, M 01146 A( I, J ) = CZERO 01147 360 CONTINUE 01148 370 CONTINUE 01149 * 01150 ELSE IF( IPACK.EQ.2 ) THEN 01151 * 01152 * 'L' -- Lower triangular, not packed 01153 * 01154 DO 390 J = 2, M 01155 DO 380 I = 1, J - 1 01156 A( I, J ) = CZERO 01157 380 CONTINUE 01158 390 CONTINUE 01159 * 01160 ELSE IF( IPACK.EQ.3 ) THEN 01161 * 01162 * 'C' -- Upper triangle packed Columnwise. 01163 * 01164 ICOL = 1 01165 IROW = 0 01166 DO 410 J = 1, M 01167 DO 400 I = 1, J 01168 IROW = IROW + 1 01169 IF( IROW.GT.LDA ) THEN 01170 IROW = 1 01171 ICOL = ICOL + 1 01172 END IF 01173 A( IROW, ICOL ) = A( I, J ) 01174 400 CONTINUE 01175 410 CONTINUE 01176 * 01177 ELSE IF( IPACK.EQ.4 ) THEN 01178 * 01179 * 'R' -- Lower triangle packed Columnwise. 01180 * 01181 ICOL = 1 01182 IROW = 0 01183 DO 430 J = 1, M 01184 DO 420 I = J, M 01185 IROW = IROW + 1 01186 IF( IROW.GT.LDA ) THEN 01187 IROW = 1 01188 ICOL = ICOL + 1 01189 END IF 01190 A( IROW, ICOL ) = A( I, J ) 01191 420 CONTINUE 01192 430 CONTINUE 01193 * 01194 ELSE IF( IPACK.GE.5 ) THEN 01195 * 01196 * 'B' -- The lower triangle is packed as a band matrix. 01197 * 'Q' -- The upper triangle is packed as a band matrix. 01198 * 'Z' -- The whole matrix is packed as a band matrix. 01199 * 01200 IF( IPACK.EQ.5 ) 01201 $ UUB = 0 01202 IF( IPACK.EQ.6 ) 01203 $ LLB = 0 01204 * 01205 DO 450 J = 1, UUB 01206 DO 440 I = MIN( J+LLB, M ), 1, -1 01207 A( I-J+UUB+1, J ) = A( I, J ) 01208 440 CONTINUE 01209 450 CONTINUE 01210 * 01211 DO 470 J = UUB + 2, N 01212 DO 460 I = J - UUB, MIN( J+LLB, M ) 01213 A( I-J+UUB+1, J ) = A( I, J ) 01214 460 CONTINUE 01215 470 CONTINUE 01216 END IF 01217 * 01218 * If packed, zero out extraneous elements. 01219 * 01220 * Symmetric/Triangular Packed -- 01221 * zero out everything after A(IROW,ICOL) 01222 * 01223 IF( IPACK.EQ.3 .OR. IPACK.EQ.4 ) THEN 01224 DO 490 JC = ICOL, M 01225 DO 480 JR = IROW + 1, LDA 01226 A( JR, JC ) = CZERO 01227 480 CONTINUE 01228 IROW = 0 01229 490 CONTINUE 01230 * 01231 ELSE IF( IPACK.GE.5 ) THEN 01232 * 01233 * Packed Band -- 01234 * 1st row is now in A( UUB+2-j, j), zero above it 01235 * m-th row is now in A( M+UUB-j,j), zero below it 01236 * last non-zero diagonal is now in A( UUB+LLB+1,j ), 01237 * zero below it, too. 01238 * 01239 IR1 = UUB + LLB + 2 01240 IR2 = UUB + M + 2 01241 DO 520 JC = 1, N 01242 DO 500 JR = 1, UUB + 1 - JC 01243 A( JR, JC ) = CZERO 01244 500 CONTINUE 01245 DO 510 JR = MAX( 1, MIN( IR1, IR2-JC ) ), LDA 01246 A( JR, JC ) = CZERO 01247 510 CONTINUE 01248 520 CONTINUE 01249 END IF 01250 END IF 01251 * 01252 RETURN 01253 * 01254 * End of ZLATMS 01255 * 01256 END