LAPACK  3.10.1 LAPACK: Linear Algebra PACKage

## ◆ slatms()

 subroutine slatms ( integer M, integer N, character DIST, integer, dimension( 4 ) ISEED, character SYM, real, dimension( * ) D, integer MODE, real COND, real DMAX, integer KL, integer KU, character PACK, real, dimension( lda, * ) A, integer LDA, real, dimension( * ) WORK, integer INFO )

SLATMS

Purpose:
```    SLATMS generates random matrices with specified singular values
(or symmetric/hermitian with specified eigenvalues)
for testing LAPACK programs.

SLATMS operates by applying the following sequence of
operations:

Set the diagonal to D, where D may be input or
computed according to MODE, COND, DMAX, and SYM
as described below.

Generate a matrix with the appropriate band structure, by one
of two methods:

Method A:
Generate a dense M x N matrix by multiplying D on the left
and the right by random unitary matrices, then:

Reduce the bandwidth according to KL and KU, using
Householder transformations.

Method B:
Convert the bandwidth-0 (i.e., diagonal) matrix to a
bandwidth-1 matrix using Givens rotations, "chasing"
out-of-band elements back, much as in QR; then
convert the bandwidth-1 to a bandwidth-2 matrix, etc.
Note that for reasonably small bandwidths (relative to
M and N) this requires less storage, as a dense matrix
is not generated.  Also, for symmetric matrices, only
one triangle is generated.

Method A is chosen if the bandwidth is a large fraction of the
order of the matrix, and LDA is at least M (so a dense
matrix can be stored.)  Method B is chosen if the bandwidth
is small (< 1/2 N for symmetric, < .3 N+M for
non-symmetric), or LDA is less than M and not less than the
bandwidth.

Pack the matrix if desired. Options specified by PACK are:
no packing
zero out upper half (if symmetric)
zero out lower half (if symmetric)
store the upper half columnwise (if symmetric or upper
triangular)
store the lower half columnwise (if symmetric or lower
triangular)
store the lower triangle in banded format (if symmetric
or lower triangular)
store the upper triangle in banded format (if symmetric
or upper triangular)
store the entire matrix in banded format
If Method B is chosen, and band format is specified, then the
matrix will be generated in the band format, so no repacking
will be necessary.```
Parameters
 [in] M ``` M is INTEGER The number of rows of A. Not modified.``` [in] N ``` N is INTEGER The number of columns of A. Not modified.``` [in] DIST ``` DIST is CHARACTER*1 On entry, DIST specifies the type of distribution to be used to generate the random eigen-/singular values. 'U' => UNIFORM( 0, 1 ) ( 'U' for uniform ) 'S' => UNIFORM( -1, 1 ) ( 'S' for symmetric ) 'N' => NORMAL( 0, 1 ) ( 'N' for normal ) Not modified.``` [in,out] ISEED ``` ISEED is INTEGER array, dimension ( 4 ) On entry ISEED specifies the seed of the random number generator. They should lie between 0 and 4095 inclusive, and ISEED(4) should be odd. The random number generator uses a linear congruential sequence limited to small integers, and so should produce machine independent random numbers. The values of ISEED are changed on exit, and can be used in the next call to SLATMS to continue the same random number sequence. Changed on exit.``` [in] SYM ``` SYM is CHARACTER*1 If SYM='S' or 'H', the generated matrix is symmetric, with eigenvalues specified by D, COND, MODE, and DMAX; they may be positive, negative, or zero. If SYM='P', the generated matrix is symmetric, with eigenvalues (= singular values) specified by D, COND, MODE, and DMAX; they will not be negative. If SYM='N', the generated matrix is nonsymmetric, with singular values specified by D, COND, MODE, and DMAX; they will not be negative. Not modified.``` [in,out] D ``` D is REAL array, dimension ( MIN( M , N ) ) This array is used to specify the singular values or eigenvalues of A (see SYM, above.) If MODE=0, then D is assumed to contain the singular/eigenvalues, otherwise they will be computed according to MODE, COND, and DMAX, and placed in D. Modified if MODE is nonzero.``` [in] MODE ``` MODE is INTEGER On entry this describes how the singular/eigenvalues are to be specified: MODE = 0 means use D as input MODE = 1 sets D(1)=1 and D(2:N)=1.0/COND MODE = 2 sets D(1:N-1)=1 and D(N)=1.0/COND MODE = 3 sets D(I)=COND**(-(I-1)/(N-1)) MODE = 4 sets D(i)=1 - (i-1)/(N-1)*(1 - 1/COND) MODE = 5 sets D to random numbers in the range ( 1/COND , 1 ) such that their logarithms are uniformly distributed. MODE = 6 set D to random numbers from same distribution as the rest of the matrix. MODE < 0 has the same meaning as ABS(MODE), except that the order of the elements of D is reversed. Thus if MODE is positive, D has entries ranging from 1 to 1/COND, if negative, from 1/COND to 1, If SYM='S' or 'H', and MODE is neither 0, 6, nor -6, then the elements of D will also be multiplied by a random sign (i.e., +1 or -1.) Not modified.``` [in] COND ``` COND is REAL On entry, this is used as described under MODE above. If used, it must be >= 1. Not modified.``` [in] DMAX ``` DMAX is REAL If MODE is neither -6, 0 nor 6, the contents of D, as computed according to MODE and COND, will be scaled by DMAX / max(abs(D(i))); thus, the maximum absolute eigen- or singular value (which is to say the norm) will be abs(DMAX). Note that DMAX need not be positive: if DMAX is negative (or zero), D will be scaled by a negative number (or zero). Not modified.``` [in] KL ``` KL is INTEGER This specifies the lower bandwidth of the matrix. For example, KL=0 implies upper triangular, KL=1 implies upper Hessenberg, and KL being at least M-1 means that the matrix has full lower bandwidth. KL must equal KU if the matrix is symmetric. Not modified.``` [in] KU ``` KU is INTEGER This specifies the upper bandwidth of the matrix. For example, KU=0 implies lower triangular, KU=1 implies lower Hessenberg, and KU being at least N-1 means that the matrix has full upper bandwidth. KL must equal KU if the matrix is symmetric. Not modified.``` [in] PACK ``` PACK is CHARACTER*1 This specifies packing of matrix as follows: 'N' => no packing 'U' => zero out all subdiagonal entries (if symmetric) 'L' => zero out all superdiagonal entries (if symmetric) 'C' => store the upper triangle columnwise (only if the matrix is symmetric or upper triangular) 'R' => store the lower triangle columnwise (only if the matrix is symmetric or lower triangular) 'B' => store the lower triangle in band storage scheme (only if matrix symmetric or lower triangular) 'Q' => store the upper triangle in band storage scheme (only if matrix symmetric or upper triangular) 'Z' => store the entire matrix in band storage scheme (pivoting can be provided for by using this option to store A in the trailing rows of the allocated storage) Using these options, the various LAPACK packed and banded storage schemes can be obtained: GB - use 'Z' PB, SB or TB - use 'B' or 'Q' PP, SP or TP - use 'C' or 'R' If two calls to SLATMS differ only in the PACK parameter, they will generate mathematically equivalent matrices. Not modified.``` [in,out] A ``` A is REAL array, dimension ( LDA, N ) On exit A is the desired test matrix. A is first generated in full (unpacked) form, and then packed, if so specified by PACK. Thus, the first M elements of the first N columns will always be modified. If PACK specifies a packed or banded storage scheme, all LDA elements of the first N columns will be modified; the elements of the array which do not correspond to elements of the generated matrix are set to zero. Modified.``` [in] LDA ``` LDA is INTEGER LDA specifies the first dimension of A as declared in the calling program. If PACK='N', 'U', 'L', 'C', or 'R', then LDA must be at least M. If PACK='B' or 'Q', then LDA must be at least MIN( KL, M-1) (which is equal to MIN(KU,N-1)). If PACK='Z', LDA must be large enough to hold the packed array: MIN( KU, N-1) + MIN( KL, M-1) + 1. Not modified.``` [out] WORK ``` WORK is REAL array, dimension ( 3*MAX( N , M ) ) Workspace. Modified.``` [out] INFO ``` INFO is INTEGER Error code. On exit, INFO will be set to one of the following values: 0 => normal return -1 => M negative or unequal to N and SYM='S', 'H', or 'P' -2 => N negative -3 => DIST illegal string -5 => SYM illegal string -7 => MODE not in range -6 to 6 -8 => COND less than 1.0, and MODE neither -6, 0 nor 6 -10 => KL negative -11 => KU negative, or SYM='S' or 'H' and KU not equal to KL -12 => PACK illegal string, or PACK='U' or 'L', and SYM='N'; or PACK='C' or 'Q' and SYM='N' and KL is not zero; or PACK='R' or 'B' and SYM='N' and KU is not zero; or PACK='U', 'L', 'C', 'R', 'B', or 'Q', and M is not N. -14 => LDA is less than M, or PACK='Z' and LDA is less than MIN(KU,N-1) + MIN(KL,M-1) + 1. 1 => Error return from SLATM1 2 => Cannot scale to DMAX (max. sing. value is 0) 3 => Error return from SLAGGE or SLAGSY```

