1*> \brief \b CLATME 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 CLATME( N, DIST, ISEED, D, MODE, COND, DMAX, 12* RSIGN, 13* UPPER, SIM, DS, MODES, CONDS, KL, KU, ANORM, 14* A, 15* LDA, WORK, INFO ) 16* 17* .. Scalar Arguments .. 18* CHARACTER DIST, RSIGN, SIM, UPPER 19* INTEGER INFO, KL, KU, LDA, MODE, MODES, N 20* REAL ANORM, COND, CONDS 21* COMPLEX DMAX 22* .. 23* .. Array Arguments .. 24* INTEGER ISEED( 4 ) 25* REAL DS( * ) 26* COMPLEX A( LDA, * ), D( * ), WORK( * ) 27* .. 28* 29* 30*> \par Purpose: 31* ============= 32*> 33*> \verbatim 34*> 35*> CLATME generates random non-symmetric square matrices with 36*> specified eigenvalues for testing LAPACK programs. 37*> 38*> CLATME operates by applying the following sequence of 39*> operations: 40*> 41*> 1. Set the diagonal to D, where D may be input or 42*> computed according to MODE, COND, DMAX, and RSIGN 43*> as described below. 44*> 45*> 2. If UPPER='T', the upper triangle of A is set to random values 46*> out of distribution DIST. 47*> 48*> 3. If SIM='T', A is multiplied on the left by a random matrix 49*> X, whose singular values are specified by DS, MODES, and 50*> CONDS, and on the right by X inverse. 51*> 52*> 4. If KL < N-1, the lower bandwidth is reduced to KL using 53*> Householder transformations. If KU < N-1, the upper 54*> bandwidth is reduced to KU. 55*> 56*> 5. If ANORM is not negative, the matrix is scaled to have 57*> maximum-element-norm ANORM. 58*> 59*> (Note: since the matrix cannot be reduced beyond Hessenberg form, 60*> no packing options are available.) 61*> \endverbatim 62* 63* Arguments: 64* ========== 65* 66*> \param[in] N 67*> \verbatim 68*> N is INTEGER 69*> The number of columns (or rows) of A. Not modified. 70*> \endverbatim 71*> 72*> \param[in] DIST 73*> \verbatim 74*> DIST is CHARACTER*1 75*> On entry, DIST specifies the type of distribution to be used 76*> to generate the random eigen-/singular values, and on the 77*> upper triangle (see UPPER). 78*> 'U' => UNIFORM( 0, 1 ) ( 'U' for uniform ) 79*> 'S' => UNIFORM( -1, 1 ) ( 'S' for symmetric ) 80*> 'N' => NORMAL( 0, 1 ) ( 'N' for normal ) 81*> 'D' => uniform on the complex disc |z| < 1. 82*> Not modified. 83*> \endverbatim 84*> 85*> \param[in,out] ISEED 86*> \verbatim 87*> ISEED is INTEGER array, dimension ( 4 ) 88*> On entry ISEED specifies the seed of the random number 89*> generator. They should lie between 0 and 4095 inclusive, 90*> and ISEED(4) should be odd. The random number generator 91*> uses a linear congruential sequence limited to small 92*> integers, and so should produce machine independent 93*> random numbers. The values of ISEED are changed on 94*> exit, and can be used in the next call to CLATME 95*> to continue the same random number sequence. 96*> Changed on exit. 97*> \endverbatim 98*> 99*> \param[in,out] D 100*> \verbatim 101*> D is COMPLEX array, dimension ( N ) 102*> This array is used to specify the eigenvalues of A. If 103*> MODE=0, then D is assumed to contain the eigenvalues 104*> otherwise they will be computed according to MODE, COND, 105*> DMAX, and RSIGN and placed in D. 106*> Modified if MODE is nonzero. 107*> \endverbatim 108*> 109*> \param[in] MODE 110*> \verbatim 111*> MODE is INTEGER 112*> On entry this describes how the eigenvalues are to 113*> be specified: 114*> MODE = 0 means use D as input 115*> MODE = 1 sets D(1)=1 and D(2:N)=1.0/COND 116*> MODE = 2 sets D(1:N-1)=1 and D(N)=1.0/COND 117*> MODE = 3 sets D(I)=COND**(-(I-1)/(N-1)) 118*> MODE = 4 sets D(i)=1 - (i-1)/(N-1)*(1 - 1/COND) 119*> MODE = 5 sets D to random numbers in the range 120*> ( 1/COND , 1 ) such that their logarithms 121*> are uniformly distributed. 122*> MODE = 6 set D to random numbers from same distribution 123*> as the rest of the matrix. 124*> MODE < 0 has the same meaning as ABS(MODE), except that 125*> the order of the elements of D is reversed. 126*> Thus if MODE is between 1 and 4, D has entries ranging 127*> from 1 to 1/COND, if between -1 and -4, D has entries 128*> ranging from 1/COND to 1, 129*> Not modified. 130*> \endverbatim 131*> 132*> \param[in] COND 133*> \verbatim 134*> COND is REAL 135*> On entry, this is used as described under MODE above. 136*> If used, it must be >= 1. Not modified. 137*> \endverbatim 138*> 139*> \param[in] DMAX 140*> \verbatim 141*> DMAX is COMPLEX 142*> If MODE is neither -6, 0 nor 6, the contents of D, as 143*> computed according to MODE and COND, will be scaled by 144*> DMAX / max(abs(D(i))). Note that DMAX need not be 145*> positive or real: if DMAX is negative or complex (or zero), 146*> D will be scaled by a negative or complex number (or zero). 147*> If RSIGN='F' then the largest (absolute) eigenvalue will be 148*> equal to DMAX. 149*> Not modified. 150*> \endverbatim 151*> 152*> \param[in] RSIGN 153*> \verbatim 154*> RSIGN is CHARACTER*1 155*> If MODE is not 0, 6, or -6, and RSIGN='T', then the 156*> elements of D, as computed according to MODE and COND, will 157*> be multiplied by a random complex number from the unit 158*> circle |z| = 1. If RSIGN='F', they will not be. RSIGN may 159*> only have the values 'T' or 'F'. 160*> Not modified. 161*> \endverbatim 162*> 163*> \param[in] UPPER 164*> \verbatim 165*> UPPER is CHARACTER*1 166*> If UPPER='T', then the elements of A above the diagonal 167*> will be set to random numbers out of DIST. If UPPER='F', 168*> they will not. UPPER may only have the values 'T' or 'F'. 169*> Not modified. 170*> \endverbatim 171*> 172*> \param[in] SIM 173*> \verbatim 174*> SIM is CHARACTER*1 175*> If SIM='T', then A will be operated on by a "similarity 176*> transform", i.e., multiplied on the left by a matrix X and 177*> on the right by X inverse. X = U S V, where U and V are 178*> random unitary matrices and S is a (diagonal) matrix of 179*> singular values specified by DS, MODES, and CONDS. If 180*> SIM='F', then A will not be transformed. 181*> Not modified. 182*> \endverbatim 183*> 184*> \param[in,out] DS 185*> \verbatim 186*> DS is REAL array, dimension ( N ) 187*> This array is used to specify the singular values of X, 188*> in the same way that D specifies the eigenvalues of A. 189*> If MODE=0, the DS contains the singular values, which 190*> may not be zero. 191*> Modified if MODE is nonzero. 192*> \endverbatim 193*> 194*> \param[in] MODES 195*> \verbatim 196*> MODES is INTEGER 197*> \endverbatim 198*> 199*> \param[in] CONDS 200*> \verbatim 201*> CONDS is REAL 202*> Similar to MODE and COND, but for specifying the diagonal 203*> of S. MODES=-6 and +6 are not allowed (since they would 204*> result in randomly ill-conditioned eigenvalues.) 205*> \endverbatim 206*> 207*> \param[in] KL 208*> \verbatim 209*> KL is INTEGER 210*> This specifies the lower bandwidth of the matrix. KL=1 211*> specifies upper Hessenberg form. If KL is at least N-1, 212*> then A will have full lower bandwidth. 213*> Not modified. 214*> \endverbatim 215*> 216*> \param[in] KU 217*> \verbatim 218*> KU is INTEGER 219*> This specifies the upper bandwidth of the matrix. KU=1 220*> specifies lower Hessenberg form. If KU is at least N-1, 221*> then A will have full upper bandwidth; if KU and KL 222*> are both at least N-1, then A will be dense. Only one of 223*> KU and KL may be less than N-1. 224*> Not modified. 225*> \endverbatim 226*> 227*> \param[in] ANORM 228*> \verbatim 229*> ANORM is REAL 230*> If ANORM is not negative, then A will be scaled by a non- 231*> negative real number to make the maximum-element-norm of A 232*> to be ANORM. 233*> Not modified. 234*> \endverbatim 235*> 236*> \param[out] A 237*> \verbatim 238*> A is COMPLEX array, dimension ( LDA, N ) 239*> On exit A is the desired test matrix. 240*> Modified. 241*> \endverbatim 242*> 243*> \param[in] LDA 244*> \verbatim 245*> LDA is INTEGER 246*> LDA specifies the first dimension of A as declared in the 247*> calling program. LDA must be at least M. 248*> Not modified. 249*> \endverbatim 250*> 251*> \param[out] WORK 252*> \verbatim 253*> WORK is COMPLEX array, dimension ( 3*N ) 254*> Workspace. 255*> Modified. 256*> \endverbatim 257*> 258*> \param[out] INFO 259*> \verbatim 260*> INFO is INTEGER 261*> Error code. On exit, INFO will be set to one of the 262*> following values: 263*> 0 => normal return 264*> -1 => N negative 265*> -2 => DIST illegal string 266*> -5 => MODE not in range -6 to 6 267*> -6 => COND less than 1.0, and MODE neither -6, 0 nor 6 268*> -9 => RSIGN is not 'T' or 'F' 269*> -10 => UPPER is not 'T' or 'F' 270*> -11 => SIM is not 'T' or 'F' 271*> -12 => MODES=0 and DS has a zero singular value. 272*> -13 => MODES is not in the range -5 to 5. 273*> -14 => MODES is nonzero and CONDS is less than 1. 274*> -15 => KL is less than 1. 275*> -16 => KU is less than 1, or KL and KU are both less than 276*> N-1. 277*> -19 => LDA is less than M. 278*> 1 => Error return from CLATM1 (computing D) 279*> 2 => Cannot scale to DMAX (max. eigenvalue is 0) 280*> 3 => Error return from SLATM1 (computing DS) 281*> 4 => Error return from CLARGE 282*> 5 => Zero singular value from SLATM1. 283*> \endverbatim 284* 285* Authors: 286* ======== 287* 288*> \author Univ. of Tennessee 289*> \author Univ. of California Berkeley 290*> \author Univ. of Colorado Denver 291*> \author NAG Ltd. 292* 293*> \ingroup complex_matgen 294* 295* ===================================================================== 296 SUBROUTINE CLATME( N, DIST, ISEED, D, MODE, COND, DMAX, 297 $ RSIGN, 298 $ UPPER, SIM, DS, MODES, CONDS, KL, KU, ANORM, 299 $ A, 300 $ LDA, WORK, INFO ) 301* 302* -- LAPACK computational routine -- 303* -- LAPACK is a software package provided by Univ. of Tennessee, -- 304* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 305* 306* .. Scalar Arguments .. 307 CHARACTER DIST, RSIGN, SIM, UPPER 308 INTEGER INFO, KL, KU, LDA, MODE, MODES, N 309 REAL ANORM, COND, CONDS 310 COMPLEX DMAX 311* .. 312* .. Array Arguments .. 313 INTEGER ISEED( 4 ) 314 REAL DS( * ) 315 COMPLEX A( LDA, * ), D( * ), WORK( * ) 316* .. 317* 318* ===================================================================== 319* 320* .. Parameters .. 321 REAL ZERO 322 PARAMETER ( ZERO = 0.0E+0 ) 323 REAL ONE 324 PARAMETER ( ONE = 1.0E+0 ) 325 COMPLEX CZERO 326 PARAMETER ( CZERO = ( 0.0E+0, 0.0E+0 ) ) 327 COMPLEX CONE 328 PARAMETER ( CONE = ( 1.0E+0, 0.0E+0 ) ) 329* .. 330* .. Local Scalars .. 331 LOGICAL BADS 332 INTEGER I, IC, ICOLS, IDIST, IINFO, IR, IROWS, IRSIGN, 333 $ ISIM, IUPPER, J, JC, JCR 334 REAL RALPHA, TEMP 335 COMPLEX ALPHA, TAU, XNORMS 336* .. 337* .. Local Arrays .. 338 REAL TEMPA( 1 ) 339* .. 340* .. External Functions .. 341 LOGICAL LSAME 342 REAL CLANGE 343 COMPLEX CLARND 344 EXTERNAL LSAME, CLANGE, CLARND 345* .. 346* .. External Subroutines .. 347 EXTERNAL CCOPY, CGEMV, CGERC, CLACGV, CLARFG, CLARGE, 348 $ CLARNV, CLATM1, CLASET, CSCAL, CSSCAL, SLATM1, 349 $ XERBLA 350* .. 351* .. Intrinsic Functions .. 352 INTRINSIC ABS, CONJG, MAX, MOD 353* .. 354* .. Executable Statements .. 355* 356* 1) Decode and Test the input parameters. 357* Initialize flags & seed. 358* 359 INFO = 0 360* 361* Quick return if possible 362* 363 IF( N.EQ.0 ) 364 $ RETURN 365* 366* Decode DIST 367* 368 IF( LSAME( DIST, 'U' ) ) THEN 369 IDIST = 1 370 ELSE IF( LSAME( DIST, 'S' ) ) THEN 371 IDIST = 2 372 ELSE IF( LSAME( DIST, 'N' ) ) THEN 373 IDIST = 3 374 ELSE IF( LSAME( DIST, 'D' ) ) THEN 375 IDIST = 4 376 ELSE 377 IDIST = -1 378 END IF 379* 380* Decode RSIGN 381* 382 IF( LSAME( RSIGN, 'T' ) ) THEN 383 IRSIGN = 1 384 ELSE IF( LSAME( RSIGN, 'F' ) ) THEN 385 IRSIGN = 0 386 ELSE 387 IRSIGN = -1 388 END IF 389* 390* Decode UPPER 391* 392 IF( LSAME( UPPER, 'T' ) ) THEN 393 IUPPER = 1 394 ELSE IF( LSAME( UPPER, 'F' ) ) THEN 395 IUPPER = 0 396 ELSE 397 IUPPER = -1 398 END IF 399* 400* Decode SIM 401* 402 IF( LSAME( SIM, 'T' ) ) THEN 403 ISIM = 1 404 ELSE IF( LSAME( SIM, 'F' ) ) THEN 405 ISIM = 0 406 ELSE 407 ISIM = -1 408 END IF 409* 410* Check DS, if MODES=0 and ISIM=1 411* 412 BADS = .FALSE. 413 IF( MODES.EQ.0 .AND. ISIM.EQ.1 ) THEN 414 DO 10 J = 1, N 415 IF( DS( J ).EQ.ZERO ) 416 $ BADS = .TRUE. 417 10 CONTINUE 418 END IF 419* 420* Set INFO if an error 421* 422 IF( N.LT.0 ) THEN 423 INFO = -1 424 ELSE IF( IDIST.EQ.-1 ) THEN 425 INFO = -2 426 ELSE IF( ABS( MODE ).GT.6 ) THEN 427 INFO = -5 428 ELSE IF( ( MODE.NE.0 .AND. ABS( MODE ).NE.6 ) .AND. COND.LT.ONE ) 429 $ THEN 430 INFO = -6 431 ELSE IF( IRSIGN.EQ.-1 ) THEN 432 INFO = -9 433 ELSE IF( IUPPER.EQ.-1 ) THEN 434 INFO = -10 435 ELSE IF( ISIM.EQ.-1 ) THEN 436 INFO = -11 437 ELSE IF( BADS ) THEN 438 INFO = -12 439 ELSE IF( ISIM.EQ.1 .AND. ABS( MODES ).GT.5 ) THEN 440 INFO = -13 441 ELSE IF( ISIM.EQ.1 .AND. MODES.NE.0 .AND. CONDS.LT.ONE ) THEN 442 INFO = -14 443 ELSE IF( KL.LT.1 ) THEN 444 INFO = -15 445 ELSE IF( KU.LT.1 .OR. ( KU.LT.N-1 .AND. KL.LT.N-1 ) ) THEN 446 INFO = -16 447 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 448 INFO = -19 449 END IF 450* 451 IF( INFO.NE.0 ) THEN 452 CALL XERBLA( 'CLATME', -INFO ) 453 RETURN 454 END IF 455* 456* Initialize random number generator 457* 458 DO 20 I = 1, 4 459 ISEED( I ) = MOD( ABS( ISEED( I ) ), 4096 ) 460 20 CONTINUE 461* 462 IF( MOD( ISEED( 4 ), 2 ).NE.1 ) 463 $ ISEED( 4 ) = ISEED( 4 ) + 1 464* 465* 2) Set up diagonal of A 466* 467* Compute D according to COND and MODE 468* 469 CALL CLATM1( MODE, COND, IRSIGN, IDIST, ISEED, D, N, IINFO ) 470 IF( IINFO.NE.0 ) THEN 471 INFO = 1 472 RETURN 473 END IF 474 IF( MODE.NE.0 .AND. ABS( MODE ).NE.6 ) THEN 475* 476* Scale by DMAX 477* 478 TEMP = ABS( D( 1 ) ) 479 DO 30 I = 2, N 480 TEMP = MAX( TEMP, ABS( D( I ) ) ) 481 30 CONTINUE 482* 483 IF( TEMP.GT.ZERO ) THEN 484 ALPHA = DMAX / TEMP 485 ELSE 486 INFO = 2 487 RETURN 488 END IF 489* 490 CALL CSCAL( N, ALPHA, D, 1 ) 491* 492 END IF 493* 494 CALL CLASET( 'Full', N, N, CZERO, CZERO, A, LDA ) 495 CALL CCOPY( N, D, 1, A, LDA+1 ) 496* 497* 3) If UPPER='T', set upper triangle of A to random numbers. 498* 499 IF( IUPPER.NE.0 ) THEN 500 DO 40 JC = 2, N 501 CALL CLARNV( IDIST, ISEED, JC-1, A( 1, JC ) ) 502 40 CONTINUE 503 END IF 504* 505* 4) If SIM='T', apply similarity transformation. 506* 507* -1 508* Transform is X A X , where X = U S V, thus 509* 510* it is U S V A V' (1/S) U' 511* 512 IF( ISIM.NE.0 ) THEN 513* 514* Compute S (singular values of the eigenvector matrix) 515* according to CONDS and MODES 516* 517 CALL SLATM1( MODES, CONDS, 0, 0, ISEED, DS, N, IINFO ) 518 IF( IINFO.NE.0 ) THEN 519 INFO = 3 520 RETURN 521 END IF 522* 523* Multiply by V and V' 524* 525 CALL CLARGE( N, A, LDA, ISEED, WORK, IINFO ) 526 IF( IINFO.NE.0 ) THEN 527 INFO = 4 528 RETURN 529 END IF 530* 531* Multiply by S and (1/S) 532* 533 DO 50 J = 1, N 534 CALL CSSCAL( N, DS( J ), A( J, 1 ), LDA ) 535 IF( DS( J ).NE.ZERO ) THEN 536 CALL CSSCAL( N, ONE / DS( J ), A( 1, J ), 1 ) 537 ELSE 538 INFO = 5 539 RETURN 540 END IF 541 50 CONTINUE 542* 543* Multiply by U and U' 544* 545 CALL CLARGE( N, A, LDA, ISEED, WORK, IINFO ) 546 IF( IINFO.NE.0 ) THEN 547 INFO = 4 548 RETURN 549 END IF 550 END IF 551* 552* 5) Reduce the bandwidth. 553* 554 IF( KL.LT.N-1 ) THEN 555* 556* Reduce bandwidth -- kill column 557* 558 DO 60 JCR = KL + 1, N - 1 559 IC = JCR - KL 560 IROWS = N + 1 - JCR 561 ICOLS = N + KL - JCR 562* 563 CALL CCOPY( IROWS, A( JCR, IC ), 1, WORK, 1 ) 564 XNORMS = WORK( 1 ) 565 CALL CLARFG( IROWS, XNORMS, WORK( 2 ), 1, TAU ) 566 TAU = CONJG( TAU ) 567 WORK( 1 ) = CONE 568 ALPHA = CLARND( 5, ISEED ) 569* 570 CALL CGEMV( 'C', IROWS, ICOLS, CONE, A( JCR, IC+1 ), LDA, 571 $ WORK, 1, CZERO, WORK( IROWS+1 ), 1 ) 572 CALL CGERC( IROWS, ICOLS, -TAU, WORK, 1, WORK( IROWS+1 ), 1, 573 $ A( JCR, IC+1 ), LDA ) 574* 575 CALL CGEMV( 'N', N, IROWS, CONE, A( 1, JCR ), LDA, WORK, 1, 576 $ CZERO, WORK( IROWS+1 ), 1 ) 577 CALL CGERC( N, IROWS, -CONJG( TAU ), WORK( IROWS+1 ), 1, 578 $ WORK, 1, A( 1, JCR ), LDA ) 579* 580 A( JCR, IC ) = XNORMS 581 CALL CLASET( 'Full', IROWS-1, 1, CZERO, CZERO, 582 $ A( JCR+1, IC ), LDA ) 583* 584 CALL CSCAL( ICOLS+1, ALPHA, A( JCR, IC ), LDA ) 585 CALL CSCAL( N, CONJG( ALPHA ), A( 1, JCR ), 1 ) 586 60 CONTINUE 587 ELSE IF( KU.LT.N-1 ) THEN 588* 589* Reduce upper bandwidth -- kill a row at a time. 590* 591 DO 70 JCR = KU + 1, N - 1 592 IR = JCR - KU 593 IROWS = N + KU - JCR 594 ICOLS = N + 1 - JCR 595* 596 CALL CCOPY( ICOLS, A( IR, JCR ), LDA, WORK, 1 ) 597 XNORMS = WORK( 1 ) 598 CALL CLARFG( ICOLS, XNORMS, WORK( 2 ), 1, TAU ) 599 TAU = CONJG( TAU ) 600 WORK( 1 ) = CONE 601 CALL CLACGV( ICOLS-1, WORK( 2 ), 1 ) 602 ALPHA = CLARND( 5, ISEED ) 603* 604 CALL CGEMV( 'N', IROWS, ICOLS, CONE, A( IR+1, JCR ), LDA, 605 $ WORK, 1, CZERO, WORK( ICOLS+1 ), 1 ) 606 CALL CGERC( IROWS, ICOLS, -TAU, WORK( ICOLS+1 ), 1, WORK, 1, 607 $ A( IR+1, JCR ), LDA ) 608* 609 CALL CGEMV( 'C', ICOLS, N, CONE, A( JCR, 1 ), LDA, WORK, 1, 610 $ CZERO, WORK( ICOLS+1 ), 1 ) 611 CALL CGERC( ICOLS, N, -CONJG( TAU ), WORK, 1, 612 $ WORK( ICOLS+1 ), 1, A( JCR, 1 ), LDA ) 613* 614 A( IR, JCR ) = XNORMS 615 CALL CLASET( 'Full', 1, ICOLS-1, CZERO, CZERO, 616 $ A( IR, JCR+1 ), LDA ) 617* 618 CALL CSCAL( IROWS+1, ALPHA, A( IR, JCR ), 1 ) 619 CALL CSCAL( N, CONJG( ALPHA ), A( JCR, 1 ), LDA ) 620 70 CONTINUE 621 END IF 622* 623* Scale the matrix to have norm ANORM 624* 625 IF( ANORM.GE.ZERO ) THEN 626 TEMP = CLANGE( 'M', N, N, A, LDA, TEMPA ) 627 IF( TEMP.GT.ZERO ) THEN 628 RALPHA = ANORM / TEMP 629 DO 80 J = 1, N 630 CALL CSSCAL( N, RALPHA, A( 1, J ), 1 ) 631 80 CONTINUE 632 END IF 633 END IF 634* 635 RETURN 636* 637* End of CLATME 638* 639 END 640