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*> \date November 2011 294* 295*> \ingroup complex_matgen 296* 297* ===================================================================== 298 SUBROUTINE CLATME( N, DIST, ISEED, D, MODE, COND, DMAX, 299 $ RSIGN, 300 $ UPPER, SIM, DS, MODES, CONDS, KL, KU, ANORM, 301 $ A, 302 $ LDA, WORK, INFO ) 303* 304* -- LAPACK computational routine (version 3.4.0) -- 305* -- LAPACK is a software package provided by Univ. of Tennessee, -- 306* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 307* November 2011 308* 309* .. Scalar Arguments .. 310 CHARACTER DIST, RSIGN, SIM, UPPER 311 INTEGER INFO, KL, KU, LDA, MODE, MODES, N 312 REAL ANORM, COND, CONDS 313 COMPLEX DMAX 314* .. 315* .. Array Arguments .. 316 INTEGER ISEED( 4 ) 317 REAL DS( * ) 318 COMPLEX A( LDA, * ), D( * ), WORK( * ) 319* .. 320* 321* ===================================================================== 322* 323* .. Parameters .. 324 REAL ZERO 325 PARAMETER ( ZERO = 0.0E+0 ) 326 REAL ONE 327 PARAMETER ( ONE = 1.0E+0 ) 328 COMPLEX CZERO 329 PARAMETER ( CZERO = ( 0.0E+0, 0.0E+0 ) ) 330 COMPLEX CONE 331 PARAMETER ( CONE = ( 1.0E+0, 0.0E+0 ) ) 332* .. 333* .. Local Scalars .. 334 LOGICAL BADS 335 INTEGER I, IC, ICOLS, IDIST, IINFO, IR, IROWS, IRSIGN, 336 $ ISIM, IUPPER, J, JC, JCR 337 REAL RALPHA, TEMP 338 COMPLEX ALPHA, TAU, XNORMS 339* .. 340* .. Local Arrays .. 341 REAL TEMPA( 1 ) 342* .. 343* .. External Functions .. 344 LOGICAL LSAME 345 REAL CLANGE 346 COMPLEX CLARND 347 EXTERNAL LSAME, CLANGE, CLARND 348* .. 349* .. External Subroutines .. 350 EXTERNAL CCOPY, CGEMV, CGERC, CLACGV, CLARFG, CLARGE, 351 $ CLARNV, CLATM1, CLASET, CSCAL, CSSCAL, SLATM1, 352 $ XERBLA 353* .. 354* .. Intrinsic Functions .. 355 INTRINSIC ABS, CONJG, MAX, MOD 356* .. 357* .. Executable Statements .. 358* 359* 1) Decode and Test the input parameters. 360* Initialize flags & seed. 361* 362 INFO = 0 363* 364* Quick return if possible 365* 366 IF( N.EQ.0 ) 367 $ RETURN 368* 369* Decode DIST 370* 371 IF( LSAME( DIST, 'U' ) ) THEN 372 IDIST = 1 373 ELSE IF( LSAME( DIST, 'S' ) ) THEN 374 IDIST = 2 375 ELSE IF( LSAME( DIST, 'N' ) ) THEN 376 IDIST = 3 377 ELSE IF( LSAME( DIST, 'D' ) ) THEN 378 IDIST = 4 379 ELSE 380 IDIST = -1 381 END IF 382* 383* Decode RSIGN 384* 385 IF( LSAME( RSIGN, 'T' ) ) THEN 386 IRSIGN = 1 387 ELSE IF( LSAME( RSIGN, 'F' ) ) THEN 388 IRSIGN = 0 389 ELSE 390 IRSIGN = -1 391 END IF 392* 393* Decode UPPER 394* 395 IF( LSAME( UPPER, 'T' ) ) THEN 396 IUPPER = 1 397 ELSE IF( LSAME( UPPER, 'F' ) ) THEN 398 IUPPER = 0 399 ELSE 400 IUPPER = -1 401 END IF 402* 403* Decode SIM 404* 405 IF( LSAME( SIM, 'T' ) ) THEN 406 ISIM = 1 407 ELSE IF( LSAME( SIM, 'F' ) ) THEN 408 ISIM = 0 409 ELSE 410 ISIM = -1 411 END IF 412* 413* Check DS, if MODES=0 and ISIM=1 414* 415 BADS = .FALSE. 416 IF( MODES.EQ.0 .AND. ISIM.EQ.1 ) THEN 417 DO 10 J = 1, N 418 IF( DS( J ).EQ.ZERO ) 419 $ BADS = .TRUE. 420 10 CONTINUE 421 END IF 422* 423* Set INFO if an error 424* 425 IF( N.LT.0 ) THEN 426 INFO = -1 427 ELSE IF( IDIST.EQ.-1 ) THEN 428 INFO = -2 429 ELSE IF( ABS( MODE ).GT.6 ) THEN 430 INFO = -5 431 ELSE IF( ( MODE.NE.0 .AND. ABS( MODE ).NE.6 ) .AND. COND.LT.ONE ) 432 $ THEN 433 INFO = -6 434 ELSE IF( IRSIGN.EQ.-1 ) THEN 435 INFO = -9 436 ELSE IF( IUPPER.EQ.-1 ) THEN 437 INFO = -10 438 ELSE IF( ISIM.EQ.-1 ) THEN 439 INFO = -11 440 ELSE IF( BADS ) THEN 441 INFO = -12 442 ELSE IF( ISIM.EQ.1 .AND. ABS( MODES ).GT.5 ) THEN 443 INFO = -13 444 ELSE IF( ISIM.EQ.1 .AND. MODES.NE.0 .AND. CONDS.LT.ONE ) THEN 445 INFO = -14 446 ELSE IF( KL.LT.1 ) THEN 447 INFO = -15 448 ELSE IF( KU.LT.1 .OR. ( KU.LT.N-1 .AND. KL.LT.N-1 ) ) THEN 449 INFO = -16 450 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 451 INFO = -19 452 END IF 453* 454 IF( INFO.NE.0 ) THEN 455 CALL XERBLA( 'CLATME', -INFO ) 456 RETURN 457 END IF 458* 459* Initialize random number generator 460* 461 DO 20 I = 1, 4 462 ISEED( I ) = MOD( ABS( ISEED( I ) ), 4096 ) 463 20 CONTINUE 464* 465 IF( MOD( ISEED( 4 ), 2 ).NE.1 ) 466 $ ISEED( 4 ) = ISEED( 4 ) + 1 467* 468* 2) Set up diagonal of A 469* 470* Compute D according to COND and MODE 471* 472 CALL CLATM1( MODE, COND, IRSIGN, IDIST, ISEED, D, N, IINFO ) 473 IF( IINFO.NE.0 ) THEN 474 INFO = 1 475 RETURN 476 END IF 477 IF( MODE.NE.0 .AND. ABS( MODE ).NE.6 ) THEN 478* 479* Scale by DMAX 480* 481 TEMP = ABS( D( 1 ) ) 482 DO 30 I = 2, N 483 TEMP = MAX( TEMP, ABS( D( I ) ) ) 484 30 CONTINUE 485* 486 IF( TEMP.GT.ZERO ) THEN 487 ALPHA = DMAX / TEMP 488 ELSE 489 INFO = 2 490 RETURN 491 END IF 492* 493 CALL CSCAL( N, ALPHA, D, 1 ) 494* 495 END IF 496* 497 CALL CLASET( 'Full', N, N, CZERO, CZERO, A, LDA ) 498 CALL CCOPY( N, D, 1, A, LDA+1 ) 499* 500* 3) If UPPER='T', set upper triangle of A to random numbers. 501* 502 IF( IUPPER.NE.0 ) THEN 503 DO 40 JC = 2, N 504 CALL CLARNV( IDIST, ISEED, JC-1, A( 1, JC ) ) 505 40 CONTINUE 506 END IF 507* 508* 4) If SIM='T', apply similarity transformation. 509* 510* -1 511* Transform is X A X , where X = U S V, thus 512* 513* it is U S V A V' (1/S) U' 514* 515 IF( ISIM.NE.0 ) THEN 516* 517* Compute S (singular values of the eigenvector matrix) 518* according to CONDS and MODES 519* 520 CALL SLATM1( MODES, CONDS, 0, 0, ISEED, DS, N, IINFO ) 521 IF( IINFO.NE.0 ) THEN 522 INFO = 3 523 RETURN 524 END IF 525* 526* Multiply by V and V' 527* 528 CALL CLARGE( N, A, LDA, ISEED, WORK, IINFO ) 529 IF( IINFO.NE.0 ) THEN 530 INFO = 4 531 RETURN 532 END IF 533* 534* Multiply by S and (1/S) 535* 536 DO 50 J = 1, N 537 CALL CSSCAL( N, DS( J ), A( J, 1 ), LDA ) 538 IF( DS( J ).NE.ZERO ) THEN 539 CALL CSSCAL( N, ONE / DS( J ), A( 1, J ), 1 ) 540 ELSE 541 INFO = 5 542 RETURN 543 END IF 544 50 CONTINUE 545* 546* Multiply by U and U' 547* 548 CALL CLARGE( N, A, LDA, ISEED, WORK, IINFO ) 549 IF( IINFO.NE.0 ) THEN 550 INFO = 4 551 RETURN 552 END IF 553 END IF 554* 555* 5) Reduce the bandwidth. 556* 557 IF( KL.LT.N-1 ) THEN 558* 559* Reduce bandwidth -- kill column 560* 561 DO 60 JCR = KL + 1, N - 1 562 IC = JCR - KL 563 IROWS = N + 1 - JCR 564 ICOLS = N + KL - JCR 565* 566 CALL CCOPY( IROWS, A( JCR, IC ), 1, WORK, 1 ) 567 XNORMS = WORK( 1 ) 568 CALL CLARFG( IROWS, XNORMS, WORK( 2 ), 1, TAU ) 569 TAU = CONJG( TAU ) 570 WORK( 1 ) = CONE 571 ALPHA = CLARND( 5, ISEED ) 572* 573 CALL CGEMV( 'C', IROWS, ICOLS, CONE, A( JCR, IC+1 ), LDA, 574 $ WORK, 1, CZERO, WORK( IROWS+1 ), 1 ) 575 CALL CGERC( IROWS, ICOLS, -TAU, WORK, 1, WORK( IROWS+1 ), 1, 576 $ A( JCR, IC+1 ), LDA ) 577* 578 CALL CGEMV( 'N', N, IROWS, CONE, A( 1, JCR ), LDA, WORK, 1, 579 $ CZERO, WORK( IROWS+1 ), 1 ) 580 CALL CGERC( N, IROWS, -CONJG( TAU ), WORK( IROWS+1 ), 1, 581 $ WORK, 1, A( 1, JCR ), LDA ) 582* 583 A( JCR, IC ) = XNORMS 584 CALL CLASET( 'Full', IROWS-1, 1, CZERO, CZERO, 585 $ A( JCR+1, IC ), LDA ) 586* 587 CALL CSCAL( ICOLS+1, ALPHA, A( JCR, IC ), LDA ) 588 CALL CSCAL( N, CONJG( ALPHA ), A( 1, JCR ), 1 ) 589 60 CONTINUE 590 ELSE IF( KU.LT.N-1 ) THEN 591* 592* Reduce upper bandwidth -- kill a row at a time. 593* 594 DO 70 JCR = KU + 1, N - 1 595 IR = JCR - KU 596 IROWS = N + KU - JCR 597 ICOLS = N + 1 - JCR 598* 599 CALL CCOPY( ICOLS, A( IR, JCR ), LDA, WORK, 1 ) 600 XNORMS = WORK( 1 ) 601 CALL CLARFG( ICOLS, XNORMS, WORK( 2 ), 1, TAU ) 602 TAU = CONJG( TAU ) 603 WORK( 1 ) = CONE 604 CALL CLACGV( ICOLS-1, WORK( 2 ), 1 ) 605 ALPHA = CLARND( 5, ISEED ) 606* 607 CALL CGEMV( 'N', IROWS, ICOLS, CONE, A( IR+1, JCR ), LDA, 608 $ WORK, 1, CZERO, WORK( ICOLS+1 ), 1 ) 609 CALL CGERC( IROWS, ICOLS, -TAU, WORK( ICOLS+1 ), 1, WORK, 1, 610 $ A( IR+1, JCR ), LDA ) 611* 612 CALL CGEMV( 'C', ICOLS, N, CONE, A( JCR, 1 ), LDA, WORK, 1, 613 $ CZERO, WORK( ICOLS+1 ), 1 ) 614 CALL CGERC( ICOLS, N, -CONJG( TAU ), WORK, 1, 615 $ WORK( ICOLS+1 ), 1, A( JCR, 1 ), LDA ) 616* 617 A( IR, JCR ) = XNORMS 618 CALL CLASET( 'Full', 1, ICOLS-1, CZERO, CZERO, 619 $ A( IR, JCR+1 ), LDA ) 620* 621 CALL CSCAL( IROWS+1, ALPHA, A( IR, JCR ), 1 ) 622 CALL CSCAL( N, CONJG( ALPHA ), A( JCR, 1 ), LDA ) 623 70 CONTINUE 624 END IF 625* 626* Scale the matrix to have norm ANORM 627* 628 IF( ANORM.GE.ZERO ) THEN 629 TEMP = CLANGE( 'M', N, N, A, LDA, TEMPA ) 630 IF( TEMP.GT.ZERO ) THEN 631 RALPHA = ANORM / TEMP 632 DO 80 J = 1, N 633 CALL CSSCAL( N, RALPHA, A( 1, J ), 1 ) 634 80 CONTINUE 635 END IF 636 END IF 637* 638 RETURN 639* 640* End of CLATME 641* 642 END 643