Definition at line 319 of file slatms.f.

321 *
322 * -- LAPACK computational routine --
323 * -- LAPACK is a software package provided by Univ. of Tennessee, --
324 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
325 *
326 * .. Scalar Arguments ..
327  CHARACTER DIST, PACK, SYM
328  INTEGER INFO, KL, KU, LDA, M, MODE, N
329  REAL COND, DMAX
330 * ..
331 * .. Array Arguments ..
332  INTEGER ISEED( 4 )
333  REAL A( LDA, * ), D( * ), WORK( * )
334 * ..
335 *
336 * =====================================================================
337 *
338 * .. Parameters ..
339  REAL ZERO
340  parameter( zero = 0.0e0 )
341  REAL ONE
342  parameter( one = 1.0e0 )
343  REAL TWOPI
344  parameter( twopi = 6.28318530717958647692528676655900576839e+0 )
345 * ..
346 * .. Local Scalars ..
347  LOGICAL GIVENS, ILEXTR, ILTEMP, TOPDWN
348  INTEGER I, IC, ICOL, IDIST, IENDCH, IINFO, IL, ILDA,
349  \$ IOFFG, IOFFST, IPACK, IPACKG, IR, IR1, IR2,
350  \$ IROW, IRSIGN, ISKEW, ISYM, ISYMPK, J, JC, JCH,
351  \$ JKL, JKU, JR, K, LLB, MINLDA, MNMIN, MR, NC,
352  \$ UUB
353  REAL ALPHA, ANGLE, C, DUMMY, EXTRA, S, TEMP
354 * ..
355 * .. External Functions ..
356  LOGICAL LSAME
357  REAL SLARND
358  EXTERNAL lsame, slarnd
359 * ..
360 * .. External Subroutines ..
361  EXTERNAL scopy, slagge, slagsy, slarot, slartg, slatm1,
362  \$ slaset, sscal, xerbla
363 * ..
364 * .. Intrinsic Functions ..
365  INTRINSIC abs, cos, max, min, mod, real, sin
366 * ..
367 * .. Executable Statements ..
368 *
369 * 1) Decode and Test the input parameters.
370 * Initialize flags & seed.
371 *
372  info = 0
373 *
374 * Quick return if possible
375 *
376  IF( m.EQ.0 .OR. n.EQ.0 )
377  \$ RETURN
378 *
379 * Decode DIST
380 *
381  IF( lsame( dist, 'U' ) ) THEN
382  idist = 1
383  ELSE IF( lsame( dist, 'S' ) ) THEN
384  idist = 2
385  ELSE IF( lsame( dist, 'N' ) ) THEN
386  idist = 3
387  ELSE
388  idist = -1
389  END IF
390 *
391 * Decode SYM
392 *
393  IF( lsame( sym, 'N' ) ) THEN
394  isym = 1
395  irsign = 0
396  ELSE IF( lsame( sym, 'P' ) ) THEN
397  isym = 2
398  irsign = 0
399  ELSE IF( lsame( sym, 'S' ) ) THEN
400  isym = 2
401  irsign = 1
402  ELSE IF( lsame( sym, 'H' ) ) THEN
403  isym = 2
404  irsign = 1
405  ELSE
406  isym = -1
407  END IF
408 *
409 * Decode PACK
410 *
411  isympk = 0
412  IF( lsame( pack, 'N' ) ) THEN
413  ipack = 0
414  ELSE IF( lsame( pack, 'U' ) ) THEN
415  ipack = 1
416  isympk = 1
417  ELSE IF( lsame( pack, 'L' ) ) THEN
418  ipack = 2
419  isympk = 1
420  ELSE IF( lsame( pack, 'C' ) ) THEN
421  ipack = 3
422  isympk = 2
423  ELSE IF( lsame( pack, 'R' ) ) THEN
424  ipack = 4
425  isympk = 3
426  ELSE IF( lsame( pack, 'B' ) ) THEN
427  ipack = 5
428  isympk = 3
429  ELSE IF( lsame( pack, 'Q' ) ) THEN
430  ipack = 6
431  isympk = 2
432  ELSE IF( lsame( pack, 'Z' ) ) THEN
433  ipack = 7
434  ELSE
435  ipack = -1
436  END IF
437 *
438 * Set certain internal parameters
439 *
440  mnmin = min( m, n )
441  llb = min( kl, m-1 )
442  uub = min( ku, n-1 )
443  mr = min( m, n+llb )
444  nc = min( n, m+uub )
445 *
446  IF( ipack.EQ.5 .OR. ipack.EQ.6 ) THEN
447  minlda = uub + 1
448  ELSE IF( ipack.EQ.7 ) THEN
449  minlda = llb + uub + 1
450  ELSE
451  minlda = m
452  END IF
453 *
454 * Use Givens rotation method if bandwidth small enough,
455 * or if LDA is too small to store the matrix unpacked.
456 *
457  givens = .false.
458  IF( isym.EQ.1 ) THEN
459  IF( real( llb+uub ).LT.0.3*real( max( 1, mr+nc ) ) )
460  \$ givens = .true.
461  ELSE
462  IF( 2*llb.LT.m )
463  \$ givens = .true.
464  END IF
465  IF( lda.LT.m .AND. lda.GE.minlda )
466  \$ givens = .true.
467 *
468 * Set INFO if an error
469 *
470  IF( m.LT.0 ) THEN
471  info = -1
472  ELSE IF( m.NE.n .AND. isym.NE.1 ) THEN
473  info = -1
474  ELSE IF( n.LT.0 ) THEN
475  info = -2
476  ELSE IF( idist.EQ.-1 ) THEN
477  info = -3
478  ELSE IF( isym.EQ.-1 ) THEN
479  info = -5
480  ELSE IF( abs( mode ).GT.6 ) THEN
481  info = -7
482  ELSE IF( ( mode.NE.0 .AND. abs( mode ).NE.6 ) .AND. cond.LT.one )
483  \$ THEN
484  info = -8
485  ELSE IF( kl.LT.0 ) THEN
486  info = -10
487  ELSE IF( ku.LT.0 .OR. ( isym.NE.1 .AND. kl.NE.ku ) ) THEN
488  info = -11
489  ELSE IF( ipack.EQ.-1 .OR. ( isympk.EQ.1 .AND. isym.EQ.1 ) .OR.
490  \$ ( isympk.EQ.2 .AND. isym.EQ.1 .AND. kl.GT.0 ) .OR.
491  \$ ( isympk.EQ.3 .AND. isym.EQ.1 .AND. ku.GT.0 ) .OR.
492  \$ ( isympk.NE.0 .AND. m.NE.n ) ) THEN
493  info = -12
494  ELSE IF( lda.LT.max( 1, minlda ) ) THEN
495  info = -14
496  END IF
497 *
498  IF( info.NE.0 ) THEN
499  CALL xerbla( 'SLATMS', -info )
500  RETURN
501  END IF
502 *
503 * Initialize random number generator
504 *
505  DO 10 i = 1, 4
506  iseed( i ) = mod( abs( iseed( i ) ), 4096 )
507  10 CONTINUE
508 *
509  IF( mod( iseed( 4 ), 2 ).NE.1 )
510  \$ iseed( 4 ) = iseed( 4 ) + 1
511 *
512 * 2) Set up D if indicated.
513 *
514 * Compute D according to COND and MODE
515 *
516  CALL slatm1( mode, cond, irsign, idist, iseed, d, mnmin, iinfo )
517  IF( iinfo.NE.0 ) THEN
518  info = 1
519  RETURN
520  END IF
521 *
522 * Choose Top-Down if D is (apparently) increasing,
523 * Bottom-Up if D is (apparently) decreasing.
524 *
525  IF( abs( d( 1 ) ).LE.abs( d( mnmin ) ) ) THEN
526  topdwn = .true.
527  ELSE
528  topdwn = .false.
529  END IF
530 *
531  IF( mode.NE.0 .AND. abs( mode ).NE.6 ) THEN
532 *
533 * Scale by DMAX
534 *
535  temp = abs( d( 1 ) )
536  DO 20 i = 2, mnmin
537  temp = max( temp, abs( d( i ) ) )
538  20 CONTINUE
539 *
540  IF( temp.GT.zero ) THEN
541  alpha = dmax / temp
542  ELSE
543  info = 2
544  RETURN
545  END IF
546 *
547  CALL sscal( mnmin, alpha, d, 1 )
548 *
549  END IF
550 *
551 * 3) Generate Banded Matrix using Givens rotations.
552 * Also the special case of UUB=LLB=0
553 *
554 * Compute Addressing constants to cover all
555 * storage formats. Whether GE, SY, GB, or SB,
556 * upper or lower triangle or both,
557 * the (i,j)-th element is in
558 * A( i - ISKEW*j + IOFFST, j )
559 *
560  IF( ipack.GT.4 ) THEN
561  ilda = lda - 1
562  iskew = 1
563  IF( ipack.GT.5 ) THEN
564  ioffst = uub + 1
565  ELSE
566  ioffst = 1
567  END IF
568  ELSE
569  ilda = lda
570  iskew = 0
571  ioffst = 0
572  END IF
573 *
574 * IPACKG is the format that the matrix is generated in. If this is
575 * different from IPACK, then the matrix must be repacked at the
576 * end. It also signals how to compute the norm, for scaling.
577 *
578  ipackg = 0
579  CALL slaset( 'Full', lda, n, zero, zero, a, lda )
580 *
581 * Diagonal Matrix -- We are done, unless it
582 * is to be stored SP/PP/TP (PACK='R' or 'C')
583 *
584  IF( llb.EQ.0 .AND. uub.EQ.0 ) THEN
585  CALL scopy( mnmin, d, 1, a( 1-iskew+ioffst, 1 ), ilda+1 )
586  IF( ipack.LE.2 .OR. ipack.GE.5 )
587  \$ ipackg = ipack
588 *
589  ELSE IF( givens ) THEN
590 *
591 * Check whether to use Givens rotations,
592 * Householder transformations, or nothing.
593 *
594  IF( isym.EQ.1 ) THEN
595 *
596 * Non-symmetric -- A = U D V
597 *
598  IF( ipack.GT.4 ) THEN
599  ipackg = ipack
600  ELSE
601  ipackg = 0
602  END IF
603 *
604  CALL scopy( mnmin, d, 1, a( 1-iskew+ioffst, 1 ), ilda+1 )
605 *
606  IF( topdwn ) THEN
607  jkl = 0
608  DO 50 jku = 1, uub
609 *
610 * Transform from bandwidth JKL, JKU-1 to JKL, JKU
611 *
612 * Last row actually rotated is M
613 * Last column actually rotated is MIN( M+JKU, N )
614 *
615  DO 40 jr = 1, min( m+jku, n ) + jkl - 1
616  extra = zero
617  angle = twopi*slarnd( 1, iseed )
618  c = cos( angle )
619  s = sin( angle )
620  icol = max( 1, jr-jkl )
621  IF( jr.LT.m ) THEN
622  il = min( n, jr+jku ) + 1 - icol
623  CALL slarot( .true., jr.GT.jkl, .false., il, c,
624  \$ s, a( jr-iskew*icol+ioffst, icol ),
625  \$ ilda, extra, dummy )
626  END IF
627 *
628 * Chase "EXTRA" back up
629 *
630  ir = jr
631  ic = icol
632  DO 30 jch = jr - jkl, 1, -jkl - jku
633  IF( ir.LT.m ) THEN
634  CALL slartg( a( ir+1-iskew*( ic+1 )+ioffst,
635  \$ ic+1 ), extra, c, s, dummy )
636  END IF
637  irow = max( 1, jch-jku )
638  il = ir + 2 - irow
639  temp = zero
640  iltemp = jch.GT.jku
641  CALL slarot( .false., iltemp, .true., il, c, -s,
642  \$ a( irow-iskew*ic+ioffst, ic ),
643  \$ ilda, temp, extra )
644  IF( iltemp ) THEN
645  CALL slartg( a( irow+1-iskew*( ic+1 )+ioffst,
646  \$ ic+1 ), temp, c, s, dummy )
647  icol = max( 1, jch-jku-jkl )
648  il = ic + 2 - icol
649  extra = zero
650  CALL slarot( .true., jch.GT.jku+jkl, .true.,
651  \$ il, c, -s, a( irow-iskew*icol+
652  \$ ioffst, icol ), ilda, extra,
653  \$ temp )
654  ic = icol
655  ir = irow
656  END IF
657  30 CONTINUE
658  40 CONTINUE
659  50 CONTINUE
660 *
661  jku = uub
662  DO 80 jkl = 1, llb
663 *
664 * Transform from bandwidth JKL-1, JKU to JKL, JKU
665 *
666  DO 70 jc = 1, min( n+jkl, m ) + jku - 1
667  extra = zero
668  angle = twopi*slarnd( 1, iseed )
669  c = cos( angle )
670  s = sin( angle )
671  irow = max( 1, jc-jku )
672  IF( jc.LT.n ) THEN
673  il = min( m, jc+jkl ) + 1 - irow
674  CALL slarot( .false., jc.GT.jku, .false., il, c,
675  \$ s, a( irow-iskew*jc+ioffst, jc ),
676  \$ ilda, extra, dummy )
677  END IF
678 *
679 * Chase "EXTRA" back up
680 *
681  ic = jc
682  ir = irow
683  DO 60 jch = jc - jku, 1, -jkl - jku
684  IF( ic.LT.n ) THEN
685  CALL slartg( a( ir+1-iskew*( ic+1 )+ioffst,
686  \$ ic+1 ), extra, c, s, dummy )
687  END IF
688  icol = max( 1, jch-jkl )
689  il = ic + 2 - icol
690  temp = zero
691  iltemp = jch.GT.jkl
692  CALL slarot( .true., iltemp, .true., il, c, -s,
693  \$ a( ir-iskew*icol+ioffst, icol ),
694  \$ ilda, temp, extra )
695  IF( iltemp ) THEN
696  CALL slartg( a( ir+1-iskew*( icol+1 )+ioffst,
697  \$ icol+1 ), temp, c, s, dummy )
698  irow = max( 1, jch-jkl-jku )
699  il = ir + 2 - irow
700  extra = zero
701  CALL slarot( .false., jch.GT.jkl+jku, .true.,
702  \$ il, c, -s, a( irow-iskew*icol+
703  \$ ioffst, icol ), ilda, extra,
704  \$ temp )
705  ic = icol
706  ir = irow
707  END IF
708  60 CONTINUE
709  70 CONTINUE
710  80 CONTINUE
711 *
712  ELSE
713 *
714 * Bottom-Up -- Start at the bottom right.
715 *
716  jkl = 0
717  DO 110 jku = 1, uub
718 *
719 * Transform from bandwidth JKL, JKU-1 to JKL, JKU
720 *
721 * First row actually rotated is M
722 * First column actually rotated is MIN( M+JKU, N )
723 *
724  iendch = min( m, n+jkl ) - 1
725  DO 100 jc = min( m+jku, n ) - 1, 1 - jkl, -1
726  extra = zero
727  angle = twopi*slarnd( 1, iseed )
728  c = cos( angle )
729  s = sin( angle )
730  irow = max( 1, jc-jku+1 )
731  IF( jc.GT.0 ) THEN
732  il = min( m, jc+jkl+1 ) + 1 - irow
733  CALL slarot( .false., .false., jc+jkl.LT.m, il,
734  \$ c, s, a( irow-iskew*jc+ioffst,
735  \$ jc ), ilda, dummy, extra )
736  END IF
737 *
738 * Chase "EXTRA" back down
739 *
740  ic = jc
741  DO 90 jch = jc + jkl, iendch, jkl + jku
742  ilextr = ic.GT.0
743  IF( ilextr ) THEN
744  CALL slartg( a( jch-iskew*ic+ioffst, ic ),
745  \$ extra, c, s, dummy )
746  END IF
747  ic = max( 1, ic )
748  icol = min( n-1, jch+jku )
749  iltemp = jch + jku.LT.n
750  temp = zero
751  CALL slarot( .true., ilextr, iltemp, icol+2-ic,
752  \$ c, s, a( jch-iskew*ic+ioffst, ic ),
753  \$ ilda, extra, temp )
754  IF( iltemp ) THEN
755  CALL slartg( a( jch-iskew*icol+ioffst,
756  \$ icol ), temp, c, s, dummy )
757  il = min( iendch, jch+jkl+jku ) + 2 - jch
758  extra = zero
759  CALL slarot( .false., .true.,
760  \$ jch+jkl+jku.LE.iendch, il, c, s,
761  \$ a( jch-iskew*icol+ioffst,
762  \$ icol ), ilda, temp, extra )
763  ic = icol
764  END IF
765  90 CONTINUE
766  100 CONTINUE
767  110 CONTINUE
768 *
769  jku = uub
770  DO 140 jkl = 1, llb
771 *
772 * Transform from bandwidth JKL-1, JKU to JKL, JKU
773 *
774 * First row actually rotated is MIN( N+JKL, M )
775 * First column actually rotated is N
776 *
777  iendch = min( n, m+jku ) - 1
778  DO 130 jr = min( n+jkl, m ) - 1, 1 - jku, -1
779  extra = zero
780  angle = twopi*slarnd( 1, iseed )
781  c = cos( angle )
782  s = sin( angle )
783  icol = max( 1, jr-jkl+1 )
784  IF( jr.GT.0 ) THEN
785  il = min( n, jr+jku+1 ) + 1 - icol
786  CALL slarot( .true., .false., jr+jku.LT.n, il,
787  \$ c, s, a( jr-iskew*icol+ioffst,
788  \$ icol ), ilda, dummy, extra )
789  END IF
790 *
791 * Chase "EXTRA" back down
792 *
793  ir = jr
794  DO 120 jch = jr + jku, iendch, jkl + jku
795  ilextr = ir.GT.0
796  IF( ilextr ) THEN
797  CALL slartg( a( ir-iskew*jch+ioffst, jch ),
798  \$ extra, c, s, dummy )
799  END IF
800  ir = max( 1, ir )
801  irow = min( m-1, jch+jkl )
802  iltemp = jch + jkl.LT.m
803  temp = zero
804  CALL slarot( .false., ilextr, iltemp, irow+2-ir,
805  \$ c, s, a( ir-iskew*jch+ioffst,
806  \$ jch ), ilda, extra, temp )
807  IF( iltemp ) THEN
808  CALL slartg( a( irow-iskew*jch+ioffst, jch ),
809  \$ temp, c, s, dummy )
810  il = min( iendch, jch+jkl+jku ) + 2 - jch
811  extra = zero
812  CALL slarot( .true., .true.,
813  \$ jch+jkl+jku.LE.iendch, il, c, s,
814  \$ a( irow-iskew*jch+ioffst, jch ),
815  \$ ilda, temp, extra )
816  ir = irow
817  END IF
818  120 CONTINUE
819  130 CONTINUE
820  140 CONTINUE
821  END IF
822 *
823  ELSE
824 *
825 * Symmetric -- A = U D U'
826 *
827  ipackg = ipack
828  ioffg = ioffst
829 *
830  IF( topdwn ) THEN
831 *
832 * Top-Down -- Generate Upper triangle only
833 *
834  IF( ipack.GE.5 ) THEN
835  ipackg = 6
836  ioffg = uub + 1
837  ELSE
838  ipackg = 1
839  END IF
840  CALL scopy( mnmin, d, 1, a( 1-iskew+ioffg, 1 ), ilda+1 )
841 *
842  DO 170 k = 1, uub
843  DO 160 jc = 1, n - 1
844  irow = max( 1, jc-k )
845  il = min( jc+1, k+2 )
846  extra = zero
847  temp = a( jc-iskew*( jc+1 )+ioffg, jc+1 )
848  angle = twopi*slarnd( 1, iseed )
849  c = cos( angle )
850  s = sin( angle )
851  CALL slarot( .false., jc.GT.k, .true., il, c, s,
852  \$ a( irow-iskew*jc+ioffg, jc ), ilda,
853  \$ extra, temp )
854  CALL slarot( .true., .true., .false.,
855  \$ min( k, n-jc )+1, c, s,
856  \$ a( ( 1-iskew )*jc+ioffg, jc ), ilda,
857  \$ temp, dummy )
858 *
859 * Chase EXTRA back up the matrix
860 *
861  icol = jc
862  DO 150 jch = jc - k, 1, -k
863  CALL slartg( a( jch+1-iskew*( icol+1 )+ioffg,
864  \$ icol+1 ), extra, c, s, dummy )
865  temp = a( jch-iskew*( jch+1 )+ioffg, jch+1 )
866  CALL slarot( .true., .true., .true., k+2, c, -s,
867  \$ a( ( 1-iskew )*jch+ioffg, jch ),
868  \$ ilda, temp, extra )
869  irow = max( 1, jch-k )
870  il = min( jch+1, k+2 )
871  extra = zero
872  CALL slarot( .false., jch.GT.k, .true., il, c,
873  \$ -s, a( irow-iskew*jch+ioffg, jch ),
874  \$ ilda, extra, temp )
875  icol = jch
876  150 CONTINUE
877  160 CONTINUE
878  170 CONTINUE
879 *
880 * If we need lower triangle, copy from upper. Note that
881 * the order of copying is chosen to work for 'q' -> 'b'
882 *
883  IF( ipack.NE.ipackg .AND. ipack.NE.3 ) THEN
884  DO 190 jc = 1, n
885  irow = ioffst - iskew*jc
886  DO 180 jr = jc, min( n, jc+uub )
887  a( jr+irow, jc ) = a( jc-iskew*jr+ioffg, jr )
888  180 CONTINUE
889  190 CONTINUE
890  IF( ipack.EQ.5 ) THEN
891  DO 210 jc = n - uub + 1, n
892  DO 200 jr = n + 2 - jc, uub + 1
893  a( jr, jc ) = zero
894  200 CONTINUE
895  210 CONTINUE
896  END IF
897  IF( ipackg.EQ.6 ) THEN
898  ipackg = ipack
899  ELSE
900  ipackg = 0
901  END IF
902  END IF
903  ELSE
904 *
905 * Bottom-Up -- Generate Lower triangle only
906 *
907  IF( ipack.GE.5 ) THEN
908  ipackg = 5
909  IF( ipack.EQ.6 )
910  \$ ioffg = 1
911  ELSE
912  ipackg = 2
913  END IF
914  CALL scopy( mnmin, d, 1, a( 1-iskew+ioffg, 1 ), ilda+1 )
915 *
916  DO 240 k = 1, uub
917  DO 230 jc = n - 1, 1, -1
918  il = min( n+1-jc, k+2 )
919  extra = zero
920  temp = a( 1+( 1-iskew )*jc+ioffg, jc )
921  angle = twopi*slarnd( 1, iseed )
922  c = cos( angle )
923  s = -sin( angle )
924  CALL slarot( .false., .true., n-jc.GT.k, il, c, s,
925  \$ a( ( 1-iskew )*jc+ioffg, jc ), ilda,
926  \$ temp, extra )
927  icol = max( 1, jc-k+1 )
928  CALL slarot( .true., .false., .true., jc+2-icol, c,
929  \$ s, a( jc-iskew*icol+ioffg, icol ),
930  \$ ilda, dummy, temp )
931 *
932 * Chase EXTRA back down the matrix
933 *
934  icol = jc
935  DO 220 jch = jc + k, n - 1, k
936  CALL slartg( a( jch-iskew*icol+ioffg, icol ),
937  \$ extra, c, s, dummy )
938  temp = a( 1+( 1-iskew )*jch+ioffg, jch )
939  CALL slarot( .true., .true., .true., k+2, c, s,
940  \$ a( jch-iskew*icol+ioffg, icol ),
941  \$ ilda, extra, temp )
942  il = min( n+1-jch, k+2 )
943  extra = zero
944  CALL slarot( .false., .true., n-jch.GT.k, il, c,
945  \$ s, a( ( 1-iskew )*jch+ioffg, jch ),
946  \$ ilda, temp, extra )
947  icol = jch
948  220 CONTINUE
949  230 CONTINUE
950  240 CONTINUE
951 *
952 * If we need upper triangle, copy from lower. Note that
953 * the order of copying is chosen to work for 'b' -> 'q'
954 *
955  IF( ipack.NE.ipackg .AND. ipack.NE.4 ) THEN
956  DO 260 jc = n, 1, -1
957  irow = ioffst - iskew*jc
958  DO 250 jr = jc, max( 1, jc-uub ), -1
959  a( jr+irow, jc ) = a( jc-iskew*jr+ioffg, jr )
960  250 CONTINUE
961  260 CONTINUE
962  IF( ipack.EQ.6 ) THEN
963  DO 280 jc = 1, uub
964  DO 270 jr = 1, uub + 1 - jc
965  a( jr, jc ) = zero
966  270 CONTINUE
967  280 CONTINUE
968  END IF
969  IF( ipackg.EQ.5 ) THEN
970  ipackg = ipack
971  ELSE
972  ipackg = 0
973  END IF
974  END IF
975  END IF
976  END IF
977 *
978  ELSE
979 *
980 * 4) Generate Banded Matrix by first
981 * Rotating by random Unitary matrices,
982 * then reducing the bandwidth using Householder
983 * transformations.
984 *
985 * Note: we should get here only if LDA .ge. N
986 *
987  IF( isym.EQ.1 ) THEN
988 *
989 * Non-symmetric -- A = U D V
990 *
991  CALL slagge( mr, nc, llb, uub, d, a, lda, iseed, work,
992  \$ iinfo )
993  ELSE
994 *
995 * Symmetric -- A = U D U'
996 *
997  CALL slagsy( m, llb, d, a, lda, iseed, work, iinfo )
998 *
999  END IF
1000  IF( iinfo.NE.0 ) THEN
1001  info = 3
1002  RETURN
1003  END IF
1004  END IF
1005 *
1006 * 5) Pack the matrix
1007 *
1008  IF( ipack.NE.ipackg ) THEN
1009  IF( ipack.EQ.1 ) THEN
1010 *
1011 * 'U' -- Upper triangular, not packed
1012 *
1013  DO 300 j = 1, m
1014  DO 290 i = j + 1, m
1015  a( i, j ) = zero
1016  290 CONTINUE
1017  300 CONTINUE
1018 *
1019  ELSE IF( ipack.EQ.2 ) THEN
1020 *
1021 * 'L' -- Lower triangular, not packed
1022 *
1023  DO 320 j = 2, m
1024  DO 310 i = 1, j - 1
1025  a( i, j ) = zero
1026  310 CONTINUE
1027  320 CONTINUE
1028 *
1029  ELSE IF( ipack.EQ.3 ) THEN
1030 *
1031 * 'C' -- Upper triangle packed Columnwise.
1032 *
1033  icol = 1
1034  irow = 0
1035  DO 340 j = 1, m
1036  DO 330 i = 1, j
1037  irow = irow + 1
1038  IF( irow.GT.lda ) THEN
1039  irow = 1
1040  icol = icol + 1
1041  END IF
1042  a( irow, icol ) = a( i, j )
1043  330 CONTINUE
1044  340 CONTINUE
1045 *
1046  ELSE IF( ipack.EQ.4 ) THEN
1047 *
1048 * 'R' -- Lower triangle packed Columnwise.
1049 *
1050  icol = 1
1051  irow = 0
1052  DO 360 j = 1, m
1053  DO 350 i = j, m
1054  irow = irow + 1
1055  IF( irow.GT.lda ) THEN
1056  irow = 1
1057  icol = icol + 1
1058  END IF
1059  a( irow, icol ) = a( i, j )
1060  350 CONTINUE
1061  360 CONTINUE
1062 *
1063  ELSE IF( ipack.GE.5 ) THEN
1064 *
1065 * 'B' -- The lower triangle is packed as a band matrix.
1066 * 'Q' -- The upper triangle is packed as a band matrix.
1067 * 'Z' -- The whole matrix is packed as a band matrix.
1068 *
1069  IF( ipack.EQ.5 )
1070  \$ uub = 0
1071  IF( ipack.EQ.6 )
1072  \$ llb = 0
1073 *
1074  DO 380 j = 1, uub
1075  DO 370 i = min( j+llb, m ), 1, -1
1076  a( i-j+uub+1, j ) = a( i, j )
1077  370 CONTINUE
1078  380 CONTINUE
1079 *
1080  DO 400 j = uub + 2, n
1081  DO 390 i = j - uub, min( j+llb, m )
1082  a( i-j+uub+1, j ) = a( i, j )
1083  390 CONTINUE
1084  400 CONTINUE
1085  END IF
1086 *
1087 * If packed, zero out extraneous elements.
1088 *
1089 * Symmetric/Triangular Packed --
1090 * zero out everything after A(IROW,ICOL)
1091 *
1092  IF( ipack.EQ.3 .OR. ipack.EQ.4 ) THEN
1093  DO 420 jc = icol, m
1094  DO 410 jr = irow + 1, lda
1095  a( jr, jc ) = zero
1096  410 CONTINUE
1097  irow = 0
1098  420 CONTINUE
1099 *
1100  ELSE IF( ipack.GE.5 ) THEN
1101 *
1102 * Packed Band --
1103 * 1st row is now in A( UUB+2-j, j), zero above it
1104 * m-th row is now in A( M+UUB-j,j), zero below it
1105 * last non-zero diagonal is now in A( UUB+LLB+1,j ),
1106 * zero below it, too.
1107 *
1108  ir1 = uub + llb + 2
1109  ir2 = uub + m + 2
1110  DO 450 jc = 1, n
1111  DO 430 jr = 1, uub + 1 - jc
1112  a( jr, jc ) = zero
1113  430 CONTINUE
1114  DO 440 jr = max( 1, min( ir1, ir2-jc ) ), lda
1115  a( jr, jc ) = zero
1116  440 CONTINUE
1117  450 CONTINUE
1118  END IF
1119  END IF
1120 *
1121  RETURN
1122 *
1123 * End of SLATMS
1124 *
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: slaset.f:110
subroutine slartg(f, g, c, s, r)
SLARTG generates a plane rotation with real cosine and real sine.
Definition: slartg.f90:113
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine slatm1(MODE, COND, IRSIGN, IDIST, ISEED, D, N, INFO)
SLATM1
Definition: slatm1.f:135
subroutine slarot(LROWS, LLEFT, LRIGHT, NL, C, S, A, LDA, XLEFT, XRIGHT)
SLAROT
Definition: slarot.f:226
subroutine slagsy(N, K, D, A, LDA, ISEED, WORK, INFO)
SLAGSY
Definition: slagsy.f:101
real function slarnd(IDIST, ISEED)
SLARND
Definition: slarnd.f:73
subroutine slagge(M, N, KL, KU, D, A, LDA, ISEED, WORK, INFO)
SLAGGE
Definition: slagge.f:113
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:82
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:79
Here is the call graph for this function: